Context Navigation

cv3p_mixing.f90 @ 5779

Last change on this file since 5779 was 5346, checked in by fhourdin, 12 months ago

Debut de replaysation de la convection profonde.

Regroupement de cvparam, cv3param et cvthermo (récemment
passés de statut de .h à module, dans un unique module
lmdz_cv_ini.f90

Property copyright set to
Name of program: LMDZ
Creation date: 1984
Version: LMDZ5
License: CeCILL version 2
Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
See the license file in the root directory
Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 25.3 KB

Line
1	SUBROUTINE cv3p_mixing(nloc, ncum, nd, na, ntra, icb, nk, inb, &
2	ph, t, rr, rs, u, v, tra, h, lv, lf, frac, qta, &
3	unk, vnk, hp, tv, tvp, ep, clw, sig, &
4	Ment, Qent, hent, uent, vent, nent, &
5	Sigij, elij, supmax, Ments, Qents, traent)
6	! **************************************************************
7	! *
8	! CV3P_MIXING : compute mixed draught properties and, *
9	! within a scaling factor, mixed draught *
10	! mass fluxes. *
11	! written by : VTJ Philips,JY Grandpeix, 21/05/2003, 09.14.15*
12	! modified by : *
13	! **************************************************************
14
15	USE yomcst2_mod_h
16	USE lmdz_cv_ini, ONLY : cpd,cpv,minorig,nl,rrv
17	USE cvflag_mod_h
18	USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
19	USE ioipsl_getin_p_mod, ONLY: getin_p
20	USE add_phys_tend_mod, ONLY: fl_cor_ebil
21
22	IMPLICIT NONE
23
24
25	!inputs:
26	INTEGER, INTENT (IN) :: ncum, nd, na
27	INTEGER, INTENT (IN) :: ntra, nloc
28	INTEGER, DIMENSION (nloc), INTENT (IN) :: icb, inb, nk
29	REAL, DIMENSION (nloc, nd), INTENT (IN) :: sig
30	REAL, DIMENSION (nloc), INTENT (IN) :: unk, vnk
31	REAL, DIMENSION (nloc, nd), INTENT (IN) :: qta
32	REAL, DIMENSION (nloc, nd+1), INTENT (IN) :: ph
33	REAL, DIMENSION (nloc, nd), INTENT (IN) :: t, rr, rs
34	REAL, DIMENSION (nloc, nd), INTENT (IN) :: u, v
35	REAL, DIMENSION (nloc, nd, ntra), INTENT (IN) :: tra ! input of convect3
36	REAL, DIMENSION (nloc, na), INTENT (IN) :: lv
37	REAL, DIMENSION (nloc, na), INTENT (IN) :: lf
38	REAL, DIMENSION (nloc, na), INTENT (IN) :: frac !ice fraction in condensate
39	REAL, DIMENSION (nloc, na), INTENT (IN) :: h !liquid water static energy of environment
40	REAL, DIMENSION (nloc, na), INTENT (IN) :: hp !liquid water static energy of air shed from adiab. asc.
41	REAL, DIMENSION (nloc, na), INTENT (IN) :: tv, tvp
42	REAL, DIMENSION (nloc, na), INTENT (IN) :: ep, clw
43
44	!outputs:
45	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: Ment, Qent
46	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: uent, vent
47	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: Sigij, elij
48	REAL, DIMENSION (nloc, na), INTENT (OUT) :: supmax ! Highest mixing fraction of mixed
49	! updraughts with the sign of (h-hp)
50	REAL, DIMENSION (nloc, nd, nd, ntra), INTENT (OUT) :: traent
51	REAL, DIMENSION (nloc, nd, nd), INTENT (OUT) :: Ments, Qents
52	REAL, DIMENSION (nloc, nd, nd), INTENT (OUT) :: hent
53	INTEGER, DIMENSION (nloc, nd), INTENT (OUT) :: nent
54
55	!local variables:
56	INTEGER i, j, k, il, im, jm
57	INTEGER num1, num2
58	REAL :: rti, bf2, anum, denom, dei, altem, cwat, stemp
59	REAL :: alt, delp, delm
60	REAL, DIMENSION (nloc) :: Qmixmax, Rmixmax, sqmrmax
61	REAL, DIMENSION (nloc) :: Qmixmin, Rmixmin, sqmrmin
62	REAL, DIMENSION (nloc) :: signhpmh
63	REAL, DIMENSION (nloc) :: Sx
64	REAL :: Scrit2
65	REAL, DIMENSION (nloc) :: Smid, Sjmin, Sjmax
66	REAL, DIMENSION (nloc) :: Sbef, sup, smin
67	REAL, DIMENSION (nloc) :: ASij, ASij_inv, smax, Scrit
68	REAL, DIMENSION (nloc, nd, nd) :: Sij
69	REAL, DIMENSION (nloc, nd) :: csum
70	REAL :: awat
71	REAL :: cpm !Mixed draught heat capacity
72	REAL :: Tm !Mixed draught temperature
73	LOGICAL, DIMENSION (nloc) :: lwork
74
75	REAL amxupcrit, df, ff
76	INTEGER nstep
77
78	INTEGER,SAVE :: igout=1
79	!$OMP THREADPRIVATE(igout)
80
81	! -- Mixing probability distribution functions
82
83	REAL Qcoef1, Qcoef2, QFF, QFFF, Qmix, Rmix, Qmix1, Rmix1, Qmix2, Rmix2, F
84
85	Qcoef1(F) = tanh(F/gammas)
86	Qcoef2(F) = (tanh(F/gammas)+gammas*log(cosh((1.-F)/gammas)/cosh(F/gammas)))
87	QFF(F) = max(min(F,1.), 0.)
88	QFFf(F) = min(QFF(F), scut)
89	Qmix1(F) = (tanh((QFF(F)-Fmax)/gammas)+Qcoef1max)/Qcoef2max
90	Rmix1(F) = (gammaslog(cosh((QFF(F)-Fmax)/gammas))+QFF(F)Qcoef1max)/Qcoef2max
91	Qmix2(F) = -log(1.-QFFf(F))/scut
92	Rmix2(F) = (QFFf(F)+(1.-QFF(F))*log(1.-QFFf(F)))/scut
93	Qmix(F) = qqa1Qmix1(F) + qqa2Qmix2(F)
94	Rmix(F) = qqa1Rmix1(F) + qqa2Rmix2(F)
95
96	INTEGER, SAVE :: ifrst
97	DATA ifrst/0/
98	!$OMP THREADPRIVATE(ifrst)
99
100
101	! =====================================================================
102	! --- INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS
103	! =====================================================================
104
105	! -- Initialize mixing PDF coefficients
106	IF (ifrst==0) THEN
107	ifrst = 1
108	Qcoef1max = Qcoef1(Fmax)
109	Qcoef2max = Qcoef2(Fmax)
110	!<jyg
111	print*, 'fmax, gammas, qqa1, qqa2, Qcoef1max, Qcoef2max ', &
112	fmax, gammas, qqa1, qqa2, Qcoef1max, Qcoef2max
113	!>jyg
114	!
115	END IF
116
117
118	! ori do 360 i=1,ncum*nlp
119	DO j = 1, nl
120	DO i = 1, ncum
121	nent(i, j) = 0
122	! in convect3, m is computed in cv3_closure
123	! ori m(i,1)=0.0
124	END DO
125	END DO
126
127	! ori do 400 k=1,nlp
128	! ori do 390 j=1,nlp
129	DO j = 1, nl
130	DO k = 1, nl
131	DO i = 1, ncum
132	Qent(i, k, j) = rr(i, j)
133	uent(i, k, j) = u(i, j)
134	vent(i, k, j) = v(i, j)
135	elij(i, k, j) = 0.0
136	hent(i, k, j) = 0.0
137	!AC! Ment(i,k,j)=0.0
138	!AC! Sij(i,k,j)=0.0
139	END DO
140	END DO
141	END DO
142
143	!AC!
144	Ment(1:ncum, 1:nd, 1:nd) = 0.0
145	Sij(1:ncum, 1:nd, 1:nd) = 0.0
146	!AC!
147	!ym
148	Sigij(1:ncum, 1:nd, 1:nd) = 0.0
149	!ym
150
151	!jyg! DO k = 1, ntra
152	!jyg! DO j = 1, nd ! instead nlp
153	!jyg! DO i = 1, nd ! instead nlp
154	!jyg! DO il = 1, ncum
155	!jyg! traent(il, i, j, k) = tra(il, j, k)
156	!jyg! END DO
157	!jyg! END DO
158	!jyg! END DO
159	!jyg! END DO
160
161	! =====================================================================
162	! --- CALCULATE ENTRAINED AIR MASS FLUX (Ment), TOTAL WATER MIXING
163	! --- RATIO (QENT), TOTAL CONDENSED WATER (elij), AND MIXING
164	! --- FRACTION (Sij)
165	! =====================================================================
166
167	DO i = minorig + 1, nl
168
169	IF (ok_entrain) THEN
170	DO j = minorig, nl
171	DO il = 1, ncum
172	IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (j>=(icb(il)-1)) &
173	.AND. (j<=inb(il))) THEN
174
175	!! rti = qnk(il) - ep(il, i)*clw(il, i)
176	rti = qta(il,i-1) - ep(il, i)*clw(il, i)
177	bf2 = 1. + lv(il, j)lv(il, j)rs(il, j)/(rrvt(il,j)t(il,j)*cpd)
178	!jyg(from aj)<
179	IF (cvflag_ice) THEN
180	! print*,cvflag_ice,'cvflag_ice dans do 700'
181	IF (t(il,j)<=263.15) THEN
182	bf2 = 1. + (lf(il,j)+lv(il,j))(lv(il,j)+frac(il,j) &
183	lf(il,j))rs(il, j)/(rrvt(il,j)t(il,j)cpd)
184	END IF
185	END IF
186	!>jyg
187	anum = h(il, j) - hp(il, i) + (cpv-cpd)t(il, j)(rti-rr(il,j))
188	denom = h(il, i) - hp(il, i) + (cpd-cpv)(rr(il,i)-rti)t(il, j)
189	dei = denom
190	IF (abs(dei)<0.01) dei = 0.01
191	Sij(il, i, j) = anum/dei
192	Sij(il, i, i) = 1.0
193	altem = Sij(il, i, j)rr(il, i) + (1.-Sij(il,i,j))rti - rs(il, j)
194	altem = altem/bf2
195	cwat = clw(il, j)*(1.-ep(il,j))
196	stemp = Sij(il, i, j)
197	IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN
198	!jyg(from aj)<
199	IF (cvflag_ice) THEN
200	anum = anum - (lv(il,j)+frac(il,j)lf(il,j))(rti-rs(il,j)-cwat*bf2)
201	denom = denom + (lv(il,j)+frac(il,j)lf(il,j))(rr(il,i)-rti)
202	ELSE
203	anum = anum - lv(il, j)(rti-rs(il,j)-cwatbf2)
204	denom = denom + lv(il, j)*(rr(il,i)-rti)
205	END IF
206	!>jyg
207	IF (abs(denom)<0.01) denom = 0.01
208	Sij(il, i, j) = anum/denom
209	altem = Sij(il, i, j)rr(il, i) + (1.-Sij(il,i,j))rti - rs(il, j)
210	altem = altem - (bf2-1.)*cwat
211	END IF
212	IF (Sij(il,i,j)>0.0) THEN
213	!!! Ment(il,i,j)=m(il,i)
214	Ment(il, i, j) = 1.
215	elij(il, i, j) = altem
216	elij(il, i, j) = amax1(0.0, elij(il,i,j))
217	nent(il, i) = nent(il, i) + 1
218	END IF
219
220	Sij(il, i, j) = amax1(0.0, Sij(il,i,j))
221	Sij(il, i, j) = amin1(1.0, Sij(il,i,j))
222	ELSE IF (j > i) THEN
223	IF (prt_level >= 10) THEN
224	print *,'cv3p_mixing i, j, Sij given by the no-precip eq. ', i, j, Sij(il,i,j)
225	ENDIF
226	END IF ! new
227	END DO
228	END DO
229	ELSE ! (ok_entrain)
230	DO il = 1,ncum
231	nent(il,i) = 0
232	ENDDO
233	ENDIF ! (ok_entrain)
234
235	!jygdebug<
236	IF (prt_level >= 10) THEN
237	print *,'cv3p_mixing i, nent(i), icb, inb ',i, nent(igout,i), icb(igout), inb(igout)
238	IF (nent(igout,i) .gt. 0) THEN
239	print *,'i,(j,Sij(i,j),j=icb-1,inb) ',i,(j,Sij(igout,i,j),j=icb(igout)-1,inb(igout))
240	ENDIF
241	ENDIF
242	!>jygdebug
243
244	! * if no air can entrain at level i assume that updraft detrains *
245	! * at that level and calculate detrained air flux and properties *
246
247
248	! @ do 170 i=icb(il),inb(il)
249
250	DO il = 1, ncum
251	IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (nent(il,i)==0)) THEN
252	! @ if(nent(il,i).eq.0)then
253	!!! Ment(il,i,i)=m(il,i)
254	Ment(il, i, i) = 1.
255	!! Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
256	Qent(il, i, i) = qta(il,i-1) - ep(il, i)*clw(il, i)
257	uent(il, i, i) = unk(il)
258	vent(il, i, i) = vnk(il)
259	IF (fl_cor_ebil .GE. 2) THEN
260	hent(il, i, i) = hp(il,i)
261	ENDIF
262	elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
263	Sij(il, i, i) = 0.0
264	END IF
265	END DO
266	END DO ! i = minorig + 1, nl
267
268	!jyg! DO j = 1, ntra
269	!jyg! DO i = minorig + 1, nl
270	!jyg! DO il = 1, ncum
271	!jyg! IF (i>=icb(il) .AND. i<=inb(il) .AND. nent(il,i)==0) THEN
272	!jyg! traent(il, i, i, j) = tra(il, nk(il), j)
273	!jyg! END IF
274	!jyg! END DO
275	!jyg! END DO
276	!jyg! END DO
277
278	DO j = minorig, nl
279	DO i = minorig, nl
280	DO il = 1, ncum
281	IF ((j>=(icb(il)-1)) .AND. (j<=inb(il)) .AND. &
282	(i>=icb(il)) .AND. (i<=inb(il))) THEN
283	Sigij(il, i, j) = Sij(il, i, j)
284	END IF
285	END DO
286	END DO
287	END DO
288	! @ enddo
289
290	! @170 continue
291
292	! =====================================================================
293	! --- NORMALIZE ENTRAINED AIR MASS FLUXES
294	! --- TO REPRESENT EQUAL PROBABILITIES OF MIXING
295	! =====================================================================
296
297	CALL zilch(csum, nloc*nd)
298
299	DO il = 1, ncum
300	lwork(il) = .FALSE.
301	END DO
302
303	! ---------------------------------------------------------------
304	DO i = minorig + 1, nl !Loop on origin level "i"
305	! ---------------------------------------------------------------
306
307	num1 = 0
308	DO il = 1, ncum
309	IF (i>=icb(il) .AND. i<=inb(il)) num1 = num1 + 1
310	END DO
311	IF (num1<=0) GO TO 789
312
313
314	!JYG1 Find maximum of SIJ for J>I, if any.
315
316	Sx(:) = 0.
317
318	DO il = 1, ncum
319	IF (i>=icb(il) .AND. i<=inb(il)) THEN
320	signhpmh(il) = sign(1., hp(il,i)-h(il,i))
321	Sbef(il) = max(0., signhpmh(il))
322	END IF
323	END DO
324
325	DO j = i + 1, nl
326	DO il = 1, ncum
327	IF (i>=icb(il) .AND. i<=inb(il) .AND. j<=inb(il)) THEN
328	IF (Sbef(il)<Sij(il,i,j)) THEN
329	Sx(il) = max(Sij(il,i,j), Sx(il))
330	END IF
331	Sbef(il) = Sij(il, i, j)
332	END IF
333	END DO
334	END DO
335
336
337	DO il = 1, ncum
338	IF (i>=icb(il) .AND. i<=inb(il)) THEN
339	lwork(il) = (nent(il,i)/=0)
340	!! rti = qnk(il) - ep(il, i)*clw(il, i)
341	rti = qta(il,i-1) - ep(il, i)*clw(il, i)
342	!jyg<
343	IF (cvflag_ice) THEN
344
345	anum = h(il, i) - hp(il, i) - (lv(il,i)+frac(il,i)lf(il,i)) &
346	(rti-rs(il,i)) + (cpv-cpd)t(il, i)(rti-rr(il,i))
347	denom = h(il, i) - hp(il, i) + (lv(il,i)+frac(il,i)lf(il,i)) &
348	(rr(il,i)-rti) + (cpd-cpv)t(il, i)(rr(il,i)-rti)
349	ELSE
350
351	anum = h(il, i) - hp(il, i) - lv(il, i)*(rti-rs(il,i)) + &
352	(cpv-cpd)t(il, i)(rti-rr(il,i))
353	denom = h(il, i) - hp(il, i) + lv(il, i)*(rr(il,i)-rti) + &
354	(cpd-cpv)t(il, i)(rr(il,i)-rti)
355	END IF
356	!>jyg
357	IF (abs(denom)<0.01) denom = 0.01
358	Scrit(il) = min(anum/denom, 1.)
359	alt = rti - rs(il, i) + Scrit(il)*(rr(il,i)-rti)
360
361	!JYG1 Find new critical value Scrit2
362	! such that : Sij > Scrit2 => mixed draught will detrain at J<I
363	! Sij < Scrit2 => mixed draught will detrain at J>I
364
365	Scrit2 = min(Scrit(il), Sx(il))*max(0., -signhpmh(il)) + &
366	Scrit(il)*max(0., signhpmh(il))
367
368	Scrit(il) = Scrit2
369
370	!JYG Correction pour la nouvelle logique; la correction pour ALT
371	! est un peu au hazard
372	IF (Scrit(il)<=0.0) Scrit(il) = 0.0
373	IF (alt<=0.0) Scrit(il) = 1.0
374
375	smax(il) = 0.0
376	ASij(il) = 0.0
377	sup(il) = 0. ! upper S-value reached by descending draughts
378	END IF
379	END DO
380
381	! ---------------------------------------------------------------
382	DO j = minorig, nl !Loop on destination level "j"
383	! ---------------------------------------------------------------
384
385	num2 = 0
386	DO il = 1, ncum
387	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
388	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
389	lwork(il)) num2 = num2 + 1
390	END DO
391	IF (num2<=0) GO TO 175
392
393	! -----------------------------------------------
394	IF (j>i) THEN
395	! -----------------------------------------------
396	DO il = 1, ncum
397	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
398	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
399	lwork(il)) THEN
400	IF (Sij(il,i,j)>0.0) THEN
401	Smid(il) = min(Sij(il,i,j), Scrit(il))
402	Sjmax(il) = Smid(il)
403	Sjmin(il) = Smid(il)
404	IF (Smid(il)<smin(il) .AND. Sij(il,i,j+1)<Smid(il)) THEN
405	smin(il) = Smid(il)
406	Sjmax(il) = min((Sij(il,i,j+1)+Sij(il,i,j))/2., Sij(il,i,j), Scrit(il))
407	Sjmin(il) = max((Sbef(il)+Sij(il,i,j))/2., Sij(il,i,j))
408	Sjmin(il) = min(Sjmin(il), Scrit(il))
409	Sbef(il) = Sij(il, i, j)
410	END IF
411	END IF
412	END IF
413	END DO
414	! -----------------------------------------------
415	ELSE IF (j==i) THEN
416	! -----------------------------------------------
417	DO il = 1, ncum
418	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
419	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
420	lwork(il)) THEN
421	IF (Sij(il,i,j)>0.0) THEN
422	Smid(il) = 1.
423	Sjmin(il) = max((Sij(il,i,j-1)+Smid(il))/2., Scrit(il))*max(0., -signhpmh(il)) + &
424	min((Sij(il,i,j+1)+Smid(il))/2., Scrit(il))*max(0., signhpmh(il))
425	Sjmin(il) = max(Sjmin(il), sup(il))
426	Sjmax(il) = 1.
427
428	! - preparation des variables Scrit, Smin et Sbef pour la partie j>i
429	Scrit(il) = min(Sjmin(il), Sjmax(il), Scrit(il))
430
431	smin(il) = 1.
432	Sbef(il) = max(0., signhpmh(il))
433	supmax(il, i) = sign(Scrit(il), -signhpmh(il))
434	END IF
435	END IF
436	END DO
437	! -----------------------------------------------
438	ELSE IF (j<i) THEN
439	! -----------------------------------------------
440	DO il = 1, ncum
441	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
442	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
443	lwork(il)) THEN
444	IF (Sij(il,i,j)>0.0) THEN
445	Smid(il) = max(Sij(il,i,j), Scrit(il))
446	Sjmax(il) = Smid(il)
447	Sjmin(il) = Smid(il)
448	IF (Smid(il)>smax(il) .AND. Sij(il,i,j+1)>Smid(il)) THEN
449	smax(il) = Smid(il)
450	Sjmax(il) = max((Sij(il,i,j+1)+Sij(il,i,j))/2., Sij(il,i,j))
451	Sjmax(il) = max(Sjmax(il), Scrit(il))
452	Sjmin(il) = min((Sbef(il)+Sij(il,i,j))/2., Sij(il,i,j))
453	Sjmin(il) = max(Sjmin(il), Scrit(il))
454	Sbef(il) = Sij(il, i, j)
455	END IF
456	IF (abs(Sjmin(il)-Sjmax(il))>1.E-10) &
457	sup(il) = max(Sjmin(il), Sjmax(il), sup(il))
458	END IF
459	END IF
460	END DO
461	! -----------------------------------------------
462	END IF
463	! -----------------------------------------------
464
465
466	DO il = 1, ncum
467	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
468	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
469	lwork(il)) THEN
470	IF (Sij(il,i,j)>0.0) THEN
471	!! rti = qnk(il) - ep(il, i)*clw(il, i)
472	rti = qta(il,i-1) - ep(il, i)*clw(il, i)
473	Qmixmax(il) = Qmix(Sjmax(il))
474	Qmixmin(il) = Qmix(Sjmin(il))
475	Rmixmax(il) = Rmix(Sjmax(il))
476	Rmixmin(il) = Rmix(Sjmin(il))
477	sqmrmax(il) = Sjmax(il)*Qmix(Sjmax(il)) - Rmix(Sjmax(il))
478	sqmrmin(il) = Sjmin(il)*Qmix(Sjmin(il)) - Rmix(Sjmin(il))
479
480	Ment(il, i, j) = abs(Qmixmax(il)-Qmixmin(il))*Ment(il, i, j)
481
482	! Sigij(i,j) is the 'true' mixing fraction of mixture Ment(i,j)
483	IF (abs(Qmixmax(il)-Qmixmin(il))>1.E-10) THEN
484	Sigij(il, i, j) = (sqmrmax(il)-sqmrmin(il))/(Qmixmax(il)-Qmixmin(il))
485	ELSE
486	Sigij(il, i, j) = 0.
487	END IF
488
489	! -- Compute Qent, uent, vent according to the true mixing fraction
490	Qent(il, i, j) = (1.-Sigij(il,i,j))rti + Sigij(il, i, j)rr(il, i)
491	uent(il, i, j) = (1.-Sigij(il,i,j))unk(il) + Sigij(il, i, j)u(il, i)
492	vent(il, i, j) = (1.-Sigij(il,i,j))vnk(il) + Sigij(il, i, j)v(il, i)
493
494	! -- Compute liquid water static energy of mixed draughts
495	! IF (j .GT. i) THEN
496	! awat=elij(il,i,j)-(1.-ep(il,j))*clw(il,j)
497	! awat=amax1(awat,0.0)
498	! ELSE
499	! awat = 0.
500	! ENDIF
501	! Hent(il,i,j) = (1.-Sigij(il,i,j))*HP(il,i)
502	! : + Sigij(il,i,j)*H(il,i)
503	! : + (LV(il,j)+(cpd-cpv)t(il,j))awat
504	!IM 301008 beg
505	hent(il, i, j) = (1.-Sigij(il,i,j))hp(il, i) + Sigij(il, i, j)h(il, i)
506
507	!jyg<
508	! elij(il, i, j) = Qent(il, i, j) - rs(il, j)
509	! elij(il, i, j) = elij(il, i, j) + &
510	! ((h(il,j)-hent(il,i,j))rs(il,j)lv(il,j) / &
511	! ((cpd(1.-Qent(il,i,j))+Qent(il,i,j)cpv)rrvt(il,j)*t(il,j)))
512	! elij(il, i, j) = elij(il, i, j) / &
513	! (1.+lv(il,j)lv(il,j)rs(il,j) / &
514	! ((cpd(1.-Qent(il,i,j))+Qent(il,i,j)cpv)rrvt(il,j)*t(il,j)))
515	!
516	! Computation of condensate amount Elij, taking into account the ice fraction frac
517	! Warning : the same saturation humidity rs is used over both liquid water and ice; this
518	! should be corrected.
519	!
520	! Heat capacity of mixed draught
521	cpm = cpd+Qent(il,i,j)*(cpv-cpd)
522	!
523	IF (cvflag_ice .and. frac(il,j) .gt. 0.) THEN
524	elij(il, i, j) = Qent(il, i, j) - rs(il, j)
525	elij(il, i, j) = elij(il, i, j) + &
526	(h(il,j)-hent(il,i,j)+(cpv-cpd)(Qent(il,i,j)-rr(il,j))t(il,j))* &
527	rs(il,j)lv(il,j) / (cpmrrvt(il,j)t(il,j))
528	elij(il, i, j) = elij(il, i, j) / &
529	(1.+(lv(il,j)+frac(il,j)lf(il,j))lv(il,j)*rs(il,j) / &
530	(cpmrrvt(il,j)*t(il,j)))
531	ELSE
532	elij(il, i, j) = Qent(il, i, j) - rs(il, j)
533	elij(il, i, j) = elij(il, i, j) + &
534	(h(il,j)-hent(il,i,j)+(cpv-cpd)(Qent(il,i,j)-rr(il,j))t(il,j))* &
535	rs(il,j)lv(il,j) / (cpmrrvt(il,j)t(il,j))
536	elij(il, i, j) = elij(il, i, j) / &
537	(1.+lv(il,j)lv(il,j)rs(il,j) / &
538	(cpmrrvt(il,j)*t(il,j)))
539	ENDIF
540	!>jyg
541	elij(il, i, j) = max(elij(il,i,j), 0.)
542
543	elij(il, i, j) = min(elij(il,i,j), Qent(il,i,j))
544
545	IF (j>i) THEN
546	awat = elij(il, i, j) - (1.-ep(il,j))*clw(il, j)
547	awat = amax1(awat, 0.0)
548	ELSE
549	awat = 0.
550	END IF
551
552	! print ,h(il,j)-hent(il,i,j),LV(il,j)rs(il,j)/(cpdrrvt(il,j)*
553	! : t(il,j))
554
555	!jyg<
556	! hent(il, i, j) = hent(il, i, j) + (lv(il,j)+(cpd-cpv)t(il,j))awat
557	! Mixed draught temperature at level j
558	IF (cvflag_ice .and. frac(il,j) .gt. 0.) THEN
559	Tm = t(il,j) + (Qent(il,i,j)-elij(il,i,j)-rs(il,j))rrvt(il,j)t(il,j)/(lv(il,j)rs(il,j))
560	hent(il, i, j) = hent(il, i, j) + (lv(il,j)+frac(il,j)lf(il,j)+(cpd-cpv)Tm)*awat
561	ELSE
562	Tm = t(il,j) + (Qent(il,i,j)-elij(il,i,j)-rs(il,j))rrvt(il,j)t(il,j)/(lv(il,j)rs(il,j))
563	hent(il, i, j) = hent(il, i, j) + (lv(il,j)+(cpd-cpv)Tm)awat
564	ENDIF
565	!>jyg
566
567	!IM 301008 end
568
569	! print *,'mix : i,j,hent(il,i,j),Sigij(il,i,j) ',
570	! : i,j,hent(il,i,j),Sigij(il,i,j)
571
572	! -- ASij is the integral of P(F) over the relevant F interval
573	ASij(il) = ASij(il) + abs(Qmixmax(il)*(1.-Sjmax(il))+Rmixmax(il) - &
574	Qmixmin(il)*(1.-Sjmin(il))-Rmixmin(il))
575
576	END IF
577	END IF
578	END DO
579	!jyg! DO k = 1, ntra
580	!jyg! DO il = 1, ncum
581	!jyg! IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. &
582	!jyg! (j>=(icb(il)-1)) .AND. (j<=inb(il)) .AND. &
583	!jyg! lwork(il)) THEN
584	!jyg! IF (Sij(il,i,j)>0.0) THEN
585	!jyg! traent(il, i, j, k) = Sigij(il, i, j)*tra(il, i, k) + &
586	!jyg! (1.-Sigij(il,i,j))*tra(il, nk(il), k)
587	!jyg! END IF
588	!jyg! END IF
589	!jyg! END DO
590	!jyg! END DO
591
592	! -- If I=J (detrainement and entrainement at the same level), then only the
593	! -- adiabatic ascent part of the mixture is considered
594	IF (i==j) THEN
595	DO il = 1, ncum
596	IF (i>=icb(il) .AND. i<=inb(il) .AND. &
597	j>=(icb(il)-1) .AND. j<=inb(il) .AND. &
598	lwork(il)) THEN
599	IF (Sij(il,i,j)>0.0) THEN
600	!! rti = qnk(il) - ep(il, i)*clw(il, i)
601	rti = qta(il,i-1) - ep(il, i)*clw(il, i)
602	!!! Ment(il,i,i) = m(il,i)abs(Qmixmax(il)(1.-Sjmax(il))
603	Ment(il, i, i) = abs(Qmixmax(il)*(1.-Sjmax(il))+Rmixmax(il) - &
604	Qmixmin(il)*(1.-Sjmin(il))-Rmixmin(il))
605	Qent(il, i, i) = rti
606	uent(il, i, i) = unk(il)
607	vent(il, i, i) = vnk(il)
608	hent(il, i, i) = hp(il, i)
609	elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
610	Sigij(il, i, i) = 0.
611	END IF
612	END IF
613	END DO
614	!jyg! DO k = 1, ntra
615	!jyg! DO il = 1, ncum
616	!jyg! IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. &
617	!jyg! (j>=(icb(il)-1)) .AND. (j<=inb(il)) .AND. &
618	!jyg! lwork(il)) THEN
619	!jyg! IF (Sij(il,i,j)>0.0) THEN
620	!jyg! traent(il, i, i, k) = tra(il, nk(il), k)
621	!jyg! END IF
622	!jyg! END IF
623	!jyg! END DO
624	!jyg! END DO
625
626	END IF
627
628	! ---------------------------------------------------------------
629	175 END DO ! End loop on destination level "j"
630	! ---------------------------------------------------------------
631
632	DO il = 1, ncum
633	IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il)) THEN
634	ASij(il) = amax1(1.0E-16, ASij(il))
635	!jyg+lluis<
636	!! ASij(il) = 1.0/ASij(il)
637	ASij_inv(il) = 1.0/ASij(il)
638	! IF the F-interval spanned by possible mixtures is less than 0.01, no mixing occurs
639	IF (ASij_inv(il) > 100.) ASij_inv(il) = 0.
640	!>jyg+lluis
641	csum(il, i) = 0.0
642	END IF
643	END DO
644
645	DO j = minorig, nl
646	DO il = 1, ncum
647	IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il) .AND. &
648	j>=(icb(il)-1) .AND. j<=inb(il)) THEN
649	!jyg Ment(il, i, j) = Ment(il, i, j)*ASij(il)
650	Ment(il, i, j) = Ment(il, i, j)*ASij_inv(il)
651	END IF
652	END DO
653	END DO
654
655	DO j = minorig, nl
656	DO il = 1, ncum
657	IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il) .AND. &
658	j>=(icb(il)-1) .AND. j<=inb(il)) THEN
659	csum(il, i) = csum(il, i) + Ment(il, i, j)
660	END IF
661	END DO
662	END DO
663
664	DO il = 1, ncum
665	IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il) .AND. csum(il,i)<1.) THEN
666	! cc : .and. csum(il,i).lt.m(il,i) ) then
667	nent(il, i) = 0
668	! cc Ment(il,i,i)=m(il,i)
669	Ment(il, i, i) = 1.
670	!! Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
671	Qent(il, i, i) = qta(il,i-1) - ep(il, i)*clw(il, i)
672	uent(il, i, i) = unk(il)
673	vent(il, i, i) = vnk(il)
674	elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
675	IF (fl_cor_ebil .GE. 2) THEN
676	hent(il, i, i) = hp(il,i)
677	Sigij(il, i, i) = 0.0
678	ELSE
679	Sij(il, i, i) = 0.0
680	ENDIF
681	END IF
682	END DO ! il
683
684	!jyg! DO j = 1, ntra
685	!jyg! DO il = 1, ncum
686	!jyg! IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il) .AND. csum(il,i)<1.) THEN
687	!jyg!! cc : .and. csum(il,i).lt.m(il,i) ) then
688	!jyg! traent(il, i, i, j) = tra(il, nk(il), j)
689	!jyg! END IF
690	!jyg! END DO
691	!jyg! END DO
692
693	! ---------------------------------------------------------------
694	789 END DO ! End loop on origin level "i"
695	! ---------------------------------------------------------------
696
697
698	RETURN
699	END SUBROUTINE cv3p_mixing
700

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

Download in other formats:

Original Format