CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rrtmgp_cloud_mp.F90
1
6 use machine, only: kind_phys
7 use radiation_tools, only: check_error_msg
9 use rrtmgp_lw_cloud_optics, only: &
10 radliq_lwr => radliq_lwrlw, radliq_upr => radliq_uprlw,&
11 radice_lwr => radice_lwrlw, radice_upr => radice_uprlw
12 use module_mp_thompson, only: calc_effectrad, nt_c_l, nt_c_o, re_qc_min, re_qc_max, &
13 re_qi_min, re_qi_max, re_qs_min, re_qs_max
16
17 real (kind_phys), parameter :: &
18 cld_limit_lower = 0.001, &
19 cld_limit_ovcst = 1.0 - 1.0e-8, &
20 reliq_def = 10.0 , &
21 reice_def = 50.0, &
22 rerain_def = 1000.0, &
23 resnow_def = 250.0, &
24 reice_min = 10.0, &
25 reice_max = 150.0
26
27 public gfs_rrtmgp_cloud_mp_init, gfs_rrtmgp_cloud_mp_run, gfs_rrtmgp_cloud_mp_finalize
28
29contains
30
37 subroutine gfs_rrtmgp_cloud_mp_run(nCol, nLev, nTracers, ncnd, i_cldliq, i_cldice, &
38 i_cldrain, i_cldsnow, i_cldgrpl, i_cldtot, i_cldliq_nc, i_cldice_nc, i_twa, kdt, &
39 imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_samf, doSWrad, doLWrad, effr_in, lmfshal, &
40 ltaerosol,mraerosol, icloud, imp_physics, imp_physics_thompson, imp_physics_gfdl, &
41 lgfdlmprad, do_mynnedmf, uni_cld, lmfdeep2, p_lev, p_lay, t_lay, qs_lay, q_lay, &
42 relhum, lsmask, xlon, xlat, dx, tv_lay, effrin_cldliq, effrin_cldice, &
43 effrin_cldrain, effrin_cldsnow, tracer, cnv_mixratio, cld_cnv_frac, qci_conv, &
44 deltaZ, deltaZc, deltaP, qc_mynn, qi_mynn, cld_pbl_frac, con_g, con_rd, con_eps, &
45 con_ttp, doGP_cldoptics_PADE, doGP_cldoptics_LUT, doGP_smearclds, &
46 cld_frac, cld_lwp, cld_reliq, &
47 cld_iwp, cld_reice, cld_swp, cld_resnow, cld_rwp, cld_rerain, precip_frac, &
48 cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, cld_cnv_reice, cld_pbl_lwp, &
49 cld_pbl_reliq, cld_pbl_iwp, cld_pbl_reice, lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
50 cldfra2d, errmsg, errflg)
51 implicit none
52
53 ! Inputs
54 integer, intent(in) :: &
55 ncol, & !< Number of horizontal grid points
56 nlev, & !< Number of vertical layers
57 ncnd, & !< Number of cloud condensation types.
58 ntracers, & !< Number of tracers from model.
59 i_cldliq, & !< Index into tracer array for cloud liquid.
60 i_cldice, & !< Index into tracer array for cloud ice.
61 i_cldrain, & !< Index into tracer array for cloud rain.
62 i_cldsnow, & !< Index into tracer array for cloud snow.
63 i_cldgrpl, & !< Index into tracer array for cloud groupel.
64 i_cldtot, & !< Index into tracer array for cloud total amount.
65 i_cldliq_nc, & !< cloud liquid number concentration.
66 i_cldice_nc, & !< cloud ice number concentration.
67 i_twa, & !< water friendly aerosol.
68 imfdeepcnv, & !< Choice of mass-flux deep convection scheme
69 imfdeepcnv_gf, & !< Flag for Grell-Freitas deep convection scheme
70 imfdeepcnv_samf, & !< Flag for scale awware mass flux convection scheme
71 kdt, & !< Current forecast iteration
72 imp_physics, & !< Choice of microphysics scheme
73 imp_physics_thompson, & !< Choice of Thompson
74 imp_physics_gfdl, & !< Choice of GFDL
75 icloud
76 logical, intent(in) :: &
77 doswrad, & !< Call SW radiation?
78 dolwrad, & !< Call LW radiation?
79 effr_in, & !< Provide hydrometeor radii from macrophysics?
80 lmfshal, & !< Flag for mass-flux shallow convection scheme used by Xu-Randall
81 ltaerosol, & !< Flag for aerosol option
82 mraerosol, & !< Flag for aerosol option
83 lgfdlmprad, & !< Flag for GFDLMP radiation interaction
84 do_mynnedmf, & !< Flag to activate MYNN-EDMF
85 uni_cld, & !< Flag for unified cloud scheme
86 lmfdeep2, & !< Flag for mass flux deep convection
87 dogp_cldoptics_lut, & !< Flag to do GP cloud-optics (LUTs)
88 dogp_cldoptics_pade, & !< (PADE approximation)
89 dogp_smearclds
90 real(kind_phys), intent(in) :: &
91 con_g, & !< Physical constant: gravitational constant
92 con_rd, & !< Physical constant: gas-constant for dry air
93 con_ttp, & !< Triple point temperature of water (K)
94 con_eps
95 real(kind_phys), dimension(:), intent(in) :: &
96 lsmask, & !< Land/Sea mask
97 xlon, & !< Longitude
98 xlat, & !< Latitude
99 dx
100 real(kind_phys), dimension(:,:), intent(in), optional :: &
101 tv_lay, & !< Virtual temperature (K)
102 t_lay, & !< Temperature (K)
103 qs_lay, & !< Saturation vapor pressure (Pa)
104 q_lay, & !< water-vapor mixing ratio (kg/kg)
105 relhum, & !< Relative humidity
106 p_lay
107 real(kind_phys), dimension(:,:), intent(in) :: &
108 cnv_mixratio
109 real(kind_phys), dimension(:,:), intent(in), optional :: &
110 qci_conv, & !< Convective cloud condesate after rainout (kg/kg)
111 deltaz, & !< Layer-thickness (m)
112 deltazc, & !< Layer-thickness, from layer centers (m)
113 deltap, & !< Layer-thickness (Pa)
114 qc_mynn, & !<
115 qi_mynn
116 real(kind_phys), dimension(:,:), intent(in), optional :: &
117 cld_pbl_frac
118 real(kind_phys), dimension(:,:), intent(inout), optional :: &
119 effrin_cldliq, & !< Effective radius for stratiform liquid cloud-particles (microns)
120 effrin_cldice, & !< Effective radius for stratiform ice cloud-particles (microns)
121 effrin_cldsnow
122 real(kind_phys), dimension(:,:), intent(in), optional :: &
123 effrin_cldrain
124 real(kind_phys), dimension(:,:), intent(in), optional :: &
125 p_lev
126 real(kind_phys), dimension(:,:,:),intent(in) :: &
127 tracer
128
129 ! Outputs
130 real(kind_phys), dimension(:), intent(inout) :: &
131 lwp_ex, & !< Total liquid water path from explicit microphysics
132 iwp_ex, & !< Total ice water path from explicit microphysics
133 lwp_fc, & !< Total liquid water path from cloud fraction scheme
134 iwp_fc
135 real(kind_phys), dimension(:), intent(out) :: &
136 cldfra2d
137 real(kind_phys), dimension(:,:),intent(inout) :: &
138 cld_frac, & !< Cloud-fraction for stratiform clouds
139 cld_lwp, & !< Water path for stratiform liquid cloud-particles
140 cld_reliq, & !< Effective radius for stratiform liquid cloud-particles
141 cld_iwp, & !< Water path for stratiform ice cloud-particles
142 cld_reice, & !< Effective radius for stratiform ice cloud-particles
143 cld_swp, & !< Water path for snow hydrometeors
144 cld_resnow, & !< Effective radius for snow hydrometeors
145 cld_rwp, & !< Water path for rain hydrometeors
146 cld_rerain
147 real(kind_phys), dimension(:,:),intent(inout), optional :: &
148 precip_frac, & !< Precipitation fraction
149 cld_cnv_frac, & !< Cloud-fraction for convective clouds
150 cld_cnv_lwp, & !< Water path for convective liquid cloud-particles
151 cld_cnv_reliq, & !< Effective radius for convective liquid cloud-particles
152 cld_cnv_iwp, & !< Water path for convective ice cloud-particles
153 cld_cnv_reice, & !< Effective radius for convective ice cloud-particles
154 cld_pbl_lwp, & !< Water path for SGS PBL liquid cloud-particles
155 cld_pbl_reliq, & !< Effective radius for SGS PBL liquid cloud-particles
156 cld_pbl_iwp, & !< Water path for SGS PBL ice cloud-particles
157 cld_pbl_reice
158 character(len=*), intent(out) :: &
159 errmsg
160 integer, intent(out) :: &
161 errflg
162
163 ! Local
164 integer :: icol, ilay
165 real(kind_phys) :: alpha0
166 real(kind_phys), dimension(nCol,nLev) :: cldcov, cldtot, cldcnv
167
168 if (.not. (doswrad .or. dolwrad)) return
169
170 ! Initialize CCPP error handling variables
171 errmsg = ''
172 errflg = 0
173
174 ! ###################################################################################
175 ! GFDL Microphysics
176 ! ("Implicit" SGS cloud-coupling to the radiation)
177 ! ###################################################################################
178 if (imp_physics == imp_physics_gfdl) then
179 ! GFDL-Lin
180 if (.not. lgfdlmprad) then
181 errflg = 1
182 errmsg = "ERROR: MP choice not available with RRTMGP"
183 return
184 ! GFDL-EMC
185 else
186
187 ! "cld_frac" is modified prior to include subgrid scale cloudiness, see
188 ! module_SGSCloud_RadPre.F90.
189 do ilay = 1, nlev
190 do icol = 1, ncol
191 !
192 ! SGS clouds present, use cloud-fraction modified to include sgs clouds.
193 !
194 if ((imfdeepcnv==imfdeepcnv_gf .or. do_mynnedmf) .and. kdt>1) then
195 ! MYNN sub-grid cloud fraction.
196 if (do_mynnedmf) then
197 ! If rain/snow present, use GFDL MP cloud-fraction...
198 if (tracer(icol,ilay,i_cldrain)>1.0e-7 .OR. tracer(icol,ilay,i_cldsnow)>1.0e-7) then
199 cld_frac(icol,ilay) = tracer(icol,ilay,i_cldtot)
200 endif
201 ! GF sub-grid cloud fraction.
202 else
203 ! If no convective cloud condensate present, use GFDL MP cloud-fraction....
204 if (qci_conv(icol,ilay) <= 0.) then
205 cld_frac(icol,ilay) = tracer(icol,ilay,i_cldtot)
206 endif
207 endif
208 !
209 ! No SGS clouds, use GFDL MP cloud-fraction...
210 !
211 else
212 cld_frac(icol,ilay) = tracer(icol,ilay,i_cldtot)
213 endif
214 enddo
215 enddo
216
217 call cloud_mp_uni(ncol, nlev, ntracers, ncnd, i_cldliq, i_cldice, i_cldrain, &
218 i_cldsnow, i_cldgrpl, i_cldtot, effr_in, kdt, lsmask, p_lev, p_lay, &
219 t_lay, tv_lay, effrin_cldliq, effrin_cldice, effrin_cldsnow, tracer, &
220 con_g, con_rd, con_ttp, cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice,&
221 cld_swp, cld_resnow, cld_rwp, cld_rerain, effrin_cldrain=effrin_cldrain)
222 end if
223 endif
224
225 ! ###################################################################################
226 ! Thompson Microphysics
227 ! ("Explicit" SGS cloud-coupling to the radiation)
228 ! ###################################################################################
229 if (imp_physics == imp_physics_thompson) then
230
231 ! MYNN-EDMF PBL clouds?
232 if(do_mynnedmf) then
233 call cloud_mp_mynn(ncol, nlev, lsmask, t_lay, p_lev, p_lay, qs_lay, relhum, &
234 qc_mynn, qi_mynn, con_ttp, con_g, &
235 cld_pbl_lwp, cld_pbl_reliq, cld_pbl_iwp, cld_pbl_reice, cld_pbl_frac)
236 endif
237
238 ! Grell-Freitas convective clouds?
239 if (imfdeepcnv == imfdeepcnv_gf) then
240 alpha0 = 100.
241 call cloud_mp_gf(ncol, nlev, lsmask, t_lay, p_lev, p_lay, qs_lay, relhum, &
242 qci_conv, con_ttp, con_g, alpha0, &
243 cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, cld_cnv_reice, cld_cnv_frac)
244 endif
245
246 ! SAMF scale & aerosol-aware mass-flux convective clouds?
247 if (imfdeepcnv == imfdeepcnv_samf) then
248 alpha0 = 200.
249 call cloud_mp_samf(ncol, nlev, t_lay, p_lev, p_lay, qs_lay, relhum, &
250 cnv_mixratio, con_ttp, con_g, alpha0, &
251 cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, cld_cnv_reice, cld_cnv_frac)
252 endif
253
254 ! Update particle size using modified mixing-ratios from Thompson.
255 call cmp_reff_thompson(nlev, ncol, i_cldliq, i_cldice, i_cldsnow, i_cldice_nc, &
256 i_cldliq_nc, i_twa, q_lay, p_lay, t_lay, tracer, con_eps, con_rd, ltaerosol,&
257 mraerosol, lsmask, effrin_cldliq, effrin_cldice, effrin_cldsnow)
258 cld_reliq = effrin_cldliq
259 cld_reice = effrin_cldice
260 cld_resnow = effrin_cldsnow
261
262 ! Thomson MP using modified Xu-Randall cloud-fraction (additionally conditioned on RH)
263 alpha0 = 2000.
264 if (lmfshal) then
265 alpha0 = 100.
266 if (lmfdeep2) alpha0 = 200.
267 endif
268 call cloud_mp_thompson(ncol, nlev, ntracers, ncnd, i_cldliq, i_cldice, i_cldrain,&
269 i_cldsnow, i_cldgrpl, p_lev, p_lay, tv_lay, t_lay, tracer, qs_lay, q_lay, &
270 relhum, con_ttp, con_g, con_rd, con_eps, alpha0, cnv_mixratio, lwp_ex, &
271 iwp_ex, lwp_fc, iwp_fc, cld_frac, cld_lwp, cld_iwp, cld_swp, cld_rwp, &
272 cond_cfrac_onrh = .true., dogp_smearclds = dogp_smearclds)
273 endif
274
275 ! Bound effective radii for RRTMGP, LUT's for cloud-optics go from
276 ! 2.5 - 21.5 microns for liquid clouds,
277 ! 10 - 180 microns for ice-clouds
278 if (dogp_cldoptics_pade .or. dogp_cldoptics_lut) then
279 where(cld_reliq .lt. radliq_lwr) cld_reliq = radliq_lwr
280 where(cld_reliq .gt. radliq_upr) cld_reliq = radliq_upr
281 where(cld_reice .lt. radice_lwr) cld_reice = radice_lwr
282 where(cld_reice .gt. radice_upr) cld_reice = radice_upr
283 if (imfdeepcnv == imfdeepcnv_samf .or. imfdeepcnv == imfdeepcnv_gf) then
284 where(cld_cnv_reliq .lt. radliq_lwr) cld_cnv_reliq = radliq_lwr
285 where(cld_cnv_reliq .gt. radliq_upr) cld_cnv_reliq = radliq_upr
286 where(cld_cnv_reice .lt. radice_lwr) cld_cnv_reice = radice_lwr
287 where(cld_cnv_reice .gt. radice_upr) cld_cnv_reice = radice_upr
288 endif
289 if (do_mynnedmf) then
290 where(cld_pbl_reliq .lt. radliq_lwr) cld_pbl_reliq = radliq_lwr
291 where(cld_pbl_reliq .gt. radliq_upr) cld_pbl_reliq = radliq_upr
292 where(cld_pbl_reice .lt. radice_lwr) cld_pbl_reice = radice_lwr
293 where(cld_pbl_reice .gt. radice_upr) cld_pbl_reice = radice_upr
294 endif
295 endif
296
297 ! Instantaneous 2D (max-in-column) cloud fraction
298 do icol = 1, ncol
299 cldfra2d(icol) = 0._kind_phys
300 do ilay = 1, nlev-1
301 cldfra2d(icol) = max(cldfra2d(icol), cld_frac(icol,ilay))
302 enddo
303 enddo
304
305 precip_frac(1:ncol,1:nlev) = cld_frac(1:ncol,1:nlev)
306
307 end subroutine gfs_rrtmgp_cloud_mp_run
308
323 subroutine cloud_mp_gf(nCol, nLev, lsmask, t_lay, p_lev, p_lay, qs_lay, relhum, &
324 qci_conv, con_ttp, con_g, alpha0, cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, &
325 cld_cnv_reice, cld_cnv_frac)
326 implicit none
327
328 ! Inputs
329 integer, intent(in) :: &
330 nCol, & !< Number of horizontal grid points
331 nLev
332 real(kind_phys), dimension(:), intent(in) :: &
333 lsmask
334 real(kind_phys), intent(in) :: &
335 con_g, & !< Physical constant: gravitational constant
336 con_ttp, & !< Triple point temperature of water (K)
337 alpha0
338 real(kind_phys), dimension(:,:),intent(in) :: &
339 t_lay, & !< Temperature at layer centers (K)
340 p_lev, & !< Pressure at layer interfaces (Pa)
341 p_lay, & !<
342 qs_lay, & !<
343 relhum, & !<
344 qci_conv
345 ! Outputs
346 real(kind_phys), dimension(:,:),intent(inout) :: &
347 cld_cnv_lwp, & !< Convective cloud liquid water path
348 cld_cnv_reliq, & !< Convective cloud liquid effective radius
349 cld_cnv_iwp, & !< Convective cloud ice water path
350 cld_cnv_reice, & !< Convective cloud ice effecive radius
351 cld_cnv_frac
352 ! Local
353 integer :: iCol, iLay
354 real(kind_phys) :: tem1, deltaP, clwc, qc, qi
355
356 tem1 = 1.0e5/con_g
357 do ilay = 1, nlev
358 do icol = 1, ncol
359 if (qci_conv(icol,ilay) > 0.) then
360 ! Partition the convective clouds by phase.
361 qc = qci_conv(icol,ilay)*( min(1., max(0., (t_lay(icol,ilay)-244.)*0.04)))
362 qi = qci_conv(icol,ilay)*(1. - min(1., max(0., (t_lay(icol,ilay)-244.)*0.04)))
363
364 ! Compute LWP/IWP
365 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))*0.01
366 cld_cnv_lwp(icol,ilay) = max(0., qc * tem1 * deltap)
367 cld_cnv_iwp(icol,ilay) = max(0., qi * tem1 * deltap)
368
369 ! Particle sizes
370 if (nint(lsmask(icol)) == 1) then !land
371 if(qc > 1.e-8) cld_cnv_reliq(icol,ilay) = 5.4
372 else
373 !eff radius cloud water (microns), from Miles et al.
374 if(qc > 1.e-8) cld_cnv_reliq(icol,ilay) = 9.6
375 endif
376 !eff radius cloud ice (microns), from Mishra et al. (2014, JGR Atmos, fig 6b)
377 if(qi > 1.e-8) cld_cnv_reice(icol,ilay) = max(173.45 + 2.14*(t_lay(icol,ilay)-273.15), 20.)
378
379 ! Xu-Randall (1996) cloud-fraction.
380 cld_cnv_frac(icol,ilay) = cld_frac_xurandall(p_lay(icol,ilay), &
381 qs_lay(icol,ilay), relhum(icol,ilay), qc+qi, alpha0)
382 endif
383 enddo
384 enddo
385 end subroutine cloud_mp_gf
386
396 subroutine cloud_mp_mynn(nCol, nLev, lsmask, t_lay, p_lev, p_lay, qs_lay, relhum, &
397 qc_mynn, qi_mynn, con_ttp, con_g, cld_pbl_lwp, cld_pbl_reliq, cld_pbl_iwp, &
398 cld_pbl_reice, cld_pbl_frac)
399 implicit none
400
401 ! Inputs
402 integer, intent(in) :: &
403 nCol, & !< Number of horizontal grid points
404 nLev
405 real(kind_phys), dimension(:), intent(in) :: &
406 lsmask
407 real(kind_phys), intent(in) :: &
408 con_g, & !< Physical constant: gravitational constant
409 con_ttp
410 real(kind_phys), dimension(:,:),intent(in) :: &
411 t_lay, & !< Temperature at layer centers (K)
412 p_lev, & !< Pressure at layer interfaces (Pa)
413 p_lay, & !<
414 qs_lay, & !<
415 relhum, & !<
416 qc_mynn, & !< Liquid cloud mixing-ratio (MYNN PBL cloud)
417 qi_mynn, & !< Ice cloud mixing-ratio (MYNN PBL cloud)
418 cld_pbl_frac
419 ! Outputs
420 real(kind_phys), dimension(:,:),intent(inout) :: &
421 cld_pbl_lwp, & !< Convective cloud liquid water path
422 cld_pbl_reliq, & !< Convective cloud liquid effective radius
423 cld_pbl_iwp, & !< Convective cloud ice water path
424 cld_pbl_reice
425
426 ! Local
427 integer :: iCol, iLay
428 real(kind_phys) :: tem1, qc, qi, deltaP
429
430 tem1 = 1.0e5/con_g
431 do ilay = 1, nlev
432 do icol = 1, ncol
433 if (cld_pbl_frac(icol,ilay) > cld_limit_lower) then
434 ! Cloud mixing-ratios (DJS asks: Why is this done?)
435 qc = qc_mynn(icol,ilay)*cld_pbl_frac(icol,ilay)
436 qi = qi_mynn(icol,ilay)*cld_pbl_frac(icol,ilay)
437
438 ! LWP/IWP
439 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))
440 cld_pbl_lwp(icol,ilay) = max(0., qc * tem1 * deltap)
441 cld_pbl_iwp(icol,ilay) = max(0., qi * tem1 * deltap)
442
443 ! Particle sizes
444 if (nint(lsmask(icol)) == 1) then
445 if(qc > 1.e-8) cld_pbl_reliq(icol,ilay) = 5.4
446 else
447 ! Cloud water (microns), from Miles et al.
448 if(qc > 1.e-8) cld_pbl_reliq(icol,ilay) = 9.6
449 endif
450 ! Cloud ice (microns), from Mishra et al. (2014, JGR Atmos, fig 6b)
451 if(qi > 1.e-8) cld_pbl_reice(icol,ilay) = max(173.45 + 2.14*(t_lay(icol,ilay)-273.15), 20.)
452 endif
453 enddo
454 enddo
455 end subroutine cloud_mp_mynn
456
457
468 subroutine cloud_mp_samf(nCol, nLev, t_lay, p_lev, p_lay, qs_lay, relhum, &
469 cnv_mixratio, con_ttp, con_g, alpha0, cld_cnv_lwp, cld_cnv_reliq, cld_cnv_iwp, &
470 cld_cnv_reice, cld_cnv_frac)
471 implicit none
472
473 ! Inputs
474 integer, intent(in) :: &
475 nCol, & !< Number of horizontal grid points
476 nLev
477 real(kind_phys), intent(in) :: &
478 con_g, & !< Physical constant: gravity (m s-2)
479 con_ttp, & !< Triple point temperature of water (K)
480 alpha0
481 real(kind_phys), dimension(:,:),intent(in) :: &
482 t_lay, & !< Temperature at layer-centers (K)
483 p_lev, & !< Pressure at layer-interfaces (Pa)
484 p_lay, & !< Presure at layer-centers (Pa)
485 qs_lay, & !< Specific-humidity at layer-centers (kg/kg)
486 relhum, & !< Relative-humidity (1)
487 cnv_mixratio
488 ! Outputs
489 real(kind_phys), dimension(:,:),intent(inout) :: &
490 cld_cnv_lwp, & !< Convective cloud liquid water path
491 cld_cnv_reliq, & !< Convective cloud liquid effective radius
492 cld_cnv_iwp, & !< Convective cloud ice water path
493 cld_cnv_reice, & !< Convective cloud ice effecive radius
494 cld_cnv_frac
495 ! Local
496 integer :: iCol, iLay
497 real(kind_phys) :: tem0, tem1, deltaP, clwc
498
499 tem0 = 1.0e5/con_g
500 do ilay = 1, nlev
501 do icol = 1, ncol
502 if (cnv_mixratio(icol,ilay) > 0._kind_phys) then
503 tem1 = min(1.0, max(0.0, (con_ttp-t_lay(icol,ilay))*0.05))
504 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))*0.01
505 clwc = max(0.0, cnv_mixratio(icol,ilay)) * tem0 * deltap
506 cld_cnv_iwp(icol,ilay) = clwc * tem1
507 cld_cnv_lwp(icol,ilay) = clwc - cld_cnv_iwp(icol,ilay)
508 cld_cnv_reliq(icol,ilay) = reliq_def
509 cld_cnv_reice(icol,ilay) = reice_def
510
511 ! Xu-Randall (1996) cloud-fraction.
512 cld_cnv_frac(icol,ilay) = cld_frac_xurandall(p_lay(icol,ilay), &
513 qs_lay(icol,ilay), relhum(icol,ilay), cnv_mixratio(icol,ilay), alpha0)
514 endif
515 enddo
516 enddo
517
518 end subroutine cloud_mp_samf
519
524 subroutine cloud_mp_uni(nCol, nLev, nTracers, ncnd, i_cldliq, i_cldice, i_cldrain, &
525 i_cldsnow, i_cldgrpl, i_cldtot, effr_in, kdt, lsmask, p_lev, p_lay, t_lay, tv_lay,&
526 effrin_cldliq, effrin_cldice, effrin_cldsnow, tracer, con_g, con_rd, con_ttp, &
527 cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, cld_swp, cld_resnow, cld_rwp, &
528 cld_rerain, effrin_cldrain)
529 implicit none
530
531 ! Inputs
532 integer, intent(in) :: &
533 nCol, & !< Number of horizontal grid points
534 nLev, & !< Number of vertical layers
535 ncnd, & !< Number of cloud condensation types.
536 nTracers, & !< Number of tracers from model.
537 i_cldliq, & !< Index into tracer array for cloud liquid.
538 i_cldice, & !< Index into tracer array for cloud ice.
539 i_cldrain, & !< Index into tracer array for cloud rain.
540 i_cldsnow, & !< Index into tracer array for cloud snow.
541 i_cldgrpl, & !< Index into tracer array for cloud groupel.
542 i_cldtot, & !< Index into tracer array for cloud total amount.
543 kdt
544 logical, intent(in) :: &
545 effr_in
546 real(kind_phys), intent(in) :: &
547 con_g, & !< Physical constant: gravitational constant
548 con_ttp, & !< Triple point temperature of water (K)
549 con_rd
550 real(kind_phys), dimension(:), intent(in) :: &
551 lsmask
552 real(kind_phys), dimension(:,:), intent(in) :: &
553 t_lay, & !< Temperature at model-layers (K)
554 tv_lay, & !< Virtual temperature (K)
555 p_lay, & !< Pressure at model-layers (Pa)
556 cld_frac, & !< Total cloud fraction
557 effrin_cldliq, & !< Effective radius for liquid cloud-particles (microns)
558 effrin_cldice, & !< Effective radius for ice cloud-particles (microns)
559 effrin_cldsnow
560 real(kind_phys), dimension(:,:), intent(in), optional :: &
561 effrin_cldrain
562 real(kind_phys), dimension(:,:), intent(in) :: &
563 p_lev
564 real(kind_phys), dimension(:,:,:),intent(in) :: &
565 tracer
566
567 ! Outputs
568 real(kind_phys), dimension(:,:),intent(inout) :: &
569 cld_lwp, & !< Cloud liquid water path
570 cld_reliq, & !< Cloud liquid effective radius
571 cld_iwp, & !< Cloud ice water path
572 cld_reice, & !< Cloud ice effecive radius
573 cld_swp, & !< Cloud snow water path
574 cld_resnow, & !< Cloud snow effective radius
575 cld_rwp, & !< Cloud rain water path
576 cld_rerain
577
578 ! Local variables
579 real(kind_phys) :: tem1,tem2,tem3,pfac,deltaP
580 real(kind_phys), dimension(nCol, nLev, min(4,ncnd)) :: cld_condensate
581 integer :: iCol,iLay,l,ncndl
582
583 ! Cloud condensate
584 cld_condensate(1:ncol,1:nlev,1) = tracer(1:ncol,1:nlev,i_cldliq) ! -liquid water
585 cld_condensate(1:ncol,1:nlev,2) = tracer(1:ncol,1:nlev,i_cldice) ! -ice water
586 if (ncnd > 2) then
587 cld_condensate(1:ncol,1:nlev,3) = tracer(1:ncol,1:nlev,i_cldrain) ! -rain water
588 cld_condensate(1:ncol,1:nlev,4) = tracer(1:ncol,1:nlev,i_cldsnow) + &! -snow + grapuel
589 tracer(1:ncol,1:nlev,i_cldgrpl)
590 endif
591
592 ! Cloud water path (g/m2)
593 tem1 = 1.0e5/con_g
594 do ilay = 1, nlev
595 do icol = 1, ncol
596 ! Compute liquid/ice condensate path from mixing ratios (kg/kg)->(g/m2)
597 if (cld_frac(icol,ilay) > cld_limit_lower) then
598 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))*0.01
599 cld_lwp(icol,ilay) = max(0., cld_condensate(icol,ilay,1) * tem1 * deltap)
600 cld_iwp(icol,ilay) = max(0., cld_condensate(icol,ilay,2) * tem1 * deltap)
601 if (ncnd > 2) then
602 cld_rwp(icol,ilay) = max(0., cld_condensate(icol,ilay,3) * tem1 * deltap)
603 cld_swp(icol,ilay) = max(0., cld_condensate(icol,ilay,4) * tem1 * deltap)
604 endif
605 endif
606 enddo
607 enddo
608
609 ! Particle size
610 do ilay = 1, nlev
611 do icol = 1, ncol
612 ! Use radii provided from the macrophysics
613 if (effr_in) then
614 cld_reliq(icol,ilay) = effrin_cldliq(icol,ilay)
615 cld_reice(icol,ilay) = max(reice_min, min(reice_max,effrin_cldice(icol,ilay)))
616 cld_resnow(icol,ilay) = effrin_cldsnow(icol,ilay)
617 if (present(effrin_cldrain)) then
618 cld_rerain(icol,ilay) = effrin_cldrain(icol,ilay)
619 else
620 cld_rerain(icol,ilay) = rerain_def
621 endif
622 else
623 ! Compute effective liquid cloud droplet radius over land.
624 if (nint(lsmask(icol)) == 1) then
625 cld_reliq(icol,ilay) = 5.0 + 5.0 * min(1.0, max(0.0, (con_ttp-t_lay(icol,ilay))*0.05))
626 endif
627 ! Compute effective ice cloud droplet radius following Heymsfield
628 ! and McFarquhar (1996) \cite heymsfield_and_mcfarquhar_1996.
629 tem2 = t_lay(icol,ilay) - con_ttp
630 if (cld_iwp(icol,ilay) > 0.0) then
631 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))*0.01
632 tem3 = (con_g/con_rd ) * cld_iwp(icol,ilay) * (0.01*p_lay(icol,ilay)) / (deltap*tv_lay(icol,ilay))
633 if (tem2 < -50.0) then
634 cld_reice(icol,ilay) = (1250.0/9.917) * tem3 ** 0.109
635 elseif (tem2 < -40.0) then
636 cld_reice(icol,ilay) = (1250.0/9.337) * tem3 ** 0.08
637 elseif (tem2 < -30.0) then
638 cld_reice(icol,ilay) = (1250.0/9.208) * tem3 ** 0.055
639 else
640 cld_reice(icol,ilay) = (1250.0/9.387) * tem3 ** 0.031
641 endif
642 cld_reice(icol,ilay) = max(10.0, min(cld_reice(icol,ilay), 150.0))
643 endif
644 endif ! effr_in
645 enddo ! nCol
646 enddo ! nLev
647
648 end subroutine cloud_mp_uni
649
661 subroutine cloud_mp_thompson(nCol, nLev, nTracers, ncnd, i_cldliq, i_cldice, i_cldrain,&
662 i_cldsnow, i_cldgrpl, p_lev, p_lay, tv_lay, t_lay, tracer, qs_lay, q_lay, relhum, &
663 con_ttp, con_g, con_rd, con_eps, alpha0, cnv_mixratio, lwp_ex, iwp_ex, lwp_fc, &
664 iwp_fc, cld_frac, cld_lwp, cld_iwp, cld_swp, cld_rwp, cond_cfrac_onRH, doGP_smearclds)
665 implicit none
666
667 ! Inputs
668 logical, intent(in), optional :: &
669 cond_cfrac_onRH, & !< If true, cloud-fracion set to unity when rh>99%
670 doGP_smearclds
671 integer, intent(in) :: &
672 nCol, & !< Number of horizontal grid points
673 nLev, & !< Number of vertical layers
674 ncnd, & !< Number of cloud condensation types.
675 nTracers, & !< Number of tracers from model.
676 i_cldliq, & !< Index into tracer array for cloud liquid amount.
677 i_cldice, & !< cloud ice amount.
678 i_cldrain, & !< cloud rain amount.
679 i_cldsnow, & !< cloud snow amount.
680 i_cldgrpl
681 real(kind_phys), intent(in) :: &
682 con_ttp, & !< Triple point temperature of water (K)
683 con_g, & !< Physical constant: gravitational constant
684 con_rd, & !< Physical constant: gas-constant for dry air
685 con_eps, & !< Physical constant: gas constant air / gas constant H2O
686 alpha0
687 real(kind_phys), dimension(:,:), intent(in) :: &
688 tv_lay, & !< Virtual temperature (K)
689 t_lay, & !< Temperature (K)
690 qs_lay, & !< Saturation vapor pressure (Pa)
691 q_lay, & !< water-vapor mixing ratio (kg/kg)
692 relhum, & !< Relative humidity
693 p_lay, & !< Pressure at model-layers (Pa)
694 cnv_mixratio
695 real(kind_phys), dimension(:,:), intent(in) :: &
696 p_lev
697 real(kind_phys), dimension(:,:,:),intent(in) :: &
698 tracer
699
700 ! In/Outs
701 real(kind_phys), dimension(:), intent(inout) :: &
702 lwp_ex, & !< total liquid water path from explicit microphysics
703 iwp_ex, & !< total ice water path from explicit microphysics
704 lwp_fc, & !< total liquid water path from cloud fraction scheme
705 iwp_fc
706 real(kind_phys), dimension(:,:), intent(inout) :: &
707 cld_frac, & !< Total cloud fraction
708 cld_lwp, & !< Cloud liquid water path
709 cld_iwp, & !< Cloud ice water path
710 cld_swp, & !< Cloud snow water path
711 cld_rwp
712
713 ! Local variables
714 real(kind_phys) :: tem1, pfac, cld_mr, deltaP, tem2
715 real(kind_phys), dimension(nCol, nLev, min(4,ncnd)) :: cld_condensate
716 integer :: iCol,iLay,l
717
718 ! Cloud condensate
719 cld_condensate(1:ncol,1:nlev,1) = tracer(1:ncol,1:nlev,i_cldliq) ! -liquid water
720 cld_condensate(1:ncol,1:nlev,2) = tracer(1:ncol,1:nlev,i_cldice) ! -ice water
721 cld_condensate(1:ncol,1:nlev,3) = tracer(1:ncol,1:nlev,i_cldrain) ! -rain hydrometeors
722 cld_condensate(1:ncol,1:nlev,4) = tracer(1:ncol,1:nlev,i_cldsnow) ! -snow hydrometeors
723
724 cld_lwp(:,:) = 0.0
725 cld_iwp(:,:) = 0.0
726 cld_rwp(:,:) = 0.0
727 cld_swp(:,:) = 0.0
728 cld_frac(:,:) = 0.0
729 tem1 = 1.0e5/con_g
730 do ilay = 1, nlev-1
731 do icol = 1, ncol
732 ! Add convective cloud to gridmean cloud?
733 if (dogp_smearclds) then
734 tem2 = min(1.0, max(0.0, (con_ttp-t_lay(icol,ilay))*0.05))
735 cld_condensate(icol,ilay,1) = cld_condensate(icol,ilay,1) + cnv_mixratio(icol,ilay)*(1._kind_phys - tem2)
736 cld_condensate(icol,ilay,2) = cld_condensate(icol,ilay,2) + cnv_mixratio(icol,ilay)*tem2
737 endif
738 ! Compute liquid/ice condensate path from mixing ratios (kg/kg)->(g/m2)
739 deltap = abs(p_lev(icol,ilay+1)-p_lev(icol,ilay))*0.01
740 cld_lwp(icol,ilay) = max(0., cld_condensate(icol,ilay,1) * tem1 * deltap)
741 cld_iwp(icol,ilay) = max(0., cld_condensate(icol,ilay,2) * tem1 * deltap)
742 cld_rwp(icol,ilay) = max(0., cld_condensate(icol,ilay,3) * tem1 * deltap)
743 cld_swp(icol,ilay) = max(0., cld_condensate(icol,ilay,4) * tem1 * deltap)
744
745 ! Xu-Randall (1996) cloud-fraction. **Additionally, Conditioned on relative-humidity**
746 if (present(cond_cfrac_onrh) .and. relhum(icol,ilay) > 0.99) then
747 cld_frac(icol,ilay) = 1._kind_phys
748 else
749 cld_mr = cld_condensate(icol,ilay,1) + cld_condensate(icol,ilay,2) + &
750 cld_condensate(icol,ilay,3) + cld_condensate(icol,ilay,4)
751 cld_frac(icol,ilay) = cld_frac_xurandall(p_lay(icol,ilay), &
752 qs_lay(icol,ilay), relhum(icol,ilay), cld_mr, alpha0)
753 endif
754 enddo
755 enddo
756
757 ! Sum the liquid water and ice paths that come from explicit micro
758 ! What portion of water and ice contents is associated with the partly cloudy boxes?
759 do icol = 1, ncol
760 lwp_ex(icol) = 0.0
761 iwp_ex(icol) = 0.0
762 lwp_fc(icol) = 0.0
763 iwp_fc(icol) = 0.0
764 do ilay = 1, nlev-1
765 lwp_ex(icol) = lwp_ex(icol) + cld_lwp(icol,ilay)
766 iwp_ex(icol) = iwp_ex(icol) + cld_iwp(icol,ilay) + cld_swp(icol,ilay)
767 if (cld_frac(icol,ilay) .ge. cld_limit_lower .and. &
768 cld_frac(icol,ilay) .lt. cld_limit_ovcst) then
769 lwp_fc(icol) = lwp_fc(icol) + cld_lwp(icol,ilay)
770 iwp_fc(icol) = iwp_fc(icol) + cld_iwp(icol,ilay) + cld_swp(icol,ilay)
771 endif
772 enddo
773 lwp_fc(icol) = lwp_fc(icol)*1.e-3
774 iwp_fc(icol) = iwp_fc(icol)*1.e-3
775 lwp_ex(icol) = lwp_ex(icol)*1.e-3
776 iwp_ex(icol) = iwp_ex(icol)*1.e-3
777 enddo
778
779 end subroutine cloud_mp_thompson
780
784 function cld_frac_xurandall(p_lay, qs_lay, relhum, cld_mr, alpha)
785 implicit none
786 ! Inputs
787 real(kind_phys), intent(in) :: &
788 p_lay, & !< Pressure (Pa)
789 qs_lay, & !< Saturation vapor-pressure (Pa)
790 relhum, & !< Relative humidity
791 cld_mr, & !< Total cloud mixing ratio
792 alpha
793
794 ! Outputs
795 real(kind_phys) :: cld_frac_xurandall
796
797 ! Locals
798 real(kind_phys) :: clwt, clwm, onemrh, tem1, tem2, tem3
799
800 ! Parameters
801 real(kind_phys) :: &
802 lambda = 0.50, & !
803 p = 0.25
804
805 clwt = 1.0e-6 * (p_lay*0.001)
806 if (cld_mr > clwt) then
807 onemrh = max(1.e-10, 1.0 - relhum)
808 tem1 = alpha / min(max((onemrh*qs_lay)**lambda,0.0001),1.0)
809 tem2 = max(min(tem1*(cld_mr - clwt), 50.0 ), 0.0 )
810 tem3 = sqrt(sqrt(relhum)) ! This assumes "p" = 0.25. Identical, but cheaper than relhum**p
811 !
812 cld_frac_xurandall = max( tem3*(1.0-exp(-tem2)), 0.0 )
813 else
814 cld_frac_xurandall = 0.0
815 endif
816
817 return
818 end function
819
822 subroutine cmp_reff_thompson(nLev, nCol, i_cldliq, i_cldice, i_cldsnow, i_cldice_nc, &
823 i_cldliq_nc, i_twa, q_lay, p_lay, t_lay, tracer, con_eps, con_rd, ltaerosol, &
824 mraerosol, lsmask, effrin_cldliq, effrin_cldice, effrin_cldsnow)
825 implicit none
826
827 ! Inputs
828 integer, intent(in) :: nLev, nCol, i_cldliq, i_cldice, i_cldsnow, i_cldice_nc, &
829 i_cldliq_nc, i_twa
830 logical, intent(in) :: ltaerosol, mraerosol
831 real(kind_phys), intent(in) :: con_eps,con_rd
832 real(kind_phys), dimension(:,:),intent(in) :: q_lay, p_lay, t_lay
833 real(kind_phys), dimension(:,:,:),intent(in) :: tracer
834 real(kind_phys), dimension(:), intent(in) :: lsmask
835
836 ! Outputs
837 real(kind_phys), dimension(:,:), intent(inout) :: effrin_cldliq, effrin_cldice, &
838 effrin_cldsnow
839
840 ! Local
841 integer :: iCol, iLay
842 real(kind_phys) :: rho, orho
843 real(kind_phys),dimension(nCol,nLev) :: qv_mp, qc_mp, qi_mp, qs_mp, ni_mp, nc_mp, &
844 nwfa, re_cloud, re_ice, re_snow
845 integer :: ilsmask
846
847 ! Prepare cloud mixing-ratios and number concentrations for calc_effectRa
848 do ilay = 1, nlev
849 do icol = 1, ncol
850 qv_mp(icol,ilay) = q_lay(icol,ilay)/(1.-q_lay(icol,ilay))
851 rho = con_eps*p_lay(icol,ilay)/(con_rd*t_lay(icol,ilay)*(qv_mp(icol,ilay)+con_eps))
852 orho = 1./rho
853 qc_mp(icol,ilay) = tracer(icol,ilay,i_cldliq) / (1.-q_lay(icol,ilay))
854 qi_mp(icol,ilay) = tracer(icol,ilay,i_cldice) / (1.-q_lay(icol,ilay))
855 qs_mp(icol,ilay) = tracer(icol,ilay,i_cldsnow) / (1.-q_lay(icol,ilay))
856 ni_mp(icol,ilay) = tracer(icol,ilay,i_cldice_nc) / (1.-q_lay(icol,ilay))
857 if (ltaerosol .or. mraerosol) then
858 nc_mp(icol,ilay) = tracer(icol,ilay,i_cldliq_nc) / (1.-q_lay(icol,ilay))
859 nwfa(icol,ilay) = tracer(icol,ilay,i_twa)
860 if (qc_mp(icol,ilay) > 1.e-12 .and. nc_mp(icol,ilay) < 100.) then
861 nc_mp(icol,ilay) = make_dropletnumber(qc_mp(icol,ilay)*rho, nwfa(icol,ilay)*rho) * orho
862 endif
863 else
864 if (nint(lsmask(icol)) == 1) then !land
865 nc_mp(icol,ilay) = nt_c_l*orho
866 else
867 nc_mp(icol,ilay) = nt_c_o*orho
868 endif
869 endif
870 if (qi_mp(icol,ilay) > 1.e-12 .and. ni_mp(icol,ilay) < 100.) then
871 ni_mp(icol,ilay) = make_icenumber(qi_mp(icol,ilay)*rho, t_lay(icol,ilay)) * orho
872 endif
873 enddo
874 enddo
875
876 ! Compute effective radii for liquid/ice/snow.
877 do icol=1,ncol
878 ilsmask = nint(lsmask(icol))
879 call calc_effectrad (t_lay(icol,:), p_lay(icol,:), qv_mp(icol,:), qc_mp(icol,:), &
880 nc_mp(icol,:), qi_mp(icol,:), ni_mp(icol,:), qs_mp(icol,:), &
881 re_cloud(icol,:), re_ice(icol,:), re_snow(icol,:), ilsmask, &
882 1, nlev )
883 do ilay = 1, nlev
884 re_cloud(icol,ilay) = max(re_qc_min, min(re_cloud(icol,ilay), re_qc_max))
885 re_ice(icol,ilay) = max(re_qi_min, min(re_ice(icol,ilay), re_qi_max))
886 re_snow(icol,ilay) = max(re_qs_min, min(re_snow(icol,ilay), re_qs_max))
887 enddo
888 enddo
889
890 ! Scale to microns.
891 do ilay = 1, nlev
892 do icol = 1, ncol
893 effrin_cldliq(icol,ilay) = re_cloud(icol,ilay)*1.e6
894 effrin_cldice(icol,ilay) = re_ice(icol,ilay)*1.e6
895 effrin_cldsnow(icol,ilay) = re_snow(icol,ilay)*1.e6
896 enddo
897 enddo
898
899 end subroutine cmp_reff_thompson
900
901end module gfs_rrtmgp_cloud_mp
subroutine calc_effectrad(t1d, p1d, qv1d, qc1d, nc1d, qi1d, ni1d, qs1d, re_qc1d, re_qi1d, re_qs1d, lsml, kts, kte)
Compute radiation effective radii of cloud water, ice, and snow. These are entirely consistent with m...
elemental real function, public make_icenumber(q_ice, temp)
Table of lookup values of radiative effective radius of ice crystals as a function of Temperature fro...
subroutine, public progcld_thompson(plyr, plvl, tlyr, qlyr, qstl, rhly, clw, xlat, xlon, slmsk, dz, delp, ntrac, ntcw, ntiw, ntrw, ntsw, ntgl, ix, nlay, nlp1, uni_cld, lmfshal, lmfdeep2, cldcov, re_cloud, re_ice, re_snow, lwp_ex, iwp_ex, lwp_fc, iwp_fc, dzlay, gridkm, top_at_1, cldtot, cldcnv, cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, cld_rwp, cld_rerain, cld_swp, cld_resnow)
This subroutine added by G. Thompson specifically to account for explicit (microphysics-produced) clo...
real(kind=kind_phys), parameter reliq_def
default liq radius to 10 micron
real(kind=kind_phys), parameter reice_def
default ice radius to 50 micron
This module computes cloud related quantities for radiation computations.