/[lmdze]/trunk/libf/phylmd/Conflx/flxasc.f90
ViewVC logotype

Annotation of /trunk/libf/phylmd/Conflx/flxasc.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 71 - (hide annotations)
Mon Jul 8 18:12:18 2013 UTC (10 years, 11 months ago) by guez
File size: 10572 byte(s)
No reason to call inidissip in ce0l.

In inidissip, set random seed to 1 beacuse PGI compiler does not
accept all zeros.

dq was computed needlessly in caladvtrac. Arguments masse and dq of
calfis not used.

Replaced real*8 by double precision.

Pass arrays with inverted order of vertical levels to conflx instead
of creating local variables for this inside conflx.

1 guez 70 module flxasc_m
2    
3 guez 62 IMPLICIT none
4    
5 guez 70 contains
6    
7     SUBROUTINE flxasc(pdtime, ptenh, pqenh, pten, pqen, pqsen, pgeo, pgeoh, &
8     pap, paph, pqte, pvervel, ldland, ldcum, ktype, klab, ptu, pqu, plu, &
9     pmfu, pmfub, pentr, pmfus, pmfuq, pmful, plude, pdmfup, kcbot, kctop, &
10     kctop0, kcum, pen_u, pde_u)
11    
12 guez 71 ! This routine does the calculations for cloud ascents for cumulus
13     ! parameterization.
14    
15     USE dimphy, ONLY: klev, klon
16 guez 70 use flxadjtq_m, only: flxadjtq
17     USE suphec_m, ONLY: rcpd, rd, retv, rg, rtt
18     USE yoecumf, ONLY: cmfcmin, cmfctop, cprcon, entrmid, lmfmid
19    
20     REAL, intent(in):: pdtime
21 guez 71 REAL, intent(in):: ptenh(klon, klev)
22     REAL, intent(in):: pqenh(klon, klev)
23     REAL, intent(in):: pten(klon, klev)
24     REAL, intent(in):: pqen(klon, klev)
25     REAL, intent(in):: pqsen(klon, klev)
26 guez 70 REAL, intent(in):: pgeo(klon, klev), pgeoh(klon, klev)
27     REAL pap(klon, klev), paph(klon, klev+1)
28     REAL pqte(klon, klev)
29     REAL pvervel(klon, klev) ! vitesse verticale en Pa/s
30 guez 71 LOGICAL ldland(klon)
31     LOGICAL, intent(inout):: ldcum(klon)
32     INTEGER, intent(inout):: ktype(klon)
33     integer klab(klon, klev)
34 guez 70 REAL ptu(klon, klev), pqu(klon, klev), plu(klon, klev)
35 guez 71 REAL pmfu(klon, klev)
36     REAL, intent(inout):: pmfub(klon)
37     real pentr(klon)
38     real pmfus(klon, klev)
39     REAL pmfuq(klon, klev), pmful(klon, klev)
40 guez 70 REAL plude(klon, klev)
41     REAL pdmfup(klon, klev)
42 guez 71 integer kcbot(klon), kctop(klon)
43 guez 70 INTEGER kctop0(klon)
44 guez 71 integer, intent(out):: kcum
45     REAL pen_u(klon, klev), pde_u(klon, klev)
46 guez 70
47 guez 71 ! Local:
48    
49 guez 70 REAL zqold(klon)
50     REAL zdland(klon)
51     LOGICAL llflag(klon)
52 guez 71 INTEGER k, i, is, icall
53 guez 70 REAL ztglace, zdphi, zqeen, zseen, zscde, zqude
54     REAL zmfusk, zmfuqk, zmfulk, zbuo, zdnoprc, zprcon, zlnew
55    
56     REAL zpbot(klon), zptop(klon), zrho(klon)
57     REAL zdprho, zentr, zpmid, zmftest, zmfmax
58     LOGICAL llo1, llo2
59    
60     REAL zwmax(klon), zzzmb
61     INTEGER klwmin(klon) ! level of maximum vertical velocity
62    
63     !----------------------------------------------------------------------
64    
65     ztglace = RTT - 13.
66    
67     ! Chercher le niveau où la vitesse verticale est maximale :
68    
69     DO i = 1, klon
70     klwmin(i) = klev
71     zwmax(i) = 0.0
72     ENDDO
73    
74     DO k = klev, 3, -1
75     DO i = 1, klon
76     IF (pvervel(i, k) < zwmax(i)) THEN
77     zwmax(i) = pvervel(i, k)
78     klwmin(i) = k
79     ENDIF
80     ENDDO
81     ENDDO
82    
83     ! Set default values:
84    
85     DO i = 1, klon
86     IF (.NOT. ldcum(i)) ktype(i)=0
87     ENDDO
88    
89     DO k=1, klev
90     DO i = 1, klon
91     plu(i, k)=0.
92     pmfu(i, k)=0.
93     pmfus(i, k)=0.
94     pmfuq(i, k)=0.
95     pmful(i, k)=0.
96     plude(i, k)=0.
97     pdmfup(i, k)=0.
98     IF(.NOT. ldcum(i).OR.ktype(i) == 3) klab(i, k)=0
99     IF(.NOT. ldcum(i).AND.paph(i, k) < 4.E4) kctop0(i)=k
100     ENDDO
101     ENDDO
102    
103     DO i = 1, klon
104     IF (ldland(i)) THEN
105     zdland(i)=3.0E4
106     zdphi=pgeoh(i, kctop0(i))-pgeoh(i, kcbot(i))
107 guez 71 IF (ptu(i, kctop0(i)) >= ztglace) zdland(i)=zdphi
108 guez 70 zdland(i)=MAX(3.0E4, zdland(i))
109     zdland(i)=MIN(5.0E4, zdland(i))
110     ENDIF
111     ENDDO
112    
113     ! Initialiser les valeurs au niveau d'ascendance
114    
115     DO i = 1, klon
116     kctop(i) = klev-1
117     IF (.NOT. ldcum(i)) THEN
118     kcbot(i) = klev-1
119     pmfub(i) = 0.
120     pqu(i, klev) = 0.
121     ENDIF
122     pmfu(i, klev) = pmfub(i)
123 guez 71 pmfus(i, klev) = pmfub(i) * (RCPD * ptu(i, klev)+pgeoh(i, klev))
124     pmfuq(i, klev) = pmfub(i) * pqu(i, klev)
125 guez 70 ENDDO
126    
127     DO i = 1, klon
128     ldcum(i) = .FALSE.
129     ENDDO
130    
131     ! Do ascent: subcloud layer (klab=1), clouds (klab=2) by doing
132     ! first dry-adiabatic ascent and then by adjusting t, q and l
133     ! accordingly in flxadjtq, then check for buoyancy and set flags
134     ! accordingly.
135    
136     DO k = klev - 1, 3, -1
137     IF (LMFMID .AND. k < klev - 1 .AND. k > klev / 2) THEN
138     DO i = 1, klon
139     IF (.NOT. ldcum(i) .AND. klab(i, k + 1) == 0 .AND. &
140     pqen(i, k) > 0.9 * pqsen(i, k)) THEN
141     ptu(i, k+1) = pten(i, k) +(pgeo(i, k)-pgeoh(i, k+1))/RCPD
142     pqu(i, k+1) = pqen(i, k)
143     plu(i, k+1) = 0.0
144     zzzmb = MAX(CMFCMIN, -pvervel(i, k)/RG)
145 guez 71 zmfmax = (paph(i, k)-paph(i, k-1))/(RG * pdtime)
146 guez 70 pmfub(i) = MIN(zzzmb, zmfmax)
147     pmfu(i, k+1) = pmfub(i)
148 guez 71 pmfus(i, k+1) = pmfub(i) * (RCPD * ptu(i, k+1)+pgeoh(i, k+1))
149     pmfuq(i, k+1) = pmfub(i) * pqu(i, k+1)
150 guez 70 pmful(i, k+1) = 0.0
151     pdmfup(i, k+1) = 0.0
152     kcbot(i) = k
153     klab(i, k+1) = 1
154     ktype(i) = 3
155     pentr(i) = ENTRMID
156     ENDIF
157     ENDDO
158     ENDIF
159    
160     is = 0
161     DO i = 1, klon
162     is = is + klab(i, k+1)
163     IF (klab(i, k+1) == 0) klab(i, k) = 0
164     llflag(i) = .FALSE.
165     IF (klab(i, k+1) > 0) llflag(i) = .TRUE.
166     ENDDO
167     IF (is == 0) cycle
168    
169     ! Calculer le taux d'entraînement et de détraînement :
170    
171     DO i = 1, klon
172     pen_u(i, k) = 0.0
173     pde_u(i, k) = 0.0
174 guez 71 zrho(i)=paph(i, k+1)/(RD * ptenh(i, k+1))
175 guez 70 zpbot(i)=paph(i, kcbot(i))
176     zptop(i)=paph(i, kctop0(i))
177     ENDDO
178    
179     DO i = 1, klon
180     IF(ldcum(i)) THEN
181 guez 71 zdprho=(paph(i, k+1)-paph(i, k))/(RG * zrho(i))
182     zentr=pentr(i) * pmfu(i, k+1) * zdprho
183 guez 70 llo1=k < kcbot(i)
184     IF(llo1) pde_u(i, k)=zentr
185 guez 71 zpmid=0.5 * (zpbot(i)+zptop(i))
186 guez 70 llo2=llo1.AND.ktype(i) == 2.AND. &
187     (zpbot(i)-paph(i, k) < 0.2E5.OR. &
188     paph(i, k) > zpmid)
189     IF(llo2) pen_u(i, k)=zentr
190     llo2=llo1.AND.(ktype(i) == 1.OR.ktype(i) == 3).AND. &
191 guez 71 (k >= MAX(klwmin(i), kctop0(i)+2).OR.pap(i, k) > zpmid)
192 guez 70 IF(llo2) pen_u(i, k)=zentr
193     llo1=pen_u(i, k) > 0..AND.(ktype(i) == 1.OR.ktype(i) == 2)
194     IF(llo1) THEN
195 guez 71 zentr=zentr * (1.+3. * (1.-MIN(1., (zpbot(i)-pap(i, k))/1.5E4)))
196     pen_u(i, k)=pen_u(i, k) * (1.+3. * (1.-MIN(1., &
197 guez 70 (zpbot(i)-pap(i, k))/1.5E4)))
198 guez 71 pde_u(i, k)=pde_u(i, k) * (1.+3. * (1.-MIN(1., &
199 guez 70 (zpbot(i)-pap(i, k))/1.5E4)))
200     ENDIF
201     IF(llo2.AND.pqenh(i, k+1) > 1.E-5) &
202 guez 71 pen_u(i, k)=zentr+MAX(pqte(i, k), 0.)/pqenh(i, k+1) * &
203     zrho(i) * zdprho
204 guez 70 ENDIF
205     end DO
206    
207     ! Do adiabatic ascent for entraining/detraining plume
208    
209     DO i = 1, klon
210     IF (llflag(i)) THEN
211     IF (k < kcbot(i)) THEN
212     zmftest = pmfu(i, k+1)+pen_u(i, k)-pde_u(i, k)
213 guez 71 zmfmax = MIN(zmftest, (paph(i, k)-paph(i, k-1))/(RG * pdtime))
214 guez 70 pen_u(i, k)=MAX(pen_u(i, k)-MAX(0.0, zmftest-zmfmax), 0.0)
215     ENDIF
216 guez 71 pde_u(i, k)=MIN(pde_u(i, k), 0.75 * pmfu(i, k+1))
217 guez 70 ! calculer le flux de masse du niveau k a partir de celui du k+1
218     pmfu(i, k)=pmfu(i, k+1)+pen_u(i, k)-pde_u(i, k)
219     ! calculer les valeurs Su, Qu et l du niveau k dans le
220     ! panache montant
221 guez 71 zqeen=pqenh(i, k+1) * pen_u(i, k)
222     zseen=(RCPD * ptenh(i, k+1)+pgeoh(i, k+1)) * pen_u(i, k)
223     zscde=(RCPD * ptu(i, k+1)+pgeoh(i, k+1)) * pde_u(i, k)
224     zqude=pqu(i, k+1) * pde_u(i, k)
225     plude(i, k)=plu(i, k+1) * pde_u(i, k)
226 guez 70 zmfusk=pmfus(i, k+1)+zseen-zscde
227     zmfuqk=pmfuq(i, k+1)+zqeen-zqude
228     zmfulk=pmful(i, k+1) -plude(i, k)
229 guez 71 plu(i, k)=zmfulk * (1./MAX(CMFCMIN, pmfu(i, k)))
230     pqu(i, k)=zmfuqk * (1./MAX(CMFCMIN, pmfu(i, k)))
231     ptu(i, k)=(zmfusk * (1./MAX(CMFCMIN, pmfu(i, k)))- &
232 guez 70 pgeoh(i, k))/RCPD
233     ptu(i, k)=MAX(100., ptu(i, k))
234     ptu(i, k)=MIN(400., ptu(i, k))
235     zqold(i)=pqu(i, k)
236     ELSE
237     zqold(i)=0.0
238     ENDIF
239     end DO
240    
241     ! Do corrections for moist ascent by adjusting t, q and l
242    
243     icall = 1
244     CALL flxadjtq(paph(1, k), ptu(1, k), pqu(1, k), llflag, icall)
245    
246     DO i = 1, klon
247     IF(llflag(i).AND.pqu(i, k).NE.zqold(i)) THEN
248     klab(i, k) = 2
249     plu(i, k) = plu(i, k)+zqold(i)-pqu(i, k)
250 guez 71 zbuo = ptu(i, k) * (1.+RETV * pqu(i, k))- &
251     ptenh(i, k) * (1.+RETV * pqenh(i, k))
252 guez 70 IF (klab(i, k+1) == 1) zbuo=zbuo+0.5
253 guez 71 IF (zbuo > 0. .AND. pmfu(i, k) >= 0.1 * pmfub(i)) THEN
254 guez 70 kctop(i) = k
255     ldcum(i) = .TRUE.
256     zdnoprc = 1.5E4
257     IF (ldland(i)) zdnoprc = zdland(i)
258     zprcon = CPRCON
259     IF ((zpbot(i)-paph(i, k)) < zdnoprc) zprcon = 0.0
260 guez 71 zlnew=plu(i, k)/(1.+zprcon * (pgeoh(i, k)-pgeoh(i, k+1)))
261     pdmfup(i, k)=MAX(0., (plu(i, k)-zlnew) * pmfu(i, k))
262 guez 70 plu(i, k)=zlnew
263     ELSE
264     klab(i, k)=0
265     pmfu(i, k)=0.
266     ENDIF
267     ENDIF
268     end DO
269     DO i = 1, klon
270     IF (llflag(i)) THEN
271 guez 71 pmful(i, k)=plu(i, k) * pmfu(i, k)
272     pmfus(i, k)=(RCPD * ptu(i, k)+pgeoh(i, k)) * pmfu(i, k)
273     pmfuq(i, k)=pqu(i, k) * pmfu(i, k)
274 guez 70 ENDIF
275     end DO
276    
277     end DO
278    
279     ! Determine convective fluxes above non-buoyancy level (note:
280     ! cloud variables like t, q and l are not affected by detrainment
281     ! and are already known from previous calculations above).
282    
283     DO i = 1, klon
284     IF (kctop(i) == klev-1) ldcum(i) = .FALSE.
285     kcbot(i) = MAX(kcbot(i), kctop(i))
286     ENDDO
287    
288     ldcum(1)=ldcum(1)
289    
290     is = 0
291     DO i = 1, klon
292     if (ldcum(i)) is = is + 1
293     ENDDO
294     kcum = is
295     IF (is /= 0) then
296     DO i = 1, klon
297     IF (ldcum(i)) THEN
298     k=kctop(i)-1
299 guez 71 pde_u(i, k)=(1.-CMFCTOP) * pmfu(i, k+1)
300     plude(i, k)=pde_u(i, k) * plu(i, k+1)
301 guez 70 pmfu(i, k)=pmfu(i, k+1)-pde_u(i, k)
302     zlnew=plu(i, k)
303 guez 71 pdmfup(i, k)=MAX(0., (plu(i, k)-zlnew) * pmfu(i, k))
304 guez 70 plu(i, k)=zlnew
305 guez 71 pmfus(i, k)=(RCPD * ptu(i, k)+pgeoh(i, k)) * pmfu(i, k)
306     pmfuq(i, k)=pqu(i, k) * pmfu(i, k)
307     pmful(i, k)=plu(i, k) * pmfu(i, k)
308 guez 70 plude(i, k-1)=pmful(i, k)
309     ENDIF
310     end DO
311     end IF
312    
313     END SUBROUTINE flxasc
314    
315     end module flxasc_m

  ViewVC Help
Powered by ViewVC 1.1.21