/[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 195 - (show annotations)
Wed May 18 17:56:44 2016 UTC (7 years, 11 months ago) by guez
File size: 19697 byte(s)
In cv30_feed, iflag1 is 0 on entry so we can simplify the test for
iflag1 = 7.

In cv30_feed, for the computation of icb, replaced sequential search
(with a useless end of loop on k) by a call to locate.

In CV30 routines, replaced len, nloc, nd, na by klon or
klev. Philosophy: no more generality than actually necessary.

Converted as many variables as possible to named constants in
cv30_param_m and downgraded pbcrit, ptcrit, dtovsh, dpbase, dttrig,
tau, delta to local objects in procedures. spfac, betad and omtrain
are useless and removed.

Instead of filling the array sigp with the constant spfac in
cv30_undilute2, just made sigp a constant in cv30_unsat.

In cv_driver, define as allocatable variables that are only
used on the range (ncum, nl).

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

  ViewVC Help
Powered by ViewVC 1.1.21