trunk/dyn3d/vlspltqs.f

module vlspltqs_m

  IMPLICIT NONE

contains

  SUBROUTINE vlspltqs(q, pente_max, masse, w, pbaru, pbarv, pdt, p, pk, teta)

    ! From LMDZ4/libf/dyn3d/vlspltqs.F, version 1.2 2005/02/24 12:16:57 fairhead

    ! Authors: P. Le Van, F. Hourdin, F. Forget, F. Codron

    ! Schéma d'advection "pseudo amont"
    ! + test sur humidité spécifique : Q advecté < Qsat aval
    ! (F. Codron, 10/99)

    ! q, pbaru, pbarv, w sont des arguments d'entree pour le sous-programme

    ! pente_max facteur de limitation des pentes: 2 en général
    ! 0 pour un schéma amont
    ! pbaru, pbarv, w flux de masse en u , v , w
    ! pdt pas de temps

    ! teta température potentielle, p pression aux interfaces,
    ! pk exner au milieu des couches nécessaire pour calculer Qsat

    USE dimens_m, ONLY : iim, llm
    use FCTTRE, only: foeew
    USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, llmp1
    USE comconst, ONLY : cpp
    use SUPHEC_M, only: rtt

    ! Arguments:

    REAL masse(ip1jmp1, llm), pente_max
    REAL, intent(in):: pbaru(ip1jmp1, llm), pbarv(ip1jm, llm)
    REAL q(ip1jmp1, llm)
    REAL w(ip1jmp1, llm)
    real, intent(in):: pdt
    REAL, intent(in):: p(ip1jmp1, llmp1)
    real, intent(in):: teta(ip1jmp1, llm)
    real, intent(in):: pk(ip1jmp1, llm)

    ! Local

    INTEGER ij, l

    REAL qsat(ip1jmp1, llm)
    REAL zm(ip1jmp1, llm)
    REAL mu(ip1jmp1, llm)
    REAL mv(ip1jm, llm)
    REAL mw(ip1jmp1, llm+1)
    REAL zq(ip1jmp1, llm)
    REAL temps1, temps2, temps3
    REAL zzpbar, zzw
    LOGICAL testcpu
    SAVE testcpu
    SAVE temps1, temps2, temps3

    REAL qmin, qmax
    DATA qmin, qmax/0., 1.e33/
    DATA testcpu/.false./
    DATA temps1, temps2, temps3/0., 0., 0./

    !--pour rapport de melange saturant--

    REAL retv, r2es, play
    logical zdelta
    REAL tempe(ip1jmp1)

    !------------------------------------------------------------------

    r2es = 380.11733
    retv = 0.6077667

    !-- Calcul de Qsat en chaque point
    !-- approximation: au milieu des couches play(l)=(p(l)+p(l+1))/2
    ! pour eviter une exponentielle.
    DO l = 1, llm
       DO ij = 1, ip1jmp1
          tempe(ij) = teta(ij, l) * pk(ij, l) /cpp
       ENDDO
       DO ij = 1, ip1jmp1
          zdelta = rtt > tempe(ij)
          play = 0.5*(p(ij, l)+p(ij, l+1))
          qsat(ij, l) = MIN(0.5, r2es* FOEEW(tempe(ij), zdelta) / play)
          qsat(ij, l) = qsat(ij, l) / (1. - retv * qsat(ij, l))
       ENDDO
    ENDDO

    zzpbar = 0.5 * pdt
    zzw = pdt
    DO l=1, llm
       DO ij = iip2, ip1jm
          mu(ij, l)=pbaru(ij, l) * zzpbar
       ENDDO
       DO ij=1, ip1jm
          mv(ij, l)=pbarv(ij, l) * zzpbar
       ENDDO
       DO ij=1, ip1jmp1
          mw(ij, l)=w(ij, l) * zzw
       ENDDO
    ENDDO

    DO ij=1, ip1jmp1
       mw(ij, llm+1)=0.
    ENDDO

    CALL SCOPY(ijp1llm, q, 1, zq, 1)
    CALL SCOPY(ijp1llm, masse, 1, zm, 1)

    ! call minmaxq(zq, qmin, qmax, 'avant vlxqs ')
    call vlxqs(zq, pente_max, zm, mu, qsat)

    ! call minmaxq(zq, qmin, qmax, 'avant vlyqs ')

    call vlyqs(zq, pente_max, zm, mv, qsat)

    ! call minmaxq(zq, qmin, qmax, 'avant vlz ')

    call vlz(zq, pente_max, zm, mw)

    ! call minmaxq(zq, qmin, qmax, 'avant vlyqs ')
    ! call minmaxq(zm, qmin, qmax, 'M avant vlyqs ')

    call vlyqs(zq, pente_max, zm, mv, qsat)

    ! call minmaxq(zq, qmin, qmax, 'avant vlxqs ')
    ! call minmaxq(zm, qmin, qmax, 'M avant vlxqs ')

    call vlxqs(zq, pente_max, zm, mu, qsat)

    ! call minmaxq(zq, qmin, qmax, 'apres vlxqs ')
    ! call minmaxq(zm, qmin, qmax, 'M apres vlxqs ')

    DO l=1, llm
       DO ij=1, ip1jmp1
          q(ij, l)=zq(ij, l)
       ENDDO
       DO ij=1, ip1jm+1, iip1
          q(ij+iim, l)=q(ij, l)
       ENDDO
    ENDDO

  END SUBROUTINE vlspltqs

