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

  ViewVC Help
Powered by ViewVC 1.1.21