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)
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
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)
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)
95 real(kind_phys),
dimension(:),
intent(in) :: &
96 lsmask, & !< Land/Sea mask
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
107 real(kind_phys),
dimension(:,:),
intent(in) :: &
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)
116 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
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)
122 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
124 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
126 real(kind_phys),
dimension(:,:,:),
intent(in) :: &
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
135 real(kind_phys),
dimension(:),
intent(out) :: &
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
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
158 character(len=*),
intent(out) :: &
160 integer,
intent(out) :: &
164 integer :: icol, ilay
165 real(kind_phys) :: alpha0
166 real(kind_phys),
dimension(nCol,nLev) :: cldcov, cldtot, cldcnv
168 if (.not. (doswrad .or. dolwrad))
return
178 if (imp_physics == imp_physics_gfdl)
then
180 if (.not. lgfdlmprad)
then
182 errmsg =
"ERROR: MP choice not available with RRTMGP"
194 if ((imfdeepcnv==imfdeepcnv_gf .or. do_mynnedmf) .and. kdt>1)
then
196 if (do_mynnedmf)
then
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)
204 if (qci_conv(icol,ilay) <= 0.)
then
205 cld_frac(icol,ilay) = tracer(icol,ilay,i_cldtot)
212 cld_frac(icol,ilay) = tracer(icol,ilay,i_cldtot)
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)
229 if (imp_physics == imp_physics_thompson)
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)
239 if (imfdeepcnv == imfdeepcnv_gf)
then
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)
247 if (imfdeepcnv == imfdeepcnv_samf)
then
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)
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
266 if (lmfdeep2) alpha0 = 200.
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)
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
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
299 cldfra2d(icol) = 0._kind_phys
301 cldfra2d(icol) = max(cldfra2d(icol), cld_frac(icol,ilay))
305 precip_frac(1:ncol,1:nlev) = cld_frac(1:ncol,1:nlev)
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)
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.
544 logical,
intent(in) :: &
546 real(kind_phys),
intent(in) :: &
547 con_g, & !< Physical constant: gravitational constant
548 con_ttp, & !< Triple point temperature of water (K)
550 real(kind_phys),
dimension(:),
intent(in) :: &
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)
560 real(kind_phys),
dimension(:,:),
intent(in),
optional :: &
562 real(kind_phys),
dimension(:,:),
intent(in) :: &
564 real(kind_phys),
dimension(:,:,:),
intent(in) :: &
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
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
584 cld_condensate(1:ncol,1:nlev,1) = tracer(1:ncol,1:nlev,i_cldliq)
585 cld_condensate(1:ncol,1:nlev,2) = tracer(1:ncol,1:nlev,i_cldice)
587 cld_condensate(1:ncol,1:nlev,3) = tracer(1:ncol,1:nlev,i_cldrain)
588 cld_condensate(1:ncol,1:nlev,4) = tracer(1:ncol,1:nlev,i_cldsnow) + &
589 tracer(1:ncol,1:nlev,i_cldgrpl)
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)
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)
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)
620 cld_rerain(icol,ilay) = rerain_def
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))
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
640 cld_reice(icol,ilay) = (1250.0/9.387) * tem3 ** 0.031
642 cld_reice(icol,ilay) = max(10.0, min(cld_reice(icol,ilay), 150.0))
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)
668 logical,
intent(in),
optional :: &
669 cond_cfrac_onRH, & !< If true, cloud-fracion set to unity when rh>99%
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.
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
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)
695 real(kind_phys),
dimension(:,:),
intent(in) :: &
697 real(kind_phys),
dimension(:,:,:),
intent(in) :: &
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
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
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
719 cld_condensate(1:ncol,1:nlev,1) = tracer(1:ncol,1:nlev,i_cldliq)
720 cld_condensate(1:ncol,1:nlev,2) = tracer(1:ncol,1:nlev,i_cldice)
721 cld_condensate(1:ncol,1:nlev,3) = tracer(1:ncol,1:nlev,i_cldrain)
722 cld_condensate(1:ncol,1:nlev,4) = tracer(1:ncol,1:nlev,i_cldsnow)
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
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)
746 if (
present(cond_cfrac_onrh) .and. relhum(icol,ilay) > 0.99)
then
747 cld_frac(icol,ilay) = 1._kind_phys
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)
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)
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
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)
828 integer,
intent(in) :: nLev, nCol, i_cldliq, i_cldice, i_cldsnow, i_cldice_nc, &
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
837 real(kind_phys),
dimension(:,:),
intent(inout) :: effrin_cldliq, effrin_cldice, &
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
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))
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
864 if (nint(lsmask(icol)) == 1)
then
865 nc_mp(icol,ilay) = nt_c_l*orho
867 nc_mp(icol,ilay) = nt_c_o*orho
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
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, &
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))
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