MODULE module_ra_cam use module_ra_cam_support use module_cam_support, only: endrun implicit none real(r8) abarl(4) real(r8) bbarl(4) real(r8) cbarl(4) real(r8) dbarl(4) real(r8) ebarl(4) real(r8) fbarl(4) save abarl, bbarl, cbarl, dbarl, ebarl, fbarl data abarl/ 2.817e-02, 2.682e-02,2.264e-02,1.281e-02/ data bbarl/ 1.305 , 1.346 ,1.454 ,1.641 / data cbarl/-5.62e-08 ,-6.94e-06 ,4.64e-04 ,0.201 / data dbarl/ 1.63e-07 , 2.35e-05 ,1.24e-03 ,7.56e-03 / data ebarl/ 0.829 , 0.794 ,0.754 ,0.826 / data fbarl/ 2.482e-03, 4.226e-03,6.560e-03,4.353e-03/ real(r8) abari(4) real(r8) bbari(4) real(r8) cbari(4) real(r8) dbari(4) real(r8) ebari(4) real(r8) fbari(4) save abari, bbari, cbari, dbari, ebari, fbari data abari/ 3.448e-03, 3.448e-03,3.448e-03,3.448e-03/ data bbari/ 2.431 , 2.431 ,2.431 ,2.431 / data cbari/ 1.00e-05 , 1.10e-04 ,1.861e-02,.46658 / data dbari/ 0.0 , 1.405e-05,8.328e-04,2.05e-05 / data ebari/ 0.7661 , 0.7730 ,0.794 ,0.9595 / data fbari/ 5.851e-04, 5.665e-04,7.267e-04,1.076e-04/ real(r8) delta real(r8) o2mmr save delta, o2mmr data delta / 0.0014257179260883 / data o2mmr / .23143 / real(r8) frcsol(nspint) real(r8) wavmin(nspint) real(r8) wavmax(nspint) real(r8) raytau(nspint) real(r8) abh2o(nspint) real(r8) abo3 (nspint) real(r8) abco2(nspint) real(r8) abo2 (nspint) real(r8) ph2o(nspint) real(r8) pco2(nspint) real(r8) po2 (nspint) real(r8) nirwgt(nspint) save frcsol ,wavmin ,wavmax ,raytau ,abh2o ,abo3 , & abco2 ,abo2 ,ph2o ,pco2 ,po2 ,nirwgt data frcsol / .001488, .001389, .001290, .001686, .002877, & .003869, .026336, .360739, .065392, .526861, & .526861, .526861, .526861, .526861, .526861, & .526861, .006239, .001834, .001834/ data nirwgt / 0.0, 0.0, 0.0, 0.0, 0.0, & 0.0, 0.0, 0.0, 0.320518, 1.0, 1.0, & 1.0, 1.0, 1.0, 1.0, 1.0, & 1.0, 1.0, 1.0 / data wavmin / .200, .245, .265, .275, .285, & .295, .305, .350, .640, .700, .701, & .701, .701, .701, .702, .702, & 2.630, 4.160, 4.160/ data wavmax / .245, .265, .275, .285, .295, & .305, .350, .640, .700, 5.000, 5.000, & 5.000, 5.000, 5.000, 5.000, 5.000, & 2.860, 4.550, 4.550/ real(r8) v_raytau_35 real(r8) v_raytau_64 real(r8) v_abo3_35 real(r8) v_abo3_64 parameter( & v_raytau_35 = 0.155208, & v_raytau_64 = 0.0392, & v_abo3_35 = 2.4058030e+01, & v_abo3_64 = 2.210e+01 & ) data raytau / 4.020, 2.180, 1.700, 1.450, 1.250, & 1.085, 0.730, v_raytau_35, v_raytau_64, & 0.02899756, 0.01356763, 0.00537341, & 0.00228515, 0.00105028, 0.00046631, & 0.00025734, & .0001, .0001, .0001/ data abh2o / .000, .000, .000, .000, .000, & .000, .000, .000, .000, & 0.00256608, 0.06310504, 0.42287445, 2.45397941, & 11.20070807, 47.66091389, 240.19010243, & .000, .000, .000/ data abo3 /5.370e+04, 13.080e+04, 9.292e+04, 4.530e+04, 1.616e+04, & 4.441e+03, 1.775e+02, v_abo3_35, v_abo3_64, .000, & .000, .000 , .000 , .000 , .000, & .000, .000 , .000 , .000 / data abco2 / .000, .000, .000, .000, .000, & .000, .000, .000, .000, .000, & .000, .000, .000, .000, .000, & .000, .094, .196, 1.963/ data abo2 / .000, .000, .000, .000, .000, & .000, .000, .000,1.11e-05,6.69e-05, & .000, .000, .000, .000, .000, & .000, .000, .000, .000/ data ph2o / .000, .000, .000, .000, .000, & .000, .000, .000, .000, .505, & .210, .120, .070, .048, .029, & .018, .000, .000, .000/ data pco2 / .000, .000, .000, .000, .000, & .000, .000, .000, .000, .000, & .000, .000, .000, .000, .000, & .000, 1.000, .640, .360/ data po2 / .000, .000, .000, .000, .000, & .000, .000, .000, 1.000, 1.000, & .000, .000, .000, .000, .000, & .000, .000, .000, .000/ real(r8) amo save amo data amo / 48.0000 / contains subroutine camrad(RTHRATENLW,RTHRATENSW, & dolw,dosw, & SWUPT,SWUPTC,SWDNT,SWDNTC, & LWUPT,LWUPTC,LWDNT,LWDNTC, & SWUPB,SWUPBC,SWDNB,SWDNBC, & LWUPB,LWUPBC,LWDNB,LWDNBC, & swcf,lwcf,olr,cemiss,taucldc,taucldi,coszr, & GSW,GLW,XLAT,XLONG, & ALBEDO,t_phy,TSK,EMISS, & QV3D,QC3D,QR3D,QI3D,QS3D,QG3D, & ALSWVISDIR,ALSWVISDIF, & ALSWNIRDIR,ALSWNIRDIF, & SWVISDIR,SWVISDIF, & SWNIRDIR,SWNIRDIF, & sf_surface_physics, & SWDDIR,SWDDIF,SWDDNI, & F_QV,F_QC,F_QR,F_QI,F_QS,F_QG, & f_ice_phy,f_rain_phy, & p_phy,p8w,z,pi_phy,rho_phy,dz8w, & CLDFRA,XLAND,XICE,SNOW, & ozmixm,pin0,levsiz,num_months, & m_psp,m_psn,aerosolcp,aerosolcn,m_hybi0, & cam_abs_dim1, cam_abs_dim2, & paerlev,naer_c, & GMT,JULDAY,JULIAN,YR,DT,XTIME,DECLIN,SOLCON, & RADT,DEGRAD,n_cldadv, & abstot_3d, absnxt_3d, emstot_3d, & doabsems, & ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte, & coszen ) USE MODULE_RA_CLWRF_SUPPORT, ONLY : read_CAMgases USE module_wrf_error USE module_state_description, ONLY : SSIBSCHEME, CLMSCHEME IMPLICIT NONE INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte LOGICAL, INTENT(IN ) :: F_QV,F_QC,F_QR,F_QI,F_QS,F_QG LOGICAL, INTENT(INout) :: doabsems LOGICAL, INTENT(IN ) :: dolw,dosw INTEGER, INTENT(IN ) :: n_cldadv INTEGER, INTENT(IN ) :: JULDAY REAL, INTENT(IN ) :: JULIAN INTEGER, INTENT(IN ) :: YR REAL, INTENT(IN ) :: DT INTEGER, INTENT(IN ) :: levsiz, num_months INTEGER, INTENT(IN ) :: paerlev, naer_c INTEGER, INTENT(IN ) :: cam_abs_dim1, cam_abs_dim2 REAL, INTENT(IN ) :: RADT,DEGRAD, & XTIME,DECLIN,SOLCON,GMT REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), & INTENT(IN ) :: P_PHY, & P8W, & Z, & pi_PHY, & rho_PHY, & dz8w, & T_PHY, & QV3D, & QC3D, & QR3D, & QI3D, & QS3D, & QG3D, & CLDFRA REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), & INTENT(INOUT) :: RTHRATENLW, & RTHRATENSW REAL, DIMENSION( ims:ime, jms:jme ), & INTENT(IN ) :: XLAT, & XLONG, & XLAND, & XICE, & SNOW, & EMISS, & TSK, & ALBEDO REAL, DIMENSION( ims:ime, levsiz, jms:jme, num_months ), & INTENT(IN ) :: OZMIXM REAL, DIMENSION(levsiz), INTENT(IN ) :: PIN0 REAL, DIMENSION(ims:ime,jms:jme), INTENT(IN ) :: m_psp,m_psn REAL, DIMENSION(paerlev), intent(in) :: m_hybi0 REAL, DIMENSION( ims:ime, paerlev, jms:jme, naer_c ), & INTENT(IN ) :: aerosolcp, aerosolcn REAL, DIMENSION( ims:ime, jms:jme ), & INTENT(INOUT) :: GSW, GLW REAL, DIMENSION( ims:ime, jms:jme ), & INTENT(IN) :: ALSWVISDIR, & ALSWVISDIF, & ALSWNIRDIR, & ALSWNIRDIF REAL, DIMENSION( ims:ime, jms:jme ), & INTENT(OUT) :: SWVISDIR, & SWVISDIF, & SWNIRDIR, & SWNIRDIF, & SWDDIR, & SWDDNI, & SWDDIF INTEGER, INTENT(IN) :: sf_surface_physics REAL, DIMENSION( ims:ime, kms:kme, cam_abs_dim2 , jms:jme ), & INTENT(INOUT) :: abstot_3d REAL, DIMENSION( ims:ime, kms:kme, cam_abs_dim1 , jms:jme ), & INTENT(INOUT) :: absnxt_3d REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), & INTENT(INOUT) :: emstot_3d REAL, DIMENSION( ims:ime, jms:jme ), OPTIONAL, INTENT(INOUT) ::& SWUPT,SWUPTC,SWDNT,SWDNTC, & LWUPT,LWUPTC,LWDNT,LWDNTC, & SWUPB,SWUPBC,SWDNB,SWDNBC, & LWUPB,LWUPBC,LWDNB,LWDNBC REAL, DIMENSION( ims:ime, jms:jme ), & INTENT(INOUT) :: swcf, & lwcf, & olr, & coszr REAL, DIMENSION( ims:ime, kms:kme, jms:jme ) , & INTENT(OUT ) :: cemiss, & taucldc, & taucldi REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), & INTENT(IN ) :: & F_ICE_PHY, & F_RAIN_PHY real, dimension(ims:ime,jms:jme), optional, intent(in) :: coszen INTEGER :: lchnk, ncol, pcols, pver, pverp, pverr, pverrp INTEGER :: pcnst, pnats, ppcnst, i, j, k, ii, kk, kk1, m, n integer :: begchunk, endchunk integer :: nyrm, nyrp real(r8) doymodel, doydatam, doydatap, deltat, fact1, fact2 REAL :: XT24, TLOCTM, HRANG, XXLAT, oldXT24 real(r8), DIMENSION( 1:ite-its+1 ) :: coszrs, landfrac, landm, snowh, icefrac, lwups real(r8), DIMENSION( 1:ite-its+1 ) :: asdir, asdif, aldir, aldif, ps real(r8), DIMENSION( 1:ite-its+1, 1:kte-kts+1 ) :: cld, pmid, lnpmid, pdel, zm, t real(r8), DIMENSION( 1:ite-its+1, 1:kte-kts+2 ) :: pint, lnpint real(r8), DIMENSION( 1:ite-its+1, 1:kte-kts+1, n_cldadv) :: q real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: cicewp real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: cliqwp real(r8), dimension( 1:ite-its+1, 0:kte-kts+1 ) :: tauxcl real(r8), dimension( 1:ite-its+1, 0:kte-kts+1 ) :: tauxci real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: emis real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: rel real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: rei real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: pmxrgn integer , dimension( 1:ite-its+1 ) :: nmxrgn real(r8), dimension( 1:ite-its+1 ) :: fsns real(r8), dimension( 1:ite-its+1 ) :: fsnt real(r8), dimension( 1:ite-its+1 ) :: flns real(r8), dimension( 1:ite-its+1 ) :: flnt real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsup real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsupc real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdn real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdnc real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdndir real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdncdir real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdndif real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fsdncdif real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: flup real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: flupc real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fldn real(r8), dimension( 1:ite-its+1, 1:kte-kts+2 ) :: fldnc real(r8), dimension( 1:ite-its+1 ) :: swcftoa real(r8), dimension( 1:ite-its+1 ) :: lwcftoa real(r8), dimension( 1:ite-its+1 ) :: olrtoa real(r8), dimension( 1:ite-its+1 ) :: sols real(r8), dimension( 1:ite-its+1 ) :: soll real(r8), dimension( 1:ite-its+1 ) :: solsd real(r8), dimension( 1:ite-its+1 ) :: solld real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: qrs real(r8), dimension( 1:ite-its+1 ) :: fsds real(r8), dimension( 1:ite-its+1 ) :: fsdsdir real(r8), dimension( 1:ite-its+1 ) :: fsdsdif real(r8), dimension( 1:ite-its+1, 1:kte-kts+1 ) :: qrl real(r8), dimension( 1:ite-its+1 ) :: flwds real(r8), dimension( 1:ite-its+1, levsiz, num_months ) :: ozmixmj real(r8), dimension( 1:ite-its+1, levsiz ) :: ozmix real(r8), dimension(levsiz) :: pin real(r8), dimension(1:ite-its+1) :: m_psjp,m_psjn real(r8), dimension( 1:ite-its+1, paerlev, naer_c ) :: aerosoljp real(r8), dimension( 1:ite-its+1, paerlev, naer_c ) :: aerosoljn real(r8), dimension(paerlev) :: m_hybi real(r8), dimension(1:ite-its+1 ) :: clat real(r8), dimension(its:ite,kts:kte+1,kts:kte+1) :: abstot real(r8), dimension(its:ite,kts:kte,4) :: absnxt real(r8), dimension(its:ite,kts:kte+1) :: emstot CHARACTER(LEN=256) :: msgstr LOGICAL, EXTERNAL :: wrf_dm_on_monitor CHARACTER(LEN=256) :: message lchnk = 1 begchunk = ims endchunk = ime ncol = ite - its + 1 pcols= ite - its + 1 pver = kte - kts + 1 pverp= pver + 1 pverr = kte - kts + 1 pverrp= pverr + 1 ppcnst = n_cldadv pnats = 0 pcnst = ppcnst-pnats if(naer_c.ne.naer_all) then WRITE( wrf_err_message , * ) 'naer_c-1 ne naer_all ', naer_c, naer_all CALL wrf_error_fatal3("",451,& wrf_err_message ) endif nyrm = yr - yrdata(1) + 1 nyrp = nyrm + 1 doymodel = yr*365. + julian doydatam = yrdata(nyrm)*365. + 1. doydatap = yrdata(nyrp)*365. + 1. deltat = doydatap - doydatam fact1 = (doydatap - doymodel)/deltat fact2 = (doymodel - doydatam)/deltat co2vmr = (co2(nyrm)*fact1 + co2(nyrp)*fact2)*1.e-06 IF ( wrf_dm_on_monitor() ) THEN WRITE(message,*)'CAM interpolated values_____ year:',yr,' julian day:',julian call wrf_debug( 100, message) WRITE(message,*)'CAM co2vmr: ',co2vmr,' n2ovmr:',n2ovmr,' ch4vmr:',ch4vmr,' cfc11:',f11vmr,& ' cfc12:',f12vmr call wrf_debug( 100, message) ENDIF co2mmr=co2vmr*mwco2/mwdry do k=1,levsiz pin(k)=pin0(k) enddo do k=1,paerlev m_hybi(k)=m_hybi0(k) enddo if(abstot_3d(its,kts,kts,jts) .eq. 0.0 .and. .not.doabsems .and. dolw)then CALL wrf_debug(0, 'camrad lw: CAUTION: re-calculating abstot, absnxt, emstot on restart') doabsems = .true. endif do j =jts,jte if (present(coszen)) then do i=its,ite ii=i-its+1 clat(ii)=XLAT(I,J)*DEGRAD coszrs(ii)=coszen(i,j) enddo else do i = its,ite ii = i - its + 1 XT24=MOD(XTIME+RADT*0.5,1440.) TLOCTM=GMT+XT24/60.+XLONG(I,J)/15. HRANG=15.*(TLOCTM-12.)*DEGRAD XXLAT=XLAT(I,J)*DEGRAD clat(ii)=xxlat coszrs(II)=SIN(XXLAT)*SIN(DECLIN)+COS(XXLAT)*COS(DECLIN)*COS(HRANG) enddo end if do k = kts,kte kk = kte - k + kts do i = its,ite ii = i - its + 1 q(ii,kk,1) = max(1.e-10,qv3d(i,k,j)/(1.+qv3d(i,k,j))) IF ( F_QI .and. F_QC .and. F_QS ) THEN q(ii,kk,ixcldliq) = max(0.,qc3d(i,k,j)/(1.+qv3d(i,k,j))) q(ii,kk,ixcldice) = max(0.,(qi3d(i,k,j)+qs3d(i,k,j))/(1.+qv3d(i,k,j))) ELSE IF ( F_QC .and. F_QI ) THEN q(ii,kk,ixcldice) = max(0.,qi3d(i,k,j)/(1.+qv3d(i,k,j))) q(ii,kk,ixcldliq) = max(0.,qc3d(i,k,j)/(1.+qv3d(i,k,j))) ELSE IF ( F_QC .and. F_QR ) THEN q(ii,kk,ixcldliq) = 0. q(ii,kk,ixcldice) = 0. if(t_phy(i,k,j).gt.273.15)q(ii,kk,ixcldliq) = max(0.,qc3d(i,k,j)/(1.+qv3d(i,k,j))) if(t_phy(i,k,j).le.273.15)q(ii,kk,ixcldice) = max(0.,qc3d(i,k,j)/(1.+qv3d(i,k,j))) ELSE q(ii,kk,ixcldliq) = 0. q(ii,kk,ixcldice) = 0. ENDIF cld(ii,kk) = CLDFRA(I,K,J) enddo enddo do i = its,ite ii = i - its + 1 landfrac(ii) = 2.-XLAND(I,J) landm(ii) = landfrac(ii) snowh(ii) = 0.001*SNOW(I,J) icefrac(ii) = XICE(I,J) enddo do m=1,num_months-1 do k=1,levsiz do i = its,ite ii = i - its + 1 ozmixmj(ii,k,m) = ozmixm(i,k,j,m+1) enddo enddo enddo do i = its,ite ii = i - its + 1 m_psjp(ii) = m_psp(i,j) m_psjn(ii) = m_psn(i,j) enddo do n=1,naer_c do k=1,paerlev do i = its,ite ii = i - its + 1 aerosoljp(ii,k,n) = aerosolcp(i,k,j,n) aerosoljn(ii,k,n) = aerosolcn(i,k,j,n) enddo enddo enddo do i = its,ite ii = i - its + 1 lwups(ii) = stebol*EMISS(I,J)*TSK(I,J)**4 enddo do k = kts,kte+1 kk = kte - k + kts + 1 do i = its,ite ii = i - its + 1 pint(ii,kk) = p8w(i,k,j) if(k.eq.kts)ps(ii)=pint(ii,kk) lnpint(ii,kk) = log(pint(ii,kk)) enddo enddo if(.not.doabsems .and. dolw)then do kk = 1,cam_abs_dim2 do kk1 = kts,kte+1 do i = its,ite abstot(i,kk1,kk) = abstot_3d(i,kk1,kk,j) enddo enddo enddo do kk = 1,cam_abs_dim1 do kk1 = kts,kte do i = its,ite absnxt(i,kk1,kk) = absnxt_3d(i,kk1,kk,j) enddo enddo enddo do kk = kts,kte+1 do i = its,ite emstot(i,kk) = emstot_3d(i,kk,j) enddo enddo endif do k = kts,kte kk = kte - k + kts do i = its,ite ii = i - its + 1 pmid(ii,kk) = p_phy(i,k,j) lnpmid(ii,kk) = log(pmid(ii,kk)) lnpint(ii,kk) = log(pint(ii,kk)) pdel(ii,kk) = pint(ii,kk+1) - pint(ii,kk) t(ii,kk) = t_phy(i,k,j) zm(ii,kk) = z(i,k,j) enddo enddo call param_cldoptics_calc(ncol, pcols, pver, pverp, pverr, pverrp, ppcnst, q, cld, landfrac, landm,icefrac, & pdel, t, ps, pmid, pint, cicewp, cliqwp, emis, rel, rei, pmxrgn, nmxrgn, snowh) SELECT CASE(sf_surface_physics) CASE (SSIBSCHEME) if (xtime .gt. 1.0) then do i = its,ite ii = i - its + 1 if (xland(i,j).lt.1.5) then asdir(ii) = ALSWVISDIR(i,j) asdif(ii) = ALSWVISDIF(i,j) aldir(ii) = ALSWNIRDIR(i,j) aldif(ii) = ALSWNIRDIF(i,j) else asdir(ii) = albedo(i,j) asdif(ii) = albedo(i,j) aldir(ii) = albedo(i,j) aldif(ii) = albedo(i,j) endif enddo else do i = its,ite ii = i - its + 1 asdir(ii) = albedo(i,j) asdif(ii) = albedo(i,j) aldir(ii) = albedo(i,j) aldif(ii) = albedo(i,j) enddo endif CASE (CLMSCHEME) if (xtime .gt. 1.0) then do i = its,ite ii = i - its + 1 if (xland(i,j).lt.1.5) then asdir(ii) = ALSWVISDIR(i,j) asdif(ii) = ALSWVISDIF(i,j) aldir(ii) = ALSWNIRDIR(i,j) aldif(ii) = ALSWNIRDIF(i,j) else asdir(ii) = albedo(i,j) asdif(ii) = albedo(i,j) aldir(ii) = albedo(i,j) aldif(ii) = albedo(i,j) endif enddo else do i = its,ite ii = i - its + 1 asdir(ii) = albedo(i,j) asdif(ii) = albedo(i,j) aldir(ii) = albedo(i,j) aldif(ii) = albedo(i,j) enddo endif CASE DEFAULT do i = its,ite ii = i - its + 1 asdir(ii) = albedo(i,j) asdif(ii) = albedo(i,j) aldir(ii) = albedo(i,j) aldif(ii) = albedo(i,j) enddo END SELECT call radctl (j,lchnk, ncol, pcols, pver, pverp, pverr, pverrp, ppcnst, pcnst, lwups, emis, pmid, & pint, lnpmid, lnpint, pdel, t, q, & cld, cicewp, cliqwp, tauxcl, tauxci, coszrs, clat, asdir, asdif, & aldir, aldif, solcon, GMT,JULDAY,JULIAN,DT,XTIME, & pin, ozmixmj, ozmix, levsiz, num_months, & m_psjp,m_psjn, aerosoljp, aerosoljn, m_hybi, paerlev, naer_c, pmxrgn, nmxrgn, & dolw, dosw, doabsems, abstot, absnxt, emstot, & fsup, fsupc, fsdn, fsdnc, & fsdndir ,fsdncdir,fsdndif ,fsdncdif, & flup, flupc, fldn, fldnc, swcftoa, lwcftoa, olrtoa, & fsns, fsnt ,flns ,flnt , & qrs, qrl, flwds, rel, rei, & sols, soll, solsd, solld, & landfrac, zm, fsds, fsdsdir, fsdsdif) do k = kts,kte kk = kte - k + kts do i = its,ite ii = i - its + 1 if(dolw)RTHRATENLW(I,K,J) = 1.e4*qrl(ii,kk)/(cpair*pi_phy(i,k,j)) if(dosw)RTHRATENSW(I,K,J) = 1.e4*qrs(ii,kk)/(cpair*pi_phy(i,k,j)) cemiss(i,k,j) = emis(ii,kk) taucldc(i,k,j) = tauxcl(ii,kk) taucldi(i,k,j) = tauxci(ii,kk) enddo enddo if(doabsems .and. dolw)then do kk = 1,cam_abs_dim2 do kk1 = kts,kte+1 do i = its,ite abstot_3d(i,kk1,kk,j) = abstot(i,kk1,kk) enddo enddo enddo do kk = 1,cam_abs_dim1 do kk1 = kts,kte do i = its,ite absnxt_3d(i,kk1,kk,j) = absnxt(i,kk1,kk) enddo enddo enddo do kk = kts,kte+1 do i = its,ite emstot_3d(i,kk,j) = emstot(i,kk) enddo enddo endif IF(PRESENT(SWUPT))THEN if(dosw)then do k = kts,kte+1 kk = kte +1 - k + kts do i = its,ite ii = i - its + 1 if(k.eq.kte+1)then swupt(i,j) = fsup(ii,kk) swuptc(i,j) = fsupc(ii,kk) swdnt(i,j) = fsdn(ii,kk) swdntc(i,j) = fsdnc(ii,kk) endif if(k.eq.kts)then swupb(i,j) = fsup(ii,kk) swupbc(i,j) = fsupc(ii,kk) swdnb(i,j) = fsdn(ii,kk) swdnbc(i,j) = fsdnc(ii,kk) endif enddo enddo endif if(dolw)then do k = kts,kte+1 kk = kte +1 - k + kts do i = its,ite ii = i - its + 1 if(k.eq.kte+1)then lwupt(i,j) = flup(ii,kk) lwuptc(i,j) = flupc(ii,kk) lwdnt(i,j) = fldn(ii,kk) lwdntc(i,j) = fldnc(ii,kk) endif if(k.eq.kts)then lwupb(i,j) = flup(ii,kk) lwupbc(i,j) = flupc(ii,kk) lwdnb(i,j) = fldn(ii,kk) lwdnbc(i,j) = fldnc(ii,kk) endif enddo enddo endif ENDIF do i = its,ite ii = i - its + 1 if(dolw)then GLW(I,J) = flwds(ii) lwcf(i,j) = lwcftoa(ii) olr(i,j) = olrtoa(ii) endif if(dosw)then GSW(I,J) = fsns(ii) swcf(i,j) = swcftoa(ii) coszr(i,j) = coszrs(ii) SWDDIR(i,j)= fsdsdir(ii) SWDDNI(i,j)= fsdsdir(ii)/coszrs(ii) SWDDIF(i,j)= fsdsdif(ii) endif enddo SELECT CASE(sf_surface_physics) CASE (SSIBSCHEME) if(dosw)then do i = its,ite ii = i - its + 1 SWVISDIR(I,J) = sols(ii) SWVISDIF(I,J) = solsd(ii) SWNIRDIR(I,J) = soll(ii) SWNIRDIF(I,J) = solld(ii) enddo endif CASE (CLMSCHEME) if(dosw)then do i = its,ite ii = i - its + 1 SWVISDIR(I,J) = sols(ii) SWVISDIF(I,J) = solsd(ii) SWNIRDIR(I,J) = soll(ii) SWNIRDIF(I,J) = solld(ii) enddo endif END SELECT enddo end subroutine camrad SUBROUTINE camradinit( & R_D,R_V,CP,G,STBOLT,EP_2,shalf,pptop, & ozmixm,pin,levsiz,XLAT,num_months, & m_psp,m_psn,m_hybi,aerosolcp,aerosolcn, & paerlev,naer_c, & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte ) USE module_wrf_error USE module_state_description IMPLICIT NONE INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte REAL, intent(in) :: pptop REAL, INTENT(IN) :: R_D,R_V,CP,G,STBOLT,EP_2 REAL, DIMENSION( kms:kme ) :: shalf INTEGER, INTENT(IN ) :: levsiz, num_months INTEGER, INTENT(IN ) :: paerlev, naer_c REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN ) :: XLAT REAL, DIMENSION( ims:ime, levsiz, jms:jme, num_months ), & INTENT(INOUT ) :: OZMIXM REAL, DIMENSION(levsiz), INTENT(INOUT ) :: PIN REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT ) :: m_psp,m_psn REAL, DIMENSION(paerlev), INTENT(INOUT ) :: m_hybi REAL, DIMENSION( ims:ime, paerlev, jms:jme, naer_c ), & INTENT(INOUT) :: aerosolcp,aerosolcn REAL(r8) :: pstd REAL(r8) :: rh2o, cpair IF ( .NOT. ALLOCATED( ksul ) ) ALLOCATE( ksul( nrh, nspint ) ) IF ( .NOT. ALLOCATED( wsul ) ) ALLOCATE( wsul( nrh, nspint ) ) IF ( .NOT. ALLOCATED( gsul ) ) ALLOCATE( gsul( nrh, nspint ) ) IF ( .NOT. ALLOCATED( ksslt ) ) ALLOCATE( ksslt( nrh, nspint ) ) IF ( .NOT. ALLOCATED( wsslt ) ) ALLOCATE( wsslt( nrh, nspint ) ) IF ( .NOT. ALLOCATED( gsslt ) ) ALLOCATE( gsslt( nrh, nspint ) ) IF ( .NOT. ALLOCATED( kcphil ) ) ALLOCATE( kcphil( nrh, nspint ) ) IF ( .NOT. ALLOCATED( wcphil ) ) ALLOCATE( wcphil( nrh, nspint ) ) IF ( .NOT. ALLOCATED( gcphil ) ) ALLOCATE( gcphil( nrh, nspint ) ) IF ( .NOT. ALLOCATED(ah2onw ) ) ALLOCATE( ah2onw(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(eh2onw ) ) ALLOCATE( eh2onw(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(ah2ow ) ) ALLOCATE( ah2ow(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(cn_ah2ow) ) ALLOCATE( cn_ah2ow(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(cn_eh2ow) ) ALLOCATE( cn_eh2ow(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(ln_ah2ow) ) ALLOCATE( ln_ah2ow(n_p, n_tp, n_u, n_te, n_rh) ) IF ( .NOT. ALLOCATED(ln_eh2ow) ) ALLOCATE( ln_eh2ow(n_p, n_tp, n_u, n_te, n_rh) ) ozncyc = .true. indirect = .true. ixcldliq = 2 ixcldice = 3 idxSUL = P_SUL idxSSLT = P_SSLT idxDUSTfirst = P_DUST1 idxOCPHO = P_OCPHO idxCARBONfirst = P_OCPHO idxBCPHO = P_BCPHO idxOCPHI = P_OCPHI idxBCPHI = P_BCPHI idxBG = P_BG idxVOLC = P_VOLC pstd = 101325.0 mwdry = 28.966 mwco2 = 44. mwh2o = 18.016 mwch4 = 16. mwn2o = 44. mwf11 = 136. mwf12 = 120. cappa = R_D/CP rair = R_D tmelt = 273.16 r_universal = 6.02214e26 * STBOLT latvap = 2.501e6 latice = 3.336e5 zvir = R_V/R_D - 1. rh2o = R_V cpair = CP epsqs = EP_2 CALL radini(G, CP, EP_2, STBOLT, pstd*10.0 ) CALL esinti(epsqs ,latvap ,latice ,rh2o ,cpair ,tmelt ) CALL oznini(ozmixm,pin,levsiz,num_months,XLAT, & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte) CALL aerosol_init(m_psp,m_psn,m_hybi,aerosolcp,aerosolcn,paerlev,naer_c,shalf,pptop, & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte) END SUBROUTINE camradinit subroutine oznint(julday,julian,dt,gmt,xtime,ozmixmj,ozmix,levsiz,num_months,pcols) IMPLICIT NONE INTEGER, INTENT(IN ) :: levsiz, num_months,pcols REAL(r8), DIMENSION( pcols, levsiz, num_months ), & INTENT(IN ) :: ozmixmj REAL, INTENT(IN ) :: XTIME,GMT INTEGER, INTENT(IN ) :: JULDAY REAL, INTENT(IN ) :: JULIAN REAL, INTENT(IN ) :: DT REAL(r8), DIMENSION( pcols, levsiz ), & INTENT(OUT ) :: ozmix REAL(r8) :: intJULIAN integer :: np1,np,nm,m,k,i integer :: IJUL integer, dimension(12) :: date_oz data date_oz/16, 45, 75, 105, 136, 166, 197, 228, 258, 289, 319, 350/ real(r8) :: cdayozp, cdayozm real(r8) :: fact1, fact2 logical :: finddate CHARACTER(LEN=256) :: msgstr intJULIAN = JULIAN + 1.0_r8 IJUL=INT(intJULIAN) intJULIAN=intJULIAN-FLOAT(IJUL) IJUL=MOD(IJUL,365) IF(IJUL.EQ.0)IJUL=365 intJULIAN=intJULIAN+IJUL np1=1 finddate=.false. do m=1,12 if(date_oz(m).gt.intjulian.and..not.finddate) then np1=m finddate=.true. endif enddo cdayozp=date_oz(np1) if(np1.gt.1) then cdayozm=date_oz(np1-1) np=np1 nm=np-1 else cdayozm=date_oz(12) np=np1 nm=12 endif call getfactors(ozncyc,np1, cdayozm, cdayozp,intjulian, & fact1, fact2) do k=1,levsiz do i=1,pcols ozmix(i,k) = ozmixmj(i,k,nm)*fact1 + ozmixmj(i,k,np)*fact2 end do end do END subroutine oznint subroutine get_aerosol(c, julday, julian, dt, gmt, xtime, m_psp, m_psn, aerosoljp, & aerosoljn, m_hybi, paerlev, naer_c, pint, pcols, pver, pverp, pverr, pverrp, AEROSOLt, scale) integer, intent(in) :: c integer, intent(in) :: paerlev, naer_c, pcols, pver, pverp, pverr, pverrp real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: scale(naer_all) REAL, INTENT(IN ) :: XTIME,GMT INTEGER, INTENT(IN ) :: JULDAY REAL, INTENT(IN ) :: JULIAN REAL, INTENT(IN ) :: DT real(r8), intent(in ) :: m_psp(pcols),m_psn(pcols) real(r8), intent(in ) :: aerosoljp(pcols,paerlev,naer_c) real(r8), intent(in ) :: aerosoljn(pcols,paerlev,naer_c) real(r8), intent(in ) :: m_hybi(paerlev) real(r8), intent(out) :: AEROSOLt(pcols, pver, naer_all) real(r8) caldayloc real(r8) fact1, fact2 integer :: nm = 1 integer :: np = 2 integer :: mo_nxt = bigint integer :: mo_prv real(r8) :: cdaym = inf real(r8) :: cdayp = inf real(r8) :: Mid(12) data Mid/16.5, 46.0, 75.5, 106.0, 136.5, 167.0, 197.5, 228.5, 259.0, 289.5, 320.0, 350.5 / integer i, k, j integer m integer lats(pcols),lons(pcols) integer ncol INTEGER IJUL REAL(r8) intJULIAN real(r8) speciesmin(naer) real(r8) AEROSOLm(pcols,pver,naer) real(r8) AEROSOLp(pcols,pver,naer) CHARACTER(LEN=256) :: msgstr intJULIAN = JULIAN + 1.0_r8 IJUL=INT(intJULIAN) intJULIAN=intJULIAN-FLOAT(IJUL) IJUL=MOD(IJUL,365) IF(IJUL.EQ.0)IJUL=365 caldayloc=intJULIAN+IJUL if (caldayloc < Mid(1)) then mo_prv = 12 mo_nxt = 1 else if (caldayloc >= Mid(12)) then mo_prv = 12 mo_nxt = 1 else do i = 2 , 12 if (caldayloc < Mid(i)) then mo_prv = i-1 mo_nxt = i exit end if end do end if cdaym = Mid(mo_prv) cdayp = Mid(mo_nxt) call getfactors (.true., mo_nxt, cdaym, cdayp, caldayloc, & fact1, fact2) ncol = pcols call vert_interpolate (m_psp, aerosoljp, m_hybi, paerlev, naer_c, pint, nm, AEROSOLm, pcols, pver, pverp, ncol, c) call vert_interpolate (m_psn, aerosoljn, m_hybi, paerlev, naer_c, pint, np, AEROSOLp, pcols, pver, pverp, ncol, c) do m=1,naer do k=1,pver do i=1,ncol AEROSOLt(i,k,m) = AEROSOLm(i,k,m)*fact1 + AEROSOLp(i,k,m)*fact2 end do end do end do call background (c, ncol, pint, pcols, pverr, pverrp, AEROSOLt(:, :, idxBG)) AEROSOLt(:,:,idxVOLC) = 0._r8 speciesmin(:) = 0. do m=1,naer do k=1,pver do i=1,ncol if (AEROSOLt(i, k, m) < speciesmin(m)) then write(6,*) 'AEROSOL_INTERPOLATE: negative mass mixing ratio, exiting' write(6,*) 'm, column, pver',m, i, k ,AEROSOLt(i, k, m) call endrun () end if end do end do end do call scale_aerosols (AEROSOLt, pcols, pver, ncol, c, scale) return end subroutine get_aerosol subroutine aerosol_indirect(ncol,lchnk,pcols,pver,ppcnst,landfrac,pmid,t,qm1,cld,zm,rel) integer, intent(in) :: ncol integer, intent(in) :: lchnk integer, intent(in) :: pcols,pver,ppcnst real(r8), intent(in) :: landfrac(pcols) real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: t(pcols,pver) real(r8), intent(in) :: qm1(pcols,pver,ppcnst) real(r8), intent(in) :: cld(pcols,pver) real(r8), intent(in) :: zm(pcols,pver) real(r8), intent(in) :: rel(pcols,pver) real(r8) locrhoair(pcols,pver) real(r8) lwcwat(pcols,pver) real(r8) sulfmix(pcols,pver) real(r8) so4mass(pcols,pver) real(r8) Aso4(pcols,pver) real(r8) Ntot(pcols,pver) real(r8) relmod(pcols,pver) real(r8) wrel(pcols,pver) real(r8) wlwc(pcols,pver) real(r8) cldfrq(pcols,pver) real(r8) locPi real(r8) Rdryair real(r8) rhowat real(r8) Acoef real(r8) rekappa real(r8) recoef real(r8) reexp real(r8) Ntotb real(r8) Cmarn real(r8) Cland real(r8) Hmarn real(r8) Hland parameter ( Cmarn = 50.0, Cland = 100.0 ) parameter ( Hmarn = 1000.0, Hland = 2000.0 ) real(r8) bgaer integer i,k logical land land(i) = nint(landfrac(i)).gt.0.5_r8 if (indirect) then sulfmix(1:ncol,1:pver) = 0._r8 locPi = 3.141592654 Rdryair = 287.04 rhowat = 1000.0 Acoef = 1.2930E14 recoef = 3.0/(4.0*locPi*rhowat) reexp = 1.0/3.0 do k=pver,1,-1 do i = 1,ncol locrhoair(i,k) = pmid(i,k)/( Rdryair*t(i,k) ) lwcwat(i,k) = ( qm1(i,k,ixcldliq)/max(0.01_r8,cld(i,k)) )* & locrhoair(i,k) so4mass(i,k) = sulfmix(i,k)*locrhoair(i,k)*0.001 Aso4(i,k) = so4mass(i,k)*Acoef if (Aso4(i,k) <= 280.0) then Aso4(i,k) = max(36.0_r8,Aso4(i,k)) Ntot(i,k) = -1.15E-3*Aso4(i,k)**2 + 0.963*Aso4(i,k)+5.30 rekappa = 0.80 else Aso4(i,k) = min(1500.0_r8,Aso4(i,k)) Ntot(i,k) = -2.10E-4*Aso4(i,k)**2 + 0.568*Aso4(i,k)-27.9 rekappa = 0.67 end if if (land(i)) then bgaer = Cland*exp(-(zm(i,k)/Hland)) Ntot(i,k) = max(bgaer,Ntot(i,k)) else bgaer = Cmarn*exp(-(zm(i,k)/Hmarn)) Ntot(i,k) = max(bgaer,Ntot(i,k)) end if if (k == pver) then Ntotb = Ntot(i,k) else Ntotb = Ntot(i,k+1) end if relmod(i,k) = (( (recoef*lwcwat(i,k))/(rekappa*Ntotb))**reexp)*10000.0 relmod(i,k) = max(4.0_r8,relmod(i,k)) relmod(i,k) = min(20.0_r8,relmod(i,k)) if (cld(i,k) >= 0.01) then cldfrq(i,k) = 1.0 else cldfrq(i,k) = 0.0 end if wrel(i,k) = relmod(i,k)*cldfrq(i,k) wlwc(i,k) = lwcwat(i,k)*cldfrq(i,k) end do end do else do k = 1, pver do i = 1, ncol relmod(i,k) = rel(i,k) end do end do endif return end subroutine aerosol_indirect subroutine aer_trn(aer_mpp, aer_trn_ttl, pcols, plev, plevp ) implicit none integer, intent(in) :: pcols, plev, plevp real(r8), intent(in) :: aer_mpp(pcols,plevp) real(r8), intent(out) :: aer_trn_ttl(pcols,plevp,plevp,bnd_nbr_LW) integer bnd_idx integer i integer k1 integer k2 real(r8) aer_pth_dlt real(r8) odap_aer_ttl if (strat_volcanic) then do bnd_idx=1,bnd_nbr_LW do i=1,pcols aer_trn_ttl(i,1,1,bnd_idx)=1.0 end do do k1=2,plevp do i=1,pcols aer_trn_ttl(i,k1,k1,bnd_idx)=1.0 aer_pth_dlt = abs(aer_mpp(i,k1) - aer_mpp(i,1)) odap_aer_ttl = abs_cff_mss_aer(bnd_idx) * aer_pth_dlt aer_trn_ttl(i,1,k1,bnd_idx) = exp(-1.66 * odap_aer_ttl) end do end do do k1=2,plev do k2=k1+1,plevp do i=1,pcols aer_trn_ttl(i,k1,k2,bnd_idx) = & aer_trn_ttl(i,1,k2,bnd_idx) / & aer_trn_ttl(i,1,k1,bnd_idx) end do end do end do do k1=2,plevp do k2=1,k1-1 do i=1,pcols aer_trn_ttl(i,k1,k2,bnd_idx)=aer_trn_ttl(i,k2,k1,bnd_idx) end do end do end do end do else aer_trn_ttl = 1.0 endif return end subroutine aer_trn subroutine aer_pth(aer_mass, aer_mpp, ncol, pcols, plev, plevp) implicit none integer, intent(in) :: pcols, plev, plevp real(r8), intent(in):: aer_mass(pcols,plev) integer, intent(in):: ncol real(r8), intent(out):: aer_mpp(pcols,plevp) integer i integer k aer_mpp(1:ncol,1) = 0._r8 do k=2,plevp aer_mpp(1:ncol,k) = aer_mpp(1:ncol,k-1) + aer_mass(1:ncol,k-1) enddo return end subroutine aer_pth subroutine radctl(j, lchnk ,ncol , pcols, pver, pverp, pverr, pverrp, ppcnst, pcnst, & lwups ,emis , & pmid ,pint ,pmln ,piln ,pdel ,t , & qm1 ,cld ,cicewp ,cliqwp ,tauxcl, tauxci, coszrs, clat, & asdir ,asdif ,aldir ,aldif ,solcon, GMT,JULDAY,JULIAN,DT,XTIME, & pin, ozmixmj, ozmix, levsiz, num_months, & m_psp, m_psn, aerosoljp, aerosoljn, m_hybi, paerlev, naer_c, pmxrgn , & nmxrgn , & dolw, dosw, doabsems, abstot, absnxt, emstot, & fsup ,fsupc ,fsdn ,fsdnc , & fsdndir ,fsdncdir,fsdndif ,fsdncdif, & flup ,flupc ,fldn ,fldnc , & swcf ,lwcf ,flut , & fsns ,fsnt ,flns ,flnt , & qrs ,qrl ,flwds ,rel ,rei , & sols ,soll ,solsd ,solld , & landfrac,zm ,fsds, fsdsdir,fsdsdif ) implicit none integer, intent(in) :: lchnk,j integer, intent(in) :: ncol integer, intent(in) :: levsiz integer, intent(in) :: num_months integer, intent(in) :: paerlev,naer_c integer, intent(in) :: pcols, pver, pverp, pverr, pverrp, ppcnst, pcnst logical, intent(in) :: dolw,dosw,doabsems integer nspint integer naer_groups parameter ( nspint = 19 ) parameter ( naer_groups = 7 ) real(r8), intent(in) :: lwups(pcols) real(r8), intent(in) :: emis(pcols,pver) real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: pmln(pcols,pver) real(r8), intent(in) :: rel(pcols,pver) real(r8), intent(in) :: rei(pcols,pver) real(r8), intent(in) :: piln(pcols,pverp) real(r8), intent(in) :: pdel(pcols,pverp) real(r8), intent(in) :: t(pcols,pver) real(r8), intent(in) :: qm1(pcols,pver,ppcnst) real(r8), intent(in) :: cld(pcols,pver) real(r8), intent(in) :: cicewp(pcols,pver) real(r8), intent(in) :: cliqwp(pcols,pver) real(r8), intent(inout) :: tauxcl(pcols,0:pver) real(r8), intent(inout) :: tauxci(pcols,0:pver) real(r8), intent(in) :: coszrs(pcols) real(r8), intent(in) :: clat(pcols) real(r8), intent(in) :: asdir(pcols) real(r8), intent(in) :: asdif(pcols) real(r8), intent(in) :: aldir(pcols) real(r8), intent(in) :: aldif(pcols) real(r8), intent(in) :: landfrac(pcols) real(r8), intent(in) :: zm(pcols,pver) real(r8), intent(in) :: pin(levsiz) real(r8), intent(in) :: ozmixmj(pcols,levsiz,num_months) real(r8), intent(inout) :: ozmix(pcols,levsiz) real, intent(in) :: solcon REAL, INTENT(IN ) :: XTIME,GMT INTEGER, INTENT(IN ) :: JULDAY REAL, INTENT(IN ) :: JULIAN REAL, INTENT(IN ) :: DT real(r8), intent(in) :: m_psp(pcols),m_psn(pcols) real(r8), intent(in) :: aerosoljp(pcols,paerlev,naer_c) real(r8), intent(in) :: aerosoljn(pcols,paerlev,naer_c) real(r8), intent(in) :: m_hybi(paerlev) real(r8), intent(inout) :: pmxrgn(pcols,pverp) integer, intent(inout) :: nmxrgn(pcols) real(r8) :: pmxrgnrf(pcols,pverp) integer :: nmxrgnrf(pcols) real(r8), intent(out) :: fsns(pcols) real(r8), intent(out) :: fsnt(pcols) real(r8), intent(out) :: flns(pcols) real(r8), intent(out) :: flnt(pcols) real(r8), intent(out) :: sols(pcols) real(r8), intent(out) :: soll(pcols) real(r8), intent(out) :: solsd(pcols) real(r8), intent(out) :: solld(pcols) real(r8), intent(out) :: qrs(pcols,pver) real(r8), intent(out) :: fsds(pcols) real(r8), intent(out) :: fsdsdir(pcols) real(r8), intent(out) :: fsdsdif(pcols) real(r8), intent(out) :: fsup(pcols,pverp) real(r8), intent(out) :: fsupc(pcols,pverp) real(r8), intent(out) :: fsdn(pcols,pverp) real(r8), intent(out) :: fsdnc(pcols,pverp) real(r8), intent(out) :: fsdndir(pcols,pverp) real(r8), intent(out) :: fsdncdir(pcols,pverp) real(r8), intent(out) :: fsdndif(pcols,pverp) real(r8), intent(out) :: fsdncdif(pcols,pverp) real(r8), intent(out) :: flup(pcols,pverp) real(r8), intent(out) :: flupc(pcols,pverp) real(r8), intent(out) :: fldn(pcols,pverp) real(r8), intent(out) :: fldnc(pcols,pverp) real(r8), intent(out) :: swcf(pcols) real(r8), intent(out) :: lwcf(pcols) real(r8), intent(out) :: flut(pcols) real(r8), intent(out) :: qrl(pcols,pver) real(r8), intent(out) :: flwds(pcols) real(r8), intent(inout) :: abstot(pcols,pverp,pverp) real(r8), intent(inout) :: absnxt(pcols,pver,4) real(r8), intent(inout) :: emstot(pcols,pverp) integer i, k integer :: in2o, ich4, if11, if12 real(r8) solin(pcols) real(r8) fsntoa(pcols) real(r8) fsntoac(pcols) real(r8) fsnirt(pcols) real(r8) fsnrtc(pcols) real(r8) fsnirtsq(pcols) real(r8) fsntc(pcols) real(r8) fsnsc(pcols) real(r8) fsdsc(pcols) real(r8) fsdscdir(pcols) real(r8) fsdscdif(pcols) real(r8) flutc(pcols) real(r8) flntc(pcols) real(r8) flnsc(pcols) real(r8) ftem(pcols,pver) real(r8) pbr(pcols,pverr) real(r8) pnm(pcols,pverrp) real(r8) o3vmr(pcols,pverr) real(r8) o3mmr(pcols,pverr) real(r8) eccf real(r8) n2o(pcols,pver) real(r8) ch4(pcols,pver) real(r8) cfc11(pcols,pver) real(r8) cfc12(pcols,pver) real(r8) rh(pcols,pverr) real(r8) lwupcgs(pcols) real(r8) esat(pcols,pverr) real(r8) qsat(pcols,pverr) real(r8) :: frc_day(pcols) real(r8) :: aertau(pcols,nspint,naer_groups) real(r8) :: aerssa(pcols,nspint,naer_groups) real(r8) :: aerasm(pcols,nspint,naer_groups) real(r8) :: aerfwd(pcols,nspint,naer_groups) real(r8) aerosol(pcols, pver, naer_all) real(r8) scales(naer_all) LOGICAL, EXTERNAL :: wrf_dm_on_monitor CHARACTER (LEN=256) :: message call oznint(julday,julian,dt,gmt,xtime,ozmixmj,ozmix,levsiz,num_months,pcols) call radozn(lchnk ,ncol & ,pcols, pver & ,pmid ,pin, levsiz, ozmix, o3vmr ) call radinp(lchnk ,ncol ,pcols, pver, pverp, & pmid ,pint ,o3vmr , pbr ,& pnm ,eccf ,o3mmr ) if (dosw) then call aqsat(t, pmid, esat, qsat, pcols, & ncol, pver, 1, pver) rh(1:ncol,1:pver) = qm1(1:ncol,1:pver,1) / qsat(1:ncol,1:pver) * & ((1.0 - epsilo) * qsat(1:ncol,1:pver) + epsilo) / & ((1.0 - epsilo) * qm1(1:ncol,1:pver,1) + epsilo) if (radforce) then pmxrgnrf = pmxrgn nmxrgnrf = nmxrgn call get_rf_scales(scales) call get_aerosol(lchnk, julday, julian, dt, gmt, xtime, m_psp, m_psn, aerosoljp, & aerosoljn, m_hybi, paerlev, naer, pint, pcols, pver, pverp, pverr, pverrp, aerosol, scales) call aerosol_indirect(ncol,lchnk,pcols,pver,ppcnst,landfrac,pmid,t,qm1,cld,zm,rel) call radcswmx(j, lchnk ,ncol ,pcols, pver, pverp, & pnm ,pbr ,qm1 ,rh ,o3mmr , & aerosol ,cld ,cicewp ,cliqwp ,rel , & rei ,tauxcl ,tauxci ,eccf ,coszrs ,scon ,solin ,solcon , & asdir ,asdif ,aldir ,aldif ,nmxrgnrf, & pmxrgnrf,qrs ,fsnt ,fsntc ,fsntoa , & fsntoac ,fsnirt ,fsnrtc ,fsnirtsq,fsns , & fsnsc ,fsdsc ,fsds ,sols ,soll , & solsd ,solld ,frc_day , & fsup ,fsupc ,fsdn ,fsdnc , & fsdndir ,fsdncdir,fsdndif ,fsdncdif, & fsdsdir ,fsdsdif ,fsdscdir,fsdscdif, & aertau ,aerssa ,aerasm ,aerfwd ) do i = 1, ncol solin(i) = solin(i)*1.e-3 fsnt(i) = fsnt(i) *1.e-3 fsns(i) = fsns(i) *1.e-3 fsntc(i) = fsntc(i)*1.e-3 fsnsc(i) = fsnsc(i)*1.e-3 end do ftem(:ncol,:pver) = qrs(:ncol,:pver)/cpair endif call get_int_scales(scales) call get_aerosol(lchnk, julday, julian, dt, gmt, xtime, m_psp, m_psn, aerosoljp, aerosoljn, & m_hybi, paerlev, naer, pint, pcols, pver, pverp, pverr, pverrp, aerosol, scales) call aerosol_indirect(ncol,lchnk,pcols,pver,ppcnst,landfrac,pmid,t,qm1,cld,zm,rel) call radcswmx(j, lchnk ,ncol ,pcols, pver, pverp, & pnm ,pbr ,qm1 ,rh ,o3mmr , & aerosol ,cld ,cicewp ,cliqwp ,rel , & rei ,tauxcl ,tauxci ,eccf ,coszrs ,scon ,solin ,solcon , & asdir ,asdif ,aldir ,aldif ,nmxrgn , & pmxrgn ,qrs ,fsnt ,fsntc ,fsntoa , & fsntoac ,fsnirt ,fsnrtc ,fsnirtsq,fsns , & fsnsc ,fsdsc ,fsds ,sols ,soll , & solsd ,solld ,frc_day , & fsup ,fsupc ,fsdn ,fsdnc , & fsdndir ,fsdncdir,fsdndif ,fsdncdif, & fsdsdir ,fsdsdif ,fsdscdir,fsdscdif, & aertau ,aerssa ,aerasm ,aerfwd ) do i=1,ncol solin(i) = solin(i)*1.e-3 fsds(i) = fsds(i)*1.e-3 fsnirt(i)= fsnirt(i)*1.e-3 fsnrtc(i)= fsnrtc(i)*1.e-3 fsnirtsq(i)= fsnirtsq(i)*1.e-3 fsnt(i) = fsnt(i) *1.e-3 fsns(i) = fsns(i) *1.e-3 fsntc(i) = fsntc(i)*1.e-3 fsnsc(i) = fsnsc(i)*1.e-3 fsdsc(i) = fsdsc(i)*1.e-3 fsntoa(i)=fsntoa(i)*1.e-3 fsntoac(i)=fsntoac(i)*1.e-3 swcf(i) = fsntoa(i) - fsntoac(i) fsdsdir(i) = fsdsdir(i)*1.e-3 fsdsdif(i) = fsdsdif(i)*1.e-3 fsdscdir(i) = fsdscdir(i)*1.e-3 fsdscdif(i) = fsdscdif(i)*1.e-3 end do ftem(:ncol,:pver) = qrs(:ncol,:pver)/cpair do k = 1, pverp do i = 1, ncol fsup(i,k) = fsup(i,k)*1.e-3 fsupc(i,k) = fsupc(i,k)*1.e-3 fsdn(i,k) = fsdn(i,k)*1.e-3 fsdnc(i,k) = fsdnc(i,k)*1.e-3 end do end do end if if (dolw) then call get_int_scales(scales) call get_aerosol(lchnk, julday, julian, dt, gmt, xtime, m_psp, m_psn, aerosoljp, aerosoljn, & m_hybi, paerlev, naer, pint, pcols, pver, pverp, pverr, pverrp, aerosol, scales) do i=1,ncol lwupcgs(i) = lwups(i) end do if (trace_gas) then call radclwmx(lchnk ,ncol ,pcols, pver, pverp , & lwupcgs ,t ,qm1(1,1,1) ,o3vmr , & pbr ,pnm ,pmln ,piln , & qm1(1,1,in2o) ,qm1(1,1,ich4) , & qm1(1,1,if11) ,qm1(1,1,if12) , & cld ,emis ,pmxrgn ,nmxrgn ,qrl , & doabsems, abstot, absnxt, emstot, & flns ,flnt ,flnsc ,flntc ,flwds , & flut ,flutc , & flup ,flupc ,fldn ,fldnc , & aerosol(:,:,idxVOLC)) else call trcmix(lchnk ,ncol ,pcols, pver, & pmid ,clat, n2o ,ch4 , & cfc11 ,cfc12 ) IF ( wrf_dm_on_monitor() ) THEN WRITE(message,*)'CLWRF post_trcmix_values. n2o:', n2o(pcols/2,pver/2), ' ch4:', & ch4(pcols/2,pver/2),' cfc11:', cfc11(pcols/2,pver/2),' cfc12:', cfc12(pcols/2,pver/2) CALL wrf_debug(1, message) END IF call radclwmx(lchnk ,ncol ,pcols, pver, pverp , & lwupcgs ,t ,qm1(1,1,1) ,o3vmr , & pbr ,pnm ,pmln ,piln , & n2o ,ch4 ,cfc11 ,cfc12 , & cld ,emis ,pmxrgn ,nmxrgn ,qrl , & doabsems, abstot, absnxt, emstot, & flns ,flnt ,flnsc ,flntc ,flwds , & flut ,flutc , & flup ,flupc ,fldn ,fldnc , & aerosol(:,:,idxVOLC)) endif do i=1,ncol flnt(i) = flnt(i)*1.e-3 flut(i) = flut(i)*1.e-3 flutc(i) = flutc(i)*1.e-3 flns(i) = flns(i)*1.e-3 flntc(i) = flntc(i)*1.e-3 flnsc(i) = flnsc(i)*1.e-3 flwds(i) = flwds(i)*1.e-3 lwcf(i) = flutc(i) - flut(i) end do do k = 1, pverp do i = 1, ncol flup(i,k) = flup(i,k)*1.e-3 flupc(i,k) = flupc(i,k)*1.e-3 fldn(i,k) = fldn(i,k)*1.e-3 fldnc(i,k) = fldnc(i,k)*1.e-3 end do end do end if return end subroutine radctl subroutine param_cldoptics_calc(ncol, pcols, pver, pverp, pverr, pverrp, ppcnst, & q, cldn, landfrac, landm,icefrac, & pdel, t, ps, pmid, pint, cicewp, cliqwp, emis, rel, rei, pmxrgn, nmxrgn, snowh ) implicit none integer, intent(in) :: ncol, pcols, pver, pverp, pverr, pverrp, ppcnst real(r8), intent(in) :: q(pcols,pver,ppcnst) real(r8), intent(in) :: cldn(pcols,pver) real(r8), intent(in) :: pdel(pcols,pver) real(r8), intent(in) :: t(pcols,pver) real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: ps(pcols) real(r8), intent(in) :: landfrac(pcols) real(r8), intent(in) :: icefrac(pcols) real(r8), intent(in) :: landm(pcols) real(r8), intent(in) :: snowh(pcols) real(r8), intent(out) :: cicewp(pcols,pver) real(r8), intent(out) :: cliqwp(pcols,pver) real(r8), intent(out) :: emis (pcols,pver) real(r8), intent(out) :: rel (pcols,pver) real(r8), intent(out) :: rei (pcols,pver) real(r8), intent(out) :: pmxrgn(pcols,pver+1) integer , intent(out) :: nmxrgn(pcols) real(r8) :: cwp (pcols,pver) real(r8) :: effcld(pcols,pver) real(r8) :: gicewp(pcols,pver) real(r8) :: gliqwp(pcols,pver) real(r8) :: gwp (pcols,pver) real(r8) :: hl (pcols) real(r8) :: tgicewp(pcols) real(r8) :: tgliqwp(pcols) real(r8) :: tgwp (pcols) real(r8) :: tpw (pcols) real(r8) :: clwpold(pcols,pver) real(r8) :: ficemr (pcols,pver) real(r8) :: rgrav integer :: i,k integer :: lchnk tgicewp(:ncol) = 0. tgliqwp(:ncol) = 0. do k=1,pver do i = 1,ncol gicewp(i,k) = q(i,k,ixcldice)*pdel(i,k)/gravmks*1000.0 gliqwp(i,k) = q(i,k,ixcldliq)*pdel(i,k)/gravmks*1000.0 cicewp(i,k) = gicewp(i,k) / max(0.01_r8,cldn(i,k)) cliqwp(i,k) = gliqwp(i,k) / max(0.01_r8,cldn(i,k)) ficemr(i,k) = q(i,k,ixcldice) / & max(1.e-10_r8,(q(i,k,ixcldice)+q(i,k,ixcldliq))) tgicewp(i) = tgicewp(i) + gicewp(i,k) tgliqwp(i) = tgliqwp(i) + gliqwp(i,k) end do end do tgwp(:ncol) = tgicewp(:ncol) + tgliqwp(:ncol) gwp(:ncol,:pver) = gicewp(:ncol,:pver) + gliqwp(:ncol,:pver) cwp(:ncol,:pver) = cicewp(:ncol,:pver) + cliqwp(:ncol,:pver) tpw(:ncol) = 0.0 rgrav = 1.0/gravmks do k=1,pver do i=1,ncol tpw(i) = tpw(i) + pdel(i,k)*q(i,k,1)*rgrav end do end do call cldefr(lchnk, ncol, pcols, pver, pverp, landfrac, t, rel, rei, ps, pmid, landm, icefrac, snowh) call cldems(lchnk, ncol, pcols, pver, pverp, cwp, ficemr, rei, emis) do k=1,pver do i=1,ncol effcld(i,k) = cldn(i,k)*emis(i,k) end do end do call cldovrlap(lchnk, ncol, pcols, pver, pverp, pint, cldn, nmxrgn, pmxrgn) end subroutine param_cldoptics_calc subroutine radabs(lchnk ,ncol ,pcols, pver, pverp, & pbr ,pnm ,co2em ,co2eml ,tplnka , & s2c ,tcg ,w ,h2otr ,plco2 , & plh2o ,co2t ,tint ,tlayr ,plol , & plos ,pmln ,piln ,ucfc11 ,ucfc12 , & un2o0 ,un2o1 ,uch4 ,uco211 ,uco212 , & uco213 ,uco221 ,uco222 ,uco223 ,uptype , & bn2o0 ,bn2o1 ,bch4 ,abplnk1 ,abplnk2 , & abstot ,absnxt ,plh2ob ,wb , & aer_mpp ,aer_trn_ttl) integer, intent(in) :: lchnk integer, intent(in) :: ncol integer, intent(in) :: pcols, pver, pverp real(r8), intent(in) :: pbr(pcols,pver) real(r8), intent(in) :: pnm(pcols,pverp) real(r8), intent(in) :: co2em(pcols,pverp) real(r8), intent(in) :: co2eml(pcols,pver) real(r8), intent(in) :: tplnka(pcols,pverp) real(r8), intent(in) :: s2c(pcols,pverp) real(r8), intent(in) :: tcg(pcols,pverp) real(r8), intent(in) :: w(pcols,pverp) real(r8), intent(in) :: h2otr(pcols,pverp) real(r8), intent(in) :: plco2(pcols,pverp) real(r8), intent(in) :: plh2o(pcols,pverp) real(r8), intent(in) :: co2t(pcols,pverp) real(r8), intent(in) :: tint(pcols,pverp) real(r8), intent(in) :: tlayr(pcols,pverp) real(r8), intent(in) :: plol(pcols,pverp) real(r8), intent(in) :: plos(pcols,pverp) real(r8), intent(in) :: pmln(pcols,pver) real(r8), intent(in) :: piln(pcols,pverp) real(r8), intent(in) :: plh2ob(nbands,pcols,pverp) real(r8), intent(in) :: wb(nbands,pcols,pverp) real(r8), intent(in) :: aer_mpp(pcols,pverp) real(r8), intent(in) :: aer_trn_ttl(pcols,pverp,pverp,bnd_nbr_LW) real(r8), intent(in) :: ucfc11(pcols,pverp) real(r8), intent(in) :: ucfc12(pcols,pverp) real(r8), intent(in) :: un2o0(pcols,pverp) real(r8), intent(in) :: un2o1(pcols,pverp) real(r8), intent(in) :: uch4(pcols,pverp) real(r8), intent(in) :: uco211(pcols,pverp) real(r8), intent(in) :: uco212(pcols,pverp) real(r8), intent(in) :: uco213(pcols,pverp) real(r8), intent(in) :: uco221(pcols,pverp) real(r8), intent(in) :: uco222(pcols,pverp) real(r8), intent(in) :: uco223(pcols,pverp) real(r8), intent(in) :: uptype(pcols,pverp) real(r8), intent(in) :: bn2o0(pcols,pverp) real(r8), intent(in) :: bn2o1(pcols,pverp) real(r8), intent(in) :: bch4(pcols,pverp) real(r8), intent(in) :: abplnk1(14,pcols,pverp) real(r8), intent(in) :: abplnk2(14,pcols,pverp) real(r8), intent(out) :: abstot(pcols,pverp,pverp) real(r8), intent(out) :: absnxt(pcols,pver,4) integer i integer k integer k1 integer k2 integer kn integer wvl real(r8) abstrc(pcols) real(r8) bplnk(14,pcols,4) real(r8) pnew(pcols) real(r8) pnewb(nbands) real(r8) u(pcols) real(r8) ub(nbands) real(r8) tbar(pcols,4) real(r8) emm(pcols,4) real(r8) o3emm(pcols,4) real(r8) o3bndi real(r8) temh2o(pcols,4) real(r8) k21 real(r8) k22 real(r8) uc1(pcols) real(r8) to3h2o(pcols) real(r8) pi real(r8) sqti(pcols) real(r8) et real(r8) et2 real(r8) et4 real(r8) omet real(r8) f1co2 real(r8) f2co2(pcols) real(r8) f3co2(pcols) real(r8) t1co2(pcols) real(r8) sqwp real(r8) f1sqwp(pcols) real(r8) oneme real(r8) alphat real(r8) wco2 real(r8) posqt real(r8) u7(pcols) real(r8) u8 real(r8) u9 real(r8) u13 real(r8) rbeta7(pcols) real(r8) rbeta8 real(r8) rbeta9 real(r8) rbeta13 real(r8) tpatha real(r8) abso(pcols,4) real(r8) dtx(pcols) real(r8) dty(pcols) real(r8) term7(pcols,2) real(r8) term8(pcols,2) real(r8) tr1 real(r8) tr10(pcols) real(r8) tr2 real(r8) tr5 real(r8) tr6 real(r8) tr9(pcols) real(r8) sqrtu(pcols) real(r8) fwk(pcols) real(r8) fwku(pcols) real(r8) to3co2(pcols) real(r8) dpnm(pcols) real(r8) pnmsq(pcols,pverp) real(r8) dw(pcols) real(r8) uinpl(pcols,4) real(r8) winpl(pcols,4) real(r8) zinpl(pcols,4) real(r8) pinpl(pcols,4) real(r8) dplh2o(pcols) real(r8) r293 real(r8) r250 real(r8) r3205 real(r8) r300 real(r8) rsslp real(r8) r2sslp real(r8) ds2c real(r8) dplos real(r8) dplol real(r8) tlocal real(r8) beta real(r8) rphat real(r8) tcrfac real(r8) tmp1 real(r8) u1 real(r8) realnu real(r8) tmp2 real(r8) u2 real(r8) rsqti real(r8) tpath real(r8) tmp3 real(r8) rdpnmsq real(r8) rdpnm real(r8) p1 real(r8) p2 real(r8) dtym10 real(r8) dplco2 real(r8) te real(r8) denom real(r8) th2o(pcols) real(r8) tco2(pcols) real(r8) to3(pcols) real(r8) trab2(pcols) real(r8) absbnd real(r8) dbvtit(pcols,pverp) real(r8) dbvtly(pcols,pver) real(r8) fa real(r8) a_star real(r8) l_star real(r8) c_star real(r8) te1 real(r8) te2 real(r8) te3 real(r8) te4 real(r8) te5 real(r8) log_u real(r8) log_uc real(r8) log_p real(r8) t_p real(r8) t_e integer iu integer iu1 integer iuc integer iuc1 integer ip integer ip1 integer itp integer itp1 integer ite integer ite1 integer irh integer irh1 real(r8) dvar real(r8) uvar real(r8) uscl real(r8) wu real(r8) wu1 real(r8) wuc real(r8) wuc1 real(r8) wp real(r8) wp1 real(r8) wtp real(r8) wtp1 real(r8) wte real(r8) wte1 real(r8) wrh real(r8) wrh1 real(r8) w_0_0_ real(r8) w_0_1_ real(r8) w_1_0_ real(r8) w_1_1_ real(r8) w_0_00 real(r8) w_0_01 real(r8) w_0_10 real(r8) w_0_11 real(r8) w_1_00 real(r8) w_1_01 real(r8) w_1_10 real(r8) w_1_11 real(r8) w00_00 real(r8) w00_01 real(r8) w00_10 real(r8) w00_11 real(r8) w01_00 real(r8) w01_01 real(r8) w01_10 real(r8) w01_11 real(r8) w10_00 real(r8) w10_01 real(r8) w10_10 real(r8) w10_11 real(r8) w11_00 real(r8) w11_01 real(r8) w11_10 real(r8) w11_11 integer ib real(r8) pch2o real(r8) fch2o real(r8) uch2o real(r8) fdif real(r8) sslp_mks real(r8) esx real(r8) qsx real(r8) pnew_mks real(r8) q_path real(r8) rh_path real(r8) omeps integer iest integer bnd_idx real(r8) aer_pth_dlt real(r8) aer_pth_ngh(pcols) real(r8) odap_aer_ttl real(r8) aer_trn_ngh(pcols,bnd_nbr_LW) real(r8) dbvt,t dbvt(t)=(-2.8911366682e-4+(2.3771251896e-6+1.1305188929e-10*t)*t)/ & (1.0+(-6.1364820707e-3+1.5550319767e-5*t)*t) do k2=1,ntoplw-1 do k1=1,ntoplw-1 abstot(:,k1,k2) = inf end do end do do k2=1,4 do k1=1,ntoplw-1 absnxt(:,k1,k2) = inf end do end do do k=ntoplw,pverp abstot(:,k,k) = inf end do do k=ntoplw,pver do i=1,ncol dbvtly(i,k) = dbvt(tlayr(i,k+1)) dbvtit(i,k) = dbvt(tint(i,k)) end do end do do i=1,ncol dbvtit(i,pverp) = dbvt(tint(i,pverp)) end do r293 = 1./293. r250 = 1./250. r3205 = 1./.3205 r300 = 1./300. rsslp = 1./sslp r2sslp = 1./(2.*sslp) fdif = 1.66 sslp_mks = sslp / 10.0 omeps = 1.0 - epsilo do k=ntoplw,pverp do i=1,ncol pnmsq(i,k) = pnm(i,k)**2 dtx(i) = tplnka(i,k) - 250. end do end do do k1=pverp,ntoplw,-1 do k2=pverp,ntoplw,-1 if (k1 == k2) cycle do i=1,ncol dplh2o(i) = plh2o(i,k1) - plh2o(i,k2) u(i) = abs(dplh2o(i)) sqrtu(i) = sqrt(u(i)) ds2c = abs(s2c(i,k1) - s2c(i,k2)) dw(i) = abs(w(i,k1) - w(i,k2)) uc1(i) = (ds2c + 1.7e-3*u(i))*(1. + 2.*ds2c)/(1. + 15.*ds2c) pch2o = ds2c pnew(i) = u(i)/dw(i) pnew_mks = pnew(i) * sslp_mks tpatha = abs(tcg(i,k1) - tcg(i,k2))/dw(i) t_p = min(max(tpatha, min_tp_h2o), max_tp_h2o) iest = floor(t_p) - min_tp_h2o esx = estblh2o(iest) + (estblh2o(iest+1)-estblh2o(iest)) * & (t_p - min_tp_h2o - iest) qsx = epsilo * esx / (pnew_mks - omeps * esx) q_path = dw(i) / abs(pnm(i,k1) - pnm(i,k2)) / rga ub(1) = abs(plh2ob(1,i,k1) - plh2ob(1,i,k2)) / psi(t_p,1) ub(2) = abs(plh2ob(2,i,k1) - plh2ob(2,i,k2)) / psi(t_p,2) pnewb(1) = ub(1) / abs(wb(1,i,k1) - wb(1,i,k2)) * phi(t_p,1) pnewb(2) = ub(2) / abs(wb(2,i,k1) - wb(2,i,k2)) * phi(t_p,2) dtx(i) = tplnka(i,k2) - 250. dty(i) = tpatha - 250. fwk(i) = fwcoef + fwc1/(1. + fwc2*u(i)) fwku(i) = fwk(i)*u(i) te1 = tplnka(i,k2) te2 = te1 * te1 te3 = te2 * te1 te4 = te3 * te1 te5 = te4 * te1 dvar = (t_p - min_tp_h2o) / dtp_h2o itp = min(max(int(aint(dvar,r8)) + 1, 1), n_tp - 1) itp1 = itp + 1 wtp = dvar - floor(dvar) wtp1 = 1.0 - wtp t_e = min(max(tplnka(i,k2)-t_p, min_te_h2o), max_te_h2o) dvar = (t_e - min_te_h2o) / dte_h2o ite = min(max(int(aint(dvar,r8)) + 1, 1), n_te - 1) ite1 = ite + 1 wte = dvar - floor(dvar) wte1 = 1.0 - wte rh_path = min(max(q_path / qsx, min_rh_h2o), max_rh_h2o) dvar = (rh_path - min_rh_h2o) / drh_h2o irh = min(max(int(aint(dvar,r8)) + 1, 1), n_rh - 1) irh1 = irh + 1 wrh = dvar - floor(dvar) wrh1 = 1.0 - wrh w_0_0_ = wtp * wte w_0_1_ = wtp * wte1 w_1_0_ = wtp1 * wte w_1_1_ = wtp1 * wte1 w_0_00 = w_0_0_ * wrh w_0_01 = w_0_0_ * wrh1 w_0_10 = w_0_1_ * wrh w_0_11 = w_0_1_ * wrh1 w_1_00 = w_1_0_ * wrh w_1_01 = w_1_0_ * wrh1 w_1_10 = w_1_1_ * wrh w_1_11 = w_1_1_ * wrh1 fch2o = fh2oself(t_p) uch2o = (pch2o * epsilo) / (q_path * fch2o) ib = 1 uvar = ub(ib) * fdif log_u = min(log10(max(uvar, min_u_h2o)), max_lu_h2o) dvar = (log_u - min_lu_h2o) / dlu_h2o iu = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iu1 = iu + 1 wu = dvar - floor(dvar) wu1 = 1.0 - wu log_p = min(log10(max(pnewb(ib), min_p_h2o)), max_lp_h2o) dvar = (log_p - min_lp_h2o) / dlp_h2o ip = min(max(int(aint(dvar,r8)) + 1, 1), n_p - 1) ip1 = ip + 1 wp = dvar - floor(dvar) wp1 = 1.0 - wp w00_00 = wp * w_0_00 w00_01 = wp * w_0_01 w00_10 = wp * w_0_10 w00_11 = wp * w_0_11 w01_00 = wp * w_1_00 w01_01 = wp * w_1_01 w01_10 = wp * w_1_10 w01_11 = wp * w_1_11 w10_00 = wp1 * w_0_00 w10_01 = wp1 * w_0_01 w10_10 = wp1 * w_0_10 w10_11 = wp1 * w_0_11 w11_00 = wp1 * w_1_00 w11_01 = wp1 * w_1_01 w11_10 = wp1 * w_1_10 w11_11 = wp1 * w_1_11 fa = fat(1,ib) + & fat(2,ib) * te1 + & fat(3,ib) * te2 + & fat(4,ib) * te3 + & fat(5,ib) * te4 + & fat(6,ib) * te5 a_star = & ah2onw(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & ah2onw(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & ah2onw(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & ah2onw(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & ah2onw(ip , itp , iu1, ite , irh ) * w11_11 * wu + & ah2onw(ip , itp , iu1, ite , irh1) * w11_10 * wu + & ah2onw(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & ah2onw(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & ah2onw(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & ah2onw(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & ah2onw(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & ah2onw(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & ah2onw(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & ah2onw(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & ah2onw(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & ah2onw(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & ah2onw(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & ah2onw(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & ah2onw(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & ah2onw(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & ah2onw(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & ah2onw(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & ah2onw(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & ah2onw(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & ah2onw(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & ah2onw(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & ah2onw(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & ah2onw(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & ah2onw(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & ah2onw(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & ah2onw(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & ah2onw(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu abso(i,ib) = min(max(fa * (1.0 - (1.0 - a_star) * & aer_trn_ttl(i,k1,k2,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o abso(i,ib) = abso(i,ib) * uscl endif ib = 2 uvar = ub(ib) * fdif log_u = min(log10(max(uvar, min_u_h2o)), max_lu_h2o) dvar = (log_u - min_lu_h2o) / dlu_h2o iu = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iu1 = iu + 1 wu = dvar - floor(dvar) wu1 = 1.0 - wu log_p = min(log10(max(pnewb(ib), min_p_h2o)), max_lp_h2o) dvar = (log_p - min_lp_h2o) / dlp_h2o ip = min(max(int(aint(dvar,r8)) + 1, 1), n_p - 1) ip1 = ip + 1 wp = dvar - floor(dvar) wp1 = 1.0 - wp w00_00 = wp * w_0_00 w00_01 = wp * w_0_01 w00_10 = wp * w_0_10 w00_11 = wp * w_0_11 w01_00 = wp * w_1_00 w01_01 = wp * w_1_01 w01_10 = wp * w_1_10 w01_11 = wp * w_1_11 w10_00 = wp1 * w_0_00 w10_01 = wp1 * w_0_01 w10_10 = wp1 * w_0_10 w10_11 = wp1 * w_0_11 w11_00 = wp1 * w_1_00 w11_01 = wp1 * w_1_01 w11_10 = wp1 * w_1_10 w11_11 = wp1 * w_1_11 log_uc = min(log10(max(uch2o * fdif, min_u_h2o)), max_lu_h2o) dvar = (log_uc - min_lu_h2o) / dlu_h2o iuc = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iuc1 = iuc + 1 wuc = dvar - floor(dvar) wuc1 = 1.0 - wuc fa = fat(1,ib) + & fat(2,ib) * te1 + & fat(3,ib) * te2 + & fat(4,ib) * te3 + & fat(5,ib) * te4 + & fat(6,ib) * te5 l_star = & ln_ah2ow(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & ln_ah2ow(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & ln_ah2ow(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & ln_ah2ow(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & ln_ah2ow(ip , itp , iu1, ite , irh ) * w11_11 * wu + & ln_ah2ow(ip , itp , iu1, ite , irh1) * w11_10 * wu + & ln_ah2ow(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & ln_ah2ow(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & ln_ah2ow(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & ln_ah2ow(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & ln_ah2ow(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & ln_ah2ow(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & ln_ah2ow(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & ln_ah2ow(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & ln_ah2ow(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & ln_ah2ow(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & ln_ah2ow(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & ln_ah2ow(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & ln_ah2ow(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & ln_ah2ow(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & ln_ah2ow(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & ln_ah2ow(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & ln_ah2ow(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & ln_ah2ow(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & ln_ah2ow(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & ln_ah2ow(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & ln_ah2ow(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & ln_ah2ow(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & ln_ah2ow(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & ln_ah2ow(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & ln_ah2ow(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & ln_ah2ow(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu c_star = & cn_ah2ow(ip , itp , iuc , ite , irh ) * w11_11 * wuc1 + & cn_ah2ow(ip , itp , iuc , ite , irh1) * w11_10 * wuc1 + & cn_ah2ow(ip , itp , iuc , ite1, irh ) * w11_01 * wuc1 + & cn_ah2ow(ip , itp , iuc , ite1, irh1) * w11_00 * wuc1 + & cn_ah2ow(ip , itp , iuc1, ite , irh ) * w11_11 * wuc + & cn_ah2ow(ip , itp , iuc1, ite , irh1) * w11_10 * wuc + & cn_ah2ow(ip , itp , iuc1, ite1, irh ) * w11_01 * wuc + & cn_ah2ow(ip , itp , iuc1, ite1, irh1) * w11_00 * wuc + & cn_ah2ow(ip , itp1, iuc , ite , irh ) * w10_11 * wuc1 + & cn_ah2ow(ip , itp1, iuc , ite , irh1) * w10_10 * wuc1 + & cn_ah2ow(ip , itp1, iuc , ite1, irh ) * w10_01 * wuc1 + & cn_ah2ow(ip , itp1, iuc , ite1, irh1) * w10_00 * wuc1 + & cn_ah2ow(ip , itp1, iuc1, ite , irh ) * w10_11 * wuc + & cn_ah2ow(ip , itp1, iuc1, ite , irh1) * w10_10 * wuc + & cn_ah2ow(ip , itp1, iuc1, ite1, irh ) * w10_01 * wuc + & cn_ah2ow(ip , itp1, iuc1, ite1, irh1) * w10_00 * wuc + & cn_ah2ow(ip1, itp , iuc , ite , irh ) * w01_11 * wuc1 + & cn_ah2ow(ip1, itp , iuc , ite , irh1) * w01_10 * wuc1 + & cn_ah2ow(ip1, itp , iuc , ite1, irh ) * w01_01 * wuc1 + & cn_ah2ow(ip1, itp , iuc , ite1, irh1) * w01_00 * wuc1 + & cn_ah2ow(ip1, itp , iuc1, ite , irh ) * w01_11 * wuc + & cn_ah2ow(ip1, itp , iuc1, ite , irh1) * w01_10 * wuc + & cn_ah2ow(ip1, itp , iuc1, ite1, irh ) * w01_01 * wuc + & cn_ah2ow(ip1, itp , iuc1, ite1, irh1) * w01_00 * wuc + & cn_ah2ow(ip1, itp1, iuc , ite , irh ) * w00_11 * wuc1 + & cn_ah2ow(ip1, itp1, iuc , ite , irh1) * w00_10 * wuc1 + & cn_ah2ow(ip1, itp1, iuc , ite1, irh ) * w00_01 * wuc1 + & cn_ah2ow(ip1, itp1, iuc , ite1, irh1) * w00_00 * wuc1 + & cn_ah2ow(ip1, itp1, iuc1, ite , irh ) * w00_11 * wuc + & cn_ah2ow(ip1, itp1, iuc1, ite , irh1) * w00_10 * wuc + & cn_ah2ow(ip1, itp1, iuc1, ite1, irh ) * w00_01 * wuc + & cn_ah2ow(ip1, itp1, iuc1, ite1, irh1) * w00_00 * wuc abso(i,ib) = min(max(fa * (1.0 - l_star * c_star * & aer_trn_ttl(i,k1,k2,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o abso(i,ib) = abso(i,ib) * uscl endif end do do i=1,ncol term7(i,1) = coefj(1,1) + coefj(2,1)*dty(i)*(1. + c16*dty(i)) term8(i,1) = coefk(1,1) + coefk(2,1)*dty(i)*(1. + c17*dty(i)) term7(i,2) = coefj(1,2) + coefj(2,2)*dty(i)*(1. + c26*dty(i)) term8(i,2) = coefk(1,2) + coefk(2,2)*dty(i)*(1. + c27*dty(i)) end do do i=1,ncol k21 = term7(i,1) + term8(i,1)/ & (1. + (c30 + c31*(dty(i)-10.)*(dty(i)-10.))*sqrtu(i)) k22 = term7(i,2) + term8(i,2)/ & (1. + (c28 + c29*(dty(i)-10.))*sqrtu(i)) tr1 = exp(-(k21*(sqrtu(i) + fc1*fwku(i)))) tr2 = exp(-(k22*(sqrtu(i) + fc1*fwku(i)))) tr1=tr1*aer_trn_ttl(i,k1,k2,idx_LW_0650_0800) tr2=tr2*aer_trn_ttl(i,k1,k2,idx_LW_0500_0650) tr5 = exp(-((coefh(1,3) + coefh(2,3)*dtx(i))*uc1(i))) tr6 = exp(-((coefh(1,4) + coefh(2,4)*dtx(i))*uc1(i))) tr9(i) = tr1*tr5 tr10(i) = tr2*tr6 th2o(i) = tr10(i) trab2(i) = 0.65*tr9(i) + 0.35*tr10(i) end do if (k2 < k1) then do i=1,ncol to3h2o(i) = h2otr(i,k1)/h2otr(i,k2) end do else do i=1,ncol to3h2o(i) = h2otr(i,k2)/h2otr(i,k1) end do end if do i=1,ncol dpnm(i) = pnm(i,k1) - pnm(i,k2) to3co2(i) = (pnm(i,k1)*co2t(i,k1) - pnm(i,k2)*co2t(i,k2))/dpnm(i) te = (to3co2(i)*r293)**.7 dplos = plos(i,k1) - plos(i,k2) dplol = plol(i,k1) - plol(i,k2) u1 = 18.29*abs(dplos)/te u2 = .5649*abs(dplos)/te rphat = dplol/dplos tlocal = tint(i,k2) tcrfac = sqrt(tlocal*r250)*te beta = r3205*(rphat + dpfo3*tcrfac) realnu = te/beta tmp1 = u1/sqrt(4. + u1*(1. + realnu)) tmp2 = u2/sqrt(4. + u2*(1. + realnu)) o3bndi = 74.*te*log(1. + tmp1 + tmp2) abso(i,3) = o3bndi*to3h2o(i)*dbvtit(i,k2) to3(i) = 1.0/(1. + 0.1*tmp1 + 0.1*tmp2) end do do i=1,ncol sqwp = sqrt(abs(plco2(i,k1) - plco2(i,k2))) et = exp(-480./to3co2(i)) sqti(i) = sqrt(to3co2(i)) rsqti = 1./sqti(i) et2 = et*et et4 = et2*et2 omet = 1. - 1.5*et2 f1co2 = 899.70*omet*(1. + 1.94774*et + 4.73486*et2)*rsqti f1sqwp(i) = f1co2*sqwp t1co2(i) = 1./(1. + (245.18*omet*sqwp*rsqti)) oneme = 1. - et2 alphat = oneme**3*rsqti pi = abs(dpnm(i)) wco2 = 2.5221*co2vmr*pi*rga u7(i) = 4.9411e4*alphat*et2*wco2 u8 = 3.9744e4*alphat*et4*wco2 u9 = 1.0447e5*alphat*et4*et2*wco2 u13 = 2.8388e3*alphat*et4*wco2 tpath = to3co2(i) tlocal = tint(i,k2) tcrfac = sqrt(tlocal*r250*tpath*r300) posqt = ((pnm(i,k2) + pnm(i,k1))*r2sslp + dpfco2*tcrfac)*rsqti rbeta7(i) = 1./(5.3228*posqt) rbeta8 = 1./(10.6576*posqt) rbeta9 = rbeta7(i) rbeta13 = rbeta9 f2co2(i) = (u7(i)/sqrt(4. + u7(i)*(1. + rbeta7(i)))) + & (u8 /sqrt(4. + u8*(1. + rbeta8))) + & (u9 /sqrt(4. + u9*(1. + rbeta9))) f3co2(i) = u13/sqrt(4. + u13*(1. + rbeta13)) end do if (k2 >= k1) then do i=1,ncol sqti(i) = sqrt(tlayr(i,k2)) end do end if do i=1,ncol tmp1 = log(1. + f1sqwp(i)) tmp2 = log(1. + f2co2(i)) tmp3 = log(1. + f3co2(i)) absbnd = (tmp1 + 2.*t1co2(i)*tmp2 + 2.*tmp3)*sqti(i) abso(i,4) = trab2(i)*co2em(i,k2)*absbnd tco2(i) = 1./(1.0+10.0*(u7(i)/sqrt(4. + u7(i)*(1. + rbeta7(i))))) end do call trcab( lchnk ,ncol ,pcols, pverp, & k1 ,k2 ,ucfc11 ,ucfc12 ,un2o0 , & un2o1 ,uch4 ,uco211 ,uco212 ,uco213 , & uco221 ,uco222 ,uco223 ,bn2o0 ,bn2o1 , & bch4 ,to3co2 ,pnm ,dw ,pnew , & s2c ,uptype ,u ,abplnk1 ,tco2 , & th2o ,to3 ,abstrc , & aer_trn_ttl) do i=1,ncol abstot(i,k1,k2) = abso(i,1) + abso(i,2) + & abso(i,3) + abso(i,4) + abstrc(i) end do end do end do do k2=pver,ntoplw,-1 do i=1,ncol tbar(i,1) = 0.5*(tint(i,k2+1) + tlayr(i,k2+1)) emm(i,1) = 0.5*(co2em(i,k2+1) + co2eml(i,k2)) tbar(i,2) = 0.5*(tlayr(i,k2+1) + tint(i,k2)) emm(i,2) = 0.5*(co2em(i,k2) + co2eml(i,k2)) tbar(i,3) = 0.5*(tbar(i,2) + tbar(i,1)) emm(i,3) = emm(i,1) tbar(i,4) = tbar(i,3) emm(i,4) = emm(i,2) o3emm(i,1) = 0.5*(dbvtit(i,k2+1) + dbvtly(i,k2)) o3emm(i,2) = 0.5*(dbvtit(i,k2) + dbvtly(i,k2)) o3emm(i,3) = o3emm(i,1) o3emm(i,4) = o3emm(i,2) temh2o(i,1) = tbar(i,1) temh2o(i,2) = tbar(i,2) temh2o(i,3) = tbar(i,1) temh2o(i,4) = tbar(i,2) dpnm(i) = pnm(i,k2+1) - pnm(i,k2) end do do wvl = 1,14 do i = 1,ncol bplnk(wvl,i,1) = 0.5*(abplnk1(wvl,i,k2+1) + abplnk2(wvl,i,k2)) bplnk(wvl,i,2) = 0.5*(abplnk1(wvl,i,k2) + abplnk2(wvl,i,k2)) bplnk(wvl,i,3) = bplnk(wvl,i,1) bplnk(wvl,i,4) = bplnk(wvl,i,2) end do end do do i=1,ncol rdpnmsq = 1./(pnmsq(i,k2+1) - pnmsq(i,k2)) rdpnm = 1./dpnm(i) p1 = .5*(pbr(i,k2) + pnm(i,k2+1)) p2 = .5*(pbr(i,k2) + pnm(i,k2 )) uinpl(i,1) = (pnmsq(i,k2+1) - p1**2)*rdpnmsq uinpl(i,2) = -(pnmsq(i,k2 ) - p2**2)*rdpnmsq uinpl(i,3) = -(pnmsq(i,k2 ) - p1**2)*rdpnmsq uinpl(i,4) = (pnmsq(i,k2+1) - p2**2)*rdpnmsq winpl(i,1) = (.5*( pnm(i,k2+1) - pbr(i,k2)))*rdpnm winpl(i,2) = (.5*(-pnm(i,k2 ) + pbr(i,k2)))*rdpnm winpl(i,3) = (.5*( pnm(i,k2+1) + pbr(i,k2)) - pnm(i,k2 ))*rdpnm winpl(i,4) = (.5*(-pnm(i,k2 ) - pbr(i,k2)) + pnm(i,k2+1))*rdpnm tmp1 = 1./(piln(i,k2+1) - piln(i,k2)) tmp2 = piln(i,k2+1) - pmln(i,k2) tmp3 = piln(i,k2 ) - pmln(i,k2) zinpl(i,1) = (.5*tmp2 )*tmp1 zinpl(i,2) = ( - .5*tmp3)*tmp1 zinpl(i,3) = (.5*tmp2 - tmp3)*tmp1 zinpl(i,4) = ( tmp2 - .5*tmp3)*tmp1 pinpl(i,1) = 0.5*(p1 + pnm(i,k2+1)) pinpl(i,2) = 0.5*(p2 + pnm(i,k2 )) pinpl(i,3) = 0.5*(p1 + pnm(i,k2 )) pinpl(i,4) = 0.5*(p2 + pnm(i,k2+1)) if(strat_volcanic) then aer_pth_ngh(i) = abs(aer_mpp(i,k2)-aer_mpp(i,k2+1)) endif end do do kn=1,4 do i=1,ncol u(i) = uinpl(i,kn)*abs(plh2o(i,k2) - plh2o(i,k2+1)) sqrtu(i) = sqrt(u(i)) dw(i) = abs(w(i,k2) - w(i,k2+1)) pnew(i) = u(i)/(winpl(i,kn)*dw(i)) pnew_mks = pnew(i) * sslp_mks t_p = min(max(tbar(i,kn), min_tp_h2o), max_tp_h2o) iest = floor(t_p) - min_tp_h2o esx = estblh2o(iest) + (estblh2o(iest+1)-estblh2o(iest)) * & (t_p - min_tp_h2o - iest) qsx = epsilo * esx / (pnew_mks - omeps * esx) q_path = dw(i) / ABS(dpnm(i)) / rga ds2c = abs(s2c(i,k2) - s2c(i,k2+1)) uc1(i) = uinpl(i,kn)*ds2c pch2o = uc1(i) uc1(i) = (uc1(i) + 1.7e-3*u(i))*(1. + 2.*uc1(i))/(1. + 15.*uc1(i)) dtx(i) = temh2o(i,kn) - 250. dty(i) = tbar(i,kn) - 250. fwk(i) = fwcoef + fwc1/(1. + fwc2*u(i)) fwku(i) = fwk(i)*u(i) if(strat_volcanic) then aer_pth_dlt=uinpl(i,kn)*aer_pth_ngh(i) do bnd_idx=1,bnd_nbr_LW odap_aer_ttl=abs_cff_mss_aer(bnd_idx) * aer_pth_dlt aer_trn_ngh(i,bnd_idx)=exp(-fdif * odap_aer_ttl) end do else aer_trn_ngh(i,:) = 1.0 endif te1 = temh2o(i,kn) te2 = te1 * te1 te3 = te2 * te1 te4 = te3 * te1 te5 = te4 * te1 uvar = u(i)*fdif log_u = min(log10(max(uvar, min_u_h2o)), max_lu_h2o) dvar = (log_u - min_lu_h2o) / dlu_h2o iu = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iu1 = iu + 1 wu = dvar - floor(dvar) wu1 = 1.0 - wu log_p = min(log10(max(pnew(i), min_p_h2o)), max_lp_h2o) dvar = (log_p - min_lp_h2o) / dlp_h2o ip = min(max(int(aint(dvar,r8)) + 1, 1), n_p - 1) ip1 = ip + 1 wp = dvar - floor(dvar) wp1 = 1.0 - wp dvar = (t_p - min_tp_h2o) / dtp_h2o itp = min(max(int(aint(dvar,r8)) + 1, 1), n_tp - 1) itp1 = itp + 1 wtp = dvar - floor(dvar) wtp1 = 1.0 - wtp t_e = min(max(temh2o(i,kn)-t_p,min_te_h2o),max_te_h2o) dvar = (t_e - min_te_h2o) / dte_h2o ite = min(max(int(aint(dvar,r8)) + 1, 1), n_te - 1) ite1 = ite + 1 wte = dvar - floor(dvar) wte1 = 1.0 - wte rh_path = min(max(q_path / qsx, min_rh_h2o), max_rh_h2o) dvar = (rh_path - min_rh_h2o) / drh_h2o irh = min(max(int(aint(dvar,r8)) + 1, 1), n_rh - 1) irh1 = irh + 1 wrh = dvar - floor(dvar) wrh1 = 1.0 - wrh w_0_0_ = wtp * wte w_0_1_ = wtp * wte1 w_1_0_ = wtp1 * wte w_1_1_ = wtp1 * wte1 w_0_00 = w_0_0_ * wrh w_0_01 = w_0_0_ * wrh1 w_0_10 = w_0_1_ * wrh w_0_11 = w_0_1_ * wrh1 w_1_00 = w_1_0_ * wrh w_1_01 = w_1_0_ * wrh1 w_1_10 = w_1_1_ * wrh w_1_11 = w_1_1_ * wrh1 w00_00 = wp * w_0_00 w00_01 = wp * w_0_01 w00_10 = wp * w_0_10 w00_11 = wp * w_0_11 w01_00 = wp * w_1_00 w01_01 = wp * w_1_01 w01_10 = wp * w_1_10 w01_11 = wp * w_1_11 w10_00 = wp1 * w_0_00 w10_01 = wp1 * w_0_01 w10_10 = wp1 * w_0_10 w10_11 = wp1 * w_0_11 w11_00 = wp1 * w_1_00 w11_01 = wp1 * w_1_01 w11_10 = wp1 * w_1_10 w11_11 = wp1 * w_1_11 ib = 1 fa = fat(1,ib) + & fat(2,ib) * te1 + & fat(3,ib) * te2 + & fat(4,ib) * te3 + & fat(5,ib) * te4 + & fat(6,ib) * te5 a_star = & ah2onw(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & ah2onw(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & ah2onw(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & ah2onw(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & ah2onw(ip , itp , iu1, ite , irh ) * w11_11 * wu + & ah2onw(ip , itp , iu1, ite , irh1) * w11_10 * wu + & ah2onw(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & ah2onw(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & ah2onw(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & ah2onw(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & ah2onw(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & ah2onw(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & ah2onw(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & ah2onw(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & ah2onw(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & ah2onw(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & ah2onw(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & ah2onw(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & ah2onw(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & ah2onw(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & ah2onw(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & ah2onw(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & ah2onw(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & ah2onw(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & ah2onw(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & ah2onw(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & ah2onw(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & ah2onw(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & ah2onw(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & ah2onw(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & ah2onw(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & ah2onw(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu abso(i,ib) = min(max(fa * (1.0 - (1.0 - a_star) * & aer_trn_ngh(i,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o abso(i,ib) = abso(i,ib) * uscl endif ib = 2 fa = fat(1,ib) + & fat(2,ib) * te1 + & fat(3,ib) * te2 + & fat(4,ib) * te3 + & fat(5,ib) * te4 + & fat(6,ib) * te5 a_star = & ah2ow(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & ah2ow(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & ah2ow(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & ah2ow(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & ah2ow(ip , itp , iu1, ite , irh ) * w11_11 * wu + & ah2ow(ip , itp , iu1, ite , irh1) * w11_10 * wu + & ah2ow(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & ah2ow(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & ah2ow(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & ah2ow(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & ah2ow(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & ah2ow(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & ah2ow(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & ah2ow(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & ah2ow(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & ah2ow(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & ah2ow(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & ah2ow(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & ah2ow(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & ah2ow(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & ah2ow(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & ah2ow(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & ah2ow(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & ah2ow(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & ah2ow(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & ah2ow(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & ah2ow(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & ah2ow(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & ah2ow(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & ah2ow(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & ah2ow(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & ah2ow(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu abso(i,ib) = min(max(fa * (1.0 - (1.0 - a_star) * & aer_trn_ngh(i,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o abso(i,ib) = abso(i,ib) * uscl endif end do do i=1,ncol term7(i,1) = coefj(1,1) + coefj(2,1)*dty(i)*(1. + c16*dty(i)) term8(i,1) = coefk(1,1) + coefk(2,1)*dty(i)*(1. + c17*dty(i)) term7(i,2) = coefj(1,2) + coefj(2,2)*dty(i)*(1. + c26*dty(i)) term8(i,2) = coefk(1,2) + coefk(2,2)*dty(i)*(1. + c27*dty(i)) end do do i=1,ncol dtym10 = dty(i) - 10. denom = 1. + (c30 + c31*dtym10*dtym10)*sqrtu(i) k21 = term7(i,1) + term8(i,1)/denom denom = 1. + (c28 + c29*dtym10 )*sqrtu(i) k22 = term7(i,2) + term8(i,2)/denom tr1 = exp(-(k21*(sqrtu(i) + fc1*fwku(i)))) tr2 = exp(-(k22*(sqrtu(i) + fc1*fwku(i)))) tr1=tr1*aer_trn_ngh(i,idx_LW_0650_0800) tr2=tr2*aer_trn_ngh(i,idx_LW_0500_0650) tr5 = exp(-((coefh(1,3) + coefh(2,3)*dtx(i))*uc1(i))) tr6 = exp(-((coefh(1,4) + coefh(2,4)*dtx(i))*uc1(i))) tr9(i) = tr1*tr5 tr10(i) = tr2*tr6 trab2(i)= 0.65*tr9(i) + 0.35*tr10(i) th2o(i) = tr10(i) end do do i=1,ncol te = (tbar(i,kn)*r293)**.7 dplos = abs(plos(i,k2+1) - plos(i,k2)) u1 = zinpl(i,kn)*18.29*dplos/te u2 = zinpl(i,kn)*.5649*dplos/te tlocal = tbar(i,kn) tcrfac = sqrt(tlocal*r250)*te beta = r3205*(pinpl(i,kn)*rsslp + dpfo3*tcrfac) realnu = te/beta tmp1 = u1/sqrt(4. + u1*(1. + realnu)) tmp2 = u2/sqrt(4. + u2*(1. + realnu)) o3bndi = 74.*te*log(1. + tmp1 + tmp2) abso(i,3) = o3bndi*o3emm(i,kn)*(h2otr(i,k2+1)/h2otr(i,k2)) to3(i) = 1.0/(1. + 0.1*tmp1 + 0.1*tmp2) end do do i=1,ncol dplco2 = plco2(i,k2+1) - plco2(i,k2) sqwp = sqrt(uinpl(i,kn)*dplco2) et = exp(-480./tbar(i,kn)) sqti(i) = sqrt(tbar(i,kn)) rsqti = 1./sqti(i) et2 = et*et et4 = et2*et2 omet = (1. - 1.5*et2) f1co2 = 899.70*omet*(1. + 1.94774*et + 4.73486*et2)*rsqti f1sqwp(i)= f1co2*sqwp t1co2(i) = 1./(1. + (245.18*omet*sqwp*rsqti)) oneme = 1. - et2 alphat = oneme**3*rsqti pi = abs(dpnm(i))*winpl(i,kn) wco2 = 2.5221*co2vmr*pi*rga u7(i) = 4.9411e4*alphat*et2*wco2 u8 = 3.9744e4*alphat*et4*wco2 u9 = 1.0447e5*alphat*et4*et2*wco2 u13 = 2.8388e3*alphat*et4*wco2 tpath = tbar(i,kn) tlocal = tbar(i,kn) tcrfac = sqrt((tlocal*r250)*(tpath*r300)) posqt = (pinpl(i,kn)*rsslp + dpfco2*tcrfac)*rsqti rbeta7(i)= 1./(5.3228*posqt) rbeta8 = 1./(10.6576*posqt) rbeta9 = rbeta7(i) rbeta13 = rbeta9 f2co2(i) = u7(i)/sqrt(4. + u7(i)*(1. + rbeta7(i))) + & u8 /sqrt(4. + u8*(1. + rbeta8)) + & u9 /sqrt(4. + u9*(1. + rbeta9)) f3co2(i) = u13/sqrt(4. + u13*(1. + rbeta13)) tmp1 = log(1. + f1sqwp(i)) tmp2 = log(1. + f2co2(i)) tmp3 = log(1. + f3co2(i)) absbnd = (tmp1 + 2.*t1co2(i)*tmp2 + 2.*tmp3)*sqti(i) abso(i,4)= trab2(i)*emm(i,kn)*absbnd tco2(i) = 1.0/(1.0+ 10.0*u7(i)/sqrt(4. + u7(i)*(1. + rbeta7(i)))) end do call trcabn(lchnk ,ncol ,pcols, pverp, & k2 ,kn ,ucfc11 ,ucfc12 ,un2o0 , & un2o1 ,uch4 ,uco211 ,uco212 ,uco213 , & uco221 ,uco222 ,uco223 ,tbar ,bplnk , & winpl ,pinpl ,tco2 ,th2o ,to3 , & uptype ,dw ,s2c ,u ,pnew , & abstrc ,uinpl , & aer_trn_ngh) do i=1,ncol absnxt(i,k2,kn) = abso(i,1) + abso(i,2) + & abso(i,3) + abso(i,4) + abstrc(i) end do end do end do return end subroutine radabs subroutine radems(lchnk ,ncol ,pcols, pver, pverp, & s2c ,tcg ,w ,tplnke ,plh2o , & pnm ,plco2 ,tint ,tint4 ,tlayr , & tlayr4 ,plol ,plos ,ucfc11 ,ucfc12 , & un2o0 ,un2o1 ,uch4 ,uco211 ,uco212 , & uco213 ,uco221 ,uco222 ,uco223 ,uptype , & bn2o0 ,bn2o1 ,bch4 ,co2em ,co2eml , & co2t ,h2otr ,abplnk1 ,abplnk2 ,emstot , & plh2ob ,wb , & aer_trn_ttl) integer, intent(in) :: lchnk integer, intent(in) :: ncol integer, intent(in) :: pcols, pver, pverp real(r8), intent(in) :: s2c(pcols,pverp) real(r8), intent(in) :: tcg(pcols,pverp) real(r8), intent(in) :: w(pcols,pverp) real(r8), intent(in) :: tplnke(pcols) real(r8), intent(in) :: plh2o(pcols,pverp) real(r8), intent(in) :: pnm(pcols,pverp) real(r8), intent(in) :: plco2(pcols,pverp) real(r8), intent(in) :: tint(pcols,pverp) real(r8), intent(in) :: tint4(pcols,pverp) real(r8), intent(in) :: tlayr(pcols,pverp) real(r8), intent(in) :: tlayr4(pcols,pverp) real(r8), intent(in) :: plol(pcols,pverp) real(r8), intent(in) :: plos(pcols,pverp) real(r8), intent(in) :: plh2ob(nbands,pcols,pverp) real(r8), intent(in) :: wb(nbands,pcols,pverp) real(r8), intent(in) :: aer_trn_ttl(pcols,pverp,pverp,bnd_nbr_LW) real(r8), intent(in) :: ucfc11(pcols,pverp) real(r8), intent(in) :: ucfc12(pcols,pverp) real(r8), intent(in) :: un2o0(pcols,pverp) real(r8), intent(in) :: un2o1(pcols,pverp) real(r8), intent(in) :: uch4(pcols,pverp) real(r8), intent(in) :: uco211(pcols,pverp) real(r8), intent(in) :: uco212(pcols,pverp) real(r8), intent(in) :: uco213(pcols,pverp) real(r8), intent(in) :: uco221(pcols,pverp) real(r8), intent(in) :: uco222(pcols,pverp) real(r8), intent(in) :: uco223(pcols,pverp) real(r8), intent(in) :: bn2o0(pcols,pverp) real(r8), intent(in) :: bn2o1(pcols,pverp) real(r8), intent(in) :: bch4(pcols,pverp) real(r8), intent(in) :: uptype(pcols,pverp) real(r8), intent(out) :: emstot(pcols,pverp) real(r8), intent(out) :: co2em(pcols,pverp) real(r8), intent(out) :: co2eml(pcols,pver) real(r8), intent(out) :: co2t(pcols,pverp) real(r8), intent(out) :: h2otr(pcols,pverp) real(r8), intent(out) :: abplnk1(14,pcols,pverp) real(r8), intent(out) :: abplnk2(14,pcols,pverp) integer i integer k integer k1 real(r8) h2oems(pcols,pverp) real(r8) tpathe real(r8) dtx(pcols) real(r8) dty(pcols) real(r8) emis(pcols,2) real(r8) term7(pcols,2) real(r8) term8(pcols,2) real(r8) tr1(pcols) real(r8) tr2(pcols) real(r8) tr3(pcols) real(r8) tr4(pcols) real(r8) tr7(pcols) real(r8) tr8(pcols) real(r8) k21(pcols) real(r8) k22(pcols) real(r8) u(pcols) real(r8) ub(nbands) real(r8) pnew real(r8) pnewb(nbands) real(r8) uc1(pcols) real(r8) fwk real(r8) troco2(pcols,pverp) real(r8) emplnk(14,pcols) real(r8) emstrc(pcols,pverp) real(r8) co2ems(pcols,pverp) real(r8) co2plk(pcols) real(r8) sum(pcols) real(r8) t1i real(r8) sqti real(r8) pi real(r8) et real(r8) et2 real(r8) et4 real(r8) omet real(r8) ex real(r8) f1co2 real(r8) f2co2 real(r8) f3co2 real(r8) t1co2 real(r8) sqwp real(r8) f1sqwp real(r8) oneme real(r8) alphat real(r8) wco2 real(r8) posqt real(r8) rbeta7 real(r8) rbeta8 real(r8) rbeta9 real(r8) rbeta13 real(r8) tpath real(r8) tmp1 real(r8) tmp2 real(r8) tmp3 real(r8) tlayr5 real(r8) rsqti real(r8) exm1sq real(r8) u7 real(r8) u8 real(r8) u9 real(r8) u13 real(r8) r250 real(r8) r300 real(r8) rsslp real(r8) o3ems(pcols,pverp) real(r8) dbvtt(pcols) real(r8) dbvt,fo3,t,ux,vx real(r8) te real(r8) u1 real(r8) u2 real(r8) phat real(r8) tlocal real(r8) tcrfac real(r8) beta real(r8) realnu real(r8) o3bndi real(r8) absbnd real(r8) tco2(pcols) real(r8) th2o(pcols) real(r8) to3(pcols) real(r8) fe real(r8) e_star real(r8) l_star real(r8) c_star real(r8) te1 real(r8) te2 real(r8) te3 real(r8) te4 real(r8) te5 real(r8) log_u real(r8) log_uc real(r8) log_p real(r8) t_p real(r8) t_e integer iu integer iu1 integer iuc integer iuc1 integer ip integer ip1 integer itp integer itp1 integer ite integer ite1 integer irh integer irh1 real(r8) dvar real(r8) uvar real(r8) uscl real(r8) wu real(r8) wu1 real(r8) wuc real(r8) wuc1 real(r8) wp real(r8) wp1 real(r8) wtp real(r8) wtp1 real(r8) wte real(r8) wte1 real(r8) wrh real(r8) wrh1 real(r8) w_0_0_ real(r8) w_0_1_ real(r8) w_1_0_ real(r8) w_1_1_ real(r8) w_0_00 real(r8) w_0_01 real(r8) w_0_10 real(r8) w_0_11 real(r8) w_1_00 real(r8) w_1_01 real(r8) w_1_10 real(r8) w_1_11 real(r8) w00_00 real(r8) w00_01 real(r8) w00_10 real(r8) w00_11 real(r8) w01_00 real(r8) w01_01 real(r8) w01_10 real(r8) w01_11 real(r8) w10_00 real(r8) w10_01 real(r8) w10_10 real(r8) w10_11 real(r8) w11_00 real(r8) w11_01 real(r8) w11_10 real(r8) w11_11 integer ib real(r8) pch2o real(r8) fch2o real(r8) uch2o real(r8) fdif real(r8) sslp_mks real(r8) esx real(r8) qsx real(r8) pnew_mks real(r8) q_path real(r8) rh_path real(r8) omeps integer iest dbvt(t)=(-2.8911366682e-4+(2.3771251896e-6+1.1305188929e-10*t)*t)/ & (1.0+(-6.1364820707e-3+1.5550319767e-5*t)*t) fo3(ux,vx)=ux/sqrt(4.+ux*(1.+vx)) r250 = 1./250. r300 = 1./300. rsslp = 1./sslp fdif = 1.66 sslp_mks = sslp / 10.0 omeps = 1.0 - epsilo do i=1,ncol ex = exp(960./tplnke(i)) co2plk(i) = 5.e8/((tplnke(i)**4)*(ex - 1.)) co2t(i,ntoplw) = tplnke(i) sum(i) = co2t(i,ntoplw)*pnm(i,ntoplw) end do k = ntoplw do k1=pverp,ntoplw+1,-1 k = k + 1 do i=1,ncol sum(i) = sum(i) + tlayr(i,k)*(pnm(i,k)-pnm(i,k-1)) ex = exp(960./tlayr(i,k1)) tlayr5 = tlayr(i,k1)*tlayr4(i,k1) co2eml(i,k1-1) = 1.2e11*ex/(tlayr5*(ex - 1.)**2) co2t(i,k) = sum(i)/pnm(i,k) end do end do do i=1,ncol dbvtt(i) = dbvt(tplnke(i)) end do call trcplk(lchnk ,ncol ,pcols, pver, pverp, & tint ,tlayr ,tplnke ,emplnk ,abplnk1 , & abplnk2 ) do k1=ntoplw,pverp do i=1,ncol u(i) = plh2o(i,k1) pnew = u(i)/w(i,k1) pnew_mks = pnew * sslp_mks uc1(i) = (s2c(i,k1) + 1.7e-3*plh2o(i,k1))*(1. + 2.*s2c(i,k1))/ & (1. + 15.*s2c(i,k1)) pch2o = s2c(i,k1) tpathe = tcg(i,k1)/w(i,k1) t_p = min(max(tpathe, min_tp_h2o), max_tp_h2o) iest = floor(t_p) - min_tp_h2o esx = estblh2o(iest) + (estblh2o(iest+1)-estblh2o(iest)) * & (t_p - min_tp_h2o - iest) qsx = epsilo * esx / (pnew_mks - omeps * esx) q_path = w(i,k1) / pnm(i,k1) / rga ub(1) = plh2ob(1,i,k1) / psi(t_p,1) ub(2) = plh2ob(2,i,k1) / psi(t_p,2) pnewb(1) = ub(1) / wb(1,i,k1) * phi(t_p,1) pnewb(2) = ub(2) / wb(2,i,k1) * phi(t_p,2) dtx(i) = tplnke(i) - 250. dty(i) = tpathe - 250. te1 = tplnke(i) te2 = te1 * te1 te3 = te2 * te1 te4 = te3 * te1 te5 = te4 * te1 dvar = (t_p - min_tp_h2o) / dtp_h2o itp = min(max(int(aint(dvar,r8)) + 1, 1), n_tp - 1) itp1 = itp + 1 wtp = dvar - floor(dvar) wtp1 = 1.0 - wtp t_e = min(max(tplnke(i) - t_p, min_te_h2o), max_te_h2o) dvar = (t_e - min_te_h2o) / dte_h2o ite = min(max(int(aint(dvar,r8)) + 1, 1), n_te - 1) ite1 = ite + 1 wte = dvar - floor(dvar) wte1 = 1.0 - wte rh_path = min(max(q_path / qsx, min_rh_h2o), max_rh_h2o) dvar = (rh_path - min_rh_h2o) / drh_h2o irh = min(max(int(aint(dvar,r8)) + 1, 1), n_rh - 1) irh1 = irh + 1 wrh = dvar - floor(dvar) wrh1 = 1.0 - wrh w_0_0_ = wtp * wte w_0_1_ = wtp * wte1 w_1_0_ = wtp1 * wte w_1_1_ = wtp1 * wte1 w_0_00 = w_0_0_ * wrh w_0_01 = w_0_0_ * wrh1 w_0_10 = w_0_1_ * wrh w_0_11 = w_0_1_ * wrh1 w_1_00 = w_1_0_ * wrh w_1_01 = w_1_0_ * wrh1 w_1_10 = w_1_1_ * wrh w_1_11 = w_1_1_ * wrh1 fch2o = fh2oself(t_p) uch2o = (pch2o * epsilo) / (q_path * fch2o) ib = 1 uvar = ub(ib) * fdif log_u = min(log10(max(uvar, min_u_h2o)), max_lu_h2o) dvar = (log_u - min_lu_h2o) / dlu_h2o iu = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iu1 = iu + 1 wu = dvar - floor(dvar) wu1 = 1.0 - wu log_p = min(log10(max(pnewb(ib), min_p_h2o)), max_lp_h2o) dvar = (log_p - min_lp_h2o) / dlp_h2o ip = min(max(int(aint(dvar,r8)) + 1, 1), n_p - 1) ip1 = ip + 1 wp = dvar - floor(dvar) wp1 = 1.0 - wp w00_00 = wp * w_0_00 w00_01 = wp * w_0_01 w00_10 = wp * w_0_10 w00_11 = wp * w_0_11 w01_00 = wp * w_1_00 w01_01 = wp * w_1_01 w01_10 = wp * w_1_10 w01_11 = wp * w_1_11 w10_00 = wp1 * w_0_00 w10_01 = wp1 * w_0_01 w10_10 = wp1 * w_0_10 w10_11 = wp1 * w_0_11 w11_00 = wp1 * w_1_00 w11_01 = wp1 * w_1_01 w11_10 = wp1 * w_1_10 w11_11 = wp1 * w_1_11 fe = fet(1,ib) + & fet(2,ib) * te1 + & fet(3,ib) * te2 + & fet(4,ib) * te3 + & fet(5,ib) * te4 + & fet(6,ib) * te5 e_star = & eh2onw(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & eh2onw(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & eh2onw(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & eh2onw(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & eh2onw(ip , itp , iu1, ite , irh ) * w11_11 * wu + & eh2onw(ip , itp , iu1, ite , irh1) * w11_10 * wu + & eh2onw(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & eh2onw(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & eh2onw(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & eh2onw(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & eh2onw(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & eh2onw(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & eh2onw(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & eh2onw(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & eh2onw(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & eh2onw(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & eh2onw(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & eh2onw(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & eh2onw(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & eh2onw(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & eh2onw(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & eh2onw(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & eh2onw(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & eh2onw(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & eh2onw(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & eh2onw(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & eh2onw(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & eh2onw(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & eh2onw(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & eh2onw(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & eh2onw(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & eh2onw(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu emis(i,ib) = min(max(fe * (1.0 - (1.0 - e_star) * & aer_trn_ttl(i,k1,1,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o emis(i,ib) = emis(i,ib) * uscl endif ib = 2 uvar = ub(ib) * fdif log_u = min(log10(max(uvar, min_u_h2o)), max_lu_h2o) dvar = (log_u - min_lu_h2o) / dlu_h2o iu = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iu1 = iu + 1 wu = dvar - floor(dvar) wu1 = 1.0 - wu log_p = min(log10(max(pnewb(ib), min_p_h2o)), max_lp_h2o) dvar = (log_p - min_lp_h2o) / dlp_h2o ip = min(max(int(aint(dvar,r8)) + 1, 1), n_p - 1) ip1 = ip + 1 wp = dvar - floor(dvar) wp1 = 1.0 - wp w00_00 = wp * w_0_00 w00_01 = wp * w_0_01 w00_10 = wp * w_0_10 w00_11 = wp * w_0_11 w01_00 = wp * w_1_00 w01_01 = wp * w_1_01 w01_10 = wp * w_1_10 w01_11 = wp * w_1_11 w10_00 = wp1 * w_0_00 w10_01 = wp1 * w_0_01 w10_10 = wp1 * w_0_10 w10_11 = wp1 * w_0_11 w11_00 = wp1 * w_1_00 w11_01 = wp1 * w_1_01 w11_10 = wp1 * w_1_10 w11_11 = wp1 * w_1_11 log_uc = min(log10(max(uch2o * fdif, min_u_h2o)), max_lu_h2o) dvar = (log_uc - min_lu_h2o) / dlu_h2o iuc = min(max(int(aint(dvar,r8)) + 1, 1), n_u - 1) iuc1 = iuc + 1 wuc = dvar - floor(dvar) wuc1 = 1.0 - wuc fe = fet(1,ib) + & fet(2,ib) * te1 + & fet(3,ib) * te2 + & fet(4,ib) * te3 + & fet(5,ib) * te4 + & fet(6,ib) * te5 l_star = & ln_eh2ow(ip , itp , iu , ite , irh ) * w11_11 * wu1 + & ln_eh2ow(ip , itp , iu , ite , irh1) * w11_10 * wu1 + & ln_eh2ow(ip , itp , iu , ite1, irh ) * w11_01 * wu1 + & ln_eh2ow(ip , itp , iu , ite1, irh1) * w11_00 * wu1 + & ln_eh2ow(ip , itp , iu1, ite , irh ) * w11_11 * wu + & ln_eh2ow(ip , itp , iu1, ite , irh1) * w11_10 * wu + & ln_eh2ow(ip , itp , iu1, ite1, irh ) * w11_01 * wu + & ln_eh2ow(ip , itp , iu1, ite1, irh1) * w11_00 * wu + & ln_eh2ow(ip , itp1, iu , ite , irh ) * w10_11 * wu1 + & ln_eh2ow(ip , itp1, iu , ite , irh1) * w10_10 * wu1 + & ln_eh2ow(ip , itp1, iu , ite1, irh ) * w10_01 * wu1 + & ln_eh2ow(ip , itp1, iu , ite1, irh1) * w10_00 * wu1 + & ln_eh2ow(ip , itp1, iu1, ite , irh ) * w10_11 * wu + & ln_eh2ow(ip , itp1, iu1, ite , irh1) * w10_10 * wu + & ln_eh2ow(ip , itp1, iu1, ite1, irh ) * w10_01 * wu + & ln_eh2ow(ip , itp1, iu1, ite1, irh1) * w10_00 * wu + & ln_eh2ow(ip1, itp , iu , ite , irh ) * w01_11 * wu1 + & ln_eh2ow(ip1, itp , iu , ite , irh1) * w01_10 * wu1 + & ln_eh2ow(ip1, itp , iu , ite1, irh ) * w01_01 * wu1 + & ln_eh2ow(ip1, itp , iu , ite1, irh1) * w01_00 * wu1 + & ln_eh2ow(ip1, itp , iu1, ite , irh ) * w01_11 * wu + & ln_eh2ow(ip1, itp , iu1, ite , irh1) * w01_10 * wu + & ln_eh2ow(ip1, itp , iu1, ite1, irh ) * w01_01 * wu + & ln_eh2ow(ip1, itp , iu1, ite1, irh1) * w01_00 * wu + & ln_eh2ow(ip1, itp1, iu , ite , irh ) * w00_11 * wu1 + & ln_eh2ow(ip1, itp1, iu , ite , irh1) * w00_10 * wu1 + & ln_eh2ow(ip1, itp1, iu , ite1, irh ) * w00_01 * wu1 + & ln_eh2ow(ip1, itp1, iu , ite1, irh1) * w00_00 * wu1 + & ln_eh2ow(ip1, itp1, iu1, ite , irh ) * w00_11 * wu + & ln_eh2ow(ip1, itp1, iu1, ite , irh1) * w00_10 * wu + & ln_eh2ow(ip1, itp1, iu1, ite1, irh ) * w00_01 * wu + & ln_eh2ow(ip1, itp1, iu1, ite1, irh1) * w00_00 * wu c_star = & cn_eh2ow(ip , itp , iuc , ite , irh ) * w11_11 * wuc1 + & cn_eh2ow(ip , itp , iuc , ite , irh1) * w11_10 * wuc1 + & cn_eh2ow(ip , itp , iuc , ite1, irh ) * w11_01 * wuc1 + & cn_eh2ow(ip , itp , iuc , ite1, irh1) * w11_00 * wuc1 + & cn_eh2ow(ip , itp , iuc1, ite , irh ) * w11_11 * wuc + & cn_eh2ow(ip , itp , iuc1, ite , irh1) * w11_10 * wuc + & cn_eh2ow(ip , itp , iuc1, ite1, irh ) * w11_01 * wuc + & cn_eh2ow(ip , itp , iuc1, ite1, irh1) * w11_00 * wuc + & cn_eh2ow(ip , itp1, iuc , ite , irh ) * w10_11 * wuc1 + & cn_eh2ow(ip , itp1, iuc , ite , irh1) * w10_10 * wuc1 + & cn_eh2ow(ip , itp1, iuc , ite1, irh ) * w10_01 * wuc1 + & cn_eh2ow(ip , itp1, iuc , ite1, irh1) * w10_00 * wuc1 + & cn_eh2ow(ip , itp1, iuc1, ite , irh ) * w10_11 * wuc + & cn_eh2ow(ip , itp1, iuc1, ite , irh1) * w10_10 * wuc + & cn_eh2ow(ip , itp1, iuc1, ite1, irh ) * w10_01 * wuc + & cn_eh2ow(ip , itp1, iuc1, ite1, irh1) * w10_00 * wuc + & cn_eh2ow(ip1, itp , iuc , ite , irh ) * w01_11 * wuc1 + & cn_eh2ow(ip1, itp , iuc , ite , irh1) * w01_10 * wuc1 + & cn_eh2ow(ip1, itp , iuc , ite1, irh ) * w01_01 * wuc1 + & cn_eh2ow(ip1, itp , iuc , ite1, irh1) * w01_00 * wuc1 + & cn_eh2ow(ip1, itp , iuc1, ite , irh ) * w01_11 * wuc + & cn_eh2ow(ip1, itp , iuc1, ite , irh1) * w01_10 * wuc + & cn_eh2ow(ip1, itp , iuc1, ite1, irh ) * w01_01 * wuc + & cn_eh2ow(ip1, itp , iuc1, ite1, irh1) * w01_00 * wuc + & cn_eh2ow(ip1, itp1, iuc , ite , irh ) * w00_11 * wuc1 + & cn_eh2ow(ip1, itp1, iuc , ite , irh1) * w00_10 * wuc1 + & cn_eh2ow(ip1, itp1, iuc , ite1, irh ) * w00_01 * wuc1 + & cn_eh2ow(ip1, itp1, iuc , ite1, irh1) * w00_00 * wuc1 + & cn_eh2ow(ip1, itp1, iuc1, ite , irh ) * w00_11 * wuc + & cn_eh2ow(ip1, itp1, iuc1, ite , irh1) * w00_10 * wuc + & cn_eh2ow(ip1, itp1, iuc1, ite1, irh ) * w00_01 * wuc + & cn_eh2ow(ip1, itp1, iuc1, ite1, irh1) * w00_00 * wuc emis(i,ib) = min(max(fe * (1.0 - l_star * c_star * & aer_trn_ttl(i,k1,1,ib)), & 0.0_r8), 1.0_r8) if (uvar < min_u_h2o) then uscl = uvar / min_u_h2o emis(i,ib) = emis(i,ib) * uscl endif h2oems(i,k1) = emis(i,1)+emis(i,2) end do do i=1,ncol term7(i,1) = coefj(1,1) + coefj(2,1)*dty(i)*(1.+c16*dty(i)) term8(i,1) = coefk(1,1) + coefk(2,1)*dty(i)*(1.+c17*dty(i)) term7(i,2) = coefj(1,2) + coefj(2,2)*dty(i)*(1.+c26*dty(i)) term8(i,2) = coefk(1,2) + coefk(2,2)*dty(i)*(1.+c27*dty(i)) end do do i=1,ncol k21(i) = term7(i,1) + term8(i,1)/ & (1. + (c30 + c31*(dty(i)-10.)*(dty(i)-10.))*sqrt(u(i))) k22(i) = term7(i,2) + term8(i,2)/ & (1. + (c28 + c29*(dty(i)-10.))*sqrt(u(i))) fwk = fwcoef + fwc1/(1.+fwc2*u(i)) tr1(i) = exp(-(k21(i)*(sqrt(u(i)) + fc1*fwk*u(i)))) tr2(i) = exp(-(k22(i)*(sqrt(u(i)) + fc1*fwk*u(i)))) tr1(i)=tr1(i)*aer_trn_ttl(i,k1,1,idx_LW_0650_0800) tr2(i)=tr2(i)*aer_trn_ttl(i,k1,1,idx_LW_0500_0650) tr3(i) = exp(-((coefh(1,1) + coefh(2,1)*dtx(i))*uc1(i))) tr4(i) = exp(-((coefh(1,2) + coefh(2,2)*dtx(i))*uc1(i))) tr7(i) = tr1(i)*tr3(i) tr8(i) = tr2(i)*tr4(i) troco2(i,k1) = 0.65*tr7(i) + 0.35*tr8(i) th2o(i) = tr8(i) end do do i=1,ncol t1i = exp(-480./co2t(i,k1)) sqti = sqrt(co2t(i,k1)) rsqti = 1./sqti et = t1i et2 = et*et et4 = et2*et2 omet = 1. - 1.5*et2 f1co2 = 899.70*omet*(1. + 1.94774*et + 4.73486*et2)*rsqti sqwp = sqrt(plco2(i,k1)) f1sqwp = f1co2*sqwp t1co2 = 1./(1. + 245.18*omet*sqwp*rsqti) oneme = 1. - et2 alphat = oneme**3*rsqti wco2 = 2.5221*co2vmr*pnm(i,k1)*rga u7 = 4.9411e4*alphat*et2*wco2 u8 = 3.9744e4*alphat*et4*wco2 u9 = 1.0447e5*alphat*et4*et2*wco2 u13 = 2.8388e3*alphat*et4*wco2 tpath = co2t(i,k1) tlocal = tplnke(i) tcrfac = sqrt((tlocal*r250)*(tpath*r300)) pi = pnm(i,k1)*rsslp + 2.*dpfco2*tcrfac posqt = pi/(2.*sqti) rbeta7 = 1./( 5.3288*posqt) rbeta8 = 1./ (10.6576*posqt) rbeta9 = rbeta7 rbeta13= rbeta9 f2co2 = (u7/sqrt(4. + u7*(1. + rbeta7))) + & (u8/sqrt(4. + u8*(1. + rbeta8))) + & (u9/sqrt(4. + u9*(1. + rbeta9))) f3co2 = u13/sqrt(4. + u13*(1. + rbeta13)) tmp1 = log(1. + f1sqwp) tmp2 = log(1. + f2co2) tmp3 = log(1. + f3co2) absbnd = (tmp1 + 2.*t1co2*tmp2 + 2.*tmp3)*sqti tco2(i)=1.0/(1.0+10.0*(u7/sqrt(4. + u7*(1. + rbeta7)))) co2ems(i,k1) = troco2(i,k1)*absbnd*co2plk(i) ex = exp(960./tint(i,k1)) exm1sq = (ex - 1.)**2 co2em(i,k1) = 1.2e11*ex/(tint(i,k1)*tint4(i,k1)*exm1sq) end do do i=1,ncol h2otr(i,k1) = exp(-12.*s2c(i,k1)) h2otr(i,k1)=h2otr(i,k1)*aer_trn_ttl(i,k1,1,idx_LW_1000_1200) te = (co2t(i,k1)/293.)**.7 u1 = 18.29*plos(i,k1)/te u2 = .5649*plos(i,k1)/te phat = plos(i,k1)/plol(i,k1) tlocal = tplnke(i) tcrfac = sqrt(tlocal*r250)*te beta = (1./.3205)*((1./phat) + (dpfo3*tcrfac)) realnu = (1./beta)*te o3bndi = 74.*te*(tplnke(i)/375.)*log(1. + fo3(u1,realnu) + fo3(u2,realnu)) o3ems(i,k1) = dbvtt(i)*h2otr(i,k1)*o3bndi to3(i)=1.0/(1. + 0.1*fo3(u1,realnu) + 0.1*fo3(u2,realnu)) end do call trcems(lchnk ,ncol ,pcols, pverp, & k1 ,co2t ,pnm ,ucfc11 ,ucfc12 , & un2o0 ,un2o1 ,bn2o0 ,bn2o1 ,uch4 , & bch4 ,uco211 ,uco212 ,uco213 ,uco221 , & uco222 ,uco223 ,uptype ,w ,s2c , & u ,emplnk ,th2o ,tco2 ,to3 , & emstrc , & aer_trn_ttl) do i=1,ncol emstot(i,k1) = h2oems(i,k1) + co2ems(i,k1) + o3ems(i,k1) & + emstrc(i,k1) end do end do return end subroutine radems subroutine radtpl(lchnk ,ncol ,pcols, pver, pverp, & tnm ,lwupcgs ,qnm ,pnm ,plco2 ,plh2o , & tplnka ,s2c ,tcg ,w ,tplnke , & tint ,tint4 ,tlayr ,tlayr4 ,pmln , & piln ,plh2ob ,wb ) integer, intent(in) :: lchnk integer, intent(in) :: ncol integer, intent(in) :: pcols, pver, pverp real(r8), intent(in) :: tnm(pcols,pver) real(r8), intent(in) :: lwupcgs(pcols) real(r8), intent(in) :: qnm(pcols,pver) real(r8), intent(in) :: pnm(pcols,pverp) real(r8), intent(in) :: pmln(pcols,pver) real(r8), intent(in) :: piln(pcols,pverp) real(r8), intent(out) :: plco2(pcols,pverp) real(r8), intent(out) :: plh2o(pcols,pverp) real(r8), intent(out) :: tplnka(pcols,pverp) real(r8), intent(out) :: s2c(pcols,pverp) real(r8), intent(out) :: tcg(pcols,pverp) real(r8), intent(out) :: w(pcols,pverp) real(r8), intent(out) :: tplnke(pcols) real(r8), intent(out) :: tint(pcols,pverp) real(r8), intent(out) :: tint4(pcols,pverp) real(r8), intent(out) :: tlayr(pcols,pverp) real(r8), intent(out) :: tlayr4(pcols,pverp) real(r8), intent(out) :: plh2ob(nbands,pcols,pverp) real(r8), intent(out) :: wb(nbands,pcols,pverp) integer i integer k integer kp1 real(r8) repsil real(r8) dy real(r8) dpnm real(r8) dpnmsq real(r8) dw real(r8) dplh2o real(r8) cpwpl repsil = 1./epsilo cpwpl = amco2/amd * 0.5/(gravit*p0) do i=1,ncol plh2o(i,ntoplw) = rgsslp*qnm(i,ntoplw)*pnm(i,ntoplw)*pnm(i,ntoplw) plco2(i,ntoplw) = co2vmr*cpwpl*pnm(i,ntoplw)*pnm(i,ntoplw) end do do k=ntoplw,pver do i=1,ncol plh2o(i,k+1) = plh2o(i,k) + rgsslp* & (pnm(i,k+1)**2 - pnm(i,k)**2)*qnm(i,k) plco2(i,k+1) = co2vmr*cpwpl*pnm(i,k+1)**2 end do end do do i=1,ncol tint4(i,pverp) = lwupcgs(i)/stebol tint(i,pverp) = sqrt(sqrt(tint4(i,pverp))) tplnka(i,ntoplw) = tnm(i,ntoplw) tint(i,ntoplw) = tplnka(i,ntoplw) tlayr4(i,ntoplw) = tplnka(i,ntoplw)**4 tint4(i,ntoplw) = tlayr4(i,ntoplw) end do do k=ntoplw+1,pver do i=1,ncol dy = (piln(i,k) - pmln(i,k))/(pmln(i,k-1) - pmln(i,k)) tint(i,k) = tnm(i,k) - dy*(tnm(i,k)-tnm(i,k-1)) tint4(i,k) = tint(i,k)**4 end do end do do k=ntoplw+1,pverp do i=1,ncol tlayr(i,k) = tnm(i,k-1) tlayr4(i,k) = tlayr(i,k)**4 tplnka(i,k) = .5*(tint(i,k) + tint(i,k-1)) end do end do do i=1,ncol tplnke(i) = tplnka(i,ntoplw) tlayr(i,ntoplw) = tint(i,ntoplw) end do do i=1,ncol tcg(i,ntoplw) = rga*qnm(i,ntoplw)*pnm(i,ntoplw)*tnm(i,ntoplw) w(i,ntoplw) = sslp * (plh2o(i,ntoplw)*2.) / pnm(i,ntoplw) wb(1,i,ntoplw) = w(i,ntoplw) * phi(tnm(i,ntoplw),1) wb(2,i,ntoplw) = w(i,ntoplw) * phi(tnm(i,ntoplw),2) plh2ob(1,i,ntoplw) = plh2o(i,ntoplw) * psi(tnm(i,ntoplw),1) plh2ob(2,i,ntoplw) = plh2o(i,ntoplw) * psi(tnm(i,ntoplw),2) s2c(i,ntoplw) = plh2o(i,ntoplw)*fh2oself(tnm(i,ntoplw))*qnm(i,ntoplw)*repsil end do do k=ntoplw,pver do i=1,ncol dpnm = pnm(i,k+1) - pnm(i,k) dpnmsq = pnm(i,k+1)**2 - pnm(i,k)**2 dw = rga*qnm(i,k)*dpnm kp1 = k+1 w(i,kp1) = w(i,k) + dw wb(1,i,kp1) = wb(1,i,k) + dw * phi(tnm(i,k),1) wb(2,i,kp1) = wb(2,i,k) + dw * phi(tnm(i,k),2) dplh2o = plh2o(i,kp1) - plh2o(i,k) plh2ob(1,i,kp1) = plh2ob(1,i,k) + dplh2o * psi(tnm(i,k),1) plh2ob(2,i,kp1) = plh2ob(2,i,k) + dplh2o * psi(tnm(i,k),2) tcg(i,kp1) = tcg(i,k) + dw*tnm(i,k) s2c(i,kp1) = s2c(i,k) + rgsslp*dpnmsq*qnm(i,k)* & fh2oself(tnm(i,k))*qnm(i,k)*repsil end do end do return end subroutine radtpl subroutine radclwmx(lchnk ,ncol ,pcols, pver, pverp, & lwupcgs ,tnm ,qnm ,o3vmr , & pmid ,pint ,pmln ,piln , & n2o ,ch4 ,cfc11 ,cfc12 , & cld ,emis ,pmxrgn ,nmxrgn ,qrl , & doabsems, abstot, absnxt, emstot, & flns ,flnt ,flnsc ,flntc ,flwds , & flut ,flutc , & flup ,flupc ,fldn ,fldnc , & aer_mass) implicit none integer pverp2,pverp3,pverp4 real(r8) cldmin parameter (cldmin = 1.0d-80) integer, intent(in) :: lchnk integer, intent(in) :: pcols, pver, pverp integer, intent(in) :: ncol integer, intent(in) :: nmxrgn(pcols) logical, intent(in) :: doabsems real(r8), intent(in) :: pmxrgn(pcols,pverp) real(r8), intent(in) :: lwupcgs(pcols) real(r8), intent(in) :: tnm(pcols,pver) real(r8), intent(in) :: qnm(pcols,pver) real(r8), intent(in) :: o3vmr(pcols,pver) real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: pmln(pcols,pver) real(r8), intent(in) :: piln(pcols,pverp) real(r8), intent(in) :: n2o(pcols,pver) real(r8), intent(in) :: ch4(pcols,pver) real(r8), intent(in) :: cfc11(pcols,pver) real(r8), intent(in) :: cfc12(pcols,pver) real(r8), intent(in) :: cld(pcols,pver) real(r8), intent(in) :: emis(pcols,pver) real(r8), intent(in) :: aer_mass(pcols,pver) real(r8), intent(out) :: qrl(pcols,pver) real(r8), intent(out) :: flns(pcols) real(r8), intent(out) :: flnt(pcols) real(r8), intent(out) :: flut(pcols) real(r8), intent(out) :: flnsc(pcols) real(r8), intent(out) :: flntc(pcols) real(r8), intent(out) :: flutc(pcols) real(r8), intent(out) :: flwds(pcols) real(r8), intent(out) :: flup(pcols,pverp) real(r8), intent(out) :: flupc(pcols,pverp) real(r8), intent(out) :: fldn(pcols,pverp) real(r8), intent(out) :: fldnc(pcols,pverp) real(r8), intent(inout) :: abstot(pcols,pverp,pverp) real(r8), intent(inout) :: absnxt(pcols,pver,4) real(r8), intent(inout) :: emstot(pcols,pverp) integer i integer ilon integer ii integer iimx integer k integer k1 integer k2 integer k3 integer km integer km1 integer km3 integer km4 integer irgn integer l integer l1 integer n real(r8) :: plco2(pcols,pverp) real(r8) :: plh2o(pcols,pverp) real(r8) tmp(pcols) real(r8) tmp2(pcols) real(r8) absbt(pcols) real(r8) plol(pcols,pverp) real(r8) plos(pcols,pverp) real(r8) aer_mpp(pcols,pverp) real(r8) co2em(pcols,pverp) real(r8) co2eml(pcols,pver) real(r8) delt(pcols) real(r8) delt1(pcols) real(r8) bk1(pcols) real(r8) bk2(pcols) real(r8) cldp(pcols,pverp) real(r8) ful(pcols,pverp) real(r8) fsul(pcols,pverp) real(r8) fdl(pcols,pverp) real(r8) fsdl(pcols,pverp) real(r8) fclb4(pcols,-1:pver) real(r8) fclt4(pcols,0:pver) real(r8) s(pcols,pverp,pverp) real(r8) tplnka(pcols,pverp) real(r8) s2c(pcols,pverp) real(r8) tcg(pcols,pverp) real(r8) w(pcols,pverp) real(r8) tplnke(pcols) real(r8) h2otr(pcols,pverp) real(r8) co2t(pcols,pverp) real(r8) tint(pcols,pverp) real(r8) tint4(pcols,pverp) real(r8) tlayr(pcols,pverp) real(r8) tlayr4(pcols,pverp) real(r8) plh2ob(nbands,pcols,pverp) real(r8) wb(nbands,pcols,pverp) real(r8) cld0 real(r8) cld1 real(r8) emx(0:pverp) real(r8) emx0 real(r8) trans real(r8) asort(pver) real(r8) atmp real(r8) maxcld(pcols) integer indx(pcols) integer indxmx(pcols,pverp) integer nrgn(pcols) integer npts integer ncolmx(pverp) integer kx1(pcols,pverp) integer kx2(pcols,0:pverp) integer kxs(0:pverp,pcols,pverp) integer nxs(pcols,pverp) integer nxsk integer ksort(0:pverp) integer ktmp real(r8) aer_trn_ttl(pcols,pverp,pverp,bnd_nbr_LW) real(r8) ucfc11(pcols,pverp) real(r8) ucfc12(pcols,pverp) real(r8) un2o0(pcols,pverp) real(r8) un2o1(pcols,pverp) real(r8) uch4(pcols,pverp) real(r8) uco211(pcols,pverp) real(r8) uco212(pcols,pverp) real(r8) uco213(pcols,pverp) real(r8) uco221(pcols,pverp) real(r8) uco222(pcols,pverp) real(r8) uco223(pcols,pverp) real(r8) bn2o0(pcols,pverp) real(r8) bn2o1(pcols,pverp) real(r8) bch4(pcols,pverp) real(r8) uptype(pcols,pverp) real(r8) abplnk1(14,pcols,pverp) real(r8) abplnk2(14,pcols,pverp) pverp2=pver+2 pverp3=pver+3 pverp4=pver+4 call aer_pth(aer_mass, aer_mpp, ncol, pcols, pver, pverp) call radtpl(lchnk ,ncol ,pcols, pver, pverp, & tnm ,lwupcgs ,qnm ,pint ,plco2 ,plh2o , & tplnka ,s2c ,tcg ,w ,tplnke , & tint ,tint4 ,tlayr ,tlayr4 ,pmln , & piln ,plh2ob ,wb ) if (doabsems) then call radoz2(lchnk, ncol, pcols, pver, pverp, o3vmr ,pint ,plol ,plos, ntoplw ) call trcpth(lchnk ,ncol ,pcols, pver, pverp, & tnm ,pint ,cfc11 ,cfc12 ,n2o , & ch4 ,qnm ,ucfc11 ,ucfc12 ,un2o0 , & un2o1 ,uch4 ,uco211 ,uco212 ,uco213 , & uco221 ,uco222 ,uco223 ,bn2o0 ,bn2o1 , & bch4 ,uptype ) call aer_trn(aer_mpp, aer_trn_ttl, pcols, pver, pverp) call radems(lchnk ,ncol ,pcols, pver, pverp, & s2c ,tcg ,w ,tplnke ,plh2o , & pint ,plco2 ,tint ,tint4 ,tlayr , & tlayr4 ,plol ,plos ,ucfc11 ,ucfc12 , & un2o0 ,un2o1 ,uch4 ,uco211 ,uco212 , & uco213 ,uco221 ,uco222 ,uco223 ,uptype , & bn2o0 ,bn2o1 ,bch4 ,co2em ,co2eml , & co2t ,h2otr ,abplnk1 ,abplnk2 ,emstot , & plh2ob ,wb , & aer_trn_ttl) call radabs(lchnk ,ncol ,pcols, pver, pverp, & pmid ,pint ,co2em ,co2eml ,tplnka , & s2c ,tcg ,w ,h2otr ,plco2 , & plh2o ,co2t ,tint ,tlayr ,plol , & plos ,pmln ,piln ,ucfc11 ,ucfc12 , & un2o0 ,un2o1 ,uch4 ,uco211 ,uco212 , & uco213 ,uco221 ,uco222 ,uco223 ,uptype , & bn2o0 ,bn2o1 ,bch4 ,abplnk1 ,abplnk2 , & abstot ,absnxt ,plh2ob ,wb , & aer_mpp ,aer_trn_ttl) end if do i=1,ncol delt(i) = tint4(i,pver) - tlayr4(i,pverp) delt1(i) = tlayr4(i,pverp) - tint4(i,pverp) s(i,pverp,pverp) = stebol*(delt1(i)*absnxt(i,pver,1) + delt (i)*absnxt(i,pver,4)) s(i,pver,pverp) = stebol*(delt (i)*absnxt(i,pver,2) + delt1(i)*absnxt(i,pver,3)) end do do k=ntoplw,pver-1 do i=1,ncol bk2(i) = (abstot(i,k,pver) + abstot(i,k,pverp))*0.5 bk1(i) = bk2(i) s(i,k,pverp) = stebol*(bk2(i)*delt(i) + bk1(i)*delt1(i)) end do end do do km=pver,ntoplw+1,-1 do i=1,ncol delt(i) = tint4(i,km-1) - tlayr4(i,km) delt1(i) = tlayr4(i,km) - tint4(i,km) end do do k=pverp,ntoplw,-1 if (k == km) then do i=1,ncol bk2(i) = absnxt(i,km-1,4) bk1(i) = absnxt(i,km-1,1) end do else if (k == km-1) then do i=1,ncol bk2(i) = absnxt(i,km-1,2) bk1(i) = absnxt(i,km-1,3) end do else do i=1,ncol bk2(i) = (abstot(i,k,km-1) + abstot(i,k,km))*0.5 bk1(i) = bk2(i) end do end if do i=1,ncol s(i,k,km) = s(i,k,km+1) + stebol*(bk2(i)*delt(i) + bk1(i)*delt1(i)) end do end do end do do i=1,ncol fsul(i,pverp) = lwupcgs(i) end do do i=1,ncol tmp(i) = fsul(i,pverp) - stebol*tint4(i,pverp) fsul(i,ntoplw) = fsul(i,pverp) - abstot(i,ntoplw,pverp)*tmp(i) + s(i,ntoplw,ntoplw+1) fsdl(i,ntoplw) = stebol*(tplnke(i)**4)*emstot(i,ntoplw) end do do k=ntoplw+1,pver do i=1,ncol fsul(i,k) = fsul(i,pverp) - abstot(i,k,pverp)*tmp(i) + s(i,k,k+1) fsdl(i,k) = stebol*(tplnke(i)**4)*emstot(i,k) - (s(i,k,ntoplw+1) - s(i,k,k+1)) end do end do do i=1,ncol absbt(i) = stebol*(tplnke(i)**4)*emstot(i,pverp) fsdl(i,pverp) = absbt(i) - s(i,pverp,ntoplw+1) end do cldp(:ncol,ntoplw:pver) = cld(:ncol,ntoplw:pver) cldp(:ncol,pverp) = 0.0 maxcld(1:ncol) = maxval(cldp(1:ncol,ntoplw:pver),dim=2) npts = 0 do i=1,ncol if (maxcld(i) < cldmin) then npts = npts + 1 indx(npts) = i end if end do do ii = 1, npts i = indx(ii) do k = ntoplw, pverp fdl(i,k) = fsdl(i,k) ful(i,k) = fsul(i,k) end do end do npts = 0 do i=1,ncol if (maxcld(i) >= cldmin) then npts = npts + 1 indx(npts) = i end if end do do ii = 1, npts i = indx(ii) fdl(i,ntoplw) = fsdl(i,ntoplw) fdl(i,pverp) = 0.0 ful(i,ntoplw) = 0.0 ful(i,pverp) = fsul(i,pverp) do k = ntoplw+1, pver fdl(i,k) = 0.0 ful(i,k) = 0.0 end do do k = ntoplw, pver fclt4(i,k-1) = stebol*tint4(i,k) fclb4(i,k-1) = stebol*tint4(i,k+1) enddo fclb4(i,ntoplw-2) = stebol*tint4(i,ntoplw) fclt4(i,pver) = stebol*tint4(i,pverp) do irgn = 0, nmxrgn(i) kx2(i,irgn) = ntoplw-1 end do nrgn(i) = 0 end do do ii = 1, npts ilon = indx(ii) do irgn = 1, nmxrgn(ilon) n = 0 if (kx2(ilon,irgn-1) < pver) then nrgn(ilon) = irgn k1 = kx2(ilon,irgn-1)+1 kx1(ilon,irgn) = k1 kx2(ilon,irgn) = 0 do k2 = pver, k1, -1 if (pmid(ilon,k2) <= pmxrgn(ilon,irgn)) then kx2(ilon,irgn) = k2 exit end if end do do k = k1, k2 if (cldp(ilon,k) >= cldmin) then n = n+1 indxmx(n,irgn) = ilon exit endif end do endif ncolmx(irgn) = n do iimx = 1, ncolmx(irgn) i = indxmx(iimx,irgn) n = 0 do k = kx1(i,irgn),kx2(i,irgn) if (cldp(i,k) >= cldmin) then n = n+1 ksort(n) = k asort(n) = 1.0-cldp(i,k) end if end do nxs(i,irgn) = n if (nxs(i,irgn) == 2) then if (asort(2) < asort(1)) then ktmp = ksort(1) ksort(1) = ksort(2) ksort(2) = ktmp atmp = asort(1) asort(1) = asort(2) asort(2) = atmp endif else if (nxs(i,irgn) >= 3) then call sortarray(nxs(i,irgn),asort,ksort(1:)) endif do l = 1, nxs(i,irgn) kxs(l,i,irgn) = ksort(l) end do end do end do do irgn = 1, nmxrgn(ilon) iimx = 1 if (ilon < indxmx(iimx,irgn) .and. irgn <= nrgn(ilon)) then k1 = kx1(ilon,irgn) do km1 = ntoplw-2, k1-2 km4 = km1+3 k2 = k1 k3 = k2+1 tmp(ilon) = s(ilon,k2,min(k3,pverp))*min(1,pverp2-k3) emx0 = (fdl(ilon,k1)-fsdl(ilon,k1))/ & ((fclb4(ilon,km1)-s(ilon,k2,km4)+tmp(ilon))- fsdl(ilon,k1)) if (emx0 >= 0.0 .and. emx0 <= 1.0) exit end do km1 = min(km1,k1-2) do k2 = kx1(ilon,irgn)+1, kx2(ilon,irgn)+1 k3 = k2+1 tmp(ilon) = s(ilon,k2,min(k3,pverp))*min(1,pverp2-k3) fdl(ilon,k2) = (1.0-emx0)*fsdl(ilon,k2) + & emx0*(fclb4(ilon,km1)-s(ilon,k2,km4)+tmp(ilon)) end do else if (ilon==indxmx(iimx,irgn) .and. iimx<=ncolmx(irgn)) then iimx = iimx+1 end if do iimx = 1, ncolmx(irgn) i = indxmx(iimx,irgn) k1 = kx1(i,irgn) do km1 = ntoplw-2,k1-2 km4 = km1+3 k2 = k1 k3 = k2 + 1 tmp(i) = s(i,k2,min(k3,pverp))*min(1,pverp2-k3) tmp2(i) = s(i,k2,min(km4,pverp))*min(1,pverp2-km4) emx0 = (fdl(i,k1)-fsdl(i,k1))/((fclb4(i,km1)-tmp2(i)+tmp(i))-fsdl(i,k1)) if (emx0 >= 0.0 .and. emx0 <= 1.0) exit end do km1 = min(km1,k1-2) ksort(0) = km1 + 1 nxsk = 0 do k = kx1(i,irgn), kx2(i,irgn) if (nxsk < nxs(i,irgn)) then nxsk = 0 do l = 1, nxs(i,irgn) k1 = kxs(l,i,irgn) if (k >= k1) then nxsk = nxsk + 1 ksort(nxsk) = k1 endif end do endif ksort(nxsk+1) = pverp do l = 1, nxsk emx(l) = emis(i,ksort(l)) end do emx(0) = emx0 cld0 = 1.0 k2 = k+1 k3 = k2+1 tmp(i) = s(i,k2,min(k3,pverp))*min(1,pverp2-k3) do l = 1, nxsk+1 cld1 = cldp(i,ksort(l))*min(1,nxsk+1-l) if (cld0 /= cld1) then fdl(i,k2) = fdl(i,k2)+(cld0-cld1)*fsdl(i,k2) do l1 = 0, l - 1 km1 = ksort(l1)-1 km4 = km1+3 tmp2(i) = s(i,k2,min(km4,pverp))* min(1,pverp2-km4) fdl(i,k2) = fdl(i,k2)+(cld0-cld1)*emx(l1)*(fclb4(i,km1)-tmp2(i)+tmp(i)- & fsdl(i,k2)) end do endif cld0 = cld1 if (l <= nxsk) then k1 = ksort(l) trans = 1.0-emis(i,k1) do l1 = 0, nxsk if (ksort(l1) < k1) then emx(l1) = emx(l1)*trans endif end do end if end do end do end do end do do irgn = nmxrgn(ilon), 1, -1 iimx = 1 if (ilon < indxmx(iimx,irgn) .and. irgn <= nrgn(ilon)) then k1 = kx2(ilon,irgn)+1 if (k1 < pverp) then do km1 = pver-1,kx2(ilon,irgn),-1 km3 = km1+2 k2 = k1 k3 = k2+1 tmp(ilon) = s(ilon,k2,min(km3,pverp))* min(1,pverp2-km3) emx0 = (ful(ilon,k1)-fsul(ilon,k1))/ & ((fclt4(ilon,km1)+s(ilon,k2,k3)-tmp(ilon))- fsul(ilon,k1)) if (emx0 >= 0.0 .and. emx0 <= 1.0) exit end do km1 = max(km1,kx2(ilon,irgn)) else km1 = k1-1 km3 = km1+2 emx0 = 1.0 endif do k2 = kx1(ilon,irgn), kx2(ilon,irgn) k3 = k2+1 tmp(ilon) = s(ilon,k2,min(km3,pverp))* min(1,pverp2-km3) ful(ilon,k2) =(1.0-emx0)*fsul(ilon,k2) + emx0* & (fclt4(ilon,km1)+s(ilon,k2,k3)-tmp(ilon)) end do else if (ilon==indxmx(iimx,irgn) .and. iimx<=ncolmx(irgn)) then iimx = iimx+1 end if do iimx = 1, ncolmx(irgn) i = indxmx(iimx,irgn) k1 = kx2(i,irgn)+1 if (k1 < pverp) then do km1 = pver-1,kx2(i,irgn),-1 km3 = km1+2 k2 = k1 k3 = k2+1 tmp(i) = s(i,k2,min(km3,pverp))*min(1,pverp2-km3) emx0 = (ful(i,k1)-fsul(i,k1))/((fclt4(i,km1)+s(i,k2,k3)-tmp(i))-fsul(i,k1)) if (emx0 >= 0.0 .and. emx0 <= 1.0) exit end do km1 = max(km1,kx2(i,irgn)) else emx0 = 1.0 km1 = k1-1 endif ksort(0) = km1 + 1 nxsk = 0 do k = kx2(i,irgn), kx1(i,irgn), -1 if (nxsk < nxs(i,irgn)) then nxsk = 0 do l = 1, nxs(i,irgn) k1 = kxs(l,i,irgn) if (k <= k1) then nxsk = nxsk + 1 ksort(nxsk) = k1 endif end do endif ksort(nxsk+1) = pverp do l = 1, nxsk emx(l) = emis(i,ksort(l)) end do emx(0) = emx0 cld0 = 1.0 k2 = k k3 = k2+1 do l = 1, nxsk+1 cld1 = cldp(i,ksort(l))*min(1,nxsk+1-l) if (cld0 /= cld1) then ful(i,k2) = ful(i,k2)+(cld0-cld1)*fsul(i,k2) do l1 = 0, l - 1 km1 = ksort(l1)-1 km3 = km1+2 tmp(i) = s(i,k2,min(km3,pverp))* min(1,pverp2-km3) ful(i,k2) = ful(i,k2)+(cld0-cld1)*emx(l1)* & (fclt4(i,km1)+s(i,k2,k3)-tmp(i)- fsul(i,k2)) end do endif cld0 = cld1 if (l <= nxsk) then k1 = ksort(l) trans = 1.0-emis(i,k1) do l1 = 0, nxsk if (ksort(l1) > k1) then emx(l1) = emx(l1)*trans endif end do end if end do end do end do end do end do do i=1,ncol flwds(i) = fdl (i,pverp ) flns(i) = ful (i,pverp ) - fdl (i,pverp ) flnsc(i) = fsul(i,pverp ) - fsdl(i,pverp ) flnt(i) = ful (i,ntoplw) - fdl (i,ntoplw) flntc(i) = fsul(i,ntoplw) - fsdl(i,ntoplw) flut(i) = ful (i,ntoplw) flutc(i) = fsul(i,ntoplw) end do do k=ntoplw,pver do i=1,ncol qrl(i,k) = (ful(i,k) - fdl(i,k) - ful(i,k+1) + fdl(i,k+1))* & 1.E-4*gravit/((pint(i,k) - pint(i,k+1))) end do end do if ( ntoplw > 1 )then qrl(:ncol,:ntoplw-1) = 0. end if do k=ntoplw,pverp do i=1,ncol flup(i,k) = ful(i,k) flupc(i,k) = fsul(i,k) fldn(i,k) = fdl(i,k) fldnc(i,k) = fsdl(i,k) end do end do if ( ntoplw > 1 )then flup(:ncol,:ntoplw-1) = 0. flupc(:ncol,:ntoplw-1) = 0. fldn(:ncol,:ntoplw-1) = 0. fldnc(:ncol,:ntoplw-1) = 0. end if return end subroutine radclwmx subroutine radcswmx(jj, lchnk ,ncol ,pcols, pver, pverp, & pint ,pmid ,h2ommr ,rh ,o3mmr , & aermmr ,cld ,cicewp ,cliqwp ,rel , & rei ,tauxcl ,tauxci ,eccf ,coszrs ,scon ,solin ,solcon, & asdir ,asdif ,aldir ,aldif ,nmxrgn , & pmxrgn ,qrs ,fsnt ,fsntc ,fsntoa , & fsntoac ,fsnirtoa,fsnrtoac,fsnrtoaq,fsns , & fsnsc ,fsdsc ,fsds ,sols ,soll , & solsd ,solld ,frc_day , & fsup ,fsupc ,fsdn ,fsdnc , & fsdndir ,fsdncdir,fsdndif ,fsdncdif, & fsdsdir ,fsdsdif ,fsdscdir,fsdscdif, & aertau ,aerssa ,aerasm ,aerfwd ) implicit none integer nspint integer naer_groups parameter ( nspint = 19 ) parameter ( naer_groups = 7 ) real(r8) cldmin parameter (cldmin = 1.0e-80_r8) real(r8) areamin parameter (areamin = 0.01_r8) real(r8) cldeps parameter (cldeps = 0.0_r8) integer nconfgmax parameter (nconfgmax = 15) integer, intent(in) :: lchnk,jj integer, intent(in) :: pcols, pver, pverp integer, intent(in) :: ncol real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: h2ommr(pcols,pver) real(r8), intent(in) :: o3mmr(pcols,pver) real(r8), intent(in) :: aermmr(pcols,pver,naer_all) real(r8), intent(in) :: rh(pcols,pver) real(r8), intent(in) :: cld(pcols,pver) real(r8), intent(in) :: cicewp(pcols,pver) real(r8), intent(in) :: cliqwp(pcols,pver) real(r8), intent(in) :: rel(pcols,pver) real(r8), intent(in) :: rei(pcols,pver) real(r8), intent(in) :: eccf real, intent(in) :: solcon real(r8), intent(in) :: coszrs(pcols) real(r8), intent(in) :: asdir(pcols) real(r8), intent(in) :: aldir(pcols) real(r8), intent(in) :: asdif(pcols) real(r8), intent(in) :: aldif(pcols) real(r8), intent(in) :: scon real(r8), intent(inout) :: pmxrgn(pcols,pverp) integer, intent(inout) :: nmxrgn(pcols) real(r8), intent(out) :: solin(pcols) real(r8), intent(out) :: qrs(pcols,pver) real(r8), intent(out) :: fsns(pcols) real(r8), intent(out) :: fsnt(pcols) real(r8), intent(out) :: fsntoa(pcols) real(r8), intent(out) :: fsds(pcols) real(r8), intent(out) :: fsnsc(pcols) real(r8), intent(out) :: fsdsc(pcols) real(r8), intent(out) :: fsdscdir(pcols) real(r8), intent(out) :: fsdscdif(pcols) real(r8), intent(out) :: fsdsdir(pcols) real(r8), intent(out) :: fsdsdif(pcols) real(r8), intent(out) :: fsntc(pcols) real(r8), intent(out) :: fsntoac(pcols) real(r8), intent(out) :: sols(pcols) real(r8), intent(out) :: soll(pcols) real(r8), intent(out) :: solsd(pcols) real(r8), intent(out) :: solld(pcols) real(r8), intent(out) :: fsnirtoa(pcols) real(r8), intent(out) :: fsnrtoac(pcols) real(r8), intent(out) :: fsnrtoaq(pcols) real(r8), intent(out) :: tauxcl(pcols,0:pver) real(r8), intent(out) :: tauxci(pcols,0:pver) real(r8), intent(out) :: fsup(pcols,pverp) real(r8), intent(out) :: fsupc(pcols,pverp) real(r8), intent(out) :: fsdn(pcols,pverp) real(r8), intent(out) :: fsdnc(pcols,pverp) real(r8), intent(out) :: fsdndir(pcols,pverp) real(r8), intent(out) :: fsdncdir(pcols,pverp) real(r8), intent(out) :: fsdndif(pcols,pverp) real(r8), intent(out) :: fsdncdif(pcols,pverp) real(r8) , intent(out) :: frc_day(pcols) real(r8) :: aertau(pcols,nspint,naer_groups) real(r8) :: aerssa(pcols,nspint,naer_groups) real(r8) :: aerasm(pcols,nspint,naer_groups) real(r8) :: aerfwd(pcols,nspint,naer_groups) real(r8) asort(pverp) real(r8) atmp real(r8) cld0 real(r8) totwgt real(r8) wgtv(nconfgmax) real(r8) wstr(pverp,pverp) real(r8) xexpt real(r8) xrdnd real(r8) xrupd real(r8) xrups real(r8) xtdnt real(r8) xwgt real(r8) yexpt real(r8) yrdnd real(r8) yrupd real(r8) ytdnd real(r8) ytupd real(r8) zexpt real(r8) zrdnd real(r8) zrupd real(r8) zrups real(r8) ztdnt logical new_term logical region_found integer ccon(0:pverp,nconfgmax) integer cstr(0:pverp,pverp) integer icond(0:pverp,nconfgmax) integer iconu(0:pverp,nconfgmax) integer iconfig integer irgn integer is0 integer is1 integer isn integer istr(pverp+1) integer istrtd(0:pverp,0:nconfgmax+1) integer istrtu(0:pverp,0:nconfgmax+1) integer j integer k1 integer k2 integer ksort(pverp) integer ktmp integer kx1(0:pverp) integer kx2(0:pverp) integer l integer l0 integer mrgn integer mstr integer n0 integer n1 integer nconfig integer nconfigm integer npasses integer nrgn integer nstr(pverp) integer nuniq integer nuniqd(0:pverp) integer nuniqu(0:pverp) integer nxs integer ptr0(nconfgmax) integer ptr1(nconfgmax) integer ptrc(nconfgmax) integer ns integer i integer k integer km1 integer kp1 integer n integer ndayc integer idayc(pcols) integer indxsl integer ksz integer krh integer kaer real(r8) wrh real(r8) :: rhtrunc real(r8) albdir(pcols,nspint) real(r8) albdif(pcols,nspint) real(r8) wgtint real(r8) solflx real(r8) sfltot real(r8) totfld(0:pver) real(r8) fswup(0:pverp) real(r8) fswdn(0:pverp) real(r8) fswupc(0:pverp) real(r8) fswdnc(0:pverp) real(r8) fswdndir(0:pverp) real(r8) fswdncdir(0:pverp) real(r8) wcl(pcols,0:pver) real(r8) gcl(pcols,0:pver) real(r8) fcl(pcols,0:pver) real(r8) wci(pcols,0:pver) real(r8) gci(pcols,0:pver) real(r8) fci(pcols,0:pver) real(r8) usul(pcols,pver) real(r8) ubg(pcols,pver) real(r8) usslt(pcols,pver) real(r8) ucphil(pcols,pver) real(r8) ucphob(pcols,pver) real(r8) ucb(pcols,pver) real(r8) uvolc(pcols,pver) real(r8) udst(ndstsz,pcols,pver) real(r8) tau_sul real(r8) tau_bg real(r8) tau_sslt real(r8) tau_cphil real(r8) tau_cphob real(r8) tau_cb real(r8) tau_volc real(r8) tau_dst(ndstsz) real(r8) tau_dst_tot real(r8) tau_tot real(r8) tau_w_sul real(r8) tau_w_bg real(r8) tau_w_sslt real(r8) tau_w_cphil real(r8) tau_w_cphob real(r8) tau_w_cb real(r8) tau_w_volc real(r8) tau_w_dst(ndstsz) real(r8) tau_w_dst_tot real(r8) tau_w_tot real(r8) tau_w_g_sul real(r8) tau_w_g_bg real(r8) tau_w_g_sslt real(r8) tau_w_g_cphil real(r8) tau_w_g_cphob real(r8) tau_w_g_cb real(r8) tau_w_g_volc real(r8) tau_w_g_dst(ndstsz) real(r8) tau_w_g_dst_tot real(r8) tau_w_g_tot real(r8) f_sul real(r8) f_bg real(r8) f_sslt real(r8) f_cphil real(r8) f_cphob real(r8) f_cb real(r8) f_volc real(r8) f_dst(ndstsz) real(r8) f_dst_tot real(r8) f_tot real(r8) tau_w_f_sul real(r8) tau_w_f_bg real(r8) tau_w_f_sslt real(r8) tau_w_f_cphil real(r8) tau_w_f_cphob real(r8) tau_w_f_cb real(r8) tau_w_f_volc real(r8) tau_w_f_dst(ndstsz) real(r8) tau_w_f_dst_tot real(r8) tau_w_f_tot real(r8) w_dst_tot real(r8) w_tot real(r8) g_dst_tot real(r8) g_tot real(r8) ksuli real(r8) ksslti real(r8) kcphili real(r8) wsuli real(r8) wsslti real(r8) wcphili real(r8) gsuli real(r8) gsslti real(r8) gcphili real(r8) tauxar(pcols,0:pver) real(r8) wa(pcols,0:pver) real(r8) ga(pcols,0:pver) real(r8) fa(pcols,0:pver) real(r8) pflx(pcols,0:pverp) real(r8) zenfac(pcols) real(r8) sqrco2 real(r8) tmp1 real(r8) tmp2 real(r8) pdel real(r8) path real(r8) ptop real(r8) ptho2 real(r8) ptho3 real(r8) pthco2 real(r8) pthh2o real(r8) h2ostr real(r8) wavmid(nspint) real(r8) trayoslp real(r8) tmp1l real(r8) tmp2l real(r8) tmp3l real(r8) tmp1i real(r8) tmp2i real(r8) tmp3i real(r8) rdenom real(r8) rdirexp real(r8) tdnmexp real(r8) psf(nspint) real(r8) uh2o(pcols,0:pver) real(r8) uo3(pcols,0:pver) real(r8) uco2(pcols,0:pver) real(r8) uo2(pcols,0:pver) real(r8) uaer(pcols,0:pver) real(r8) uth2o(pcols) real(r8) uto3(pcols) real(r8) utco2(pcols) real(r8) uto2(pcols) real(r8) rdir(nspint,pcols,0:pver) real(r8) rdif(nspint,pcols,0:pver) real(r8) tdir(nspint,pcols,0:pver) real(r8) tdif(nspint,pcols,0:pver) real(r8) explay(nspint,pcols,0:pver) real(r8) rdirc(nspint,pcols,0:pver) real(r8) rdifc(nspint,pcols,0:pver) real(r8) tdirc(nspint,pcols,0:pver) real(r8) tdifc(nspint,pcols,0:pver) real(r8) explayc(nspint,pcols,0:pver) real(r8) flxdiv real(r8) exptdn(0:pverp,nconfgmax) real(r8) rdndif(0:pverp,nconfgmax) real(r8) rupdif(0:pverp,nconfgmax) real(r8) rupdir(0:pverp,nconfgmax) real(r8) tdntot(0:pverp,nconfgmax) real(r8) exptdnc(0:pverp) real(r8) rdndifc(0:pverp) real(r8) rupdifc(0:pverp) real(r8) rupdirc(0:pverp) real(r8) tdntotc(0:pverp) real(r8) fluxup(0:pverp) real(r8) fluxdn(0:pverp) real(r8) fluxdndir(0:pverp) real(r8) wexptdn real(r8) abarli real(r8) bbarli real(r8) cbarli real(r8) dbarli real(r8) ebarli real(r8) fbarli real(r8) abarii real(r8) bbarii real(r8) cbarii real(r8) dbarii real(r8) ebarii real(r8) fbarii do i=1, ncol fsds(i) = 0.0_r8 fsnirtoa(i) = 0.0_r8 fsnrtoac(i) = 0.0_r8 fsnrtoaq(i) = 0.0_r8 fsns(i) = 0.0_r8 fsnsc(i) = 0.0_r8 fsdsc(i) = 0.0_r8 fsdscdir(i) = 0.0_r8 fsdscdif(i) = 0.0_r8 fsdsdir(i) = 0.0_r8 fsdsdif(i) = 0.0_r8 fsnt(i) = 0.0_r8 fsntc(i) = 0.0_r8 fsntoa(i) = 0.0_r8 fsntoac(i) = 0.0_r8 solin(i) = 0.0_r8 sols(i) = 0.0_r8 soll(i) = 0.0_r8 solsd(i) = 0.0_r8 solld(i) = 0.0_r8 do k=1,pverp fsup(i,k) = 0.0_r8 fsupc(i,k) = 0.0_r8 fsdn(i,k) = 0.0_r8 fsdnc(i,k) = 0.0_r8 fsdndir(i,k) = 0.0_r8 fsdncdir(i,k) = 0.0_r8 fsdndif(i,k) = 0.0_r8 fsdncdif(i,k) = 0.0_r8 tauxcl(i,k-1) = 0.0_r8 tauxci(i,k-1) = 0.0_r8 end do do k=1, pver qrs(i,k) = 0.0_r8 end do do kaer = 1, naer_groups do ns = 1, nspint frc_day(i) = 0.0_r8 aertau(i,ns,kaer) = 0.0_r8 aerssa(i,ns,kaer) = 0.0_r8 aerasm(i,ns,kaer) = 0.0_r8 aerfwd(i,ns,kaer) = 0.0_r8 end do end do end do ndayc = 0 do i=1,ncol if (coszrs(i) > 0.0_r8) then ndayc = ndayc + 1 idayc(ndayc) = i end if end do if (ndayc == 0) return tmp1 = 0.5_r8/(gravit*sslp) tmp2 = delta/gravit sqrco2 = sqrt(co2mmr) do n=1,ndayc i=idayc(n) solin(i) = solcon*coszrs(i)*1000. pflx(i,0) = 0._r8 do k=1,pverp pflx(i,k) = pint(i,k) end do ptop = pflx(i,1) ptho2 = o2mmr * ptop / gravit ptho3 = o3mmr(i,1) * ptop / gravit pthco2 = sqrco2 * (ptop / gravit) h2ostr = sqrt( 1._r8 / h2ommr(i,1) ) zenfac(i) = sqrt(coszrs(i)) pthh2o = ptop**2*tmp1 + (ptop*rga)* & (h2ostr*zenfac(i)*delta) uh2o(i,0) = h2ommr(i,1)*pthh2o uco2(i,0) = zenfac(i)*pthco2 uo2 (i,0) = zenfac(i)*ptho2 uo3 (i,0) = ptho3 uaer(i,0) = 0.0_r8 do k=1,pver pdel = pflx(i,k+1) - pflx(i,k) path = pdel / gravit ptho2 = o2mmr * path ptho3 = o3mmr(i,k) * path pthco2 = sqrco2 * path h2ostr = sqrt(1.0_r8/h2ommr(i,k)) pthh2o = (pflx(i,k+1)**2 - pflx(i,k)**2)*tmp1 + pdel*h2ostr*zenfac(i)*tmp2 uh2o(i,k) = h2ommr(i,k)*pthh2o uco2(i,k) = zenfac(i)*pthco2 uo2 (i,k) = zenfac(i)*ptho2 uo3 (i,k) = ptho3 usul(i,k) = aermmr(i,k,idxSUL) * path ubg(i,k) = aermmr(i,k,idxBG) * path usslt(i,k) = aermmr(i,k,idxSSLT) * path if (usslt(i,k) .lt. 0.0) then usslt(i,k) = 0.0 end if ucphil(i,k) = aermmr(i,k,idxOCPHI) * path ucphob(i,k) = aermmr(i,k,idxOCPHO) * path ucb(i,k) = ( aermmr(i,k,idxBCPHO) + aermmr(i,k,idxBCPHI) ) * path uvolc(i,k) = aermmr(i,k,idxVOLC) do ksz = 1, ndstsz udst(ksz,i,k) = aermmr(i,k,idxDUSTfirst-1+ksz) * path end do end do uth2o(i) = 0.0_r8 uto3(i) = 0.0_r8 utco2(i) = 0.0_r8 uto2(i) = 0.0_r8 do k=1,pver uth2o(i) = uth2o(i) + uh2o(i,k) uto3(i) = uto3(i) + uo3(i,k) utco2(i) = utco2(i) + uco2(i,k) uto2(i) = uto2(i) + uo2(i,k) end do tauxcl(i,0) = 0._r8 wcl(i,0) = 0.999999_r8 gcl(i,0) = 0.85_r8 fcl(i,0) = 0.725_r8 tauxci(i,0) = 0._r8 wci(i,0) = 0.999999_r8 gci(i,0) = 0.85_r8 fci(i,0) = 0.725_r8 tauxar(i,0) = 0._r8 wa(i,0) = 0.925_r8 ga(i,0) = 0.850_r8 fa(i,0) = 0.7225_r8 end do do ns=1,nspint if(wavmax(ns) <= 0.7_r8) then indxsl = 1 else if(wavmin(ns) == 0.700_r8) then indxsl = 2 else if(wavmin(ns) == 0.701_r8) then indxsl = 3 else if(wavmin(ns) == 0.702_r8 .or. wavmin(ns) > 2.38_r8) then indxsl = 4 end if abarli = abarl(indxsl) bbarli = bbarl(indxsl) cbarli = cbarl(indxsl) dbarli = dbarl(indxsl) ebarli = ebarl(indxsl) fbarli = fbarl(indxsl) abarii = abari(indxsl) bbarii = bbari(indxsl) cbarii = cbari(indxsl) dbarii = dbari(indxsl) ebarii = ebari(indxsl) fbarii = fbari(indxsl) psf(ns) = 1.0_r8 if(ph2o(ns)/=0._r8) psf(ns) = psf(ns)*ph2o(ns) if(pco2(ns)/=0._r8) psf(ns) = psf(ns)*pco2(ns) if(po2 (ns)/=0._r8) psf(ns) = psf(ns)*po2 (ns) do n=1,ndayc i=idayc(n) frc_day(i) = 1.0_r8 do kaer = 1, naer_groups aertau(i,ns,kaer) = 0.0 aerssa(i,ns,kaer) = 0.0 aerasm(i,ns,kaer) = 0.0 aerfwd(i,ns,kaer) = 0.0 end do do k=1,pver tmp1l = abarli + bbarli/rel(i,k) tmp2l = 1._r8 - cbarli - dbarli*rel(i,k) tmp3l = fbarli*rel(i,k) tmp1i = abarii + bbarii/rei(i,k) tmp2i = 1._r8 - cbarii - dbarii*rei(i,k) tmp3i = fbarii*rei(i,k) if (cld(i,k) >= cldmin .and. cld(i,k) >= cldeps) then tauxcl(i,k) = cliqwp(i,k)*tmp1l tauxci(i,k) = cicewp(i,k)*tmp1i else tauxcl(i,k) = 0.0 tauxci(i,k) = 0.0 endif wcl(i,k) = min(tmp2l,.999999_r8) gcl(i,k) = ebarli + tmp3l fcl(i,k) = gcl(i,k)*gcl(i,k) wci(i,k) = min(tmp2i,.999999_r8) gci(i,k) = ebarii + tmp3i fci(i,k) = gci(i,k)*gci(i,k) rhtrunc = rh(i,k) rhtrunc = min(rh(i,k),1._r8) krh = min(floor( rhtrunc * nrh ) + 1, nrh - 1) wrh = rhtrunc * nrh - krh ksuli = ksul(krh + 1, ns) * (wrh + 1) - ksul(krh, ns) * wrh ksslti = ksslt(krh + 1, ns) * (wrh + 1) - ksslt(krh, ns) * wrh kcphili = kcphil(krh + 1, ns) * (wrh + 1) - kcphil(krh, ns) * wrh wsuli = wsul(krh + 1, ns) * (wrh + 1) - wsul(krh, ns) * wrh wsslti = wsslt(krh + 1, ns) * (wrh + 1) - wsslt(krh, ns) * wrh wcphili = wcphil(krh + 1, ns) * (wrh + 1) - wcphil(krh, ns) * wrh gsuli = gsul(krh + 1, ns) * (wrh + 1) - gsul(krh, ns) * wrh gsslti = gsslt(krh + 1, ns) * (wrh + 1) - gsslt(krh, ns) * wrh gcphili = gcphil(krh + 1, ns) * (wrh + 1) - gcphil(krh, ns) * wrh tau_sul = 1.e4 * ksuli * usul(i,k) tau_sslt = 1.e4 * ksslti * usslt(i,k) tau_cphil = 1.e4 * kcphili * ucphil(i,k) tau_cphob = 1.e4 * kcphob(ns) * ucphob(i,k) tau_cb = 1.e4 * kcb(ns) * ucb(i,k) tau_volc = 1.e3 * kvolc(ns) * uvolc(i,k) tau_dst(:) = 1.e4 * kdst(:,ns) * udst(:,i,k) tau_bg = 1.e4 * kbg(ns) * ubg(i,k) tau_w_sul = tau_sul * wsuli tau_w_sslt = tau_sslt * wsslti tau_w_cphil = tau_cphil * wcphili tau_w_cphob = tau_cphob * wcphob(ns) tau_w_cb = tau_cb * wcb(ns) tau_w_volc = tau_volc * wvolc(ns) tau_w_dst(:) = tau_dst(:) * wdst(:,ns) tau_w_bg = tau_bg * wbg(ns) tau_w_g_sul = tau_w_sul * gsuli tau_w_g_sslt = tau_w_sslt * gsslti tau_w_g_cphil = tau_w_cphil * gcphili tau_w_g_cphob = tau_w_cphob * gcphob(ns) tau_w_g_cb = tau_w_cb * gcb(ns) tau_w_g_volc = tau_w_volc * gvolc(ns) tau_w_g_dst(:) = tau_w_dst(:) * gdst(:,ns) tau_w_g_bg = tau_w_bg * gbg(ns) f_sul = gsuli * gsuli f_sslt = gsslti * gsslti f_cphil = gcphili * gcphili f_cphob = gcphob(ns) * gcphob(ns) f_cb = gcb(ns) * gcb(ns) f_volc = gvolc(ns) * gvolc(ns) f_dst(:) = gdst(:,ns) * gdst(:,ns) f_bg = gbg(ns) * gbg(ns) tau_w_f_sul = tau_w_sul * f_sul tau_w_f_bg = tau_w_bg * f_bg tau_w_f_sslt = tau_w_sslt * f_sslt tau_w_f_cphil = tau_w_cphil * f_cphil tau_w_f_cphob = tau_w_cphob * f_cphob tau_w_f_cb = tau_w_cb * f_cb tau_w_f_volc = tau_w_volc * f_volc tau_w_f_dst(:) = tau_w_dst(:) * f_dst(:) tau_dst_tot = sum(tau_dst) tau_w_dst_tot = sum(tau_w_dst) tau_w_g_dst_tot = sum(tau_w_g_dst) tau_w_f_dst_tot = sum(tau_w_f_dst) if (tau_dst_tot .gt. 0.0) then w_dst_tot = tau_w_dst_tot / tau_dst_tot else w_dst_tot = 0.0 endif if (tau_w_dst_tot .gt. 0.0) then g_dst_tot = tau_w_g_dst_tot / tau_w_dst_tot f_dst_tot = tau_w_f_dst_tot / tau_w_dst_tot else g_dst_tot = 0.0 f_dst_tot = 0.0 endif tau_tot = tau_sul + tau_sslt & + tau_cphil + tau_cphob + tau_cb + tau_dst_tot tau_tot = tau_tot + tau_bg + tau_volc tau_w_tot = tau_w_sul + tau_w_sslt & + tau_w_cphil + tau_w_cphob + tau_w_cb + tau_w_dst_tot tau_w_tot = tau_w_tot + tau_w_bg + tau_w_volc tau_w_g_tot = tau_w_g_sul + tau_w_g_sslt & + tau_w_g_cphil + tau_w_g_cphob + tau_w_g_cb + tau_w_g_dst_tot tau_w_g_tot = tau_w_g_tot + tau_w_g_bg + tau_w_g_volc tau_w_f_tot = tau_w_f_sul + tau_w_f_sslt & + tau_w_f_cphil + tau_w_f_cphob + tau_w_f_cb + tau_w_f_dst_tot tau_w_f_tot = tau_w_f_tot + tau_w_f_bg + tau_w_f_volc if (tau_tot .gt. 0.0) then w_tot = tau_w_tot / tau_tot else w_tot = 0.0 endif if (tau_w_tot .gt. 0.0) then g_tot = tau_w_g_tot / tau_w_tot f_tot = tau_w_f_tot / tau_w_tot else g_tot = 0.0 f_tot = 0.0 endif tauxar(i,k) = tau_tot wa(i,k) = min(w_tot, 0.999999_r8) if (g_tot.gt.1._r8) write(6,*) "g_tot > 1" if (g_tot.lt.-1._r8) write(6,*) "g_tot < -1" ga(i,k) = g_tot if (f_tot.gt.1._r8) write(6,*)"f_tot > 1" if (f_tot.lt.0._r8) write(6,*)"f_tot < 0" fa(i,k) = f_tot aertau(i,ns,1) = aertau(i,ns,1) + tau_sul aertau(i,ns,2) = aertau(i,ns,2) + tau_sslt aertau(i,ns,3) = aertau(i,ns,3) + tau_cphil + tau_cphob + tau_cb aertau(i,ns,4) = aertau(i,ns,4) + tau_dst_tot aertau(i,ns,5) = aertau(i,ns,5) + tau_bg aertau(i,ns,6) = aertau(i,ns,6) + tau_volc aertau(i,ns,7) = aertau(i,ns,7) + tau_tot aerssa(i,ns,1) = aerssa(i,ns,1) + tau_w_sul aerssa(i,ns,2) = aerssa(i,ns,2) + tau_w_sslt aerssa(i,ns,3) = aerssa(i,ns,3) + tau_w_cphil + tau_w_cphob + tau_w_cb aerssa(i,ns,4) = aerssa(i,ns,4) + tau_w_dst_tot aerssa(i,ns,5) = aerssa(i,ns,5) + tau_w_bg aerssa(i,ns,6) = aerssa(i,ns,6) + tau_w_volc aerssa(i,ns,7) = aerssa(i,ns,7) + tau_w_tot aerasm(i,ns,1) = aerasm(i,ns,1) + tau_w_g_sul aerasm(i,ns,2) = aerasm(i,ns,2) + tau_w_g_sslt aerasm(i,ns,3) = aerasm(i,ns,3) + tau_w_g_cphil + tau_w_g_cphob + tau_w_g_cb aerasm(i,ns,4) = aerasm(i,ns,4) + tau_w_g_dst_tot aerasm(i,ns,5) = aerasm(i,ns,5) + tau_w_g_bg aerasm(i,ns,6) = aerasm(i,ns,6) + tau_w_g_volc aerasm(i,ns,7) = aerasm(i,ns,7) + tau_w_g_tot aerfwd(i,ns,1) = aerfwd(i,ns,1) + tau_w_f_sul aerfwd(i,ns,2) = aerfwd(i,ns,2) + tau_w_f_sslt aerfwd(i,ns,3) = aerfwd(i,ns,3) + tau_w_f_cphil + tau_w_f_cphob + tau_w_f_cb aerfwd(i,ns,4) = aerfwd(i,ns,4) + tau_w_f_dst_tot aerfwd(i,ns,5) = aerfwd(i,ns,5) + tau_w_f_bg aerfwd(i,ns,6) = aerfwd(i,ns,6) + tau_w_f_volc aerfwd(i,ns,7) = aerfwd(i,ns,7) + tau_w_f_tot end do do kaer = 1, naer_groups if (aerssa(i,ns,kaer) .gt. 0.0) then aerasm(i,ns,kaer) = aerasm(i,ns,kaer) / aerssa(i,ns,kaer) aerfwd(i,ns,kaer) = aerfwd(i,ns,kaer) / aerssa(i,ns,kaer) else aerasm(i,ns,kaer) = 0.0_r8 aerfwd(i,ns,kaer) = 0.0_r8 end if if (aertau(i,ns,kaer) .gt. 0.0) then aerssa(i,ns,kaer) = aerssa(i,ns,kaer) / aertau(i,ns,kaer) else aerssa(i,ns,kaer) = 0.0_r8 end if end do end do wavmid(ns) = 0.5_r8*(wavmin(ns) + wavmax(ns)) if (wavmid(ns) < 0.7_r8 ) then do n=1,ndayc i=idayc(n) albdir(i,ns) = asdir(i) albdif(i,ns) = asdif(i) end do else do n=1,ndayc i=idayc(n) albdir(i,ns) = aldir(i) albdif(i,ns) = aldif(i) end do end if trayoslp = raytau(ns)/sslp call raddedmx(pver, pverp, pcols, coszrs ,ndayc ,idayc , & abh2o(ns),abo3(ns) ,abco2(ns),abo2(ns) , & uh2o ,uo3 ,uco2 ,uo2 , & trayoslp ,pflx ,ns , & tauxcl ,wcl ,gcl ,fcl , & tauxci ,wci ,gci ,fci , & tauxar ,wa ,ga ,fa , & rdir ,rdif ,tdir ,tdif ,explay , & rdirc ,rdifc ,tdirc ,tdifc ,explayc ) end do do n=1,ndayc i=idayc(n) npasses = 0 do do irgn = 0, nmxrgn(i) kx2(irgn) = 0 end do mrgn = 0 do irgn = 1, nmxrgn(i) region_found = .false. if (kx2(irgn-1) < pver) then k1 = kx2(irgn-1)+1 kx1(irgn) = k1 kx2(irgn) = k1-1 do k2 = pver, k1, -1 if (pmid(i,k2) <= pmxrgn(i,irgn)) then kx2(irgn) = k2 mrgn = mrgn+1 region_found = .true. exit end if end do else exit endif if (region_found) then nxs = 0 if (cldeps > 0) then do k = k1,k2 if (cld(i,k) >= cldmin .and. cld(i,k) >= cldeps) then nxs = nxs+1 ksort(nxs) = k asort(nxs) = 1.0_r8-(floor(cld(i,k)/cldeps)*cldeps) end if end do else do k = k1,k2 if (cld(i,k) >= cldmin) then nxs = nxs+1 ksort(nxs) = k asort(nxs) = 1.0_r8-cld(i,k) end if end do endif if (nxs == 2) then if (asort(2) < asort(1)) then ktmp = ksort(1) ksort(1) = ksort(2) ksort(2) = ktmp atmp = asort(1) asort(1) = asort(2) asort(2) = atmp endif else if (nxs >= 3) then call sortarray(nxs,asort,ksort) endif cstr(k1:k2,1:nxs+1) = 0 mstr = 1 cld0 = 0.0_r8 do l = 1, nxs if (asort(l) /= cld0) then wstr(mstr,mrgn) = asort(l) - cld0 cld0 = asort(l) mstr = mstr + 1 endif cstr(ksort(l),mstr:nxs+1) = 1 end do nstr(mrgn) = mstr wstr(mstr,mrgn) = 1.0_r8 - cld0 endif end do nrgn = mrgn cstr(0,1:nstr(1)) = 0 nconfigm = product(nstr(1: nrgn)) istr(1: nrgn) = 1 nconfig = 0 totwgt = 0.0_r8 new_term = .true. do iconfig = 1, nconfigm xwgt = 1.0_r8 do mrgn = 1, nrgn xwgt = xwgt * wstr(istr(mrgn),mrgn) end do if (xwgt >= areamin) then nconfig = nconfig + 1 if (nconfig <= nconfgmax) then j = nconfig ptrc(nconfig) = nconfig else nconfig = nconfgmax if (new_term) then j = findvalue(1,nconfig,wgtv,ptrc) endif if (wgtv(j) < xwgt) then totwgt = totwgt - wgtv(j) new_term = .true. else new_term = .false. endif endif if (new_term) then wgtv(j) = xwgt totwgt = totwgt + xwgt do mrgn = 1, nrgn ccon(kx1(mrgn):kx2(mrgn),j) = cstr(kx1(mrgn):kx2(mrgn),istr(mrgn)) end do endif endif mrgn = nrgn istr(mrgn) = istr(mrgn) + 1 do while (istr(mrgn) > nstr(mrgn) .and. mrgn > 1) istr(mrgn) = 1 mrgn = mrgn - 1 istr(mrgn) = istr(mrgn) + 1 end do end do if (totwgt > 0.) then exit else npasses = npasses + 1 if (npasses >= 2 ) then write(6,*)'RADCSWMX: Maximum overlap of column ','failed' call endrun endif nmxrgn(i)=1 pmxrgn(i,1)=1.0e30 end if end do ccon(0,:) = 0 ccon(pverp,:) = 0 nuniqd(0) = 1 nuniqu(pverp) = 1 istrtd(0,1) = 1 istrtu(pverp,1) = 1 do j = 1, nconfig icond(0,j)=j iconu(pverp,j)=j end do istrtd(0,2) = nconfig+1 istrtu(pverp,2) = nconfig+1 do k = 1, pverp km1 = k-1 nuniq = 0 istrtd(k,1) = 1 do l0 = 1, nuniqd(km1) is0 = istrtd(km1,l0) is1 = istrtd(km1,l0+1)-1 n0 = 0 n1 = 0 do isn = is0, is1 j = icond(km1,isn) if (ccon(k,j) == 0) then n0 = n0 + 1 ptr0(n0) = j endif if (ccon(k,j) == 1) then n1 = n1 + 1 ptr1(n1) = j endif end do if (n0 > 0) then nuniq = nuniq + 1 istrtd(k,nuniq+1) = istrtd(k,nuniq)+n0 icond(k,istrtd(k,nuniq):istrtd(k,nuniq+1)-1) = ptr0(1:n0) endif if (n1 > 0) then nuniq = nuniq + 1 istrtd(k,nuniq+1) = istrtd(k,nuniq)+n1 icond(k,istrtd(k,nuniq):istrtd(k,nuniq+1)-1) = ptr1(1:n1) endif end do nuniqd(k) = nuniq end do do k = pver, 0, -1 kp1 = k+1 nuniq = 0 istrtu(k,1) = 1 do l0 = 1, nuniqu(kp1) is0 = istrtu(kp1,l0) is1 = istrtu(kp1,l0+1)-1 n0 = 0 n1 = 0 do isn = is0, is1 j = iconu(kp1,isn) if (ccon(k,j) == 0) then n0 = n0 + 1 ptr0(n0) = j endif if (ccon(k,j) == 1) then n1 = n1 + 1 ptr1(n1) = j endif end do if (n0 > 0) then nuniq = nuniq + 1 istrtu(k,nuniq+1) = istrtu(k,nuniq)+n0 iconu(k,istrtu(k,nuniq):istrtu(k,nuniq+1)-1) = ptr0(1:n0) endif if (n1 > 0) then nuniq = nuniq + 1 istrtu(k,nuniq+1) = istrtu(k,nuniq)+n1 iconu(k,istrtu(k,nuniq):istrtu(k,nuniq+1)-1) = ptr1(1:n1) endif end do nuniqu(k) = nuniq end do do k=0,pver totfld(k) = 0.0_r8 fswup (k) = 0.0_r8 fswdn (k) = 0.0_r8 fswupc (k) = 0.0_r8 fswdnc (k) = 0.0_r8 fswdndir(k) = 0.0_r8 fswdncdir(k)= 0.0_r8 end do sfltot = 0.0_r8 fswup (pverp) = 0.0_r8 fswdn (pverp) = 0.0_r8 fswupc (pverp) = 0.0_r8 fswdnc (pverp) = 0.0_r8 fswdndir(pverp) = 0.0_r8 fswdncdir(pverp)= 0.0_r8 do ns = 1,nspint wgtint = nirwgt(ns) rdndif(0,1:nconfig) = 0.0_r8 exptdn(0,1:nconfig) = 1.0_r8 tdntot(0,1:nconfig) = 1.0_r8 do k = 1, pverp km1 = k - 1 do l0 = 1, nuniqd(km1) is0 = istrtd(km1,l0) is1 = istrtd(km1,l0+1)-1 j = icond(km1,is0) xexpt = exptdn(km1,j) xrdnd = rdndif(km1,j) tdnmexp = tdntot(km1,j) - xexpt if (ccon(km1,j) == 1) then ytdnd = tdif(ns,i,km1) yrdnd = rdif(ns,i,km1) rdenom = 1._r8/(1._r8-yrdnd*xrdnd) rdirexp = rdir(ns,i,km1)*xexpt zexpt = xexpt * explay(ns,i,km1) zrdnd = yrdnd + xrdnd*(ytdnd**2)*rdenom ztdnt = xexpt*tdir(ns,i,km1) + ytdnd*(tdnmexp + xrdnd*rdirexp)*rdenom else ytdnd = tdifc(ns,i,km1) yrdnd = rdifc(ns,i,km1) rdenom = 1._r8/(1._r8-yrdnd*xrdnd) rdirexp = rdirc(ns,i,km1)*xexpt zexpt = xexpt * explayc(ns,i,km1) zrdnd = yrdnd + xrdnd*(ytdnd**2)*rdenom ztdnt = xexpt*tdirc(ns,i,km1) + ytdnd* & (tdnmexp + xrdnd*rdirexp)*rdenom endif do isn = is0, is1 j = icond(km1,isn) exptdn(k,j) = zexpt rdndif(k,j) = zrdnd tdntot(k,j) = ztdnt end do end do end do rupdir(pverp,1:nconfig) = albdir(i,ns) rupdif(pverp,1:nconfig) = albdif(i,ns) do k = pver, 0, -1 do l0 = 1, nuniqu(k) is0 = istrtu(k,l0) is1 = istrtu(k,l0+1)-1 j = iconu(k,is0) xrupd = rupdif(k+1,j) xrups = rupdir(k+1,j) if (ccon(k,j) == 1) then yexpt = explay(ns,i,k) yrupd = rdif(ns,i,k) ytupd = tdif(ns,i,k) rdenom = 1._r8/( 1._r8 - yrupd*xrupd) tdnmexp = (tdir(ns,i,k)-yexpt) rdirexp = xrups*yexpt zrupd = yrupd + xrupd*(ytupd**2)*rdenom zrups = rdir(ns,i,k) + ytupd*(rdirexp + xrupd*tdnmexp)*rdenom else yexpt = explayc(ns,i,k) yrupd = rdifc(ns,i,k) ytupd = tdifc(ns,i,k) rdenom = 1._r8/( 1._r8 - yrupd*xrupd) tdnmexp = (tdirc(ns,i,k)-yexpt) rdirexp = xrups*yexpt zrupd = yrupd + xrupd*(ytupd**2)*rdenom zrups = rdirc(ns,i,k) + ytupd*(rdirexp + xrupd*tdnmexp)*rdenom endif do isn = is0, is1 j = iconu(k,isn) rupdif(k,j) = zrupd rupdir(k,j) = zrups end do end do end do do k = 0,pverp fluxup(k)=0.0_r8 fluxdn(k)=0.0_r8 fluxdndir(k) = 0.0_r8 do iconfig = 1, nconfig xwgt = wgtv(iconfig) xexpt = exptdn(k,iconfig) xtdnt = tdntot(k,iconfig) xrdnd = rdndif(k,iconfig) xrupd = rupdif(k,iconfig) xrups = rupdir(k,iconfig) rdenom = 1._r8/(1._r8 - xrdnd * xrupd) fluxup(k) = fluxup(k) + xwgt * & ((xexpt * xrups + (xtdnt - xexpt) * xrupd) * rdenom) fluxdn(k) = fluxdn(k) + xwgt * & (xexpt + (xtdnt - xexpt + xexpt * xrups * xrdnd) * rdenom) fluxdndir(k) = fluxdndir(k) + xwgt * xexpt end do fluxup(k)=fluxup(k) / totwgt fluxdn(k)=fluxdn(k) / totwgt fluxdndir(k)=fluxdndir(k) / totwgt end do wexptdn = 0.0_r8 do iconfig = 1, nconfig wexptdn = wexptdn + wgtv(iconfig) * exptdn(pverp,iconfig) end do wexptdn = wexptdn / totwgt solflx = solin(i)*frcsol(ns)*psf(ns) fsnt(i) = fsnt(i) + solflx*(fluxdn(1) - fluxup(1)) fsntoa(i)= fsntoa(i) + solflx*(fluxdn(0) - fluxup(0)) fsns(i) = fsns(i) + solflx*(fluxdn(pverp)-fluxup(pverp)) sfltot = sfltot + solflx fswup(0) = fswup(0) + solflx*fluxup(0) fswdn(0) = fswdn(0) + solflx*fluxdn(0) fswdndir(0) = fswdndir(0) + solflx*fluxdndir(0) if (wavmid(ns) < 0.7_r8) then sols(i) = sols(i) + wexptdn*solflx*0.001_r8 solsd(i) = solsd(i)+(fluxdn(pverp)-wexptdn)*solflx*0.001_r8 else soll(i) = soll(i) + wexptdn*solflx*0.001_r8 solld(i) = solld(i)+(fluxdn(pverp)-wexptdn)*solflx*0.001_r8 fsnrtoaq(i) = fsnrtoaq(i) + solflx*(fluxdn(0) - fluxup(0)) end if fsnirtoa(i) = fsnirtoa(i) + wgtint*solflx*(fluxdn(0) - fluxup(0)) do k=0,pver kp1 = k+1 flxdiv = (fluxdn(k ) - fluxdn(kp1)) + (fluxup(kp1) - fluxup(k )) totfld(k) = totfld(k) + solflx*flxdiv fswdn(kp1) = fswdn(kp1) + solflx*fluxdn(kp1) fswup(kp1) = fswup(kp1) + solflx*fluxup(kp1) fswdndir(kp1) = fswdndir(kp1) + solflx*fluxdndir(kp1) end do exptdnc(0) = 1.0_r8 rdndifc(0) = 0.0_r8 tdntotc(0) = 1.0_r8 rupdirc(pverp) = albdir(i,ns) rupdifc(pverp) = albdif(i,ns) do k = 1, pverp km1 = k - 1 xexpt = exptdnc(km1) xrdnd = rdndifc(km1) yrdnd = rdifc(ns,i,km1) ytdnd = tdifc(ns,i,km1) exptdnc(k) = xexpt*explayc(ns,i,km1) rdenom = 1._r8/(1._r8 - yrdnd*xrdnd) rdirexp = rdirc(ns,i,km1)*xexpt tdnmexp = tdntotc(km1) - xexpt tdntotc(k) = xexpt*tdirc(ns,i,km1) + ytdnd*(tdnmexp + xrdnd*rdirexp)* & rdenom rdndifc(k) = yrdnd + xrdnd*(ytdnd**2)*rdenom end do do k=pver,0,-1 xrupd = rupdifc(k+1) yexpt = explayc(ns,i,k) yrupd = rdifc(ns,i,k) ytupd = tdifc(ns,i,k) rdenom = 1._r8/( 1._r8 - yrupd*xrupd) rupdirc(k) = rdirc(ns,i,k) + ytupd*(rupdirc(k+1)*yexpt + & xrupd*(tdirc(ns,i,k)-yexpt))*rdenom rupdifc(k) = yrupd + xrupd*ytupd**2*rdenom end do do k=0,1 rdenom = 1._r8/(1._r8 - rdndifc(k)*rupdifc(k)) fluxup(k) = (exptdnc(k)*rupdirc(k) + (tdntotc(k)-exptdnc(k))*rupdifc(k))* & rdenom fluxdn(k) = exptdnc(k) + & (tdntotc(k) - exptdnc(k) + exptdnc(k)*rupdirc(k)*rdndifc(k))* & rdenom fswupc(k) = fswupc(k) + solflx*fluxup(k) fswdnc(k) = fswdnc(k) + solflx*fluxdn(k) fswdncdir(k) = fswdncdir(k) + solflx*exptdnc(k) end do do k=2,pverp rdenom = 1._r8/(1._r8 - rdndifc(k)*rupdifc(k)) fluxup(k) = (exptdnc(k)*rupdirc(k) + (tdntotc(k)-exptdnc(k))*rupdifc(k))* & rdenom fluxdn(k) = exptdnc(k) + (tdntotc(k) - exptdnc(k) + & exptdnc(k)*rupdirc(k)*rdndifc(k))*rdenom fswupc(k) = fswupc(k) + solflx*fluxup(k) fswdnc(k) = fswdnc(k) + solflx*fluxdn(k) fswdncdir(k) = fswdncdir(k) + solflx*exptdnc(k) end do fsntc(i) = fsntc(i)+solflx*(fluxdn(1)-fluxup(1)) fsntoac(i) = fsntoac(i)+solflx*(fluxdn(0)-fluxup(0)) fsnsc(i) = fsnsc(i)+solflx*(fluxdn(pverp)-fluxup(pverp)) fsdsc(i) = fsdsc(i)+solflx*(fluxdn(pverp)) fsnrtoac(i) = fsnrtoac(i)+wgtint*solflx*(fluxdn(0)-fluxup(0)) end do do k=1,pver qrs(i,k) = -1.E-4*gravit*totfld(k)/(pint(i,k) - pint(i,k+1)) end do do k=1,pverp fsup(i,k) = fswup(k) fsupc(i,k) = fswupc(k) fsdn(i,k) = fswdn(k) fsdnc(i,k) = fswdnc(k) fsdndir(i,k) = fswdndir(k) fsdncdir(i,k) = fswdncdir(k) fsdndif(i,k) = fswdn(k)-fswdndir(k) fsdncdif(i,k) = fswdnc(k)-fswdncdir(k) end do fsds(i) = fswdn(pverp) fsdscdir(i) = fsdncdir(i,pverp) fsdscdif(i) = fsdncdif(i,pverp) fsdsdir(i) = fsdndir(i,pverp) fsdsdif(i) = fsdndif(i,pverp) end do return end subroutine radcswmx subroutine raddedmx(pver, pverp, pcols, coszrs ,ndayc ,idayc ,abh2o , & abo3 ,abco2 ,abo2 ,uh2o ,uo3 , & uco2 ,uo2 ,trayoslp,pflx ,ns , & tauxcl ,wcl ,gcl ,fcl ,tauxci , & wci ,gci ,fci ,tauxar ,wa , & ga ,fa ,rdir ,rdif ,tdir , & tdif ,explay ,rdirc ,rdifc ,tdirc , & tdifc ,explayc ) implicit none integer nspint parameter ( nspint = 19 ) real(r8) trmin real(r8) wray real(r8) gray real(r8) fray parameter (trmin = 1.e-3) parameter (wray = 0.999999) parameter (gray = 0.0) parameter (fray = 0.1) integer, intent(in) :: pver, pverp, pcols real(r8), intent(in) :: coszrs(pcols) real(r8), intent(in) :: trayoslp real(r8), intent(in) :: pflx(pcols,0:pverp) real(r8), intent(in) :: abh2o real(r8), intent(in) :: abo3 real(r8), intent(in) :: abco2 real(r8), intent(in) :: abo2 real(r8), intent(in) :: uh2o(pcols,0:pver) real(r8), intent(in) :: uo3(pcols,0:pver) real(r8), intent(in) :: uco2(pcols,0:pver) real(r8), intent(in) :: uo2(pcols,0:pver) real(r8), intent(in) :: tauxcl(pcols,0:pver) real(r8), intent(in) :: wcl(pcols,0:pver) real(r8), intent(in) :: gcl(pcols,0:pver) real(r8), intent(in) :: fcl(pcols,0:pver) real(r8), intent(in) :: tauxci(pcols,0:pver) real(r8), intent(in) :: wci(pcols,0:pver) real(r8), intent(in) :: gci(pcols,0:pver) real(r8), intent(in) :: fci(pcols,0:pver) real(r8), intent(in) :: tauxar(pcols,0:pver) real(r8), intent(in) :: wa(pcols,0:pver) real(r8), intent(in) :: ga(pcols,0:pver) real(r8), intent(in) :: fa(pcols,0:pver) integer, intent(in) :: ndayc integer, intent(in) :: idayc(pcols) integer, intent(in) :: ns real(r8), intent(inout) :: rdir(nspint,pcols,0:pver) real(r8), intent(inout) :: rdif(nspint,pcols,0:pver) real(r8), intent(inout) :: tdir(nspint,pcols,0:pver) real(r8), intent(inout) :: tdif(nspint,pcols,0:pver) real(r8), intent(inout) :: explay(nspint,pcols,0:pver) real(r8), intent(inout) :: rdirc(nspint,pcols,0:pver) real(r8), intent(inout) :: rdifc(nspint,pcols,0:pver) real(r8), intent(inout) :: tdirc(nspint,pcols,0:pver) real(r8), intent(inout) :: tdifc(nspint,pcols,0:pver) real(r8), intent(inout) :: explayc(nspint,pcols,0:pver) integer i integer k integer nn real(r8) taugab(pcols) real(r8) tauray(pcols) real(r8) taucsc real(r8) tautot real(r8) wtot real(r8) gtot real(r8) ftot real(r8) wtau real(r8) wt real(r8) ts real(r8) ws real(r8) gs real(r8) alpha real(r8) gamma real(r8) el real(r8) taus real(r8) omgs real(r8) asys real(r8) u real(r8) n real(r8) lm real(r8) ne real(r8) w real(r8) uu real(r8) g real(r8) e real(r8) f real(r8) t real(r8) et real(r8) alp real(r8) gam real(r8) ue real(r8) arg real(r8) extins real(r8) amg real(r8) apg alpha(w,uu,g,e) = .75_r8*w*uu*((1._r8 + g*(1._r8-w))/(1._r8 - e*e*uu*uu)) gamma(w,uu,g,e) = .50_r8*w*((3._r8*g*(1._r8-w)*uu*uu + 1._r8)/(1._r8-e*e*uu*uu)) el(w,g) = sqrt(3._r8*(1._r8-w)*(1._r8 - w*g)) taus(w,f,t) = (1._r8 - w*f)*t omgs(w,f) = (1._r8 - f)*w/(1._r8 - w*f) asys(g,f) = (g - f)/(1._r8 - f) u(w,g,e) = 1.5_r8*(1._r8 - w*g)/e n(uu,et) = ((uu+1._r8)*(uu+1._r8)/et ) - ((uu-1._r8)*(uu-1._r8)*et) do k=0,pver do nn=1,ndayc i=idayc(nn) tauray(i) = trayoslp*(pflx(i,k+1)-pflx(i,k)) taugab(i) = abh2o*uh2o(i,k) + abo3*uo3(i,k) + abco2*uco2(i,k) + abo2*uo2(i,k) tautot = tauxcl(i,k) + tauxci(i,k) + tauray(i) + taugab(i) + tauxar(i,k) taucsc = tauxcl(i,k)*wcl(i,k) + tauxci(i,k)*wci(i,k) + tauxar(i,k)*wa(i,k) wtau = wray*tauray(i) wt = wtau + taucsc wtot = wt/tautot gtot = (wtau*gray + gcl(i,k)*wcl(i,k)*tauxcl(i,k) & + gci(i,k)*wci(i,k)*tauxci(i,k) + ga(i,k) *wa(i,k) *tauxar(i,k))/wt ftot = (wtau*fray + fcl(i,k)*wcl(i,k)*tauxcl(i,k) & + fci(i,k)*wci(i,k)*tauxci(i,k) + fa(i,k) *wa(i,k) *tauxar(i,k))/wt ts = taus(wtot,ftot,tautot) ws = omgs(wtot,ftot) gs = asys(gtot,ftot) lm = el(ws,gs) alp = alpha(ws,coszrs(i),gs,lm) gam = gamma(ws,coszrs(i),gs,lm) ue = u(ws,gs,lm) arg = min(lm*ts,25._r8) extins = exp(-arg) ne = n(ue,extins) rdif(ns,i,k) = (ue+1._r8)*(ue-1._r8)*(1._r8/extins - extins)/ne tdif(ns,i,k) = 4._r8*ue/ne arg = min(ts/coszrs(i),25._r8) explay(ns,i,k) = exp(-arg) apg = alp + gam amg = alp - gam rdir(ns,i,k) = amg*(tdif(ns,i,k)*explay(ns,i,k)-1._r8) + apg*rdif(ns,i,k) tdir(ns,i,k) = apg*tdif(ns,i,k) + (amg*rdif(ns,i,k)-(apg-1._r8))*explay(ns,i,k) rdir(ns,i,k) = max(rdir(ns,i,k),0.0_r8) tdir(ns,i,k) = max(tdir(ns,i,k),0.0_r8) rdif(ns,i,k) = max(rdif(ns,i,k),0.0_r8) tdif(ns,i,k) = max(tdif(ns,i,k),0.0_r8) if (tauxcl(i,k) == 0.0_r8 .and. tauxci(i,k) == 0.0_r8) then rdirc(ns,i,k) = rdir(ns,i,k) tdirc(ns,i,k) = tdir(ns,i,k) rdifc(ns,i,k) = rdif(ns,i,k) tdifc(ns,i,k) = tdif(ns,i,k) explayc(ns,i,k) = explay(ns,i,k) else tautot = tauray(i) + taugab(i) + tauxar(i,k) taucsc = tauxar(i,k)*wa(i,k) wt = wtau + taucsc wtot = wt/tautot gtot = (wtau*gray + ga(i,k)*wa(i,k)*tauxar(i,k))/wt ftot = (wtau*fray + fa(i,k)*wa(i,k)*tauxar(i,k))/wt ts = taus(wtot,ftot,tautot) ws = omgs(wtot,ftot) gs = asys(gtot,ftot) lm = el(ws,gs) alp = alpha(ws,coszrs(i),gs,lm) gam = gamma(ws,coszrs(i),gs,lm) ue = u(ws,gs,lm) arg = min(lm*ts,25._r8) extins = exp(-arg) ne = n(ue,extins) rdifc(ns,i,k) = (ue+1._r8)*(ue-1._r8)*(1._r8/extins - extins)/ne tdifc(ns,i,k) = 4._r8*ue/ne arg = min(ts/coszrs(i),25._r8) explayc(ns,i,k) = exp(-arg) apg = alp + gam amg = alp - gam rdirc(ns,i,k) = amg*(tdifc(ns,i,k)*explayc(ns,i,k)-1._r8)+ & apg*rdifc(ns,i,k) tdirc(ns,i,k) = apg*tdifc(ns,i,k) + (amg*rdifc(ns,i,k) - (apg-1._r8))* & explayc(ns,i,k) rdirc(ns,i,k) = max(rdirc(ns,i,k),0.0_r8) tdirc(ns,i,k) = max(tdirc(ns,i,k),0.0_r8) rdifc(ns,i,k) = max(rdifc(ns,i,k),0.0_r8) tdifc(ns,i,k) = max(tdifc(ns,i,k),0.0_r8) end if end do end do return end subroutine raddedmx subroutine radinp(lchnk ,ncol , pcols, pver, pverp, & pmid ,pint ,o3vmr , pmidrd ,& pintrd ,eccf ,o3mmr ) implicit none integer, intent(in) :: lchnk integer, intent(in) :: pcols, pver, pverp integer, intent(in) :: ncol real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) real(r8), intent(in) :: o3vmr(pcols,pver) real(r8), intent(out) :: pmidrd(pcols,pver) real(r8), intent(out) :: pintrd(pcols,pverp) real(r8), intent(out) :: eccf real(r8), intent(out) :: o3mmr(pcols,pver) integer i integer k real(r8) :: calday real(r8) vmmr real(r8) delta eccf = 1. do k=1,pver do i=1,ncol pmidrd(i,k) = pmid(i,k)*10.0 pintrd(i,k) = pint(i,k)*10.0 end do end do do i=1,ncol pintrd(i,pverp) = pint(i,pverp)*10.0 end do vmmr = amo/amd do k=1,pver do i=1,ncol o3mmr(i,k) = vmmr*o3vmr(i,k) end do end do return end subroutine radinp subroutine radoz2(lchnk ,ncol ,pcols, pver, pverp, o3vmr ,pint ,plol ,plos, ntoplw ) implicit none integer, intent(in) :: lchnk integer, intent(in) :: ncol integer, intent(in) :: pcols, pver, pverp real(r8), intent(in) :: o3vmr(pcols,pver) real(r8), intent(in) :: pint(pcols,pverp) integer, intent(in) :: ntoplw real(r8), intent(out) :: plol(pcols,pverp) real(r8), intent(out) :: plos(pcols,pverp) integer i integer k do i=1,ncol plos(i,ntoplw) = 0.1 *cplos*o3vmr(i,ntoplw)*pint(i,ntoplw) plol(i,ntoplw) = 0.01*cplol*o3vmr(i,ntoplw)*pint(i,ntoplw)*pint(i,ntoplw) end do do k=ntoplw+1,pverp do i=1,ncol plos(i,k) = plos(i,k-1) + 0.1*cplos*o3vmr(i,k-1)*(pint(i,k) - pint(i,k-1)) plol(i,k) = plol(i,k-1) + 0.01*cplol*o3vmr(i,k-1)* & (pint(i,k)*pint(i,k) - pint(i,k-1)*pint(i,k-1)) end do end do return end subroutine radoz2 subroutine radozn (lchnk, ncol, pcols, pver,pmid, pin, levsiz, ozmix, o3vmr) implicit none integer, intent(in) :: lchnk integer, intent(in) :: pcols, pver integer, intent(in) :: ncol integer, intent(in) :: levsiz real(r8), intent(in) :: pmid(pcols,pver) real(r8), intent(in) :: pin(levsiz) real(r8), intent(in) :: ozmix(pcols,levsiz) real(r8), intent(out) :: o3vmr(pcols,pver) integer i integer k, kk, kkstart integer kupper(pcols) integer kount integer lats(pcols) integer lons(pcols) real(r8) dpu real(r8) dpl do i=1,ncol kupper(i) = 1 end do do k=1,pver kkstart = levsiz do i=1,ncol kkstart = min0(kkstart,kupper(i)) end do kount = 0 do kk=kkstart,levsiz-1 do i=1,ncol if (pin(kk).lt.pmid(i,k) .and. pmid(i,k).le.pin(kk+1)) then kupper(i) = kk kount = kount + 1 end if end do if (kount.eq.ncol) then do i=1,ncol dpu = pmid(i,k) - pin(kupper(i)) dpl = pin(kupper(i)+1) - pmid(i,k) o3vmr(i,k) = (ozmix(i,kupper(i))*dpl + & ozmix(i,kupper(i)+1)*dpu)/(dpl + dpu) end do goto 35 end if end do do i=1,ncol if (pmid(i,k) .lt. pin(1)) then o3vmr(i,k) = ozmix(i,1)*pmid(i,k)/pin(1) else if (pmid(i,k) .gt. pin(levsiz)) then o3vmr(i,k) = ozmix(i,levsiz) else dpu = pmid(i,k) - pin(kupper(i)) dpl = pin(kupper(i)+1) - pmid(i,k) o3vmr(i,k) = (ozmix(i,kupper(i))*dpl + & ozmix(i,kupper(i)+1)*dpu)/(dpl + dpu) end if end do if (kount.gt.ncol) then call endrun ('RADOZN: Bad ozone data: non-monotonicity suspected') end if 35 continue end do return end subroutine radozn end MODULE module_ra_cam