Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

cv3p_mixing.f90

Last change on this file was 5703, checked in by yann meurdesoif, 2 months ago

Convection GPU porting : change some loop name indices to avoid GPUmorphosis confusion (solve in future)

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.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, &
11	ph, t, rr, rs, u, v, tra, h, lv, lf, frac, qta, &
12	unk, vnk, hp, tv, tvp, ep, clw, sig, &
13	Ment, Qent, hent, uent, vent, nent, &
14	Sigij, elij, supmax, Ments, Qents, traent)
15	! **************************************************************
16	! *
17	! CV3P_MIXING : compute mixed draught properties and, *
18	! within a scaling factor, mixed draught *
19	! mass fluxes. *
20	! written by : VTJ Philips,JY Grandpeix, 21/05/2003, 09.14.15*
21	! modified by : *
22	! **************************************************************
23
24	USE yomcst2_mod_h, ONLY : Fmax, gammas, scut, qqa1, qqa2
25	USE lmdz_cv_ini, ONLY : cpd,cpv,minorig,nl,rrv
26	USE cvflag_mod_h
27	USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
28	USE ioipsl_getin_p_mod, ONLY: getin_p
29	USE add_phys_tend_mod, ONLY: fl_cor_ebil
30
31	IMPLICIT NONE
32
33
34	!inputs:
35	INTEGER, INTENT (IN) :: ncum, nd, na
36	INTEGER, INTENT (IN) :: ntra, nloc
37	INTEGER, DIMENSION (nloc), INTENT (IN) :: icb, inb, nk
38	REAL, DIMENSION (nloc, nd), INTENT (IN) :: sig
39	REAL, DIMENSION (nloc), INTENT (IN) :: unk, vnk
40	REAL, DIMENSION (nloc, nd), INTENT (IN) :: qta
41	REAL, DIMENSION (nloc, nd+1), INTENT (IN) :: ph
42	REAL, DIMENSION (nloc, nd), INTENT (IN) :: t, rr, rs
43	REAL, DIMENSION (nloc, nd), INTENT (IN) :: u, v
44	REAL, DIMENSION (nloc, nd, ntra), INTENT (IN) :: tra ! input of convect3
45	REAL, DIMENSION (nloc, na), INTENT (IN) :: lv
46	REAL, DIMENSION (nloc, na), INTENT (IN) :: lf
47	REAL, DIMENSION (nloc, na), INTENT (IN) :: frac !ice fraction in condensate
48	REAL, DIMENSION (nloc, na), INTENT (IN) :: h !liquid water static energy of environment
49	REAL, DIMENSION (nloc, na), INTENT (IN) :: hp !liquid water static energy of air shed from adiab. asc.
50	REAL, DIMENSION (nloc, na), INTENT (IN) :: tv, tvp
51	REAL, DIMENSION (nloc, na), INTENT (IN) :: ep, clw
52
53	!outputs:
54	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: Ment, Qent
55	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: uent, vent
56	REAL, DIMENSION (nloc, na, na), INTENT (OUT) :: Sigij, elij
57	REAL, DIMENSION (nloc, na), INTENT (OUT) :: supmax ! Highest mixing fraction of mixed
58	! updraughts with the sign of (h-hp)
59	REAL, DIMENSION (nloc, nd, nd, ntra), INTENT (OUT) :: traent
60	REAL, DIMENSION (nloc, nd, nd), INTENT (OUT) :: Ments, Qents
61	REAL, DIMENSION (nloc, nd, nd), INTENT (OUT) :: hent
62	INTEGER, DIMENSION (nloc, nd), INTENT (OUT) :: nent
63
64	!local variables:
65	INTEGER i, j, k, il, im, jm
66	INTEGER num1, num2
67	REAL :: rti, bf2, anum, denom, dei, altem, cwat, stemp
68	REAL :: alt, delp, delm
69	REAL, DIMENSION (nloc) :: Qmixmax, Rmixmax, sqmrmax
70	REAL, DIMENSION (nloc) :: Qmixmin, Rmixmin, sqmrmin
71	REAL, DIMENSION (nloc) :: signhpmh
72	REAL, DIMENSION (nloc) :: Sx
73	REAL :: Scrit2
74	REAL, DIMENSION (nloc) :: Smid, Sjmin, Sjmax
75	REAL, DIMENSION (nloc) :: Sbef, sup, smin
76	REAL, DIMENSION (nloc) :: ASij, ASij_inv, smax, Scrit
77	REAL, DIMENSION (nloc, nd, nd) :: Sij
78	REAL, DIMENSION (nloc, nd) :: csum
79	REAL :: awat
80	REAL :: cpm !Mixed draught heat capacity
81	REAL :: Tm !Mixed draught temperature
82	LOGICAL, DIMENSION (nloc) :: lwork
83
84	REAL amxupcrit, df, ff
85	INTEGER nstep
86
87	INTEGER,PARAMETER :: igout=1
88
89	! -- Mixing probability distribution functions
90
91	REAL Qcoef1, Qcoef2, QFF, QFFF, Qmix, Rmix, Qmix1, Rmix1, Qmix2, Rmix2, F
92	REAL :: Qcoef1max,Qcoef2max !ym WARNING
93	!ym redefine local variable instead use module variable
94	!ym to check when refactoring deep convection
95	!ym => eliminate "first" SAVE variable
96	!ym probably all these folowing lines will be removed
97	Qcoef1(F) = tanh(F/gammas)
98	Qcoef2(F) = (tanh(F/gammas)+gammas*log(cosh((1.-F)/gammas)/cosh(F/gammas)))
99	QFF(F) = max(min(F,1.), 0.)
100	QFFf(F) = min(QFF(F), scut)
101	Qmix1(F) = (tanh((QFF(F)-Fmax)/gammas)+Qcoef1max)/Qcoef2max
102	Rmix1(F) = (gammaslog(cosh((QFF(F)-Fmax)/gammas))+QFF(F)Qcoef1max)/Qcoef2max
103	Qmix2(F) = -log(1.-QFFf(F))/scut
104	Rmix2(F) = (QFFf(F)+(1.-QFF(F))*log(1.-QFFf(F)))/scut
105	Qmix(F) = qqa1Qmix1(F) + qqa2Qmix2(F)
106	Rmix(F) = qqa1Rmix1(F) + qqa2Rmix2(F)
107
108
109
110	! =====================================================================
111	! --- INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS
112	! =====================================================================
113
114	Qcoef1max = Qcoef1(Fmax)
115	Qcoef2max = Qcoef2(Fmax)
116
117	! ori do 360 i=1,ncum*nlp
118	DO j = 1, nl
119	DO il = 1, ncum
120	nent(il, j) = 0
121	! in convect3, m is computed in cv3_closure
122	! ori m(i,1)=0.0
123	END DO
124	END DO
125
126	! ori do 400 k=1,nlp
127	! ori do 390 j=1,nlp
128	DO j = 1, nl
129	DO k = 1, nl
130	DO il = 1, ncum
131	Qent(il, k, j) = rr(il, j)
132	uent(il, k, j) = u(il, j)
133	vent(il, k, j) = v(il, j)
134	elij(il, k, j) = 0.0
135	hent(il, k, j) = 0.0
136	!AC! Ment(i,k,j)=0.0
137	!AC! Sij(i,k,j)=0.0
138	END DO
139	END DO
140	END DO
141
142	!AC!
143	Ment(1:ncum, 1:nd, 1:nd) = 0.0
144	Sij(1:ncum, 1:nd, 1:nd) = 0.0
145	!AC!
146	!ym
147	Sigij(1:ncum, 1:nd, 1:nd) = 0.0
148	!ym
149
150	!jyg! DO k = 1, ntra
151	!jyg! DO j = 1, nd ! instead nlp
152	!jyg! DO i = 1, nd ! instead nlp
153	!jyg! DO il = 1, ncum
154	!jyg! traent(il, i, j, k) = tra(il, j, k)
155	!jyg! END DO
156	!jyg! END DO
157	!jyg! END DO
158	!jyg! END DO
159
160	! =====================================================================
161	! --- CALCULATE ENTRAINED AIR MASS FLUX (Ment), TOTAL WATER MIXING
162	! --- RATIO (QENT), TOTAL CONDENSED WATER (elij), AND MIXING
163	! --- FRACTION (Sij)
164	! =====================================================================
165
166	DO i = minorig + 1, nl
167
168	IF (ok_entrain) THEN
169	DO j = minorig, nl
170	DO il = 1, ncum
171	IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (j>=(icb(il)-1)) &
172	.AND. (j<=inb(il))) THEN
173
174	!! rti = qnk(il) - ep(il, i)*clw(il, i)
175	rti = qta(il,i-1) - ep(il, i)*clw(il, i)
176	bf2 = 1. + lv(il, j)lv(il, j)rs(il, j)/(rrvt(il,j)t(il,j)*cpd)
177	!jyg(from aj)<
178	IF (cvflag_ice) THEN
179	! print*,cvflag_ice,'cvflag_ice dans do 700'
180	IF (t(il,j)<=263.15) THEN
181	bf2 = 1. + (lf(il,j)+lv(il,j))(lv(il,j)+frac(il,j) &
182	lf(il,j))rs(il, j)/(rrvt(il,j)t(il,j)cpd)
183	END IF
184	END IF
185	!>jyg
186	anum = h(il, j) - hp(il, i) + (cpv-cpd)t(il, j)(rti-rr(il,j))
187	denom = h(il, i) - hp(il, i) + (cpd-cpv)(rr(il,i)-rti)t(il, j)
188	dei = denom
189	IF (abs(dei)<0.01) dei = 0.01
190	Sij(il, i, j) = anum/dei
191	Sij(il, i, i) = 1.0
192	altem = Sij(il, i, j)rr(il, i) + (1.-Sij(il,i,j))rti - rs(il, j)
193	altem = altem/bf2
194	cwat = clw(il, j)*(1.-ep(il,j))
195	stemp = Sij(il, i, j)
196	IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN
197	!jyg(from aj)<
198	IF (cvflag_ice) THEN
199	anum = anum - (lv(il,j)+frac(il,j)lf(il,j))(rti-rs(il,j)-cwat*bf2)
200	denom = denom + (lv(il,j)+frac(il,j)lf(il,j))(rr(il,i)-rti)
201	ELSE
202	anum = anum - lv(il, j)(rti-rs(il,j)-cwatbf2)
203	denom = denom + lv(il, j)*(rr(il,i)-rti)
204	END IF
205	!>jyg
206	IF (abs(denom)<0.01) denom = 0.01
207	Sij(il, i, j) = anum/denom
208	altem = Sij(il, i, j)rr(il, i) + (1.-Sij(il,i,j))rti - rs(il, j)
209	altem = altem - (bf2-1.)*cwat
210	END IF
211	IF (Sij(il,i,j)>0.0) THEN
212	!!! Ment(il,i,j)=m(il,i)
213	Ment(il, i, j) = 1.
214	elij(il, i, j) = altem
215	elij(il, i, j) = amax1(0.0, elij(il,i,j))
216	nent(il, i) = nent(il, i) + 1
217	END IF
218
219	Sij(il, i, j) = amax1(0.0, Sij(il,i,j))
220	Sij(il, i, j) = amin1(1.0, Sij(il,i,j))
221	ELSE IF (j > i) THEN
222	IF (prt_level >= 10) THEN
223	print *,'cv3p_mixing i, j, Sij given by the no-precip eq. ', i, j, Sij(il,i,j)
224	ENDIF
225	END IF ! new
226	END DO
227	END DO
228	ELSE ! (ok_entrain)
229	DO il = 1,ncum
230	nent(il,i) = 0
231	ENDDO
232	ENDIF ! (ok_entrain)
233
234	!jygdebug<
235	IF (prt_level >= 10) THEN
236	print *,'cv3p_mixing i, nent(i), icb, inb ',i, nent(igout,i), icb(igout), inb(igout)
237	IF (nent(igout,i) .gt. 0) THEN
238	print *,'i,(j,Sij(i,j),j=icb-1,inb) ',i,(j,Sij(igout,i,j),j=icb(igout)-1,inb(igout))
239	ENDIF
240	ENDIF
241	!>jygdebug
242
243	! * if no air can entrain at level i assume that updraft detrains *
244	! * at that level and calculate detrained air flux and properties *
245
246
247	! @ do 170 i=icb(il),inb(il)
248
249	DO il = 1, ncum
250	IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (nent(il,i)==0)) THEN
251	! @ if(nent(il,i).eq.0)then
252	!!! Ment(il,i,i)=m(il,i)
253	Ment(il, i, i) = 1.
254	!! Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
255	Qent(il, i, i) = qta(il,i-1) - ep(il, i)*clw(il, i)
256	uent(il, i, i) = unk(il)
257	vent(il, i, i) = vnk(il)
258	IF (fl_cor_ebil .GE. 2) THEN
259	hent(il, i, i) = hp(il,i)
260	ENDIF
261	elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
262	Sij(il, i, i) = 0.0
263	END IF
264	END DO
265	END DO ! i = minorig + 1, nl
266
267	!jyg! DO j = 1, ntra
268	!jyg! DO i = minorig + 1, nl
269	!jyg! DO il = 1, ncum
270	!jyg! IF (i>=icb(il) .AND. i<=inb(il) .AND. nent(il,i)==0) THEN
271	!jyg! traent(il, i, i, j) = tra(il, nk(il), j)
272	!jyg! END IF
273	!jyg! END DO
274	!jyg! END DO
275	!jyg! END DO
276
277	DO j = minorig, nl
278	DO i = minorig, nl
279	DO il = 1, ncum
280	IF ((j>=(icb(il)-1)) .AND. (j<=inb(il)) .AND. &
281	(i>=icb(il)) .AND. (i<=inb(il))) THEN
282	Sigij(il, i, j) = Sij(il, i, j)
283	END IF
284	END DO
285	END DO
286	END DO
287	! @ enddo
288
289	! @170 continue
290
291	! =====================================================================
292	! --- NORMALIZE ENTRAINED AIR MASS FLUXES
293	! --- TO REPRESENT EQUAL PROBABILITIES OF MIXING
294	! =====================================================================
295
296	csum(:,:) = 0.
297
298	DO il = 1, ncum
299	lwork(il) = .FALSE.
300	END DO
301
302	! ---------------------------------------------------------------
303	DO i = minorig + 1, nl !Loop on origin level "i"
304	! ---------------------------------------------------------------
305
306	num1 = 0
307	DO il = 1, ncum
308	IF (i>=icb(il) .AND. i<=inb(il)) num1 = num1 + 1
309	END DO
310	!ym IF (num1<=0) GO TO 789
311	IF (num1<=0) CYCLE
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	!ym IF (num2<=0) GO TO 175
392	IF (num2<=0) CYCLE
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	END DO !ym label 175 ! 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	END DO !ym label 789 ! End loop on origin level "i"
695
696	! ---------------------------------------------------------------
697
698
699	RETURN
700	END SUBROUTINE cv3p_mixing
701
702	END MODULE cv3p_mixing_mod

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

Download in other formats:

Original Format