Context Navigation

cv3p_mixing.F90 @ 2272

Last change on this file since 2272 was 2258, checked in by Laurent Fairhead, 10 years ago
Merged trunk changes 2216:2237 into testing branch
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: 20.4 KB

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

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

Download in other formats:

Original Format