11 use wv_saturation,
only : epsqs, ttrice, hlatv, hlatf, pcf, rgasv
18 use mapl_constantsmod,
only: mapl_tice , mapl_cp , mapl_grav ,
19 & mapl_alhs , mapl_alhl , mapl_alhf , mapl_rgas , mapl_h2omw,
20 & mapl_airmw, mapl_rvap , mapl_pi , mapl_r8 , mapl_r4
21 use mapl_basemod,
only: mapl_undef
24 use physcons, mapl_tice => con_t0c, mapl_grav => con_g,
25 & mapl_cp => con_cp, mapl_alhl => con_hvap,
26 & mapl_alhf => con_hfus, mapl_pi => con_pi,
27 & mapl_rgas => con_rd, mapl_rvap => con_rv
46 integer :: nsmax, disable_rad, icefrpwr
47 &, fr_ls_wat, fr_ls_ice, fr_an_wat, fr_an_ice
89 real :: turnrhcrit_upper
90 real :: min_ri, max_ri, min_rl, max_rl, ri_anv
93 real,
parameter :: t_ice_max = mapl_tice
94 real,
parameter :: rho_w = 1.0e3
95 real,
parameter :: min_cld_frac = 1.0e-8
96 real,
parameter :: mapl_alhs = mapl_alhl+mapl_alhf
98 real,
parameter :: alhlbcp = mapl_alhl/mapl_cp
99 &, alhfbcp = mapl_alhf/mapl_cp
100 &, alhsbcp = alhlbcp+alhfbcp
104 real :: omeps, trinv, t_ice_denom
115 &, QLW_LS_dev, QLW_AN_dev, QIW_LS_dev &
116 &, QIW_AN_dev, ANVFRC_dev, CLDFRC_dev &
117 &, PHYSPARAMS, SCLMFDFR &
120 &, CNV_NDROP_dev, CNV_NICE_dev, SCICE_dev &
121 &, NCPL_dev, NCPI_dev, PFRZ_dev &
122 &, lprnt, ipr, rhc, pdfflag, qc_min )
125 integer,
intent(in ) :: irun, lm, pdfflag
126 real,
intent(in ) :: dt, alf_fac, qc_min(2)
127 real,
intent(in ),
dimension(IRUN, LM) :: pp_dev
128 real,
intent(in ),
dimension(IRUN,0:LM) :: ppe_dev
131 real,
intent(in ),
dimension(IRUN, LM) :: qlwdtr_dev
136 real,
intent(in ),
dimension(IRUN, LM) :: rhc
137 real,
intent(inout),
dimension(IRUN, LM) :: th_dev
138 real,
intent(inout),
dimension(IRUN, LM) :: q_dev
139 real,
intent(inout),
dimension(IRUN, LM) :: qlw_ls_dev
140 real,
intent(inout),
dimension(IRUN, LM) :: qlw_an_dev
141 real,
intent(inout),
dimension(IRUN, LM) :: qiw_ls_dev
142 real,
intent(inout),
dimension(IRUN, LM) :: qiw_an_dev
143 real,
intent(inout),
dimension(IRUN, LM) :: anvfrc_dev
144 real,
intent(inout),
dimension(IRUN, LM) :: cldfrc_dev
148 real,
intent(in ),
dimension(58 ) :: physparams
149 real,
intent(in ) :: sclmfdfr
162 real,
intent( out),
dimension(IRUN, LM) :: alpht_dev
167 real,
intent(inout),
dimension(IRUN, LM) :: cnv_fice_dev
168 real,
intent(inout),
dimension(IRUN, LM) :: cnv_ndrop_dev
169 real,
intent(inout),
dimension(IRUN, LM) :: cnv_nice_dev
170 real,
intent(inout),
dimension(IRUN, LM) :: scice_dev
171 real,
intent(inout),
dimension(IRUN, LM) :: ncpl_dev
172 real,
intent(inout),
dimension(IRUN, LM) :: ncpi_dev
173 real,
intent(out),
dimension(IRUN, LM) :: pfrz_dev
188 integer :: i , j , k , l
192 real :: mass, imass, totfrc, temp, alpha
193 &, dti, tx1, tend, fqi
205 logical :: use_autoconv_timescale
207 real,
parameter :: rl_cub = 1.0e-15, ri_cub = 6.4e-14
214 cnv_beta = physparams(1)
215 anv_beta = physparams(2)
216 ls_beta = physparams(3)
219 lwcrit = physparams(6)
220 c_acc = physparams(7)
221 c_ev_r = physparams(8)
222 c_ev_s = physparams(56)
223 cldvol2frc = physparams(9)
224 rhsup_ice = physparams(10)
225 shr_evap_fac = physparams(11)
226 min_cld_water = physparams(12)
227 cld_evp_eff = physparams(13)
228 nsmax = int( physparams(14) )
229 ls_sdqv2 = physparams(15)
230 ls_sdqv3 = physparams(16)
231 ls_sdqvt1 = physparams(17)
232 anv_sdqv2 = physparams(18)
233 anv_sdqv3 = physparams(19)
234 anv_sdqvt1 = physparams(20)
235 anv_to_ls = physparams(21)
236 n_warm = physparams(22)
237 n_ice = physparams(23)
238 n_anvil = physparams(24)
239 n_pbl = physparams(25)
240 disable_rad = int( physparams(26) )
241 anv_icefall_c = physparams(28)
242 ls_icefall_c = physparams(29)
243 revap_off_p = physparams(30)
244 cnvenvfc = physparams(31)
245 wrhodep = physparams(32)
246 t_ice_all = physparams(33) + mapl_tice
247 cnviceparam = physparams(34)
248 icefrpwr = int( physparams(35) + .001 )
249 cnvddrfc = physparams(36)
250 anvddrfc = physparams(37)
251 lsddrfc = physparams(38)
255 turnrhcrit = physparams(45) * 0.001
257 fr_ls_wat = int( physparams(47) )
258 fr_ls_ice = int( physparams(48) )
259 fr_an_wat = int( physparams(49) )
260 fr_an_ice = int( physparams(50) )
261 min_rl = physparams(51)
262 min_ri = physparams(52)
263 max_rl = physparams(53)
264 max_ri = physparams(54)
265 ri_anv = physparams(55)
269 turnrhcrit_upper = physparams(58) * 0.001
271 use_autoconv_timescale = .false.
273 t_ice_denom = 1.0 / (t_ice_max-t_ice_all)
275 run_loop:
DO i = 1, irun
340 mass = (ppe_dev(i,k) - ppe_dev(i,k-1)) * (100./mapl_grav)
353 if (temp < t_ice_all)
then
355 elseif (temp > t_ice_max)
then
358 fqi = cnv_fice_dev(i,k)
360 tx1 = (1.0-fqi)*qlwdtr_dev(i,k)
361 if (tx1 > 0.0 .and. cnv_ndrop_dev(i,k) <= 0.0)
then
362 cnv_ndrop_dev(i,k) = tx1 / ( 1.333 * mapl_pi *rl_cub*997.0)
365 tx1 = fqi*qlwdtr_dev(i,k)
366 if (tx1 > 0.0 .and. cnv_nice_dev(i,k) <= 0.0)
then
367 cnv_nice_dev(i,k) = tx1 / ( 1.333 * mapl_pi *ri_cub*500.0)
371 ncpl_dev(i,k) = max(ncpl_dev(i,k)+cnv_ndrop_dev(i,k)*tx1,0.0)
372 ncpi_dev(i,k) = max(ncpi_dev(i,k)+cnv_nice_dev(i,k)*tx1,0.0)
415 alpha = max(1.0e-4, 1.0-rhc(i,k)) * alf_fac
416 alpht_dev(i,k) = alpha
423 call pfreezing (alpha , pp_dev(i,k) , temp , q_dev(i,k),
424 & qlw_ls_dev(i,k), qlw_an_dev(i,k),
425 & qiw_ls_dev(i,k), qiw_an_dev(i,k),
426 & scice_dev(i,k) , cldfrc_dev(i,k),
427 & anvfrc_dev(i,k), pfrz_dev(i,k), pdfflag)
560 & qiw_ls_dev(i,k), cldfrc_dev(i,k), qlw_an_dev(i,k),
561 & qiw_an_dev(i,k), anvfrc_dev(i,k), ncpl_dev(i, k),
562 & ncpi_dev(i, k), qc_min)
584 subroutine pdf_spread (K, LM, PP, ALPHA, ALPHT_DIAG, FRLAND, rhc)
586 integer,
intent(in) :: k,lm
587 real,
intent(in) :: PP, FRLAND, rhc
589 real,
intent(out) :: ALPHA, ALPHT_DIAG
592 real,
parameter :: slope = 0.02, slope_up = 0.02
594 real :: aux1, aux2, maxalpha
599 aux1 = min(max((pp - turnrhcrit)/slope, -20.0), 20.0)
600 aux2 = min(max((turnrhcrit_upper - pp)/slope_up, -20.0), 20.0)
602 if (frland > 0.05)
then
604 aux1 = 1.0 / (1.0+exp(aux1+aux1))
606 aux1 = 2.0 / (1.0+exp(aux1+aux1))
609 aux2 = 1.0 / (1.0+exp(aux2))
611 alpha = max(1.0e-4, min(0.3, maxalpha*aux1*aux2))
624 real,
intent(in) :: qc_min(2)
625 real,
intent(inout) :: te,qv,qlc,cf,qla,af,qic,qia, nl, ni
629 real,
parameter :: nmin = 1.0, cfmin = 1.0e-5
630 &, ri_cub = 6.4e-14, rl_cub = 1.0e-15
633 if (af <= cfmin)
then
635 te = te - alhlbcp*qla - alhsbcp*qia
640 if ( cf <= cfmin)
then
642 te = te - alhlbcp*qlc - alhsbcp*qic
649 if (qlc <= qc_min(1))
then
651 te = te - alhlbcp*qlc
655 if (qic <= qc_min(2))
then
657 te = te - alhsbcp*qic
661 if (qla <= qc_min(1))
then
663 te = te - alhlbcp*qla
667 if (qia <= qc_min(2))
then
669 te = te - alhsbcp*qia
673 if (qla+qia <= qc_min(1))
then
675 te = te - alhlbcp*qla - alhsbcp*qia
681 if (qlc+qic <= qc_min(1))
then
683 te = te - alhlbcp*qlc - alhsbcp*qic
689 if (qla+qlc <= qc_min(1))
then
691 elseif (nl <= nmin)
then
692 nl = max((qla+qlc)/( fourb3 * mapl_pi *rl_cub*997.0), nmin)
695 if (qia+qic <= qc_min(2))
then
697 elseif (ni <= nmin)
then
698 ni = max((qia+qic)/( fourb3 * mapl_pi *ri_cub*500.0), nmin)
710 & PDFSHAPE, PL, QV, QCl, QAl, &
711 & QCi, QAi, TE, CF, AF, &
715 integer,
intent(in) :: irun, lm, pdfshape
716 real,
intent(in) :: dt, qc_min(2)
717 real,
intent(in),
dimension(irun,lm) :: alpha, pl
719 real,
intent(inout),
dimension(irun,lm) :: te, qv, qcl, qci &
720 &, CF, QAl, QAi, AF, NI, NL, SCICE
723 real :: cfo, pl100, qt, dq, qsx, qcx, qc, qa &
724 &, QX, QSLIQ, QSICE, CFALL, DQx, FQA, tem
726 real :: esl, esi, esn
732 if (qv(i,k) > 1.0e-6)
then
733 qc = qcl(i,k) + qci(i,k)
734 qa = qal(i,k) + qai(i,k)
736 if(qc <= 0.) cf(i,k) = 0.
737 if(qa <= 0.) af(i,k) = 0.
740 cfall = af(i,k) + cf(i,k)
750 esl = min(
fpvsl(te(i,k)), pl100)
751 qsliq = min(epsqs*esl/(pl100-omeps*esl), 1.)
752 esi = min(
fpvsi(te(i,k)), pl100)
753 qsice = min(epsqs*esi/(pl100-omeps*esi), 1.)
755 qsx = ( (qcl(i,k)+qal(i,k))*qsliq
756 & + (qci(i,k)+qai(i,k))*qsice ) *tem
759 esn = min(fpvs(te(i,k)), pl100)
760 qsx = min(epsqs*esn/(pl100-omeps*esn), 1.)
765 qx = qt - qsx*scice(i,k)
769 if (qcx > qc_min(1))
then
771 cfo = 1.0 / (1.0 + sqrt(1.0-qx/qcx) )
781 cfall = min(1.0, max(cfo, 0.0))
783 af(i,k) = cfall * fqa
784 cf(i,k) = cfall - af(i,k)
788 call hystpdf( dt, alpha(i,k), pdfshape, qc_min, pl(i,k)
789 &, qv(i,k), qcl(i,k), qal(i,k), qci(i,k)
790 &, qai(i,k), te(i,k), cf(i,k), af(i,k)
791 &, scice(i,k), ni(i,k), nl(i,k))
795 if(qcl(i,k)+qci(i,k) <= 0.0) cf(i,k) = 0.0
796 if(qal(i,k)+qai(i,k) <= 0.0) af(i,k) = 0.0
811 subroutine hystpdf( DT, ALPHA, PDFSHAPE, qc_min, PL, QV, QCl, QAl&
812 &, QCi, QAi, TE, CF, AF, SCICE, NI, NL)
815 real,
intent(in) :: DT, ALPHA, PL, qc_min(2)
816 integer,
intent(in) :: pdfshape
817 real,
intent(inout) :: TE, QV, QCl, QCi, CF, QAl, QAi, AF, &
820 integer,
parameter :: nmax=10
822 real :: QCO, QVO, CFO, QAO, TAU
823 real :: QT, QMX, QMN, DQ, QVtop, sigmaqt1, sigmaqt2, qsnx
826 real :: QSx, DQSx, QS, DQs
827 &, tep, qsp, cfp, qvp, qcp
828 &, ten, qsn, cfn, qvn, qcn
831 real :: QCx, QVx, CFx, QAx, QC, QA, fQi
832 &, dqai, dqal, dqci, dqcl
835 real :: QX, QSLIQ, QSICE, DQx, pl100, tmpARR
836 &, alhx, dqcall, esn, desdt, tc, hltalt, tterm
855 if (te <= t_ice_all)
then
857 elseif (te >= t_ice_max)
then
860 fqi = (1.0 - (te-t_ice_all)*t_ice_denom) ** icefrpwr
865 esn = min(fpvs(te), pl100)
866 qsx = min(epsqs*esn/(pl100-omeps*esn), 1.)
870 if (te < mapl_tice)
then
871 hltalt = hlatv + hlatf * min(-tc*trinv,1.0)
873 hltalt = hlatv - 2369.0*tc
875 if (tc >= -ttrice .and. tc < 0.0)
then
876 tterm = pcf(1) + tc*(pcf(2) + tc*(pcf(3)
877 & + tc*(pcf(4) + tc*pcf(5))))
881 desdt = hltalt*esn/(rgasv*te*te) + tterm*trinv
882 dqsx = qsx*pl100*desdt/(esn*(pl100-omeps*esn))
888 tmparr = 1.0 / (1.0-af)
895 qvx = (qv - qsx*af) * tmparr
921 if(pdfshape < 2)
then
925 elseif(pdfshape == 2)
then
928 elseif(pdfshape == 4)
then
929 sigmaqt1 = max(alpha/sqrt(3.0), 0.001)
931 write(0,*)
' Aborting : invalid pdfshape=',pdfshape
938 call pdffrac(pdfshape,qt,sigmaqt1,sigmaqt2,qsnx,cfn)
939 call pdfcondensate(pdfshape,qt,sigmaqt1,sigmaqt2,qsnx,qcn,cfn)
954 alhx = (1.0-fqi)*alhlbcp + fqi*alhsbcp
956 if(pdfshape == 1)
then
957 qcn = qcp + (qcn- qcp)
958 & / (1.0 - (cfn*(alpha-1.0) - qcn/qsn) *dqs*alhx)
959 elseif(pdfshape == 2)
then
960 if (n < nmax) qcn = qcp + ( qcn - qcp ) * 0.5
963 qvn = qvp - (qcn - qcp)
964 ten = tep + alhx * ((qcn-qcp)*(1.0-af) + (qao-qax)*af)
966 if (abs(ten-tep) < 0.00001)
exit
968 if (ten <= t_ice_all)
then
970 elseif (ten >= t_ice_max)
then
973 fqi = (1.0 - (te-t_ice_all)*t_ice_denom) ** icefrpwr
979 esn = min(fpvs(ten), pl100)
980 qsn = min(epsqs*esn/(pl100-omeps*esn), 1.0)
984 if (ten < mapl_tice)
then
985 hltalt = hlatv + hlatf * min(-tc*trinv,1.0)
987 hltalt = hlatv - 2369.0*tc
989 if (tc >= -ttrice .and. tc < 0.0)
then
990 tterm = pcf(1) + tc*(pcf(2) + tc*(pcf(3)
991 & + tc*(pcf(4) + tc*pcf(5))))
995 desdt = hltalt*esn/(rgasv*ten*ten) + tterm*trinv
996 dqs = qsn*pl100*desdt / (esn*(pl100-omeps*esn))
1016 qco = qco * (1.0-af)
1029 qao = max(qt-qsx, 0.)
1043 dqcl = max(qcx, -qcl)
1044 dqci = max(qcx-dqcl, -qci)
1056 dqal = max(qax, -qal)
1057 dqai = max(qax-dqal, -qai)
1064 if ( af < 1.e-5 )
then
1069 if ( cf < 1.e-5 )
then
1080 qv = qv - ( dqai+dqci+dqal+dqcl)
1083 te = te + (alhlbcp * (dqai+dqci+dqal+dqcl)
1084 & + alhfbcp * (dqai+dqci))
1088 if ( qao <= 0. )
then
1090 te = te - alhsbcp*qai - alhlbcp*qal
1104 subroutine pdffrac (flag,qtmean,sigmaqt1,sigmaqt2,qstar,clfrac)
1109 real :: qtmean, sigmaqt1, sigmaqt2, qstar, clfrac
1111 real :: qtmode, qtmin, qtmax, qtmedian, aux
1114 aux = qtmean + sigmaqt1 - qstar
1118 if(sigmaqt1 > 0.0)
then
1119 clfrac = min(0.5*aux/sigmaqt1, 1.0)
1124 elseif(flag == 2)
then
1125 qtmode = qtmean + (sigmaqt1-sigmaqt2)/3.
1126 qtmin = min(qtmode-sigmaqt1,0.)
1127 qtmax = qtmode + sigmaqt2
1128 if(qtmax < qstar)
then
1130 elseif ( (qtmode <= qstar).and.(qstar < qtmax) )
then
1131 clfrac = (qtmax-qstar)*(qtmax-qstar) /
1132 & ((qtmax-qtmin)*(qtmax-qtmode))
1133 elseif ( (qtmin <= qstar).and.(qstar < qtmode) )
then
1134 clfrac = 1. - ((qstar-qtmin)*(qstar-qtmin)
1135 & /( (qtmax-qtmin)*(qtmode-qtmin)))
1136 elseif ( qstar <= qtmin )
then
1139 elseif(flag == 4)
then
1140 if (qtmean > 1.0e-20)
then
1141 qtmedian = qtmean*exp(-0.5*sigmaqt1*sigmaqt1)
1142 aux = log(qtmedian/qstar)/sqrt(2.0)/sigmaqt1
1143 aux = min(max(aux, -20.0), 20.0)
1144 clfrac = 0.5*(1.0+
erf_app(aux))
1158 & qstar4,condensate4, clfrac4)
1163 real qtmean4, sigmaqt14, sigmaqt24, qstar4, condensate4, clfrac4
1165 real *8 :: qtmode, qtmin, qtmax, constA, constB, cloudf
1166 &, term1, term2, term3
1167 &, qtmean, sigmaqt1, sigmaqt2, qstar, condensate
1168 &, qtmedian, aux, clfrac, tx1
1170 qtmean = dble(qtmean4)
1171 sigmaqt1 = dble(sigmaqt14)
1172 sigmaqt2 = dble(sigmaqt24)
1173 qstar = dble(qstar4)
1174 clfrac = dble(clfrac4)
1177 if(qtmean+sigmaqt1 < qstar)
then
1179 elseif(qstar > qtmean-sigmaqt1)
then
1180 if(sigmaqt1 > 0.d0)
then
1181 tx1 = min(qtmean+sigmaqt1-qstar, 2.d0*sigmaqt1)
1182 condensate = tx1*tx1 / (4.d0*sigmaqt1)
1184 condensate = qtmean - qstar
1187 condensate = qtmean - qstar
1189 elseif(flag == 2)
then
1190 qtmode = qtmean + (sigmaqt1-sigmaqt2)/3.d0
1191 qtmin = min(qtmode-sigmaqt1,0.d0)
1192 qtmax = qtmode + sigmaqt2
1193 if ( qtmax < qstar )
then
1195 elseif ( qtmode <= qstar .and. qstar < qtmax )
then
1196 constb = 2.d0 / ( (qtmax - qtmin)*(qtmax-qtmode) )
1197 cloudf = (qtmax-qstar)*(qtmax-qstar) * 0.5d0 * constb
1198 term1 = (qstar*qstar*qstar)/3.d0
1199 term2 = (qtmax*qstar*qstar)/2.d0
1200 term3 = (qtmax*qtmax*qtmax)/6.d0
1201 condensate = constb * (term1-term2+term3) - qstar*cloudf
1202 elseif ( qtmin <= qstar .and. qstar < qtmode )
then
1203 consta = 2.d0 / ((qtmax-qtmin)*(qtmode-qtmin))
1204 cloudf = 1.d0 - (qstar-qtmin)*(qstar-qtmin)*0.5d0*consta
1205 term1 = qstar*qstar*qstar/3.d0
1206 term2 = qtmin*qstar*qstar/2.d0
1207 term3 = qtmin*qtmin*qtmin/6.d0
1208 condensate = qtmean - (consta*(term1-term2+term3))
1210 elseif ( qstar <= qtmin )
then
1211 condensate = qtmean - qstar
1214 elseif(flag == 4)
then
1216 if (qtmean > 1.0e-20)
then
1217 aux = 0.5*sigmaqt1*sigmaqt1
1218 qtmedian = qtmean*exp(-aux)
1219 aux = (aux + log(qtmedian/qstar))/(sqrt(2.0)*sigmaqt1)
1220 aux = min(max(aux, -20.0), 20.0)
1221 aux = 1.0 + dble(
erf_app(sngl(aux)))
1222 condensate = (0.5*qtmean*aux/qstar - clfrac) * qstar
1223 condensate= min(condensate, qtmean)
1229 condensate4 = real(condensate)
1237 subroutine cnvsrc( DT, ICEPARAM, SCLMFDFR, MASS, iMASS, PL, &
1238 & TE, QV, DCF, DMF, QLA, QIA, CF, AF, QS, &
1239 & NL, NI, CNVFICE, CNVNDROP, CNVNICE)
1241 real,
intent(in) :: DT, ICEPARAM, SCLMFDFR, MASS, iMASS, QS &
1243 real,
intent(inout) :: TE, AF,QV, QLA, QIA, NI, NL &
1244 &, CNVFICE, CNVNDROP, CNVNICE
1246 real :: TEND,QVx,QCA,fQi
1248 integer,
parameter :: STRATEGY = 3
1249 real,
parameter :: RL_cub = 1.0e-15, ri_cub = 6.4e-14
1254 if (te < t_ice_all)
then
1256 elseif (te > t_ice_max)
then
1269 if ( ( (1.0-fqi)*dcf > 0.0) .and. (cnvndrop <= 0.0))
then
1270 cnvndrop = (1.0-fqi)*dcf/( 1.333 * mapl_pi *rl_cub*997.0)
1273 if ((fqi*dcf > 0.0) .and. (cnvnice <= 0.0))
then
1274 cnvnice = fqi*dcf/( 1.333 * mapl_pi *ri_cub*500.0)
1277 nl = max(nl + cnvndrop*imass*dt, 0.0)
1278 ni = max(ni + cnvnice*imass*dt, 0.0)
1284 tend = dmf*imass * sclmfdfr
1285 af = min(af + tend*dt, 1.0)
1287 if ( af < 1.0 )
then
1288 qvx = ( qv - qs * af )/(1.-af)
1317 subroutine precip3( K,LM , DT , FRLAND , RHCR3 , QPl , QPi , &
1318 & QCl , QCi , TE , QV , mass , imass , PL , dZE , QDDF3 , AA , BB ,&
1319 & AREA , RAIN , SNOW , PFl_above , PFi_above , EVAP_DD_above, &
1320 & SUBL_DD_above, REVAP_DIAG , RSUBL_DIAG , ACRLL_DIAG , &
1321 & ACRIL_DIAG , PFL_DIAG , PFI_DIAG , VFALLRN , VFALLSN , FRZ_DIAG ,&
1322 & ENVFC,DDRFC, AF, CF, PCBL,i )
1325 integer,
intent(in) :: K,LM,i
1327 real,
intent(in ) :: DT
1329 real,
intent(inout) :: QV,QPl,QPi,QCl,QCi,TE
1331 real,
intent(in ) :: mass,imass
1332 real,
intent(in ) :: PL
1333 real,
intent(in ) :: AA,BB
1334 real,
intent(in ) :: RHCR3
1335 real,
intent(in ) :: dZE
1336 real,
intent(in ) :: QDDF3
1337 real,
intent( out) :: RAIN,SNOW
1338 real,
intent(in ) :: AREA
1339 real,
intent(in ) :: FRLAND
1341 real,
intent(inout) :: PFl_above, PFi_above
1342 real,
intent(inout) :: EVAP_DD_above, SUBL_DD_above
1344 real,
intent( out) :: REVAP_DIAG
1345 real,
intent( out) :: RSUBL_DIAG
1346 real,
intent( out) :: ACRLL_DIAG,ACRIL_DIAG
1347 real,
intent( out) :: PFL_DIAG, PFI_DIAG
1348 real,
intent(inout) :: FRZ_DIAG
1349 real,
intent( out) :: VFALLSN, VFALLRN
1351 real,
intent(in ) :: ENVFC,DDRFC
1353 real,
intent(in ) :: AF,CF, PCBL
1356 real :: PFi,PFl,QS,dQS,ENVFRAC
1357 real,
save :: TKo,QKo,QSTKo,DQSTKo,RH_BOX,T_ED,QPlKo,QPiKo
1358 real :: Ifactor,RAINRAT0,SNOWRAT0
1359 real :: FALLRN,FALLSN,VEsn,VErn,NRAIN,NSNOW,Efactor
1361 real :: TinLAYERrn,DIAMrn,DROPRAD
1362 real :: TinLAYERsn,DIAMsn,FLAKRAD,pl100
1364 real :: EVAP,SUBL,ACCR,MLTFRZ,EVAPx,SUBLx
1365 real :: EVAP_DD,SUBL_DD,DDFRACT
1368 real :: tmpARR, CFR, aux
1370 real,
parameter :: TRMV_L = 1.0
1372 real :: TAU_FRZ, TAU_MLT, QSICE, DQSI
1374 integer :: NS, NSMX, itr,L
1376 logical,
parameter :: taneff = .true.
1378 real,
parameter :: B_SUB = 1.00
1379 real esl,esi, esn,desdt,weight,tc,hlatsb,hlatvp,hltalt,tterm,
1386 aux = min(max((pl- pcbl)/10.0, -20.0), 20.0)
1387 aux = 1.0/(1.0+exp(-aux))
1388 envfrac = envfc + (1.0-envfc)*aux
1389 envfrac = min(envfrac,1.)
1395 if ( cfr < 0.99)
then
1396 tmparr = 1./(1.-cfr)
1402 IF ( area > 0. )
THEN
1403 ifactor = 1./ ( area )
1408 ifactor = max( ifactor, 1.)
1420 esn=min(fpvs(te),pl100)
1421 qs= min(epsqs*esn/(pl100-omeps*esn),1.)
1427 esi=min(
fpvsi(tko),pl100)
1428 qsice= min(epsqs*esi/(pl100-omeps*esi),1.)
1429 hltalt = hlatv + hlatf
1430 desdt = hltalt*esi/(rgasv*tko*tko)
1431 if (qsice < 1.0)
then
1432 gam = hltalt*qsice*pl100*desdt/(mapl_cp*esi
1433 & * (pl100 - omeps*esi))
1437 dqsi = (mapl_cp/hltalt)*gam
1452 qpl = qpl + pfl_above * imass
1454 qpi = qpi + pfi_above * imass
1458 accr = b_sub * c_acc * ( qpl*mass ) *qcl
1460 accr = min( accr , qcl )
1465 acrll_diag = accr / dt
1468 accr = b_sub * c_acc * ( qpi*mass ) *qcl
1470 accr = min( accr , qcl )
1475 te = te + alhfbcp*accr
1477 acril_diag = accr / dt
1479 rainrat0 = ifactor*qpl*mass/dt
1480 snowrat0 = ifactor*qpi*mass/dt
1482 call marshpalmq2(rainrat0,pl,diamrn,nrain,fallrn,vern)
1483 call marshpalmq2(snowrat0,pl,diamsn,nsnow,fallsn,vesn)
1485 IF ( frland < 0.1 )
THEN
1492 tinlayerrn = dze / ( max(fallrn,0.)+0.01 )
1493 tinlayersn = dze / ( max(fallsn,0.)+0.01 )
1498 IF ( (te > mapl_tice ) .and. (te <= mapl_tice+5. ) )
THEN
1499 mltfrz = tinlayersn * qpi *( te - mapl_tice ) / tau_frz
1500 mltfrz = min( qpi , mltfrz )
1501 te = te - alhfbcp*mltfrz
1505 frz_diag = frz_diag - mltfrz / dt
1508 IF ( te > mapl_tice+5. )
THEN
1510 te = te - alhfbcp*mltfrz
1514 frz_diag = frz_diag - mltfrz / dt
1517 if ( k >= lm-1 )
THEN
1518 IF ( te > mapl_tice+0. )
THEN
1520 te = te - alhfbcp*mltfrz
1525 frz_diag = frz_diag - mltfrz / dt
1529 IF ( te <= mapl_tice )
THEN
1530 te = te + alhfbcp*qpl
1535 frz_diag = frz_diag + mltfrz / dt
1547 esn=min(fpvs(tko),pl100)
1548 qstko= min(epsqs*esn/(pl100-omeps*esn),1.)
1549 tc = tko - mapl_tice
1550 lflg = (tc >= -ttrice .and. tc < 0.)
1551 weight = min(-tc*trinv,1.0)
1552 hlatsb = hlatv + weight*hlatf
1553 hlatvp = hlatv - 2369.0*tc
1554 if (tko < mapl_tice)
then
1560 tterm = pcf(1) + tc*(pcf(2) + tc*(pcf(3)
1561 & +tc*(pcf(4) + tc*pcf(5))))
1565 desdt = hltalt*esn/(rgasv*tko*tko) + tterm*trinv
1566 dqstko=(epsqs + omeps*qstko)/(pl100 - omeps*esn)*desdt
1571 esi=min(
fpvsi(tko),pl100)
1572 qsice= min(epsqs*esi/(pl100-omeps*esi),1.)
1573 hltalt = hlatv + hlatf
1574 desdt = hltalt*esi/(rgasv*tko*tko)
1575 if (qsice < 1.0)
then
1576 gam = hltalt*qsice*pl100*desdt/(mapl_cp*esi
1577 & * (pl100 - omeps*esi))
1581 dqsi = (mapl_cp/hltalt)*gam
1583 qstko = max( qstko , 1.0e-7 )
1584 qsice = max( qsice , 1.0e-7 )
1586 if (tmparr > 0.0)
then
1587 qko =(qko -qstko*cfr)*tmparr
1593 IF ( rh_box < rhcr3 )
THEN
1594 efactor = rho_w * ( aa + bb ) / (rhcr3 - rh_box )
1603 if ( ( rh_box < rhcr3 ) .AND. ( diamrn > 0.00 ) .AND.
1604 & ( pl > 100.) .AND. ( pl < revap_off_p ) )
then
1605 droprad = 0.5*diamrn
1606 t_ed = efactor * droprad**2
1607 t_ed = t_ed * ( 1.0 + dqstko*alhlbcp )
1609 evap = qpl*(1.0 - exp( -c_ev_r * vern * landseaf *envfrac*
1610 & tinlayerrn / t_ed ) )
1617 if (tmparr > 0.0)
then
1618 qko = (qko -qsice*cfr)*tmparr
1623 IF ( rh_box < rhcr3 )
THEN
1624 efactor = 0.5*rho_w * ( aa+bb ) / (rhcr3-rh_box )
1630 if ( ( rh_box < rhcr3 ) .AND. ( diamsn > 0.00 ) .AND.
1631 & ( pl > 100. ) .AND. ( pl < revap_off_p ) )
then
1632 flakrad = 0.5*diamsn
1633 t_ed = efactor * flakrad**2
1634 t_ed = t_ed * ( 1.0 + dqsi*alhsbcp)
1635 subl = qpi*(1.0 - exp( -c_ev_s * vesn * landseaf * envfrac
1636 & * tinlayersn / t_ed ) )
1646 evap = (evap+evapx) /2.0
1647 subl = (subl+sublx) /2.0
1650 evap = evap*(1.-cfr)
1651 subl = subl*(1.-cfr)
1653 subl = min(qpi, max(subl,0.))
1654 evap = min(qpl, max(evap,0.))
1657 qko =qko + evap + subl
1658 tko = tko - evap * alhlbcp - subl * alhsbcp
1665 evap_dd = evap_dd_above + ddfract*evap*mass
1666 evap = evap - ddfract*evap
1667 subl_dd = subl_dd_above + ddfract*subl*mass
1668 subl = subl - ddfract*subl
1671 qv = qv + evap + subl
1672 te = te - evap * alhlbcp - subl * alhsbcp
1674 revap_diag = evap / dt
1675 rsubl_diag = subl / dt
1686 evap = qddf3*evap_dd/mass
1687 subl = qddf3*subl_dd/mass
1689 subl = min(qpi, max(subl,0.))
1690 evap = min(qpl, max(evap,0.))
1691 qv = qv + evap + subl
1692 te = te - evap * alhlbcp - subl * alhsbcp
1693 revap_diag = revap_diag + evap / dt
1694 rsubl_diag = rsubl_diag + subl / dt
1707 evap_dd_above = evap_dd
1708 subl_dd_above = subl_dd
1717 real,
intent(in ) :: RAIN,PR
1718 real,
intent(out) :: DIAM3,NTOTAL,W,VE
1720 real :: RAIN_DAY,LAMBDA,A,B,SLOPR,DIAM1
1722 real,
parameter :: N0 = 0.08
1726 real :: RX(8) , D3X(8)
1730 rx = (/ 0. , 5. , 20. , 80. , 320. , 1280., 4*1280., 16*1280. /)
1731 d3x = (/ 0.019, 0.032, 0.043, 0.057, 0.076, 0.102, 0.137 ,
1734 rain_day = rain * 3600. *24.
1736 IF ( rain_day <= 0.00 )
THEN
1744 IF ( (rain_day <= rx(iqd+1)) .AND. (rain_day > rx(iqd)))
THEN
1745 slopr =( d3x(iqd+1)-d3x(iqd) ) / ( rx(iqd+1)-rx(iqd))
1746 diam3 = d3x(iqd) + (rain_day-rx(iqd))*slopr
1750 IF ( rain_day >= rx(8) )
THEN
1754 ntotal = 0.019*diam3
1756 diam3 = 0.664 * diam3
1758 w = (2483.8 * diam3 + 80.)*sqrt(1000./pr)
1760 ve = max( 0.99*w/100. , 1.000 )
1767 ntotal = ntotal*1.0e6
1775 real,
intent(in ) :: TEMP,Q_SAT
1776 real,
intent(in ) :: PR
1777 real,
intent(out) :: AA,BB
1781 real,
parameter :: EPSILON = 0.622
1782 real,
parameter :: K_COND = 2.4e-2
1783 real,
parameter :: DIFFU = 2.2e-5
1785 e_sat = 100.* pr * q_sat /( (epsilon) + (1.0-(epsilon))*q_sat )
1787 aa = (
get_alhx3(temp)**2 ) / ( k_cond*mapl_rvap*(temp**2) )
1790 bb = mapl_rvap*temp / ( diffu*(1000./pr)*e_sat )
1799 real,
intent(in) :: te,pl,nn,qcl
1805 rho = 100.*pl / (mapl_rgas*te )
1807 radius = muu/(nn*rho_w*(4./3.)*mapl_pi)
1808 radius = radius**(1./3.)
1817 real,
intent(in) :: temp
1821 if ( temp <= t_ice_all )
then
1823 else if ( (temp > t_ice_all) .AND. (temp <= t_ice_max) )
then
1824 icefrct = 1.00 - ( temp - t_ice_all ) / ( t_ice_max -
1827 icefrct = min(icefrct,1.00)
1828 icefrct = max(icefrct,0.00)
1830 icefrct = icefrct**icefrpwr
1838 real,
intent(in) :: t
1845 if ( t < t_ice_all )
then
1853 if ( (t <= t_x) .and. (t >= t_ice_all) )
then
1854 alhx3 = mapl_alhs + (mapl_alhl-mapl_alhs)
1855 & *( t - t_ice_all )/( t_x - t_ice_all )
1866 real,
intent(in) :: t
1867 real,
intent(in),
optional :: t_trans
1868 real,
intent(in),
optional :: t_freez
1872 if (
present( t_trans ))
then
1877 if (
present( t_freez ))
then
1888 if ( t <= t_x .and. t >= t_f )
then
1889 icefrac = 1.00 - ( t - t_f ) /( t_x - t_f )
1902 & QLLS , QLCN , CF , AF , NL , NI , DQALL , FQI )
1904 real ,
intent(in ) :: DTIME, PL, TE
1905 real ,
intent(inout ) :: DQALL
1906 real ,
intent(in) :: QV, QLLS, QLCN, QICN, QILS
1907 real ,
intent(in) :: CF, AF, NL, NI
1908 real,
intent (out) :: FQI
1909 real :: DC, TEFF,QCm,DEP, QC, QS, RHCR, DQSL, DQSI, QI, TC, DIFF,
1910 & denair, denice, aux, dcf, qtot, lhcorr, ql, dqi, dql, qvinc,
1911 & qsliq, cfall, new_qi, new_ql, qsice, fqi_0, qs_0, dqs_0, fqa,
1913 real esl,esi, esn,desdt,weight,hlatsb,hlatvp,hltalt,tterm,
1924 if (qtot .gt. 0.0) fqa = (qicn+qils)/qtot
1928 cfall= min(cf+af, 1.0)
1934 if (te .ge. t_ice_max)
then
1936 elseif(te .le. t_ice_all)
then
1942 if (qils .le. 0.0)
return
1947 esl=min(fpvsl(te),pl100)
1948 qsliq= min(epsqs*esl/(pl100-omeps*esl),1.)
1952 esi=min(fpvsi(te),pl100)
1953 qsice= min(epsqs*esi/(pl100-omeps*esi),1.)
1954 hltalt = hlatv + hlatf
1955 desdt = hltalt*esi/(rgasv*te*te)
1956 if (qsice < 1.0)
then
1957 gam = hltalt*qsice*pl100*desdt/(mapl_cp*esi
1958 & * (pl100 - omeps*esi))
1962 dqsi = (mapl_cp/hltalt)*gam
1967 qvinc =min(qvinc, qsliq)
1969 diff=(0.211*1013.25/(pl+0.1))*(((te+0.1)/273.0)**1.94)*1e-4
1970 denair = pl*100.0/mapl_rgas/te
1971 denice = 1000.0*(0.9167 - 1.75e-4*tc -5.0e-7*tc*tc)
1972 lhcorr = ( 1.0 + dqsi*alhsbcp)
1974 if ((nix .gt. 1.0) .and. (qils .gt. 1.0e-10))
then
1975 dc = max((qils/(nix*denice*mapl_pi))**(0.333), 20.0e-6)
1981 teff= nix*denair*2.0*mapl_pi*diff*dc/lhcorr
1984 if ((teff .gt. 0.0) .and. (qils .gt. 1.0e-14))
then
1985 aux =max(min(dtime*teff, 20.0), 0.0)
1986 dep=(qvinc-qsice)*(1.0-exp(-aux))/dtime
1989 dep=max(dep, -qils/dtime)
1995 if (dqall .ge. 0.0)
then
1997 if (dep .gt. 0.0)
then
1998 dqi = min(dep, dqall + qlls/dtime)
2007 if (dqall .lt. 0.0)
then
2008 dql = max(dqall, -qlls/dtime)
2009 dqi = max(dqall - dql, -qils/dtime)
2012 if (dqall .ne. 0.0) fqi=max(min(dqi/dqall, 1.0), 0.0)
2023 subroutine pfreezing ( ALPHA , PL , TE , QV , QCl , QAl , QCi , &
2024 & QAi , SC_ICE , CF , AF , PF, pdfflag)
2026 integer,
intent(in) :: pdfflag
2027 real ,
intent(in) :: PL, ALPHA, QV, SC_ICE, AF, TE, &
2028 & QCl, QCi, QAl, QAi, CF
2029 real ,
intent(out) :: PF
2031 real :: qt, QCx, QSn, tmpARR, CFALL, QVx, CFio, QA, QAx, QC, QI, &
2032 & QL, DQSx, sigmaqt1, sigmaqt2, qsnx, esl, esi,pl100
2040 if ( cfall >= 1.0 )
then
2048 esi = min(fpvsi(te),pl100)
2049 qsn = max(min(epsqs*esi/(pl100-omeps*esi), 1.0), 1.0e-9)
2052 if ( cfall < 0.99 )
then
2053 tmparr = 1./(1.0-cfall)
2057 qvx = ( qv - qsn*cfall )*tmparr
2066 if(pdfflag < 2)
then
2067 sigmaqt1 = max(alpha, 0.1) * qsn
2068 sigmaqt2 = max(alpha, 0.1) * qsn
2069 elseif(pdfflag == 2)
then
2074 sigmaqt1 = alpha * qsn
2075 sigmaqt2 = alpha * qsn
2076 elseif(pdfflag == 4)
then
2077 sigmaqt1 = max(alpha/sqrt(3.0), 0.001)
2079 write(0,*)
' Aborting : invalid pdfflag=',pdfflag
2083 call pdffrac(pdfflag,qt,sigmaqt1,sigmaqt2,qsn,cfio)
2085 pf = min(max(cfio*(1.0-cfall), 0.0), 0.999)
2099 integer,
intent(in) :: im, lm
2100 real ,
intent(inout),
dimension(:,:) :: te,qcl,qci, qal, qai, &
2103 real ,
dimension(im,lm) :: fqi,dqil, dqmax, qltot, qitot, fqal,
2110 where (qitot+qltot > 0.0)
2111 fqa= (qai+qal)/(qitot+qltot)
2118 where( te <= t_ice_all )
2119 dqmax = (t_ice_all - te)/(alhsbcp-alhlbcp)
2120 dqil = min(qltot , dqmax)
2123 where ((dqil <= dqmax) .and. (dqil > 0.0))
2127 where ((dqil > dqmax) .and. (dqil > 0.0))
2128 dnil = nl*dqmax/dqil
2131 dqil = max( 0., dqil )
2135 dqil = min(qltot,dqil)
2136 qitot = qitot + dqil
2137 qltot = qltot - dqil
2140 te = te + (alhsbcp-alhlbcp)*dqil
2147 where( te > t_ice_max )
2148 dqmax = (te-t_ice_max) / (alhsbcp-alhlbcp)
2149 dqil = min(qitot, dqmax)
2152 where ((dqil .le. dqmax) .and. (dqil .gt. 0.0))
2155 where ((dqil .gt. dqmax) .and. (dqil .gt. 0.0))
2156 dnil = ni*dqmax/dqil
2158 dqil = max( 0., dqil )
2162 dqil = min(qitot,dqil)
2163 qitot = qitot - dqil
2164 qltot = qltot + dqil
2168 te = te - (alhsbcp-alhlbcp)*dqil
2170 qci = qitot*(1.0-fqa)
2172 qcl = qltot*(1.0-fqa)
2182 & SMAXL, SMAXI, WSUB, CCN01, CCN04, CCN1, NHET_NUC, NLIM_NUC, SO4, &
2183 & ORG, BCARBON, DUST, SEASALT, NCPL_VOL, NCPI_VOL, NRAIN, NSNOW, &
2184 & CDNC_NUC, INC_NUC, SAT_RAT, QSTOT, QRTOT, CLDREFFS, CLDREFFR, &
2185 & DQVDT_micro,DQIDT_micro, DQLDT_micro, DTDT_micro, RL_MASK, &
2186 & RI_MASK, KAPPA, SC_ICE, CFICE, CFLIQ, RHICE, RHLIQ, RAD_CF, &
2187 & RAD_QL, RAD_QI, RAD_QS, RAD_QR, RAD_QV, CLDREFFI, CLDREFFL, &
2188 & NHET_IMM, NHET_DEP, NHET_DHF, DUST_IMM, DUST_DEP, DUST_DHF, SCF, &
2189 & SCF_ALL, SIGW_GW, SIGW_CNV, SIGW_TURB, SIGW_RC, RHCmicro, &
2190 & DNHET_IMM, NONDUST_IMM, NONDUST_DEP, BERG, BERGSO, MELT, &
2191 & DNHET_CT, DTDT_macro, QCRES, DT_RASP, FRZPP_LS, SNOWMELT_LS, &
2192 & QIRES, AUTICE, PFRZ, DNCNUC, DNCHMSPLIT, DNCSUBL, DNCAUTICE, &
2193 & DNCACRIS, DNDCCN, DNDACRLS, DNDEVAPC, DNDACRLR, DNDAUTLIQ)
2198 real ,
pointer ,
dimension(:,:,:) :: smaxl,smaxi, wsub, ccn01, &
2199 & CCN04, CCN1, NHET_NUC, NLIM_NUC, SO4, ORG, BCARBON, DUST, &
2200 & SEASALT, NCPL_VOL, NCPI_VOL, NRAIN, NSNOW, CDNC_NUC, INC_NUC, &
2201 & SAT_RAT, QSTOT, QRTOT, CLDREFFS, CLDREFFR, DQVDT_micro, &
2202 &DQIDT_micro, DQLDT_micro, DTDT_micro, RL_MASK, RI_MASK, KAPPA, &
2203 & SC_ICE, CFICE, CFLIQ, RHICE, RHLIQ, ALPH, RAD_CF, RAD_QL, RAD_QI,&
2204 & RAD_QS, RAD_QR, RAD_QV, CLDREFFI, CLDREFFL, NHET_IMM, NHET_DEP, &
2205 & NHET_DHF, DUST_IMM, DUST_DEP, DUST_DHF, SCF, SCF_ALL, SIGW_GW, &
2206 & SIGW_CNV, SIGW_TURB, SIGW_RC, RHCmicro, DNHET_IMM, NONDUST_IMM, &
2207 & NONDUST_DEP, BERG, BERGSO, MELT, DNHET_CT, DTDT_macro, QCRES, &
2208 & DT_RASP, FRZPP_LS, SNOWMELT_LS, QIRES, AUTICE, PFRZ, DNCNUC, &
2209 & DNCHMSPLIT, DNCSUBL, DNCAUTICE, DNCACRIS, DNDCCN, DNDACRLS, &
2210 & DNDEVAPC, DNDACRLR, DNDAUTLIQ
2216 IF(
ASSOCIATED(smaxl) ) smaxl = 0.
2217 IF(
ASSOCIATED(smaxi) ) smaxi = 0.
2218 IF(
ASSOCIATED(wsub) ) wsub = 0.
2219 IF(
ASSOCIATED(ccn01) ) ccn01 = 0.
2220 IF(
ASSOCIATED(ccn04) ) ccn04 = 0.
2221 IF(
ASSOCIATED(ccn1) ) ccn1 = 0.
2222 IF(
ASSOCIATED(nhet_nuc) ) nhet_nuc = 0.
2223 IF(
ASSOCIATED(nlim_nuc) ) nlim_nuc = 0.
2224 IF(
ASSOCIATED(so4) ) so4 = 0.
2225 IF(
ASSOCIATED(org) ) org = 0.
2226 IF(
ASSOCIATED(bcarbon) ) bcarbon = 0.
2227 IF(
ASSOCIATED(dust) ) dust = 0.
2228 IF(
ASSOCIATED(seasalt) ) seasalt = 0.
2229 IF(
ASSOCIATED(ncpl_vol) ) ncpl_vol = 0.
2230 IF(
ASSOCIATED(ncpi_vol) ) ncpi_vol = 0.
2232 IF(
ASSOCIATED(nrain) ) nrain = 0.
2233 IF(
ASSOCIATED(nsnow) ) nsnow = 0.
2234 IF(
ASSOCIATED(cdnc_nuc) ) cdnc_nuc = 0.
2235 IF(
ASSOCIATED(inc_nuc) ) inc_nuc = 0.
2236 IF(
ASSOCIATED(sat_rat) ) sat_rat = 0.
2237 IF(
ASSOCIATED(qstot) ) qstot = 0.
2238 IF(
ASSOCIATED(qrtot) ) qrtot = 0.
2240 IF(
ASSOCIATED(dqvdt_micro) ) dqvdt_micro = 0.
2241 IF(
ASSOCIATED(dqidt_micro) ) dqidt_micro = 0.
2242 IF(
ASSOCIATED(dqldt_micro) ) dqldt_micro = 0.
2243 IF(
ASSOCIATED(dtdt_micro) ) dtdt_micro = 0.
2244 IF(
ASSOCIATED(dtdt_macro) ) dtdt_macro = 0.
2246 IF(
ASSOCIATED(rl_mask) ) rl_mask = 0.
2247 IF(
ASSOCIATED(ri_mask) ) ri_mask = 0.
2248 IF(
ASSOCIATED(kappa) ) kappa = 0.
2249 IF(
ASSOCIATED(sc_ice)) sc_ice = 0.
2250 IF(
ASSOCIATED(rhice) ) rhice = 0.
2251 IF(
ASSOCIATED(rhliq) ) rhliq = 0.
2252 IF(
ASSOCIATED(cfice) ) cfice = 0.
2253 IF(
ASSOCIATED(cfliq) ) cfliq = 0.
2254 IF(
ASSOCIATED(alph) ) alph = 0.
2257 IF(
ASSOCIATED(rad_cf) ) rad_cf = 0.
2258 IF(
ASSOCIATED(rad_ql) ) rad_ql = 0.
2259 IF(
ASSOCIATED(rad_qi) ) rad_qi = 0.
2260 IF(
ASSOCIATED(rad_qs) ) rad_qs = 0.
2261 IF(
ASSOCIATED(rad_qr) ) rad_qr = 0.
2262 IF(
ASSOCIATED(rad_qv) ) rad_qv = 0.
2263 IF(
ASSOCIATED(cldreffi) ) cldreffi = 0.
2264 IF(
ASSOCIATED(cldreffl) ) cldreffl = 0.
2265 IF(
ASSOCIATED(cldreffs) ) cldreffs = 0.
2266 IF(
ASSOCIATED(cldreffr) ) cldreffr = 0.
2268 IF(
ASSOCIATED(nhet_imm) ) nhet_imm = 0.
2269 IF(
ASSOCIATED(nhet_dep) ) nhet_dep = 0.
2270 IF(
ASSOCIATED(nhet_dhf) ) nhet_dhf = 0.
2271 IF(
ASSOCIATED(dust_imm) ) dust_imm = 0.
2272 IF(
ASSOCIATED(dust_dep) ) dust_dep = 0.
2273 IF(
ASSOCIATED(dust_dhf) ) dust_dhf = 0.
2274 IF(
ASSOCIATED(nondust_imm) ) nondust_imm = 0.
2275 IF(
ASSOCIATED(nondust_dep) ) nondust_dep = 0.
2278 IF(
ASSOCIATED(scf) ) scf = 0.
2279 IF(
ASSOCIATED(scf_all) ) scf_all = 0.
2280 IF(
ASSOCIATED(sigw_gw) ) sigw_gw = 0.
2281 IF(
ASSOCIATED(sigw_cnv) ) sigw_cnv = 0.
2282 IF(
ASSOCIATED(sigw_turb) ) sigw_turb = 0.
2283 IF(
ASSOCIATED(sigw_rc) ) sigw_rc = 0.
2284 IF(
ASSOCIATED(rhcmicro) ) rhcmicro = 0.
2285 IF(
ASSOCIATED(dnhet_imm) ) dnhet_imm = 0.
2286 IF(
ASSOCIATED(berg) ) berg = 0.
2287 IF(
ASSOCIATED(bergso)) bergso = 0.
2288 IF(
ASSOCIATED(melt) ) melt = 0.
2289 IF(
ASSOCIATED(dnhet_ct) ) dnhet_ct = 0.
2290 IF(
ASSOCIATED(dt_rasp) ) dt_rasp = 0.
2292 IF(
ASSOCIATED(qcres) ) qcres = 0.
2293 IF(
ASSOCIATED(qires) ) qires = 0.
2294 IF(
ASSOCIATED(autice) ) autice = 0.
2295 IF(
ASSOCIATED(frzpp_ls) ) frzpp_ls = 0.
2296 IF(
ASSOCIATED(snowmelt_ls) ) snowmelt_ls = 0.
2297 IF(
ASSOCIATED(pfrz) ) pfrz = 0.
2299 IF(
ASSOCIATED(dncnuc) ) dncnuc = 0.
2300 IF(
ASSOCIATED(dncsubl) ) dncsubl = 0.
2301 IF(
ASSOCIATED(dnchmsplit) ) dnchmsplit = 0.
2302 IF(
ASSOCIATED(dncautice) ) dncautice = 0.
2303 IF(
ASSOCIATED(dncacris) ) dncacris = 0.
2304 IF(
ASSOCIATED(dndccn) ) dndccn = 0.
2305 IF(
ASSOCIATED(dndacrls) ) dndacrls = 0.
2306 IF(
ASSOCIATED(dndacrlr) ) dndacrlr = 0.
2307 IF(
ASSOCIATED(dndevapc) ) dndevapc = 0.
2308 IF(
ASSOCIATED(dndautliq) ) dndautliq = 0.
2321 real*8:: aa(4), axx, y
2322 DATA aa /0.278393d0,0.230389d0,0.000972d0,0.078108d0/
2325 axx = 1.d0 + y*(aa(1)+y*(aa(2)+y*(aa(3)+y*aa(4))))
2328 axx = 1.d0 - (1.d0/axx)
real function erf_app(x)
REAL FUNCTION erf (overwrites previous versions) THIS SUBROUTINE CALCULATES THE ERROR FUNCTION USING ...
subroutine pdffrac(flag, qtmean, sigmaqt1, sigmaqt2, qstar, clfrac)
This subroutine calculates.
subroutine, public fix_up_clouds_2m(qv, te, qlc, qic, cf, qla, qia, af, nl, ni, qc_min)
subroutine pdfcondensate(flag, qtmean4, sigmaqt14, sigmaqt24, qstar4, condensate4, clfrac4)
This subroutine calculates.
real function icefrac(t, t_trans, t_freez)
subroutine micro_aa_bb_3(temp, pr, q_sat, aa, bb)
This subroutine calculates.
subroutine, public update_cld(irun, lm, dt, alpha, qc_min, pdfshape, pl, qv, qcl, qal, qci, qai, te, cf, af, scice, ni, nl)
This subroutine calculates the cloud fraction that would correspond to the current condensate.
subroutine, public cloud_ptr_stubs(smaxl, smaxi, wsub, ccn01, ccn04, ccn1, nhet_nuc, nlim_nuc, so4, org, bcarbon, dust, seasalt, ncpl_vol, ncpi_vol, nrain, nsnow, cdnc_nuc, inc_nuc, sat_rat, qstot, qrtot, cldreffs, cldreffr, dqvdt_micro, dqidt_micro, dqldt_micro, dtdt_micro, rl_mask, ri_mask, kappa, sc_ice, cfice, cfliq, rhice, rhliq, rad_cf, rad_ql, rad_qi, rad_qs, rad_qr, rad_qv, cldreffi, cldreffl, nhet_imm, nhet_dep, nhet_dhf, dust_imm, dust_dep, dust_dhf, scf, scf_all, sigw_gw, sigw_cnv, sigw_turb, sigw_rc, rhcmicro, dnhet_imm, nondust_imm, nondust_dep, berg, bergso, melt, dnhet_ct, dtdt_macro, qcres, dt_rasp, frzpp_ls, snowmelt_ls, qires, autice, pfrz, dncnuc, dnchmsplit, dncsubl, dncautice, dncacris, dndccn, dndacrls, dndevapc, dndacrlr, dndautliq)
real function get_alhx3(t)
subroutine, public macro_cloud(irun, lm, dt, alf_fac, pp_dev, ppe_dev, qlwdtr_dev, th_dev, q_dev, qlw_ls_dev, qlw_an_dev, qiw_ls_dev, qiw_an_dev, anvfrc_dev, cldfrc_dev, physparams, sclmfdfr, alpht_dev, cnv_fice_dev, cnv_ndrop_dev, cnv_nice_dev, scice_dev, ncpl_dev, ncpi_dev, pfrz_dev, lprnt, ipr, rhc, pdfflag, qc_min)
This subroutine is the cloud macrophysics scheme in MG micriphysics.
subroutine marshpalmq2(rain, pr, diam3, ntotal, w, ve)
This subroutine calculates.
subroutine pdf_spread(k, lm, pp, alpha, alpht_diag, frland, rhc)
subroutine hystpdf(dt, alpha, pdfshape, qc_min, pl, qv, qcl, qal, qci, qai, te, cf, af, scice, ni, nl)
This subroutine calculates.
subroutine bergeron_iter(dtime, pl, te, qv, qils, qicn,
This subroutine calculates.
subroutine cnvsrc(dt, iceparam, sclmfdfr, mass, imass, pl, te, qv, dcf, dmf, qla, qia, cf, af, qs, nl, ni, cnvfice, cnvndrop, cnvnice)
This subroutine calculates.
subroutine, public meltfrz_inst(im, lm, te, qcl, qal, qci, qai, nl, ni)
This subroutine calculates instantaneous freezing of condensate.
subroutine pfreezing(alpha, pl, te, qv, qcl, qal, qci, qai, sc_ice, cf, af, pf, pdfflag)
This calculates the probability of finding a supersaturated parcel in the grid cell....
real function ice_fraction(temp)
real function ldradius3(pl, te, qcl, nn)
subroutine precip3(k, lm, dt, frland, rhcr3, qpl, qpi, qcl, qci, te, qv, mass, imass, pl, dze, qddf3, aa, bb, area, rain, snow, pfl_above, pfi_above, evap_dd_above, subl_dd_above, revap_diag, rsubl_diag, acrll_diag, acril_diag, pfl_diag, pfi_diag, vfallrn, vfallsn, frz_diag, envfc, ddrfc, af, cf, pcbl, i)
This subroutine calculates.
This module provides an Application Program Interface (API) for computing basic thermodynamic physics...
This module contains some of the most frequently used math and physics constants for GCM models.