end module vlspltqs_m
1	module vlspltqs_m
2
3	IMPLICIT NONE
4
5	contains
6
7	SUBROUTINE vlspltqs(q, pente_max, masse, w, pbaru, pbarv, pdt, p, pk, teta)
8
9	! From LMDZ4/libf/dyn3d/vlspltqs.F, version 1.2 2005/02/24 12:16:57 fairhead
10
11	! Authors: P. Le Van, F. Hourdin, F. Forget, F. Codron
12
13	! Schéma d'advection "pseudo amont"
14	! + test sur humidité spécifique : Q advecté < Qsat aval
15	! (F. Codron, 10/99)
16
17	! q, pbaru, pbarv, w sont des arguments d'entree pour le sous-programme
18
19	! pente_max facteur de limitation des pentes: 2 en général
20	! 0 pour un schéma amont
21	! pbaru, pbarv, w flux de masse en u , v , w
22	! pdt pas de temps
23
24	! teta température potentielle, p pression aux interfaces,
25	! pk exner au milieu des couches nécessaire pour calculer Qsat
26
27	USE dimens_m, ONLY : iim, llm
28	use FCTTRE, only: foeew
29	USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, llmp1
30	USE comconst, ONLY : cpp
31	use SUPHEC_M, only: rtt
32
33	! Arguments:
34
35	REAL masse(ip1jmp1, llm), pente_max
36	REAL, intent(in):: pbaru(ip1jmp1, llm), pbarv(ip1jm, llm)
37	REAL q(ip1jmp1, llm)
38	REAL w(ip1jmp1, llm)
39	real, intent(in):: pdt
40	REAL, intent(in):: p(ip1jmp1, llmp1)
41	real, intent(in):: teta(ip1jmp1, llm)
42	real, intent(in):: pk(ip1jmp1, llm)
43
44	! Local
45
46	INTEGER ij, l
47
48	REAL qsat(ip1jmp1, llm)
49	REAL zm(ip1jmp1, llm)
50	REAL mu(ip1jmp1, llm)
51	REAL mv(ip1jm, llm)
52	REAL mw(ip1jmp1, llm+1)
53	REAL zq(ip1jmp1, llm)
54	REAL temps1, temps2, temps3
55	REAL zzpbar, zzw
56	LOGICAL testcpu
57	SAVE testcpu
58	SAVE temps1, temps2, temps3
59
60	REAL qmin, qmax
61	DATA qmin, qmax/0., 1.e33/
62	DATA testcpu/.false./
63	DATA temps1, temps2, temps3/0., 0., 0./
64
65	!--pour rapport de melange saturant--
66
67	REAL retv, r2es, play
68	logical zdelta
69	REAL tempe(ip1jmp1)
70
71	!------------------------------------------------------------------
72
73	r2es = 380.11733
74	retv = 0.6077667
75
76	!-- Calcul de Qsat en chaque point
77	!-- approximation: au milieu des couches play(l)=(p(l)+p(l+1))/2
78	! pour eviter une exponentielle.
79	DO l = 1, llm
80	DO ij = 1, ip1jmp1
81	tempe(ij) = teta(ij, l) * pk(ij, l) /cpp
82	ENDDO
83	DO ij = 1, ip1jmp1
84	zdelta = rtt > tempe(ij)
85	play = 0.5*(p(ij, l)+p(ij, l+1))
86	qsat(ij, l) = MIN(0.5, r2es* FOEEW(tempe(ij), zdelta) / play)
87	qsat(ij, l) = qsat(ij, l) / (1. - retv * qsat(ij, l))
88	ENDDO
89	ENDDO
90
91	zzpbar = 0.5 * pdt
92	zzw = pdt
93	DO l=1, llm
94	DO ij = iip2, ip1jm
95	mu(ij, l)=pbaru(ij, l) * zzpbar
96	ENDDO
97	DO ij=1, ip1jm
98	mv(ij, l)=pbarv(ij, l) * zzpbar
99	ENDDO
100	DO ij=1, ip1jmp1
101	mw(ij, l)=w(ij, l) * zzw
102	ENDDO
103	ENDDO
104
105	DO ij=1, ip1jmp1
106	mw(ij, llm+1)=0.
107	ENDDO
108
109	CALL SCOPY(ijp1llm, q, 1, zq, 1)
110	CALL SCOPY(ijp1llm, masse, 1, zm, 1)
111
112	! call minmaxq(zq, qmin, qmax, 'avant vlxqs ')
113	call vlxqs(zq, pente_max, zm, mu, qsat)
114
115	! call minmaxq(zq, qmin, qmax, 'avant vlyqs ')
116
117	call vlyqs(zq, pente_max, zm, mv, qsat)
118
119	! call minmaxq(zq, qmin, qmax, 'avant vlz ')
120
121	call vlz(zq, pente_max, zm, mw)
122
123	! call minmaxq(zq, qmin, qmax, 'avant vlyqs ')
124	! call minmaxq(zm, qmin, qmax, 'M avant vlyqs ')
125
126	call vlyqs(zq, pente_max, zm, mv, qsat)
127
128	! call minmaxq(zq, qmin, qmax, 'avant vlxqs ')
129	! call minmaxq(zm, qmin, qmax, 'M avant vlxqs ')
130
131	call vlxqs(zq, pente_max, zm, mu, qsat)
132
133	! call minmaxq(zq, qmin, qmax, 'apres vlxqs ')
134	! call minmaxq(zm, qmin, qmax, 'M apres vlxqs ')
135
136	DO l=1, llm
137	DO ij=1, ip1jmp1
138	q(ij, l)=zq(ij, l)
139	ENDDO
140	DO ij=1, ip1jm+1, iip1
141	q(ij+iim, l)=q(ij, l)
142	ENDDO
143	ENDDO
144
145	END SUBROUTINE vlspltqs
146
147	end module vlspltqs_m