Context Navigation

cv3p2_closure.F90 @ 2400

Last change on this file since 2400 was 2398, checked in by jyg, 9 years ago

Introduction of two new flags (ok_conv_stop
[Def=F], ok_intermittent [Def=F]) and one new
parameter (tau_stop [Def=15000]) in
conv_param.data:

. ok_conv_stop=T => convective mass fluxes are

set to zero if there is no trigger for a time
lapse longer than tau_stop.

. ok_intermittent=T => intermittent convection is

represented; the change concerns the upper
bound on the cloud base mass flux in
cv3p2_closure.F90.

The upper bound on Alp in physiq.F90 should

also be changed: this is still to be done.

File size: 23.0 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format