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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 201 - (hide annotations)
Mon Jun 6 17:42:15 2016 UTC (7 years, 11 months ago) by guez
File size: 19375 byte(s) 
Removed intermediary objects of cv_thermo_m, access suphec_m
directly. Procedure cv_thermo disappeared, all objects are named
constants.

In cv_driver and below, limited extents of arrays to what is needed.

lv, cpn and th in cv30_compress were set at level nl + 1 but lv1, cpn1
and th1 are not defined at this level. This did not lead to an error
because values at nl + 1 were not used.

Removed test on ok_sync in phystokenc because it is not read at run
time. Printing min and max of output NetCDF variables is heavy and
archaic.

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