module module_cu_tiedtke use module_model_constants, only:rd=>r_d, rv=>r_v, & & cpd=>cp, alv=>xlv, als=>xls, alf=>xlf, g implicit none real :: rcpd,vtmpc1,t000,hgfr,rhoh2o,tmelt, & c1es,c2es,c3les,c3ies,c4les,c4ies,c5les,c5ies,zrg real :: entrpen,entrscv,entrmid,entrdd,cmfctop,rhm,rhc, & cmfcmax,cmfcmin,cmfdeps,cprcon,crirh,zbuo0, & fdbk,ztau real :: cevapcu1, cevapcu2, zdnoprc parameter( & rcpd=1.0/cpd, & rhoh2o=1.0e03, & tmelt=273.16, & t000= 273.15, & hgfr= 233.15, & zrg=1.0/g, & c1es=610.78, & c2es=c1es*rd/rv, & c3les=17.269, & c4les=35.86, & c5les=c3les*(tmelt-c4les), & c3ies=21.875, & c4ies=7.66, & c5ies=c3ies*(tmelt-c4ies), & vtmpc1=rv/rd-1.0, & cevapcu1=1.93e-6*261.0*0.5/g, & cevapcu2=1.e3/(38.3*0.293) ) parameter(entrpen=1.0e-4) parameter(entrscv=1.2e-3) parameter(entrmid=1.0e-4) parameter(entrdd =2.0e-4) parameter(cmfctop=0.30) parameter(cmfcmax=1.0) parameter(cmfcmin=1.e-10) parameter(cmfdeps=0.30) parameter(cprcon = 1.1e-3/g) parameter(zdnoprc = 1.5e4) parameter(rhc=0.80,rhm=1.0,zbuo0=0.50) parameter(crirh=0.70,fdbk = 1.0,ztau = 2400.0) logical :: lmfpen,lmfmid,lmfscv,lmfdd,lmfdudv parameter(lmfpen=.true.,lmfmid=.true.,lmfscv=.true.,lmfdd=.true.,lmfdudv=.true.) contains subroutine cu_tiedtke( & dt,itimestep,stepcu & ,raincv,pratec,qfx,znu & ,u3d,v3d,w,t3d,qv3d,qc3d,qi3d,pi3d,rho3d & ,qvften,qvpblten & ,dz8w,pcps,p8w,xland,cu_act_flag & ,ids,ide, jds,jde, kds,kde & ,ims,ime, jms,jme, kms,kme & ,its,ite, jts,jte, kts,kte & ,rthcuten,rqvcuten,rqccuten,rqicuten & ,rucuten, rvcuten & ,f_qv ,f_qc ,f_qr ,f_qi ,f_qs & ) implicit none integer, intent(in) :: ids,ide, jds,jde, kds,kde, & ims,ime, jms,jme, kms,kme, & its,ite, jts,jte, kts,kte, & itimestep, & stepcu real, intent(in) :: & dt real, dimension(ims:ime, jms:jme), intent(in) :: & xland real, dimension(ims:ime, jms:jme), intent(inout) :: & raincv, pratec logical, dimension(ims:ime,jms:jme), intent(inout) :: & cu_act_flag real, dimension(ims:ime, kms:kme, jms:jme), intent(in) :: & dz8w, & p8w, & pcps, & pi3d, & qc3d, & qvften, & qvpblten, & qi3d, & qv3d, & rho3d, & t3d, & u3d, & v3d, & w real, dimension(ims:ime, kms:kme, jms:jme), & optional, intent(inout) :: & rqccuten, & rqicuten, & rqvcuten, & rthcuten, & rucuten, & rvcuten logical, optional :: & f_qv & ,f_qc & ,f_qr & ,f_qi & ,f_qs real, dimension(ims:ime, jms:jme) :: & qfx real :: & delt, & rdelt real , dimension(its:ite) :: & rcs, & rn, & evap integer , dimension(its:ite) :: slimsk real , dimension(its:ite, kts:kte+1) :: & prsi real , dimension(its:ite, kts:kte) :: & del, & dot, & phil, & prsl, & q1, & q2, & q3, & q1b, & q1bl, & q11, & q12, & t1, & u1, & v1, & zi, & zl, & omg, & ght integer, dimension(its:ite) :: & kbot, & ktop integer :: & i, & im, & j, & k, & km, & kp, & kx logical :: run_param , doing_adapt_dt , decided integer,dimension( its:ite ) :: ktype real, dimension( kts:kte ) :: sig1 real, dimension( kms:kme ) :: znu integer :: zz do j=jts,jte do i=its,ite cu_act_flag(i,j)=.true. enddo enddo im=ite-its+1 kx=kte-kts+1 delt=dt*stepcu rdelt=1./delt do j=jts,jte do i=its,ite zi(i,kts)=0.0 enddo do k=kts+1,kte km=k-1 do i=its,ite zi(i,k)=zi(i,km)+dz8w(i,km,j) enddo enddo do k=kts+1,kte km=k-1 do i=its,ite zl(i,km)=(zi(i,k)+zi(i,km))*0.5 enddo enddo do i=its,ite zl(i,kte)=2.*zi(i,kte)-zl(i,kte-1) enddo do i=its,ite slimsk(i)=int(abs(xland(i,j)-2.)) enddo do k=kts,kte kp=k+1 do i=its,ite dot(i,k)=-0.5*g*rho3d(i,k,j)*(w(i,k,j)+w(i,kp,j)) enddo enddo do k=kts,kte zz = kte+1-k do i=its,ite u1(i,zz)=u3d(i,k,j) v1(i,zz)=v3d(i,k,j) t1(i,zz)=t3d(i,k,j) q1(i,zz)= qv3d(i,k,j) if(itimestep == 1) then q1b(i,zz)=0. q1bl(i,zz)=0. else q1b(i,zz)=qvften(i,k,j) q1bl(i,zz)=qvpblten(i,k,j) endif q2(i,zz)=qc3d(i,k,j) q3(i,zz)=qi3d(i,k,j) omg(i,zz)=dot(i,k) ght(i,zz)=zl(i,k) prsl(i,zz) = pcps(i,k,j) enddo enddo do k=kts,kte+1 zz = kte+2-k do i=its,ite prsi(i,zz) = p8w(i,k,j) enddo enddo do k=kts,kte zz = kte+1-k sig1(zz) = znu(k) enddo do i=its,ite evap(i) = qfx(i,j) enddo call tiecnv(u1,v1,t1,q1,q2,q3,q1b,q1bl,ght,omg,prsl,prsi,evap, & rn,slimsk,ktype,im,kx,kx+1,sig1,delt) do i=its,ite raincv(i,j)=rn(i)/stepcu pratec(i,j)=rn(i)/(stepcu * dt) enddo do k=kts,kte zz = kte+1-k do i=its,ite rthcuten(i,k,j)=(t1(i,zz)-t3d(i,k,j))/pi3d(i,k,j)*rdelt rqvcuten(i,k,j)=(q1(i,zz)-qv3d(i,k,j))*rdelt rucuten(i,k,j) =(u1(i,zz)-u3d(i,k,j))*rdelt rvcuten(i,k,j) =(v1(i,zz)-v3d(i,k,j))*rdelt enddo enddo if(present(rqccuten))then if ( f_qc ) then do k=kts,kte zz = kte+1-k do i=its,ite rqccuten(i,k,j)=(q2(i,zz)-qc3d(i,k,j))*rdelt enddo enddo endif endif if(present(rqicuten))then if ( f_qi ) then do k=kts,kte zz = kte+1-k do i=its,ite rqicuten(i,k,j)=(q3(i,zz)-qi3d(i,k,j))*rdelt enddo enddo endif endif enddo end subroutine cu_tiedtke subroutine tiedtkeinit(rthcuten,rqvcuten,rqccuten,rqicuten, & rucuten,rvcuten, & restart,p_qc,p_qi,p_first_scalar, & allowed_to_read, & ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte) implicit none logical , intent(in) :: allowed_to_read,restart integer , intent(in) :: ids, ide, jds, jde, kds, kde, & ims, ime, jms, jme, kms, kme, & its, ite, jts, jte, kts, kte integer , intent(in) :: p_first_scalar, p_qi, p_qc real, dimension( ims:ime , kms:kme , jms:jme ) , intent(out) :: & rthcuten, & rqvcuten, & rqccuten, & rqicuten, & rucuten,rvcuten integer :: i, j, k, itf, jtf, ktf jtf=min0(jte,jde-1) ktf=min0(kte,kde-1) itf=min0(ite,ide-1) if(.not.restart)then do j=jts,jtf do k=kts,ktf do i=its,itf rthcuten(i,k,j)=0. rqvcuten(i,k,j)=0. rucuten(i,k,j)=0. rvcuten(i,k,j)=0. enddo enddo enddo if (p_qc .ge. p_first_scalar) then do j=jts,jtf do k=kts,ktf do i=its,itf rqccuten(i,k,j)=0. enddo enddo enddo endif if (p_qi .ge. p_first_scalar) then do j=jts,jtf do k=kts,ktf do i=its,itf rqicuten(i,k,j)=0. enddo enddo enddo endif endif end subroutine tiedtkeinit subroutine tiecnv(pu,pv,pt,pqv,pqc,pqi,pqvf,pqvbl,poz,pomg, & pap,paph,evap,zprecc,lndj,ktype,lq,km,km1,sig1,dt) implicit none real pu(lq,km),pv(lq,km),pt(lq,km),pqv(lq,km),pqvf(lq,km) real poz(lq,km),pomg(lq,km),evap(lq),zprecc(lq),pqvbl(lq,km) real pum1(lq,km), pvm1(lq,km), & ptte(lq,km), pqte(lq,km), pvom(lq,km), pvol(lq,km), & pverv(lq,km), pgeo(lq,km), pap(lq,km), paph(lq,km1) real pqhfl(lq), zqq(lq,km), paprc(lq), paprs(lq), & prsfc(lq), pssfc(lq), paprsm(lq), pcte(lq,km) real ztp1(lq,km), zqp1(lq,km), ztu(lq,km), zqu(lq,km), & zlu(lq,km), zlude(lq,km), zmfu(lq,km), zmfd(lq,km), & zqsat(lq,km), pqc(lq,km), pqi(lq,km), zrain(lq) real sig(km1),sig1(km) integer icbot(lq), ictop(lq), ktype(lq), lndj(lq) real dt logical locum(lq) real psheat,psrain,psevap,psmelt,psdiss,tt real ztmst,ztpp1,fliq,fice,ztc,zalf integer i,j,k,lq,lp,km,km1 ztmst=dt psheat=0.0 psrain=0.0 psevap=0.0 psmelt=0.0 psdiss=0.0 do j=1,lq zrain(j)=0.0 locum(j)=.false. prsfc(j)=0.0 pssfc(j)=0.0 paprc(j)=0.0 paprs(j)=0.0 paprsm(j)=0.0 pqhfl(j)=evap(j) end do do k=1,km do j=1,lq ptte(j,k)=0.0 pcte(j,k)=0.0 pvom(j,k)=0.0 pvol(j,k)=0.0 ztp1(j,k)=pt(j,k) zqp1(j,k)=pqv(j,k)/(1.0+pqv(j,k)) pum1(j,k)=pu(j,k) pvm1(j,k)=pv(j,k) pverv(j,k)=pomg(j,k) pgeo(j,k)=g*poz(j,k) tt=ztp1(j,k) zqsat(j,k)=tlucua(tt)/pap(j,k) zqsat(j,k)=min(0.5,zqsat(j,k)) zqsat(j,k)=zqsat(j,k)/(1.-vtmpc1*zqsat(j,k)) pqte(j,k)=pqvf(j,k)+pqvbl(j,k) zqq(j,k)=pqte(j,k) end do end do call cumastr_new & (lq, km, km1, km-1, ztp1, & zqp1, pum1, pvm1, pverv, zqsat, & pqhfl, ztmst, pap, paph, pgeo, & ptte, pqte, pvom, pvol, prsfc, & pssfc, paprc, paprsm, paprs, locum, & ktype, icbot, ictop, ztu, zqu, & zlu, zlude, zmfu, zmfd, zrain, & psrain, psevap, psheat, psdiss, psmelt, & pcte, sig1, lndj) if(fdbk.ge.1.0e-9) then do k=1,km do j=1,lq if(pcte(j,k).gt.0.0) then ztpp1=pt(j,k)+ptte(j,k)*ztmst if(ztpp1.ge.t000) then fliq=1.0 zalf=0.0 else if(ztpp1.le.hgfr) then fliq=0.0 zalf=alf else ztc=ztpp1-t000 fliq=0.0059+0.9941*exp(-0.003102*ztc*ztc) zalf=alf endif fice=1.0-fliq pqc(j,k)=pqc(j,k)+fliq*pcte(j,k)*ztmst pqi(j,k)=pqi(j,k)+fice*pcte(j,k)*ztmst ptte(j,k)=ptte(j,k)-zalf*rcpd*fliq*pcte(j,k) endif end do end do endif do k=1,km do j=1,lq pt(j,k)=ztp1(j,k)+ptte(j,k)*ztmst zqp1(j,k)=zqp1(j,k)+(pqte(j,k)-zqq(j,k))*ztmst pqv(j,k)=zqp1(j,k)/(1.0-zqp1(j,k)) end do end do do j=1,lq zprecc(j)=amax1(0.0,(prsfc(j)+pssfc(j))*ztmst) end do if (lmfdudv) then do k=1,km do j=1,lq pu(j,k)=pu(j,k)+pvom(j,k)*ztmst pv(j,k)=pv(j,k)+pvol(j,k)*ztmst end do end do endif return end subroutine tiecnv subroutine cumastr_new & (klon, klev, klevp1, klevm1, pten, & pqen, puen, pven, pverv, pqsen, & pqhfl, ztmst, pap, paph, pgeo, & ptte, pqte, pvom, pvol, prsfc, & pssfc, paprc, paprsm, paprs, ldcum, & ktype, kcbot, kctop, ptu, pqu, & plu, plude, pmfu, pmfd, prain, & psrain, psevap, psheat, psdiss, psmelt,& pcte, sig1, lndj) implicit none integer klon, klev, klevp1 integer klevm1 real ztmst real psrain, psevap, psheat, psdiss, psmelt, zcons2 integer jk,jl,ikb real zqumqe, zdqmin, zmfmax, zalvdcp, zqalv real zhsat, zgam, zzz, zhhat, zbi, zro, zdz, zdhdz, zdepth real zfac, zrh, zpbmpt, dept, zht, zeps integer icum, itopm2 real pten(klon,klev), pqen(klon,klev), & puen(klon,klev), pven(klon,klev), & ptte(klon,klev), pqte(klon,klev), & pvom(klon,klev), pvol(klon,klev), & pqsen(klon,klev), pgeo(klon,klev), & pap(klon,klev), paph(klon,klevp1),& pverv(klon,klev), pqhfl(klon) real ptu(klon,klev), pqu(klon,klev), & plu(klon,klev), plude(klon,klev), & pmfu(klon,klev), pmfd(klon,klev), & paprc(klon), paprs(klon), & paprsm(klon), prain(klon), & prsfc(klon), pssfc(klon) real ztenh(klon,klev), zqenh(klon,klev),& zgeoh(klon,klev), zqsenh(klon,klev),& ztd(klon,klev), zqd(klon,klev), & zmfus(klon,klev), zmfds(klon,klev), & zmfuq(klon,klev), zmfdq(klon,klev), & zdmfup(klon,klev), zdmfdp(klon,klev),& zmful(klon,klev), zrfl(klon), & zuu(klon,klev), zvu(klon,klev), & zud(klon,klev), zvd(klon,klev) real zentr(klon), zhcbase(klon), & zmfub(klon), zmfub1(klon), & zdqpbl(klon), zdqcv(klon) real zsfl(klon), zdpmel(klon,klev), & pcte(klon,klev), zcape(klon), & zheat(klon), zhhatt(klon,klev), & zhmin(klon), zrelh(klon) real sig1(klev) integer ilab(klon,klev), idtop(klon), & ictop0(klon), ilwmin(klon) integer kcbot(klon), kctop(klon), & ktype(klon), ihmin(klon), & ktop0, lndj(klon) logical ldcum(klon) logical loddraf(klon), llo1 zcons2=1./(g*ztmst) call cuini & (klon, klev, klevp1, klevm1, pten, & pqen, pqsen, puen, pven, pverv, & pgeo, paph, zgeoh, ztenh, zqenh, & zqsenh, ilwmin, ptu, pqu, ztd, & zqd, zuu, zvu, zud, zvd, & pmfu, pmfd, zmfus, zmfds, zmfuq, & zmfdq, zdmfup, zdmfdp, zdpmel, plu, & plude, ilab) call cubase & (klon, klev, klevp1, klevm1, ztenh, & zqenh, zgeoh, paph, ptu, pqu, & plu, puen, pven, zuu, zvu, & ldcum, kcbot, ilab) jk=1 do jl=1,klon zdqcv(jl) =pqte(jl,jk)*(paph(jl,jk+1)-paph(jl,jk)) zdqpbl(jl)=0.0 idtop(jl)=0 end do do jk=2,klev do jl=1,klon zdqcv(jl)=zdqcv(jl)+pqte(jl,jk)*(paph(jl,jk+1)-paph(jl,jk)) if(jk.ge.kcbot(jl)) zdqpbl(jl)=zdqpbl(jl)+pqte(jl,jk) & *(paph(jl,jk+1)-paph(jl,jk)) end do end do do jl=1,klon ktype(jl)=0 if(zdqcv(jl).gt.max(0.,1.1*pqhfl(jl)*g)) then ktype(jl)=1 else ktype(jl)=2 endif ikb=kcbot(jl) zqumqe=pqu(jl,ikb)+plu(jl,ikb)-zqenh(jl,ikb) zdqmin=max(0.01*zqenh(jl,ikb),1.e-10) if(zdqpbl(jl).gt.0..and.zqumqe.gt.zdqmin.and.ldcum(jl)) then zmfub(jl)=zdqpbl(jl)/(g*max(zqumqe,zdqmin)) else zmfub(jl)=0.01 ldcum(jl)=.false. endif zmfmax=(paph(jl,ikb)-paph(jl,ikb-1))*zcons2 zmfub(jl)=min(zmfub(jl),zmfmax) ikb=kcbot(jl) zhcbase(jl)=cpd*ptu(jl,ikb)+zgeoh(jl,ikb)+alv*pqu(jl,ikb) ictop0(jl)=kcbot(jl)-1 end do zalvdcp=alv/cpd zqalv=1./alv do jk=klevm1,3,-1 do jl=1,klon zhsat=cpd*ztenh(jl,jk)+zgeoh(jl,jk)+alv*zqsenh(jl,jk) zgam=c5les*zalvdcp*zqsenh(jl,jk)/ & ((1.-vtmpc1*zqsenh(jl,jk))*(ztenh(jl,jk)-c4les)**2) zzz=cpd*ztenh(jl,jk)*0.608 zhhat=zhsat-(zzz+zgam*zzz)/(1.+zgam*zzz*zqalv)* & max(zqsenh(jl,jk)-zqenh(jl,jk),0.) zhhatt(jl,jk)=zhhat if(jk.lt.ictop0(jl).and.zhcbase(jl).gt.zhhat) ictop0(jl)=jk end do end do do jl=1,klon jk=kcbot(jl) zhsat=cpd*ztenh(jl,jk)+zgeoh(jl,jk)+alv*zqsenh(jl,jk) zgam=c5les*zalvdcp*zqsenh(jl,jk)/ & ((1.-vtmpc1*zqsenh(jl,jk))*(ztenh(jl,jk)-c4les)**2) zzz=cpd*ztenh(jl,jk)*0.608 zhhat=zhsat-(zzz+zgam*zzz)/(1.+zgam*zzz*zqalv)* & max(zqsenh(jl,jk)-zqenh(jl,jk),0.) zhhatt(jl,jk)=zhhat end do do jl = 1, klon zhmin(jl) = 0. if( ldcum(jl).and.ktype(jl).eq.1 ) then ihmin(jl) = kcbot(jl) else ihmin(jl) = -1 end if end do zbi = 1./(25.*g) do jk = klev, 1, -1 do jl = 1, klon llo1 = ldcum(jl).and.ktype(jl).eq.1.and.ihmin(jl).eq.kcbot(jl) if (llo1.and.jk.lt.kcbot(jl).and.jk.ge.ictop0(jl)) then ikb = kcbot(jl) zro = rd*ztenh(jl,jk)/(g*paph(jl,jk)) zdz = (paph(jl,jk)-paph(jl,jk-1))*zro zdhdz=(cpd*(pten(jl,jk-1)-pten(jl,jk))+alv*(pqen(jl,jk-1)- & pqen(jl,jk))+(pgeo(jl,jk-1)-pgeo(jl,jk)))*g/(pgeo(jl, & jk-1)-pgeo(jl,jk)) zdepth = zgeoh(jl,jk) - zgeoh(jl,ikb) zfac = sqrt(1.+zdepth*zbi) zhmin(jl) = zhmin(jl) + zdhdz*zfac*zdz zrh = -alv*(zqsenh(jl,jk)-zqenh(jl,jk))*zfac if (zhmin(jl).gt.zrh) ihmin(jl) = jk end if end do end do do jl = 1, klon if (ldcum(jl).and.ktype(jl).eq.1) then if (ihmin(jl).lt.ictop0(jl)) ihmin(jl) = ictop0(jl) end if if(ktype(jl).eq.1) then zentr(jl)=entrpen else zentr(jl)=entrscv endif if(lndj(jl).eq.1) zentr(jl)=zentr(jl)*1.1 end do call cuasc_new & (klon, klev, klevp1, klevm1, ztenh, & zqenh, puen, pven, pten, pqen, & pqsen, pgeo, zgeoh, pap, paph, & pqte, pverv, ilwmin, ldcum, zhcbase, & ktype, ilab, ptu, pqu, plu, & zuu, zvu, pmfu, zmfub, zentr, & zmfus, zmfuq, zmful, plude, zdmfup, & kcbot, kctop, ictop0, icum, ztmst, & ihmin, zhhatt, zqsenh) if(icum.eq.0) return do jl=1,klon zpbmpt=paph(jl,kcbot(jl))-paph(jl,kctop(jl)) if(ldcum(jl)) ictop0(jl)=kctop(jl) if(ldcum(jl).and.ktype(jl).eq.1.and.zpbmpt.lt.zdnoprc) ktype(jl)=2 if(ktype(jl).eq.2) then zentr(jl)=entrscv if(lndj(jl).eq.1) zentr(jl)=zentr(jl)*1.1 endif zrfl(jl)=zdmfup(jl,1) end do do jk=2,klev do jl=1,klon zrfl(jl)=zrfl(jl)+zdmfup(jl,jk) end do end do if(lmfdd) then call cudlfs & (klon, klev, klevp1, ztenh, zqenh, & puen, pven, zgeoh, paph, ptu, & pqu, zuu, zvu, ldcum, kcbot, & kctop, zmfub, zrfl, ztd, zqd, & zud, zvd, pmfd, zmfds, zmfdq, & zdmfdp, idtop, loddraf) call cuddraf & (klon, klev, klevp1, ztenh, zqenh, & puen, pven, zgeoh, paph, zrfl, & loddraf, ztd, zqd, zud, zvd, & pmfd, zmfds, zmfdq, zdmfdp) end if do jl=1,klon zheat(jl)=0.0 zcape(jl)=0.0 zrelh(jl)=0.0 zmfub1(jl)=zmfub(jl) end do do jl=1,klon if(ldcum(jl).and.ktype(jl).eq.1) then ktop0=max(12,kctop(jl)) ikb = kcbot(jl) do jk=2,klev if(jk.le.kcbot(jl).and.jk.gt.kctop(jl)) then zro=paph(jl,jk)/(rd*ztenh(jl,jk)) zdz=(paph(jl,jk)-paph(jl,jk-1))/(g*zro) zheat(jl)=zheat(jl)+((pten(jl,jk-1)-pten(jl,jk) & +g*zdz/cpd)/ztenh(jl,jk)+0.608*(pqen(jl,jk-1)- & pqen(jl,jk)))*(pmfu(jl,jk)+pmfd(jl,jk))*g/zro zcape(jl)=zcape(jl)+g*((ptu(jl,jk)*(1.+.608*pqu(jl,jk) & -plu(jl,jk)))/(ztenh(jl,jk)*(1.+.608*zqenh(jl,jk))) & -1.0)*zdz endif if(jk.le.kcbot(jl).and.jk.gt.ktop0) then dept=(paph(jl,jk+1)-paph(jl,jk))/(paph(jl,ikb+1)- & paph(jl,ktop0+1)) zrelh(jl)=zrelh(jl)+dept*pqen(jl,jk)/pqsen(jl,jk) endif enddo if(zrelh(jl).ge.crirh) then zht=max(0.0,(zcape(jl)-0.0))/(ztau*zheat(jl)) zmfub1(jl)=max(zmfub(jl)*zht,0.01) zmfmax=(paph(jl,ikb)-paph(jl,ikb-1))*zcons2 zmfub1(jl)=min(zmfub1(jl),zmfmax) else zmfub1(jl)=0.01 zmfub(jl)=0.01 ldcum(jl)=.false. endif endif end do do jl=1,klon if(ktype(jl).ne.1) then ikb=kcbot(jl) if(pmfd(jl,ikb).lt.0.0.and.loddraf(jl)) then zeps=cmfdeps else zeps=0. endif zqumqe=pqu(jl,ikb)+plu(jl,ikb)- & zeps*zqd(jl,ikb)-(1.-zeps)*zqenh(jl,ikb) zdqmin=max(0.01*zqenh(jl,ikb),1.e-10) zmfmax=(paph(jl,ikb)-paph(jl,ikb-1))*zcons2 if(zdqpbl(jl).gt.0..and.zqumqe.gt.zdqmin.and.ldcum(jl) & .and.zmfub(jl).lt.zmfmax) then zmfub1(jl)=zdqpbl(jl)/(g*max(zqumqe,zdqmin)) else zmfub1(jl)=zmfub(jl) endif llo1=(ktype(jl).eq.2).and.abs(zmfub1(jl) & -zmfub(jl)).lt.0.2*zmfub(jl) if(.not.llo1) zmfub1(jl)=zmfub(jl) zmfub1(jl)=min(zmfub1(jl),zmfmax) end if end do do jk=1,klev do jl=1,klon if(ldcum(jl)) then zfac=zmfub1(jl)/max(zmfub(jl),1.e-10) pmfd(jl,jk)=pmfd(jl,jk)*zfac zmfds(jl,jk)=zmfds(jl,jk)*zfac zmfdq(jl,jk)=zmfdq(jl,jk)*zfac zdmfdp(jl,jk)=zdmfdp(jl,jk)*zfac else pmfd(jl,jk)=0.0 zmfds(jl,jk)=0.0 zmfdq(jl,jk)=0.0 zdmfdp(jl,jk)=0.0 endif end do end do do jl=1,klon if(ldcum(jl)) then zmfub(jl)=zmfub1(jl) else zmfub(jl)=0.0 endif end do call cuasc_new & (klon, klev, klevp1, klevm1, ztenh, & zqenh, puen, pven, pten, pqen, & pqsen, pgeo, zgeoh, pap, paph, & pqte, pverv, ilwmin, ldcum, zhcbase,& ktype, ilab, ptu, pqu, plu, & zuu, zvu, pmfu, zmfub, zentr, & zmfus, zmfuq, zmful, plude, zdmfup, & kcbot, kctop, ictop0, icum, ztmst, & ihmin, zhhatt, zqsenh) call cuflx & (klon, klev, klevp1, pqen, pqsen, & ztenh, zqenh, paph, zgeoh, kcbot, & kctop, idtop, ktype, loddraf, ldcum, & pmfu, pmfd, zmfus, zmfds, zmfuq, & zmfdq, zmful, plude, zdmfup, zdmfdp, & zrfl, prain, pten, zsfl, zdpmel, & itopm2, ztmst, sig1) call cudtdq & (klon, klev, klevp1, itopm2, paph, & ldcum, pten, ptte, pqte, zmfus, & zmfds, zmfuq, zmfdq, zmful, zdmfup, & zdmfdp, ztmst, zdpmel, prain, zrfl, & zsfl, psrain, psevap, psheat, psmelt, & prsfc, pssfc, paprc, paprsm, paprs, & pqen, pqsen, plude, pcte) if(lmfdudv) then call cududv & (klon, klev, klevp1, itopm2, ktype, & kcbot, paph, ldcum, puen, pven, & pvom, pvol, zuu, zud, zvu, & zvd, pmfu, pmfd, psdiss) end if return end subroutine cumastr_new subroutine cuini & (klon, klev, klevp1, klevm1, pten, & pqen, pqsen, puen, pven, pverv, & pgeo, paph, pgeoh, ptenh, pqenh, & pqsenh, klwmin, ptu, pqu, ptd, & pqd, puu, pvu, pud, pvd, & pmfu, pmfd, pmfus, pmfds, pmfuq, & pmfdq, pdmfup, pdmfdp, pdpmel, plu, & plude, klab) implicit none integer klon, klev, klevp1 integer klevm1 integer jk,jl,ik, icall real zdp, zzs real pten(klon,klev), pqen(klon,klev), & puen(klon,klev), pven(klon,klev), & pqsen(klon,klev), pverv(klon,klev), & pgeo(klon,klev), pgeoh(klon,klev), & paph(klon,klevp1), ptenh(klon,klev), & pqenh(klon,klev), pqsenh(klon,klev) real ptu(klon,klev), pqu(klon,klev), & ptd(klon,klev), pqd(klon,klev), & puu(klon,klev), pud(klon,klev), & pvu(klon,klev), pvd(klon,klev), & pmfu(klon,klev), pmfd(klon,klev), & pmfus(klon,klev), pmfds(klon,klev), & pmfuq(klon,klev), pmfdq(klon,klev), & pdmfup(klon,klev), pdmfdp(klon,klev), & plu(klon,klev), plude(klon,klev) real zwmax(klon), zph(klon), & pdpmel(klon,klev) integer klab(klon,klev), klwmin(klon) logical loflag(klon) zdp=0.5 do jk=2,klev do jl=1,klon pgeoh(jl,jk)=pgeo(jl,jk)+(pgeo(jl,jk-1)-pgeo(jl,jk))*zdp ptenh(jl,jk)=(max(cpd*pten(jl,jk-1)+pgeo(jl,jk-1), & cpd*pten(jl,jk)+pgeo(jl,jk))-pgeoh(jl,jk))*rcpd pqsenh(jl,jk)=pqsen(jl,jk-1) zph(jl)=paph(jl,jk) loflag(jl)=.true. end do ik=jk icall=0 call cuadjtq(klon,klev,ik,zph,ptenh,pqsenh,loflag,icall) do jl=1,klon pqenh(jl,jk)=min(pqen(jl,jk-1),pqsen(jl,jk-1)) & +(pqsenh(jl,jk)-pqsen(jl,jk-1)) pqenh(jl,jk)=max(pqenh(jl,jk),0.) end do end do do jl=1,klon ptenh(jl,klev)=(cpd*pten(jl,klev)+pgeo(jl,klev)- & pgeoh(jl,klev))*rcpd pqenh(jl,klev)=pqen(jl,klev) ptenh(jl,1)=pten(jl,1) pqenh(jl,1)=pqen(jl,1) pgeoh(jl,1)=pgeo(jl,1) klwmin(jl)=klev zwmax(jl)=0. end do do jk=klevm1,2,-1 do jl=1,klon zzs=max(cpd*ptenh(jl,jk)+pgeoh(jl,jk), & cpd*ptenh(jl,jk+1)+pgeoh(jl,jk+1)) ptenh(jl,jk)=(zzs-pgeoh(jl,jk))*rcpd end do end do do jk=klev,3,-1 do jl=1,klon if(pverv(jl,jk).lt.zwmax(jl)) then zwmax(jl)=pverv(jl,jk) klwmin(jl)=jk end if end do end do do jk=1,klev ik=jk-1 if(jk.eq.1) ik=1 do jl=1,klon ptu(jl,jk)=ptenh(jl,jk) ptd(jl,jk)=ptenh(jl,jk) pqu(jl,jk)=pqenh(jl,jk) pqd(jl,jk)=pqenh(jl,jk) plu(jl,jk)=0. puu(jl,jk)=puen(jl,ik) pud(jl,jk)=puen(jl,ik) pvu(jl,jk)=pven(jl,ik) pvd(jl,jk)=pven(jl,ik) pmfu(jl,jk)=0. pmfd(jl,jk)=0. pmfus(jl,jk)=0. pmfds(jl,jk)=0. pmfuq(jl,jk)=0. pmfdq(jl,jk)=0. pdmfup(jl,jk)=0. pdmfdp(jl,jk)=0. pdpmel(jl,jk)=0. plude(jl,jk)=0. klab(jl,jk)=0 end do end do return end subroutine cuini subroutine cubase & (klon, klev, klevp1, klevm1, ptenh, & pqenh, pgeoh, paph, ptu, pqu, & plu, puen, pven, puu, pvu, & ldcum, kcbot, klab) implicit none integer klon, klev, klevp1 integer klevm1 integer jl,jk,is,ik,icall,ikb real zbuo,zz real ptenh(klon,klev), pqenh(klon,klev), & pgeoh(klon,klev), paph(klon,klevp1) real ptu(klon,klev), pqu(klon,klev), & plu(klon,klev) real puen(klon,klev), pven(klon,klev), & puu(klon,klev), pvu(klon,klev) real zqold(klon,klev), zph(klon) integer klab(klon,klev), kcbot(klon) logical ldcum(klon), loflag(klon) logical ldbase(klon) logical llo1 do jl=1,klon klab(jl,klev)=1 kcbot(jl)=klevm1 ldcum(jl)=.false. ldbase(jl)=.false. puu(jl,klev)=puen(jl,klev)*(paph(jl,klevp1)-paph(jl,klev)) pvu(jl,klev)=pven(jl,klev)*(paph(jl,klevp1)-paph(jl,klev)) end do do jk=1,klev do jl=1,klon zqold(jl,jk)=0.0 end do end do do jk=klevm1,2,-1 is=0 do jl=1,klon if(klab(jl,jk+1).eq.1 .or.(ldcum(jl).and.kcbot(jl).eq.jk+1)) then is=is+1 loflag(jl)=.true. else loflag(jl)=.false. endif zph(jl)=paph(jl,jk) end do if(is.eq.0) cycle if(lmfdudv) then do jl=1,klon if(.not.ldbase(jl)) then puu(jl,klev)=puu(jl,klev)+ & puen(jl,jk)*(paph(jl,jk+1)-paph(jl,jk)) pvu(jl,klev)=pvu(jl,klev)+ & pven(jl,jk)*(paph(jl,jk+1)-paph(jl,jk)) endif enddo endif do jl=1,klon if(loflag(jl)) then pqu(jl,jk)=pqu(jl,jk+1) ptu(jl,jk)=(cpd*ptu(jl,jk+1)+pgeoh(jl,jk+1) & -pgeoh(jl,jk))*rcpd zqold(jl,jk)=pqu(jl,jk) end if end do ik=jk icall=1 call cuadjtq(klon,klev,ik,zph,ptu,pqu,loflag,icall) do jl=1,klon if(loflag(jl)) then if(pqu(jl,jk).eq.zqold(jl,jk)) then zbuo=ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk))- & ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))+zbuo0 if(zbuo.gt.0.) klab(jl,jk)=1 else klab(jl,jk)=2 plu(jl,jk)=plu(jl,jk)+zqold(jl,jk)-pqu(jl,jk) zbuo=ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk))- & ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))+zbuo0 llo1=zbuo.gt.0..and.klab(jl,jk+1).eq.1 if(llo1) then kcbot(jl)=jk ldcum(jl)=.true. ldbase(jl)=.true. end if end if end if end do end do if(lmfdudv) then do jl=1,klon if(ldcum(jl)) then ikb=kcbot(jl) zz=1./(paph(jl,klevp1)-paph(jl,ikb)) puu(jl,klev)=puu(jl,klev)*zz pvu(jl,klev)=pvu(jl,klev)*zz else puu(jl,klev)=puen(jl,klevm1) pvu(jl,klev)=pven(jl,klevm1) end if end do end if return end subroutine cubase subroutine cuasc_new & (klon, klev, klevp1, klevm1, ptenh, & pqenh, puen, pven, pten, pqen, & pqsen, pgeo, pgeoh, pap, paph, & pqte, pverv, klwmin, ldcum, phcbase,& ktype, klab, ptu, pqu, plu, & puu, pvu, pmfu, pmfub, pentr, & pmfus, pmfuq, pmful, plude, pdmfup, & kcbot, kctop, kctop0, kcum, ztmst, & khmin, phhatt, pqsenh) implicit none integer klon, klev, klevp1 integer klevm1,kcum real ztmst,zcons2,zdz,zdrodz integer jl,jk,ikb,ik,is,ikt,icall real zmfmax,zfac,zmftest,zdprho,zmse,znevn,zodmax real zqeen,zseen,zscde,zga,zdt,zscod real zqude,zqcod, zmfusk, zmfuqk,zmfulk real zbuo, zprcon, zlnew, zz, zdmfeu, zdmfdu real zbuoyz,zzdmf real ptenh(klon,klev), pqenh(klon,klev), & puen(klon,klev), pven(klon,klev), & pten(klon,klev), pqen(klon,klev), & pgeo(klon,klev), pgeoh(klon,klev), & pap(klon,klev), paph(klon,klevp1), & pqsen(klon,klev), pqte(klon,klev), & pverv(klon,klev), pqsenh(klon,klev) real ptu(klon,klev), pqu(klon,klev), & puu(klon,klev), pvu(klon,klev), & pmfu(klon,klev), zph(klon), & pmfub(klon), pentr(klon), & pmfus(klon,klev), pmfuq(klon,klev), & plu(klon,klev), plude(klon,klev), & pmful(klon,klev), pdmfup(klon,klev) real zdmfen(klon), zdmfde(klon), & zmfuu(klon), zmfuv(klon), & zpbase(klon), zqold(klon), & phhatt(klon,klev), zodetr(klon,klev), & zoentr(klon,klev), zbuoy(klon) real phcbase(klon) integer klwmin(klon), ktype(klon), & klab(klon,klev), kcbot(klon), & kctop(klon), kctop0(klon), & khmin(klon) logical ldcum(klon), loflag(klon) zcons2=1./(g*ztmst) do jl=1,klon zmfuu(jl)=0. zmfuv(jl)=0. zbuoy(jl)=0. if(.not.ldcum(jl)) ktype(jl)=0 end do do jk=1,klev do jl=1,klon plu(jl,jk)=0. pmfu(jl,jk)=0. pmfus(jl,jk)=0. pmfuq(jl,jk)=0. pmful(jl,jk)=0. plude(jl,jk)=0. pdmfup(jl,jk)=0. zoentr(jl,jk)=0. zodetr(jl,jk)=0. if(.not.ldcum(jl).or.ktype(jl).eq.3) klab(jl,jk)=0 if(.not.ldcum(jl).and.paph(jl,jk).lt.4.e4) kctop0(jl)=jk end do end do do jl=1,klon kctop(jl)=klevm1 if(.not.ldcum(jl)) then kcbot(jl)=klevm1 pmfub(jl)=0. pqu(jl,klev)=0. end if pmfu(jl,klev)=pmfub(jl) pmfus(jl,klev)=pmfub(jl)*(cpd*ptu(jl,klev)+pgeoh(jl,klev)) pmfuq(jl,klev)=pmfub(jl)*pqu(jl,klev) if(lmfdudv) then zmfuu(jl)=pmfub(jl)*puu(jl,klev) zmfuv(jl)=pmfub(jl)*pvu(jl,klev) end if end do do jl=1,klon ldcum(jl)=.false. if (ktype(jl).eq.1) then ikb = kcbot(jl) zbuoy(jl)=g*((ptu(jl,ikb)-ptenh(jl,ikb))/ptenh(jl,ikb)+ & 0.608*(pqu(jl,ikb)-pqenh(jl,ikb))) if (zbuoy(jl).gt.0.) then zdz = (pgeo(jl,ikb-1)-pgeo(jl,ikb))*zrg zdrodz = -log(pten(jl,ikb-1)/pten(jl,ikb))/zdz - & g/(rd*ptenh(jl,ikb)) zoentr(jl,ikb-1)=zbuoy(jl)*0.5/(1.+zbuoy(jl)*zdz) & +zdrodz zoentr(jl,ikb-1) = min(zoentr(jl,ikb-1),1.e-3) zoentr(jl,ikb-1) = max(zoentr(jl,ikb-1),0.) end if end if end do do jk=klevm1,2,-1 ik=jk if(lmfmid.and.ik.lt.klevm1.and.ik.gt.klev-13) then call cubasmc & (klon, klev, klevm1, ik, pten, & pqen, pqsen, puen, pven, pverv, & pgeo, pgeoh, ldcum, ktype, klab, & pmfu, pmfub, pentr, kcbot, ptu, & pqu, plu, puu, pvu, pmfus, & pmfuq, pmful, pdmfup, zmfuu, zmfuv) endif is=0 do jl=1,klon zqold(jl)=0.0 is=is+klab(jl,jk+1) if(klab(jl,jk+1).eq.0) klab(jl,jk)=0 loflag(jl)=klab(jl,jk+1).gt.0 zph(jl)=paph(jl,jk) if(ktype(jl).eq.3.and.jk.eq.kcbot(jl)) then zmfmax=(paph(jl,jk)-paph(jl,jk-1))*zcons2 if(pmfub(jl).gt.zmfmax) then zfac=zmfmax/pmfub(jl) pmfu(jl,jk+1)=pmfu(jl,jk+1)*zfac pmfus(jl,jk+1)=pmfus(jl,jk+1)*zfac pmfuq(jl,jk+1)=pmfuq(jl,jk+1)*zfac zmfuu(jl)=zmfuu(jl)*zfac zmfuv(jl)=zmfuv(jl)*zfac pmfub(jl)=zmfmax end if end if end do if(is.eq.0) cycle ik=jk call cuentr_new & (klon, klev, klevp1, ik, ptenh,& paph, pap, pgeoh, klwmin, ldcum,& ktype, kcbot, kctop0, zpbase, pmfu, & pentr, zdmfen, zdmfde, zodetr, khmin) do jl=1,klon if(loflag(jl)) then if(jk.lt.kcbot(jl)) then zmftest=pmfu(jl,jk+1)+zdmfen(jl)-zdmfde(jl) zmfmax=min(zmftest,(paph(jl,jk)-paph(jl,jk-1))*zcons2) zdmfen(jl)=max(zdmfen(jl)-max(zmftest-zmfmax,0.),0.) end if zdmfde(jl)=min(zdmfde(jl),0.75*pmfu(jl,jk+1)) pmfu(jl,jk)=pmfu(jl,jk+1)+zdmfen(jl)-zdmfde(jl) if (jk.lt.kcbot(jl)) then zdprho = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg zoentr(jl,jk) = zoentr(jl,jk)*zdprho*pmfu(jl,jk+1) zmftest = pmfu(jl,jk) + zoentr(jl,jk)-zodetr(jl,jk) zmfmax = min(zmftest,(paph(jl,jk)-paph(jl,jk-1))*zcons2) zoentr(jl,jk) = max(zoentr(jl,jk)-max(zmftest-zmfmax,0.),0.) end if if (ktype(jl).eq.1.and.jk.lt.kcbot(jl).and.jk.le.khmin(jl)) then zmse = cpd*ptu(jl,jk+1) + alv*pqu(jl,jk+1) + pgeoh(jl,jk+1) ikt = kctop0(jl) znevn=(pgeoh(jl,ikt)-pgeoh(jl,jk+1))*(zmse-phhatt(jl, & jk+1))*zrg if (znevn.le.0.) znevn = 1. zdprho = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg zodmax = ((phcbase(jl)-zmse)/znevn)*zdprho*pmfu(jl,jk+1) zodmax = max(zodmax,0.) zodetr(jl,jk) = min(zodetr(jl,jk),zodmax) end if zodetr(jl,jk) = min(zodetr(jl,jk),0.75*pmfu(jl,jk)) pmfu(jl,jk) = pmfu(jl,jk) + zoentr(jl,jk) - zodetr(jl,jk) zqeen=pqenh(jl,jk+1)*zdmfen(jl) zqeen=zqeen + pqenh(jl,jk+1)*zoentr(jl,jk) zseen=(cpd*ptenh(jl,jk+1)+pgeoh(jl,jk+1))*zdmfen(jl) zseen=zseen+(cpd*ptenh(jl,jk+1)+pgeoh(jl,jk+1))* & zoentr(jl,jk) zscde=(cpd*ptu(jl,jk+1)+pgeoh(jl,jk+1))*zdmfde(jl) zga = alv*pqsenh(jl,jk+1)/(rv*(ptenh(jl,jk+1)**2)) zdt = (plu(jl,jk+1)-0.608*(pqsenh(jl,jk+1)-pqenh(jl, & jk+1)))/(1./ptenh(jl,jk+1)+0.608*zga) zscod = cpd*ptenh(jl,jk+1) + pgeoh(jl,jk+1) + cpd*zdt zscde = zscde + zodetr(jl,jk)*zscod zqude = pqu(jl,jk+1)*zdmfde(jl) zqcod = pqsenh(jl,jk+1) + zga*zdt zqude = zqude + zodetr(jl,jk)*zqcod plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl) plude(jl,jk) = plude(jl,jk)+plu(jl,jk+1)*zodetr(jl,jk) zmfusk = pmfus(jl,jk+1) + zseen - zscde zmfuqk = pmfuq(jl,jk+1) + zqeen - zqude zmfulk = pmful(jl,jk+1) - plude(jl,jk) plu(jl,jk) = zmfulk*(1./max(cmfcmin,pmfu(jl,jk))) pqu(jl,jk) = zmfuqk*(1./max(cmfcmin,pmfu(jl,jk))) ptu(jl,jk)=(zmfusk*(1./max(cmfcmin,pmfu(jl,jk)))- & pgeoh(jl,jk))*rcpd ptu(jl,jk) = max(100.,ptu(jl,jk)) ptu(jl,jk) = min(400.,ptu(jl,jk)) zqold(jl) = pqu(jl,jk) end if end do ik=jk icall=1 call cuadjtq(klon,klev,ik,zph,ptu,pqu,loflag,icall) do jl=1,klon if(loflag(jl).and.pqu(jl,jk).ne.zqold(jl)) then klab(jl,jk)=2 plu(jl,jk)=plu(jl,jk)+zqold(jl)-pqu(jl,jk) zbuo=ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk)-plu(jl,jk))- & ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk)) if(klab(jl,jk+1).eq.1) zbuo=zbuo+zbuo0 if(zbuo.gt.0..and.pmfu(jl,jk).gt.0.01*pmfub(jl).and. & jk.ge.kctop0(jl)) then kctop(jl)=jk ldcum(jl)=.true. if(zpbase(jl)-paph(jl,jk).ge.zdnoprc) then zprcon=cprcon else zprcon=0. endif zlnew=plu(jl,jk)/(1.+zprcon*(pgeoh(jl,jk)-pgeoh(jl,jk+1))) pdmfup(jl,jk)=max(0.,(plu(jl,jk)-zlnew)*pmfu(jl,jk)) plu(jl,jk)=zlnew else klab(jl,jk)=0 pmfu(jl,jk)=0. end if end if if(loflag(jl)) then pmful(jl,jk)=plu(jl,jk)*pmfu(jl,jk) pmfus(jl,jk)=(cpd*ptu(jl,jk)+pgeoh(jl,jk))*pmfu(jl,jk) pmfuq(jl,jk)=pqu(jl,jk)*pmfu(jl,jk) end if end do if(lmfdudv) then do jl=1,klon zdmfen(jl) = zdmfen(jl) + zoentr(jl,jk) zdmfde(jl) = zdmfde(jl) + zodetr(jl,jk) if(loflag(jl)) then if(ktype(jl).eq.1.or.ktype(jl).eq.3) then if(zdmfen(jl).le.1.e-20) then zz=3. else zz=2. endif else if(zdmfen(jl).le.1.0e-20) then zz=1. else zz=0. endif end if zdmfeu=zdmfen(jl)+zz*zdmfde(jl) zdmfdu=zdmfde(jl)+zz*zdmfde(jl) zdmfdu=min(zdmfdu,0.75*pmfu(jl,jk+1)) zmfuu(jl)=zmfuu(jl)+ & zdmfeu*puen(jl,jk)-zdmfdu*puu(jl,jk+1) zmfuv(jl)=zmfuv(jl)+ & zdmfeu*pven(jl,jk)-zdmfdu*pvu(jl,jk+1) if(pmfu(jl,jk).gt.0.) then puu(jl,jk)=zmfuu(jl)*(1./pmfu(jl,jk)) pvu(jl,jk)=zmfuv(jl)*(1./pmfu(jl,jk)) end if end if end do end if do jl = 1, klon if (loflag(jl).and.ktype(jl).eq.1) then zbuoyz=g*((ptu(jl,jk)-ptenh(jl,jk))/ptenh(jl,jk)+ & 0.608*(pqu(jl,jk)-pqenh(jl,jk))-plu(jl,jk)) zbuoyz = max(zbuoyz,0.0) zdz = (pgeo(jl,jk-1)-pgeo(jl,jk))*zrg zdrodz = -log(pten(jl,jk-1)/pten(jl,jk))/zdz - & g/(rd*ptenh(jl,jk)) zbuoy(jl) = zbuoy(jl) + zbuoyz*zdz zoentr(jl,jk-1) = zbuoyz*0.5/(1.+zbuoy(jl))+zdrodz zoentr(jl,jk-1) = min(zoentr(jl,jk-1),1.e-3) zoentr(jl,jk-1) = max(zoentr(jl,jk-1),0.) end if end do end do do jl=1,klon if(kctop(jl).eq.klevm1) ldcum(jl)=.false. kcbot(jl)=max(kcbot(jl),kctop(jl)) end do is=0 do jl=1,klon if(ldcum(jl)) then is=is+1 endif end do kcum=is if(is.eq.0) return do jl=1,klon if(ldcum(jl)) then jk=kctop(jl)-1 zzdmf=cmfctop zdmfde(jl)=(1.-zzdmf)*pmfu(jl,jk+1) plude(jl,jk)=zdmfde(jl)*plu(jl,jk+1) pmfu(jl,jk)=pmfu(jl,jk+1)-zdmfde(jl) pmfus(jl,jk)=(cpd*ptu(jl,jk)+pgeoh(jl,jk))*pmfu(jl,jk) pmfuq(jl,jk)=pqu(jl,jk)*pmfu(jl,jk) pmful(jl,jk)=plu(jl,jk)*pmfu(jl,jk) plude(jl,jk-1)=pmful(jl,jk) pdmfup(jl,jk)=0. end if end do if(lmfdudv) then do jl=1,klon if(ldcum(jl)) then jk=kctop(jl)-1 puu(jl,jk)=puu(jl,jk+1) pvu(jl,jk)=pvu(jl,jk+1) end if end do end if return end subroutine cuasc_new subroutine cudlfs & (klon, klev, klevp1, ptenh, pqenh, & puen, pven, pgeoh, paph, ptu, & pqu, puu, pvu, ldcum, kcbot, & kctop, pmfub, prfl, ptd, pqd, & pud, pvd, pmfd, pmfds, pmfdq, & pdmfdp, kdtop, lddraf) implicit none integer klon, klev, klevp1 integer jl,ke,jk,is,ik,icall real zttest, zqtest, zbuo, zmftop real ptenh(klon,klev), pqenh(klon,klev), & puen(klon,klev), pven(klon,klev), & pgeoh(klon,klev), paph(klon,klevp1), & ptu(klon,klev), pqu(klon,klev), & puu(klon,klev), pvu(klon,klev), & pmfub(klon), prfl(klon) real ptd(klon,klev), pqd(klon,klev), & pud(klon,klev), pvd(klon,klev), & pmfd(klon,klev), pmfds(klon,klev), & pmfdq(klon,klev), pdmfdp(klon,klev) real ztenwb(klon,klev), zqenwb(klon,klev), & zcond(klon), zph(klon) integer kcbot(klon), kctop(klon), & kdtop(klon) logical ldcum(klon), llo2(klon), & lddraf(klon) do jl=1,klon lddraf(jl)=.false. kdtop(jl)=klevp1 end do if(.not.lmfdd) return ke=klev-3 do jk=3,ke is=0 do jl=1,klon ztenwb(jl,jk)=ptenh(jl,jk) zqenwb(jl,jk)=pqenh(jl,jk) zph(jl)=paph(jl,jk) llo2(jl)=ldcum(jl).and.prfl(jl).gt.0..and..not.lddraf(jl).and. & (jk.lt.kcbot(jl).and.jk.gt.kctop(jl)) if(llo2(jl))then is=is+1 endif end do if(is.eq.0) cycle ik=jk icall=2 call cuadjtq(klon,klev,ik,zph,ztenwb,zqenwb,llo2,icall) do jl=1,klon if(llo2(jl)) then zttest=0.5*(ptu(jl,jk)+ztenwb(jl,jk)) zqtest=0.5*(pqu(jl,jk)+zqenwb(jl,jk)) zbuo=zttest*(1.+vtmpc1*zqtest)- & ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk)) zcond(jl)=pqenh(jl,jk)-zqenwb(jl,jk) zmftop=-cmfdeps*pmfub(jl) if(zbuo.lt.0..and.prfl(jl).gt.10.*zmftop*zcond(jl)) then kdtop(jl)=jk lddraf(jl)=.true. ptd(jl,jk)=zttest pqd(jl,jk)=zqtest pmfd(jl,jk)=zmftop pmfds(jl,jk)=pmfd(jl,jk)*(cpd*ptd(jl,jk)+pgeoh(jl,jk)) pmfdq(jl,jk)=pmfd(jl,jk)*pqd(jl,jk) pdmfdp(jl,jk-1)=-0.5*pmfd(jl,jk)*zcond(jl) prfl(jl)=prfl(jl)+pdmfdp(jl,jk-1) end if end if end do if(lmfdudv) then do jl=1,klon if(pmfd(jl,jk).lt.0.) then pud(jl,jk)=0.5*(puu(jl,jk)+puen(jl,jk-1)) pvd(jl,jk)=0.5*(pvu(jl,jk)+pven(jl,jk-1)) end if end do end if end do return end subroutine cudlfs subroutine cuddraf & (klon, klev, klevp1, ptenh, pqenh, & puen, pven, pgeoh, paph, prfl, & lddraf, ptd, pqd, pud, pvd, & pmfd, pmfds, pmfdq, pdmfdp) implicit none integer klon, klev, klevp1 integer jk,is,jl,itopde, ik, icall real zentr,zseen, zqeen, zsdde, zqdde,zmfdsk, zmfdqk real zbuo, zdmfdp, zmfduk, zmfdvk real ptenh(klon,klev), pqenh(klon,klev), & puen(klon,klev), pven(klon,klev), & pgeoh(klon,klev), paph(klon,klevp1) real ptd(klon,klev), pqd(klon,klev), & pud(klon,klev), pvd(klon,klev), & pmfd(klon,klev), pmfds(klon,klev), & pmfdq(klon,klev), pdmfdp(klon,klev), & prfl(klon) real zdmfen(klon), zdmfde(klon), & zcond(klon), zph(klon) logical lddraf(klon), llo2(klon) do jk=3,klev is=0 do jl=1,klon zph(jl)=paph(jl,jk) llo2(jl)=lddraf(jl).and.pmfd(jl,jk-1).lt.0. if(llo2(jl)) then is=is+1 endif end do if(is.eq.0) cycle do jl=1,klon if(llo2(jl)) then zentr=entrdd*pmfd(jl,jk-1)*rd*ptenh(jl,jk-1)/ & (g*paph(jl,jk-1))*(paph(jl,jk)-paph(jl,jk-1)) zdmfen(jl)=zentr zdmfde(jl)=zentr end if end do itopde=klev-2 if(jk.gt.itopde) then do jl=1,klon if(llo2(jl)) then zdmfen(jl)=0. zdmfde(jl)=pmfd(jl,itopde)* & (paph(jl,jk)-paph(jl,jk-1))/ & (paph(jl,klevp1)-paph(jl,itopde)) end if end do end if do jl=1,klon if(llo2(jl)) then pmfd(jl,jk)=pmfd(jl,jk-1)+zdmfen(jl)-zdmfde(jl) zseen=(cpd*ptenh(jl,jk-1)+pgeoh(jl,jk-1))*zdmfen(jl) zqeen=pqenh(jl,jk-1)*zdmfen(jl) zsdde=(cpd*ptd(jl,jk-1)+pgeoh(jl,jk-1))*zdmfde(jl) zqdde=pqd(jl,jk-1)*zdmfde(jl) zmfdsk=pmfds(jl,jk-1)+zseen-zsdde zmfdqk=pmfdq(jl,jk-1)+zqeen-zqdde pqd(jl,jk)=zmfdqk*(1./min(-cmfcmin,pmfd(jl,jk))) ptd(jl,jk)=(zmfdsk*(1./min(-cmfcmin,pmfd(jl,jk)))- & pgeoh(jl,jk))*rcpd ptd(jl,jk)=min(400.,ptd(jl,jk)) ptd(jl,jk)=max(100.,ptd(jl,jk)) zcond(jl)=pqd(jl,jk) end if end do ik=jk icall=2 call cuadjtq(klon,klev,ik,zph,ptd,pqd,llo2,icall) do jl=1,klon if(llo2(jl)) then zcond(jl)=zcond(jl)-pqd(jl,jk) zbuo=ptd(jl,jk)*(1.+vtmpc1*pqd(jl,jk))- & ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk)) if(zbuo.ge.0..or.prfl(jl).le.(pmfd(jl,jk)*zcond(jl))) then pmfd(jl,jk)=0. endif pmfds(jl,jk)=(cpd*ptd(jl,jk)+pgeoh(jl,jk))*pmfd(jl,jk) pmfdq(jl,jk)=pqd(jl,jk)*pmfd(jl,jk) zdmfdp=-pmfd(jl,jk)*zcond(jl) pdmfdp(jl,jk-1)=zdmfdp prfl(jl)=prfl(jl)+zdmfdp end if end do if(lmfdudv) then do jl=1,klon if(llo2(jl).and.pmfd(jl,jk).lt.0.) then zmfduk=pmfd(jl,jk-1)*pud(jl,jk-1)+ & zdmfen(jl)*puen(jl,jk-1)-zdmfde(jl)*pud(jl,jk-1) zmfdvk=pmfd(jl,jk-1)*pvd(jl,jk-1)+ & zdmfen(jl)*pven(jl,jk-1)-zdmfde(jl)*pvd(jl,jk-1) pud(jl,jk)=zmfduk*(1./min(-cmfcmin,pmfd(jl,jk))) pvd(jl,jk)=zmfdvk*(1./min(-cmfcmin,pmfd(jl,jk))) end if end do end if end do return end subroutine cuddraf subroutine cuflx & (klon, klev, klevp1, pqen, pqsen, & ptenh, pqenh, paph, pgeoh, kcbot, & kctop, kdtop, ktype, lddraf, ldcum, & pmfu, pmfd, pmfus, pmfds, pmfuq, & pmfdq, pmful, plude, pdmfup, pdmfdp, & prfl, prain, pten, psfl, pdpmel, & ktopm2, ztmst, sig1) implicit none integer klon, klev, klevp1 integer ktopm2, itop, jl, jk, ikb real ztmst, zcons1, zcons2, zcucov, ztmelp2 real zzp, zfac, zsnmlt, zrfl, cevapcu, zrnew real zrmin, zrfln, zdrfl, zdpevap real pqen(klon,klev), pqsen(klon,klev), & ptenh(klon,klev), pqenh(klon,klev), & paph(klon,klevp1), pgeoh(klon,klev) real pmfu(klon,klev), pmfd(klon,klev), & pmfus(klon,klev), pmfds(klon,klev), & pmfuq(klon,klev), pmfdq(klon,klev), & pdmfup(klon,klev), pdmfdp(klon,klev), & pmful(klon,klev), plude(klon,klev), & prfl(klon), prain(klon) real pten(klon,klev), pdpmel(klon,klev), & psfl(klon), zpsubcl(klon) real sig1(klev) integer kcbot(klon), kctop(klon), & kdtop(klon), ktype(klon) logical lddraf(klon), ldcum(klon) zcons1=cpd/(alf*g*ztmst) zcons2=1./(g*ztmst) zcucov=0.05 ztmelp2=tmelt+2. itop=klev do jl=1,klon prfl(jl)=0. psfl(jl)=0. prain(jl)=0. if(.not.lmfscv.and.ktype(jl).eq.2)then ldcum(jl)=.false. lddraf(jl)=.false. endif itop=min(itop,kctop(jl)) if(.not.ldcum(jl).or.kdtop(jl).lt.kctop(jl)) lddraf(jl)=.false. if(.not.ldcum(jl)) ktype(jl)=0 end do ktopm2=itop-2 do jk=ktopm2,klev do jl=1,klon if(ldcum(jl).and.jk.ge.kctop(jl)-1) then pmfus(jl,jk)=pmfus(jl,jk)-pmfu(jl,jk)* & (cpd*ptenh(jl,jk)+pgeoh(jl,jk)) pmfuq(jl,jk)=pmfuq(jl,jk)-pmfu(jl,jk)*pqenh(jl,jk) if(lddraf(jl).and.jk.ge.kdtop(jl)) then pmfds(jl,jk)=pmfds(jl,jk)-pmfd(jl,jk)* & (cpd*ptenh(jl,jk)+pgeoh(jl,jk)) pmfdq(jl,jk)=pmfdq(jl,jk)-pmfd(jl,jk)*pqenh(jl,jk) else pmfd(jl,jk)=0. pmfds(jl,jk)=0. pmfdq(jl,jk)=0. pdmfdp(jl,jk-1)=0. end if else pmfu(jl,jk)=0. pmfd(jl,jk)=0. pmfus(jl,jk)=0. pmfds(jl,jk)=0. pmfuq(jl,jk)=0. pmfdq(jl,jk)=0. pmful(jl,jk)=0. pdmfup(jl,jk-1)=0. pdmfdp(jl,jk-1)=0. plude(jl,jk-1)=0. end if end do end do do jk=ktopm2,klev do jl=1,klon if(ldcum(jl).and.jk.gt.kcbot(jl)) then ikb=kcbot(jl) zzp=((paph(jl,klevp1)-paph(jl,jk))/ & (paph(jl,klevp1)-paph(jl,ikb))) if(ktype(jl).eq.3) then zzp=zzp**2 endif pmfu(jl,jk)=pmfu(jl,ikb)*zzp pmfus(jl,jk)=pmfus(jl,ikb)*zzp pmfuq(jl,jk)=pmfuq(jl,ikb)*zzp pmful(jl,jk)=pmful(jl,ikb)*zzp end if if(ldcum(jl)) then prain(jl)=prain(jl)+pdmfup(jl,jk) if(pten(jl,jk).gt.tmelt) then prfl(jl)=prfl(jl)+pdmfup(jl,jk)+pdmfdp(jl,jk) if(psfl(jl).gt.0..and.pten(jl,jk).gt.ztmelp2) then zfac=zcons1*(paph(jl,jk+1)-paph(jl,jk)) zsnmlt=min(psfl(jl),zfac*(pten(jl,jk)-ztmelp2)) pdpmel(jl,jk)=zsnmlt psfl(jl)=psfl(jl)-zsnmlt prfl(jl)=prfl(jl)+zsnmlt end if else psfl(jl)=psfl(jl)+pdmfup(jl,jk)+pdmfdp(jl,jk) end if end if end do end do do jl=1,klon prfl(jl)=max(prfl(jl),0.) psfl(jl)=max(psfl(jl),0.) zpsubcl(jl)=prfl(jl)+psfl(jl) end do do jk=ktopm2,klev do jl=1,klon if(ldcum(jl).and.jk.ge.kcbot(jl).and. & zpsubcl(jl).gt.1.e-20) then zrfl=zpsubcl(jl) cevapcu=cevapcu1*sqrt(cevapcu2*sqrt(sig1(jk))) zrnew=(max(0.,sqrt(zrfl/zcucov)- & cevapcu*(paph(jl,jk+1)-paph(jl,jk))* & max(0.,pqsen(jl,jk)-pqen(jl,jk))))**2*zcucov zrmin=zrfl-zcucov*max(0.,0.8*pqsen(jl,jk)-pqen(jl,jk)) & *zcons2*(paph(jl,jk+1)-paph(jl,jk)) zrnew=max(zrnew,zrmin) zrfln=max(zrnew,0.) zdrfl=min(0.,zrfln-zrfl) pdmfup(jl,jk)=pdmfup(jl,jk)+zdrfl zpsubcl(jl)=zrfln end if end do end do do jl=1,klon zdpevap=zpsubcl(jl)-(prfl(jl)+psfl(jl)) prfl(jl)=prfl(jl)+zdpevap*prfl(jl)* & (1./max(1.e-20,prfl(jl)+psfl(jl))) psfl(jl)=psfl(jl)+zdpevap*psfl(jl)* & (1./max(1.e-20,prfl(jl)+psfl(jl))) end do return end subroutine cuflx subroutine cudtdq & (klon, klev, klevp1, ktopm2, paph, & ldcum, pten, ptte, pqte, pmfus, & pmfds, pmfuq, pmfdq, pmful, pdmfup, & pdmfdp, ztmst, pdpmel, prain, prfl, & psfl, psrain, psevap, psheat, psmelt, & prsfc, pssfc, paprc, paprsm, paprs, & pqen, pqsen, plude, pcte) implicit none integer klon, klev, klevp1 integer ktopm2,jl, jk real ztmst, psrain, psevap, psheat, psmelt, zdiagt, zdiagw real zalv, rhk, rhcoe, pldfd, zdtdt, zdqdt real ptte(klon,klev), pqte(klon,klev), & pten(klon,klev), plude(klon,klev), & pgeo(klon,klev), paph(klon,klevp1), & paprc(klon), paprs(klon), & paprsm(klon), pcte(klon,klev), & prsfc(klon), pssfc(klon) real pmfus(klon,klev), pmfds(klon,klev), & pmfuq(klon,klev), pmfdq(klon,klev), & pmful(klon,klev), pqsen(klon,klev), & pdmfup(klon,klev), pdmfdp(klon,klev),& prfl(klon), prain(klon), & pqen(klon,klev) real pdpmel(klon,klev), psfl(klon) real zsheat(klon), zmelt(klon) logical ldcum(klon) zdiagt=ztmst zdiagw=zdiagt/rhoh2o do jl=1,klon zmelt(jl)=0. zsheat(jl)=0. end do do jk=ktopm2,klev if(jk.lt.klev) then do jl=1,klon if(ldcum(jl)) then if(pten(jl,jk).gt.tmelt) then zalv=alv else zalv=als endif rhk=min(1.0,pqen(jl,jk)/pqsen(jl,jk)) rhcoe=max(0.0,(rhk-rhc)/(rhm-rhc)) pldfd=max(0.0,rhcoe*fdbk*plude(jl,jk)) zdtdt=(g/(paph(jl,jk+1)-paph(jl,jk)))*rcpd* & (pmfus(jl,jk+1)-pmfus(jl,jk)+ & pmfds(jl,jk+1)-pmfds(jl,jk)-alf*pdpmel(jl,jk) & -zalv*(pmful(jl,jk+1)-pmful(jl,jk)-pldfd- & (pdmfup(jl,jk)+pdmfdp(jl,jk)))) ptte(jl,jk)=ptte(jl,jk)+zdtdt zdqdt=(g/(paph(jl,jk+1)-paph(jl,jk)))*& (pmfuq(jl,jk+1)-pmfuq(jl,jk)+ & pmfdq(jl,jk+1)-pmfdq(jl,jk)+ & pmful(jl,jk+1)-pmful(jl,jk)-pldfd- & (pdmfup(jl,jk)+pdmfdp(jl,jk))) pqte(jl,jk)=pqte(jl,jk)+zdqdt pcte(jl,jk)=(g/(paph(jl,jk+1)-paph(jl,jk)))*pldfd zsheat(jl)=zsheat(jl)+zalv*(pdmfup(jl,jk)+pdmfdp(jl,jk)) zmelt(jl)=zmelt(jl)+pdpmel(jl,jk) end if end do else do jl=1,klon if(ldcum(jl)) then if(pten(jl,jk).gt.tmelt) then zalv=alv else zalv=als endif rhk=min(1.0,pqen(jl,jk)/pqsen(jl,jk)) rhcoe=max(0.0,(rhk-rhc)/(rhm-rhc)) pldfd=max(0.0,rhcoe*fdbk*plude(jl,jk)) zdtdt=-(g/(paph(jl,jk+1)-paph(jl,jk)))*rcpd* & (pmfus(jl,jk)+pmfds(jl,jk)+alf*pdpmel(jl,jk)-zalv* & (pmful(jl,jk)+pdmfup(jl,jk)+pdmfdp(jl,jk)+pldfd)) ptte(jl,jk)=ptte(jl,jk)+zdtdt zdqdt=-(g/(paph(jl,jk+1)-paph(jl,jk)))* & (pmfuq(jl,jk)+pmfdq(jl,jk)+pldfd+ & (pmful(jl,jk)+pdmfup(jl,jk)+pdmfdp(jl,jk))) pqte(jl,jk)=pqte(jl,jk)+zdqdt pcte(jl,jk)=(g/(paph(jl,jk+1)-paph(jl,jk)))*pldfd zsheat(jl)=zsheat(jl)+zalv*(pdmfup(jl,jk)+pdmfdp(jl,jk)) zmelt(jl)=zmelt(jl)+pdpmel(jl,jk) end if end do end if end do do jl=1,klon prsfc(jl)=prfl(jl) pssfc(jl)=psfl(jl) paprc(jl)=paprc(jl)+zdiagw*(prfl(jl)+psfl(jl)) paprs(jl)=paprsm(jl)+zdiagw*psfl(jl) psheat=psheat+zsheat(jl) psrain=psrain+prain(jl) psevap=psevap-(prfl(jl)+psfl(jl)) psmelt=psmelt+zmelt(jl) end do psevap=psevap+psrain return end subroutine cudtdq subroutine cududv & (klon, klev, klevp1, ktopm2, ktype, & kcbot, paph, ldcum, puen, pven, & pvom, pvol, puu, pud, pvu, & pvd, pmfu, pmfd, psdiss) implicit none integer klon, klev, klevp1 integer ktopm2, jk, ik, jl, ikb real psdiss,zzp, zdudt ,zdvdt, zsum real puen(klon,klev), pven(klon,klev), & pvol(klon,klev), pvom(klon,klev), & paph(klon,klevp1) real puu(klon,klev), pud(klon,klev), & pvu(klon,klev), pvd(klon,klev), & pmfu(klon,klev), pmfd(klon,klev) real zmfuu(klon,klev), zmfdu(klon,klev), & zmfuv(klon,klev), zmfdv(klon,klev), & zdiss(klon) integer ktype(klon), kcbot(klon) logical ldcum(klon) do jk=ktopm2,klev ik=jk-1 do jl=1,klon if(ldcum(jl)) then zmfuu(jl,jk)=pmfu(jl,jk)*(puu(jl,jk)-puen(jl,ik)) zmfuv(jl,jk)=pmfu(jl,jk)*(pvu(jl,jk)-pven(jl,ik)) zmfdu(jl,jk)=pmfd(jl,jk)*(pud(jl,jk)-puen(jl,ik)) zmfdv(jl,jk)=pmfd(jl,jk)*(pvd(jl,jk)-pven(jl,ik)) end if end do end do do jk=ktopm2,klev do jl=1,klon if(ldcum(jl).and.jk.gt.kcbot(jl)) then ikb=kcbot(jl) zzp=((paph(jl,klevp1)-paph(jl,jk))/ & (paph(jl,klevp1)-paph(jl,ikb))) if(ktype(jl).eq.3) then zzp=zzp**2 endif zmfuu(jl,jk)=zmfuu(jl,ikb)*zzp zmfuv(jl,jk)=zmfuv(jl,ikb)*zzp zmfdu(jl,jk)=zmfdu(jl,ikb)*zzp zmfdv(jl,jk)=zmfdv(jl,ikb)*zzp end if end do end do do jl=1,klon zdiss(jl)=0. end do do jk=ktopm2,klev if(jk.lt.klev) then do jl=1,klon if(ldcum(jl)) then zdudt=(g/(paph(jl,jk+1)-paph(jl,jk)))* & (zmfuu(jl,jk+1)-zmfuu(jl,jk)+ & zmfdu(jl,jk+1)-zmfdu(jl,jk)) zdvdt=(g/(paph(jl,jk+1)-paph(jl,jk)))* & (zmfuv(jl,jk+1)-zmfuv(jl,jk)+ & zmfdv(jl,jk+1)-zmfdv(jl,jk)) zdiss(jl)=zdiss(jl)+ & puen(jl,jk)*(zmfuu(jl,jk+1)-zmfuu(jl,jk)+ & zmfdu(jl,jk+1)-zmfdu(jl,jk))+ & pven(jl,jk)*(zmfuv(jl,jk+1)-zmfuv(jl,jk)+ & zmfdv(jl,jk+1)-zmfdv(jl,jk)) pvom(jl,jk)=pvom(jl,jk)+zdudt pvol(jl,jk)=pvol(jl,jk)+zdvdt end if end do else do jl=1,klon if(ldcum(jl)) then zdudt=-(g/(paph(jl,jk+1)-paph(jl,jk)))* & (zmfuu(jl,jk)+zmfdu(jl,jk)) zdvdt=-(g/(paph(jl,jk+1)-paph(jl,jk)))* & (zmfuv(jl,jk)+zmfdv(jl,jk)) zdiss(jl)=zdiss(jl)- & (puen(jl,jk)*(zmfuu(jl,jk)+zmfdu(jl,jk))+ & pven(jl,jk)*(zmfuv(jl,jk)+zmfdv(jl,jk))) pvom(jl,jk)=pvom(jl,jk)+zdudt pvol(jl,jk)=pvol(jl,jk)+zdvdt end if end do end if end do zsum=ssum(klon,zdiss(1),1) psdiss=psdiss+zsum return end subroutine cududv subroutine cubasmc & (klon, klev, klevm1, kk, pten, & pqen, pqsen, puen, pven, pverv, & pgeo, pgeoh, ldcum, ktype, klab, & pmfu, pmfub, pentr, kcbot, ptu, & pqu, plu, puu, pvu, pmfus, & pmfuq, pmful, pdmfup, pmfuu, pmfuv) implicit none integer klon, klev, klevp1 integer klevm1,kk, jl real zzzmb real pten(klon,klev), pqen(klon,klev), & puen(klon,klev), pven(klon,klev), & pqsen(klon,klev), pverv(klon,klev), & pgeo(klon,klev), pgeoh(klon,klev) real ptu(klon,klev), pqu(klon,klev), & puu(klon,klev), pvu(klon,klev), & plu(klon,klev), pmfu(klon,klev), & pmfub(klon), pentr(klon), & pmfus(klon,klev), pmfuq(klon,klev), & pmful(klon,klev), pdmfup(klon,klev), & pmfuu(klon), pmfuv(klon) integer ktype(klon), kcbot(klon), & klab(klon,klev) logical ldcum(klon) do jl=1,klon if( .not. ldcum(jl).and.klab(jl,kk+1).eq.0.0.and. & pqen(jl,kk).gt.0.80*pqsen(jl,kk)) then ptu(jl,kk+1)=(cpd*pten(jl,kk)+pgeo(jl,kk)-pgeoh(jl,kk+1)) & *rcpd pqu(jl,kk+1)=pqen(jl,kk) plu(jl,kk+1)=0. zzzmb=max(cmfcmin,-pverv(jl,kk)/g) zzzmb=min(zzzmb,cmfcmax) pmfub(jl)=zzzmb pmfu(jl,kk+1)=pmfub(jl) pmfus(jl,kk+1)=pmfub(jl)*(cpd*ptu(jl,kk+1)+pgeoh(jl,kk+1)) pmfuq(jl,kk+1)=pmfub(jl)*pqu(jl,kk+1) pmful(jl,kk+1)=0. pdmfup(jl,kk+1)=0. kcbot(jl)=kk klab(jl,kk+1)=1 ktype(jl)=3 pentr(jl)=entrmid if(lmfdudv) then puu(jl,kk+1)=puen(jl,kk) pvu(jl,kk+1)=pven(jl,kk) pmfuu(jl)=pmfub(jl)*puu(jl,kk+1) pmfuv(jl)=pmfub(jl)*pvu(jl,kk+1) end if end if end do return end subroutine cubasmc subroutine cuadjtq(klon,klev,kk,pp,pt,pq,ldflag,kcall) implicit none integer klon, klev integer kk, kcall, isum, jl real zqsat, zcor, zcond1, tt real pt(klon,klev), pq(klon,klev), & zcond(klon), zqp(klon), & pp(klon) logical ldflag(klon) if (kcall.eq.1 ) then isum=0 do jl=1,klon zcond(jl)=0. if(ldflag(jl)) then zqp(jl)=1./pp(jl) tt=pt(jl,kk) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond(jl)=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) zcond(jl)=max(zcond(jl),0.) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond(jl) pq(jl,kk)=pq(jl,kk)-zcond(jl) if(zcond(jl).ne.0.0) isum=isum+1 end if end do if(isum.eq.0) return do jl=1,klon if(ldflag(jl).and.zcond(jl).ne.0.) then tt=pt(jl,kk) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond1=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond1 pq(jl,kk)=pq(jl,kk)-zcond1 end if end do end if if(kcall.eq.2) then isum=0 do jl=1,klon zcond(jl)=0. if(ldflag(jl)) then tt=pt(jl,kk) zqp(jl)=1./pp(jl) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond(jl)=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) zcond(jl)=min(zcond(jl),0.) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond(jl) pq(jl,kk)=pq(jl,kk)-zcond(jl) if(zcond(jl).ne.0.0) isum=isum+1 end if end do if(isum.eq.0) return do jl=1,klon if(ldflag(jl).and.zcond(jl).ne.0.) then tt=pt(jl,kk) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond1=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond1 pq(jl,kk)=pq(jl,kk)-zcond1 end if end do end if if(kcall.eq.0) then isum=0 do jl=1,klon tt=pt(jl,kk) zqp(jl)=1./pp(jl) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond(jl)=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond(jl) pq(jl,kk)=pq(jl,kk)-zcond(jl) if(zcond(jl).ne.0.0) isum=isum+1 end do if(isum.eq.0) return do jl=1,klon tt=pt(jl,kk) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond1=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond1 pq(jl,kk)=pq(jl,kk)-zcond1 end do end if if(kcall.eq.4) then do jl=1,klon tt=pt(jl,kk) zqp(jl)=1./pp(jl) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond(jl)=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond(jl) pq(jl,kk)=pq(jl,kk)-zcond(jl) end do do jl=1,klon tt=pt(jl,kk) zqsat=tlucua(tt)*zqp(jl) zqsat=min(0.5,zqsat) zcor=1./(1.-vtmpc1*zqsat) zqsat=zqsat*zcor zcond1=(pq(jl,kk)-zqsat)/(1.+zqsat*zcor*tlucub(tt)) pt(jl,kk)=pt(jl,kk)+tlucuc(tt)*zcond1 pq(jl,kk)=pq(jl,kk)-zcond1 end do end if return end subroutine cuadjtq subroutine cuentr_new & (klon, klev, klevp1, kk, ptenh, & paph, pap, pgeoh, klwmin, ldcum, & ktype, kcbot, kctop0, zpbase, pmfu, & pentr, zdmfen, zdmfde, zodetr, khmin) implicit none integer klon, klev, klevp1 integer kk, jl, iklwmin,ikb, ikt, ikh real zrrho, zdprho, zpmid, zentr, zzmzk, ztmzk, arg, zorgde real ptenh(klon,klev), & pap(klon,klev), paph(klon,klevp1), & pmfu(klon,klev), pgeoh(klon,klev), & pentr(klon), zpbase(klon), & zdmfen(klon), zdmfde(klon), & zodetr(klon,klev) integer klwmin(klon), ktype(klon), & kcbot(klon), kctop0(klon), & khmin(klon) logical ldcum(klon),llo1,llo2 do jl = 1, klon zpbase(jl) = paph(jl,kcbot(jl)) zrrho = (rd*ptenh(jl,kk+1))/paph(jl,kk+1) zdprho = (paph(jl,kk+1)-paph(jl,kk))*zrg zpmid = 0.5*(zpbase(jl)+paph(jl,kctop0(jl))) zentr = pentr(jl)*pmfu(jl,kk+1)*zdprho*zrrho llo1 = kk.lt.kcbot(jl).and.ldcum(jl) if(llo1) then zdmfde(jl) = zentr else zdmfde(jl) = 0.0 endif llo2 = llo1.and.ktype(jl).eq.2.and.((zpbase(jl)-paph(jl,kk)) & .lt.zdnoprc.or.paph(jl,kk).gt.zpmid) if(llo2) then zdmfen(jl) = zentr else zdmfen(jl) = 0.0 endif iklwmin = max(klwmin(jl),kctop0(jl)+2) llo2 = llo1.and.ktype(jl).eq.3.and.(kk.ge.iklwmin.or.pap(jl,kk) & .gt.zpmid) if (llo2) zdmfen(jl) = zentr llo2 = llo1.and.ktype(jl).eq.1 if (llo2) zdmfen(jl) = zentr ikb = kcbot(jl) zodetr(jl,kk) = 0. if (llo2.and.kk.le.khmin(jl).and.kk.ge.kctop0(jl)) then ikt = kctop0(jl) ikh = khmin(jl) if (ikh.gt.ikt) then zzmzk = -(pgeoh(jl,ikh)-pgeoh(jl,kk))*zrg ztmzk = -(pgeoh(jl,ikh)-pgeoh(jl,ikt))*zrg arg = 3.1415*(zzmzk/ztmzk)*0.5 zorgde = tan(arg)*3.1415*0.5/ztmzk zdprho = (paph(jl,kk+1)-paph(jl,kk))*(zrg*zrrho) zodetr(jl,kk) = min(zorgde,1.e-3)*pmfu(jl,kk+1)*zdprho end if end if enddo return end subroutine cuentr_new real function ssum ( n, x, ix ) implicit none real x(*) real zsum integer n, ix, jx, jl jx = 1 zsum = 0.0 do jl = 1, n zsum = zsum + x(jx) jx = jx + ix enddo ssum=zsum return end function ssum real function tlucua(tt) implicit none real zcvm3,zcvm4,tt if(tt-tmelt.gt.0.) then zcvm3=c3les zcvm4=c4les else zcvm3=c3ies zcvm4=c4ies end if tlucua=c2es*exp(zcvm3*(tt-tmelt)*(1./(tt-zcvm4))) return end function tlucua real function tlucub(tt) implicit none real z5alvcp,z5alscp,zcvm4,zcvm5,tt z5alvcp=c5les*alv/cpd z5alscp=c5ies*als/cpd if(tt-tmelt.gt.0.) then zcvm4=c4les zcvm5=z5alvcp else zcvm4=c4ies zcvm5=z5alscp end if tlucub=zcvm5*(1./(tt-zcvm4))**2 return end function tlucub real function tlucuc(tt) implicit none real zalvdcp,zalsdcp,tt,zldcp zalvdcp=alv/cpd zalsdcp=als/cpd if(tt-tmelt.gt.0.) then zldcp=zalvdcp else zldcp=zalsdcp end if tlucuc=zldcp return end function tlucuc end module module_cu_tiedtke