GCC Code Coverage Report


Directory: ./
File: phys/cv3p1_closure.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 237 287 82.6%
Branches: 272 358 76.0%

Line Branch Exec Source
1
2 ! $Id: cv3p1_closure.F90 3671 2020-04-29 13:48:22Z jyg $
3
4 480 SUBROUTINE cv3p1_closure(nloc, ncum, nd, icb, inb, pbase, plcl, p, ph, tv, &
5 240 tvp, buoy, supmax, ok_inhib, ale, alp, omega,sig, w0, ptop2, cape, cin, m, &
6 iflag, coef, plim1, plim2, asupmax, supmax0, asupmaxmin, cbmf, plfc, &
7 wbeff)
8
9
10 ! **************************************************************
11 ! *
12 ! CV3P1_CLOSURE *
13 ! Ale & Alp Closure of Convect3 *
14 ! *
15 ! written by : Kerry Emanuel *
16 ! vectorization: S. Bony *
17 ! modified by : Jean-Yves Grandpeix, 18/06/2003, 19.32.10 *
18 ! Julie Frohwirth, 14/10/2005 17.44.22 *
19 ! **************************************************************
20
21 USE print_control_mod, ONLY: prt_level, lunout
22 IMPLICIT NONE
23
24 include "cvthermo.h"
25 include "cv3param.h"
26 include "YOMCST2.h"
27 include "YOMCST.h"
28 include "conema3.h"
29
30 ! input:
31 INTEGER, INTENT (IN) :: ncum, nd, nloc
32 INTEGER, DIMENSION (nloc), INTENT (IN) :: icb, inb
33 REAL, DIMENSION (nloc), INTENT (IN) :: pbase, plcl
34 REAL, DIMENSION (nloc, nd), INTENT (IN) :: p
35 REAL, DIMENSION (nloc, nd+1), INTENT (IN) :: ph
36 REAL, DIMENSION (nloc, nd), INTENT (IN) :: tv, tvp, buoy
37 REAL, DIMENSION (nloc, nd), INTENT (IN) :: supmax
38 LOGICAL, INTENT (IN) :: ok_inhib ! enable convection inhibition by dryness
39 REAL, DIMENSION (nloc), INTENT (IN) :: ale, alp
40 REAL, DIMENSION (nloc, nd), INTENT (IN) :: omega
41
42 ! input/output:
43 INTEGER, DIMENSION (nloc), INTENT (INOUT) :: iflag
44 REAL, DIMENSION (nloc, nd), INTENT (INOUT) :: sig, w0
45 REAL, DIMENSION (nloc), INTENT (INOUT) :: ptop2
46
47 ! output:
48 REAL, DIMENSION (nloc), INTENT (OUT) :: cape, cin
49 REAL, DIMENSION (nloc, nd), INTENT (OUT) :: m
50 REAL, DIMENSION (nloc), INTENT (OUT) :: plim1, plim2
51 REAL, DIMENSION (nloc, nd), INTENT (OUT) :: asupmax
52 REAL, DIMENSION (nloc), INTENT (OUT) :: supmax0
53 REAL, DIMENSION (nloc), INTENT (OUT) :: asupmaxmin
54 REAL, DIMENSION (nloc), INTENT (OUT) :: cbmf, plfc
55 REAL, DIMENSION (nloc), INTENT (OUT) :: wbeff
56
57 ! local variables:
58 480 INTEGER il, i, j, k, icbmax, i0(nloc), klfc(nloc)
59 REAL deltap, fac, w, amu
60 REAL rhodp, dz
61 REAL pbmxup
62 480 REAL dtmin(nloc, nd), sigold(nloc, nd)
63 480 REAL coefmix(nloc, nd)
64 480 REAL pzero(nloc), ptop2old(nloc)
65 480 REAL cina(nloc), cinb(nloc)
66 INTEGER ibeg(nloc)
67 480 INTEGER nsupmax(nloc)
68 REAL supcrit, temp(nloc, nd)
69 480 REAL p1(nloc), pmin(nloc)
70 480 REAL asupmax0(nloc)
71 480 LOGICAL ok(nloc)
72 480 REAL siglim(nloc, nd), wlim(nloc, nd), mlim(nloc, nd)
73 480 REAL wb2(nloc)
74 480 REAL cbmflim(nloc), cbmf1(nloc), cbmfmax(nloc)
75 480 REAL cbmflast(nloc)
76 REAL coef(nloc)
77 REAL xp(nloc), xq(nloc), xr(nloc), discr(nloc), b3(nloc), b4(nloc)
78 REAL theta(nloc), bb(nloc)
79 REAL term1, term2, term3
80 480 REAL alp2(nloc) ! Alp with offset
81 !CR: variables for new erosion of adiabiatic ascent
82 480 REAL mad(nloc, nd), me(nloc, nd), betalim(nloc, nd), beta_coef(nloc, nd)
83 480 REAL med(nloc, nd), md(nloc,nd)
84 !jyg<
85 ! coef_peel is now in the common cv3_param
86 !! REAL coef_peel
87 !! PARAMETER (coef_peel=0.25)
88 !>jyg
89
90 REAL sigmax
91 PARAMETER (sigmax=0.1)
92
93 CHARACTER (LEN=20) :: modname = 'cv3p1_closure'
94 CHARACTER (LEN=80) :: abort_message
95
96 ! print *,' -> cv3p1_closure, Ale ',ale(1)
97
98
99 ! -------------------------------------------------------
100 ! -- Initialization
101 ! -------------------------------------------------------
102
103
104
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
105 127097 alp2(il) = max(alp(il), 1.E-5)
106 ! IM
107 127337 alp2(il) = max(alp(il), 1.E-12)
108 END DO
109
110 pbmxup = 50. ! PBMXUP+PBCRIT = cloud depth above which mixed updraughts
111 ! exist (if any)
112
113
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param nloc ncum nd icb inb nl', nloc, &
114 ncum, nd, icb(nloc), inb(nloc), nl
115
2/2
✓ Branch 0 taken 9360 times.
✓ Branch 1 taken 240 times.
9600 DO k = 1, nd !jyg: initialization up to nd
116
2/2
✓ Branch 0 taken 4956783 times.
✓ Branch 1 taken 9360 times.
4966383 DO il = 1, ncum
117 4966143 m(il, k) = 0.0
118 END DO
119 END DO
120
121 !CR: initializations for erosion of adiabatic ascent
122
2/2
✓ Branch 0 taken 9360 times.
✓ Branch 1 taken 240 times.
9600 DO k = 1,nd !jyg: initialization up to nd
123
2/2
✓ Branch 0 taken 4956783 times.
✓ Branch 1 taken 9360 times.
4966383 DO il = 1, ncum
124 4956783 mad(il,k)=0.
125 4956783 me(il,k)=0.
126 4956783 betalim(il,k)=1.
127 4966143 wlim(il,k)=0.
128 ENDDO
129 ENDDO
130
131 ! -------------------------------------------------------
132 ! -- Reset sig(i) and w0(i) for i>inb and i<icb
133 ! -------------------------------------------------------
134
135 ! update sig and w0 above LNB:
136
137
2/2
✓ Branch 0 taken 6240 times.
✓ Branch 1 taken 240 times.
6480 DO k = 1, nl - 1
138
2/2
✓ Branch 0 taken 3304522 times.
✓ Branch 1 taken 6240 times.
3311002 DO il = 1, ncum
139
3/4
✓ Branch 0 taken 3304522 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1639367 times.
✓ Branch 3 taken 1665155 times.
3310762 IF ((inb(il)<(nl-1)) .AND. (k>=(inb(il)+1))) THEN
140 sig(il, k) = beta*sig(il, k) + 2.*alpha*buoy(il, inb(il))*abs(buoy(il &
141 1639367 ,inb(il)))
142 1639367 sig(il, k) = amax1(sig(il,k), 0.0)
143 1639367 w0(il, k) = beta*w0(il, k)
144 END IF
145 END DO
146 END DO
147
148 ! if(prt.level.GE.20) print*,'cv3p1_param apres 100'
149 ! compute icbmax:
150
151 icbmax = 2
152
2/2
✓ Branch 0 taken 240 times.
✓ Branch 1 taken 127097 times.
127337 DO il = 1, ncum
153 127337 icbmax = max(icbmax, icb(il))
154 END DO
155 ! if(prt.level.GE.20) print*,'cv3p1_param apres 200'
156
157 ! update sig and w0 below cloud base:
158
159
2/2
✓ Branch 0 taken 2482 times.
✓ Branch 1 taken 240 times.
2722 DO k = 1, icbmax
160
2/2
✓ Branch 0 taken 1314288 times.
✓ Branch 1 taken 2482 times.
1317010 DO il = 1, ncum
161
2/2
✓ Branch 0 taken 647548 times.
✓ Branch 1 taken 666740 times.
1316770 IF (k<=icb(il)) THEN
162 sig(il, k) = beta*sig(il, k) - 2.*alpha*buoy(il, icb(il))*buoy(il, &
163 647548 icb(il))
164 647548 sig(il, k) = amax1(sig(il,k), 0.0)
165 647548 w0(il, k) = beta*w0(il, k)
166 END IF
167 END DO
168 END DO
169
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 300'
170 ! -------------------------------------------------------------
171 ! -- Reset fractional areas of updrafts and w0 at initial time
172 ! -- and after 10 time steps of no convection
173 ! -------------------------------------------------------------
174
175
2/2
✓ Branch 0 taken 6240 times.
✓ Branch 1 taken 240 times.
6480 DO k = 1, nl - 1
176
2/2
✓ Branch 0 taken 3304522 times.
✓ Branch 1 taken 6240 times.
3311002 DO il = 1, ncum
177
3/4
✓ Branch 0 taken 3304522 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2148588 times.
✓ Branch 3 taken 1155934 times.
3310762 IF (sig(il,nd)<1.5 .OR. sig(il,nd)>12.0) THEN
178 2148588 sig(il, k) = 0.0
179 2148588 w0(il, k) = 0.0
180 END IF
181 END DO
182 END DO
183
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 400'
184
185 ! -------------------------------------------------------------
186 ! jyg1
187 ! -- Calculate adiabatic ascent top pressure (ptop)
188 ! -------------------------------------------------------------
189
190
191 ! c 1. Start at first level where precipitations form
192
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
193 127337 pzero(il) = plcl(il) - pbcrit
194 END DO
195
196 ! c 2. Add offset
197
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
198 127337 pzero(il) = pzero(il) - pbmxup
199 END DO
200
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
201 127337 ptop2old(il) = ptop2(il)
202 END DO
203
204
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
205 ! CR:c est quoi ce 300??
206 127337 p1(il) = pzero(il) - 300.
207 END DO
208
209 ! compute asupmax=abs(supmax) up to lnm+1
210
211
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
212 127097 ok(il) = .TRUE.
213 127337 nsupmax(il) = inb(il)
214 END DO
215
216
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
217
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
218
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (i>icb(il) .AND. i<=inb(il)) THEN
219
6/6
✓ Branch 0 taken 574178 times.
✓ Branch 1 taken 443429 times.
✓ Branch 2 taken 529976 times.
✓ Branch 3 taken 44202 times.
✓ Branch 4 taken 82667 times.
✓ Branch 5 taken 447309 times.
1017607 IF (p(il,i)<=pzero(il) .AND. supmax(il,i)<0 .AND. ok(il)) THEN
220 82667 nsupmax(il) = i
221 82667 ok(il) = .FALSE.
222 END IF ! end IF (P(i) ... )
223 END IF ! end IF (icb+1 le i le inb)
224 END DO
225 END DO
226
227
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 2.'
228
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
229
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
230 3438099 asupmax(il, i) = abs(supmax(il,i))
231 END DO
232 END DO
233
234
235
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
236 127097 asupmaxmin(il) = 10.
237 127097 pmin(il) = 100.
238 ! IM ??
239 127337 asupmax0(il) = 0.
240 END DO
241
242 ! c 3. Compute in which level is Pzero
243
244 ! IM bug i0 = 18
245
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
246 127337 i0(il) = nl
247 END DO
248
249
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
250
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
251
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (i>icb(il) .AND. i<=inb(il)) THEN
252
4/4
✓ Branch 0 taken 574178 times.
✓ Branch 1 taken 443429 times.
✓ Branch 2 taken 292791 times.
✓ Branch 3 taken 281387 times.
1017607 IF (p(il,i)<=pzero(il) .AND. p(il,i)>=p1(il)) THEN
253
3/4
✓ Branch 0 taken 292791 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 96676 times.
✓ Branch 3 taken 196115 times.
292791 IF (pzero(il)>p(il,i) .AND. pzero(il)<p(il,i-1)) THEN
254 96676 i0(il) = i
255 END IF
256 END IF
257 END IF
258 END DO
259 END DO
260
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 3.'
261
262 ! c 4. Compute asupmax at Pzero
263
264
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
265
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
266
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (i>icb(il) .AND. i<=inb(il)) THEN
267
4/4
✓ Branch 0 taken 574178 times.
✓ Branch 1 taken 443429 times.
✓ Branch 2 taken 292791 times.
✓ Branch 3 taken 281387 times.
1017607 IF (p(il,i)<=pzero(il) .AND. p(il,i)>=p1(il)) THEN
268 asupmax0(il) = ((pzero(il)-p(il,i0(il)-1))*asupmax(il,i0(il))-( &
269 pzero(il)-p(il,i0(il)))*asupmax(il,i0(il)-1))/(p(il,i0(il))-p(il, &
270 292791 i0(il)-1))
271 END IF
272 END IF
273 END DO
274 END DO
275
276
277
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
278
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
279
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 3431619 times.
3438099 IF (p(il,i)==pzero(il)) THEN
280 asupmax(i, il) = asupmax0(il)
281 END IF
282 END DO
283 END DO
284
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 4.'
285
286 ! c 5. Compute asupmaxmin, minimum of asupmax
287
288
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO i = 1, nl
289
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
290
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (i>icb(il) .AND. i<=inb(il)) THEN
291
4/4
✓ Branch 0 taken 574178 times.
✓ Branch 1 taken 443429 times.
✓ Branch 2 taken 292791 times.
✓ Branch 3 taken 281387 times.
1017607 IF (p(il,i)<=pzero(il) .AND. p(il,i)>=p1(il)) THEN
292
2/2
✓ Branch 0 taken 206816 times.
✓ Branch 1 taken 85975 times.
292791 IF (asupmax(il,i)<asupmaxmin(il)) THEN
293 206816 asupmaxmin(il) = asupmax(il, i)
294 206816 pmin(il) = p(il, i)
295 END IF
296 END IF
297 END IF
298 END DO
299 END DO
300
301
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
302 ! IM
303
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 127097 times.
127097 IF (prt_level>=20) THEN
304 PRINT *, 'cv3p1_closure il asupmax0 asupmaxmin', il, asupmax0(il), &
305 asupmaxmin(il), pzero(il), pmin(il)
306 END IF
307
2/2
✓ Branch 0 taken 39161 times.
✓ Branch 1 taken 87936 times.
127337 IF (asupmax0(il)<asupmaxmin(il)) THEN
308 39161 asupmaxmin(il) = asupmax0(il)
309 39161 pmin(il) = pzero(il)
310 END IF
311 END DO
312
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 5.'
313
314
315 ! Compute Supmax at Pzero
316
317
1/2
✓ Branch 0 taken 2160 times.
✗ Branch 1 not taken.
2160 DO i = 1, nl
318
2/2
✓ Branch 0 taken 1017805 times.
✓ Branch 1 taken 1920 times.
1019725 DO il = 1, ncum
319
4/4
✓ Branch 0 taken 375870 times.
✓ Branch 1 taken 641935 times.
✓ Branch 2 taken 354145 times.
✓ Branch 3 taken 21725 times.
1019725 IF (i>icb(il) .AND. i<=inb(il)) THEN
320
2/2
✓ Branch 0 taken 240 times.
✓ Branch 1 taken 353905 times.
354145 IF (p(il,i)<=pzero(il)) THEN
321 supmax0(il) = ((p(il,i)-pzero(il))*asupmax(il,i-1)-(p(il, &
322 240 i-1)-pzero(il))*asupmax(il,i))/(p(il,i)-p(il,i-1))
323 240 GO TO 425
324 END IF ! end IF (P(i) ... )
325 END IF ! end IF (icb+1 le i le inb)
326 END DO
327 END DO
328
329 425 CONTINUE
330
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 425.'
331
332 ! c 6. Calculate ptop2
333
334
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
335
2/2
✓ Branch 0 taken 74413 times.
✓ Branch 1 taken 52684 times.
127097 IF (asupmaxmin(il)<supcrit1) THEN
336 74413 ptop2(il) = pmin(il)
337 END IF
338
339
4/4
✓ Branch 0 taken 52684 times.
✓ Branch 1 taken 74413 times.
✓ Branch 2 taken 4780 times.
✓ Branch 3 taken 47904 times.
127097 IF (asupmaxmin(il)>supcrit1 .AND. asupmaxmin(il)<supcrit2) THEN
340 4780 ptop2(il) = ptop2old(il)
341 END IF
342
343
2/2
✓ Branch 0 taken 47904 times.
✓ Branch 1 taken 79193 times.
127337 IF (asupmaxmin(il)>supcrit2) THEN
344 47904 ptop2(il) = ph(il, inb(il))
345 END IF
346 END DO
347
348
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 6.'
349
350 ! c 7. Compute multiplying factor for adiabatic updraught mass flux
351
352
353
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (ok_inhib) THEN
354
355 DO i = 1, nl
356 DO il = 1, ncum
357 IF (i<=nl) THEN
358 coefmix(il, i) = (min(ptop2(il),ph(il,i))-ph(il,i))/(ph(il,i+1)-ph( &
359 il,i))
360 coefmix(il, i) = min(coefmix(il,i), 1.)
361 END IF
362 END DO
363 END DO
364
365
366 ELSE ! when inhibition is not taken into account, coefmix=1
367
368
369
370
2/2
✓ Branch 0 taken 240 times.
✓ Branch 1 taken 6480 times.
6720 DO i = 1, nl
371
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
372 6480 IF (i<=nl) THEN
373 3431619 coefmix(il, i) = 1.
374 END IF
375 END DO
376 END DO
377
378 END IF ! ok_inhib
379
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 7.'
380 ! -------------------------------------------------------------------
381 ! -------------------------------------------------------------------
382
383
384 ! jyg2
385
386 ! ==========================================================================
387
388
389 ! -------------------------------------------------------------
390 ! -- Calculate convective inhibition (CIN)
391 ! -------------------------------------------------------------
392
393 ! do i=1,nloc
394 ! print*,'avant cine p',pbase(i),plcl(i)
395 ! enddo
396 ! do j=1,nd
397 ! do i=1,nloc
398 ! print*,'avant cine t',tv(i),tvp(i)
399 ! enddo
400 ! enddo
401 CALL cv3_cine(nloc, ncum, nd, icb, inb, pbase, plcl, p, ph, tv, tvp, cina, &
402 240 cinb, plfc)
403
404
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
405 127337 cin(il) = cina(il) + cinb(il)
406 END DO
407
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres cv3_cine'
408 ! -------------------------------------------------------------
409 ! --Update buoyancies to account for Ale
410 ! -------------------------------------------------------------
411
412 CALL cv3_buoy(nloc, ncum, nd, icb, inb, pbase, plcl, p, ph, ale, cin, tv, &
413 240 tvp, buoy)
414
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres cv3_buoy'
415
416 ! -------------------------------------------------------------
417 ! -- Calculate convective available potential energy (cape),
418 ! -- vertical velocity (w), fractional area covered by
419 ! -- undilute updraft (sig), and updraft mass flux (m)
420 ! -------------------------------------------------------------
421
422
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
423 127337 cape(il) = 0.0
424 END DO
425
426 ! compute dtmin (minimum buoyancy between ICB and given level k):
427
428
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
429
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
430 3438099 dtmin(il, k) = 100.0
431 END DO
432 END DO
433
434
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
435
2/2
✓ Branch 0 taken 174960 times.
✓ Branch 1 taken 6480 times.
181680 DO j = minorig, nl
436
2/2
✓ Branch 0 taken 92653713 times.
✓ Branch 1 taken 174960 times.
92835153 DO il = 1, ncum
437
8/8
✓ Branch 0 taken 75169917 times.
✓ Branch 1 taken 17483796 times.
✓ Branch 2 taken 27475389 times.
✓ Branch 3 taken 47694528 times.
✓ Branch 4 taken 23518601 times.
✓ Branch 5 taken 3956788 times.
✓ Branch 6 taken 6033943 times.
✓ Branch 7 taken 17484658 times.
92653713 IF ((k>=(icb(il)+1)) .AND. (k<=inb(il)) .AND. (j>=icb(il)) .AND. (j<= &
438 174960 (k-1))) THEN
439 6033943 dtmin(il, k) = amin1(dtmin(il,k), buoy(il,j))
440 END IF
441 END DO
442 END DO
443 END DO
444
445 ! the interval on which cape is computed starts at pbase :
446
447
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
448
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
449
450
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF ((k>=(icb(il)+1)) .AND. (k<=inb(il))) THEN
451
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1017607 times.
1017607 IF (iflag_mix_adiab.eq.1) THEN
452 !CR:computation of cape from LCL: keep flag or to modify in all cases?
453 deltap = min(plcl(il), ph(il,k-1)) - min(plcl(il), ph(il,k))
454 ELSE
455 1017607 deltap = min(pbase(il), ph(il,k-1)) - min(pbase(il), ph(il,k))
456 ENDIF
457 1017607 cape(il) = cape(il) + rrd*buoy(il, k-1)*deltap/p(il, k-1)
458 1017607 cape(il) = amax1(0.0, cape(il))
459 1017607 sigold(il, k) = sig(il, k)
460
461
462 ! jyg Coefficient coefmix limits convection to levels where a
463 ! sufficient
464 ! fraction of mixed draughts are ascending.
465 1017607 siglim(il, k) = coefmix(il, k)*alpha1*dtmin(il, k)*abs(dtmin(il,k))
466 1017607 siglim(il, k) = amax1(siglim(il,k), 0.0)
467 1017607 siglim(il, k) = amin1(siglim(il,k), 0.01)
468 ! c fac=AMIN1(((dtcrit-dtmin(il,k))/dtcrit),1.0)
469 fac = 1.
470 1017607 wlim(il, k) = fac*sqrt(cape(il))
471 1017607 amu = siglim(il, k)*wlim(il, k)
472 1017607 rhodp = 0.007*p(il, k)*(ph(il,k)-ph(il,k+1))/tv(il, k)
473 1017607 mlim(il, k) = amu*rhodp
474 ! print*, 'siglim ', k,siglim(1,k)
475 END IF
476
477 END DO
478 END DO
479
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 600'
480
481
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
482 ! IM beg
483
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 127097 times.
127097 IF (prt_level>=20) THEN
484 PRINT *, 'cv3p1_closure il icb mlim ph ph+1 ph+2', il, icb(il), &
485 mlim(il, icb(il)+1), ph(il, icb(il)), ph(il, icb(il)+1), &
486 ph(il, icb(il)+2)
487 END IF
488
489
2/2
✓ Branch 0 taken 119370 times.
✓ Branch 1 taken 7727 times.
127337 IF (icb(il)+1<=inb(il)) THEN
490 ! IM end
491 mlim(il, icb(il)) = 0.5*mlim(il, icb(il)+1)*(ph(il,icb(il))-ph(il,icb( &
492 119370 il)+1))/(ph(il,icb(il)+1)-ph(il,icb(il)+2))
493 ! IM beg
494 END IF !(icb(il.le.inb(il))) then
495 ! IM end
496 END DO
497
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres 700'
498
499 ! jyg1
500 ! ------------------------------------------------------------------------
501 ! c Correct mass fluxes so that power used to overcome CIN does not
502 ! c exceed Power Available for Lifting (PAL).
503 ! ------------------------------------------------------------------------
504
505
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
506 127097 cbmflim(il) = 0.
507 127337 cbmf(il) = 0.
508 END DO
509
510 ! c 1. Compute cloud base mass flux of elementary system (Cbmf0=Cbmflim)
511
512
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
513
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
514 ! old IF (k .ge. icb(il) .and. k .le. inb(il)) THEN
515 ! IM IF (k .ge. icb(il)+1 .and. k .le. inb(il)) THEN
516 IF (k>=icb(il) .AND. k<=inb(il) & !cor jyg
517
6/6
✓ Branch 0 taken 2911168 times.
✓ Branch 1 taken 520451 times.
✓ Branch 2 taken 1144704 times.
✓ Branch 3 taken 1766464 times.
✓ Branch 4 taken 1136977 times.
✓ Branch 5 taken 7727 times.
3438099 .AND. icb(il)+1<=inb(il)) THEN !cor jyg
518 1136977 cbmflim(il) = cbmflim(il) + mlim(il, k)
519 END IF
520 END DO
521 END DO
522
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres cbmflim'
523
524 ! 1.5 Compute cloud base mass flux given by Alp closure (Cbmf1), maximum
525 ! allowed mass flux (Cbmfmax) and final target mass flux (Cbmf)
526 ! Cbmf is set to zero if Cbmflim (the mass flux of elementary cloud)
527 ! is exceedingly small.
528
529
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
530 127337 wb2(il) = sqrt(2.*max(ale(il)+cin(il),0.))
531 END DO
532
533
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
534
2/2
✓ Branch 0 taken 8259 times.
✓ Branch 1 taken 118838 times.
127337 IF (plfc(il)<100.) THEN
535 ! This is an irealistic value for plfc => no calculation of wbeff
536 8259 wbeff(il) = 100.1
537 ELSE
538 ! Calculate wbeff
539
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 118838 times.
118838 IF (NINT(flag_wb)==0) THEN
540 wbeff(il) = wbmax
541
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 118838 times.
118838 ELSE IF (NINT(flag_wb)==1) THEN
542 wbeff(il) = wbmax/(1.+500./(ph(il,1)-plfc(il)))
543
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 118838 times.
118838 ELSE IF (NINT(flag_wb)==2) THEN
544 wbeff(il) = wbmax*(0.01*(ph(il,1)-plfc(il)))**2
545 ELSE ! Option provisoire ou le iflag_wb/10 est considere comme une vitesse
546 118838 wbeff(il) = flag_wb*0.01+wbmax/(1.+500./(ph(il,1)-plfc(il)))
547 END IF
548 END IF
549 END DO
550
551 !CR:Compute k at plfc
552
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il=1,ncum
553 127337 klfc(il)=nl
554 ENDDO
555
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k=1,nl
556
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il=1,ncum
557
4/4
✓ Branch 0 taken 874742 times.
✓ Branch 1 taken 2556877 times.
✓ Branch 2 taken 118838 times.
✓ Branch 3 taken 755904 times.
3438099 if ((plfc(il).lt.ph(il,k)).and.(plfc(il).ge.ph(il,k+1))) then
558 118838 klfc(il)=k
559 endif
560 ENDDO
561 ENDDO
562 !RC
563
564
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
565 ! jyg Modification du coef de wb*wb pour conformite avec papier Wake
566 ! c cbmf1(il) = alp2(il)/(0.5*wb*wb-Cin(il))
567 127097 cbmf1(il) = alp2(il)/(2.*wbeff(il)*wbeff(il)-cin(il))
568 !CR: Add large-scale component to the mass-flux
569 !encore connu sous le nom "Experience du tube de dentifrice"
570
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 127097 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
127097 if ((coef_clos_ls.gt.0.).and.(plfc(il).gt.0.)) then
571 cbmf1(il) = cbmf1(il) - coef_clos_ls*min(0.,1./RG*omega(il,klfc(il)))
572 endif
573 !RC
574
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 127097 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
127097 IF (cbmf1(il)==0 .AND. alp2(il)/=0.) THEN
575 WRITE (lunout, *) 'cv3p1_closure cbmf1=0 and alp NE 0 il alp2 alp cin ' &
576 , il, alp2(il), alp(il), cin(il)
577 abort_message = ''
578 CALL abort_physic(modname, abort_message, 1)
579 END IF
580 127337 cbmfmax(il) = sigmax*wb2(il)*100.*p(il, icb(il))/(rrd*tv(il,icb(il)))
581 END DO
582
583
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
584
2/2
✓ Branch 0 taken 78205 times.
✓ Branch 1 taken 48892 times.
127337 IF (cbmflim(il)>1.E-6) THEN
585 ! ATTENTION TEST CR
586 ! if (cbmfmax(il).lt.1.e-12) then
587 78205 cbmf(il) = min(cbmf1(il), cbmfmax(il))
588 ! else
589 ! cbmf(il) = cbmf1(il)
590 ! endif
591 ! print*,'cbmf',cbmf1(il),cbmfmax(il)
592 END IF
593 END DO
594
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres cbmflim_testCR'
595
596 ! c 2. Compute coefficient and apply correction
597
598
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
599 127337 coef(il) = (cbmf(il)+1.E-10)/(cbmflim(il)+1.E-10)
600 END DO
601
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres coef_plantePLUS'
602
603
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
604
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
605
4/4
✓ Branch 0 taken 2784071 times.
✓ Branch 1 taken 647548 times.
✓ Branch 2 taken 1017607 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (k>=icb(il)+1 .AND. k<=inb(il)) THEN
606 amu = beta*sig(il, k)*w0(il, k) + (1.-beta)*coef(il)*siglim(il, k)* &
607 1017607 wlim(il, k)
608 w0(il, k) = wlim(il, k)
609 1017607 w0(il, k) = max(w0(il,k), 1.E-10)
610 1017607 sig(il, k) = amu/w0(il, k)
611 1017607 sig(il, k) = min(sig(il,k), 1.)
612 ! c amu = 0.5*(SIG(il,k)+sigold(il,k))*W0(il,k)
613 1017607 m(il, k) = amu*0.007*p(il, k)*(ph(il,k)-ph(il,k+1))/tv(il, k)
614 END IF
615 END DO
616 END DO
617 ! jyg2
618
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
619 127097 w0(il, icb(il)) = 0.5*w0(il, icb(il)+1)
620 m(il, icb(il)) = 0.5*m(il, icb(il)+1)*(ph(il,icb(il))-ph(il,icb(il)+1))/ &
621 127097 (ph(il,icb(il)+1)-ph(il,icb(il)+2))
622 127097 sig(il, icb(il)) = sig(il, icb(il)+1)
623 127337 sig(il, icb(il)-1) = sig(il, icb(il))
624 END DO
625
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres w0_sig_M'
626
627 !CR: new erosion of adiabatic ascent: modification of m
628 !computation of the sum of ascending fluxes
629
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (iflag_mix_adiab.eq.1) THEN
630
631 !Verification sum(me)=sum(m)
632 DO k = 1,nd !jyg: initialization up to nd
633 DO il = 1, ncum
634 md(il,k)=0.
635 med(il,k)=0.
636 ENDDO
637 ENDDO
638
639 DO k = nl,1,-1
640 DO il = 1, ncum
641 md(il,k)=md(il,k+1)+m(il,k+1)
642 ENDDO
643 ENDDO
644
645 DO k = nl,1,-1
646 DO il = 1, ncum
647 IF ((k>=(icb(il))) .AND. (k<=inb(il))) THEN
648 mad(il,k)=mad(il,k+1)+m(il,k+1)
649 ENDIF
650 ! print*,"mad",il,k,mad(il,k)
651 ENDDO
652 ENDDO
653
654 !CR: erosion of each adiabatic ascent during its ascent
655
656 !Computation of erosion coefficient beta_coef
657 DO k = 1, nl
658 DO il = 1, ncum
659 IF ((k>=(icb(il)+1)) .AND. (k<=inb(il)) .AND. (mlim(il,k).gt.0.)) THEN
660 ! print*,"beta_coef",il,k,icb(il),inb(il),buoy(il,k),tv(il,k),wlim(il,k),wlim(il,k+1)
661 beta_coef(il,k)=RG*coef_peel*buoy(il,k)/tv(il,k)/((wlim(il,k)+wlim(il,k+1))/2.)**2
662 ELSE
663 beta_coef(il,k)=0.
664 ENDIF
665 ENDDO
666 ENDDO
667
668 ! print*,"apres beta_coef"
669
670 DO k = 1, nl
671 DO il = 1, ncum
672
673 IF ((k>=(icb(il)+1)) .AND. (k<=inb(il))) THEN
674
675 ! print*,"dz",il,k,tv(il, k-1)
676 dz = (ph(il,k-1)-ph(il,k))/(p(il, k-1)/(rrd*tv(il, k-1))*RG)
677 betalim(il,k)=betalim(il,k-1)*exp(-1.*beta_coef(il,k-1)*dz)
678 ! betalim(il,k)=betalim(il,k-1)*exp(-RG*coef_peel*buoy(il,k-1)/tv(il,k-1)/5.**2*dz)
679 ! print*,"me",il,k,mlim(il,k),buoy(il,k),wlim(il,k),mad(il,k)
680 dz = (ph(il,k)-ph(il,k+1))/(p(il, k)/(rrd*tv(il, k))*RG)
681 ! me(il,k)=betalim(il,k)*(m(il,k)+RG*coef_peel*buoy(il,k)/tv(il,k)/((wlim(il,k)+wlim(il,k+1))/2.)**2*dz*mad(il,k))
682 me(il,k)=betalim(il,k)*(m(il,k)+beta_coef(il,k)*dz*mad(il,k))
683 ! print*,"B/w2",il,k,RG*coef_peel*buoy(il,k)/tv(il,k)/((wlim(il,k)+wlim(il,k+1))/2.)**2*dz
684
685 END IF
686
687 !Modification of m
688 m(il,k)=me(il,k)
689 END DO
690 END DO
691
692 ! DO il = 1, ncum
693 ! dz = (ph(il,icb(il))-ph(il,icb(il)+1))/(p(il, icb(il))/(rrd*tv(il, icb(il)))*RG)
694 ! m(il,icb(il))=m(il,icb(il))+RG*coef_peel*buoy(il,icb(il))/tv(il,icb(il)) &
695 ! /((wlim(il,icb(il))+wlim(il,icb(il)+1))/2.)**2*dz*mad(il,icb(il))
696 ! print*,"wlim(icb)",icb(il),wlim(il,icb(il)),m(il,icb(il))
697 ! ENDDO
698
699 !Verification sum(me)=sum(m)
700 DO k = nl,1,-1
701 DO il = 1, ncum
702 med(il,k)=med(il,k+1)+m(il,k+1)
703 ! print*,"somme(me),somme(m)",il,k,icb(il),med(il,k),md(il,k),me(il,k),m(il,k),wlim(il,k)
704 ENDDO
705 ENDDO
706
707
708 ENDIF !(iflag_mix_adiab)
709 !RC
710
711
712
713 ! c 3. Compute final cloud base mass flux and set iflag to 3 if
714 ! c cloud base mass flux is exceedingly small and is decreasing (i.e. if
715 ! c the final mass flux (cbmflast) is greater than the target mass flux
716 ! c (cbmf)).
717
718
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
719 127337 cbmflast(il) = 0.
720 END DO
721
722
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
723
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
724
4/4
✓ Branch 0 taken 2911168 times.
✓ Branch 1 taken 520451 times.
✓ Branch 2 taken 1144704 times.
✓ Branch 3 taken 1766464 times.
3438099 IF (k>=icb(il) .AND. k<=inb(il)) THEN
725 !IMpropo?? IF ((k.ge.(icb(il)+1)).and.(k.le.inb(il))) THEN
726 1144704 cbmflast(il) = cbmflast(il) + m(il, k)
727 END IF
728 END DO
729 END DO
730
731
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 DO il = 1, ncum
732
4/4
✓ Branch 0 taken 83289 times.
✓ Branch 1 taken 43808 times.
✓ Branch 2 taken 83288 times.
✓ Branch 3 taken 1 times.
127337 IF (cbmflast(il)<1.E-6 .AND. cbmflast(il)>=cbmf(il)) THEN
733 83288 iflag(il) = 3
734 END IF
735 END DO
736
737
2/2
✓ Branch 0 taken 6480 times.
✓ Branch 1 taken 240 times.
6720 DO k = 1, nl
738
2/2
✓ Branch 0 taken 3431619 times.
✓ Branch 1 taken 6480 times.
3438339 DO il = 1, ncum
739
2/2
✓ Branch 0 taken 2248776 times.
✓ Branch 1 taken 1182843 times.
3438099 IF (iflag(il)>=3) THEN
740 2248776 m(il, k) = 0.
741 2248776 sig(il, k) = 0.
742 2248776 w0(il, k) = 0.
743 END IF
744 END DO
745 END DO
746
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param apres iflag'
747
748 ! c 4. Introduce a correcting factor for coef, in order to obtain an
749 ! effective
750 ! c sigdz larger in the present case (using cv3p1_closure) than in the
751 ! old
752 ! c closure (using cv3_closure).
753 IF (1==0) THEN
754 DO il = 1, ncum
755 ! c coef(il) = 2.*coef(il)
756 coef(il) = 5.*coef(il)
757 END DO
758 ! version CVS du ..2008
759 ELSE
760
1/2
✓ Branch 0 taken 240 times.
✗ Branch 1 not taken.
240 IF (iflag_cvl_sigd==0) THEN
761 ! test pour verifier qu on fait la meme chose qu avant: sid constant
762
2/2
✓ Branch 0 taken 127097 times.
✓ Branch 1 taken 240 times.
127337 coef(1:ncum) = 1.
763 ELSE
764 coef(1:ncum) = min(2.*coef(1:ncum), 5.)
765 coef(1:ncum) = max(2.*coef(1:ncum), 0.2)
766 END IF
767 END IF
768
769
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 240 times.
240 IF (prt_level>=20) PRINT *, 'cv3p1_param FIN'
770 240 RETURN
771 END SUBROUTINE cv3p1_closure
772
773
774