Context Navigation

cv3p_mixing.F90 @ 2260

Last change on this file since 2260 was 2226, checked in by Laurent Fairhead, 10 years ago
correction sur le niveau de detrainement
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
[2007]	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
[2007]	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
[2007]	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
[2007]	48	!local variables:
[1992]	49	INTEGER i, j, k, il, im, jm
	50	INTEGER num1, num2
[2007]	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
[2226]	60	!jyg REAL, DIMENSION (nloc) :: ASij, smax, Scrit
	61	REAL, DIMENSION (nloc) :: ASij, ASij_inv, smax, Scrit
[2007]	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
[2007]	70	! -- Mixing probability distribution functions
[1650]	71
[2007]	72	REAL Qcoef1, Qcoef2, QFF, QFFF, Qmix, Rmix, Qmix1, Rmix1, Qmix2, Rmix2, F
[879]	73
[2007]	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/
[2007]	87	!$OMP THREADPRIVATE(ifrst)
[879]	88
	89
[2007]	90	! =====================================================================
	91	! --- INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS
	92	! =====================================================================
[879]	93
[2007]	94	! -- Initialize mixing PDF coefficients
[1992]	95	IF (ifrst==0) THEN
	96	ifrst = 1
[2007]	97	Qcoef1max = Qcoef1(Fmax)
	98	Qcoef2max = Qcoef2(Fmax)
[879]	99
[1992]	100	END IF
[879]	101
	102
[2007]	103	! ori do 360 i=1,ncum*nlp
[1992]	104	DO j = 1, nl
	105	DO i = 1, ncum
	106	nent(i, j) = 0
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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))
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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))
[2007]	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
[2007]	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
[2007]	231	! @ enddo
[1041]	232
[2007]	233	! @170 continue
[1041]	234
[2007]	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
[2007]	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
[2007]	257	!JYG1 Find maximum of SIJ for J>I, if any.
[879]	258
[2007]	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))
[2007]	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
[2007]	271	IF (Sbef(il)<Sij(il,i,j)) THEN
	272	Sx(il) = max(Sij(il,i,j), Sx(il))
[1992]	273	END IF
[2007]	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)
[2007]	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
[2007]	289	Scrit(il) = min(anum/denom, 1.)
	290	alt = rti - rs(il, i) + Scrit(il)*(rr(il,i)-rti)
[1992]	291
[2007]	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
[2007]	296	Scrit2 = min(Scrit(il), Sx(il))*max(0., -signhpmh(il)) + &
	297	Scrit(il)*max(0., signhpmh(il))
[1992]	298
[2007]	299	Scrit(il) = Scrit2
[1992]	300
[2007]	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
[2007]	307	ASij(il) = 0.0
	308	sup(il) = 0. ! upper S-value reached by descending draughts
[1992]	309	END IF
	310	END DO
	311
[2007]	312	! ---------------------------------------------------------------
	313	DO j = minorig, nl !Loop on destination level "j"
	314	! ---------------------------------------------------------------
[1992]	315
	316	num2 = 0
	317	DO il = 1, ncum
[2007]	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
[2007]	324	! -----------------------------------------------
[1992]	325	IF (j>i) THEN
[2007]	326	! -----------------------------------------------
[1992]	327	DO il = 1, ncum
[2007]	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
[2007]	345	! -----------------------------------------------
[1992]	346	ELSE IF (j==i) THEN
[2007]	347	! -----------------------------------------------
[1992]	348	DO il = 1, ncum
[2007]	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
[2007]	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.
[2007]	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
[2007]	368	! -----------------------------------------------
[1992]	369	ELSE IF (j<i) THEN
[2007]	370	! -----------------------------------------------
[1992]	371	DO il = 1, ncum
[2007]	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
[2007]	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
[2007]	392	! -----------------------------------------------
[1992]	393	END IF
[2007]	394	! -----------------------------------------------
[1992]	395
	396
	397	DO il = 1, ncum
[2007]	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)
[2007]	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
[2007]	410	Ment(il, i, j) = abs(Qmixmax(il)-Qmixmin(il))*Ment(il, i, j)
[1992]	411
[2007]	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
[2007]	416	Sigij(il, i, j) = 0.
[1992]	417	END IF
	418
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	456	! print ,h(il,j)-hent(il,i,j),LV(il,j)rs(il,j)/(cpdrrvt(il,j)*
	457	! : t(il,j))
[1992]	458
[2007]	459	hent(il, i, j) = hent(il, i, j) + (lv(il,j)+(cpd-cpv)t(il,j))awat
	460	!IM 301008 end
[1992]	461
[2007]	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
[2007]	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
[2007]	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
[2007]	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
[2007]	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)
[2007]	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))
[2007]	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
[2007]	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
[2007]	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
[2007]	526	ASij(il) = amax1(1.0E-16, ASij(il))
[2226]	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
[2007]	539	IF (i>=icb(il) .AND. i<=inb(il) .AND. lwork(il) .AND. &
	540	j>=(icb(il)-1) .AND. j<=inb(il)) THEN
[2226]	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
[2007]	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
[2007]	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
[2007]	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))
[2007]	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
[2007]	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
[2007]	579	! ---------------------------------------------------------------
	580	789 END DO ! End loop on origin level "i"
	581	! ---------------------------------------------------------------
[1992]	582
[2007]	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