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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 70 - (show annotations)
Mon Jun 24 15:39:52 2013 UTC (10 years, 11 months ago) by guez
File size: 10327 byte(s)
In procedure, "addfi" access directly the module variable "dtphys"
instead of going through an argument.

In "conflx", do not create a local variable for temperature with
reversed order of vertical levels. Instead, give an actual argument
with reversed order in "physiq".

Changed names of variables "rmd" and "rmv" from module "suphec_m" to
"md" and "mv".

In "hgardfou", print only the first temperature out of range found.

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

  ViewVC Help
Powered by ViewVC 1.1.21