/[lmdze]/trunk/Sources/phylmd/CV30_routines/cv30_mixing.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/CV30_routines/cv30_mixing.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 178 - (show annotations)
Fri Mar 11 18:47:26 2016 UTC (8 years, 2 months ago) by guez
Original Path: trunk/Sources/phylmd/CV3_routines/cv3_mixing.f
File size: 10527 byte(s)
Moved variables date0, deltat, datasz_max, ncvar_ids, point, buff_pos,
buffer, regular from module histcom_var to modules where they are
defined.

Removed procedure ioipslmpp, useless for a sequential program.

Added argument datasz_max to histwrite_real (to avoid circular
dependency with histwrite).

Removed useless variables and computations everywhere.

Changed real litteral constants from default kind to double precision
in lwb, lwu, lwvn, sw1s, swtt, swtt1, swu.

Removed unused arguments: paer of sw, sw1s, sw2s, swclr; pcldsw of
sw1s, sw2s; pdsig, prayl of swr; co2_ppm of clmain, clqh; tsol of
transp_lay; nsrf of screenp; kcrit and kknu of gwstress; pstd of
orosetup.

Added output of relative humidity.

1 module cv3_mixing_m
2
3 implicit none
4
5 contains
6
7 SUBROUTINE cv3_mixing(nloc, ncum, nd, na, icb, nk, inb, t, rr, rs, u, v, h, &
8 lv, hp, ep, clw, m, sig, ment, qent, uent, vent, nent, sij, elij, &
9 ments, qents)
10 use cv3_param_m
11 use cvthermo
12
13 !---------------------------------------------------------------------
14 ! a faire:
15 ! - changer rr(il, 1) -> qnk(il)
16 ! - vectorisation de la partie normalisation des flux (do 789)
17 !---------------------------------------------------------------------
18
19 ! inputs:
20 integer, intent(in):: ncum, nd, na, nloc
21 integer icb(nloc), inb(nloc), nk(nloc)
22 real sig(nloc, nd)
23 real t(nloc, nd), rr(nloc, nd), rs(nloc, nd)
24 real u(nloc, nd), v(nloc, nd)
25 real lv(nloc, na), h(nloc, na), hp(nloc, na)
26 real ep(nloc, na), clw(nloc, na)
27 real m(nloc, na) ! input of convect3
28
29 ! outputs:
30 real ment(nloc, na, na), qent(nloc, na, na)
31 real uent(nloc, na, na), vent(nloc, na, na)
32 real sij(nloc, na, na), elij(nloc, na, na)
33 real ments(nloc, nd, nd), qents(nloc, nd, nd)
34 integer nent(nloc, nd)
35
36 ! local variables:
37 integer i, j, k, il, im, jm
38 integer num1, num2
39 real rti, bf2, anum, denom, dei, altem, cwat, stemp, qp
40 real alt, smid, sjmin, sjmax, delp, delm
41 real asij(nloc), smax(nloc), scrit(nloc)
42 real asum(nloc, nd), bsum(nloc, nd), csum(nloc, nd)
43 real wgh
44 real zm(nloc, na)
45 logical lwork(nloc)
46
47 ! --- INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS
48
49 do j=1, nl
50 do i=1, ncum
51 nent(i, j)=0
52 ! in convect3, m is computed in cv3_closure
53 ! ori m(i, 1)=0.0
54 end do
55 end do
56
57 do j=1, nl
58 do k=1, nl
59 do i=1, ncum
60 qent(i, k, j)=rr(i, j)
61 uent(i, k, j)=u(i, j)
62 vent(i, k, j)=v(i, j)
63 elij(i, k, j)=0.0
64 !ym ment(i, k, j)=0.0
65 !ym sij(i, k, j)=0.0
66 end do
67 end do
68 end do
69
70 !ym
71 ment(1:ncum, 1:nd, 1:nd)=0.0
72 sij(1:ncum, 1:nd, 1:nd)=0.0
73
74 zm(:, :)=0.
75
76 ! --- CALCULATE ENTRAINED AIR MASS FLUX (ment), TOTAL WATER MIXING
77 ! --- RATIO (QENT), TOTAL CONDENSED WATER (elij), AND MIXING
78 ! --- FRACTION (sij)
79
80 do i=minorig+1, nl
81
82 do j=minorig, nl
83 do il=1, ncum
84 if((i >= icb(il)).and.(i <= inb(il)).and. &
85 (j >= (icb(il)-1)).and.(j <= inb(il)))then
86
87 rti=rr(il, 1)-ep(il, i)*clw(il, i)
88 bf2=1.+lv(il, j)*lv(il, j)*rs(il, j)/(rrv*t(il, j)*t(il, j)*cpd)
89 anum=h(il, j)-hp(il, i)+(cpv-cpd)*t(il, j)*(rti-rr(il, j))
90 denom=h(il, i)-hp(il, i)+(cpd-cpv)*(rr(il, i)-rti)*t(il, j)
91 dei=denom
92 if(abs(dei) < 0.01)dei=0.01
93 sij(il, i, j)=anum/dei
94 sij(il, i, i)=1.0
95 altem=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti-rs(il, j)
96 altem=altem/bf2
97 cwat=clw(il, j)*(1.-ep(il, j))
98 stemp=sij(il, i, j)
99 if((stemp < 0.0.or.stemp > 1.0.or.altem > cwat) &
100 .and.j > i)then
101 anum=anum-lv(il, j)*(rti-rs(il, j)-cwat*bf2)
102 denom=denom+lv(il, j)*(rr(il, i)-rti)
103 if(abs(denom) < 0.01)denom=0.01
104 sij(il, i, j)=anum/denom
105 altem=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti-rs(il, j)
106 altem=altem-(bf2-1.)*cwat
107 end if
108 if(sij(il, i, j) > 0.0.and.sij(il, i, j) < 0.95)then
109 qent(il, i, j)=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti
110 uent(il, i, j)=sij(il, i, j)*u(il, i)+(1.-sij(il, i, j))*u(il, nk(il))
111 vent(il, i, j)=sij(il, i, j)*v(il, i)+(1.-sij(il, i, j))*v(il, nk(il))
112 elij(il, i, j)=altem
113 elij(il, i, j)=amax1(0.0, elij(il, i, j))
114 ment(il, i, j)=m(il, i)/(1.-sij(il, i, j))
115 nent(il, i)=nent(il, i)+1
116 end if
117 sij(il, i, j)=amax1(0.0, sij(il, i, j))
118 sij(il, i, j)=amin1(1.0, sij(il, i, j))
119 endif ! new
120 end do
121 end do
122
123 ! *** if no air can entrain at level i assume that updraft detrains ***
124 ! *** at that level and calculate detrained air flux and properties ***
125
126 !@ do 170 i=icb(il), inb(il)
127
128 do il=1, ncum
129 if ((i >= icb(il)).and.(i <= inb(il)).and.(nent(il, i) == 0)) then
130 !@ if(nent(il, i) == 0)then
131 ment(il, i, i)=m(il, i)
132 qent(il, i, i)=rr(il, nk(il))-ep(il, i)*clw(il, i)
133 uent(il, i, i)=u(il, nk(il))
134 vent(il, i, i)=v(il, nk(il))
135 elij(il, i, i)=clw(il, i)
136 !MAF sij(il, i, i)=1.0
137 sij(il, i, i)=0.0
138 end if
139 end do
140 end do
141
142 ! --- NORMALIZE ENTRAINED AIR MASS FLUXES
143 ! --- TO REPRESENT EQUAL PROBABILITIES OF MIXING
144
145 call zilch(asum, nloc*nd)
146 call zilch(csum, nloc*nd)
147 call zilch(csum, nloc*nd)
148
149 do il=1, ncum
150 lwork(il) = .FALSE.
151 enddo
152
153 DO i=minorig+1, nl
154
155 num1=0
156 do il=1, ncum
157 if (i >= icb(il) .and. i <= inb(il)) num1=num1+1
158 enddo
159 if (num1 <= 0) cycle
160
161 do il=1, ncum
162 if (i >= icb(il) .and. i <= inb(il)) then
163 lwork(il)=(nent(il, i) /= 0)
164 qp=rr(il, 1)-ep(il, i)*clw(il, i)
165 anum=h(il, i)-hp(il, i)-lv(il, i)*(qp-rs(il, i)) &
166 +(cpv-cpd)*t(il, i)*(qp-rr(il, i))
167 denom=h(il, i)-hp(il, i)+lv(il, i)*(rr(il, i)-qp) &
168 +(cpd-cpv)*t(il, i)*(rr(il, i)-qp)
169 if(abs(denom) < 0.01)denom=0.01
170 scrit(il)=anum/denom
171 alt=qp-rs(il, i)+scrit(il)*(rr(il, i)-qp)
172 if(scrit(il) <= 0.0.or.alt <= 0.0)scrit(il)=1.0
173 smax(il)=0.0
174 asij(il)=0.0
175 endif
176 end do
177
178 do j=nl, minorig, -1
179
180 num2=0
181 do il=1, ncum
182 if (i >= icb(il) .and. i <= inb(il) .and. &
183 j >= (icb(il)-1) .and. j <= inb(il) &
184 .and. lwork(il)) num2=num2+1
185 enddo
186 if (num2 <= 0) cycle
187
188 do il=1, ncum
189 if (i >= icb(il) .and. i <= inb(il) .and. &
190 j >= (icb(il)-1) .and. j <= inb(il) &
191 .and. lwork(il)) then
192
193 if(sij(il, i, j) > 1.0e-16.and.sij(il, i, j) < 0.95)then
194 wgh=1.0
195 if(j > i)then
196 sjmax=amax1(sij(il, i, j+1), smax(il))
197 sjmax=amin1(sjmax, scrit(il))
198 smax(il)=amax1(sij(il, i, j), smax(il))
199 sjmin=amax1(sij(il, i, j-1), smax(il))
200 sjmin=amin1(sjmin, scrit(il))
201 if(sij(il, i, j) < (smax(il)-1.0e-16))wgh=0.0
202 smid=amin1(sij(il, i, j), scrit(il))
203 else
204 sjmax=amax1(sij(il, i, j+1), scrit(il))
205 smid=amax1(sij(il, i, j), scrit(il))
206 sjmin=0.0
207 if(j > 1)sjmin=sij(il, i, j-1)
208 sjmin=amax1(sjmin, scrit(il))
209 endif
210 delp=abs(sjmax-smid)
211 delm=abs(sjmin-smid)
212 asij(il)=asij(il)+wgh*(delp+delm)
213 ment(il, i, j)=ment(il, i, j)*(delp+delm)*wgh
214 endif
215 endif
216 end do
217
218 end do
219
220 do il=1, ncum
221 if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
222 asij(il)=amax1(1.0e-16, asij(il))
223 asij(il)=1.0/asij(il)
224 asum(il, i)=0.0
225 bsum(il, i)=0.0
226 csum(il, i)=0.0
227 endif
228 enddo
229
230 do j=minorig, nl
231 do il=1, ncum
232 if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
233 .and. j >= (icb(il)-1) .and. j <= inb(il)) then
234 ment(il, i, j)=ment(il, i, j)*asij(il)
235 endif
236 enddo
237 end do
238
239 do j=minorig, nl
240 do il=1, ncum
241 if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
242 .and. j >= (icb(il)-1) .and. j <= inb(il)) then
243 asum(il, i)=asum(il, i)+ment(il, i, j)
244 ment(il, i, j)=ment(il, i, j)*sig(il, j)
245 bsum(il, i)=bsum(il, i)+ment(il, i, j)
246 endif
247 enddo
248 end do
249
250 do il=1, ncum
251 if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
252 bsum(il, i)=amax1(bsum(il, i), 1.0e-16)
253 bsum(il, i)=1.0/bsum(il, i)
254 endif
255 enddo
256
257 do j=minorig, nl
258 do il=1, ncum
259 if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
260 .and. j >= (icb(il)-1) .and. j <= inb(il)) then
261 ment(il, i, j)=ment(il, i, j)*asum(il, i)*bsum(il, i)
262 endif
263 enddo
264 end do
265
266 do j=minorig, nl
267 do il=1, ncum
268 if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
269 .and. j >= (icb(il)-1) .and. j <= inb(il)) then
270 csum(il, i)=csum(il, i)+ment(il, i, j)
271 endif
272 enddo
273 end do
274
275 do il=1, ncum
276 if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
277 .and. csum(il, i) < m(il, i)) then
278 nent(il, i)=0
279 ment(il, i, i)=m(il, i)
280 qent(il, i, i)=rr(il, 1)-ep(il, i)*clw(il, i)
281 uent(il, i, i)=u(il, nk(il))
282 vent(il, i, i)=v(il, nk(il))
283 elij(il, i, i)=clw(il, i)
284 !MAF sij(il, i, i)=1.0
285 sij(il, i, i)=0.0
286 endif
287 enddo ! il
288
289 end DO
290
291 ! MAF: renormalisation de MENT
292 do jm=1, nd
293 do im=1, nd
294 do il=1, ncum
295 zm(il, im)=zm(il, im)+(1.-sij(il, im, jm))*ment(il, im, jm)
296 end do
297 end do
298 end do
299
300 do jm=1, nd
301 do im=1, nd
302 do il=1, ncum
303 if(zm(il, im) /= 0.) then
304 ment(il, im, jm)=ment(il, im, jm)*m(il, im)/zm(il, im)
305 endif
306 end do
307 end do
308 end do
309
310 do jm=1, nd
311 do im=1, nd
312 do il=1, ncum
313 qents(il, im, jm)=qent(il, im, jm)
314 ments(il, im, jm)=ment(il, im, jm)
315 end do
316 enddo
317 enddo
318
319 end SUBROUTINE cv3_mixing
320
321 end module cv3_mixing_m

  ViewVC Help
Powered by ViewVC 1.1.21