subroutine precpd_ad( km, dt, del, sl, ps, rhc, q_in, & cwm_in, t_in, q_out, cwm_out, t_out, rn_out, q_in_ad, cwm_in_ad, & t_in_ad, q_out_ad, cwm_out_ad, t_out_ad, rn_out_ad, & adjoint) !$$$ subprogram documentation block ! . . . . ! subprogram: precpd_ad compute precipitation processes from suspended cloud water/ice ! prgmmr: treadon org: np23 date: 2003-12-18 ! ! abstract: This routine parameterizes the effect of precipitation processes from ! suspended cloud water and ice on temperature, moisture, and cloud ! condensate. For more details regarding the algorithm, please refer ! to ! Zhao and Carr (1997), Monthly Weather Review (August) ! Sundqvist et al., (1989) Monthly Weather review. (August) ! ! ! program history log: ! 1995-01-01 zhao - original routine ! 1998-10-10 moorthi,pan - modified and rewritten for GFS application ! 2003-12-18 treadon - add adjoint code ! 2004-06-14 treadon - reformat documenation ! 2004-08-04 treadon - add only on use declarations; add intent in/out ! 2006-04-12 treadon - change del and sl from 1d to 2d arrays ! 2008-04-24 safford - rm unused vars ! ! input argument list: ! km - number of vertical levels ! dt - time step in seconds ! del - pressure layer thickness (bottom to top) ! sl - sigma level ! ps - surface pressure (centibars) ! rhc - critical relative humidity threshold ! q_in - specific humidity ! cwm_in - condensate mixing ratio ! t_in - temperature ! q_out_ad - specific humidity perturbation ! cwm_out_ad - cloud condensate mixing perturbation ! t_out_ad - temperature perturbation ! rn_out_ad - precipitation perturbation ! adjoint - logical flag (.false.=forward model only, .true.=forward and ajoint) ! ! output argument list: ! q_out - q following grid scale precipitation ! cwm_out - cloud condensate mixing ratio following precipitation ! t_out - temperature following precipitation ! rn_out - precipitation over one time step, dt ! q_in_ad - change in rain rate with respect to moisture ! cwm_in_ad - change in rain rate with respect to cloud condensate mixing ratio ! t_in_ad - change in rain rate with respect to temperature ! ! remarks: ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ ! !*************************************************************** !============================================== ! all entries are defined explicitly !============================================== use kinds, only: r_kind,i_kind use constants, only: cmr,cws,half,epsm1,h300,rrow,two,& hfus,hvap,hsub,rcp,h1000,ke2,zero,one,ttp,eps,epsq,& grav,climit implicit none !============================================== ! define arguments !============================================== logical ,intent(in ):: adjoint integer(i_kind),intent(in ):: km real(r_kind) ,intent(inout):: cwm_in_ad(km) real(r_kind) ,intent(inout):: cwm_out_ad(km) real(r_kind) ,intent(inout):: q_in_ad(km) real(r_kind) ,intent(inout):: q_out_ad(km) real(r_kind) ,intent(inout):: rn_out_ad real(r_kind) ,intent(inout):: t_in_ad(km) real(r_kind) ,intent(inout):: t_out_ad(km) real(r_kind) ,intent(in ):: cwm_in(km) real(r_kind) ,intent( out):: cwm_out(km) real(r_kind) ,intent(in ):: del(km) real(r_kind) ,intent(in ):: dt real(r_kind) ,intent(in ):: ps real(r_kind) ,intent(in ):: q_in(km) real(r_kind) ,intent( out):: q_out(km) real(r_kind) ,intent(in ):: rhc(km) real(r_kind) ,intent( out):: rn_out real(r_kind) ,intent(in ):: sl(km) real(r_kind) ,intent(in ):: t_in(km) real(r_kind) ,intent( out):: t_out(km) !============================================== ! define local variables !============================================== real(r_kind) aa2 real(r_kind) amaxcmr_ad real(r_kind) amaxcms_ad real(r_kind) amaxps_ad real(r_kind) amaxrq_ad real(r_kind) ccr_ad real(r_kind) cwmin_ad real(r_kind) cwmk_ad real(r_kind) cwmks_ad real(r_kind) erk_ad real(r_kind) err_ad real(r_kind) ers_ad real(r_kind) es_ad real(r_kind) expf_ad real(r_kind) fi_ad real(r_kind) ppr0_ad real(r_kind) ppr1_ad real(r_kind) ppr2_ad real(r_kind) pps0_ad real(r_kind) pps1_ad real(r_kind) pps2_ad real(r_kind) praut_ad real(r_kind) praut0_ad real(r_kind) precrk0_ad real(r_kind) precrk1_ad real(r_kind) precrl1k_ad real(r_kind) precrl_0_ad real(r_kind) precrl_1_ad real(r_kind) precrl_2_ad real(r_kind) precrl_3_ad real(r_kind) precsk0_ad real(r_kind) precsk1_ad real(r_kind) precsl1k_ad real(r_kind) precsl_0_ad real(r_kind) precsl_1_ad real(r_kind) precsl_2_ad real(r_kind) precsl_3_ad real(r_kind) psaci_ad real(r_kind) psaut_ad real(r_kind) psm1_ad real(r_kind) psm2_ad real(r_kind) psm3_ad real(r_kind) q2_ad real(r_kind) qin_ad real(r_kind) qk_ad real(r_kind) qs_ad real(r_kind) qs0_ad real(r_kind) rprs_ad real(r_kind) rq_ad real(r_kind) rqkll_ad real(r_kind) t2_ad real(r_kind) tem1_ad real(r_kind) tem2_ad real(r_kind) tem3_ad real(r_kind) tem4_ad real(r_kind) tem5_ad real(r_kind) term_ad real(r_kind) term1_ad real(r_kind) term2r_ad real(r_kind) term2s_ad real(r_kind) term3r_ad real(r_kind) term3s_ad real(r_kind) term4_ad real(r_kind) term5_ad real(r_kind) term6_ad real(r_kind) tin_ad real(r_kind) tmt0_ad real(r_kind) tmt0k_ad real(r_kind) ww0_ad real(r_kind) ww1_ad real(r_kind) ww2_ad real(r_kind) wwn_ad real(r_kind) amaxcmr real(r_kind) amaxcms real(r_kind) amaxps real(r_kind) amaxrq real(r_kind) c00 real(r_kind) ccr logical comput,comput0 real(r_kind) conde real(r_kind) condt real(r_kind) const real(r_kind) cr real(r_kind) crs1 real(r_kind) crs2 real(r_kind) csm1 real(r_kind) cwmin real(r_kind) cwmk real(r_kind) cwmks real(r_kind) dtcp real(r_kind) dum1 real(r_kind) dum2 real(r_kind) erk real(r_kind) err real(r_kind) ers real(r_kind) es real(r_kind) expf real(r_kind) fi integer(i_kind) iwl integer(i_kind) iwl1,iwl1k(km+1) integer(i_kind) k real(r_kind) ppr0 real(r_kind) ppr1 real(r_kind) ppr2 real(r_kind) pps0 real(r_kind) pps1 real(r_kind) pps2 real(r_kind) praut real(r_kind) praut0 real(r_kind) precrk0 real(r_kind) precrk1 real(r_kind) precrl1k,precrl1ki(km+1) real(r_kind) precrl_0 real(r_kind) precrl_1 real(r_kind) precrl_2 real(r_kind) precrl_3 real(r_kind) precsk0 real(r_kind) precsk1 real(r_kind) precsl1k,precsl1ki(km+1) real(r_kind) precsl_0 real(r_kind) precsl_1 real(r_kind) precsl_2 real(r_kind) precsl_3 real(r_kind) pres real(r_kind) psaci real(r_kind) psaut real(r_kind) psm1 real(r_kind) psm2 real(r_kind) psm3 real(r_kind) q2 real(r_kind) qin real(r_kind) qk real(r_kind) qs real(r_kind) qs0 real(r_kind) rconde real(r_kind) rdt real(r_kind) rhci(km) real(r_kind) rke2 real(r_kind) rprs real(r_kind) rq real(r_kind) rqkll real(r_kind) t2 real(r_kind) tem real(r_kind) tem1 real(r_kind) tem2 real(r_kind) tem3 real(r_kind) tem4 real(r_kind) tem5 real(r_kind) term real(r_kind) term1 real(r_kind) term2r real(r_kind) term2s real(r_kind) term3r real(r_kind) term3s real(r_kind) term4 real(r_kind) term5 real(r_kind) term6 real(r_kind) tin real(r_kind) tmt0 real(r_kind) tmt0k real(r_kind) wmin real(r_kind) ww0 real(r_kind) ww1 real(r_kind) ww2 real(r_kind) wwn real(r_kind) zaodt real(r_kind) onem10,onem20 !---------------------------------------------- ! RESET LOCAL ADJOINT VARIABLES !---------------------------------------------- onem10 = 1.0e-10_r_kind onem20 = 1.0e-20_r_kind amaxcmr_ad = zero amaxcms_ad = zero amaxps_ad = zero amaxrq_ad = zero ccr_ad = zero cwmin_ad = zero cwmk_ad = zero cwmks_ad = zero erk_ad = zero err_ad = zero ers_ad = zero es_ad = zero expf_ad = zero fi_ad = zero ppr0_ad = zero ppr1_ad = zero ppr2_ad = zero pps0_ad = zero pps1_ad = zero pps2_ad = zero praut_ad = zero praut0_ad = zero precrk0_ad = zero precrk1_ad = zero precrl1k_ad = zero precrl_0_ad = zero precrl_1_ad = zero precrl_2_ad = zero precrl_3_ad = zero precsk0_ad = zero precsk1_ad = zero precsl1k_ad = zero precsl_0_ad = zero precsl_1_ad = zero precsl_2_ad = zero precsl_3_ad = zero psaci_ad = zero psaut_ad = zero psm1_ad = zero psm2_ad = zero psm3_ad = zero q2_ad = zero qin_ad = zero qk_ad = zero qs_ad = zero qs0_ad = zero rprs_ad = zero rq_ad = zero rqkll_ad = zero t2_ad = zero tem1_ad = zero tem2_ad = zero tem3_ad = zero tem4_ad = zero tem5_ad = zero term_ad = zero term1_ad = zero term2r_ad = zero term2s_ad = zero term3r_ad = zero term3s_ad = zero term4_ad = zero term5_ad = zero term6_ad = zero tin_ad = zero tmt0_ad = zero tmt0k_ad = zero ww0_ad = zero ww1_ad = zero ww2_ad = zero wwn_ad = zero !---------------------------------------------- ! ROUTINE BODY !---------------------------------------------- !---------------------------------------------- ! FUNCTION AND TAPE COMPUTATIONS !---------------------------------------------- rdt = one/dt zaodt = 800._r_kind*rdt csm1 = 5.e-8_r_kind*zaodt crs1 = 5.e-6_r_kind*zaodt crs2 = 6.666e-10_r_kind*zaodt cr = 0.0005_r_kind*zaodt aa2 = 0.00125_r_kind*zaodt rke2 = ke2*sqrt(rdt) dtcp = dt*rcp c00 = 0.0001_r_kind*dt do k = 1, km rhci(k) = one/rhc(k) end do !new comput0 = .false. do k = 1, km tem = 0.00001_r_kind*sl(k)*ps*0.01_r_kind if (cwm_in(k) > tem) then comput0 = .true. endif end do !new iwl1 = 0 precrl1k = zero precsl1k = zero rn_out = zero const = ps*(h1000*dt/grav) !new if (comput0) then !n do k = km, 1, -1 wmin = 0.00001_r_kind*sl(k)*ps*0.01_r_kind precrl_0 = precrl1k precsl_0 = precsl1k ! Save forward model values for use in adjoint model precrl1ki(k) = precrl_0 precsl1ki(k) = precsl_0 iwl1k(k) = iwl1 iwl = 0 tin = t_in(k) qin = q_in(k) cwmin = cwm_in(k) pres = ps*sl(k) if (precrl_0 > zero) then precrk0 = precrl_0 else precrk0 = zero endif if (precsl_0 > zero) then precsk0 = precsl_0 else precsk0 = zero endif if (cwmin > climit) then wwn = cwmin else wwn = climit endif if (wwn > climit .or. precrk0+precsk0 > zero) then comput = .true. else comput = .false. endif if (comput) then conde = const*del(k) condt = conde*rdt rconde = one/conde if (qin > epsq) then qk = qin else qk = epsq endif tmt0 = tin-ttp call fpvsx_ad(tin,es,dum1,dum2,.false.) qs0 = eps*es/(pres+epsm1*es) if (qs0 > epsq) then qs = qs0 else qs = epsq endif if (tmt0 < (-15._r_kind)) then fi = qk-rhc(k)*qs if (fi > zero .or. wwn > climit) then iwl = 1 else iwl = 0 endif else if (tmt0 >= zero) then iwl = 0 else iwl = 0 if (iwl1 == 1 .and. wwn > climit) then iwl = 1 endif endif if (qs <= onem10) then rq = zero else rq = qk/qs endif if (rq < rhc(k)) then ccr = zero else if (rq >= one) then ccr = one else if (rq < one) then rqkll = rq else rqkll = one endif ccr = one-sqrt((one-rqkll)/(one-rhc(k))) endif if (ccr > zero) then ww0 = cwmin if (ww0 > zero) then cwmk = ww0 else cwmk = zero endif if (iwl == 1) then term1 = cwmk-wmin if (term1 > zero) then amaxcms = term1 else amaxcms = zero endif expf = dt*exp(0.025_r_kind*tmt0) term2s = 0.0004_r_kind*expf*amaxcms if (term2s < cwmk) then psaut = term2s else psaut = cwmk endif ww1 = ww0-psaut if (ww1 > zero) then cwmks = ww1 else cwmks = zero endif term3s = aa2*expf*precsl_0*cwmks if (term3s < cwmks) then psaci = term3s else psaci = cwmks endif ww2 = ww1-psaci precsl_1 = precsl_0+(ww0-ww2)*condt precrl_1 = precrl_0 else amaxcmr = cwmk tem1 = precsl_0+precrl_0 term2r = 268._r_kind-tin if (term2r > zero) then term3r = term2r else term3r = zero endif if (term3r < 20._r_kind) then term4 = term3r else term4 = 20._r_kind endif tem2 = term4 tem3 = (one+h300*sqrt(tem1*rdt))*(one+half*sqrt(tem2)) if (ccr > 0.01_r_kind) then term5 = ccr else term5 = 0.01_r_kind endif tem4 = amaxcmr*cmr*tem3/term5 term6 = tem4*tem4 if (term6 < 50._r_kind) then tem5 = term6 else tem5 = 50._r_kind endif praut = c00*tem3*amaxcmr*(one-exp(-tem5)) if (praut < cwmk) then praut0 = praut else praut0 = cwmk endif ww2 = ww0-praut0 precrl_1 = precrl_0+(ww0-ww2)*condt precsl_1 = precsl_0 endif else ww2 = cwmin precrl_1 = precrl_0 precsl_1 = precsl_0 endif if (tmt0 > (-30._r_kind)) then tmt0k = tmt0 else tmt0k = -30._r_kind endif if (precrl_1 > zero) then precrk1 = precrl_1 else precrk1 = zero endif if (precsl_1 > zero) then precsk1 = precsl_1 else precsk1 = zero endif term = rhc(k)-rq if (term > zero) then amaxrq = term*conde else amaxrq = zero endif ppr0 = rke2*amaxrq*sqrt(precrk1) if (tmt0 >= zero) then pps0 = zero else pps0 = (crs1+crs2*tmt0k)*amaxrq*precsk1*rhci(k) endif erk = precrk1+precsk1 if (rq >= onem10) then erk = amaxrq*qk*rdt/rq endif if (ppr0+pps0 > abs(erk)) then rprs = erk/(precrk1+precsk1) ppr1 = precrk1*rprs pps1 = precsk1*rprs else ppr1 = ppr0 pps1 = pps0 endif if (ppr1 < precrk1) then ppr2 = ppr1 else ppr2 = precrk1 endif if (pps1 < precsk1) then pps2 = pps1 else pps2 = precsk1 endif err = ppr2*rconde ers = pps2*rconde precrl_2 = precrl_1-ppr2 precsl_2 = precsl_1-pps2 if (tmt0 > zero) then if (precsl_2 > zero) then amaxps = precsl_2 else amaxps = zero endif psm1 = csm1*tmt0*tmt0*amaxps if (ww2 > zero) then psm2 = cws*cr*ww2*amaxps else psm2 = zero endif ppr0 = (psm1+psm2)*conde if (ppr0 > amaxps) then ppr1 = amaxps psm3 = amaxps*rconde else ppr1 = ppr0 psm3 = psm1 endif precrl_3 = precrl_2+ppr1 precsl_3 = precsl_2-ppr1 else psm3 = zero precrl_3 = precrl_2 precsl_3 = precsl_2 endif t2 = tin-dtcp*(hvap*err+hsub*ers+hfus*psm3) q2 = qin+dt*(err+ers) else t2 = tin q2 = qin precrl_3 = precrl_0 precsl_3 = precsl_0 ww2 = cwmin endif iwl1 = iwl if (precrl_3 > zero) then precrl1k = precrl_3 else precrl1k = zero endif if (precsl_3 > zero) then precsl1k = precsl_3 else precsl1k = zero endif if (ww2 < zero) then q_out(k) = q2+ww2 t_out(k) = t2-hvap*rcp*ww2 cwm_out(k) = zero else q_out(k) = q2 t_out(k) = t2 cwm_out(k) = ww2 endif end do rn_out = (precrl1k+precsl1k)*rrow !new - no compute case else do k = km, 1, -1 if (cwm_in(k) < zero) then q_out(k) = q_in(k)+cwm_in(k) t_out(k) = t_in(k)-hvap*rcp*cwm_in(k) cwm_out(k) = zero else q_out(k) = q_in(k) t_out(k) = t_in(k) cwm_out(k) = cwm_in(k) endif end do rn_out = zero endif if (.not.adjoint) return !---------------------------------------------- ! ADJOINT COMPUTATIONS !---------------------------------------------- if (comput0) then const = ps*(h1000*dt/grav) precrl1k_ad = precrl1k_ad+rn_out_ad*rrow precsl1k_ad = precsl1k_ad+rn_out_ad*rrow rn_out_ad = zero do k = 1, km ! Load saved forward model values precrl1k = precrl1ki(k) precsl1k = precsl1ki(k) iwl1 = iwl1k(k) ! Recompute forward model values wmin = 0.00001_r_kind*sl(k)*ps*0.01_r_kind precrl_0 = precrl1k precsl_0 = precsl1k tin = t_in(k) qin = q_in(k) cwmin = cwm_in(k) pres = ps*sl(k) if (precrl_0 > zero) then precrk0 = precrl_0 else precrk0 = zero endif if (precsl_0 > zero) then precsk0 = precsl_0 else precsk0 = zero endif if (cwmin > climit) then wwn = cwmin else wwn = climit endif if (wwn > climit .or. precrk0+precsk0 > zero) then comput = .true. else comput = .false. endif if (comput) then conde = const*del(k) condt = conde*rdt if (qin > epsq) then qk = qin else qk = epsq endif tmt0 = tin-ttp call fpvsx_ad(tin,es,dum1,dum2,.false.) qs0 = eps*es/(pres+epsm1*es) if (qs0 > epsq) then qs = qs0 else qs = epsq endif if (tmt0 < (-15._r_kind)) then fi = qk-rhc(k)*qs if (fi > zero .or. wwn > climit) then iwl = 1 else iwl = 0 endif else if (tmt0 >= zero) then iwl = 0 else iwl = 0 if (iwl1 == 1 .and. wwn > climit) then iwl = 1 endif endif if (qs <= onem10) then rq = zero else rq = qk/qs endif if (rq < rhc(k)) then ccr = zero else if (rq >= one) then ccr = one else if (rq < one) then rqkll = rq else rqkll = one endif ccr = one-sqrt((one-rqkll)/(one-rhc(k))) endif if (ccr > zero) then ww0 = cwmin if (ww0 > zero) then cwmk = ww0 else cwmk = zero endif if (iwl == 1) then term1 = cwmk-wmin if (term1 > zero) then amaxcms = term1 else amaxcms = zero endif expf = dt*exp(0.025_r_kind*tmt0) term2s = 0.0004_r_kind*expf*amaxcms if (term2s < cwmk) then psaut = term2s else psaut = cwmk endif ww1 = ww0-psaut if (ww1 > zero) then cwmks = ww1 else cwmks = zero endif term3s = aa2*expf*precsl_0*cwmks if (term3s < cwmks) then psaci = term3s else psaci = cwmks endif ww2 = ww1-psaci precsl_1 = precsl_0+(ww0-ww2)*condt precrl_1 = precrl_0 else amaxcmr = cwmk tem1 = precsl_0+precrl_0 term2r = 268._r_kind-tin if (term2r > zero) then term3r = term2r else term3r = zero endif if (term3r < 20._r_kind) then term4 = term3r else term4 = 20._r_kind endif tem2 = term4 tem3 = (one+h300*sqrt(tem1*rdt))*(one+half*sqrt(tem2)) if (ccr > 0.01_r_kind) then term5 = ccr else term5 = 0.01_r_kind endif tem4 = amaxcmr*cmr*tem3/term5 term6 = tem4*tem4 if (term6 < 50._r_kind) then tem5 = term6 else tem5 = 50._r_kind endif praut = c00*tem3*amaxcmr*(one-exp(-tem5)) if (praut < cwmk) then praut0 = praut else praut0 = cwmk endif ww2 = ww0-praut0 precrl_1 = precrl_0+(ww0-ww2)*condt precsl_1 = precsl_0 endif else ww2 = cwmin precrl_1 = precrl_0 precsl_1 = precsl_0 endif if (tmt0 > (-30._r_kind)) then tmt0k = tmt0 else tmt0k = -30._r_kind endif if (precrl_1 > zero) then precrk1 = precrl_1 else precrk1 = zero endif if (precsl_1 > zero) then precsk1 = precsl_1 else precsk1 = zero endif term = rhc(k)-rq if (term > zero) then amaxrq = term*conde else amaxrq = zero endif ! ppr0 = rke2*amaxrq*sqrt(precrk1) if (precrk1>zero) then ppr0 = rke2*amaxrq*sqrt(precrk1) else ppr0 = zero endif if (tmt0 >= zero) then pps0 = zero else pps0 = (crs1+crs2*tmt0k)*amaxrq*precsk1*rhci(k) endif erk = precrk1+precsk1 if (rq >= onem10) then erk = amaxrq*qk*rdt/rq endif if (ppr0+pps0 > abs(erk)) then rprs = erk/(precrk1+precsk1) ppr1 = precrk1*rprs pps1 = precsk1*rprs else ppr1 = ppr0 pps1 = pps0 endif if (ppr1 < precrk1) then ppr2 = ppr1 else ppr2 = precrk1 endif if (pps1 < precsk1) then pps2 = pps1 else pps2 = precsk1 endif precrl_2 = precrl_1-ppr2 precsl_2 = precsl_1-pps2 if (tmt0 > zero) then if (precsl_2 > zero) then amaxps = precsl_2 else amaxps = zero endif psm1 = csm1*tmt0*tmt0*amaxps if (ww2 > zero) then psm2 = cws*cr*ww2*amaxps else psm2 = zero endif ppr0 = (psm1+psm2)*conde if (ppr0 > amaxps) then ppr1 = amaxps else ppr1 = ppr0 endif precrl_3 = precrl_2+ppr1 precsl_3 = precsl_2-ppr1 else precrl_3 = precrl_2 precsl_3 = precsl_2 endif else precrl_3 = precrl_0 precsl_3 = precsl_0 ww2 = cwmin endif ! Adjoint model if (ww2 < zero) then cwm_out_ad(k) = zero t2_ad = t2_ad+t_out_ad(k) ww2_ad = ww2_ad-t_out_ad(k)*hvap*rcp t_out_ad(k) = zero q2_ad = q2_ad+q_out_ad(k) ww2_ad = ww2_ad+q_out_ad(k) q_out_ad(k) = zero else ww2_ad = ww2_ad+cwm_out_ad(k) cwm_out_ad(k) = zero t2_ad = t2_ad+t_out_ad(k) t_out_ad(k) = zero q2_ad = q2_ad+q_out_ad(k) q_out_ad(k) = zero endif if (precsl_3 > zero) then precsl_3_ad = precsl_3_ad+precsl1k_ad precsl1k_ad = zero else precsl1k_ad = zero endif if (precrl_3 > zero) then precrl_3_ad = precrl_3_ad+precrl1k_ad precrl1k_ad = zero else precrl1k_ad = zero endif if (comput) then rconde = one/conde err_ad = err_ad+q2_ad*dt ers_ad = ers_ad+q2_ad*dt qin_ad = qin_ad+q2_ad q2_ad = zero err_ad = err_ad-t2_ad*dtcp*hvap ers_ad = ers_ad-t2_ad*dtcp*hsub psm3_ad = psm3_ad-t2_ad*dtcp*hfus tin_ad = tin_ad+t2_ad t2_ad = zero if (tmt0 > zero) then ppr1_ad = ppr1_ad-precsl_3_ad precsl_2_ad = precsl_2_ad+precsl_3_ad precsl_3_ad = zero ppr1_ad = ppr1_ad+precrl_3_ad precrl_2_ad = precrl_2_ad+precrl_3_ad precrl_3_ad = zero if (ppr0 > amaxps) then amaxps_ad = amaxps_ad+psm3_ad*rconde psm3_ad = zero amaxps_ad = amaxps_ad+ppr1_ad ppr1_ad = zero else psm1_ad = psm1_ad+psm3_ad psm3_ad = zero ppr0_ad = ppr0_ad+ppr1_ad ppr1_ad = zero endif psm1_ad = psm1_ad+ppr0_ad*conde psm2_ad = psm2_ad+ppr0_ad*conde ppr0_ad = zero if (ww2 > zero) then amaxps_ad = amaxps_ad+psm2_ad*cws*cr*ww2 ww2_ad = ww2_ad+psm2_ad*cws*cr*amaxps psm2_ad = zero else psm2_ad = zero endif amaxps_ad = amaxps_ad+psm1_ad*csm1*tmt0*tmt0 tmt0_ad = tmt0_ad+2*psm1_ad*csm1*tmt0*amaxps psm1_ad = zero if (precsl_2 > zero) then precsl_2_ad = precsl_2_ad+amaxps_ad amaxps_ad = zero else amaxps_ad = zero endif else precsl_2_ad = precsl_2_ad+precsl_3_ad precsl_3_ad = zero precrl_2_ad = precrl_2_ad+precrl_3_ad precrl_3_ad = zero psm3_ad = zero endif pps2_ad = pps2_ad-precsl_2_ad precsl_1_ad = precsl_1_ad+precsl_2_ad precsl_2_ad = zero ppr2_ad = ppr2_ad-precrl_2_ad precrl_1_ad = precrl_1_ad+precrl_2_ad precrl_2_ad = zero pps2_ad = pps2_ad+ers_ad*rconde ers_ad = zero ppr2_ad = ppr2_ad+err_ad*rconde err_ad = zero if (pps1 < precsk1) then pps1_ad = pps1_ad+pps2_ad pps2_ad = zero else precsk1_ad = precsk1_ad+pps2_ad pps2_ad = zero endif ppr0 = rke2*amaxrq*sqrt(precrk1) if (ppr0+pps0 > abs(erk)) then ppr1 = precrk1*rprs else ppr1 = ppr0 endif if (ppr1 < precrk1) then ppr1_ad = ppr1_ad+ppr2_ad ppr2_ad = zero else precrk1_ad = precrk1_ad+ppr2_ad ppr2_ad = zero endif ppr0 = rke2*amaxrq*sqrt(precrk1) if (ppr0+pps0 > abs(erk)) then precsk1_ad = precsk1_ad+pps1_ad*rprs rprs_ad = rprs_ad+pps1_ad*precsk1 pps1_ad = zero precrk1_ad = precrk1_ad+ppr1_ad*rprs rprs_ad = rprs_ad+ppr1_ad*precrk1 ppr1_ad = zero erk_ad = erk_ad+rprs_ad/(precrk1+precsk1) precrk1_ad = precrk1_ad-rprs_ad*(erk/((precrk1+precsk1)* & (precrk1+precsk1))) precsk1_ad = precsk1_ad-rprs_ad*(erk/((precrk1+precsk1)* & (precrk1+precsk1))) rprs_ad = zero else pps0_ad = pps0_ad+pps1_ad pps1_ad = zero ppr0_ad = ppr0_ad+ppr1_ad ppr1_ad = zero endif if (rq >= onem10) then amaxrq_ad = amaxrq_ad+erk_ad*(qk*rdt/rq) qk_ad = qk_ad+erk_ad*(amaxrq*rdt/rq) rq_ad = rq_ad-erk_ad*(amaxrq*qk*rdt/(rq*rq)) erk_ad = zero endif precrk1_ad = precrk1_ad+erk_ad precsk1_ad = precsk1_ad+erk_ad erk_ad = zero if (tmt0 >= zero) then pps0_ad = zero else amaxrq_ad = amaxrq_ad+pps0_ad*(crs1+crs2*tmt0k)*precsk1* & rhci(k) precsk1_ad = precsk1_ad+pps0_ad*(crs1+crs2*tmt0k)*amaxrq* & rhci(k) tmt0k_ad = tmt0k_ad+pps0_ad*crs2*amaxrq*precsk1*rhci(k) pps0_ad = zero endif ! amaxrq_ad = amaxrq_ad+ppr0_ad*rke2*sqrt(precrk1) ! precrk1_ad = precrk1_ad+ppr0_ad*rke2*amaxrq*one/(two* & ! sqrt(precrk1)) ! ppr0_ad = zero if (precrk1>zero) then amaxrq_ad = amaxrq_ad+ppr0_ad*rke2*sqrt(precrk1) precrk1_ad = precrk1_ad+ppr0_ad*rke2*amaxrq*one/(two* & sqrt(precrk1)) ppr0_ad = zero else ppr0_ad = zero endif if (term > zero) then term_ad = term_ad+amaxrq_ad*conde amaxrq_ad = zero else amaxrq_ad = zero endif rq_ad = rq_ad-term_ad term_ad = zero if (precsl_1 > zero) then precsl_1_ad = precsl_1_ad+precsk1_ad precsk1_ad = zero else precsk1_ad = zero endif if (precrl_1 > zero) then precrl_1_ad = precrl_1_ad+precrk1_ad precrk1_ad = zero else precrk1_ad = zero endif if (tmt0 > (-30._r_kind)) then tmt0_ad = tmt0_ad+tmt0k_ad tmt0k_ad = zero else tmt0k_ad = zero endif if (ccr > zero) then if (iwl == 1) then precrl_0_ad = precrl_0_ad+precrl_1_ad precrl_1_ad = zero precsl_0_ad = precsl_0_ad+precsl_1_ad ww0_ad = ww0_ad+precsl_1_ad*condt ww2_ad = ww2_ad-precsl_1_ad*condt precsl_1_ad = zero psaci_ad = psaci_ad-ww2_ad ww1_ad = ww1_ad+ww2_ad ww2_ad = zero if (term3s < cwmks) then term3s_ad = term3s_ad+psaci_ad psaci_ad = zero else cwmks_ad = cwmks_ad+psaci_ad psaci_ad = zero endif cwmks_ad = cwmks_ad+term3s_ad*aa2*expf*precsl_0 expf_ad = expf_ad+term3s_ad*aa2*precsl_0*cwmks precsl_0_ad = precsl_0_ad+term3s_ad*aa2*expf*cwmks term3s_ad = zero if (ww1 > zero) then ww1_ad = ww1_ad+cwmks_ad cwmks_ad = zero else cwmks_ad = zero endif psaut_ad = psaut_ad-ww1_ad ww0_ad = ww0_ad+ww1_ad ww1_ad = zero if (term2s < cwmk) then term2s_ad = term2s_ad+psaut_ad psaut_ad = zero else cwmk_ad = cwmk_ad+psaut_ad psaut_ad = zero endif amaxcms_ad = amaxcms_ad+0.0004_r_kind*term2s_ad*expf expf_ad = expf_ad+0.0004_r_kind*term2s_ad*amaxcms term2s_ad = zero tmt0_ad = tmt0_ad+0.025_r_kind*expf_ad*dt*exp(0.025_r_kind*tmt0) expf_ad = zero if (term1 > zero) then term1_ad = term1_ad+amaxcms_ad amaxcms_ad = zero else amaxcms_ad = zero endif cwmk_ad = cwmk_ad+term1_ad term1_ad = zero else precsl_0_ad = precsl_0_ad+precsl_1_ad precsl_1_ad = zero precrl_0_ad = precrl_0_ad+precrl_1_ad ww0_ad = ww0_ad+precrl_1_ad*condt ww2_ad = ww2_ad-precrl_1_ad*condt precrl_1_ad = zero praut0_ad = praut0_ad-ww2_ad ww0_ad = ww0_ad+ww2_ad ww2_ad = zero if (praut < cwmk) then praut_ad = praut_ad+praut0_ad praut0_ad = zero else cwmk_ad = cwmk_ad+praut0_ad praut0_ad = zero endif amaxcmr_ad = amaxcmr_ad+praut_ad*c00*tem3*(one-exp(-tem5)) tem3_ad = tem3_ad+praut_ad*c00*amaxcmr*(one-exp(-tem5)) tem5_ad = tem5_ad+praut_ad*c00*tem3*amaxcmr*exp(-tem5) praut_ad = zero if (term6 < 50._r_kind) then term6_ad = term6_ad+tem5_ad tem5_ad = zero else tem5_ad = zero endif tem4_ad = tem4_ad+2*term6_ad*tem4 term6_ad = zero amaxcmr_ad = amaxcmr_ad+tem4_ad*(cmr*tem3/term5) tem3_ad = tem3_ad+tem4_ad*(amaxcmr*cmr/term5) term5_ad = term5_ad-tem4_ad*(amaxcmr*cmr*tem3/(term5*term5)) tem4_ad = zero if (ccr > 0.01_r_kind) then ccr_ad = ccr_ad+term5_ad term5_ad = zero else term5_ad = zero endif !old - code below can result in divide by zero ! tem1_ad = tem1_ad+tem3_ad*h300*one/(two*sqrt(tem1*rdt))* & ! rdt*(one+half*sqrt(tem2)) ! tem2_ad = tem2_ad+tem3_ad*(one+h300*sqrt(tem1*rdt))*half* & ! (one/(two*sqrt(tem2))) if (abs(tem1)>onem20) then tem1_ad = tem1_ad+tem3_ad*h300*one/ & (two*sqrt(tem1*rdt))*rdt*(one+half*sqrt(tem2)) else tem1_ad = tem1_ad+tem3_ad*h300*one/ & (two*sqrt(onem20*rdt))*rdt*(one+half*sqrt(tem2)) endif if (abs(tem2)>onem20) then tem2_ad = tem2_ad+tem3_ad*(one+h300*sqrt(tem1*rdt))* & half*(one/(two*sqrt(tem2))) else tem2_ad = tem2_ad+tem3_ad*(one+h300*sqrt(tem1*rdt))* & half*(one/(two*sqrt(onem20))) endif !new tem3_ad = zero term4_ad = term4_ad+tem2_ad tem2_ad = zero if (term3r < 20._r_kind) then term3r_ad = term3r_ad+term4_ad term4_ad = zero else term4_ad = zero endif if (term2r > zero) then term2r_ad = term2r_ad+term3r_ad term3r_ad = zero else term3r_ad = zero endif tin_ad = tin_ad-term2r_ad term2r_ad = zero precrl_0_ad = precrl_0_ad+tem1_ad precsl_0_ad = precsl_0_ad+tem1_ad tem1_ad = zero cwmk_ad = cwmk_ad+amaxcmr_ad amaxcmr_ad = zero endif if (ww0 > zero) then ww0_ad = ww0_ad+cwmk_ad cwmk_ad = zero else cwmk_ad = zero endif cwmin_ad = cwmin_ad+ww0_ad ww0_ad = zero else precsl_0_ad = precsl_0_ad+precsl_1_ad precsl_1_ad = zero precrl_0_ad = precrl_0_ad+precrl_1_ad precrl_1_ad = zero cwmin_ad = cwmin_ad+ww2_ad ww2_ad = zero endif if (rq < rhc(k)) then ccr_ad = zero else if (rq >= one) then ccr_ad = zero else rqkll_ad = rqkll_ad+ccr_ad*(one/(two*sqrt((one-rqkll)/ & (one-rhc(k))))/(one-rhc(k))) ccr_ad = zero if (rq < one) then rq_ad = rq_ad+rqkll_ad rqkll_ad = zero else rqkll_ad = zero endif endif if (qs <= onem10) then rq_ad = zero else qk_ad = qk_ad+rq_ad/qs qs_ad = qs_ad-rq_ad*(qk/(qs*qs)) rq_ad = zero endif if (tmt0 < (-15._r_kind)) then qk_ad = qk_ad+fi_ad qs_ad = qs_ad-fi_ad*rhc(k) fi_ad = zero endif if (qs0 > epsq) then qs0_ad = qs0_ad+qs_ad qs_ad = zero else qs_ad = zero endif es_ad = es_ad+qs0_ad*(eps/(pres+epsm1*es)-eps*es*epsm1/ & ((pres+epsm1*es)*(pres+epsm1*es))) qs0_ad = zero call fpvsx_ad( tin,es,tin_ad,es_ad,adjoint ) es_ad = zero tin_ad = tin_ad+tmt0_ad tmt0_ad = zero if (qin > epsq) then qin_ad = qin_ad+qk_ad qk_ad = zero else qk_ad = zero endif else cwmin_ad = cwmin_ad+ww2_ad ww2_ad = zero precsl_0_ad = precsl_0_ad+precsl_3_ad precsl_3_ad = zero precrl_0_ad = precrl_0_ad+precrl_3_ad precrl_3_ad = zero qin_ad = qin_ad+q2_ad q2_ad = zero tin_ad = tin_ad+t2_ad t2_ad = zero endif if (cwmin > climit) then cwmin_ad = cwmin_ad+wwn_ad wwn_ad = zero else wwn_ad = zero endif if (precsl_0 > zero) then precsl_0_ad = precsl_0_ad+precsk0_ad precsk0_ad = zero else precsk0_ad = zero endif if (precrl_0 > zero) then precrl_0_ad = precrl_0_ad+precrk0_ad precrk0_ad = zero else precrk0_ad = zero endif cwm_in_ad(k) = cwm_in_ad(k)+cwmin_ad cwmin_ad = zero q_in_ad(k) = q_in_ad(k)+qin_ad qin_ad = zero t_in_ad(k) = t_in_ad(k)+tin_ad tin_ad = zero precsl1k_ad = precsl1k_ad+precsl_0_ad precsl_0_ad = zero precrl1k_ad = precrl1k_ad+precrl_0_ad precrl_0_ad = zero end do else rn_out_ad = zero do k = km, 1, -1 if (cwm_in(k) < zero) then cwm_out_ad(k) = zero cwm_in_ad(k) = cwm_in_ad(k)-t_out_ad(k)*hvap*rcp t_in_ad(k) = t_in_ad(k)+t_out_ad(k) t_out_ad(k) = zero cwm_in_ad(k) = cwm_in_ad(k)+q_out_ad(k) q_in_ad(k) = q_in_ad(k)+q_out_ad(k) q_out_ad(k) = zero else cwm_in_ad(k) = cwm_in_ad(k)+cwm_out_ad(k) cwm_out_ad(k) = zero t_in_ad(k) = t_in_ad(k)+t_out_ad(k) t_out_ad(k) = zero q_in_ad(k) = q_in_ad(k)+q_out_ad(k) q_out_ad(k) = zero endif end do endif rn_out_ad = zero precsl1k_ad = zero precrl1k_ad = zero return end subroutine precpd_ad