21 ltp, imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, me, ncnd, ntrac, &
22 num_p3d, npdf3d, xr_cnvcld, &
23 ncnvcld3d,ntqv, ntcw,ntiw, ntlnc, ntinc, ntrnc, ntsnc, ntccn, top_at_1,&
24 ntrw, ntsw, ntgl, nthl, ntwa, ntoz, ntsmoke, ntdust, ntcoarsepm, &
25 ntclamt, nleffr, nieffr, nseffr, lndp_type, kdt, &
26 ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, &
27 ntss3, ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm, &
28 imp_physics,imp_physics_nssl, nssl_ccn_on, nssl_invertccn, &
29 imp_physics_thompson, imp_physics_gfdl, imp_physics_zhao_carr, &
30 imp_physics_zhao_carr_pdf, imp_physics_mg, imp_physics_wsm6, &
31 imp_physics_fer_hires, iovr, iovr_rand, iovr_maxrand, iovr_max, &
32 iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, idcor_hogan, &
33 idcor_oreopoulos, dcorr_con, julian, yearlen, lndp_var_list, lsswr, &
34 lslwr, ltaerosol, mraerosol, lgfdlmprad, uni_cld, effr_in, do_mynnedmf,&
35 lmfshal, lcnorm, lmfdeep2, lcrick, fhswr, fhlwr, solhr, sup, con_eps, &
36 epsm1, fvirt, rog, rocp, con_rd, xlat_d, xlat, xlon, coslat, sinlat, &
37 tsfc, slmsk, prsi, prsl, prslk, tgrs, sfc_wts, mg_cld, effrr_in, &
38 pert_clds, sppt_wts, sppt_amp, cnvw_in, cnvc_in, qgrs, aer_nm, dx, &
39 icloud, iaermdl, iaerflg, con_pi, con_g, con_ttp, con_thgni, si, & !inputs from here and above
40 coszen, coszdg, effrl_inout, effri_inout, effrs_inout, &
41 clouds1, clouds2, clouds3, clouds4, clouds5, qci_conv, & !in/out from here and above
42 kd, kt, kb, mtopa, mbota, raddt, tsfg, tsfa, de_lgth, alb1d, delp, dz, & !output from here and below
43 plvl, plyr, tlvl, tlyr, qlyr, olyr, gasvmr_co2, gasvmr_n2o, gasvmr_ch4,&
44 gasvmr_o2, gasvmr_co, gasvmr_cfc11, gasvmr_cfc12, gasvmr_cfc22, &
45 gasvmr_ccl4, gasvmr_cfc113, aerodp,ext550, clouds6, clouds7, clouds8, &
46 clouds9, cldsa, cldfra, cldfra2d, lwp_ex,iwp_ex, lwp_fc,iwp_fc, &
47 faersw1, faersw2, faersw3, faerlw1, faerlw2, faerlw3, alpha, rrfs_sd, &
48 aero_dir_fdb, fdb_coef, spp_wts_rad, spp_rad, ico2, ozphys, &
53 use radcons,
only: itsfc, qmin, qme5, qme6, epsq, prsmin
77 re_qc_min, re_qc_max, &
78 re_qi_min, re_qi_max, &
88 integer,
intent(in) :: im, levs, lm, lmk, lmp, ltp, &
89 n_var_lndp, imfdeepcnv, &
90 imfdeepcnv_gf, imfdeepcnv_c3, &
92 num_p3d, npdf3d, ncnvcld3d, ntqv, &
93 ntcw, ntiw, ntlnc, ntinc, &
95 ntrw, ntsw, ntgl, nthl, ntwa, ntoz, &
96 ntsmoke, ntdust, ntcoarsepm, &
97 ntclamt, nleffr, nieffr, nseffr, &
100 imp_physics_thompson, &
102 imp_physics_zhao_carr, &
103 imp_physics_zhao_carr_pdf, &
104 imp_physics_mg, imp_physics_wsm6, &
106 imp_physics_fer_hires, &
107 yearlen, icloud, iaermdl, iaerflg
109 integer,
intent(in) :: &
123 integer,
intent(in) :: ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, ntss3, &
124 ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm
126 character(len=3),
dimension(:),
intent(in),
optional :: lndp_var_list
128 logical,
intent(in) :: lsswr, lslwr, ltaerosol, lgfdlmprad, &
129 uni_cld, effr_in, do_mynnedmf, &
130 lmfshal, lmfdeep2, pert_clds, lcrick,&
131 lcnorm, top_at_1, lextop, mraerosol
132 logical,
intent(in) :: rrfs_sd, aero_dir_fdb, xr_cnvcld
134 logical,
intent(in) :: nssl_ccn_on, nssl_invertccn
135 integer,
intent(in) :: spp_rad
136 real(kind_phys),
intent(in),
optional :: spp_wts_rad(:,:)
138 real(kind=kind_phys),
intent(in) :: fhswr, fhlwr, solhr, sup, julian, sppt_amp, dcorr_con
139 real(kind=kind_phys),
intent(in) :: con_eps, epsm1, fvirt, rog, rocp, con_rd, con_pi, con_g, con_ttp, con_thgni
141 real(kind=kind_phys),
dimension(:),
intent(in) :: xlat_d, xlat, xlon, &
142 coslat, sinlat, tsfc, &
145 real(kind=kind_phys),
dimension(:,:),
intent(in) :: prsi, prsl, prslk, &
147 real(kind=kind_phys),
dimension(:,:),
intent(in),
optional :: &
152 real(kind=kind_phys),
dimension(:,:,:),
intent(in) :: qgrs
153 real(kind=kind_phys),
dimension(:,:,:),
intent(inout) :: aer_nm
155 real(kind=kind_phys),
dimension(:),
intent(inout) :: coszen, coszdg
157 real(kind=kind_phys),
dimension(:,:),
intent(inout),
optional :: effrl_inout, &
160 real(kind=kind_phys),
dimension(:,:),
intent(inout) :: clouds1, &
163 real(kind=kind_phys),
dimension(:,:),
intent(in),
optional :: qci_conv
164 real(kind=kind_phys),
dimension(:),
intent(in) :: fdb_coef
165 real(kind=kind_phys),
dimension(:),
intent(out) :: lwp_ex,iwp_ex, &
168 integer,
intent(out) :: kd, kt, kb
170 integer,
dimension(:,:),
intent(out) :: mbota, mtopa
172 real(kind=kind_phys),
intent(out) :: raddt
174 real(kind=kind_phys),
dimension(:),
intent(out) :: tsfg, tsfa
175 real(kind=kind_phys),
dimension(:),
intent(out) :: de_lgth, &
178 real(kind=kind_phys),
dimension(:,:),
intent(out) :: delp, dz, &
182 real(kind=kind_phys),
dimension(:,:),
intent(out) :: plvl, tlvl
184 real(kind=kind_phys),
dimension(:,:),
intent(out) :: gasvmr_co2, &
194 real(kind=kind_phys),
dimension(:,:),
intent(out) :: aerodp
195 real(kind=kind_phys),
dimension(:,:),
intent(out) :: ext550
196 real(kind=kind_phys),
dimension(:,:),
intent(out) :: clouds6, &
201 real(kind=kind_phys),
dimension(:),
intent(out) :: cldfra2d
202 real(kind=kind_phys),
dimension(:,:),
intent(out) :: cldsa
204 real(kind=kind_phys),
dimension(:,:,:),
intent(out) :: faersw1,&
208 real(kind=kind_phys),
dimension(:,:,:),
intent(out) :: faerlw1,&
211 real(kind=kind_phys),
dimension(:,:),
intent(out) :: alpha
212 character(len=*),
intent(out) :: errmsg
213 integer,
intent(out) :: errflg
218 integer :: i, j, k, k1, k2, lsk, lv, n, itop, ibtc, lp1, lla, llb, lya,lyb
220 real(kind=kind_phys) :: es, qs, delt, tem0d, pfac
221 real(kind=kind_phys),
dimension(im) :: gridkm
223 real(kind=kind_phys),
dimension(im) :: cvt1, cvb1, tem1d, tskn, xland
225 real(kind=kind_phys),
dimension(im,lm+LTP) :: &
226 htswc, htlwc, gcice, grain, grime, htsw0, htlw0, &
227 rhly, tvly,qstl, vvel, clw, ciw, prslk1, tem2da, &
228 dzb, hzb, cldcov, deltaq, cnvc, cnvw, &
229 effrl, effri, effrr, effrs, rho, orho, plyrpa
232 real(kind=kind_phys),
dimension(im,lm+LTP) :: &
233 qv_mp, qc_mp, qi_mp, qs_mp, &
235 real (kind=kind_phys),
dimension(lm) :: cldfra1d, qv1d, &
236 & qc1d, qi1d, qs1d, dz1d, p1d, t1d
239 real(kind=kind_phys),
dimension(im,lm+LTP+1) :: tem2db, hz
241 real(kind=kind_phys),
dimension(im,lm+LTP,min(4,ncnd)) :: ccnd
242 real(kind=kind_phys),
dimension(im,lm+LTP,2:ntrac) :: tracer1
243 real(kind=kind_phys),
dimension(im,lm+LTP) :: &
244 & cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, &
245 & cld_rwp, cld_rerain, cld_swp, cld_resnow
246 real(kind=kind_phys),
dimension(im,lm+LTP,NF_VGAS) :: gasvmr
247 real(kind=kind_phys),
dimension(im,lm+LTP,NBDSW,NF_AESW) :: faersw
248 real(kind=kind_phys),
dimension(im,lm+LTP,NBDLW,NF_AELW) :: faerlw
251 real(kind=kind_phys),
dimension(im) :: cldp1d
252 real (kind=kind_phys) :: alpha0,beta0,m,s,cldtmp,tmp_wt,cdfz
253 real (kind=kind_phys) :: max_relh
260 integer :: ids, ide, jds, jde, kds, kde, &
261 ims, ime, jms, jme, kms, kme, &
262 its, ite, jts, jte, kts, kte
264 real(kind=kind_phys) :: qvs
272 if (.not. (lsswr .or. lslwr))
return
279 if (imp_physics == imp_physics_thompson)
then
286 gridkm(i) = dx(i)*0.001
297 if (.not. top_at_1)
then
316 if (.not. top_at_1)
then
325 raddt = min(fhswr, fhlwr)
332 if ( itsfc == 0 )
then
349 if (top_at_1 .and. lm < levs) lsk = levs - lm
356 plvl(i,k1+kb) = prsi(i,k2+kb) * 0.01
357 plyr(i,k1) = prsl(i,k2) * 0.01
358 tlyr(i,k1) = tgrs(i,k2)
359 prslk1(i,k1) = prslk(i,k2)
360 rho(i,k1) = prsl(i,k2)/(con_rd*tlyr(i,k1))
361 orho(i,k1) = 1.0/rho(i,k1)
364 es = min( prsl(i,k2), fpvs( tgrs(i,k2) ) )
365 qs = max( qmin, con_eps * es / (prsl(i,k2) + epsm1*es) )
366 rhly(i,k1) = max( 0.0, min( 1.0, max(qmin, qgrs(i,k2,ntqv))/qs ) )
376 tracer1(:,k1,j) = max(0.0, qgrs(:,k2,j))
385 plvl(i,k2) = 0.01 * prsi(i,1+kb)
386 plyr(i,k1) = 0.5 * (plvl(i,k2+1) + plvl(i,k2))
387 prslk1(i,k1) = (plyr(i,k1)*0.001) ** rocp
392 plvl(i,k1) = prsi(i,1) * 0.01
399 plvl(i,k1+1) = 0.01 * prsi(i,levs+1)
400 plyr(i,k1) = 0.5 * (plvl(i,k1+1) + plvl(i,k1))
401 prslk1(i,k1) = (plyr(i,k1)*0.001) ** rocp
406 plvl(i,k1) = prsi(i,lp1) * 0.01
414 if ( plvl(i,lla) <= prsmin ) plvl(i,lla) = 2.0*prsmin
415 plyr(i,lyb) = 0.5 * plvl(i,lla)
416 tlyr(i,lyb) = tlyr(i,lya)
417 prslk1(i,lyb) = (plyr(i,lyb)*0.001) ** rocp
418 rho(i,lyb) = plyr(i,lyb) *100.0/(con_rd*tlyr(i,lyb))
419 orho(i,lyb) = 1.0/rho(i,lyb)
420 rhly(i,lyb) = rhly(i,lya)
421 qstl(i,lyb) = qstl(i,lya)
425 tracer1(:,lyb,:) = tracer1(:,lya,:)
434 olyr(i,k) = max( qmin, tracer1(i,k,ntoz) )
438 call ozphys%run_o3clim(xlat, prslk1, con_pi, olyr)
443 call coszmn (xlon,sinlat,coslat,solhr,im,me, &
462 call getgases (plvl, xlon, xlat, im, lmk, ico2, top_at_1,&
468 gasvmr_co2(i,k) = gasvmr(i,k,1)
469 gasvmr_n2o(i,k) = gasvmr(i,k,2)
470 gasvmr_ch4(i,k) = gasvmr(i,k,3)
471 gasvmr_o2(i,k) = gasvmr(i,k,4)
472 gasvmr_co(i,k) = gasvmr(i,k,5)
473 gasvmr_cfc11(i,k) = gasvmr(i,k,6)
474 gasvmr_cfc12(i,k) = gasvmr(i,k,7)
475 gasvmr_cfc22(i,k) = gasvmr(i,k,8)
476 gasvmr_ccl4(i,k) = gasvmr(i,k,9)
477 gasvmr_cfc113(i,k) = gasvmr(i,k,10)
484 tem2da(i,k) = log( plyr(i,k) )
485 tem2db(i,k) = log( plvl(i,k) )
492 tem2da(i,1) = log( plyr(i,1) )
493 tem2db(i,1) = log( max(prsmin, plvl(i,1)) )
494 tem2db(i,lmp) = log( plvl(i,lmp) )
495 tsfa(i) = tlyr(i,lmk)
496 tlvl(i,1) = tlyr(i,1)
497 tlvl(i,lmp) = tskn(i)
503 qlyr(i,k1) = max( tem1d(i), qgrs(i,k,ntqv) )
504 tem1d(i) = min( qme5, qlyr(i,k1) )
505 tvly(i,k1) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1))
506 delp(i,k1) = plvl(i,k1+1) - plvl(i,k1)
512 qlyr(i,lyb) = qlyr(i,lya)
513 tvly(i,lyb) = tvly(i,lya)
514 delp(i,lyb) = plvl(i,lla) - plvl(i,llb)
520 tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k)) &
521 & * (tem2db(i,k) - tem2da(i,k)) &
522 & / (tem2da(i,k-1) - tem2da(i,k))
535 dz(i,k) = tem0d * (tem2db(i,k+1) - tem2db(i,k)) * tvly(i,k)
540 hz(i,k) = hz(i,k+1) + dz(i,k)
544 pfac = (tem2db(i,k+1) - tem2da(i,k)) / (tem2db(i,k+1) - tem2db(i,k))
545 hzb(i,k) = hz(i,k+1) + pfac * (hz(i,k) - hz(i,k+1))
549 dzb(i,k) = hzb(i,k) - hzb(i,k+1)
551 dzb(i,lmk) = hzb(i,lmk) - hz(i,lmp)
558 tem2da(i,1) = log( plyr(i,1) )
559 tem2db(i,1) = log( plvl(i,1) )
560 tem2db(i,lmp) = log( max(prsmin, plvl(i,lmp)) )
563 tlvl(i,lmp) = tlyr(i,lmk)
568 qlyr(i,k) = max( tem1d(i), qgrs(i,k,ntqv) )
569 tem1d(i) = min( qme5, qlyr(i,k) )
570 tvly(i,k) = tgrs(i,k) * (1.0 + fvirt*qlyr(i,k))
571 delp(i,k) = plvl(i,k) - plvl(i,k+1)
577 qlyr(i,lyb) = qlyr(i,lya)
578 tvly(i,lyb) = tvly(i,lya)
579 delp(i,lyb) = plvl(i,lla) - plvl(i,llb)
585 tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k)) &
586 & * (tem2db(i,k+1) - tem2da(i,k)) &
587 & / (tem2da(i,k+1) - tem2da(i,k))
600 dz(i,k) = tem0d * (tem2db(i,k) - tem2db(i,k+1)) * tvly(i,k)
605 hz(i,k+1) = hz(i,k) + dz(i,k)
609 pfac = (tem2db(i,k) - tem2da(i,k)) / (tem2db(i,k) - tem2db(i,k+1))
610 hzb(i,k) = hz(i,k) + pfac * (hz(i,k+1) - hz(i,k))
614 dzb(i,k) = hzb(i,k) - hzb(i,k-1)
616 dzb(i,1) = hzb(i,1) - hz(i,1)
625 if (ntchm>0 .and. iaermdl==2)
then
628 aer_nm(i,k,1) = qgrs(i,k,ntdu1)*1.e-9_kind_phys
629 aer_nm(i,k,2) = qgrs(i,k,ntdu2)*1.e-9_kind_phys
630 aer_nm(i,k,3) = qgrs(i,k,ntdu3)*1.e-9_kind_phys
631 aer_nm(i,k,4) = qgrs(i,k,ntdu4)*1.e-9_kind_phys
632 aer_nm(i,k,5) = qgrs(i,k,ntdu5)*1.e-9_kind_phys
633 aer_nm(i,k,6) = qgrs(i,k,ntss1)*1.e-9_kind_phys
634 aer_nm(i,k,7) = qgrs(i,k,ntss2)*1.e-9_kind_phys
635 aer_nm(i,k,8) = qgrs(i,k,ntss3)*1.e-9_kind_phys
636 aer_nm(i,k,9) = qgrs(i,k,ntss4)*1.e-9_kind_phys
637 aer_nm(i,k,10) = qgrs(i,k,ntss5)*1.e-9_kind_phys
638 aer_nm(i,k,11) = qgrs(i,k,ntsu)*1.e-9_kind_phys
639 aer_nm(i,k,12) = qgrs(i,k,ntbcb)*1.e-9_kind_phys
640 aer_nm(i,k,13) = qgrs(i,k,ntbcl)*1.e-9_kind_phys
641 aer_nm(i,k,14) = qgrs(i,k,ntocb)*1.e-9_kind_phys
642 aer_nm(i,k,15) = qgrs(i,k,ntocl)*1.e-9_kind_phys
648 if (rrfs_sd .and. aero_dir_fdb)
then
651 aer_nm(i,k,1 )=aer_nm(i,k,1 )+ qgrs(i,k,ntdust)*fdb_coef(1)*1.e-9
652 aer_nm(i,k,2 )=aer_nm(i,k,2 )+(qgrs(i,k,ntdust)*fdb_coef(2) &
653 +qgrs(i,k,ntcoarsepm)*fdb_coef(3))*1.e-9
654 aer_nm(i,k,3 )=aer_nm(i,k,3 )+qgrs(i,k,ntcoarsepm)*fdb_coef(4)*1.e-9
655 aer_nm(i,k,4 )=aer_nm(i,k,4 )+qgrs(i,k,ntcoarsepm)*fdb_coef(5)*1.e-9
656 aer_nm(i,k,12)=aer_nm(i,k,12)+qgrs(i,k,ntsmoke)*fdb_coef(6)*1.e-9
657 aer_nm(i,k,14)=aer_nm(i,k,14)+qgrs(i,k,ntsmoke)*fdb_coef(7)*1.e-9
665 call setaer (plvl, plyr, prslk1, tvly, rhly, slmsk, &
666 tracer1, aer_nm, xlon, xlat, im, lmk, lmp,&
667 lsswr, lslwr, iaermdl, iaerflg, top_at_1, con_pi, &
668 con_rd, con_g, faersw, faerlw, aerodp, ext550, errflg, errmsg)
675 faersw1(i,k,j) = faersw(i,k,j,1)
676 faersw2(i,k,j) = faersw(i,k,j,2)
677 faersw3(i,k,j) = faersw(i,k,j,3)
686 faerlw1(i,k,j) = faerlw(i,k,j,1)
687 faerlw2(i,k,j) = faerlw(i,k,j,2)
688 faerlw3(i,k,j) = faerlw(i,k,j,3)
703 ccnd(i,k,1) = tracer1(i,k,ntcw)
706 elseif (ncnd == 2)
then
709 ccnd(i,k,1) = tracer1(i,k,ntcw)
710 ccnd(i,k,2) = tracer1(i,k,ntiw)
713 elseif (ncnd == 4)
then
716 ccnd(i,k,1) = tracer1(i,k,ntcw)
717 ccnd(i,k,2) = tracer1(i,k,ntiw)
718 ccnd(i,k,3) = tracer1(i,k,ntrw)
719 ccnd(i,k,4) = tracer1(i,k,ntsw)
722 elseif (ncnd == 5 .or. ncnd == 6)
then
725 ccnd(i,k,1) = tracer1(i,k,ntcw)
726 ccnd(i,k,2) = tracer1(i,k,ntiw)
727 ccnd(i,k,3) = tracer1(i,k,ntrw)
728 if (imp_physics == imp_physics_fer_hires )
then
731 IF ( ncnd == 5 )
THEN
732 ccnd(i,k,4) = tracer1(i,k,ntsw) + tracer1(i,k,ntgl)
733 ELSEIF ( ncnd == 6 )
THEN
734 ccnd(i,k,4) = tracer1(i,k,ntsw) + tracer1(i,k,ntgl) + tracer1(i,k,nthl)
740 if_thompson:
if (imp_physics == imp_physics_thompson .and. (ltaerosol .or. mraerosol))
then
744 qv_mp(i,k) = qvs/(1.-qvs)
745 rho(i,k) = con_eps*plyr(i,k)*100./(con_rd*tlyr(i,k)*(qv_mp(i,k)+con_eps))
746 orho(i,k) = 1.0/rho(i,k)
747 qc_mp(i,k) = tracer1(i,k,ntcw)/(1.-qvs)
748 qi_mp(i,k) = tracer1(i,k,ntiw)/(1.-qvs)
749 qs_mp(i,k) = tracer1(i,k,ntsw)/(1.-qvs)
750 nc_mp(i,k) = tracer1(i,k,ntlnc)/(1.-qvs)
751 ni_mp(i,k) = tracer1(i,k,ntinc)/(1.-qvs)
752 nwfa(i,k) = tracer1(i,k,ntwa)
755 elseif (imp_physics == imp_physics_thompson)
then
759 qv_mp(i,k) = qvs/(1.-qvs)
760 rho(i,k) = con_eps*plyr(i,k)*100./(con_rd*tlyr(i,k)*(qv_mp(i,k)+con_eps))
761 orho(i,k) = 1.0/rho(i,k)
762 qc_mp(i,k) = tracer1(i,k,ntcw)/(1.-qvs)
763 qi_mp(i,k) = tracer1(i,k,ntiw)/(1.-qvs)
764 qs_mp(i,k) = tracer1(i,k,ntsw)/(1.-qvs)
765 if(nint(slmsk(i)) == 1)
then
766 nc_mp(i,k) = nt_c_l*orho(i,k)
768 nc_mp(i,k) = nt_c_o*orho(i,k)
770 ni_mp(i,k) = tracer1(i,k,ntinc)/(1.-qvs)
778 if (ccnd(i,k,n) < epsq) ccnd(i,k,n) = 0.0
782 if (imp_physics == imp_physics_gfdl )
then
783 if (.not. lgfdlmprad)
then
793 ccnd(:,:,1) = tracer1(:,1:lmk,ntcw)
794 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntrw)
795 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntiw)
796 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntsw)
797 ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:lmk,ntgl)
801 if (ccnd(i,k,1) < epsq ) ccnd(i,k,1) = 0.0
811 cldcov(i,k1) = mg_cld(i,k)
812 effrl(i,k1) = effrl_inout(i,k)
813 effri(i,k1) = effri_inout(i,k)
814 effrr(i,k1) = effrr_in(i,k)
815 effrs(i,k1) = effrs_inout(i,k)
822 cldcov(i,k1) = mg_cld(i,k)
826 elseif (imp_physics == imp_physics_gfdl)
then
827 if ((imfdeepcnv==imfdeepcnv_gf .or. imfdeepcnv==imfdeepcnv_c3) .and. kdt>1)
then
831 if (qci_conv(i,k)>0.)
then
833 cldcov(i,k1) = clouds1(i,k1)
835 cldcov(i,k1) = tracer1(i,k1,ntclamt)
841 cldcov(1:im,1+kd:lm+kd) = tracer1(1:im,1:lm,ntclamt)
847 effrl(i,k1) = effrl_inout(i,k)
848 effri(i,k1) = effri_inout(i,k)
849 effrr(i,k1) = effrr_in(i,k)
850 effrs(i,k1) = effrs_inout(i,k)
865 elseif (imp_physics == imp_physics_nssl )
then
871 effrl(i,k1) = effrl_inout(i,k)
872 effri(i,k1) = effri_inout(i,k)
873 effrr(i,k1) = effrr_in(i,k)
874 effrs(i,k1) = effrs_inout(i,k)
881 elseif (imp_physics == imp_physics_thompson)
then
888 if ((ltaerosol .or. mraerosol) .and. qc_mp(i,k)>1.e-12 .and. nc_mp(i,k)<100.)
then
891 if (qi_mp(i,k)>1.e-12 .and. ni_mp(i,k)<100.)
then
892 ni_mp(i,k) =
make_icenumber(qi_mp(i,k)*rho(i,k), tlyr(i,k)) * orho(i,k)
898 islmsk = nint(slmsk(i))
903 call calc_effectrad (tlyr(i,:), plyr(i,:)*100., qv_mp(i,:), qc_mp(i,:), &
904 nc_mp(i,:), qi_mp(i,:), ni_mp(i,:), qs_mp(i,:), &
905 effrl(i,:), effri(i,:), effrs(i,:), islmsk, 1, lm )
908 effrl(i,k) = max(re_qc_min, min(effrl(i,k), re_qc_max))*1.e6
909 effri(i,k) = max(re_qi_min, min(effri(i,k), re_qi_max))*1.e6
910 effrs(i,k) = max(re_qs_min, min(effrs(i,k), re_qs_max))*1.e6
912 effrl(i,lmk) = re_qc_min*1.e6
913 effri(i,lmk) = re_qi_min*1.e6
914 effrs(i,lmk) = re_qs_min*1.e6
921 effrl_inout(i,k) = effrl(i,k1)
922 effri_inout(i,k) = effri(i,k1)
923 effrs_inout(i,k) = effrs(i,k1)
937 if ((num_p3d == 4) .and. (npdf3d == 3))
then
944 cnvw(i,k1) = cnvw_in(i,k)
945 cnvc(i,k1) = cnvc_in(i,k)
948 elseif ((npdf3d == 0) .and. (ncnvcld3d == 1))
then
953 cnvw(i,k1) = cnvw_in(i,k)
967 if (imp_physics == imp_physics_zhao_carr)
then
968 ccnd(1:im,1:lmk,1) = ccnd(1:im,1:lmk,1) + cnvw(1:im,1:lmk)
973 & ( plyr, plvl, tlyr, tvly, qlyr, qstl, rhly, &
974 & ccnd, ncndl, cnvw, cnvc, tracer1, &
975 & xlat, xlon, slmsk, dz, delp, im, lm, lmk, lmp, &
976 & deltaq, sup, dcorr_con, me, icloud, kdt, &
977 & ntrac, ntcw, ntiw, ntrw, ntsw, ntgl, ntclamt, &
978 & imp_physics, imp_physics_nssl, imp_physics_fer_hires, &
979 & imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6, &
980 & imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, &
981 & imp_physics_mg, iovr, iovr_rand, iovr_maxrand, iovr_max, &
982 & iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, &
983 & idcor_hogan, idcor_oreopoulos, lcrick, lcnorm, &
984 & imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, do_mynnedmf, &
985 & lgfdlmprad, xr_cnvcld, &
986 & uni_cld, lmfshal, lmfdeep2, cldcov, clouds1, &
987 & effrl, effri, effrr, effrs, effr_in, &
988 & effrl_inout, effri_inout, effrs_inout, &
989 & lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
990 & dzb, xlat_d, julian, yearlen, gridkm, top_at_1, si, &
991 & con_ttp, con_pi, con_g, con_rd, con_thgni, &
992 & cld_frac, cld_lwp, cld_reliq, cld_iwp, cld_reice, &
993 & cld_rwp, cld_rerain, cld_swp, cld_resnow, &
994 & cldsa, mtopa, mbota, de_lgth, alpha &
1002 tmp_wt= -1*log( ( 2.0 / ( sppt_wts(i,38) ) ) - 1 )
1010 if (m<0.99 .AND. m > 0.01)
then
1011 s = sppt_amp*m*(1.-m)
1012 alpha0 = m*m*(1.-m)/(s*s)-m
1013 beta0 = alpha0*(1.-m)/m
1016 call ppfbet(cldp1d(i),alpha0,beta0,iflag,cldtmp)
1017 cld_frac(i,k) = cldtmp
1026 clouds1(i,k) = cld_frac(i,k)
1027 clouds2(i,k) = cld_lwp(i,k)
1028 clouds3(i,k) = cld_reliq(i,k)
1029 clouds4(i,k) = cld_iwp(i,k)
1030 clouds5(i,k) = cld_reice(i,k)
1031 clouds6(i,k) = cld_rwp(i,k)
1032 clouds7(i,k) = cld_rerain(i,k)
1033 clouds8(i,k) = cld_swp(i,k)
1034 clouds9(i,k) = cld_resnow(i,k)
1035 cldfra(i,k) = cld_frac(i,k)
1041 cldfra2d(i) = max(cldfra2d(i), cldfra(i,k))
1045 if ( spp_rad == 1 )
then
1049 clouds3(i,k) = clouds3(i,k) - spp_wts_rad(i,k) * clouds3(i,k)
1050 clouds5(i,k) = clouds5(i,k) - spp_wts_rad(i,k) * clouds5(i,k)
1051 clouds9(i,k) = clouds9(i,k) - spp_wts_rad(i,k) * clouds9(i,k)
1055 clouds3(i,k) = clouds3(i,k) - spp_wts_rad(i,levs) * clouds3(i,k)
1056 clouds5(i,k) = clouds5(i,k) - spp_wts_rad(i,levs) * clouds5(i,k)
1057 clouds9(i,k) = clouds9(i,k) - spp_wts_rad(i,levs) * clouds9(i,k)
1069 if (lndp_type==1)
then
1071 if (lndp_var_list(k) ==
'alb')
then
1073 call cdfnor(sfc_wts(i,k),alb1d(i))
subroutine, public gfs_rrtmg_pre_run(im, levs, lm, lmk, lmp, n_var_lndp, lextop, ltp, imfdeepcnv, imfdeepcnv_gf, imfdeepcnv_c3, me, ncnd, ntrac, num_p3d, npdf3d, xr_cnvcld, ncnvcld3d, ntqv, ntcw, ntiw, ntlnc, ntinc, ntrnc, ntsnc, ntccn, top_at_1, ntrw, ntsw, ntgl, nthl, ntwa, ntoz, ntsmoke, ntdust, ntcoarsepm, ntclamt, nleffr, nieffr, nseffr, lndp_type, kdt, ntdu1, ntdu2, ntdu3, ntdu4, ntdu5, ntss1, ntss2, ntss3, ntss4, ntss5, ntsu, ntbcb, ntbcl, ntocb, ntocl, ntchm, imp_physics, imp_physics_nssl, nssl_ccn_on, nssl_invertccn, imp_physics_thompson, imp_physics_gfdl, imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, imp_physics_mg, imp_physics_wsm6, imp_physics_fer_hires, iovr, iovr_rand, iovr_maxrand, iovr_max, iovr_dcorr, iovr_exp, iovr_exprand, idcor, idcor_con, idcor_hogan, idcor_oreopoulos, dcorr_con, julian, yearlen, lndp_var_list, lsswr, lslwr, ltaerosol, mraerosol, lgfdlmprad, uni_cld, effr_in, do_mynnedmf, lmfshal, lcnorm, lmfdeep2, lcrick, fhswr, fhlwr, solhr, sup, con_eps, epsm1, fvirt, rog, rocp, con_rd, xlat_d, xlat, xlon, coslat, sinlat, tsfc, slmsk, prsi, prsl, prslk, tgrs, sfc_wts, mg_cld, effrr_in, pert_clds, sppt_wts, sppt_amp, cnvw_in, cnvc_in, qgrs, aer_nm, dx, icloud, iaermdl, iaerflg, con_pi, con_g, con_ttp, con_thgni, si, coszen, coszdg, effrl_inout, effri_inout, effrs_inout, clouds1, clouds2, clouds3, clouds4, clouds5, qci_conv, kd, kt, kb, mtopa, mbota, raddt, tsfg, tsfa, de_lgth, alb1d, delp, dz, plvl, plyr, tlvl, tlyr, qlyr, olyr, gasvmr_co2, gasvmr_n2o, gasvmr_ch4, gasvmr_o2, gasvmr_co, gasvmr_cfc11, gasvmr_cfc12, gasvmr_cfc22, gasvmr_ccl4, gasvmr_cfc113, aerodp, ext550, clouds6, clouds7, clouds8, clouds9, cldsa, cldfra, cldfra2d, lwp_ex, iwp_ex, lwp_fc, iwp_fc, faersw1, faersw2, faersw3, faerlw1, faerlw2, faerlw3, alpha, rrfs_sd, aero_dir_fdb, fdb_coef, spp_wts_rad, spp_rad, ico2, ozphys, errmsg, errflg)