146 &, prsl_i, prsi_i, phil, phii &
147 &, omega_i, QLLS_i, QLCN_i, QILS_i, QICN_i&
148 &, lwheat_i, swheat_i, w_upi, cf_upi &
149 &, FRLAND, ZPBL, CNV_MFD_i &
150 &, CNV_DQLDT_i, CLCN_i, u_i, v_i &
152 &, TAUOROX, TAUOROY, CNV_FICE_i &
153 &, CNV_NDROP_i,CNV_NICE_i, q_io, lwm_o &
154 &, qi_o, t_io, rn_o, sr_o &
155 &, ncpl_io, ncpi_io, fprcp, rnw_io, snw_io&
156 &, qgl_io, ncpr_io, ncps_io, ncgl_io &
157 &, CLLS_io, KCBL, rainmin &
158 &, CLDREFFL, CLDREFFI, CLDREFFR, CLDREFFS &
159 &, CLDREFFG, ntrcaer, aerfld_i &
160 &, naai_i, npccn_i, iccn &
162 &, alf_fac, qc_min, pdfflag &
163 &, kdt, xlat, xlon, rhc_i, &
190 integer,
parameter :: kp = kind_phys
191 real(kind=kind_phys),
intent(in ) :: rainmin
193 integer,
parameter :: ncolmicro = 1
194 integer,
intent(in) :: im, lm, kdt, fprcp, pdfflag, iccn, ntrcaer
195 logical,
intent(in) :: flipv, skip_macro
196 real (kind=kind_phys),
intent(in):: dt_i, alf_fac, qc_min(:)
198 real (kind=kind_phys),
dimension(:,:),
intent(in) :: &
199 & prsl_i,u_i,v_i,phil, omega_i, qlls_i,qils_i, &
201 real (kind=kind_phys),
dimension(:,0:),
intent(in):: prsi_i, phii
202 real (kind=kind_phys),
dimension(:,:),
intent(in),
optional :: &
203 & cnv_dqldt_i, clcn_i, qlcn_i, qicn_i, &
204 & cnv_mfd_i, cf_upi, cnv_fice_i, cnv_ndrop_i, &
207 real (kind=kind_phys),
dimension(:,:),
intent(in) :: &
208 & rhc_i, naai_i, npccn_i
209 real (kind=kind_phys),
dimension(:,:,:),
intent(in) :: &
211 real (kind=kind_phys),
dimension(:),
intent(in):: taugwx, &
212 & taugwy, tauorox, tauoroy, frland,zpbl,xlat,xlon
217 real (kind=kind_phys),
dimension(:,:),
intent(out) :: lwm_o, qi_o
218 real (kind=kind_phys),
dimension(:,:),
intent(out),
optional :: &
219 cldreffl, cldreffi, cldreffr, cldreffs, cldreffg
220 real (kind=kind_phys),
dimension(:),
intent(out) :: rn_o, sr_o
221 character(len=*),
intent(out) :: errmsg
222 integer,
intent(out) :: errflg
226 integer,
dimension(:),
intent(inout):: kcbl
227 real (kind=kind_phys),
dimension(:,:),
intent(inout):: q_io, t_io, &
229 real (kind=kind_phys),
dimension(:,:),
intent(inout),
optional :: &
230 rnw_io, snw_io, ncpr_io, ncps_io, qgl_io, ncgl_io, ncpl_io, &
237 integer kcldtopcvn,i,k,ll, kbmin, naux, nbincontactdust,l
238 integer,
dimension(im) :: kct
239 real (kind=kind_phys) t_ice_all, use_av_v,bkgtau,lccirrus, &
240 & npre_frac, nct, wct, fcn, ksa1, tauxr8, dt_moist, dt_kp, tem, &
241 & tmaxll, usurf,lts_up, lts_low, min_exp, fracover, c2_gw, est3
243 real(kind=kind_phys),
allocatable,
dimension(:,:) :: &
244 & cnv_mfd, cnv_fice,cnv_ndrop,cnv_nice
247 real(kind=kind_phys),
dimension(IM,LM)::ncpl,ncpi,omega,sc_ice, &
248 & rad_cf, radheat,q1,u1,v1, plo, zlo, temp, &
249 & qlls, qlcn, qils,qicn, cnv_cvw,cnv_updf, &
252 & nhet_nuc, nlim_nuc, cdnc_nuc,inc_nuc, dnhet_imm, &
253 & nhet_imm,nhet_dep,nhet_dhf,dust_imm,dust_dep, dust_dhf,wsub, &
254 & sigw_gw,sigw_cnv,sigw_turb, &
256 & rnw,snw,ncpr,ncps,qgl,ncgl, &
259 & fqa,ql_tot,qi_tot,blk_l,rhc
261 real(kind=kind_phys) :: qcntot, qtot
263 real(kind=kind_phys),
dimension(IM,LM):: cnv_dqldt, clcn, clls
272 real(kind=kind_phys),
allocatable,
dimension(:,:) :: &
292 real(kind=kind_phys),
dimension(IM,LM+1) :: zet
293 real(kind=kind_phys),
dimension(IM,0:LM) :: ple, kh
298 real(kind=kind_phys),
dimension(LM) :: rhdfdar8, rhu00r8, &
299 & ttendr8,qtendr8, cwtendr8,npre8, npccninr8,ter8, &
300 & plevr8,ndropr8,qir8,qcr8,wparc_turb,qvr8, nir8,ncr8, &
301 & nimmr8,nsootr8,rnsootr8,omegr8,qrr8,qsr8,nrr8,nsr8, &
304 real(kind=kind_phys),
dimension(1:LM,10) :: rndstr8,naconr8
308 real(kind=kind_phys),
dimension(IM) :: ls_snr, ls_prc2
311 real(kind=kind_phys),
dimension (LM) :: uwind_gw,vwind_gw, &
312 & tm_gw, pm_gw, nm_gw, h_gw, rho_gw, khaux, qcaux, &
313 & dummyw , wparc_cgw, cfaux, dpre8, &
314 & wparc_ls,wparc_gw, swparc,smaxliq,smaxicer8,nheticer8, &
315 & nhet_immr8,dnhet_immr8,nhet_depr8,nhet_dhfr8,sc_icer8, &
316 & dust_immr8,dust_depr8,dust_dhfr8,nlimicer8,cldfr8,liqcldfr8, &
317 & icecldfr8,cldor8, pdelr8, &
318 & rpdelr8,lc_turb,zmr8,ficer8,rate1ord_cw2pr, tlatr8, qvlatr8, &
319 & qctendr8, qitendr8, nctendr8, nitendr8, effcr8, effc_fnr8, &
320 & effir8, nevaprr8, evapsnowr8, prainr8, &
321 & prodsnowr8, cmeoutr8, deffir8, pgamradr8, lamcradr8,qsoutr8, &
322 & qroutr8,droutr8, qcsevapr8,qisevapr8, qvresr8, &
323 & cmeioutr8, dsoutr8, qcsinksum_rate1ord,qrtend,nrtend, &
324 & qstend, nstend, alphar8, rhr8, &
326 & qgtend, ngtend, qgoutr8, ngoutr8, dgoutr8
328 real(kind=kind_phys),
dimension(1) :: prectr8, precir8
330 real(kind=kind_phys),
dimension (LM) :: vtrmcr8,vtrmir8, &
331 & qcsedtenr8,qisedtenr8, praor8,prcor8,mnucccor8, mnucctor8, &
332 & msacwior8,psacwsor8, bergsor8,bergor8,meltor8, homoor8, &
334 & prcior8, praior8,qiresor8, mnuccror8,pracsor8, meltsdtr8, &
336 & ncalr8, ncair8, mnuccdor8, nnucctor8, nsoutr8, nroutr8, &
337 & nnuccdor8, nnucccor8,naair8, &
338 & nsacwior8, nsubior8, nprcior8, npraior8, npccnor8, npsacwsor8, &
339 & nsubcor8, npraor8, nprc1or8, tlatauxr8,pfrz_inc_kp,sadice, &
340 & sadsnow, am_evp_st, reff_rain, reff_snow, &
341 & umr,ums,qrsedten,qssedten,refl,arefl,areflz,frefl,csrfl, &
342 & acsrfl,fcsrfl,rercld,qrout2,qsout2,nrout2,nsout2,drout2, &
343 & dsout2,freqs,freqr,nfice,qcrat,prer_evap, &
345 & reff_grau, umg, qgsedtenr8, mnuccrior8, &
346 & pracgr8, psacwgr8, pgsacwr8, pgracsr8, prdgr8, qmultgr8,&
347 & qmultrgr8, psacrr8, npracgr8, nscngr8, ngracsr8, nmultgr8,&
348 & nmultrgr8, npsacwgr8, qgout2, ngout2, dgout2, freqg
350 real(kind=kind_phys),
dimension (0:LM) :: pi_gw, rhoi_gw, &
353 real(kind=kind_phys),
dimension(LM+1) :: pintr8, kkvhr8
354 real(kind=kind_phys),
dimension(2:LM+1) :: lflx, iflx, rflx, &
359 real (kind=kind_phys),
parameter :: disp_liu=1.0_kp &
364 &, ninstr8 = 0.1e6_kp &
365 &, ncnstr8 = 100.0e6_kp
367 real(kind=kind_phys):: k_gw, maxkh, tausurf_gw, overscale, tx1, rh1_kp
368 real(kind=kind_phys):: t_ice_denom
370 integer,
dimension(1) :: lev_sed_strt
371 real(kind=kind_phys),
parameter :: sig_sed_strt=0.05_kp
373 real(kind=kind_phys),
dimension(3) :: ccn_diag
374 real(kind=kind_phys),
dimension(58) :: cloudparams
376 integer,
parameter :: ccn_param=2, in_param=5
378 real(kind=kind_phys),
parameter ::fdust_drop=1.0_kp, fsoot_drop=0.1_kp &
379 &, sigma_nuc_kp=0.28_kp,sclmfdfr=0.03_kp
382 type (
aerprops),
dimension (IM,LM) :: aeroprops
383 type (
aerprops) :: aeroaux, aeroaux_b
384 real,
allocatable,
dimension(:,:,:) :: aermassmix
386 logical :: use_average_v, ltrue, lprint, lprnt
394 data use_av_v/1.0_kp/, bkgtau/0.015_kp/, lccirrus/500.0_kp/, npre_frac/1.0_kp/, &
395 & tmaxll/296.0_kp/, fracover/1.0_kp/, lts_low/12.0_kp/, lts_up/24.0_kp/, &
399 & 10.0_kp, 4.0_kp , 4.0_kp , 1.0_kp , 2.e-3_kp, 8.e-4_kp, 2.0_kp , 1.0_kp , -1.0_kp &
400 &, 0.0_kp , 1.3_kp , 1.0e-9_kp, 3.3e-4_kp, 20.0_kp , 4.8_kp , 4.8_kp , 230.0_kp , 1.0_kp &
401 &, 1.0_kp , 230.0_kp, 14400._kp, 50.0_kp , 0.01_kp , 0.1_kp , 200.0_kp, 0.0_kp , 0.0_kp &
402 &, 0.5_kp , 0.5_kp , 2000.0_kp, 0.8_kp , 0.5_kp , -40.0_kp, 1.0_kp , 4.0_kp , 0.0_kp &
403 &, 0.0_kp , 0.0_kp , 1.0e-3_kp, 8.0e-4_kp, 1.0_kp , 0.95_kp , 1.0_kp , 0.0_kp , 900.0_kp&
406 &, 1.0_kp , 1.0_kp , 1.0_kp , 0.0_kp , 0.0_kp , 1.e-5_kp, 2.e-5_kp, 2.1e-5_kp, 4.e-5_kp&
408 &, 3e-5_kp, 0.1_kp , 4.0_kp , 150.0_kp/
427 omega(i,k) = omega_i(i,ll)
428 ncpl(i,k) = ncpl_io(i,ll)
429 ncpi(i,k) = ncpi_io(i,ll)
430 rnw(i,k) = rnw_io(i,ll)
431 snw(i,k) = snw_io(i,ll)
432 qgl(i,k) = qgl_io(i,ll)
433 ncpr(i,k) = ncpr_io(i,ll)
434 ncps(i,k) = ncps_io(i,ll)
435 ncgl(i,k) = ncgl_io(i,ll)
437 qlls(i,k) = qlls_i(i,ll)-qlcn_i(i,ll)
438 qlcn(i,k) = qlcn_i(i,ll)
439 qils(i,k) = qils_i(i,ll)-qicn_i(i,ll)
440 qicn(i,k) = qicn_i(i,ll)
441 cnv_cvw(i,k) = w_upi(i,ll)
442 cnv_updf(i,k) = cf_upi(i,ll)
443 cnv_dqldt(i,k) = cnv_dqldt_i(i,ll)
444 clcn(i,k) = clcn_i(i,ll)
445 clls(i,k) = max(clls_io(i,ll)-clcn_i(i,ll),zero)
446 plo(i,k) = prsl_i(i,ll)*0.01_kp
447 zlo(i,k) = phil(i,ll) * onebg
448 temp(i,k) = t_io(i,ll)
449 radheat(i,k) = lwheat_i(i,ll) + swheat_i(i,ll)
450 rhc(i,k) = rhc_i(i,ll)
452 cdnc_nuc(i,k) = npccn_i(i,ll)
453 inc_nuc(i,k) = naai_i(i,ll)
461 ple(i,k) = prsi_i(i,ll) * 0.01_kp
462 zet(i,k+1) = phii(i,ll) * onebg
465 if (.not. skip_macro)
then
467 allocate(cnv_mfd(im,lm), cnv_fice(im,lm) &
468 &, cnv_ndrop(im,lm), cnv_nice(im,lm))
472 cnv_mfd(i,k) = cnv_mfd_i(i,ll)
474 cnv_fice(i,k) = cnv_fice_i(i,ll)
475 cnv_ndrop(i,k) = cnv_ndrop_i(i,ll)
476 cnv_nice(i,k) = cnv_nice_i(i,ll)
487 omega(i,k) = omega_i(i,k)
488 ncpl(i,k) = ncpl_io(i,k)
489 ncpi(i,k) = ncpi_io(i,k)
490 rnw(i,k) = rnw_io(i,k)
491 snw(i,k) = snw_io(i,k)
492 qgl(i,k) = qgl_io(i,k)
493 ncpr(i,k) = ncpr_io(i,k)
494 ncps(i,k) = ncps_io(i,k)
495 ncgl(i,k) = ncgl_io(i,k)
497 qlls(i,k) = qlls_i(i,k)-qlcn_i(i,k)
498 qlcn(i,k) = qlcn_i(i,k)
499 qils(i,k) = qils_i(i,k)-qicn_i(i,k)
500 qicn(i,k) = qicn_i(i,k)
501 cnv_cvw(i,k) = w_upi(i,k)
502 cnv_updf(i,k) = cf_upi(i,k)
503 cnv_dqldt(i,k) = cnv_dqldt_i(i,k)
504 clcn(i,k) = clcn_i(i,k)
505 clls(i,k) = max(clls_io(i,k)-clcn_i(i,k),zero)
506 plo(i,k) = prsl_i(i,k)*0.01_kp
507 zlo(i,k) = phil(i,k) * onebg
508 temp(i,k) = t_io(i,k)
509 radheat(i,k) = lwheat_i(i,k) + swheat_i(i,k)
510 rhc(i,k) = rhc_i(i,k)
512 cdnc_nuc(i,k) = npccn_i(i,k)
513 inc_nuc(i,k) = naai_i(i,k)
520 ple(i,k) = prsi_i(i,k) * 0.01_kp
521 zet(i,k+1) = phii(i,k) * onebg
524 if (.not. skip_macro)
then
526 allocate(cnv_mfd(im,lm), cnv_fice(im,lm) &
527 &, cnv_ndrop(im,lm), cnv_nice(im,lm))
530 cnv_mfd(i,k) = cnv_mfd_i(i,k)
532 cnv_fice(i,k) = cnv_fice_i(i,k)
533 cnv_ndrop(i,k) = cnv_ndrop_i(i,k)
534 cnv_nice(i,k) = cnv_nice_i(i,k)
553 & qils(i,k), clls(i,k), qlcn(i,k), &
554 & qicn(i,k), clcn(i,k), ncpl(i,k), &
556 if (rnw(i,k) <= qc_min(1))
then
559 elseif (ncpr(i,k) <= nmin)
then
560 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
562 if (snw(i,k) <= qc_min(2))
then
565 elseif (ncps(i,k) <= nmin)
then
566 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
568 if (qgl(i,k) <= qc_min(2))
then
571 elseif (ncgl(i,k) <= nmin)
then
572 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
580 kcbl(i) = max(lm-kcbl(i),10)
586 If ((cnv_dqldt(i,k) <= 1.0e-9_kp) .and. &
587 & (cnv_dqldt(i,k+1) > 1.0e-9_kp))
then
666 tx1 =
half * (temp(i,l+1) + temp(i,l))
667 kh(i,l) = 3.55e-7_kp*tx1**2.5_kp*(rgas*0.01_kp) / ple(i,l)
672 kh(i,lm) = kh(i,lm-1)
676 blk_l(i,l) =
one / (
one/max(0.15_kp*zpbl(i),0.4_kp*zlo(i,lm-1))&
677 & +
one/(zlo(i,l)*0.4_kp) )
680 ncpl(i,l) = max( ncpl(i,l), zero)
681 ncpi(i,l) = max( ncpi(i,l), zero)
682 rad_cf(i,l) = max(zero, min(clls(i,l)+clcn(i,l),
one))
691 t_ice_all = cloudparams(33) + tice
692 t_ice_denom =
one / (tice - t_ice_all)
696 rhdfdar8(l) = 1.e-8_kp
707 rndstr8(l,k) = 2.0e-7_kp
738 allocate(aermassmix(im,lm,15))
740 aermassmix(:,:,1:ntrcaer) = aerfld_i(:,:,1:ntrcaer)
742 aermassmix(:,:,1:5) = 1.0e-6_kp
743 aermassmix(:,:,6:15) = 2.0e-14_kp
747 deallocate(aermassmix)
749 use_average_v = .false.
750 if (use_av_v > zero)
then
751 use_average_v = .true.
754 k_gw = (pi+pi) / lccirrus
761 tausurf_gw = min(
half*sqrt(tauorox(i)*tauorox(i) &
762 & + tauoroy(i)*tauoroy(i)), 10.0_kp)
765 uwind_gw(k) = min(
half*sqrt( u1(i,k)*u1(i,k) &
766 & + v1(i,k)*v1(i,k)), 50.0_kp)
771 pm_gw(k) = 100.0_kp*plo(i,k)
775 rho_gw(k) = pm_gw(k) /(rgas*tm_gw(k))
778 plevr8(k) = 100.0_kp*plo(i,k)
779 ndropr8(k) = ncpl(i,k)
780 qir8(k) = qils(i,k) + qicn(i,k)
781 qcr8(k) = qlls(i,k) + qlcn(i,k)
789 if (rad_cf(i,k) > 0.01_kp .and. qir8(k) > zero)
then
790 npre8(k) = npre_frac*ncpi(i,k)
795 omegr8(k) = omega(i,k)
796 lc_turb(k) = max(blk_l(i,k), 50.0_kp)
799 if (npre8(k) > zero .and. qir8(k) > zero)
then
800 dpre8(k) = ( qir8(k)/(6.0_kp*npre8(k)*900.0_kp*pi))**(
one/3.0_kp)
804 wparc_ls(k) = -omegr8(k) / (rho_gw(k)*grav) &
805 & + cpbg * radheat(i,k)
809 pi_gw(k) = 100.0_kp*ple(i,k)
821 call gw_prof (1, lm, 1, tm_gw, pm_gw, pi_gw, rhoi_gw, ni_gw, &
822 & ti_gw, nm_gw, q1(i,:))
825 nm_gw(k) = max(nm_gw(k), 0.005_kp)
826 h_gw(k) = k_gw*rho_gw(k)*uwind_gw(k)*nm_gw(k)
827 if (h_gw(k) > zero)
then
828 h_gw(k) = sqrt(2.0_kp*tausurf_gw/h_gw(k))
831 wparc_gw(k) = k_gw*uwind_gw(k)*h_gw(k)*0.133_kp
838 if (kcldtopcvn > 20)
then
841 nct = nm_gw(kcldtopcvn)
842 wct = max(cnv_cvw(i,kcldtopcvn), zero)
844 fcn = maxval(cnv_updf(i,kcldtopcvn:lm))
847 c2_gw = (nm_gw(k) + nct) / nct
848 wparc_cgw(k) = sqrt(ksa1*fcn*fcn*12.56_kp* &
849 & 1.806_kp*c2_gw*c2_gw)*wct*0.133_kp
855 dummyw(k) = 0.133_kp*k_gw*uwind_gw(k)/nm_gw(k)
859 if (wparc_cgw(k)+wparc_gw(k) > dummyw(k))
then
872 kbmin = min(kbmin, lm-1) - 4
874 wparc_turb(k) = kh(i,k) / lc_turb(k)
878 if (frland(i) < 0.1_kp .and. zpbl(i) < 800.0_kp .and. &
879 & temp(i,lm) < 298.0_kp .and. temp(i,lm) > 274.0_kp)
then
881 dummyw(k) = max(min((zet(i,k+1)-zpbl(i))*0.01_kp, 10.0_kp),-10.0_kp)
882 dummyw(k) =
one / (
one+exp(dummyw(k)))
884 maxkh = max(maxval(kh(i,kbmin:lm-1)*nm_gw(kbmin:lm-1)/ &
887 wparc_turb(k) = (
one-dummyw(k))*wparc_turb(k) &
893 wparc_turb(kbmin:lm) = max(wparc_turb(kbmin:lm), 0.2_kp)
900 swparc(k) = sqrt(wparc_gw(k) * wparc_gw(k) &
901 & + wparc_turb(k) * wparc_turb(k) &
902 & + wparc_cgw(k) * wparc_cgw(k))
911 if (plevr8(k) > 70.0_kp)
then
913 ccn_diag(1) = 0.001_kp
914 ccn_diag(2) = 0.004_kp
915 ccn_diag(3) = 0.01_kp
917 if (k > 2 .and. k <= lm-2)
then
918 tauxr8 = (ter8(k-1) + ter8(k+1) + ter8(k)) *
oneb3
924 aeroaux = aeroprops(i, k)
930 pfrz_inc_kp(k) = zero
939 & wparc_ls(k), aeroaux, npre8(k), dpre8(k), ccn_diag, &
940 & ndropr8(k), npccninr8(k), smaxliq(k), &
942 & naair8(k), smaxicer8(k), nheticer8(k), nhet_immr8(k), &
943 & dnhet_immr8(k), nhet_depr8(k), nhet_dhfr8(k), &
944 & sc_icer8(k), dust_immr8(k), dust_depr8(k), &
945 & dust_dhfr8(k), nlimicer8(k), use_average_v, &
946 & ccn_param, in_param, fdust_drop, &
947 & fsoot_drop,pfrz_inc_kp(k),sigma_nuc_kp, rh1_kp, &
952 if (npccninr8(k) < 1.0e-12_kp) npccninr8(k) = zero
972 dnhet_immr8(k) = zero
983 nhet_nuc(i,k) = nheticer8(k) * 1.0e-6_kp
984 nlim_nuc(i,k) = nlimicer8(k) * 1.0e-6_kp
985 sc_ice(i,k) = min(max(sc_icer8(k),
one),2.0_kp)
993 if(temp(i,k) < t_ice_all)
then
995 sc_ice(i,k) = max(sc_ice(i,k), 1.5_kp)
996 elseif(temp(i,k) > tice)
then
997 sc_ice(i,k) = rhc(i,k)
1001 tx1 = max(sc_ice(i,k), 1.5_kp)
1002 sc_ice(i,k) = ((tice-temp(i,k))*tx1 + (temp(i,k)-t_ice_all)*rhc(i,k)) &
1007 cdnc_nuc(i,k) = npccninr8(k)
1008 inc_nuc(i,k) = naair8(k)
1010 nhet_imm(i,k) = max(nhet_immr8(k), zero)
1011 dnhet_imm(i,k) = max(dnhet_immr8(k), zero)
1012 nhet_dep(i,k) = nhet_depr8(k) * 1.0e-6_kp
1013 nhet_dhf(i,k) = nhet_dhfr8(k) * 1.0e-6_kp
1014 dust_imm(i,k) = max(dust_immr8(k), zero)*1.0e-6_kp
1015 dust_dep(i,k) = max(dust_depr8(k), zero)*1.0e-6_kp
1016 dust_dhf(i,k) = max(dust_dhfr8(k), zero)*1.0e-6_kp
1017 wsub(i,k) = wparc_ls(k) + swparc(k)*0.8_kp
1018 sigw_gw(i,k) = wparc_gw(k) * wparc_gw(k)
1019 sigw_cnv(i,k) = wparc_cgw(k) * wparc_cgw(k)
1020 sigw_turb(i,k) = wparc_turb(k) * wparc_turb(k)
1054 if (.not. skip_macro)
then
1067 & alpht_x(im,lm), pfrz(im,lm))
1099 call macro_cloud (im, lm, dt_moist, alf_fac, plo, ple, &
1104 & temp, q1, qlls, qlcn, qils, qicn, &
1108 & cloudparams, sclmfdfr, &
1116 & cnv_fice, cnv_ndrop, cnv_nice, &
1117 & sc_ice, ncpl, ncpi, pfrz, &
1118 & lprnt, ipr, rhc, pdfflag, qc_min)
1130 if (cnv_mfd(i,k) > 1.0e-6_kp)
then
1131 tx1 =
one / cnv_mfd(i,k)
1132 cnv_ndrop(i,k) = cnv_ndrop(i,k) * tx1
1133 cnv_nice(i,k) = cnv_nice(i,k) * tx1
1135 cnv_ndrop(i,k) = zero
1136 cnv_nice(i,k) = zero
1139 rad_cf(i,k) = min(clls(i,k)+clcn(i,k),
one)
1142 if (pfrz(i,k) > zero)
then
1143 inc_nuc(i,k) = inc_nuc(i,k) * pfrz(i,k)
1144 nhet_nuc(i,k) = nhet_nuc(i,k) * pfrz(i,k)
1147 nhet_nuc(i,k) = zero
1169 call meltfrz_inst(im, lm, temp, qlls, qlcn, qils, qicn, ncpl, ncpi)
1200 qcntot = qlcn(i,k) + qicn(i,k)
1201 ql_tot(i,k) = qlcn(i,k) + qlls(i,k)
1202 qi_tot(i,k) = qicn(i,k) + qils(i,k)
1204 if (ql_tot(i,k) < zero)
then
1205 q1(i,k) = q1(i,k) + ql_tot(i,k)
1206 temp(i,k) = temp(i,k) - lvbcp*ql_tot(i,k)
1209 if (qi_tot(i,k) < zero)
then
1210 q1(i,k) = q1(i,k) + qi_tot(i,k)
1211 temp(i,k) = temp(i,k) - lsbcp*qi_tot(i,k)
1214 qtot = ql_tot(i,k) + qi_tot(i,k)
1215 if (qtot > zero)
then
1216 fqa(i,k) = min(max(qcntot / qtot, zero),
one)
1239 rndstr8(k,l) = 2.0e-7_kp
1248 tx1 = min(clls(i,k) + clcn(i,k),
one)
1249 if (tx1 > zero)
then
1250 cldfr8(k) = min(max(tx1, 0.00001_kp),
one)
1255 if (temp(i,k) > tice)
then
1256 liqcldfr8(k) = cldfr8(k)
1258 elseif (temp(i,k) <= t_ice_all)
then
1260 icecldfr8(k) = cldfr8(k)
1262 icecldfr8(k) = cldfr8(k) * (tice - temp(i,k))/(tice-t_ice_all)
1263 liqcldfr8(k) = cldfr8(k) - icecldfr8(k)
1267 cldor8(k) = cldfr8(k)
1271 qcr8(k) = ql_tot(i,k)
1272 qir8(k) = qi_tot(i,k)
1273 ncr8(k) = max(ncpl(i,k), zero)
1274 nir8(k) = max(ncpi(i,k), zero)
1278 nrr8(k) = max(ncpr(i,k), zero)
1279 nsr8(k) = max(ncps(i,k), zero)
1280 ngr8(k) = max(ncgl(i,k), zero)
1283 naair8(k) = inc_nuc(i,k)
1284 npccninr8(k) = cdnc_nuc(i,k)
1286 if (cldfr8(k) >= 0.001_kp)
then
1287 nimmr8(k) = min(dnhet_imm(i,k),ncr8(k)/(cldfr8(k)*dt_moist))
1294 aeroaux = aeroprops(i, k)
1300 naux = aeroaux_b%nmods
1301 if (nbincontactdust < naux)
then
1302 nbincontactdust = naux
1304 naconr8(k, 1:naux) = aeroaux_b%num(1:naux)
1305 rndstr8(k, 1:naux) = aeroaux_b%dpg(1:naux) *
half
1314 pdelr8(k) = (ple(i,k) - ple(i,k-1)) * 100.0_kp
1315 rpdelr8(k) =
one / pdelr8(k)
1316 plevr8(k) = 100.0_kp * plo(i,k)
1318 ficer8(k) = qir8(k) / (qcr8(k)+qir8(k) + 1.0e-10_kp)
1319 omegr8(k) = wsub(i,k)
1326 pintr8(k) = ple(i,k-1) * 100.0_kp
1327 kkvhr8(k) = kh(i,k-1)
1331 tx1 =
one / pintr8(lm+1)
1333 if (plevr8(k)*tx1 < sig_sed_strt)
then
1337 lev_sed_strt(1) = max(lm/6, min(lm/3,lev_sed_strt(1)))
1361 if (fprcp <= 0)
then
1367 nsootr8(k) = sum(aeroaux_b%num)
1368 naux = aeroaux_b%nmods
1369 rnsootr8(k) = sum(aeroaux_b%dpg(1:naux))/naux
1373 & dt_kp, ter8, ttendr8, &
1374 & ncolmicro, lm , qvr8, &
1375 & qtendr8, cwtendr8, qcr8, qir8, ncr8, nir8, &
1376 & abs(fprcp), qrr8, qsr8, nrr8, nsr8, &
1377 & plevr8, pdelr8, cldfr8, liqcldfr8, &
1378 & icecldfr8, cldor8, pintr8, &
1379 & rpdelr8, zmr8, rate1ord_cw2pr, &
1380 & naair8, npccninr8, &
1381 & rndstr8, naconr8, rhdfdar8, rhu00r8, ficer8, &
1382 & tlatr8, qvlatr8, qctendr8, &
1383 & qitendr8, nctendr8, nitendr8, effcr8, &
1384 & effc_fnr8, effir8, prectr8, precir8, &
1385 & nevaprr8, evapsnowr8, &
1386 & prainr8, prodsnowr8, cmeoutr8, &
1387 & deffir8, pgamradr8, lamcradr8, &
1388 & qsoutr8, dsoutr8, qroutr8, droutr8, &
1389 & qcsevapr8, qisevapr8, qvresr8, &
1390 & cmeioutr8, vtrmcr8, vtrmir8, &
1391 & qcsedtenr8, qisedtenr8, praor8, prcor8, &
1392 & mnucccor8, mnucctor8, msacwior8, &
1393 & psacwsor8, bergsor8, bergor8, &
1394 & meltor8, homoor8, qcresor8, prcior8, &
1395 & praior8, qiresor8, mnuccror8, &
1396 & pracsor8, meltsdtr8, frzrdtr8, ncalr8, &
1397 & ncair8, mnuccdor8, &
1398 & nnucctor8, nsoutr8, nroutr8, &
1399 & ncnstr8, ninstr8, nimmr8, disp_liu, &
1400 & nsootr8, rnsootr8, ui_scale, dcrit, &
1401 & nnuccdor8, nnucccor8, &
1402 & nsacwior8, nsubior8, nprcior8, &
1403 & npraior8, npccnor8, npsacwsor8, &
1404 & nsubcor8, npraor8, nprc1or8, &
1405 & tlatauxr8, nbincontactdust, &
1406 & lprnt, xlat(i), xlon(i), rhr8)
1411 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1412 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1416 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1417 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1418 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1421 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*
onebcp
1423 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k) * dt_kp, zero)
1424 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k) * dt_kp, zero)
1430 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1431 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1432 cldreffr(i,k) = max(droutr8(k)*0.5_kp*1.0e6_kp, 150.0_kp)
1433 cldreffs(i,k) = max(0.192_kp*dsoutr8(k)*0.5_kp*1.0e6_kp, 250.0_kp)
1437 elseif (fprcp == 1)
then
1444 lprint = lprnt .and. i == ipr
1459 call micro_mg_tend2_0 ( &
1460 & ncolmicro, lm, dt_kp, &
1469 & cldfr8, liqcldfr8, icecldfr8, rhr8, &
1470 & qcsinksum_rate1ord, &
1471 & naair8, npccninr8, &
1472 & rndstr8, naconr8, &
1473 & tlatr8, qvlatr8, &
1474 & qctendr8, qitendr8, &
1475 & nctendr8, nitendr8, &
1478 & effcr8, effc_fnr8, effir8, &
1479 & sadice, sadsnow, &
1480 & prectr8, precir8, &
1481 & nevaprr8, evapsnowr8, &
1483 & prainr8, prodsnowr8, &
1484 & cmeoutr8, deffir8, &
1485 & pgamradr8, lamcradr8, &
1486 & qsoutr8, dsoutr8, &
1488 & rflx,
sflx, qroutr8, &
1489 & reff_rain, reff_snow, &
1490 & qcsevapr8, qisevapr8, qvresr8, &
1491 & cmeioutr8, vtrmcr8, vtrmir8, &
1493 & qcsedtenr8, qisedtenr8, &
1494 & qrsedten, qssedten, &
1496 & mnucccor8, mnucctor8, msacwior8, &
1497 & psacwsor8, bergsor8, bergor8, &
1498 & meltor8, homoor8, &
1499 & qcresor8, prcior8, praior8, &
1500 & qiresor8, mnuccror8, pracsor8, &
1501 & meltsdtr8, frzrdtr8, mnuccdor8, &
1502 & nroutr8, nsoutr8, &
1503 & refl, arefl, areflz, &
1504 & frefl, csrfl, acsrfl, &
1512 & prer_evap, xlat(i), xlon(i), lprint, iccn, &
1515 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1516 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1518 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1519 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1520 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1521 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*
onebcp
1522 rnw(i,k) = rnw(i,k) + qrtend(k)*dt_kp
1523 snw(i,k) = snw(i,k) + qstend(k)*dt_kp
1525 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k)*dt_kp, zero)
1526 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k)*dt_kp, zero)
1527 ncpr(i,k) = max(ncpr(i,k) + nrtend(k)*dt_kp, zero)
1528 ncps(i,k) = max(ncps(i,k) + nstend(k)*dt_kp, zero)
1530 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1531 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1532 cldreffr(i,k) = max(reff_rain(k), 150.0_kp)
1533 cldreffs(i,k) = max(reff_snow(k), 250.0_kp)
1544 cldreffl(i,k) = 10.0_kp
1545 cldreffi(i,k) = 50.0_kp
1546 cldreffr(i,k) = 1000.0_kp
1547 cldreffs(i,k) = 250.0_kp
1556 lprint = lprnt .and. i == ipr
1572 call micro_mg_tend3_0 ( &
1573 & ncolmicro, lm, dt_kp, &
1583 & cldfr8, liqcldfr8, icecldfr8, rhr8, &
1584 & qcsinksum_rate1ord, &
1585 & naair8, npccninr8, &
1586 & rndstr8, naconr8, &
1587 & tlatr8, qvlatr8, &
1588 & qctendr8, qitendr8, &
1589 & nctendr8, nitendr8, &
1595 & effcr8, effc_fnr8, effir8, &
1596 & sadice, sadsnow, &
1597 & prectr8, precir8, &
1598 & nevaprr8, evapsnowr8, &
1600 & prainr8, prodsnowr8, &
1601 & cmeoutr8, deffir8, &
1602 & pgamradr8, lamcradr8, &
1603 & qsoutr8, dsoutr8, &
1605 & qgoutr8, ngoutr8, dgoutr8, &
1607 & lflx, iflx, gflx, &
1609 & rflx,
sflx, qroutr8, &
1611 & reff_rain, reff_snow, reff_grau, &
1613 & qcsevapr8, qisevapr8, qvresr8, &
1614 & cmeioutr8, vtrmcr8, vtrmir8, &
1617 & umg, qgsedtenr8, &
1619 & qcsedtenr8, qisedtenr8, &
1620 & qrsedten, qssedten, &
1622 & mnucccor8, mnucctor8, msacwior8, &
1623 & psacwsor8, bergsor8, bergor8, &
1624 & meltor8, homoor8, &
1625 & qcresor8, prcior8, praior8, &
1627 & qiresor8, mnuccror8, mnuccrior8, pracsor8, &
1629 & meltsdtr8, frzrdtr8, mnuccdor8, &
1631 & pracgr8, psacwgr8, pgsacwr8, &
1632 & pgracsr8, prdgr8, &
1633 & qmultgr8, qmultrgr8, psacrr8, &
1634 & npracgr8, nscngr8, ngracsr8, &
1635 & nmultgr8, nmultrgr8, npsacwgr8, &
1637 & nroutr8, nsoutr8, &
1638 & refl, arefl, areflz, &
1639 & frefl, csrfl, acsrfl, &
1646 & qgout2, ngout2, dgout2, freqg, &
1649 & prer_evap, xlat(i), xlon(i), lprint, iccn, &
1652 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1653 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1655 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1656 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1657 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1658 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*
onebcp
1659 rnw(i,k) = rnw(i,k) + qrtend(k)*dt_kp
1660 snw(i,k) = snw(i,k) + qstend(k)*dt_kp
1661 qgl(i,k) = qgl(i,k) + qgtend(k)*dt_kp
1663 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k)*dt_kp, zero)
1664 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k)*dt_kp, zero)
1665 ncpr(i,k) = max(ncpr(i,k) + nrtend(k)*dt_kp, zero)
1666 ncps(i,k) = max(ncps(i,k) + nstend(k)*dt_kp, zero)
1667 ncgl(i,k) = max(ncgl(i,k) + ngtend(k)*dt_kp, zero)
1669 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1670 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1671 cldreffr(i,k) = max(reff_rain(k), 150.0_kp)
1672 cldreffs(i,k) = max(reff_snow(k), 250.0_kp)
1673 cldreffg(i,k) = max(reff_grau(k), 250.0_kp)
1684 cldreffl(i,k) = 10.0_kp
1685 cldreffi(i,k) = 50.0_kp
1686 cldreffr(i,k) = 1000.0_kp
1687 cldreffs(i,k) = 250.0_kp
1688 cldreffg(i,k) = 250.0_kp
1698 if (skip_macro)
then
1701 qlcn(i,k) = ql_tot(i,k) * fqa(i,k)
1702 qlls(i,k) = ql_tot(i,k) - qlcn(i,k)
1703 qicn(i,k) = qi_tot(i,k) * fqa(i,k)
1704 qils(i,k) = qi_tot(i,k) - qicn(i,k)
1707 & qils(i,k), clls(i,k), qlcn(i,k), &
1708 & qicn(i,k), clcn(i,k), ncpl(i,k), &
1709 & ncpi(i,k), qc_min)
1711 ql_tot(i,k) = qlls(i,k) + qlcn(i,k)
1712 qi_tot(i,k) = qils(i,k) + qicn(i,k)
1713 if (rnw(i,k) <= qc_min(1))
then
1716 elseif (ncpr(i,k) <= nmin)
then
1717 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
1719 if (snw(i,k) <= qc_min(2))
then
1722 elseif (ncps(i,k) <= nmin)
then
1723 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1725 if (qgl(i,k) <= qc_min(2))
then
1728 elseif (ncgl(i,k) <= nmin)
then
1729 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1736 qlcn(i,k) = ql_tot(i,k) * fqa(i,k)
1737 qlls(i,k) = ql_tot(i,k) - qlcn(i,k)
1738 qicn(i,k) = qi_tot(i,k) * fqa(i,k)
1739 qils(i,k) = qi_tot(i,k) - qicn(i,k)
1744 call update_cld(im, lm, dt_moist, alpht_x, qc_min &
1745 &, pdfflag, plo, q1, qlls, qlcn &
1746 &, qils, qicn, temp, clls, clcn &
1747 &, sc_ice, ncpi, ncpl)
1753 ql_tot(i,k) = qlls(i,k) + qlcn(i,k)
1754 qi_tot(i,k) = qils(i,k) + qicn(i,k)
1756 if (rnw(i,k) <= qc_min(1))
then
1759 elseif (ncpr(i,k) <= nmin)
then
1760 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
1762 if (snw(i,k) <= qc_min(2))
then
1765 elseif (ncps(i,k) <= nmin)
then
1766 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1768 if (qgl(i,k) <= qc_min(2))
then
1771 elseif (ncgl(i,k) <= nmin)
then
1772 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1776 deallocate(cnv_mfd,cnv_fice,cnv_ndrop,cnv_nice)
1786 if (qi_tot(i,k) <= zero) ncpi(i,k) = zero
1787 if (ql_tot(i,k) <= zero) ncpl(i,k) = zero
1801 t_io(i,k) = temp(i,ll)
1802 q_io(i,k) = q1(i,ll)
1803 ncpi_io(i,k) = ncpi(i,ll)
1804 ncpl_io(i,k) = ncpl(i,ll)
1805 rnw_io(i,k) = rnw(i,ll)
1806 snw_io(i,k) = snw(i,ll)
1807 qgl_io(i,k) = qgl(i,ll)
1808 ncpr_io(i,k) = ncpr(i,ll)
1809 ncps_io(i,k) = ncps(i,ll)
1810 ncgl_io(i,k) = ncgl(i,ll)
1811 lwm_o(i,k) = ql_tot(i,ll)
1812 qi_o(i,k) = qi_tot(i,ll)
1815 if (skip_macro)
then
1819 clls_io(i,k) = max(zero, min(clls(i,ll)+clcn(i,ll),
one))
1826 clls_io(i,k) = clls(i,ll)
1833 t_io(i,k) = temp(i,k)
1835 ncpi_io(i,k) = ncpi(i,k)
1836 ncpl_io(i,k) = ncpl(i,k)
1837 rnw_io(i,k) = rnw(i,k)
1838 snw_io(i,k) = snw(i,k)
1839 qgl_io(i,k) = qgl(i,k)
1840 ncpr_io(i,k) = ncpr(i,k)
1841 ncps_io(i,k) = ncps(i,k)
1842 ncgl_io(i,k) = ncgl(i,k)
1843 lwm_o(i,k) = ql_tot(i,k)
1844 qi_o(i,k) = qi_tot(i,k)
1847 if (skip_macro)
then
1850 clls_io(i,k) = max(zero, min(clls(i,k)+clcn(i,k),
one))
1856 clls_io(i,k) = clls(i,k)
1863 tx1 = ls_prc2(i) + ls_snr(i)
1864 rn_o(i) = tx1 * dt_i * 0.001_kp
1866 if (rn_o(i) < rainmin)
then
1869 sr_o(i) = max(zero, min(
one, ls_snr(i)/tx1))
1873 if (
allocated(alpht_x))
deallocate (alpht_x)