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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 186 - (show annotations)
Mon Mar 21 15:36:26 2016 UTC (8 years, 1 month ago) by guez
File size: 24718 byte(s)
Removed variables nlm and nlp of module cv30_param_m. We do not
believe much in the benefit of these intermediary variables so we go
for clarity.

Removed variable noff of module cv30_param_m. Never used anywhere
else. Just set the value of nl explicitly in cv30_param.

Removed argument nd of cv30_param. Only called with nd = klev.

Replaced calls to zilch by array assignments. There was a strange
double call to zilch with the same arguments in cv30_mixing.

Removed procedure cv_flag. Just set the value of variable cvflag_grav
of module cvflag at declaration.

1 module cv30_yield_m
2
3 implicit none
4
5 contains
6
7 SUBROUTINE cv30_yield(nloc,ncum,nd,na &
8 ,icb,inb,delt &
9 ,t,rr,u,v,gz,p,ph,h,hp,lv,cpn,th &
10 ,ep,clw,m,tp,mp,rp,up,vp &
11 ,wt,water,evap,b &
12 ,ment,qent,uent,vent,nent,elij,sig &
13 ,tv,tvp &
14 ,iflag,precip,VPrecip,ft,fr,fu,fv &
15 ,upwd,dnwd,dnwd0,ma,mike,tls,tps,qcondc,wd)
16 use conema3_m
17 use cv30_param_m
18 use cvthermo
19 use cvflag
20
21 ! inputs:
22 integer, intent(in):: ncum,nd,na,nloc
23 integer icb(nloc), inb(nloc)
24 real, intent(in):: delt
25 real t(nloc,nd), rr(nloc,nd), u(nloc,nd), v(nloc,nd)
26 real sig(nloc,nd)
27 real gz(nloc,na), ph(nloc,nd+1), h(nloc,na), hp(nloc,na)
28 real th(nloc,na), p(nloc,nd), tp(nloc,na)
29 real lv(nloc,na), cpn(nloc,na), ep(nloc,na), clw(nloc,na)
30 real m(nloc,na), mp(nloc,na), rp(nloc,na), up(nloc,na)
31 real vp(nloc,na), wt(nloc,nd)
32 real water(nloc,na), evap(nloc,na), b(nloc,na)
33 real ment(nloc,na,na), qent(nloc,na,na), uent(nloc,na,na)
34 !ym real vent(nloc,na,na), nent(nloc,na), elij(nloc,na,na)
35 real vent(nloc,na,na), elij(nloc,na,na)
36 integer nent(nloc,na)
37 real tv(nloc,nd), tvp(nloc,nd)
38
39 ! input/output:
40 integer iflag(nloc)
41
42 ! outputs:
43 real precip(nloc)
44 real VPrecip(nloc,nd+1)
45 real ft(nloc,nd), fr(nloc,nd), fu(nloc,nd), fv(nloc,nd)
46 real upwd(nloc,nd), dnwd(nloc,nd), ma(nloc,nd)
47 real dnwd0(nloc,nd), mike(nloc,nd)
48 real tls(nloc,nd), tps(nloc,nd)
49 real qcondc(nloc,nd) ! cld
50 real wd(nloc) ! gust
51
52 ! local variables:
53 integer i,k,il,n,j,num1
54 real rat, awat, delti
55 real ax, bx, cx, dx
56 real cpinv, rdcp, dpinv
57 real lvcp(nloc,na)
58 real am(nloc), work(nloc), ad(nloc), amp1(nloc)
59 !!! real up1(nloc), dn1(nloc)
60 real up1(nloc,nd,nd), dn1(nloc,nd,nd)
61 real asum(nloc), bsum(nloc), csum(nloc), dsum(nloc)
62 real qcond(nloc,nd), nqcond(nloc,nd), wa(nloc,nd) ! cld
63 real siga(nloc,nd), sax(nloc,nd), mac(nloc,nd) ! cld
64
65
66 !-------------------------------------------------------------
67
68 ! initialization:
69
70 delti = 1.0/delt
71
72 do il=1,ncum
73 precip(il)=0.0
74 wd(il)=0.0 ! gust
75 VPrecip(il,nd+1)=0.
76 enddo
77
78 do i=1,nd
79 do il=1,ncum
80 VPrecip(il,i)=0.0
81 ft(il,i)=0.0
82 fr(il,i)=0.0
83 fu(il,i)=0.0
84 fv(il,i)=0.0
85 qcondc(il,i)=0.0 ! cld
86 qcond(il,i)=0.0 ! cld
87 nqcond(il,i)=0.0 ! cld
88 enddo
89 enddo
90
91
92 do i=1,nl
93 do il=1,ncum
94 lvcp(il,i)=lv(il,i)/cpn(il,i)
95 enddo
96 enddo
97
98
99 !
100 ! *** calculate surface precipitation in mm/day ***
101 !
102 do il=1,ncum
103 if(ep(il,inb(il)).ge.0.0001)then
104 if (cvflag_grav) then
105 precip(il)=wt(il,1)*sigd*water(il,1)*86400.*1000./(rowl*grav)
106 else
107 precip(il)=wt(il,1)*sigd*water(il,1)*8640.
108 endif
109 endif
110 enddo
111
112 ! *** CALCULATE VERTICAL PROFILE OF PRECIPITATIONs IN kg/m2/s ===
113 !
114 ! MAF rajout pour lessivage
115 do k=1,nl
116 do il=1,ncum
117 if (k.le.inb(il)) then
118 if (cvflag_grav) then
119 VPrecip(il,k) = wt(il,k)*sigd*water(il,k)/grav
120 else
121 VPrecip(il,k) = wt(il,k)*sigd*water(il,k)/10.
122 endif
123 endif
124 end do
125 end do
126 !
127 !
128 !
129 ! *** calculate tendencies of lowest level potential temperature ***
130 ! *** and mixing ratio ***
131 !
132 do il=1,ncum
133 work(il)=1.0/(ph(il,1)-ph(il,2))
134 am(il)=0.0
135 enddo
136
137 do k=2,nl
138 do il=1,ncum
139 if (k.le.inb(il)) then
140 am(il)=am(il)+m(il,k)
141 endif
142 enddo
143 enddo
144
145 do il=1,ncum
146
147 ! convect3 if((0.1*dpinv*am).ge.delti)iflag(il)=4
148 if (cvflag_grav) then
149 if((0.01*grav*work(il)*am(il)).ge.delti)iflag(il)=1!consist vect
150 ft(il,1)=0.01*grav*work(il)*am(il)*(t(il,2)-t(il,1) &
151 +(gz(il,2)-gz(il,1))/cpn(il,1))
152 else
153 if((0.1*work(il)*am(il)).ge.delti)iflag(il)=1 !consistency vect
154 ft(il,1)=0.1*work(il)*am(il)*(t(il,2)-t(il,1) &
155 +(gz(il,2)-gz(il,1))/cpn(il,1))
156 endif
157
158 ft(il,1)=ft(il,1)-0.5*lvcp(il,1)*sigd*(evap(il,1)+evap(il,2))
159
160 if (cvflag_grav) then
161 ft(il,1)=ft(il,1)-0.009*grav*sigd*mp(il,2) &
162 *t(il,1)*b(il,1)*work(il)
163 else
164 ft(il,1)=ft(il,1)-0.09*sigd*mp(il,2)*t(il,1)*b(il,1)*work(il)
165 endif
166
167 ft(il,1)=ft(il,1)+0.01*sigd*wt(il,1)*(cl-cpd)*water(il,2)*(t(il,2) &
168 -t(il,1))*work(il)/cpn(il,1)
169
170 if (cvflag_grav) then
171 !jyg1 Correction pour mieux conserver l'eau (conformite avec CONVECT4.3)
172 ! (sb: pour l'instant, on ne fait que le chgt concernant grav, pas evap)
173 fr(il,1)=0.01*grav*mp(il,2)*(rp(il,2)-rr(il,1))*work(il) &
174 +sigd*0.5*(evap(il,1)+evap(il,2))
175 !+tard : +sigd*evap(il,1)
176
177 fr(il,1)=fr(il,1)+0.01*grav*am(il)*(rr(il,2)-rr(il,1))*work(il)
178
179 fu(il,1)=fu(il,1)+0.01*grav*work(il)*(mp(il,2)*(up(il,2)-u(il,1)) &
180 +am(il)*(u(il,2)-u(il,1)))
181 fv(il,1)=fv(il,1)+0.01*grav*work(il)*(mp(il,2)*(vp(il,2)-v(il,1)) &
182 +am(il)*(v(il,2)-v(il,1)))
183 else ! cvflag_grav
184 fr(il,1)=0.1*mp(il,2)*(rp(il,2)-rr(il,1))*work(il) &
185 +sigd*0.5*(evap(il,1)+evap(il,2))
186 fr(il,1)=fr(il,1)+0.1*am(il)*(rr(il,2)-rr(il,1))*work(il)
187 fu(il,1)=fu(il,1)+0.1*work(il)*(mp(il,2)*(up(il,2)-u(il,1)) &
188 +am(il)*(u(il,2)-u(il,1)))
189 fv(il,1)=fv(il,1)+0.1*work(il)*(mp(il,2)*(vp(il,2)-v(il,1)) &
190 +am(il)*(v(il,2)-v(il,1)))
191 endif ! cvflag_grav
192
193 enddo ! il
194
195 do j=2,nl
196 do il=1,ncum
197 if (j.le.inb(il)) then
198 if (cvflag_grav) then
199 fr(il,1)=fr(il,1) &
200 +0.01*grav*work(il)*ment(il,j,1)*(qent(il,j,1)-rr(il,1))
201 fu(il,1)=fu(il,1) &
202 +0.01*grav*work(il)*ment(il,j,1)*(uent(il,j,1)-u(il,1))
203 fv(il,1)=fv(il,1) &
204 +0.01*grav*work(il)*ment(il,j,1)*(vent(il,j,1)-v(il,1))
205 else ! cvflag_grav
206 fr(il,1)=fr(il,1) &
207 +0.1*work(il)*ment(il,j,1)*(qent(il,j,1)-rr(il,1))
208 fu(il,1)=fu(il,1) &
209 +0.1*work(il)*ment(il,j,1)*(uent(il,j,1)-u(il,1))
210 fv(il,1)=fv(il,1) &
211 +0.1*work(il)*ment(il,j,1)*(vent(il,j,1)-v(il,1))
212 endif ! cvflag_grav
213 endif ! j
214 enddo
215 enddo
216
217 !
218 ! *** calculate tendencies of potential temperature and mixing ratio ***
219 ! *** at levels above the lowest level ***
220 !
221 ! *** first find the net saturated updraft and downdraft mass fluxes ***
222 ! *** through each level ***
223 !
224
225 do i=2,nl+1 ! newvecto: mettre nl au lieu nl+1?
226
227 num1=0
228 do il=1,ncum
229 if(i.le.inb(il))num1=num1+1
230 enddo
231 if(num1.le.0) cycle
232
233 amp1(:ncum) = 0.
234 ad(:ncum) = 0.
235
236 do k=i+1,nl+1
237 do il=1,ncum
238 if (i.le.inb(il) .and. k.le.(inb(il)+1)) then
239 amp1(il)=amp1(il)+m(il,k)
240 endif
241 end do
242 end do
243
244 do k=1,i
245 do j=i+1,nl+1
246 do il=1,ncum
247 if (i.le.inb(il) .and. j.le.(inb(il)+1)) then
248 amp1(il)=amp1(il)+ment(il,k,j)
249 endif
250 end do
251 end do
252 end do
253
254 do k=1,i-1
255 do j=i,nl+1 ! newvecto: nl au lieu nl+1?
256 do il=1,ncum
257 if (i.le.inb(il) .and. j.le.inb(il)) then
258 ad(il)=ad(il)+ment(il,j,k)
259 endif
260 end do
261 end do
262 end do
263
264 do il=1,ncum
265 if (i.le.inb(il)) then
266 dpinv=1.0/(ph(il,i)-ph(il,i+1))
267 cpinv=1.0/cpn(il,i)
268
269 ! convect3 if((0.1*dpinv*amp1).ge.delti)iflag(il)=4
270 if (cvflag_grav) then
271 if((0.01*grav*dpinv*amp1(il)).ge.delti)iflag(il)=1 ! vecto
272 else
273 if((0.1*dpinv*amp1(il)).ge.delti)iflag(il)=1 ! vecto
274 endif
275
276 if (cvflag_grav) then
277 ft(il,i)=0.01*grav*dpinv*(amp1(il)*(t(il,i+1)-t(il,i) &
278 +(gz(il,i+1)-gz(il,i))*cpinv) &
279 -ad(il)*(t(il,i)-t(il,i-1)+(gz(il,i)-gz(il,i-1))*cpinv)) &
280 -0.5*sigd*lvcp(il,i)*(evap(il,i)+evap(il,i+1))
281 rat=cpn(il,i-1)*cpinv
282 ft(il,i)=ft(il,i)-0.009*grav*sigd*(mp(il,i+1)*t(il,i)*b(il,i) &
283 -mp(il,i)*t(il,i-1)*rat*b(il,i-1))*dpinv
284 ft(il,i)=ft(il,i)+0.01*grav*dpinv*ment(il,i,i)*(hp(il,i)-h(il,i) &
285 +t(il,i)*(cpv-cpd)*(rr(il,i)-qent(il,i,i)))*cpinv
286 else ! cvflag_grav
287 ft(il,i)=0.1*dpinv*(amp1(il)*(t(il,i+1)-t(il,i) &
288 +(gz(il,i+1)-gz(il,i))*cpinv) &
289 -ad(il)*(t(il,i)-t(il,i-1)+(gz(il,i)-gz(il,i-1))*cpinv)) &
290 -0.5*sigd*lvcp(il,i)*(evap(il,i)+evap(il,i+1))
291 rat=cpn(il,i-1)*cpinv
292 ft(il,i)=ft(il,i)-0.09*sigd*(mp(il,i+1)*t(il,i)*b(il,i) &
293 -mp(il,i)*t(il,i-1)*rat*b(il,i-1))*dpinv
294 ft(il,i)=ft(il,i)+0.1*dpinv*ment(il,i,i)*(hp(il,i)-h(il,i) &
295 +t(il,i)*(cpv-cpd)*(rr(il,i)-qent(il,i,i)))*cpinv
296 endif ! cvflag_grav
297
298
299 ft(il,i)=ft(il,i)+0.01*sigd*wt(il,i)*(cl-cpd)*water(il,i+1) &
300 *(t(il,i+1)-t(il,i))*dpinv*cpinv
301
302 if (cvflag_grav) then
303 fr(il,i)=0.01*grav*dpinv*(amp1(il)*(rr(il,i+1)-rr(il,i)) &
304 -ad(il)*(rr(il,i)-rr(il,i-1)))
305 fu(il,i)=fu(il,i)+0.01*grav*dpinv*(amp1(il)*(u(il,i+1)-u(il,i)) &
306 -ad(il)*(u(il,i)-u(il,i-1)))
307 fv(il,i)=fv(il,i)+0.01*grav*dpinv*(amp1(il)*(v(il,i+1)-v(il,i)) &
308 -ad(il)*(v(il,i)-v(il,i-1)))
309 else ! cvflag_grav
310 fr(il,i)=0.1*dpinv*(amp1(il)*(rr(il,i+1)-rr(il,i)) &
311 -ad(il)*(rr(il,i)-rr(il,i-1)))
312 fu(il,i)=fu(il,i)+0.1*dpinv*(amp1(il)*(u(il,i+1)-u(il,i)) &
313 -ad(il)*(u(il,i)-u(il,i-1)))
314 fv(il,i)=fv(il,i)+0.1*dpinv*(amp1(il)*(v(il,i+1)-v(il,i)) &
315 -ad(il)*(v(il,i)-v(il,i-1)))
316 endif ! cvflag_grav
317
318 endif ! i
319 end do
320
321 do k=1,i-1
322 do il=1,ncum
323 if (i.le.inb(il)) then
324 dpinv=1.0/(ph(il,i)-ph(il,i+1))
325 cpinv=1.0/cpn(il,i)
326
327 awat=elij(il,k,i)-(1.-ep(il,i))*clw(il,i)
328 awat=amax1(awat,0.0)
329
330 if (cvflag_grav) then
331 fr(il,i)=fr(il,i) &
332 +0.01*grav*dpinv*ment(il,k,i)*(qent(il,k,i)-awat-rr(il,i))
333 fu(il,i)=fu(il,i) &
334 +0.01*grav*dpinv*ment(il,k,i)*(uent(il,k,i)-u(il,i))
335 fv(il,i)=fv(il,i) &
336 +0.01*grav*dpinv*ment(il,k,i)*(vent(il,k,i)-v(il,i))
337 else ! cvflag_grav
338 fr(il,i)=fr(il,i) &
339 +0.1*dpinv*ment(il,k,i)*(qent(il,k,i)-awat-rr(il,i))
340 fu(il,i)=fu(il,i) &
341 +0.01*grav*dpinv*ment(il,k,i)*(uent(il,k,i)-u(il,i))
342 fv(il,i)=fv(il,i) &
343 +0.1*dpinv*ment(il,k,i)*(vent(il,k,i)-v(il,i))
344 endif ! cvflag_grav
345
346 ! (saturated updrafts resulting from mixing) ! cld
347 qcond(il,i)=qcond(il,i)+(elij(il,k,i)-awat) ! cld
348 nqcond(il,i)=nqcond(il,i)+1. ! cld
349 endif ! i
350 end do
351 end do
352
353 do k=i,nl+1
354 do il=1,ncum
355 if (i.le.inb(il) .and. k.le.inb(il)) then
356 dpinv=1.0/(ph(il,i)-ph(il,i+1))
357 cpinv=1.0/cpn(il,i)
358
359 if (cvflag_grav) then
360 fr(il,i)=fr(il,i) &
361 +0.01*grav*dpinv*ment(il,k,i)*(qent(il,k,i)-rr(il,i))
362 fu(il,i)=fu(il,i) &
363 +0.01*grav*dpinv*ment(il,k,i)*(uent(il,k,i)-u(il,i))
364 fv(il,i)=fv(il,i) &
365 +0.01*grav*dpinv*ment(il,k,i)*(vent(il,k,i)-v(il,i))
366 else ! cvflag_grav
367 fr(il,i)=fr(il,i) &
368 +0.1*dpinv*ment(il,k,i)*(qent(il,k,i)-rr(il,i))
369 fu(il,i)=fu(il,i) &
370 +0.1*dpinv*ment(il,k,i)*(uent(il,k,i)-u(il,i))
371 fv(il,i)=fv(il,i) &
372 +0.1*dpinv*ment(il,k,i)*(vent(il,k,i)-v(il,i))
373 endif ! cvflag_grav
374 endif ! i and k
375 end do
376 end do
377
378 do il=1,ncum
379 if (i.le.inb(il)) then
380 dpinv=1.0/(ph(il,i)-ph(il,i+1))
381 cpinv=1.0/cpn(il,i)
382
383 if (cvflag_grav) then
384 ! sb: on ne fait pas encore la correction permettant de mieux
385 ! conserver l'eau:
386 fr(il,i)=fr(il,i)+0.5*sigd*(evap(il,i)+evap(il,i+1)) &
387 +0.01*grav*(mp(il,i+1)*(rp(il,i+1)-rr(il,i))-mp(il,i) &
388 *(rp(il,i)-rr(il,i-1)))*dpinv
389
390 fu(il,i)=fu(il,i)+0.01*grav*(mp(il,i+1)*(up(il,i+1)-u(il,i)) &
391 -mp(il,i)*(up(il,i)-u(il,i-1)))*dpinv
392 fv(il,i)=fv(il,i)+0.01*grav*(mp(il,i+1)*(vp(il,i+1)-v(il,i)) &
393 -mp(il,i)*(vp(il,i)-v(il,i-1)))*dpinv
394 else ! cvflag_grav
395 fr(il,i)=fr(il,i)+0.5*sigd*(evap(il,i)+evap(il,i+1)) &
396 +0.1*(mp(il,i+1)*(rp(il,i+1)-rr(il,i))-mp(il,i) &
397 *(rp(il,i)-rr(il,i-1)))*dpinv
398 fu(il,i)=fu(il,i)+0.1*(mp(il,i+1)*(up(il,i+1)-u(il,i)) &
399 -mp(il,i)*(up(il,i)-u(il,i-1)))*dpinv
400 fv(il,i)=fv(il,i)+0.1*(mp(il,i+1)*(vp(il,i+1)-v(il,i)) &
401 -mp(il,i)*(vp(il,i)-v(il,i-1)))*dpinv
402 endif ! cvflag_grav
403
404 endif ! i
405 end do
406
407 ! sb: interface with the cloud parameterization: ! cld
408
409 do k=i+1,nl
410 do il=1,ncum
411 if (k.le.inb(il) .and. i.le.inb(il)) then ! cld
412 ! (saturated downdrafts resulting from mixing) ! cld
413 qcond(il,i)=qcond(il,i)+elij(il,k,i) ! cld
414 nqcond(il,i)=nqcond(il,i)+1. ! cld
415 endif ! cld
416 enddo ! cld
417 enddo ! cld
418
419 ! (particular case: no detraining level is found) ! cld
420 do il=1,ncum ! cld
421 if (i.le.inb(il) .and. nent(il,i).eq.0) then ! cld
422 qcond(il,i)=qcond(il,i)+(1.-ep(il,i))*clw(il,i) ! cld
423 nqcond(il,i)=nqcond(il,i)+1. ! cld
424 endif ! cld
425 enddo ! cld
426
427 do il=1,ncum ! cld
428 if (i.le.inb(il) .and. nqcond(il,i).ne.0.) then ! cld
429 qcond(il,i)=qcond(il,i)/nqcond(il,i) ! cld
430 endif ! cld
431 enddo
432
433 end do
434
435
436 ! *** move the detrainment at level inb down to level inb-1 ***
437 ! *** in such a way as to preserve the vertically ***
438 ! *** integrated enthalpy and water tendencies ***
439 !
440 do il=1,ncum
441
442 ax=0.1*ment(il,inb(il),inb(il))*(hp(il,inb(il))-h(il,inb(il)) &
443 +t(il,inb(il))*(cpv-cpd) &
444 *(rr(il,inb(il))-qent(il,inb(il),inb(il)))) &
445 /(cpn(il,inb(il))*(ph(il,inb(il))-ph(il,inb(il)+1)))
446 ft(il,inb(il))=ft(il,inb(il))-ax
447 ft(il,inb(il)-1)=ft(il,inb(il)-1)+ax*cpn(il,inb(il)) &
448 *(ph(il,inb(il))-ph(il,inb(il)+1))/(cpn(il,inb(il)-1) &
449 *(ph(il,inb(il)-1)-ph(il,inb(il))))
450
451 bx=0.1*ment(il,inb(il),inb(il))*(qent(il,inb(il),inb(il)) &
452 -rr(il,inb(il)))/(ph(il,inb(il))-ph(il,inb(il)+1))
453 fr(il,inb(il))=fr(il,inb(il))-bx
454 fr(il,inb(il)-1)=fr(il,inb(il)-1) &
455 +bx*(ph(il,inb(il))-ph(il,inb(il)+1)) &
456 /(ph(il,inb(il)-1)-ph(il,inb(il)))
457
458 cx=0.1*ment(il,inb(il),inb(il))*(uent(il,inb(il),inb(il)) &
459 -u(il,inb(il)))/(ph(il,inb(il))-ph(il,inb(il)+1))
460 fu(il,inb(il))=fu(il,inb(il))-cx
461 fu(il,inb(il)-1)=fu(il,inb(il)-1) &
462 +cx*(ph(il,inb(il))-ph(il,inb(il)+1)) &
463 /(ph(il,inb(il)-1)-ph(il,inb(il)))
464
465 dx=0.1*ment(il,inb(il),inb(il))*(vent(il,inb(il),inb(il)) &
466 -v(il,inb(il)))/(ph(il,inb(il))-ph(il,inb(il)+1))
467 fv(il,inb(il))=fv(il,inb(il))-dx
468 fv(il,inb(il)-1)=fv(il,inb(il)-1) &
469 +dx*(ph(il,inb(il))-ph(il,inb(il)+1)) &
470 /(ph(il,inb(il)-1)-ph(il,inb(il)))
471
472 end do
473
474 !
475 ! *** homoginize tendencies below cloud base ***
476 !
477 !
478 do il=1,ncum
479 asum(il)=0.0
480 bsum(il)=0.0
481 csum(il)=0.0
482 dsum(il)=0.0
483 enddo
484
485 do i=1,nl
486 do il=1,ncum
487 if (i.le.(icb(il)-1)) then
488 asum(il)=asum(il)+ft(il,i)*(ph(il,i)-ph(il,i+1))
489 bsum(il)=bsum(il)+fr(il,i)*(lv(il,i)+(cl-cpd)*(t(il,i)-t(il,1))) &
490 *(ph(il,i)-ph(il,i+1))
491 csum(il)=csum(il)+(lv(il,i)+(cl-cpd)*(t(il,i)-t(il,1))) &
492 *(ph(il,i)-ph(il,i+1))
493 dsum(il)=dsum(il)+t(il,i)*(ph(il,i)-ph(il,i+1))/th(il,i)
494 endif
495 enddo
496 enddo
497
498 !!!! do 700 i=1,icb(il)-1
499 do i=1,nl
500 do il=1,ncum
501 if (i.le.(icb(il)-1)) then
502 ft(il,i)=asum(il)*t(il,i)/(th(il,i)*dsum(il))
503 fr(il,i)=bsum(il)/csum(il)
504 endif
505 enddo
506 enddo
507
508 !
509 ! *** reset counter and return ***
510 !
511 do il=1,ncum
512 sig(il,nd)=2.0
513 enddo
514
515
516 do i=1,nd
517 do il=1,ncum
518 upwd(il,i)=0.0
519 dnwd(il,i)=0.0
520 enddo
521 enddo
522
523 do i=1,nl
524 do il=1,ncum
525 dnwd0(il,i)=-mp(il,i)
526 enddo
527 enddo
528 do i=nl+1,nd
529 do il=1,ncum
530 dnwd0(il,i)=0.
531 enddo
532 enddo
533
534
535 do i=1,nl
536 do il=1,ncum
537 if (i.ge.icb(il) .and. i.le.inb(il)) then
538 upwd(il,i)=0.0
539 dnwd(il,i)=0.0
540 endif
541 enddo
542 enddo
543
544 do i=1,nl
545 do k=1,nl
546 do il=1,ncum
547 up1(il,k,i)=0.0
548 dn1(il,k,i)=0.0
549 enddo
550 enddo
551 enddo
552
553 do i=1,nl
554 do k=i,nl
555 do n=1,i-1
556 do il=1,ncum
557 if (i.ge.icb(il).and.i.le.inb(il).and.k.le.inb(il)) then
558 up1(il,k,i)=up1(il,k,i)+ment(il,n,k)
559 dn1(il,k,i)=dn1(il,k,i)-ment(il,k,n)
560 endif
561 enddo
562 enddo
563 enddo
564 enddo
565
566 do i=2,nl
567 do k=i,nl
568 do il=1,ncum
569 !test if (i.ge.icb(il).and.i.le.inb(il).and.k.le.inb(il)) then
570 if (i.le.inb(il).and.k.le.inb(il)) then
571 upwd(il,i)=upwd(il,i)+m(il,k)+up1(il,k,i)
572 dnwd(il,i)=dnwd(il,i)+dn1(il,k,i)
573 endif
574 enddo
575 enddo
576 enddo
577
578
579 !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
580 ! determination de la variation de flux ascendant entre
581 ! deux niveau non dilue mike
582 !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
583
584 do i=1,nl
585 do il=1,ncum
586 mike(il,i)=m(il,i)
587 enddo
588 enddo
589
590 do i=nl+1,nd
591 do il=1,ncum
592 mike(il,i)=0.
593 enddo
594 enddo
595
596 do i=1,nd
597 do il=1,ncum
598 ma(il,i)=0
599 enddo
600 enddo
601
602 do i=1,nl
603 do j=i,nl
604 do il=1,ncum
605 ma(il,i)=ma(il,i)+m(il,j)
606 enddo
607 enddo
608 enddo
609
610 do i=nl+1,nd
611 do il=1,ncum
612 ma(il,i)=0.
613 enddo
614 enddo
615
616 do i=1,nl
617 do il=1,ncum
618 if (i.le.(icb(il)-1)) then
619 ma(il,i)=0
620 endif
621 enddo
622 enddo
623
624 !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
625 ! icb represente de niveau ou se trouve la
626 ! base du nuage , et inb le top du nuage
627 !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
628
629 do i=1,nd
630 DO il=1,ncum
631 rdcp=(rrd*(1.-rr(il,i))-rr(il,i)*rrv) &
632 /(cpd*(1.-rr(il,i))+rr(il,i)*cpv)
633 tls(il,i)=t(il,i)*(1000.0/p(il,i))**rdcp
634 tps(il,i)=tp(il,i)
635 end DO
636 enddo
637
638 !
639 ! *** diagnose the in-cloud mixing ratio *** ! cld
640 ! *** of condensed water *** ! cld
641 ! ! cld
642
643 do i=1,nd ! cld
644 do il=1,ncum ! cld
645 mac(il,i)=0.0 ! cld
646 wa(il,i)=0.0 ! cld
647 siga(il,i)=0.0 ! cld
648 sax(il,i)=0.0 ! cld
649 enddo ! cld
650 enddo ! cld
651
652 do i=minorig, nl ! cld
653 do k=i+1,nl+1 ! cld
654 do il=1,ncum ! cld
655 if (i.le.inb(il) .and. k.le.(inb(il)+1)) then ! cld
656 mac(il,i)=mac(il,i)+m(il,k) ! cld
657 endif ! cld
658 enddo ! cld
659 enddo ! cld
660 enddo ! cld
661
662 do i=1,nl ! cld
663 do j=1,i ! cld
664 do il=1,ncum ! cld
665 if (i.ge.icb(il) .and. i.le.(inb(il)-1) &
666 .and. j.ge.icb(il) ) then ! cld
667 sax(il,i)=sax(il,i)+rrd*(tvp(il,j)-tv(il,j)) &
668 *(ph(il,j)-ph(il,j+1))/p(il,j) ! cld
669 endif ! cld
670 enddo ! cld
671 enddo ! cld
672 enddo ! cld
673
674 do i=1,nl ! cld
675 do il=1,ncum ! cld
676 if (i.ge.icb(il) .and. i.le.(inb(il)-1) &
677 .and. sax(il,i).gt.0.0 ) then ! cld
678 wa(il,i)=sqrt(2.*sax(il,i)) ! cld
679 endif ! cld
680 enddo ! cld
681 enddo ! cld
682
683 do i=1,nl ! cld
684 do il=1,ncum ! cld
685 if (wa(il,i).gt.0.0) &
686 siga(il,i)=mac(il,i)/wa(il,i) &
687 *rrd*tvp(il,i)/p(il,i)/100./delta ! cld
688 siga(il,i) = min(siga(il,i),1.0) ! cld
689 !IM cf. FH
690 if (iflag_clw.eq.0) then
691 qcondc(il,i)=siga(il,i)*clw(il,i)*(1.-ep(il,i)) &
692 + (1.-siga(il,i))*qcond(il,i) ! cld
693 else if (iflag_clw.eq.1) then
694 qcondc(il,i)=qcond(il,i) ! cld
695 endif
696
697 enddo ! cld
698 enddo ! cld
699
700 end SUBROUTINE cv30_yield
701
702 end module cv30_yield_m

  ViewVC Help
Powered by ViewVC 1.1.21