Context Navigation

cv3p_mixing.f90 @ 5843

Last change on this file since 5843 was 5843, checked in by jyg, 8 weeks ago
Getting rid of "ments" and "qents" arrays within cva_driver. Several comments to be cleared later.
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: 24.6 KB

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

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

Download in other formats:

Original Format