93 & grav,cp,hvap,rv,fv,t0c,rgas,cvap,cliq,eps,epsm1, &
94 & im,km,jcap,delt,delp,prslp,psp,phil,qlc,qli, &
95 & q1,t1,u1,v1,cldwrk,rn,kbot,ktop,kcnv,islimsk, &
96 & dot,ncloud,ud_mf,dd_mf,dt_mf,cnvw,cnvc, &
97 & qlcn,qicn,w_upi,cf_upi,cnv_mfd, &
98 & cnv_dqldt,clcn,cnv_fice,cnv_ndrop,cnv_nice,mp_phys, &
99 & mp_phys_mg,clam,c0,c1,betal,betas,evfact,evfactl,pgcon, &
112 real(kind=kind_phys),
intent(in) :: grav, cp, hvap, rv, fv, t0c, &
113 & rgas, cvap, cliq, eps, epsm1
114 integer,
intent(in) :: im, km, jcap, ncloud, &
115 & mp_phys, mp_phys_mg
116 integer,
intent(inout) :: kbot(:), ktop(:), kcnv(:)
117 integer,
intent(in) :: islimsk(:)
118 real(kind=kind_phys),
intent(in) :: delt, clam, c0, c1, pgcon
119 real(kind=kind_phys),
intent(in) :: betal, betas, evfact, evfactl
120 real(kind=kind_phys),
intent(in) :: psp(:), delp(:,:), &
121 & prslp(:,:), dot(:,:), &
123 real(kind=kind_phys),
intent(inout) :: &
124 & qlc(:,:), qli(:,:), &
125 & q1(:,:), t1(:,:), &
126 & u1(:,:), v1(:,:), &
127 & cnvw(:,:), cnvc(:,:)
128 real(kind=kind_phys),
intent(out) :: cldwrk(:), rn(:), &
129 & dt_mf(:,:), dd_mf(:,:)
130 real(kind=kind_phys),
intent(out),
optional :: ud_mf(:,:)
131 real(kind=kind_phys),
intent(inout),
optional :: &
132 & qlcn(:,:), qicn(:,:), &
133 & w_upi(:,:), cnv_mfd(:,:), &
134 & cnv_dqldt(:,:), clcn(:,:), &
135 & cnv_fice(:,:), cnv_ndrop(:,:),&
136 & cnv_nice(:,:), cf_upi(:,:)
137 character(len=*),
intent(out) :: errmsg
138 integer,
intent(out) :: errflg
142 integer i, indx, jmn, k, kk, km1
145 real(kind=kind_phys) cxlamu, xlamde, xlamdd
148 real(kind=kind_phys) adw, aup, aafac,
153 & dq, dqsdp, dqsdt, dt,
154 & dt2, dtmax, dtmin, dv1h,
155 & dv1q, dv2h, dv2q, dv1u,
156 & dv1v, dv2u, dv2v, dv3q,
158 & dz, dz1, e1, edtmax,
159 & edtmaxl, edtmaxs, el2orc, elocp,
160 & es, etah, cthk, dthk,
162 & fact2, factor, fjcap, fkm,
165 & rain, rfact, shear, tem1,
167 & val2, w1, w1l, w1s,
170 & w4s, xdby, xpw, xpwd,
174 integer kb(im), kbcon(im), kbcon1(im),
175 & ktcon(im), ktcon1(im),
176 & jmin(im), lmin(im), kbmax(im),
179 real(kind=kind_phys) ps(im), del(im,km), prsl(im,km)
181 real(kind=kind_phys) aa1(im), acrt(im), acrtfct(im),
182 & delhbar(im), delq(im), delq2(im),
183 & delqbar(im), delqev(im), deltbar(im),
184 & deltv(im), dtconv(im), edt(im),
185 & edto(im), edtx(im), fld(im),
186 & hcdo(im,km), hmax(im), hmin(im),
187 & ucdo(im,km), vcdo(im,km),aa2(im),
188 & pbcdif(im), pdot(im), po(im,km),
189 & pwavo(im), pwevo(im), xlamud(im),
190 & qcdo(im,km), qcond(im), qevap(im),
191 & rntot(im), vshear(im), xaa0(im),
193 & xmb(im), xmbmax(im), xpwav(im),
194 & xpwev(im), delubar(im),delvbar(im)
196 real(kind=kind_phys) cincr, cincrmax, cincrmin
205 parameter(cthk=150.,cincrmax=180.,cincrmin=120.,dthk=25.)
207 real(kind=kind_phys) pfld(im,km), to(im,km), qo(im,km),
208 & uo(im,km), vo(im,km), qeso(im,km)
211 real(kind=kind_phys) qlko_ktcon(im), dellal(im,km), tvo(im,km),
212 & dbyo(im,km), zo(im,km), xlamue(im,km),
213 & fent1(im,km), fent2(im,km), frh(im,km),
214 & heo(im,km), heso(im,km),
215 & qrcd(im,km), dellah(im,km), dellaq(im,km),
216 & dellau(im,km), dellav(im,km), hcko(im,km),
217 & ucko(im,km), vcko(im,km), qcko(im,km),
218 & eta(im,km), etad(im,km), zi(im,km),
219 & qrcko(im,km), qrcdo(im,km),
220 & pwo(im,km), pwdo(im,km),
221 & tx1(im), sumx(im), cnvwt(im,km)
224 logical totflg, cnvflg(im), flg(im)
226 real(kind=kind_phys) pcrit(15), acritt(15), acrit(15)
228 data pcrit/850.,800.,750.,700.,650.,600.,550.,500.,450.,400.,
229 & 350.,300.,250.,200.,150./
230 data acritt/.0633,.0445,.0553,.0664,.075,.1082,.1521,.2216,
231 & .3151,.3677,.41,.5255,.7663,1.1686,1.6851/
235 real(kind=kind_phys) tf, tcr, tcrf
236 parameter(tf=233.16, tcr=263.16, tcrf=1.0/(tcr-tf))
243 el2orc = hvap*hvap/(rv*cp)
245 fact1 = (cvap-cliq)/rv
246 fact2 = hvap/rv-fact1*t0c
309 if(mp_phys == mp_phys_mg)
then
327 acrit(k) = acritt(k) * (975. - pcrit(k))
331 dtmin = max(dt2, val )
333 dtmax = max(dt2, val )
354 fjcap = (float(jcap) / 126.) ** 2
356 fjcap = max(fjcap,val)
357 fkm = (float(km) / 28.) ** 2
381 if (prsl(i,k)*tx1(i) .gt. 0.04) kmax(i) = k + 1
382 if (prsl(i,k)*tx1(i) .gt. 0.45) kbmax(i) = k + 1
383 if (prsl(i,k)*tx1(i) .gt. 0.70) kbm(i) = k + 1
387 kmax(i) = min(km,kmax(i))
388 kbmax(i) = min(kbmax(i),kmax(i))
389 kbm(i) = min(kbm(i),kmax(i))
398 zo(i,k) = phil(i,k) / g
404 zi(i,k) = 0.5*(zo(i,k)+zo(i,k+1))
405 xlamue(i,k) = clam / zi(i,k)
415 if (k .le. kmax(i))
then
416 pfld(i,k) = prsl(i,k) * 10.0
458 if (k .le. kmax(i))
then
459 qeso(i,k) = 0.01 * fpvs(to(i,k))
460 qeso(i,k) = eps * qeso(i,k) / (pfld(i,k) + epsm1*qeso(i,k))
462 qeso(i,k) = max(qeso(i,k), val1)
464 qo(i,k) = max(qo(i,k), val2 )
476 if (k .le. kmax(i))
then
478 tem = phil(i,k) + cp * to(i,k)
479 heo(i,k) = tem + hvap * qo(i,k)
480 heso(i,k) = tem + hvap * qeso(i,k)
498 if (k .le. kbm(i))
then
499 if(heo(i,k).gt.hmax(i))
then
510 if (k .le. kmax(i)-1)
then
511 dz = .5 * (zo(i,k+1) - zo(i,k))
512 dp = .5 * (pfld(i,k+1) - pfld(i,k))
513 es = 0.01 * fpvs(to(i,k+1))
514 pprime = pfld(i,k+1) + epsm1 * es
515 qs = eps * es / pprime
516 dqsdp = - qs / pprime
517 desdt = es * (fact1 / to(i,k+1) + fact2 / (to(i,k+1)**2))
518 dqsdt = qs * pfld(i,k+1) * desdt / (es * pprime)
519 gamma = el2orc * qeso(i,k+1) / (to(i,k+1)**2)
520 dt = (g * dz + hvap * dqsdp * dp) / (cp * (1. + gamma))
521 dq = dqsdt * dt + dqsdp * dp
522 to(i,k) = to(i,k+1) + dt
523 qo(i,k) = qo(i,k+1) + dq
524 po(i,k) = .5 * (pfld(i,k) + pfld(i,k+1))
532 if (k .le. kmax(i)-1)
then
533 qeso(i,k) = 0.01 * fpvs(to(i,k))
534 qeso(i,k) = eps * qeso(i,k) / (po(i,k) + epsm1*qeso(i,k))
536 qeso(i,k) = max(qeso(i,k), val1)
538 qo(i,k) = max(qo(i,k), val2 )
540 frh(i,k) = 1. - min(qo(i,k)/qeso(i,k), 1.)
541 heo(i,k) = .5 * g * (zo(i,k) + zo(i,k+1)) +
542 & cp * to(i,k) + hvap * qo(i,k)
543 heso(i,k) = .5 * g * (zo(i,k) + zo(i,k+1)) +
544 & cp * to(i,k) + hvap * qeso(i,k)
545 uo(i,k) = .5 * (uo(i,k) + uo(i,k+1))
546 vo(i,k) = .5 * (vo(i,k) + vo(i,k+1))
560 if (flg(i).and.k.le.kbmax(i))
then
561 if(k.gt.kb(i).and.heo(i,kb(i)).gt.heso(i,k))
then
571 if(kbcon(i).eq.kmax(i)) cnvflg(i) = .false.
576 totflg = totflg .and. (.not. cnvflg(i))
588 pdot(i) = 0.01 * dot(i,kbcon(i))
593 if(islimsk(i) == 1)
then
604 if(pdot(i).le.w4)
then
605 tem = (pdot(i) - w4) / (w3 - w4)
606 elseif(pdot(i).ge.-w4)
then
607 tem = - (pdot(i) + w4) / (w4 - w3)
616 tem1= .5*(cincrmax-cincrmin)
617 cincr = cincrmax - tem * tem1
618 pbcdif(i) = pfld(i,kb(i)) - pfld(i,kbcon(i))
619 if(pbcdif(i).gt.cincr)
then
627 totflg = totflg .and. (.not. cnvflg(i))
643 & (k.gt.kbcon(i).and.k.lt.kmax(i)))
then
644 xlamue(i,k) = xlamue(i,kbcon(i))
655 xlamud(i) = xlamue(i,kbcon(i))
665 & (k.gt.kbcon(i).and.k.lt.kmax(i)))
then
666 tem = qeso(i,k)/qeso(i,kbcon(i))
680 & (k.ge.kbcon(i).and.k.lt.kmax(i)))
then
681 tem = cxlamu * frh(i,k) * fent2(i,k)
682 xlamue(i,k) = xlamue(i,k)*fent1(i,k) + tem
697 if(k.lt.kbcon(i).and.k.ge.kb(i))
then
698 dz = zi(i,k+1) - zi(i,k)
699 ptem = 0.5*(xlamue(i,k)+xlamue(i,k+1))-xlamud(i)
700 eta(i,k) = eta(i,k+1) / (1. + ptem * dz)
711 if(k.gt.kbcon(i).and.k.lt.kmax(i))
then
712 dz = zi(i,k) - zi(i,k-1)
713 ptem = 0.5*(xlamue(i,k)+xlamue(i,k-1))-xlamud(i)
714 eta(i,k) = eta(i,k-1) * (1 + ptem * dz)
726 hcko(i,indx) = heo(i,indx)
727 ucko(i,indx) = uo(i,indx)
728 vcko(i,indx) = vo(i,indx)
739 if(k.gt.kb(i).and.k.lt.kmax(i))
then
740 dz = zi(i,k) - zi(i,k-1)
741 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1)) * dz
742 tem1 = 0.5 * xlamud(i) * dz
743 factor = 1. + tem - tem1
744 ptem = 0.5 * tem + pgcon
745 ptem1= 0.5 * tem - pgcon
746 hcko(i,k) = ((1.-tem1)*hcko(i,k-1)+tem*0.5*
747 & (heo(i,k)+heo(i,k-1)))/factor
748 ucko(i,k) = ((1.-tem1)*ucko(i,k-1)+ptem*uo(i,k)
749 & +ptem1*uo(i,k-1))/factor
750 vcko(i,k) = ((1.-tem1)*vcko(i,k-1)+ptem*vo(i,k)
751 & +ptem1*vo(i,k-1))/factor
752 dbyo(i,k) = hcko(i,k) - heso(i,k)
768 if (flg(i).and.k.lt.kmax(i))
then
769 if(k.ge.kbcon(i).and.dbyo(i,k).gt.0.)
then
778 if(kbcon1(i).eq.kmax(i)) cnvflg(i) = .false.
783 tem = pfld(i,kbcon(i)) - pfld(i,kbcon1(i))
792 totflg = totflg .and. (.not. cnvflg(i))
806 if (flg(i).and.k .lt. kmax(i))
then
807 if(k.gt.kbcon1(i).and.dbyo(i,k).lt.0.)
then
817 tem = pfld(i,kbcon(i))-pfld(i,ktcon(i))
818 if(tem.lt.cthk) cnvflg(i) = .false.
824 totflg = totflg .and. (.not. cnvflg(i))
834 hmin(i) = heo(i,kbcon1(i))
841 if (cnvflg(i) .and. k .le. kbmax(i))
then
842 if(k.gt.kbcon1(i).and.heo(i,k).lt.hmin(i))
then
854 jmin(i) = min(lmin(i),ktcon(i)-1)
855 jmin(i) = max(jmin(i),kbcon1(i)+1)
856 if(jmin(i).ge.ktcon(i)) cnvflg(i) = .false.
868 dp = 1000. * del(i,k)
869 xmbmax(i) = dp / (g * dt2)
882 qcko(i,kb(i)) = qo(i,kb(i))
883 qrcko(i,kb(i)) = qo(i,kb(i))
891 if(k.gt.kb(i).and.k.lt.ktcon(i))
then
892 dz = zi(i,k) - zi(i,k-1)
893 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
895 & + gamma * dbyo(i,k) / (hvap * (1. + gamma))
897 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1)) * dz
898 tem1 = 0.5 * xlamud(i) * dz
899 factor = 1. + tem - tem1
900 qcko(i,k) = ((1.-tem1)*qcko(i,k-1)+tem*0.5*
901 & (qo(i,k)+qo(i,k-1)))/factor
902 qrcko(i,k) = qcko(i,k)
904 dq = eta(i,k) * (qcko(i,k) - qrch)
910 if(k.ge.kbcon(i).and.dq.gt.0.)
then
911 etah = .5 * (eta(i,k) + eta(i,k-1))
912 dp = 1000. * del(i,k)
913 if(ncloud.gt.0..and.k.gt.jmin(i))
then
914 qlk = dq / (eta(i,k) + etah * (c0 + c1) * dz)
915 dellal(i,k) = etah * c1 * dz * qlk * g / dp
917 qlk = dq / (eta(i,k) + etah * c0 * dz)
919 aa1(i) = aa1(i) - dz * g * qlk
920 qcko(i,k) = qlk + qrch
921 pwo(i,k) = etah * c0 * dz * qlk
922 pwavo(i) = pwavo(i) + pwo(i,k)
924 cnvwt(i,k) = etah * qlk * g / dp
948 if(k.ge.kbcon(i).and.k.lt.ktcon(i))
then
949 dz1 = zo(i,k+1) - zo(i,k)
950 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
951 rfact = 1. + delta * cp * gamma
954 & dz1 * (g / (cp * to(i,k)))
955 & * dbyo(i,k) / (1. + gamma)
960 & max(val,(qeso(i,k) - qo(i,k)))
967 if(cnvflg(i).and.aa1(i).le.0.) cnvflg(i) = .false.
972 totflg = totflg .and. (.not. cnvflg(i))
984 aa2(i) = aafac * aa1(i)
995 if(k.ge.ktcon(i).and.k.lt.kmax(i))
then
996 dz1 = zo(i,k+1) - zo(i,k)
997 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
998 rfact = 1. + delta * cp * gamma
1001 & dz1 * (g / (cp * to(i,k)))
1002 & * dbyo(i,k) / (1. + gamma)
1005 tem = 0.5 * (zi(i,ktcon(i))-zi(i,kbcon(i)))
1006 tem1 = zi(i,k)-zi(i,ktcon(i))
1007 if(aa2(i) < 0. .or. tem1 >= tem)
then
1023 if(k.ge.ktcon(i).and.k.lt.ktcon1(i))
then
1024 dz = zi(i,k) - zi(i,k-1)
1025 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
1027 & + gamma * dbyo(i,k) / (hvap * (1. + gamma))
1029 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1)) * dz
1030 tem1 = 0.5 * xlamud(i) * dz
1031 factor = 1. + tem - tem1
1032 qcko(i,k) = ((1.-tem1)*qcko(i,k-1)+tem*0.5*
1033 & (qo(i,k)+qo(i,k-1)))/factor
1034 qrcko(i,k) = qcko(i,k)
1036 dq = eta(i,k) * (qcko(i,k) - qrch)
1041 etah = .5 * (eta(i,k) + eta(i,k-1))
1042 dp = 1000. * del(i,k)
1043 if(ncloud.gt.0.)
then
1044 qlk = dq / (eta(i,k) + etah * (c0 + c1) * dz)
1045 dellal(i,k) = etah * c1 * dz * qlk * g / dp
1047 qlk = dq / (eta(i,k) + etah * c0 * dz)
1049 qcko(i,k) = qlk + qrch
1050 pwo(i,k) = etah * c0 * dz * qlk
1051 pwavo(i) = pwavo(i) + pwo(i,k)
1053 cnvwt(i,k) = etah * qlk * g / dp
1066 ktcon(i) = ktcon1(i)
1074 if(ncloud.gt.0)
then
1081 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
1083 & + gamma * dbyo(i,k) / (hvap * (1. + gamma))
1084 dq = qcko(i,k) - qrch
1118 if(k.gt.kb(i).and.k.le.ktcon(i))
then
1119 shear= sqrt((uo(i,k)-uo(i,k-1)) ** 2
1120 & + (vo(i,k)-vo(i,k-1)) ** 2)
1121 vshear(i) = vshear(i) + shear
1128 vshear(i) = 1.e3 * vshear(i) / (zi(i,ktcon(i))-zi(i,kb(i)))
1129 e1=1.591-.639*vshear(i)
1130 & +.0953*(vshear(i)**2)-.00496*(vshear(i)**3)
1133 edt(i) = min(edt(i),val)
1135 edt(i) = max(edt(i),val)
1155 if(cnvflg(i).and.k.ge.1.and.k.lt.kbcon(i))
then
1156 dz = zi(i,k+1) - zi(i,k)
1157 sumx(i) = sumx(i) + dz
1163 if(islimsk(i) == 1) beta = betal
1165 dz = (sumx(i)+zi(i,1))/float(kbcon(i))
1166 tem = 1./float(kbcon(i))
1167 xlamd(i) = (1.-beta**tem)/dz
1176 if (cnvflg(i) .and. k .le. kmax(i)-1)
then
1177 if(k.lt.jmin(i).and.k.ge.kbcon(i))
then
1178 dz = zi(i,k+1) - zi(i,k)
1179 ptem = xlamdd - xlamde
1180 etad(i,k) = etad(i,k+1) * (1. - ptem * dz)
1181 else if(k.lt.kbcon(i))
then
1182 dz = zi(i,k+1) - zi(i,k)
1183 ptem = xlamd(i) + xlamdd - xlamde
1184 etad(i,k) = etad(i,k+1) * (1. - ptem * dz)
1196 hcdo(i,jmn) = heo(i,jmn)
1197 qcdo(i,jmn) = qo(i,jmn)
1198 qrcdo(i,jmn)= qo(i,jmn)
1199 ucdo(i,jmn) = uo(i,jmn)
1200 vcdo(i,jmn) = vo(i,jmn)
1208 if (cnvflg(i) .and. k.lt.jmin(i))
then
1209 dz = zi(i,k+1) - zi(i,k)
1210 if(k.ge.kbcon(i))
then
1212 tem1 = 0.5 * xlamdd * dz
1215 tem1 = 0.5 * (xlamd(i)+xlamdd) * dz
1217 factor = 1. + tem - tem1
1218 ptem = 0.5 * tem - pgcon
1219 ptem1= 0.5 * tem + pgcon
1220 hcdo(i,k) = ((1.-tem1)*hcdo(i,k+1)+tem*0.5*
1221 & (heo(i,k)+heo(i,k+1)))/factor
1222 ucdo(i,k) = ((1.-tem1)*ucdo(i,k+1)+ptem*uo(i,k+1)
1223 & +ptem1*uo(i,k))/factor
1224 vcdo(i,k) = ((1.-tem1)*vcdo(i,k+1)+ptem*vo(i,k+1)
1225 & +ptem1*vo(i,k))/factor
1226 dbyo(i,k) = hcdo(i,k) - heso(i,k)
1234 if (cnvflg(i).and.k.lt.jmin(i))
then
1235 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
1236 qrcdo(i,k) = qeso(i,k)+
1237 & (1./hvap)*(gamma/(1.+gamma))*dbyo(i,k)
1240 dz = zi(i,k+1) - zi(i,k)
1241 if(k.ge.kbcon(i))
then
1243 tem1 = 0.5 * xlamdd * dz
1246 tem1 = 0.5 * (xlamd(i)+xlamdd) * dz
1248 factor = 1. + tem - tem1
1249 qcdo(i,k) = ((1.-tem1)*qrcdo(i,k+1)+tem*0.5*
1250 & (qo(i,k)+qo(i,k+1)))/factor
1257 pwdo(i,k) = etad(i,k) * (qcdo(i,k) - qrcdo(i,k))
1258 pwevo(i) = pwevo(i) + pwdo(i,k)
1270 if(islimsk(i) == 0) edtmax = edtmaxs
1272 if(pwevo(i).lt.0.)
then
1273 edto(i) = -edto(i) * pwavo(i) / pwevo(i)
1274 edto(i) = min(edto(i),edtmax)
1286 if (cnvflg(i) .and. k .lt. jmin(i))
then
1287 gamma = el2orc * qeso(i,k) / to(i,k)**2
1292 dz=-1.*(zo(i,k+1)-zo(i,k))
1293 aa1(i)=aa1(i)+edto(i)*dz*(g/(cp*dt))*((dhh-dh)/(1.+dg))
1294 & *(1.+delta*cp*dg*dt/hvap)
1296 aa1(i)=aa1(i)+edto(i)*
1297 & dz*g*delta*max(val,(qeso(i,k)-qo(i,k)))
1303 if(cnvflg(i).and.aa1(i).le.0.)
then
1310 totflg = totflg .and. (.not. cnvflg(i))
1321 if(cnvflg(i) .and. k .le. kmax(i))
then
1331 dp = 1000. * del(i,1)
1332 dellah(i,1) = edto(i) * etad(i,1) * (hcdo(i,1)
1333 & - heo(i,1)) * g / dp
1334 dellaq(i,1) = edto(i) * etad(i,1) * (qrcdo(i,1)
1335 & - qo(i,1)) * g / dp
1336 dellau(i,1) = edto(i) * etad(i,1) * (ucdo(i,1)
1337 & - uo(i,1)) * g / dp
1338 dellav(i,1) = edto(i) * etad(i,1) * (vcdo(i,1)
1339 & - vo(i,1)) * g / dp
1347 if (cnvflg(i).and.k.lt.ktcon(i))
then
1349 if(k.le.kb(i)) aup = 0.
1351 if(k.gt.jmin(i)) adw = 0.
1352 dp = 1000. * del(i,k)
1353 dz = zi(i,k) - zi(i,k-1)
1356 dv2h = .5 * (heo(i,k) + heo(i,k-1))
1359 dv2q = .5 * (qo(i,k) + qo(i,k-1))
1362 dv2u = .5 * (uo(i,k) + uo(i,k-1))
1365 dv2v = .5 * (vo(i,k) + vo(i,k-1))
1368 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1))
1371 if(k.le.kbcon(i))
then
1373 ptem1 = xlamd(i)+xlamdd
1379 dellah(i,k) = dellah(i,k) +
1380 & ((aup*eta(i,k)-adw*edto(i)*etad(i,k))*dv1h
1381 & - (aup*eta(i,k-1)-adw*edto(i)*etad(i,k-1))*dv3h
1382 & - (aup*tem*eta(i,k-1)+adw*edto(i)*ptem*etad(i,k))*dv2h*dz
1383 & + aup*tem1*eta(i,k-1)*.5*(hcko(i,k)+hcko(i,k-1))*dz
1384 & + adw*edto(i)*ptem1*etad(i,k)*.5*(hcdo(i,k)+hcdo(i,k-1))*dz
1387 dellaq(i,k) = dellaq(i,k) +
1388 & ((aup*eta(i,k)-adw*edto(i)*etad(i,k))*dv1q
1389 & - (aup*eta(i,k-1)-adw*edto(i)*etad(i,k-1))*dv3q
1390 & - (aup*tem*eta(i,k-1)+adw*edto(i)*ptem*etad(i,k))*dv2q*dz
1391 & + aup*tem1*eta(i,k-1)*.5*(qrcko(i,k)+qcko(i,k-1))*dz
1392 & + adw*edto(i)*ptem1*etad(i,k)*.5*(qrcdo(i,k)+qcdo(i,k-1))*dz
1395 dellau(i,k) = dellau(i,k) +
1396 & ((aup*eta(i,k)-adw*edto(i)*etad(i,k))*dv1u
1397 & - (aup*eta(i,k-1)-adw*edto(i)*etad(i,k-1))*dv3u
1398 & - (aup*tem*eta(i,k-1)+adw*edto(i)*ptem*etad(i,k))*dv2u*dz
1399 & + aup*tem1*eta(i,k-1)*.5*(ucko(i,k)+ucko(i,k-1))*dz
1400 & + adw*edto(i)*ptem1*etad(i,k)*.5*(ucdo(i,k)+ucdo(i,k-1))*dz
1401 & - pgcon*(aup*eta(i,k-1)-adw*edto(i)*etad(i,k))*(dv1u-dv3u)
1404 dellav(i,k) = dellav(i,k) +
1405 & ((aup*eta(i,k)-adw*edto(i)*etad(i,k))*dv1v
1406 & - (aup*eta(i,k-1)-adw*edto(i)*etad(i,k-1))*dv3v
1407 & - (aup*tem*eta(i,k-1)+adw*edto(i)*ptem*etad(i,k))*dv2v*dz
1408 & + aup*tem1*eta(i,k-1)*.5*(vcko(i,k)+vcko(i,k-1))*dz
1409 & + adw*edto(i)*ptem1*etad(i,k)*.5*(vcdo(i,k)+vcdo(i,k-1))*dz
1410 & - pgcon*(aup*eta(i,k-1)-adw*edto(i)*etad(i,k))*(dv1v-dv3v)
1422 dp = 1000. * del(i,indx)
1423 dv1h = heo(i,indx-1)
1424 dellah(i,indx) = eta(i,indx-1) *
1425 & (hcko(i,indx-1) - dv1h) * g / dp
1427 dellaq(i,indx) = eta(i,indx-1) *
1428 & (qcko(i,indx-1) - dv1q) * g / dp
1430 dellau(i,indx) = eta(i,indx-1) *
1431 & (ucko(i,indx-1) - dv1u) * g / dp
1433 dellav(i,indx) = eta(i,indx-1) *
1434 & (vcko(i,indx-1) - dv1v) * g / dp
1438 dellal(i,indx) = eta(i,indx-1) *
1439 & qlko_ktcon(i) * g / dp
1448 if (cnvflg(i).and.k .le. kmax(i))
then
1449 if(k.gt.ktcon(i))
then
1453 if(k.le.ktcon(i))
then
1454 qo(i,k) = dellaq(i,k) * mbdt + q1(i,k)
1455 dellat = (dellah(i,k) - hvap * dellaq(i,k)) / cp
1456 to(i,k) = dellat * mbdt + t1(i,k)
1458 qo(i,k) = max(qo(i,k), val )
1479 if(cnvflg(i) .and. k .le. kmax(i))
then
1480 qeso(i,k) = 0.01 * fpvs(to(i,k))
1481 qeso(i,k) = eps * qeso(i,k) / (pfld(i,k)+epsm1*qeso(i,k))
1483 qeso(i,k) = max(qeso(i,k), val )
1494 if(cnvflg(i) .and. k .le. kmax(i)-1)
then
1495 dz = .5 * (zo(i,k+1) - zo(i,k))
1496 dp = .5 * (pfld(i,k+1) - pfld(i,k))
1497 es = 0.01 * fpvs(to(i,k+1))
1498 pprime = pfld(i,k+1) + epsm1 * es
1499 qs = eps * es / pprime
1500 dqsdp = - qs / pprime
1501 desdt = es * (fact1 / to(i,k+1) + fact2 / (to(i,k+1)**2))
1502 dqsdt = qs * pfld(i,k+1) * desdt / (es * pprime)
1503 gamma = el2orc * qeso(i,k+1) / (to(i,k+1)**2)
1504 dt = (g * dz + hvap * dqsdp * dp) / (cp * (1. + gamma))
1505 dq = dqsdt * dt + dqsdp * dp
1506 to(i,k) = to(i,k+1) + dt
1507 qo(i,k) = qo(i,k+1) + dq
1508 po(i,k) = .5 * (pfld(i,k) + pfld(i,k+1))
1514 if(cnvflg(i) .and. k .le. kmax(i)-1)
then
1515 qeso(i,k) = 0.01 * fpvs(to(i,k))
1516 qeso(i,k) = eps * qeso(i,k) / (po(i,k) + epsm1 * qeso(i,k))
1518 qeso(i,k) = max(qeso(i,k), val1)
1520 qo(i,k) = max(qo(i,k), val2 )
1522 heo(i,k) = .5 * g * (zo(i,k) + zo(i,k+1)) +
1523 & cp * to(i,k) + hvap * qo(i,k)
1524 heso(i,k) = .5 * g * (zo(i,k) + zo(i,k+1)) +
1525 & cp * to(i,k) + hvap * qeso(i,k)
1532 heo(i,k) = g * zo(i,k) + cp * to(i,k) + hvap * qo(i,k)
1533 heso(i,k) = g * zo(i,k) + cp * to(i,k) + hvap * qeso(i,k)
1553 hcko(i,indx) = heo(i,indx)
1554 qcko(i,indx) = qo(i,indx)
1560 if(k.gt.kb(i).and.k.le.ktcon(i))
then
1561 dz = zi(i,k) - zi(i,k-1)
1562 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1)) * dz
1563 tem1 = 0.5 * xlamud(i) * dz
1564 factor = 1. + tem - tem1
1565 hcko(i,k) = ((1.-tem1)*hcko(i,k-1)+tem*0.5*
1566 & (heo(i,k)+heo(i,k-1)))/factor
1574 if(k.gt.kb(i).and.k.lt.ktcon(i))
then
1575 dz = zi(i,k) - zi(i,k-1)
1576 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
1577 xdby = hcko(i,k) - heso(i,k)
1579 & + gamma * xdby / (hvap * (1. + gamma))
1581 tem = 0.5 * (xlamue(i,k)+xlamue(i,k-1)) * dz
1582 tem1 = 0.5 * xlamud(i) * dz
1583 factor = 1. + tem - tem1
1584 qcko(i,k) = ((1.-tem1)*qcko(i,k-1)+tem*0.5*
1585 & (qo(i,k)+qo(i,k-1)))/factor
1587 dq = eta(i,k) * (qcko(i,k) - xqrch)
1589 if(k.ge.kbcon(i).and.dq.gt.0.)
then
1590 etah = .5 * (eta(i,k) + eta(i,k-1))
1591 if(ncloud.gt.0..and.k.gt.jmin(i))
then
1592 qlk = dq / (eta(i,k) + etah * (c0 + c1) * dz)
1594 qlk = dq / (eta(i,k) + etah * c0 * dz)
1596 if(k.lt.ktcon1(i))
then
1597 xaa0(i) = xaa0(i) - dz * g * qlk
1599 qcko(i,k) = qlk + xqrch
1600 xpw = etah * c0 * dz * qlk
1601 xpwav(i) = xpwav(i) + xpw
1604 if(k.ge.kbcon(i).and.k.lt.ktcon1(i))
then
1605 dz1 = zo(i,k+1) - zo(i,k)
1606 gamma = el2orc * qeso(i,k) / (to(i,k)**2)
1607 rfact = 1. + delta * cp * gamma
1610 & + dz1 * (g / (cp * to(i,k)))
1611 & * xdby / (1. + gamma)
1616 & max(val,(qeso(i,k) - qo(i,k)))
1630 hcdo(i,jmn) = heo(i,jmn)
1631 qcdo(i,jmn) = qo(i,jmn)
1632 qrcd(i,jmn) = qo(i,jmn)
1639 if (cnvflg(i) .and. k.lt.jmin(i))
then
1640 dz = zi(i,k+1) - zi(i,k)
1641 if(k.ge.kbcon(i))
then
1643 tem1 = 0.5 * xlamdd * dz
1646 tem1 = 0.5 * (xlamd(i)+xlamdd) * dz
1648 factor = 1. + tem - tem1
1649 hcdo(i,k) = ((1.-tem1)*hcdo(i,k+1)+tem*0.5*
1650 & (heo(i,k)+heo(i,k+1)))/factor
1657 if (cnvflg(i) .and. k .lt. jmin(i))
then
1660 gamma = el2orc * dq / dt**2
1661 dh = hcdo(i,k) - heso(i,k)
1662 qrcd(i,k)=dq+(1./hvap)*(gamma/(1.+gamma))*dh
1665 dz = zi(i,k+1) - zi(i,k)
1666 if(k.ge.kbcon(i))
then
1668 tem1 = 0.5 * xlamdd * dz
1671 tem1 = 0.5 * (xlamd(i)+xlamdd) * dz
1673 factor = 1. + tem - tem1
1674 qcdo(i,k) = ((1.-tem1)*qrcd(i,k+1)+tem*0.5*
1675 & (qo(i,k)+qo(i,k+1)))/factor
1682 xpwd = etad(i,k) * (qcdo(i,k) - qrcd(i,k))
1683 xpwev(i) = xpwev(i) + xpwd
1690 if(islimsk(i) == 0) edtmax = edtmaxs
1692 if(xpwev(i).ge.0.)
then
1695 edtx(i) = -edtx(i) * xpwav(i) / xpwev(i)
1696 edtx(i) = min(edtx(i),edtmax)
1707 if (cnvflg(i) .and. k.lt.jmin(i))
then
1708 gamma = el2orc * qeso(i,k) / to(i,k)**2
1713 dz=-1.*(zo(i,k+1)-zo(i,k))
1714 xaa0(i)=xaa0(i)+edtx(i)*dz*(g/(cp*dt))*((dhh-dh)/(1.+dg))
1715 & *(1.+delta*cp*dg*dt/hvap)
1717 xaa0(i)=xaa0(i)+edtx(i)*
1718 & dz*g*delta*max(val,(qeso(i,k)-qo(i,k)))
1729 if(pfld(i,ktcon(i)).lt.pcrit(15))
then
1730 acrt(i)=acrit(15)*(975.-pfld(i,ktcon(i)))
1732 else if(pfld(i,ktcon(i)).gt.pcrit(1))
then
1735 k = int((850. - pfld(i,ktcon(i)))/50.) + 2
1738 acrt(i)=acrit(k)+(acrit(k-1)-acrit(k))*
1739 & (pfld(i,ktcon(i))-pcrit(k))/(pcrit(k-1)-pcrit(k))
1746 if(islimsk(i) == 1)
then
1760 if(pdot(i).le.w4)
then
1761 acrtfct(i) = (pdot(i) - w4) / (w3 - w4)
1762 elseif(pdot(i).ge.-w4)
then
1763 acrtfct(i) = - (pdot(i) + w4) / (w4 - w3)
1768 acrtfct(i) = max(acrtfct(i),val1)
1770 acrtfct(i) = min(acrtfct(i),val2)
1771 acrtfct(i) = 1. - acrtfct(i)
1782 dtconv(i) = dt2 + max((1800. - dt2),0.) *
1783 & (pdot(i) - w2) / (w1 - w2)
1786 dtconv(i) = max(dtconv(i),dtmin)
1787 dtconv(i) = min(dtconv(i),dtmax)
1801 fld(i)=(aa1(i)-acrt(i)* acrtfct(i))/dtconv(i)
1802 if(fld(i).le.0.) cnvflg(i) = .false.
1811 xk(i) = (xaa0(i) - aa1(i)) / mbdt
1812 if(xk(i).ge.0.) cnvflg(i) = .false.
1822 xmb(i) = -fld(i) / xk(i)
1823 xmb(i) = min(xmb(i),xmbmax(i))
1830 totflg = totflg .and. (.not. cnvflg(i))
1840 if (cnvflg(i) .and. k .le. kmax(i))
then
1845 qeso(i,k) = 0.01 * fpvs(t1(i,k))
1846 qeso(i,k) = eps * qeso(i,k) / (pfld(i,k) + epsm1*qeso(i,k))
1848 qeso(i,k) = max(qeso(i,k), val )
1872 if (cnvflg(i) .and. k .le. kmax(i))
then
1873 if(k.le.ktcon(i))
then
1874 dellat = (dellah(i,k) - hvap * dellaq(i,k)) / cp
1875 t1(i,k) = t1(i,k) + dellat * xmb(i) * dt2
1876 q1(i,k) = q1(i,k) + dellaq(i,k) * xmb(i) * dt2
1880 u1(i,k) = u1(i,k) + dellau(i,k) * xmb(i) * dt2
1881 v1(i,k) = v1(i,k) + dellav(i,k) * xmb(i) * dt2
1882 dp = 1000. * del(i,k)
1883 delhbar(i) = delhbar(i) + dellah(i,k)*xmb(i)*dp/g
1884 delqbar(i) = delqbar(i) + dellaq(i,k)*xmb(i)*dp/g
1885 deltbar(i) = deltbar(i) + dellat*xmb(i)*dp/g
1886 delubar(i) = delubar(i) + dellau(i,k)*xmb(i)*dp/g
1887 delvbar(i) = delvbar(i) + dellav(i,k)*xmb(i)*dp/g
1895 if (cnvflg(i) .and. k .le. kmax(i))
then
1896 if(k.le.ktcon(i))
then
1897 qeso(i,k) = 0.01 * fpvs(t1(i,k))
1898 qeso(i,k) = eps * qeso(i,k)/(pfld(i,k) + epsm1*qeso(i,k))
1900 qeso(i,k) = max(qeso(i,k), val )
1915 if (cnvflg(i) .and. k .le. kmax(i))
then
1916 if(k.lt.ktcon(i))
then
1918 if(k.le.kb(i)) aup = 0.
1920 if(k.ge.jmin(i)) adw = 0.
1921 rain = aup * pwo(i,k) + adw * edto(i) * pwdo(i,k)
1922 rntot(i) = rntot(i) + rain * xmb(i) * .001 * dt2
1932 if (k .le. kmax(i))
then
1936 if(cnvflg(i).and.k.lt.ktcon(i))
then
1938 if(k.le.kb(i)) aup = 0.
1940 if(k.ge.jmin(i)) adw = 0.
1941 rain = aup * pwo(i,k) + adw * edto(i) * pwdo(i,k)
1942 rn(i) = rn(i) + rain * xmb(i) * .001 * dt2
1944 if(flg(i).and.k.lt.ktcon(i))
then
1945 evef = edt(i) * evfact
1946 if(islimsk(i) == 1) evef=edt(i) * evfactl
1949 qcond(i) = evef * (q1(i,k) - qeso(i,k))
1950 & / (1. + el2orc * qeso(i,k) / t1(i,k)**2)
1951 dp = 1000. * del(i,k)
1952 if(rn(i).gt.0..and.qcond(i).lt.0.)
then
1953 qevap(i) = -qcond(i) * (1.-exp(-.32*sqrt(dt2*rn(i))))
1954 qevap(i) = min(qevap(i), rn(i)*1000.*g/dp)
1955 delq2(i) = delqev(i) + .001 * qevap(i) * dp / g
1957 if(rn(i).gt.0..and.qcond(i).lt.0..and.
1958 & delq2(i).gt.rntot(i))
then
1959 qevap(i) = 1000.* g * (rntot(i) - delqev(i)) / dp
1962 if(rn(i).gt.0..and.qevap(i).gt.0.)
then
1963 q1(i,k) = q1(i,k) + qevap(i)
1964 t1(i,k) = t1(i,k) - elocp * qevap(i)
1965 rn(i) = rn(i) - .001 * qevap(i) * dp / g
1966 deltv(i) = - elocp*qevap(i)/dt2
1967 delq(i) = + qevap(i)/dt2
1968 delqev(i) = delqev(i) + .001*dp*qevap(i)/g
1970 dellaq(i,k) = dellaq(i,k) + delq(i) / xmb(i)
1971 delqbar(i) = delqbar(i) + delq(i)*dp/g
1972 deltbar(i) = deltbar(i) + deltv(i)*dp/g
2000 if(rn(i).lt.0..and..not.flg(i)) rn(i) = 0.
2001 if(rn(i).le.0.)
then
2017 if (cnvflg(i) .and. rn(i).gt.0.)
then
2018 if (k.ge.kbcon(i).and.k.lt.ktcon(i))
then
2019 cnvw(i,k) = cnvwt(i,k) * xmb(i) * dt2
2030 if (cnvflg(i) .and. rn(i).gt.0.)
then
2031 if (k.ge.kbcon(i).and.k.lt.ktcon(i))
then
2032 cnvc(i,k) = 0.04 * log(1. + 675. * eta(i,k) * xmb(i))
2033 cnvc(i,k) = min(cnvc(i,k), 0.6)
2034 cnvc(i,k) = max(cnvc(i,k), 0.0)
2044 if (ncloud.gt.0)
then
2048 if (cnvflg(i) .and. rn(i).gt.0.)
then
2049 if (k.gt.kb(i).and.k.le.ktcon(i))
then
2050 tem = dellal(i,k) * xmb(i) * dt2
2051 tem1 = max(0.0, min(1.0, (tcr-t1(i,k))*tcrf))
2052 if (qlc(i,k) .gt. -999.0)
then
2053 qli(i,k) = qli(i,k) + tem * tem1
2054 qlc(i,k) = qlc(i,k) + tem *(1.0-tem1)
2056 qli(i,k) = qli(i,k) + tem
2068 if(cnvflg(i).and.rn(i).le.0.)
then
2069 if (k .le. kmax(i))
then
2084 if(cnvflg(i).and.rn(i).gt.0.)
then
2085 if(k.ge.kb(i) .and. k.lt.ktop(i))
then
2086 ud_mf(i,k) = eta(i,k) * xmb(i) * dt2
2093 if(cnvflg(i).and.rn(i).gt.0.)
then
2095 dt_mf(i,k) = ud_mf(i,k)
2101 if(cnvflg(i).and.rn(i).gt.0.)
then
2102 if(k.ge.1 .and. k.le.jmin(i))
then
2103 dd_mf(i,k) = edto(i) * etad(i,k) * xmb(i) * dt2
2109 if(mp_phys == mp_phys_mg)
then
2112 qlcn(i,k) = qlc(i,k)
2113 qicn(i,k) = qli(i,k)
2114 cf_upi(i,k) = cnvc(i,k)
2115 w_upi(i,k) = ud_mf(i,k)*t1(i,k)*rgas /
2116 & (dt2*max(cf_upi(i,k),1.e-12)*prslp(i,k))
2117 cnv_mfd(i,k) = ud_mf(i,k)/dt2
2118 clcn(i,k) = cnvc(i,k)
2119 cnv_fice(i,k) = qicn(i,k)
2120 & / max(1.e-10,qlcn(i,k)+qicn(i,k))