1 | MODULE module_cu_camuwshcu_driver |
---|
2 | USE shr_kind_mod, only: r8 => shr_kind_r8 |
---|
3 | |
---|
4 | ! Roughly based on convect_shallow_tend in convect_shallow.F90 from CAM |
---|
5 | ! but tailored for the UW shallow cumulus scheme. |
---|
6 | |
---|
7 | IMPLICIT NONE |
---|
8 | |
---|
9 | PRIVATE !Default to private |
---|
10 | PUBLIC :: & !Public entities |
---|
11 | camuwshcu_driver |
---|
12 | |
---|
13 | CONTAINS |
---|
14 | |
---|
15 | !------------------------------------------------------------------------ |
---|
16 | SUBROUTINE camuwshcu_driver( & |
---|
17 | ids,ide, jds,jde, kds,kde & |
---|
18 | ,ims,ime, jms,jme, kms,kme & |
---|
19 | ,its,ite, jts,jte, kts,kte & |
---|
20 | ,num_moist, dt & |
---|
21 | ,p, p8w, pi_phy, z, z_at_w, dz8w & |
---|
22 | ,t_phy, u_phy, v_phy & |
---|
23 | ,moist, qv, qc, qi & |
---|
24 | ,pblh_in, tke_pbl, cldfra, cldfra_old, cldfrash & |
---|
25 | ,cush_inout, rainsh, pratesh, snowsh, icwmrsh & |
---|
26 | ,cmfmc, cmfmc2_inout, rprdsh_inout, cbmf_inout & |
---|
27 | ,cmfsl, cmflq, dlf, evapcsh_inout & |
---|
28 | ,rliq, rliq2_inout, cubot, cutop & |
---|
29 | ,rushten, rvshten, rthshten & |
---|
30 | ,rqvshten, rqcshten, rqrshten & |
---|
31 | ,rqishten, rqsshten, rqgshten & |
---|
32 | ,ht & |
---|
33 | ) |
---|
34 | ! This routine is based on convect_shallow_tend in CAM. It handles the |
---|
35 | ! mapping of variables from the WRF to the CAM framework for the UW |
---|
36 | ! shallow convective parameterization. |
---|
37 | ! |
---|
38 | ! Author: William.Gustafson@pnl.gov, Jan. 2010 |
---|
39 | !------------------------------------------------------------------------ |
---|
40 | USE module_state_description, only: param_first_scalar, & |
---|
41 | p_qc, p_qr, p_qi, p_qs, p_qg |
---|
42 | USE module_cam_support, only: pcols, pver |
---|
43 | USE constituents, only: cnst_get_ind |
---|
44 | USE physconst, only: cpair, gravit, latvap |
---|
45 | USE uwshcu, only: compute_uwshcu_inv |
---|
46 | USE wv_saturation, only: fqsatd |
---|
47 | |
---|
48 | ! Subroutine arguments... |
---|
49 | INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, & |
---|
50 | ims,ime, jms,jme, kms,kme, & |
---|
51 | its,ite, jts,jte, kts,kte, & |
---|
52 | num_moist |
---|
53 | |
---|
54 | REAL, DIMENSION( ims:ime, kms:kme, jms:jme, num_moist ), INTENT(IN) :: & |
---|
55 | moist !moist tracer array |
---|
56 | |
---|
57 | REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), INTENT(IN) :: & |
---|
58 | cldfra, & !cloud fraction |
---|
59 | cldfra_old, & !previous time step cloud fraction |
---|
60 | dz8w, & !height between layer interface (m) |
---|
61 | p, & !pressure at mid-level (Pa) |
---|
62 | p8w, & !pressure at level interface (Pa) |
---|
63 | pi_phy, & !exner function, (p0/p)^(R/cpair) (none) |
---|
64 | qv, & !water vapor mixing ratio (kg/kg-dry air) |
---|
65 | t_phy, & !temperature (K) |
---|
66 | tke_pbl, & !turbulent kinetic energy from PBL (m2/s2) |
---|
67 | u_phy, & !zonal wind component on T points (m/s) |
---|
68 | v_phy, & !meridional wind component on T points (m/s) |
---|
69 | z, & !height above sea level at mid-level (m) |
---|
70 | z_at_w !height above sea level at interface (m) |
---|
71 | |
---|
72 | REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), INTENT(IN), OPTIONAL :: & |
---|
73 | qc, & !cloud droplet mixing ratio (kg/kg-dry air) |
---|
74 | qi !cloud ice crystal mixing ratio (kg/kg-dry air) |
---|
75 | |
---|
76 | REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN) :: & |
---|
77 | pblh_in, & !height of PBL (m) |
---|
78 | ht !Terrain height (m) |
---|
79 | |
---|
80 | REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), INTENT(INOUT) :: & |
---|
81 | cldfrash, & !shallow convective cloud fraction |
---|
82 | cmfmc, & !deep+shalow cloud fraction (already contains deep part from ZM) |
---|
83 | cmfmc2_inout, & !shallow cloud fraction |
---|
84 | cmflq, & !convective flux of total water in energy unit (~units) |
---|
85 | cmfsl, & !convective flux of liquid water static energy (~units) |
---|
86 | dlf, & !dq/dt due to export of cloud water (input=deep from ZM, output=deep+shallow) |
---|
87 | evapcsh_inout, & !output array for evaporation of shallow convection precipitation (kg/kg/s) |
---|
88 | icwmrsh, & !shallow cumulus in-cloud water mixing ratio (kg/m2) |
---|
89 | rprdsh_inout, & !dq/dt due to deep(~?) & shallow convective rainout (~units?) |
---|
90 | rushten, & !UNcoupled zonal wind tend from shallow Cu scheme (m/s2) |
---|
91 | rvshten, & !UNcoupled meridional wind tend from shallow Cu scheme (m/s2) |
---|
92 | rthshten, & !UNcoupled potential temperature tendendcy from shallow cu scheme (K/s) |
---|
93 | rqvshten, & !UNcoupled water vapor mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
94 | rqcshten, & !UNcoupled clod droplet mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
95 | rqrshten, & !UNcoupled raindrop mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
96 | rqishten, & !UNcoupled ice crystal mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
97 | rqsshten, & !UNcoupled snow mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
98 | rqgshten !UNcoupled graupel mixing ratio tend from shallow Cu scheme (kg/kg/s) |
---|
99 | |
---|
100 | REAL, DIMENSION( ims:ime, jms:jme ), INTENT(INOUT) :: & |
---|
101 | cbmf_inout, & !cloud base mass flux (kg/m2/s) |
---|
102 | cubot, & !level number of base of convection |
---|
103 | cutop, & !level number of top of convection |
---|
104 | cush_inout, & !convective scale height (~units?) |
---|
105 | pratesh, & !time-step shallow cumulus precip rate at surface (mm/s) |
---|
106 | rainsh, & !time-step shallow cumulus precip (rain+snow) at surface (mm) |
---|
107 | rliq, & !vertically-integrated reserved cloud condensate (m/s) |
---|
108 | rliq2_inout, & !vertically-integrated reserved cloud condensate for shallow (m/s) |
---|
109 | snowsh !accumulated convective snow rate at surface for shallow Cu (m/s) ~are these the units we should use for WRF? |
---|
110 | |
---|
111 | REAL, INTENT(IN) :: & |
---|
112 | dt !time step (s) |
---|
113 | |
---|
114 | ! Local variables... |
---|
115 | !Variables dimensioned for input to CAM routines |
---|
116 | REAL(r8), DIMENSION(pcols, kte, 5) :: & ! dimensioned by 5(=ncnst) associated with CAM constituents (cnst_name size) |
---|
117 | moist8, & !tracer array for CAM routines |
---|
118 | tnd_tracer !tracer tendency |
---|
119 | |
---|
120 | REAL(r8), DIMENSION(pcols, kte+1) :: & |
---|
121 | pint8, & !pressure at layer interface (Pa) |
---|
122 | zi8, & !height above the ground at interfaces (m) |
---|
123 | tke8, & !turbulent kinetic energy at level interfaces (m2/s2) |
---|
124 | slflx, & !convective liquid water static energy flux (~units?) |
---|
125 | qtflx !convective total water flux (~units?) |
---|
126 | |
---|
127 | |
---|
128 | |
---|
129 | REAL(r8), DIMENSION(pcols, kte) :: & |
---|
130 | cld8, & !cloud fraction |
---|
131 | cldold8, & !previous time step cloud fraction ~should this be just the convective part? |
---|
132 | cmfdqs, & !convective snow production (~units?) |
---|
133 | cmfmc2, & !cloud fraction |
---|
134 | evapcsh, & !evaporation of shallow convection precipitation >= 0. (kg/kg/s) |
---|
135 | iccmr_uw, & !in-cloud cumulus LWC+IWC (kg/m2) |
---|
136 | icwmr_uw, & !in-cloud cumulus LWC (kg/m2) |
---|
137 | icimr_uw, & !in-cloud cumulus IWC (kg/m2) |
---|
138 | pdel8, & !pressure difference between layer interfaces (Pa) |
---|
139 | pdeldry8, & !pressure difference between layer interfaces for dry atm (Pa) |
---|
140 | pmid8, & !pressure at layer middle (Pa) |
---|
141 | qc2, & !dq/dt due to export of cloud water |
---|
142 | qh8, & !specific humidity (kg/kg-moist air) |
---|
143 | qc8, & !cloud liquid water (~units?) |
---|
144 | qi8, & !cloud ice (~units?) |
---|
145 | qhtnd, & !specific humidity tendency (kg/kg/s) |
---|
146 | qctnd, & !cloud water mixing ratio tendency |
---|
147 | qitnd, & !cloud ice mixing ratio tendency |
---|
148 | rprdsh, & !dq/dt due to deep(~?) & shallow convective rainout (~units?) |
---|
149 | s8, & !dry static energy (J/kg) |
---|
150 | shfrc, & !shallow cloud fraction |
---|
151 | stnd, & !heating rate (dry static energy tendency, W/kg) |
---|
152 | t8, & !temperature (K) |
---|
153 | u8, & !environment zonal wind (m/s) |
---|
154 | utnd, & !zonal wind tendency (m/s2) |
---|
155 | v8, & !environment meridional wind (m/s) |
---|
156 | vtnd, & !meridional wind tendency (m/s2) |
---|
157 | zm8 !height between interfaces (m) |
---|
158 | |
---|
159 | REAL(r8), DIMENSION(pcols) :: & |
---|
160 | cbmf, & !cloud base mass flux (kg/m2/s) |
---|
161 | cnb2, & !bottom level of convective activity |
---|
162 | cnt2, & !top level of convective activity |
---|
163 | cush, & !convective scale height (~units?) |
---|
164 | pblh, & !pblh height (m) |
---|
165 | precc, & !convective precip (rain+snow) at surface for shallow Cu (m/s) |
---|
166 | rliq2, & !vertically-integrated reserved cloud condensate for shallow (m/s) |
---|
167 | snow !convective snow rate at surface (m/s) |
---|
168 | |
---|
169 | !Other local vars |
---|
170 | REAL(r8) :: ztodt !2 delta-t (s) |
---|
171 | INTEGER :: i, j, k, kflip, m, mp1 |
---|
172 | INTEGER :: cnb, cnt !index of cloud base and top in CAM world (indices decrease with height) |
---|
173 | INTEGER :: lchnk !chunk identifier, used to map 2-D to 1-D arrays in WRF |
---|
174 | INTEGER :: ncnst !number of tracers |
---|
175 | INTEGER :: ncol !number of atmospheric columns in chunk |
---|
176 | CHARACTER(LEN=250) :: msg |
---|
177 | ! |
---|
178 | ! Initialize... |
---|
179 | ! |
---|
180 | ncol = 1 !chunk size in WRF is 1 since we loop over all columns in a tile |
---|
181 | ! ncnst = num_moist-1 !currently we only handle the moist array for tracers |
---|
182 | ncnst = 5 ! this is associated with moist8 array |
---|
183 | ztodt = 2.*dt |
---|
184 | ! |
---|
185 | ! Map variables to inputs for zm_convr and call it... |
---|
186 | ! Loop over the points in the tile and treat them each as a CAM chunk. |
---|
187 | ! |
---|
188 | ij_loops : do j = jts,jte |
---|
189 | do i = its,ite |
---|
190 | lchnk = (j-jts)*(ite-its+1) + (i-its+1) !1-D index location from the 2-D tile |
---|
191 | |
---|
192 | !Flip variables on the layer interfaces |
---|
193 | do k = kts,kte+1 |
---|
194 | kflip = kte-k+2 |
---|
195 | |
---|
196 | pint8(1,kflip) = p8w(i,k,j) |
---|
197 | zi8(1,kflip) = z_at_w(i,k,j) - ht(i,j) ! height above the ground at interfaces |
---|
198 | end do |
---|
199 | |
---|
200 | !Flip variables on the layer middles |
---|
201 | do k = kts,kte |
---|
202 | kflip = kte-k+1 |
---|
203 | |
---|
204 | cld8(1,kflip) = cldfra(i,k,j) |
---|
205 | cldold8(1,kflip) = cldfra_old(i,k,j) |
---|
206 | pdel8(1,kflip) = p8w(i,k,j) - p8w(i,k+1,j) |
---|
207 | pmid8(1,kflip) = p(i,k,j) |
---|
208 | qh8(1,kflip) = max( qv(i,k,j)/(1. + qv(i,k,j)), 1e-30 ) !values of 0 cause a crash in entropy |
---|
209 | if( present(qc) ) then |
---|
210 | qc8(1,kflip) = qc(i,k,j)/(1. + qv(i,k,j)) !Convert to moist mix ratio |
---|
211 | else |
---|
212 | qc8(1,kflip) = 0. |
---|
213 | end if |
---|
214 | if( present(qi) ) then |
---|
215 | qi8(1,kflip) = qi(i,k,j)/(1. + qv(i,k,j)) !Used in convtran, ditto for conversion |
---|
216 | else |
---|
217 | qi8(1,kflip) = 0. |
---|
218 | end if |
---|
219 | pdeldry8(1,kflip)= pdel8(1,kflip)*(1._r8 - qh8(1,kflip)) |
---|
220 | t8(1,kflip) = t_phy(i,k,j) |
---|
221 | s8(1,kflip) = cpair*t8(1,kflip) + gravit*(z(i,k,j)-ht(i,j)) |
---|
222 | u8(1,kflip) = u_phy(i,k,j) |
---|
223 | v8(1,kflip) = v_phy(i,k,j) |
---|
224 | zm8(1,kflip) = dz8w(i,k,j) |
---|
225 | end do |
---|
226 | |
---|
227 | !TKE in CAM is on the interfaces but in WRF it is on the layer |
---|
228 | !middle. We will interpolate the TKE from the to the interfaces |
---|
229 | !and then just use the lowest TKE for the surface and the highest |
---|
230 | !TKE at the top. |
---|
231 | tke8(1,kte+1) = tke_pbl(i,1,j) !surface |
---|
232 | tke8(1,1) = tke_pbl(i,kte,j) !model top interface |
---|
233 | do k = kts,kte-1 |
---|
234 | kflip = kte-k+1 |
---|
235 | tke8(1,kflip) = 0.5*(tke_pbl(i,k,j) + tke_pbl(i,k+1,j)) |
---|
236 | end do |
---|
237 | |
---|
238 | !Flip the tracer array - |
---|
239 | !shift tracer dimension down one to remove "blank" index and |
---|
240 | !convert to wet instead of dry mixing ratios. |
---|
241 | do k = kts,kte |
---|
242 | kflip = kte-k+1 |
---|
243 | !!$ do m = 1,ncnst |
---|
244 | !!$ moist8(1,kflip,m) = moist(i,k,j,m+1)/(1. + qv(i,k,j)) |
---|
245 | !!$ end do |
---|
246 | |
---|
247 | !~For now, send replicate part of the tracer array and send zeros for the |
---|
248 | ! rest since CAM treats condensate diagnostically compared to WRF that is |
---|
249 | ! prognostic. This avoids issues with hard-wired assumptions in the UW |
---|
250 | ! ShCu scheme. This should be looked at again when more time is available. |
---|
251 | ! Set to zero for most then overwrite as needed... |
---|
252 | moist8(1,kflip,1:ncnst) = 0. |
---|
253 | |
---|
254 | moist8(1,kflip,1) = qv(i,k,j)/(1. + qv(i,k,j)) |
---|
255 | |
---|
256 | call cnst_get_ind( 'CLDLIQ', m ) |
---|
257 | moist8(1,kflip,m) = qc(i,k,j)/(1. + qv(i,k,j)) |
---|
258 | |
---|
259 | call cnst_get_ind( 'CLDICE', m ) |
---|
260 | moist8(1,kflip,m) = qi(i,k,j)/(1. + qv(i,k,j)) |
---|
261 | |
---|
262 | call cnst_get_ind( 'NUMLIQ', m ) |
---|
263 | moist8(1,kflip,m) = 0. |
---|
264 | |
---|
265 | call cnst_get_ind( 'NUMICE', m ) |
---|
266 | moist8(1,kflip,m) = 0. |
---|
267 | end do |
---|
268 | |
---|
269 | !Some remapping to get arrays to pass into the routine |
---|
270 | pblh(1) = pblh_in(i,j) |
---|
271 | cush(1) = cush_inout(i,j) |
---|
272 | ! |
---|
273 | ! Main guts of the routine... |
---|
274 | ! This is a bit inefficient because we are flippling the arrays and they |
---|
275 | ! will then get flipped back again by compute_uwshcu_inv. We are doing |
---|
276 | ! this to preserve the CAM code as much as possible for maintenance. |
---|
277 | ! |
---|
278 | call compute_uwshcu_inv( & |
---|
279 | pcols, pver, ncol, ncnst, ztodt, & |
---|
280 | pint8, zi8, pmid8, zm8, pdel8, & |
---|
281 | u8, v8, qh8, qc8, qi8, & |
---|
282 | t8, s8, moist8, & |
---|
283 | tke8, cld8, cldold8, pblh, cush, & |
---|
284 | cmfmc2, slflx, qtflx, & |
---|
285 | qhtnd, qctnd, qitnd, & |
---|
286 | stnd, utnd, vtnd, tnd_tracer, & |
---|
287 | rprdsh, cmfdqs, precc, snow, & |
---|
288 | evapcsh, shfrc, iccmr_UW, icwmr_UW, & |
---|
289 | icimr_UW, cbmf, qc2, rliq2, & |
---|
290 | cnt2, cnb2, fqsatd, lchnk, pdeldry8 ) |
---|
291 | ! |
---|
292 | ! Map output into WRF-dimensioned arrays... |
---|
293 | ! |
---|
294 | cush_inout(i,j) = cush(1) |
---|
295 | |
---|
296 | do k = kts,kte |
---|
297 | kflip = kte-k+1 |
---|
298 | |
---|
299 | !Add shallow reserved cloud condensate to deep reserved cloud condensate |
---|
300 | ! dlf (kg/kg/s, qc in CAM), rliq done below |
---|
301 | dlf(i,k,j) = dlf(i,k,j) + qc2(1,kflip) |
---|
302 | |
---|
303 | evapcsh_inout(i,k,j)= evapcsh(1,kflip) |
---|
304 | icwmrsh(i,k,j) = icwmr_uw(1,kflip) |
---|
305 | |
---|
306 | rprdsh(1,kflip) = rprdsh(1,kflip) + cmfdqs(1,kflip) |
---|
307 | rprdsh_inout(i,k,j) = rprdsh(1,kflip) |
---|
308 | !Not doing rprdtot for now since not yet used by other CAM routines in WRF |
---|
309 | |
---|
310 | !Tendencies of winds, potential temperature, and moisture |
---|
311 | !fields treated specifically by UW scheme |
---|
312 | rushten(i,k,j) = utnd(1,kflip) |
---|
313 | rvshten(i,k,j) = vtnd(1,kflip) |
---|
314 | rthshten(i,k,j) = stnd(1,kflip)/cpair/pi_phy(i,k,j) |
---|
315 | rqvshten(i,k,j) = qhtnd(1,kflip)/(1. - qv(i,k,j)) |
---|
316 | if( p_qc >= param_first_scalar ) & |
---|
317 | rqcshten(i,k,j) = qctnd(1,kflip)/(1. - qv(i,k,j)) |
---|
318 | if( p_qi >= param_first_scalar ) & |
---|
319 | rqishten(i,k,j) = qitnd(1,kflip)/(1. - qv(i,k,j)) |
---|
320 | |
---|
321 | !~Turn off tendencies for most condensates since CAM treats them diagnostically. |
---|
322 | ! !Tendencies of tracers except qv,qc,qi |
---|
323 | !!~need to make sure qg tendency is propagated through to application |
---|
324 | ! do m = 4,ncnst |
---|
325 | ! mp1 = m+1 !shift to p_ value for the tracer |
---|
326 | ! if( mp1==p_qr ) then |
---|
327 | ! rqrshten(i,k,j) = tnd_tracer(1,kflip,m)/(1. - qv(i,k,j)) |
---|
328 | ! else if( mp1==p_qs ) then |
---|
329 | ! rqsshten(i,k,j) = tnd_tracer(1,kflip,m)/(1. - qv(i,k,j)) |
---|
330 | ! else if( mp1==p_qg ) then |
---|
331 | ! rqgshten(i,k,j) = tnd_tracer(1,kflip,m)/(1. - qv(i,k,j)) |
---|
332 | ! else |
---|
333 | ! write(msg,'(a,i3)') "WARNING: UW shallow Cu cannot handle tracer ",m |
---|
334 | ! call wrf_debug(100, msg) |
---|
335 | ! end if |
---|
336 | ! end do |
---|
337 | |
---|
338 | !Combine shallow and deep cumulus updraft mass flux |
---|
339 | cmfmc2_inout(i,k,j) = cmfmc2(1,kflip) |
---|
340 | cmfmc(i,k,j) = cmfmc(i,k,j) + cmfmc2(1,kflip) |
---|
341 | |
---|
342 | end do !k-loop to kte |
---|
343 | |
---|
344 | do k = kts,kte+1 |
---|
345 | kflip = kte-k+2 |
---|
346 | |
---|
347 | !Convective fluxes of 'sl' and 'qt' in energy unit |
---|
348 | cmfsl(i,k,j) = slflx(1,kflip) |
---|
349 | cmflq(i,k,j) = qtflx(1,kflip)*latvap |
---|
350 | end do !k-loop to kte+1 |
---|
351 | |
---|
352 | !Calculate fractional occurance of shallow convection |
---|
353 | !~Not doing this since it would require adding time averaging ability across output times |
---|
354 | |
---|
355 | !Rain rate for shallow convection |
---|
356 | ! rainsh(i,j) = precc(1)*1e3 |
---|
357 | pratesh(i,j) = precc(1)*1e3/dt !~this will need changing for adaptive time steps and cudt |
---|
358 | |
---|
359 | !Get indices of convection top and bottom based on deep+shallow |
---|
360 | !Note: cnt2 and cnb2 have indices decreasing with height, but |
---|
361 | ! cutop and cubot have indicies increasing with height |
---|
362 | kflip = kte - cutop(i,j) + 1 |
---|
363 | cnt = kflip |
---|
364 | if( cnt2(1) < kflip ) cnt = cnt2(1) |
---|
365 | cutop(i,j) = kte - cnt + 1 |
---|
366 | |
---|
367 | kflip = kte - cubot(i,j) + 1 |
---|
368 | cnb = kflip |
---|
369 | if( cnb2(1) > kflip ) cnb = cnb2(1) |
---|
370 | cubot(i,j) = kte - cnb + 1 |
---|
371 | |
---|
372 | !Add shallow reserved cloud condensate to deep reserved cloud condensate |
---|
373 | !dlf done above, rliq (m/s) |
---|
374 | rliq2_inout(i,j) = rliq2(1) |
---|
375 | rliq(i,j) = rliq(i,j) + rliq2(1) |
---|
376 | |
---|
377 | end do |
---|
378 | end do ij_loops |
---|
379 | END SUBROUTINE camuwshcu_driver |
---|
380 | |
---|
381 | END MODULE module_cu_camuwshcu_driver |
---|