GCC Code Coverage Report


Directory: ./
File: phys/conema3.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 0 137 0.0%
Branches: 0 156 0.0%

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