1 | ! |
---|
2 | ! $Header$ |
---|
3 | ! |
---|
4 | SUBROUTINE conema3 (dtime,paprs,pplay,t,q,u,v,tra,ntra, |
---|
5 | . work1,work2,d_t,d_q,d_u,d_v,d_tra, |
---|
6 | . rain, snow, kbas, ktop, |
---|
7 | . upwd,dnwd,dnwdbis,bas,top,Ma,cape,tvp,rflag, |
---|
8 | . pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr, |
---|
9 | . qcond_incld) |
---|
10 | |
---|
11 | IMPLICIT none |
---|
12 | c====================================================================== |
---|
13 | c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930818 |
---|
14 | c Objet: schema de convection de Emanuel (1991) interface |
---|
15 | c Mai 1998: Interface modifiee pour implementation dans LMDZ |
---|
16 | c====================================================================== |
---|
17 | c Arguments: |
---|
18 | c dtime---input-R-pas d'integration (s) |
---|
19 | c paprs---input-R-pression inter-couches (Pa) |
---|
20 | c pplay---input-R-pression au milieu des couches (Pa) |
---|
21 | c t-------input-R-temperature (K) |
---|
22 | c q-------input-R-humidite specifique (kg/kg) |
---|
23 | c u-------input-R-vitesse du vent zonal (m/s) |
---|
24 | c v-------input-R-vitesse duvent meridien (m/s) |
---|
25 | c tra-----input-R-tableau de rapport de melange des traceurs |
---|
26 | c work*: input et output: deux variables de travail, |
---|
27 | c on peut les mettre a 0 au debut |
---|
28 | c |
---|
29 | C d_t-----output-R-increment de la temperature |
---|
30 | c d_q-----output-R-increment de la vapeur d'eau |
---|
31 | c d_u-----output-R-increment de la vitesse zonale |
---|
32 | c d_v-----output-R-increment de la vitesse meridienne |
---|
33 | c d_tra---output-R-increment du contenu en traceurs |
---|
34 | c rain----output-R-la pluie (mm/s) |
---|
35 | c snow----output-R-la neige (mm/s) |
---|
36 | c kbas----output-R-bas du nuage (integer) |
---|
37 | c ktop----output-R-haut du nuage (integer) |
---|
38 | c upwd----output-R-saturated updraft mass flux (kg/m**2/s) |
---|
39 | c dnwd----output-R-saturated downdraft mass flux (kg/m**2/s) |
---|
40 | c dnwdbis-output-R-unsaturated downdraft mass flux (kg/m**2/s) |
---|
41 | c bas-----output-R-bas du nuage (real) |
---|
42 | c top-----output-R-haut du nuage (real) |
---|
43 | c Ma------output-R-flux ascendant non dilue (kg/m**2/s) |
---|
44 | c cape----output-R-CAPE |
---|
45 | c tvp-----output-R-virtual temperature of the lifted parcel |
---|
46 | c rflag---output-R-flag sur le fonctionnement de convect |
---|
47 | c pbase---output-R-pression a la base du nuage (Pa) |
---|
48 | c bbase---output-R-buoyancy a la base du nuage (K) |
---|
49 | c dtvpdt1-output-R-derivative of parcel virtual temp wrt T1 |
---|
50 | c dtvpdq1-output-R-derivative of parcel virtual temp wrt Q1 |
---|
51 | c dplcldt-output-R-derivative of the PCP pressure wrt T1 |
---|
52 | c dplcldr-output-R-derivative of the PCP pressure wrt Q1 |
---|
53 | c====================================================================== |
---|
54 | c |
---|
55 | #include "dimensions.h" |
---|
56 | #include "dimphy.h" |
---|
57 | #include "conema3.h" |
---|
58 | INTEGER i, l,m,itra |
---|
59 | INTEGER ntra,ntrac !number of tracers; if no tracer transport |
---|
60 | ! is needed, set ntra = 1 (or 0) |
---|
61 | PARAMETER (ntrac=nqmx-2) |
---|
62 | REAL dtime |
---|
63 | c |
---|
64 | REAL d_t2(klon,klev), d_q2(klon,klev) ! sbl |
---|
65 | REAL d_u2(klon,klev), d_v2(klon,klev) ! sbl |
---|
66 | REAL em_d_t2(klev), em_d_q2(klev) ! sbl |
---|
67 | REAL em_d_u2(klev), em_d_v2(klev) ! sbl |
---|
68 | c |
---|
69 | REAL paprs(klon,klev+1), pplay(klon,klev) |
---|
70 | REAL t(klon,klev), q(klon,klev), d_t(klon,klev), d_q(klon,klev) |
---|
71 | REAL u(klon,klev), v(klon,klev), tra(klon,klev,ntra) |
---|
72 | REAL d_u(klon,klev), d_v(klon,klev), d_tra(klon,klev,ntra) |
---|
73 | REAL work1(klon,klev), work2(klon,klev) |
---|
74 | REAL upwd(klon,klev), dnwd(klon,klev), dnwdbis(klon,klev) |
---|
75 | REAL rain(klon) |
---|
76 | REAL snow(klon) |
---|
77 | REAL cape(klon), tvp(klon,klev), rflag(klon) |
---|
78 | REAL pbase(klon), bbase(klon) |
---|
79 | REAL dtvpdt1(klon,klev), dtvpdq1(klon,klev) |
---|
80 | REAL dplcldt(klon), dplcldr(klon) |
---|
81 | INTEGER kbas(klon), ktop(klon) |
---|
82 | |
---|
83 | REAL wd(klon) |
---|
84 | REAL qcond_incld(klon,klev) |
---|
85 | c |
---|
86 | REAL em_t(klev) |
---|
87 | REAL em_q(klev) |
---|
88 | REAL em_qs(klev) |
---|
89 | REAL em_u(klev), em_v(klev), em_tra(klev,ntrac) |
---|
90 | REAL em_ph(klev+1), em_p(klev) |
---|
91 | REAL em_work1(klev), em_work2(klev) |
---|
92 | REAL em_precip, em_d_t(klev), em_d_q(klev) |
---|
93 | REAL em_d_u(klev), em_d_v(klev), em_d_tra(klev,ntrac) |
---|
94 | REAL em_upwd(klev), em_dnwd(klev), em_dnwdbis(klev) |
---|
95 | REAL em_dtvpdt1(klev), em_dtvpdq1(klev) |
---|
96 | REAL em_dplcldt, em_dplcldr |
---|
97 | SAVE em_t,em_q, em_qs, em_ph, em_p, em_work1, em_work2 |
---|
98 | SAVE em_u,em_v, em_tra |
---|
99 | SAVE em_d_u,em_d_v, em_d_tra |
---|
100 | SAVE em_precip, em_d_t, em_d_q, em_upwd, em_dnwd, em_dnwdbis |
---|
101 | INTEGER em_bas, em_top |
---|
102 | SAVE em_bas, em_top |
---|
103 | |
---|
104 | REAL em_wd |
---|
105 | REAL em_qcond(klev) |
---|
106 | REAL em_qcondc(klev) |
---|
107 | c |
---|
108 | REAL zx_t, zx_qs, zdelta, zcor |
---|
109 | INTEGER iflag |
---|
110 | REAL sigsum |
---|
111 | ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
112 | c VARIABLES A SORTIR |
---|
113 | cccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
114 | |
---|
115 | REAL emmip(klev) !variation de flux ascnon dilue i et i+1 |
---|
116 | SAVE emmip |
---|
117 | real emMke(klev) |
---|
118 | save emMke |
---|
119 | real top |
---|
120 | real bas |
---|
121 | real emMa(klev) |
---|
122 | save emMa |
---|
123 | real Ma(klon,klev) |
---|
124 | real Ment(klev,klev) |
---|
125 | real Qent(klev,klev) |
---|
126 | real TPS(klev),TLS(klev) |
---|
127 | real SIJ(klev,klev) |
---|
128 | real em_CAPE, em_TVP(klev) |
---|
129 | real em_pbase, em_bbase |
---|
130 | integer iw,j,k,ix,iy |
---|
131 | |
---|
132 | c -- sb: pour schema nuages: |
---|
133 | |
---|
134 | integer iflagcon |
---|
135 | integer em_ifc(klev) |
---|
136 | |
---|
137 | real em_pradj |
---|
138 | real em_cldf(klev), em_cldq(klev) |
---|
139 | real em_ftadj(klev), em_fradj(klev) |
---|
140 | |
---|
141 | integer ifc(klon,klev) |
---|
142 | real pradj(klon) |
---|
143 | real cldf(klon,klev), cldq(klon,klev) |
---|
144 | real ftadj(klon,klev), fqadj(klon,klev) |
---|
145 | |
---|
146 | c sb -- |
---|
147 | |
---|
148 | ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
149 | c |
---|
150 | #include "YOMCST.h" |
---|
151 | #include "YOETHF.h" |
---|
152 | #include "FCTTRE.h" |
---|
153 | |
---|
154 | qcond_incld(:,:) = 0. |
---|
155 | c |
---|
156 | c$$$ print*,'debut conema' |
---|
157 | |
---|
158 | DO 999 i = 1, klon |
---|
159 | DO l = 1, klev+1 |
---|
160 | em_ph(l) = paprs(i,l) / 100.0 |
---|
161 | ENDDO |
---|
162 | c |
---|
163 | DO l = 1, klev |
---|
164 | em_p(l) = pplay(i,l) / 100.0 |
---|
165 | em_t(l) = t(i,l) |
---|
166 | em_q(l) = q(i,l) |
---|
167 | em_u(l) = u(i,l) |
---|
168 | em_v(l) = v(i,l) |
---|
169 | do itra = 1, ntra |
---|
170 | em_tra(l,itra) = tra(i,l,itra) |
---|
171 | enddo |
---|
172 | c$$$ print*,'em_t',em_t |
---|
173 | c$$$ print*,'em_q',em_q |
---|
174 | c$$$ print*,'em_qs',em_qs |
---|
175 | c$$$ print*,'em_u',em_u |
---|
176 | c$$$ print*,'em_v',em_v |
---|
177 | c$$$ print*,'em_tra',em_tra |
---|
178 | c$$$ print*,'em_p',em_p |
---|
179 | |
---|
180 | |
---|
181 | c |
---|
182 | zx_t = em_t(l) |
---|
183 | zdelta=MAX(0.,SIGN(1.,rtt-zx_t)) |
---|
184 | zx_qs= r2es * FOEEW(zx_t,zdelta)/em_p(l)/100.0 |
---|
185 | zx_qs=MIN(0.5,zx_qs) |
---|
186 | c$$$ print*,'zx_qs',zx_qs |
---|
187 | zcor=1./(1.-retv*zx_qs) |
---|
188 | zx_qs=zx_qs*zcor |
---|
189 | em_qs(l) = zx_qs |
---|
190 | c$$$ print*,'em_qs',em_qs |
---|
191 | c |
---|
192 | em_work1(l) = work1(i,l) |
---|
193 | em_work2(l) = work2(i,l) |
---|
194 | emMke(l)=0 |
---|
195 | c emMa(l)=0 |
---|
196 | c Ma(i,l)=0 |
---|
197 | |
---|
198 | em_dtvpdt1(l) = 0. |
---|
199 | em_dtvpdq1(l) = 0. |
---|
200 | dtvpdt1(i,l) = 0. |
---|
201 | dtvpdq1(i,l) = 0. |
---|
202 | ENDDO |
---|
203 | c |
---|
204 | em_dplcldt = 0. |
---|
205 | em_dplcldr = 0. |
---|
206 | rain(i) = 0.0 |
---|
207 | snow(i) = 0.0 |
---|
208 | kbas(i) = 1 |
---|
209 | ktop(i) = 1 |
---|
210 | c ajout SB: |
---|
211 | bas = 1 |
---|
212 | top = 1 |
---|
213 | |
---|
214 | |
---|
215 | c sb3d write(*,1792) (em_work1(m),m=1,klev) |
---|
216 | 1792 format('sig avant convect ',/,10(1X,E13.5)) |
---|
217 | c |
---|
218 | c sb d write(*,1793) (em_work2(m),m=1,klev) |
---|
219 | 1793 format('w avant convect ',/,10(1X,E13.5)) |
---|
220 | |
---|
221 | c$$$ print*,'avant convect' |
---|
222 | ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
223 | c |
---|
224 | |
---|
225 | c print*,'avant convect i=',i |
---|
226 | CALL convect3(dtime,epmax,ok_adj_ema, |
---|
227 | . em_t, em_q, em_qs,em_u ,em_v , |
---|
228 | . em_tra, em_p, em_ph, |
---|
229 | . klev, klev+1, klev-1,ntra, dtime, iflag, |
---|
230 | . em_d_t, em_d_q,em_d_u,em_d_v, |
---|
231 | . em_d_tra, em_precip, |
---|
232 | . em_bas, em_top,em_upwd, em_dnwd, em_dnwdbis, |
---|
233 | . em_work1, em_work2,emmip,emMke,emMa,Ment, |
---|
234 | . Qent,TPS,TLS,SIJ,em_CAPE,em_TVP,em_pbase,em_bbase, |
---|
235 | . em_dtvpdt1,em_dtvpdq1,em_dplcldt,em_dplcldr, ! sbl |
---|
236 | . em_d_t2,em_d_q2,em_d_u2,em_d_v2,em_wd,em_qcond,em_qcondc)!sbl |
---|
237 | c print*,'apres convect ' |
---|
238 | c |
---|
239 | ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
240 | c |
---|
241 | c -- sb: Appel schema statistique de nuages couple a la convection |
---|
242 | c (Bony et Emanuel 2001): |
---|
243 | |
---|
244 | c -- creer cvthermo.h qui contiendra les cstes thermo de LMDZ: |
---|
245 | |
---|
246 | iflagcon = 3 |
---|
247 | c CALL cv_thermo(iflagcon) |
---|
248 | |
---|
249 | c -- appel schema de nuages: |
---|
250 | |
---|
251 | c CALL CLOUDS_SUB_LS(klev,em_q,em_qs,em_t |
---|
252 | c i ,em_p,em_ph,dtime,em_qcondc |
---|
253 | c o ,em_cldf,em_cldq,em_pradj,em_ftadj,em_fradj,em_ifc) |
---|
254 | |
---|
255 | do k = 1, klev |
---|
256 | cldf(i,k) = em_cldf(k) ! cloud fraction (0-1) |
---|
257 | cldq(i,k) = em_cldq(k) ! in-cloud water content (kg/kg) |
---|
258 | ftadj(i,k) = em_ftadj(k) ! (dT/dt)_{LS adj} (K/s) |
---|
259 | fqadj(i,k) = em_fradj(k) ! (dq/dt)_{LS adj} (kg/kg/s) |
---|
260 | ifc(i,k) = em_ifc(k) ! flag convergence clouds_gno (1 ou 2) |
---|
261 | enddo |
---|
262 | pradj(i) = em_pradj ! precip from LS supersat adj (mm/day) |
---|
263 | |
---|
264 | c sb -- |
---|
265 | c |
---|
266 | c SB: |
---|
267 | if (iflag.ne.1 .and. iflag.ne.4) then |
---|
268 | em_CAPE = 0. |
---|
269 | do l = 1, klev |
---|
270 | em_upwd(l) = 0. |
---|
271 | em_dnwd(l) = 0. |
---|
272 | em_dnwdbis(l) = 0. |
---|
273 | emMa(l) = 0. |
---|
274 | em_TVP(l) = 0. |
---|
275 | enddo |
---|
276 | endif |
---|
277 | c fin SB |
---|
278 | c |
---|
279 | c If sig has been set to zero, then set Ma to zero |
---|
280 | c |
---|
281 | sigsum = 0. |
---|
282 | do k = 1,klev |
---|
283 | sigsum = sigsum + em_work1(k) |
---|
284 | enddo |
---|
285 | if (sigsum .eq. 0.0) then |
---|
286 | do k = 1,klev |
---|
287 | emMa(k) = 0. |
---|
288 | enddo |
---|
289 | endif |
---|
290 | c |
---|
291 | c sb3d print*,'i, iflag=',i,iflag |
---|
292 | c |
---|
293 | ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
294 | c |
---|
295 | c SORTIE DES ICB ET INB |
---|
296 | c en fait inb et icb correspondent au niveau ou se trouve |
---|
297 | c le nuage,le numero d'interface |
---|
298 | cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
---|
299 | |
---|
300 | c modif SB: |
---|
301 | if (iflag.EQ.1 .or. iflag.EQ.4) then |
---|
302 | top=em_top |
---|
303 | bas=em_bas |
---|
304 | kbas(i) = em_bas |
---|
305 | ktop(i) = em_top |
---|
306 | endif |
---|
307 | |
---|
308 | pbase(i) = em_pbase |
---|
309 | bbase(i) = em_bbase |
---|
310 | rain(i) = em_precip/ 86400.0 |
---|
311 | snow(i) = 0.0 |
---|
312 | cape(i) = em_CAPE |
---|
313 | wd(i) = em_wd |
---|
314 | rflag(i) = float(iflag) |
---|
315 | c SB kbas(i) = em_bas |
---|
316 | c SB ktop(i) = em_top |
---|
317 | dplcldt(i) = em_dplcldt |
---|
318 | dplcldr(i) = em_dplcldr |
---|
319 | DO l = 1, klev |
---|
320 | d_t2(i,l) = dtime * em_d_t2(l) |
---|
321 | d_q2(i,l) = dtime * em_d_q2(l) |
---|
322 | d_u2(i,l) = dtime * em_d_u2(l) |
---|
323 | d_v2(i,l) = dtime * em_d_v2(l) |
---|
324 | |
---|
325 | d_t(i,l) = dtime * em_d_t(l) |
---|
326 | d_q(i,l) = dtime * em_d_q(l) |
---|
327 | d_u(i,l) = dtime * em_d_u(l) |
---|
328 | d_v(i,l) = dtime * em_d_v(l) |
---|
329 | do itra = 1, ntra |
---|
330 | d_tra(i,l,itra) = dtime * em_d_tra(l,itra) |
---|
331 | enddo |
---|
332 | upwd(i,l) = em_upwd(l) |
---|
333 | dnwd(i,l) = em_dnwd(l) |
---|
334 | dnwdbis(i,l) = em_dnwdbis(l) |
---|
335 | work1(i,l) = em_work1(l) |
---|
336 | work2(i,l) = em_work2(l) |
---|
337 | Ma(i,l)=emMa(l) |
---|
338 | tvp(i,l)=em_TVP(l) |
---|
339 | dtvpdt1(i,l) = em_dtvpdt1(l) |
---|
340 | dtvpdq1(i,l) = em_dtvpdq1(l) |
---|
341 | |
---|
342 | if (iflag_clw.eq.0) then |
---|
343 | qcond_incld(i,l) = em_qcondc(l) |
---|
344 | else if (iflag_clw.eq.1) then |
---|
345 | qcond_incld(i,l) = em_qcond(l) |
---|
346 | endif |
---|
347 | ENDDO |
---|
348 | 999 CONTINUE |
---|
349 | |
---|
350 | c On calcule une eau liquide diagnostique en fonction de la |
---|
351 | c precip. |
---|
352 | if ( iflag_clw.eq.2 ) then |
---|
353 | do l=1,klev |
---|
354 | do i=1,klon |
---|
355 | if (ktop(i)-kbas(i).gt.0.and. |
---|
356 | s l.ge.kbas(i).and.l.le.ktop(i)) then |
---|
357 | qcond_incld(i,l)=rain(i)*8.e4 |
---|
358 | c s *(pplay(i,l )-paprs(i,ktop(i)+1)) |
---|
359 | s /(pplay(i,kbas(i))-pplay(i,ktop(i))) |
---|
360 | c s **2 |
---|
361 | else |
---|
362 | qcond_incld(i,l)=0. |
---|
363 | endif |
---|
364 | enddo |
---|
365 | print*,'l=',l,', qcond_incld=',qcond_incld(1,l) |
---|
366 | enddo |
---|
367 | endif |
---|
368 | |
---|
369 | |
---|
370 | RETURN |
---|
371 | END |
---|
372 | |
---|