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) = (gammas*log(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) = qqa1*Qmix1(F) + qqa2*Qmix2(F) |
---|
94 | Rmix(F) = qqa1*Rmix1(F) + qqa2*Rmix2(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)/(rrv*t(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)/(rrv*t(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)-cwat*bf2) |
---|
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)*rrv*t(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)*rrv*t(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) / (cpm*rrv*t(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 | (cpm*rrv*t(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) / (cpm*rrv*t(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 | (cpm*rrv*t(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)/(cpd*rrv*t(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))*rrv*t(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))*rrv*t(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 | |
---|