160 subroutine cu_ntiedtke_run(pu,pv,pt,pqv,tdi,qvdi,pqvf,ptf,clw,poz,pzz,prsl,prsi,pomg, &
161 evap,hfx,zprecc,lmask,lq,km,dt,dx,kbot,ktop,kcnv, &
162 ktrac,ud_mf,dd_mf,dt_mf,cnvw,cnvc,errmsg,errflg)
169 integer,
intent(in) :: lq, km, ktrac
170 real(kind=kind_phys),
intent(in ) :: dt
171 integer,
dimension( : ),
intent(in) :: lmask
172 real(kind=kind_phys),
dimension( : ),
intent(in ) :: evap, hfx, dx
173 real(kind=kind_phys),
dimension( :, : ),
intent(inout) :: pu, pv, pt, pqv
174 real(kind=kind_phys),
dimension( :, :),
intent(in ) :: tdi, qvdi, poz, prsl, pomg
175 real(kind=kind_phys),
dimension( :, :),
intent(in ),
optional :: pqvf, ptf
176 real(kind=kind_phys),
dimension( :, : ),
intent(in ) :: pzz, prsi
177 real(kind=kind_phys),
dimension( :, :, : ),
intent(inout ) :: clw
179 integer,
dimension( : ),
intent(out) :: kbot, ktop, kcnv
180 real(kind=kind_phys),
dimension( : ),
intent(out) :: zprecc
181 real(kind=kind_phys),
dimension (:, :),
intent(out),
optional :: ud_mf
182 real(kind=kind_phys),
dimension (:, :),
intent(out) :: dd_mf, dt_mf, cnvw, cnvc
185 character(len=*),
intent(out) :: errmsg
186 integer,
intent(out) :: errflg
189 real(kind=kind_phys) pum1(lq,km), pvm1(lq,km), ztt(lq,km), &
190 & ptte(lq,km), pqte(lq,km), pvom(lq,km), pvol(lq,km), &
191 & pverv(lq,km), pgeo(lq,km), pap(lq,km), paph(lq,km+1)
192 real(kind=kind_phys) pqhfl(lq), zqq(lq,km), &
193 & prsfc(lq), pssfc(lq), pcte(lq,km), &
194 & phhfl(lq), pgeoh(lq,km+1)
195 real(kind=kind_phys) ztp1(lq,km), zqp1(lq,km), ztu(lq,km), zqu(lq,km),&
196 & zlu(lq,km), zlude(lq,km), zmfu(lq,km), zmfd(lq,km), zmfude_rate(lq,km),&
197 & zqsat(lq,km), zrain(lq)
198 real(kind=kind_phys),
allocatable :: pcen(:,:,:),ptenc(:,:,:)
200 integer icbot(lq), ictop(lq), ktype(lq), lndj(lq)
203 real(kind=kind_phys) ztmst,fliq,fice,ztc,zalf,tt
204 integer i,j,k,k1,n,km1,ktracer
205 real(kind=kind_phys) ztpp1
206 real(kind=kind_phys) zew,zqs,zcor
224 pgeoh(j,km1)=pzz(j,1)
225 paph(j,km1)=prsi(j,1)
226 if(lmask(j).eq.1)
then
245 pverv(j,k1)=pomg(j,k)
247 pgeoh(j,k1)=pzz(j,k+1)
249 paph(j,k1)=prsi(j,k+1)
254 zcor = 1./(1.-vtmpc1*zqs)
256 pqte(j,k1)=pqvf(j,k)+(pqv(j,k)-qvdi(j,k))/ztmst
257 zqq(j,k1) =pqte(j,k1)
258 ptte(j,k1)=ptf(j,k)+(pt(j,k)-tdi(j,k))/ztmst
259 ztt(j,k1) =ptte(j,k1)
270 allocate(pcen(lq,km,ktracer))
271 allocate(ptenc(lq,km,ktracer))
276 pcen(j,k1,n) = clw(j,k,n+2)
283 allocate(pcen(lq,km,ktracer))
284 allocate(ptenc(lq,km,ktracer))
312 & (lq, km, km1, km-1, ztp1, &
313 & zqp1, pum1, pvm1, pverv, zqsat,&
314 & pqhfl, ztmst, pap, paph, pgeo, &
315 & ptte, pqte, pvom, pvol, prsfc,&
316 & pssfc, locum, ktracer, pcen, ptenc,&
317 & ktype, icbot, ictop, ztu, zqu, &
318 & zlu, zlude, zmfu, zmfd, zrain,&
319 & pcte, phhfl, lndj, pgeoh, zmfude_rate, dx)
326 if(pcte(j,k1).gt.0.)
then
327 fliq=foealfa(ztp1(j,k1))
329 clw(j,k,2)=clw(j,k,2)+fliq*pcte(j,k1)*ztmst
330 clw(j,k,1)=clw(j,k,1)+fice*pcte(j,k1)*ztmst
338 pt(j,k) = ztp1(j,k1)+(ptte(j,k1)-ztt(j,k1))*ztmst
339 pqv(j,k)= zqp1(j,k1)+(pqte(j,k1)-zqq(j,k1))*ztmst
340 ud_mf(j,k)= zmfu(j,k1)*ztmst
341 dd_mf(j,k)= -zmfd(j,k1)*ztmst
342 dt_mf(j,k)= zmfude_rate(j,k1)*ztmst
343 cnvw(j,k) = zlude(j,k1)*ztmst*g/(prsi(j,k)-prsi(j,k+1))
344 cnvc(j,k) = 0.04 * log(1. + 675. * ud_mf(j,k))
345 cnvc(j,k) = min(cnvc(j,k), 0.6)
346 cnvc(j,k) = max(cnvc(j,k), 0.0)
351 zprecc(j)=amax1(0.0,(prsfc(j)+pssfc(j))*ztmst*0.001)
352 kbot(j) = km-icbot(j)+1
353 ktop(j) = km-ictop(j)+1
354 if(ktype(j).eq.1 .or. ktype(j).eq.3)
then
365 pu(j,k)=pu(j,k)+pvom(j,k1)*ztmst
366 pv(j,k)=pv(j,k)+pvol(j,k1)*ztmst
397 subroutine cumastrn &
398 & (klon, klev, klevp1, klevm1, pten, &
399 & pqen, puen, pven, pverv, pqsen,&
400 & pqhfl, ztmst, pap, paph, pgeo, &
401 & ptte, pqte, pvom, pvol, prsfc,&
402 & pssfc, ldcum, ktrac, pcen, ptenc,&
403 & ktype, kcbot, kctop, ptu, pqu,&
404 & plu, plude, pmfu, pmfd, prain,&
405 & pcte, phhfl, lndj, zgeoh, pmfude_rate, dx)
465 integer klev,klon,ktrac,klevp1,klevm1
466 real(kind=kind_phys) pten(klon,klev), pqen(klon,klev),&
467 & puen(klon,klev), pven(klon,klev),&
468 & ptte(klon,klev), pqte(klon,klev),&
469 & pvom(klon,klev), pvol(klon,klev),&
470 & pqsen(klon,klev), pgeo(klon,klev),&
471 & pap(klon,klev), paph(klon,klevp1),&
472 & pverv(klon,klev), pqhfl(klon),&
474 real(kind=kind_phys) ptu(klon,klev), pqu(klon,klev),&
475 & plu(klon,klev), plude(klon,klev),&
476 & pmfu(klon,klev), pmfd(klon,klev),&
478 & prsfc(klon), pssfc(klon)
479 real(kind=kind_phys) ztenh(klon,klev), zqenh(klon,klev),&
480 & zgeoh(klon,klevp1), zqsenh(klon,klev),&
481 & ztd(klon,klev), zqd(klon,klev),&
482 & zmfus(klon,klev), zmfds(klon,klev),&
483 & zmfuq(klon,klev), zmfdq(klon,klev),&
484 & zdmfup(klon,klev), zdmfdp(klon,klev),&
485 & zmful(klon,klev), zrfl(klon),&
486 & zuu(klon,klev), zvu(klon,klev),&
487 & zud(klon,klev), zvd(klon,klev),&
489 real(kind=kind_phys) pmflxr(klon,klevp1), pmflxs(klon,klevp1)
490 real(kind=kind_phys) zhcbase(klon),&
491 & zmfub(klon), zmfub1(klon),&
493 real(kind=kind_phys) zsfl(klon), zdpmel(klon,klev),&
494 & pcte(klon,klev), zcape(klon),&
495 & zcape1(klon), zcape2(klon),&
496 & ztauc(klon), ztaubl(klon),&
498 real(kind=kind_phys) pcen(klon,klev,ktrac), ptenc(klon,klev,ktrac)
499 real(kind=kind_phys) wup(klon), zdqcv(klon)
500 real(kind=kind_phys) wbase(klon), zmfuub(klon)
501 real(kind=kind_phys) upbl(klon)
502 real(kind=kind_phys) dx(klon)
503 real(kind=kind_phys) pmfude_rate(klon,klev), pmfdde_rate(klon,klev)
504 real(kind=kind_phys) zmfuus(klon,klev), zmfdus(klon,klev)
505 real(kind=kind_phys) zmfudr(klon,klev), zmfddr(klon,klev)
506 real(kind=kind_phys) zuv2(klon,klev),ztenu(klon,klev),ztenv(klon,klev)
507 real(kind=kind_phys) zmfuvb(klon),zsum12(klon),zsum22(klon)
508 integer ilab(klon,klev), idtop(klon),&
509 & ictop0(klon), ilwmin(klon)
511 integer kcbot(klon), kctop(klon),&
512 & ktype(klon), lndj(klon)
513 logical ldcum(klon), lldcum(klon)
514 logical loddraf(klon), llddraf3(klon), llo1, llo2(klon)
517 real(kind=kind_phys) zcons,zcons2,zqumqe,zdqmin,zdh,zmfmax
518 real(kind=kind_phys) zalfaw,zalv,zqalv,zc5ldcp,zc4les,zhsat,zgam,zzz,zhhat
519 real(kind=kind_phys) zpbmpt,zro,zdz,zdp,zeps,zfac,wspeed
521 integer ikb,ikt,icum,itopm2
522 real(kind=kind_phys) ztmst,ztau,zerate,zderate,zmfa
523 real(kind=kind_phys) zmfs(klon),pmean(klev),zlon
524 real(kind=kind_phys) zduten,zdvten,ztdis,pgf_u,pgf_v
532 do jk = klev , 1 , -1
533 pmean(jk) = sum(pap(:,jk))/zlon
536 do jk = klev , 3 , -1
537 if ( pmean(jk)/pmean(klev)*1.013250e5 > 650.e2 ) p650 = jk
544 & (klon, klev, klevp1, klevm1, pten, &
545 & pqen, pqsen, puen, pven, pverv,&
546 & pgeo, paph, zgeoh, ztenh, zqenh,&
547 & zqsenh, ilwmin, ptu, pqu, ztd, &
548 & zqd, zuu, zvu, zud, zvd, &
549 & pmfu, pmfd, zmfus, zmfds, zmfuq,&
550 & zmfdq, zdmfup, zdmfdp, zdpmel, plu, &
560 & ( klon, klev, klevp1, klevm1, pqen,&
561 & ztenh, zqenh, zqsenh, zgeoh, paph,&
562 & phhfl, pqhfl, pgeo, pqsen, pap,&
563 & pten, lndj, ptu, pqu, ilab,&
564 & ldcum, kcbot, ictop0, ktype, wbase, plu, kdpl)
576 if(jk.ge.kcbot(jl) .and. ldcum(jl))
then
577 zdhpbl(jl)=zdhpbl(jl)+(alv*pqte(jl,jk)+cpd*ptte(jl,jk))&
578 & *(paph(jl,jk+1)-paph(jl,jk))
579 if(lndj(jl) .eq. 0)
then
580 wspeed = sqrt(puen(jl,jk)**2 + pven(jl,jk)**2)
581 upbl(jl) = upbl(jl) + wspeed*(paph(jl,jk+1)-paph(jl,jk))
590 zmfmax = (paph(jl,ikb)-paph(jl,ikb-1))*zcons2
591 if(ktype(jl) == 1)
then
592 zmfub(jl)= 0.1*zmfmax
593 else if ( ktype(jl) == 2 )
then
594 zqumqe = pqu(jl,ikb) + plu(jl,ikb) - zqenh(jl,ikb)
595 zdqmin = max(0.01*zqenh(jl,ikb),1.e-10)
596 zdh = cpd*(ptu(jl,ikb)-ztenh(jl,ikb)) + alv*zqumqe
597 zdh = g*max(zdh,1.e5*zdqmin)
598 if ( zdhpbl(jl) > 0. )
then
599 zmfub(jl) = zdhpbl(jl)/zdh
600 zmfub(jl) = min(zmfub(jl),zmfmax)
602 zmfub(jl) = 0.1*zmfmax
616 & (klon, klev, klevp1, klevm1, ztenh,&
617 & zqenh, puen, pven, pten, pqen,&
618 & pqsen, pgeo, zgeoh, pap, paph,&
619 & pqte, pverv, ilwmin, ldcum, zhcbase,&
620 & ktype, ilab, ptu, pqu, plu,&
621 & zuu, zvu, pmfu, zmfub,&
622 & zmfus, zmfuq, zmful, plude, zdmfup,&
623 & kcbot, kctop, ictop0, icum, ztmst,&
624 & zqsenh, zlglac, lndj, wup, wbase, kdpl, pmfude_rate )
630 if ( ldcum(jl) )
then
633 zpbmpt = paph(jl,ikb) - paph(jl,itopm2)
634 if ( ktype(jl) == 1 .and. zpbmpt < zdnoprc ) ktype(jl) = 2
635 if ( ktype(jl) == 2 .and. zpbmpt >= zdnoprc ) ktype(jl) = 1
636 ictop0(jl) = kctop(jl)
638 zrfl(jl)=zdmfup(jl,1)
643 zrfl(jl)=zrfl(jl)+zdmfup(jl,jk)
665 & kcbot, kctop, lndj, ldcum, &
666 & ztenh, zqenh, puen, pven, &
667 & pten, pqsen, pgeo, &
668 & zgeoh, paph, ptu, pqu, plu, &
669 & zuu, zvu, zmfub, zrfl, &
670 & ztd, zqd, zud, zvd, &
671 & pmfd, zmfds, zmfdq, zdmfdp, &
676 & ( klon, klev, loddraf, &
677 & ztenh, zqenh, puen, pven, &
678 & pgeo, zgeoh, paph, zrfl, &
679 & ztd, zqd, zud, zvd, pmfu, &
680 & pmfd, zmfds, zmfdq, zdmfdp, pmfdde_rate )
691 if(ldcum(jl) .and. ktype(jl) .eq. 1)
then
700 ztauc(jl) = (zgeoh(jl,ikt)-zgeoh(jl,ikb)) / &
701 ((2.+ min(15.0,wup(jl)))*g)
702 if(lndj(jl) .eq. 0)
then
703 upbl(jl) = 2.+ upbl(jl)/(paph(jl,klev+1)-paph(jl,ikb))
704 ztaubl(jl) = (zgeoh(jl,ikb)-zgeoh(jl,klev+1))/(g*upbl(jl))
705 ztaubl(jl) = min(300., ztaubl(jl))
707 ztaubl(jl) = ztauc(jl)
714 llo1 = ldcum(jl) .and. ktype(jl) .eq. 1
715 if ( llo1 .and. jk <= kcbot(jl) .and. jk > kctop(jl) )
then
717 zdz = pgeo(jl,jk-1)-pgeo(jl,jk)
718 zdp = pap(jl,jk)-pap(jl,jk-1)
719 zheat(jl) = zheat(jl) + ((pten(jl,jk-1)-pten(jl,jk)+zdz*rcpd) / &
720 ztenh(jl,jk)+vtmpc1*(pqen(jl,jk-1)-pqen(jl,jk))) * &
721 (g*(pmfu(jl,jk)+pmfd(jl,jk)))
722 zcape1(jl) = zcape1(jl) + ((ptu(jl,jk)-ztenh(jl,jk))/ztenh(jl,jk) + &
723 vtmpc1*(pqu(jl,jk)-zqenh(jl,jk))-plu(jl,jk))*zdp
726 if ( llo1 .and. jk >= kcbot(jl) )
then
727 if((paph(jl,klev+1)-paph(jl,kdpl(jl)))<50.e2)
then
728 zdp = paph(jl,jk+1)-paph(jl,jk)
729 zcape2(jl) = zcape2(jl) + ztaubl(jl)* &
730 ((1.+vtmpc1*pqen(jl,jk))*ptte(jl,jk)+vtmpc1*pten(jl,jk)*pqte(jl,jk))*zdp
737 if(ldcum(jl).and.ktype(jl).eq.1)
then
740 ztau = ztauc(jl) * (1.+1.33e-5*dx(jl))
741 ztau = max(ztmst,ztau)
742 ztau = max(720.,ztau)
743 ztau = min(10800.,ztau)
745 zcape2(jl)= max(0.,zcape2(jl))
746 zcape(jl) = max(0.,min(zcape1(jl)-zcape2(jl),5000.))
748 zcape(jl) = max(0.,min(zcape1(jl),5000.))
750 zheat(jl) = max(1.e-4,zheat(jl))
751 zmfub1(jl) = (zcape(jl)*zmfub(jl))/(zheat(jl)*ztau)
752 zmfub1(jl) = max(zmfub1(jl),0.001)
753 zmfmax=(paph(jl,ikb)-paph(jl,ikb-1))*zcons2
754 zmfub1(jl)=min(zmfub1(jl),zmfmax)
762 if(ldcum(jl) .and. ktype(jl) .eq. 2)
then
764 if(pmfd(jl,ikb).lt.0.0 .and. loddraf(jl))
then
765 zeps=-pmfd(jl,ikb)/max(zmfub(jl),cmfcmin)
769 zqumqe=pqu(jl,ikb)+plu(jl,ikb)- &
770 & zeps*zqd(jl,ikb)-(1.-zeps)*zqenh(jl,ikb)
771 zdqmin=max(0.01*zqenh(jl,ikb),cmfcmin)
772 zmfmax=(paph(jl,ikb)-paph(jl,ikb-1))*zcons2
774 zdh=cpd*(ptu(jl,ikb)-zeps*ztd(jl,ikb)- &
775 & (1.-zeps)*ztenh(jl,ikb))+alv*zqumqe
776 zdh=g*max(zdh,1.e5*zdqmin)
777 if(zdhpbl(jl).gt.0.)
then
778 zmfub1(jl)=zdhpbl(jl)/zdh
780 zmfub1(jl) = zmfub(jl)
782 zmfub1(jl) = min(zmfub1(jl),zmfmax)
787 if(ldcum(jl) .and. ktype(jl) .eq. 3 )
then
788 zmfub1(jl) = zmfub(jl)
798 zfac=zmfub1(jl)/max(zmfub(jl),cmfcmin)
799 pmfd(jl,jk)=pmfd(jl,jk)*zfac
800 zmfds(jl,jk)=zmfds(jl,jk)*zfac
801 zmfdq(jl,jk)=zmfdq(jl,jk)*zfac
802 zdmfdp(jl,jk)=zdmfdp(jl,jk)*zfac
803 pmfdde_rate(jl,jk) = pmfdde_rate(jl,jk)*zfac
811 if ( ldcum(jl) ) zmfs(jl) = zmfub1(jl)/max(cmfcmin,zmfub(jl))
815 if ( ldcum(jl) .and. jk >= kctop(jl)-1 )
then
818 zdz = ((paph(jl,klev+1)-paph(jl,jk))/(paph(jl,klev+1)-paph(jl,ikb)))
819 pmfu(jl,jk) = pmfu(jl,ikb)*zdz
821 zmfmax = (paph(jl,jk)-paph(jl,jk-1))*zcons2
822 if ( pmfu(jl,jk)*zmfs(jl) > zmfmax )
then
823 zmfs(jl) = min(zmfs(jl),zmfmax/pmfu(jl,jk))
830 if ( ldcum(jl) .and. jk <= kcbot(jl) .and. jk >= kctop(jl)-1 )
then
831 pmfu(jl,jk) = pmfu(jl,jk)*zmfs(jl)
832 zmfus(jl,jk) = zmfus(jl,jk)*zmfs(jl)
833 zmfuq(jl,jk) = zmfuq(jl,jk)*zmfs(jl)
834 zmful(jl,jk) = zmful(jl,jk)*zmfs(jl)
835 zdmfup(jl,jk) = zdmfup(jl,jk)*zmfs(jl)
836 plude(jl,jk) = plude(jl,jk)*zmfs(jl)
837 pmfude_rate(jl,jk) = pmfude_rate(jl,jk)*zmfs(jl)
845 if ( ktype(jl) == 2 .and. &
846 kcbot(jl) == kctop(jl) .and. kcbot(jl) >= klev-1 )
then
852 if ( .not. lmfscv .or. .not. lmfpen )
then
855 if ( (.not. lmfscv .and. ktype(jl) == 2) .or. &
856 (.not. lmfpen .and. ktype(jl) == 1) )
then
866 if ( loddraf(jl) .and. idtop(jl) <= kctop(jl) )
then
867 idtop(jl) = kctop(jl) + 1
872 if ( loddraf(jl) )
then
873 if ( jk < idtop(jl) )
then
877 pmfdde_rate(jl,jk) = 0.
879 else if ( jk == idtop(jl) )
then
880 pmfdde_rate(jl,jk) = 0.
891 & ( klon, klev, ztmst &
892 & , pten, pqen, pqsen, ztenh, zqenh &
893 & , paph, pap, zgeoh, lndj, ldcum &
894 & , kcbot, kctop, idtop, itopm2 &
896 & , pmfu, pmfd, zmfus, zmfds &
897 & , zmfuq, zmfdq, zmful, plude &
898 & , zdmfup, zdmfdp, zdpmel, zlglac &
899 & , prain, pmfdde_rate, pmflxr, pmflxs )
908 if ( loddraf(jl) .and. jk >= idtop(jl)-1 )
then
909 zmfmax = pmfu(jl,jk)*0.98
910 if ( pmfd(jl,jk)+zmfmax+1.e-15 < 0. )
then
911 zmfs(jl) = min(zmfs(jl),-zmfmax/pmfd(jl,jk))
919 if ( zmfs(jl) < 1. .and. jk >= idtop(jl)-1 )
then
920 pmfd(jl,jk) = pmfd(jl,jk)*zmfs(jl)
921 zmfds(jl,jk) = zmfds(jl,jk)*zmfs(jl)
922 zmfdq(jl,jk) = zmfdq(jl,jk)*zmfs(jl)
923 pmfdde_rate(jl,jk) = pmfdde_rate(jl,jk)*zmfs(jl)
924 zmfuub(jl) = zmfuub(jl) - (1.-zmfs(jl))*zdmfdp(jl,jk)
925 pmflxr(jl,jk+1) = pmflxr(jl,jk+1) + zmfuub(jl)
926 zdmfdp(jl,jk) = zdmfdp(jl,jk)*zmfs(jl)
933 if ( loddraf(jl) .and. jk >= idtop(jl)-1 )
then
934 zerate = -pmfd(jl,jk) + pmfd(jl,jk-1) + pmfdde_rate(jl,jk)
935 if ( zerate < 0. )
then
936 pmfdde_rate(jl,jk) = pmfdde_rate(jl,jk) - zerate
939 if ( ldcum(jl) .and. jk >= kctop(jl)-1 )
then
940 zerate = pmfu(jl,jk) - pmfu(jl,jk+1) + pmfude_rate(jl,jk)
941 if ( zerate < 0. )
then
942 pmfude_rate(jl,jk) = pmfude_rate(jl,jk) - zerate
944 zdmfup(jl,jk) = pmflxr(jl,jk+1) + pmflxs(jl,jk+1) - &
945 pmflxr(jl,jk) - pmflxs(jl,jk)
953 if ( loddraf(jl) )
then
956 if ( zmfdq(jl,jk) < 0.3*zmfdq(jl,ik) )
then
957 zmfdq(jl,jk) = 0.3*zmfdq(jl,ik)
966 if ( ldcum(jl) .and. jk >= kctop(jl)-1 .and. jk < kcbot(jl) )
then
967 zdz = ztmst*g/(paph(jl,jk+1)-paph(jl,jk))
968 zmfa = zmfuq(jl,jk+1) + zmfdq(jl,jk+1) - &
969 zmfuq(jl,jk) - zmfdq(jl,jk) + &
970 zmful(jl,jk+1) - zmful(jl,jk) + zdmfup(jl,jk)
971 zmfa = (zmfa-plude(jl,jk))*zdz
972 if ( pqen(jl,jk)+zmfa < 0. )
then
973 plude(jl,jk) = plude(jl,jk) + 2.*(pqen(jl,jk)+zmfa)/zdz
975 if ( plude(jl,jk) < 0. ) plude(jl,jk) = 0.
977 if ( .not. ldcum(jl) ) pmfude_rate(jl,jk) = 0.
978 if ( abs(pmfd(jl,jk-1)) < 1.0e-20 ) pmfdde_rate(jl,jk) = 0.
983 prsfc(jl) = pmflxr(jl,klev+1)
984 pssfc(jl) = pmflxs(jl,klev+1)
990 call cudtdqn(klon,klev,itopm2,kctop,idtop,ldcum,loddraf, &
991 ztmst,paph,zgeoh,pgeo,pten,ztenh,pqen,zqenh,pqsen, &
992 zlglac,plude,pmfu,pmfd,zmfus,zmfds,zmfuq,zmfdq,zmful, &
993 zdmfup,zdmfdp,zdpmel,ptte,pqte,pcte)
998 do jk = klev-1 , 2 , -1
1001 if ( ldcum(jl) )
then
1002 if ( jk == kcbot(jl) .and. ktype(jl) < 3 )
then
1004 zuu(jl,jk) = puen(jl,ikb-1)
1005 zvu(jl,jk) = pven(jl,ikb-1)
1006 else if ( jk == kcbot(jl) .and. ktype(jl) == 3 )
then
1007 zuu(jl,jk) = puen(jl,jk-1)
1008 zvu(jl,jk) = pven(jl,jk-1)
1010 if ( jk < kcbot(jl) .and. jk >= kctop(jl) )
then
1011 if(momtrans .eq. 1)
then
1013 if ( ktype(jl) == 1 .or. ktype(jl) == 3 ) zfac = 2.
1014 if ( ktype(jl) == 1 .and. jk <= kctop(jl)+2 ) zfac = 3.
1015 zerate = pmfu(jl,jk) - pmfu(jl,ik) + &
1016 (1.+zfac)*pmfude_rate(jl,jk)
1017 zderate = (1.+zfac)*pmfude_rate(jl,jk)
1018 zmfa = 1./max(cmfcmin,pmfu(jl,jk))
1019 zuu(jl,jk) = (zuu(jl,ik)*pmfu(jl,ik) + &
1020 zerate*puen(jl,jk)-zderate*zuu(jl,ik))*zmfa
1021 zvu(jl,jk) = (zvu(jl,ik)*pmfu(jl,ik) + &
1022 zerate*pven(jl,jk)-zderate*zvu(jl,ik))*zmfa
1024 if(ktype(jl) == 1 .or. ktype(jl) == 3)
then
1025 pgf_u = -0.7*0.5*(pmfu(jl,ik)*(puen(jl,ik)-puen(jl,jk))+&
1026 pmfu(jl,jk)*(puen(jl,jk)-puen(jl,jk-1)))
1027 pgf_v = -0.7*0.5*(pmfu(jl,ik)*(pven(jl,ik)-pven(jl,jk))+&
1028 pmfu(jl,jk)*(pven(jl,jk)-pven(jl,jk-1)))
1033 zerate = pmfu(jl,jk) - pmfu(jl,ik) + pmfude_rate(jl,jk)
1034 zderate = pmfude_rate(jl,jk)
1035 zmfa = 1./max(cmfcmin,pmfu(jl,jk))
1036 zuu(jl,jk) = (zuu(jl,ik)*pmfu(jl,ik) + &
1037 zerate*puen(jl,jk)-zderate*zuu(jl,ik)+pgf_u)*zmfa
1038 zvu(jl,jk) = (zvu(jl,ik)*pmfu(jl,ik) + &
1039 zerate*pven(jl,jk)-zderate*zvu(jl,ik)+pgf_v)*zmfa
1050 if ( ldcum(jl) )
then
1051 if ( jk == idtop(jl) )
then
1052 zud(jl,jk) = 0.5*(zuu(jl,jk)+puen(jl,ik))
1053 zvd(jl,jk) = 0.5*(zvu(jl,jk)+pven(jl,ik))
1054 else if ( jk > idtop(jl) )
then
1055 zerate = -pmfd(jl,jk) + pmfd(jl,ik) + pmfdde_rate(jl,jk)
1056 zmfa = 1./min(-cmfcmin,pmfd(jl,jk))
1057 zud(jl,jk) = (zud(jl,ik)*pmfd(jl,ik) - &
1058 zerate*puen(jl,ik)+pmfdde_rate(jl,jk)*zud(jl,ik))*zmfa
1059 zvd(jl,jk) = (zvd(jl,ik)*pmfd(jl,ik) - &
1060 zerate*pven(jl,ik)+pmfdde_rate(jl,jk)*zvd(jl,ik))*zmfa
1072 if ( ldcum(jl) .and. jk >= kctop(jl)-1 )
then
1073 zmfmax = (paph(jl,jk)-paph(jl,jk-1))*zcons
1074 if ( pmfu(jl,jk) > zmfmax .and. jk >= kctop(jl) )
then
1075 zmfs(jl) = min(zmfs(jl),zmfmax/pmfu(jl,jk))
1082 zmfuus(jl,jk) = pmfu(jl,jk)
1083 zmfdus(jl,jk) = pmfd(jl,jk)
1084 if ( ldcum(jl) .and. jk >= kctop(jl)-1 )
then
1085 zmfuus(jl,jk) = pmfu(jl,jk)*zmfs(jl)
1086 zmfdus(jl,jk) = pmfd(jl,jk)*zmfs(jl)
1094 ztenu(jl,jk) = pvom(jl,jk)
1095 ztenv(jl,jk) = pvol(jl,jk)
1099 call cududvn(klon,klev,itopm2,ktype,kcbot,kctop, &
1100 ldcum,ztmst,paph,puen,pven,zmfuus,zmfdus,zuu, &
1101 zud,zvu,zvd,pvom,pvol)
1111 if ( ldcum(jl) .and. jk >= kctop(jl)-1 )
then
1112 zdz = (paph(jl,jk+1)-paph(jl,jk))
1113 zduten = pvom(jl,jk) - ztenu(jl,jk)
1114 zdvten = pvol(jl,jk) - ztenv(jl,jk)
1115 zuv2(jl,jk) = sqrt(zduten**2+zdvten**2)
1116 zsum22(jl) = zsum22(jl) + zuv2(jl,jk)*zdz
1117 zsum12(jl) = zsum12(jl) - &
1118 (puen(jl,jk)*zduten+pven(jl,jk)*zdvten)*zdz
1124 if ( ldcum(jl) .and. jk>=kctop(jl)-1 )
then
1125 ztdis = rcpd*zsum12(jl)*zuv2(jl,jk)/max(1.e-15,zsum22(jl))
1126 ptte(jl,jk) = ptte(jl,jk) + ztdis
1137 if ( .not. lmfscv .or. .not. lmfpen )
then
1140 if ( llo2(jl) .and. jk >= kctop(jl)-1 )
then
1141 ptu(jl,jk) = pten(jl,jk)
1142 pqu(jl,jk) = pqen(jl,jk)
1144 pmfude_rate(jl,jk) = 0.
1145 pmfdde_rate(jl,jk) = 0.
1150 if ( llo2(jl) )
then
1151 kctop(jl) = klev - 1
1152 kcbot(jl) = klev - 1
1161 if ( ktrac > 0 )
then
1164 if ( ldcum(jl) .and. ktype(jl) /= 3 .and. &
1165 kcbot(jl)-kctop(jl) >= 1 )
then
1167 llddraf3(jl) = loddraf(jl)
1169 lldcum(jl) = .false.
1170 llddraf3(jl) = .false.
1177 if ( lldcum(jl) .and. jk >= kctop(jl) )
then
1178 zmfmax = (paph(jl,jk)-paph(jl,jk-1))*0.8*zcons
1179 if ( pmfu(jl,jk) > zmfmax )
then
1180 zmfs(jl) = min(zmfs(jl),zmfmax/pmfu(jl,jk))
1188 if ( lldcum(jl) .and. jk >= kctop(jl)-1 )
then
1189 zmfuus(jl,jk) = pmfu(jl,jk)*zmfs(jl)
1190 zmfudr(jl,jk) = pmfude_rate(jl,jk)*zmfs(jl)
1195 if ( llddraf3(jl) .and. jk >= idtop(jl)-1 )
then
1196 zmfdus(jl,jk) = pmfd(jl,jk)*zmfs(jl)
1197 zmfddr(jl,jk) = pmfdde_rate(jl,jk)*zmfs(jl)
1205 call cuctracer(klon,klev,ktrac,kctop,idtop, &
1206 lldcum,llddraf3,ztmst,paph,zmfuus,zmfdus, &
1207 zmfudr,zmfddr,pcen,ptenc)
1344 subroutine cutypen &
1345 & ( klon, klev, klevp1, klevm1, pqen,&
1346 & ptenh, pqenh, pqsenh, pgeoh, paph,&
1347 & hfx, qfx, pgeo, pqsen, pap,&
1348 & pten, lndj, cutu, cuqu, culab,&
1349 & ldcum, cubot, cutop, ktype, wbase, culu, kdpl )
1393 integer klon, klev, klevp1, klevm1
1394 real(kind=kind_phys) ptenh(klon,klev), pqenh(klon,klev),&
1395 & pqsen(klon,klev), pqsenh(klon,klev),&
1396 & pgeoh(klon,klevp1), paph(klon,klevp1),&
1397 & pap(klon,klev), pqen(klon,klev)
1398 real(kind=kind_phys) pten(klon,klev)
1399 real(kind=kind_phys) ptu(klon,klev),pqu(klon,klev),plu(klon,klev)
1400 real(kind=kind_phys) pgeo(klon,klev)
1401 integer klab(klon,klev)
1402 integer kctop(klon),kcbot(klon)
1404 real(kind=kind_phys) qfx(klon),hfx(klon)
1405 real(kind=kind_phys) zph(klon)
1407 logical loflag(klon), deepflag(klon), resetflag(klon)
1410 real(kind=kind_phys) cutu(klon,klev), cuqu(klon,klev), culu(klon,klev)
1411 integer culab(klon,klev)
1412 real(kind=kind_phys) wbase(klon)
1413 integer ktype(klon),cubot(klon),cutop(klon),kdpl(klon)
1417 real(kind=kind_phys) zqold(klon)
1418 real(kind=kind_phys) rho, part1, part2, root, conw, deltt, deltq
1419 real(kind=kind_phys) eta(klon),dz(klon),coef(klon)
1420 real(kind=kind_phys) dhen(klon,klev), dh(klon,klev)
1421 real(kind=kind_phys) plude(klon,klev)
1422 real(kind=kind_phys) kup(klon,klev)
1423 real(kind=kind_phys) vptu(klon,klev),vten(klon,klev)
1424 real(kind=kind_phys) zbuo(klon,klev),abuoy(klon,klev)
1426 real(kind=kind_phys) zz,zdken,zdq
1427 real(kind=kind_phys) fscale,crirh1,pp
1428 real(kind=kind_phys) atop1,atop2,abot
1429 real(kind=kind_phys) tmix,zmix,qmix,pmix
1430 real(kind=kind_phys) zlglac,dp
1431 integer nk,is,ikb,ikt
1433 real(kind=kind_phys) zqsu,zcor,zdp,zesdp,zalfaw,zfacw,zfaci,zfac,zdsdp,zdqsdt,zdtdp
1434 real(kind=kind_phys) zpdifftop, zpdiffbot
1435 integer zcbase(klon), itoppacel(klon)
1436 integer jl,jk,ik,icall,levels
1437 logical needreset, lldcum(klon)
1455 plu(jl,jk)=culu(jl,jk)
1456 ptu(jl,jk)=cutu(jl,jk)
1457 pqu(jl,jk)=cuqu(jl,jk)
1458 klab(jl,jk)=culab(jl,jk)
1471 lldcum(jl) = .false.
1479 if(jk .eq. klevm1)
then
1482 & (rd*(pten(jl,klev)*(1.+vtmpc1*pqen(jl,klev))))
1483 hfx(jl) = hfx(jl)*rho*cpd
1484 qfx(jl) = qfx(jl)*rho
1485 part1 = 1.5*0.4*pgeo(jl,klev)/ &
1486 & (rho*pten(jl,klev))
1487 part2 = -hfx(jl)*rcpd-vtmpc1*pten(jl,klev)*qfx(jl)
1488 root = 0.001-part1*part2
1489 if(part2 .lt. 0.)
then
1490 conw = 1.2*(root)**t13
1491 deltt = max(1.5*hfx(jl)/(rho*cpd*conw),0.)
1492 deltq = max(1.5*qfx(jl)/(rho*conw),0.)
1493 kup(jl,klev) = 0.5*(conw**2)
1494 pqu(jl,klev)= pqenh(jl,klev) + deltq
1495 dhen(jl,klev)= pgeoh(jl,klev) + ptenh(jl,klev)*cpd
1496 dh(jl,klev) = dhen(jl,klev) + deltt*cpd
1497 ptu(jl,klev) = (dh(jl,klev)-pgeoh(jl,klev))*rcpd
1498 vptu(jl,klev)=ptu(jl,klev)*(1.+vtmpc1*pqu(jl,klev))
1499 vten(jl,klev)=ptenh(jl,klev)*(1.+vtmpc1*pqenh(jl,klev))
1500 zbuo(jl,klev)=(vptu(jl,klev)-vten(jl,klev))/vten(jl,klev)
1503 loflag(jl) = .false.
1519 eta(jl) = 0.55/(pgeo(jl,jk)*zrg)+1.e-4
1520 dz(jl) = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
1521 coef(jl)= 0.5*eta(jl)*dz(jl)
1522 dhen(jl,jk) = pgeoh(jl,jk) + cpd*ptenh(jl,jk)
1523 dh(jl,jk) = (coef(jl)*(dhen(jl,jk+1)+dhen(jl,jk))&
1524 & +(1.-coef(jl))*dh(jl,jk+1))/(1.+coef(jl))
1525 pqu(jl,jk) =(coef(jl)*(pqenh(jl,jk+1)+pqenh(jl,jk))&
1526 & +(1.-coef(jl))*pqu(jl,jk+1))/(1.+coef(jl))
1527 ptu(jl,jk) = (dh(jl,jk)-pgeoh(jl,jk))*rcpd
1528 zqold(jl) = pqu(jl,jk)
1535 call cuadjtqn(klon,klev,ik,zph,ptu,pqu,loflag,icall)
1537 if( loflag(jl) )
then
1538 zdq = max((zqold(jl) - pqu(jl,jk)),0.)
1539 plu(jl,jk) = plu(jl,jk+1) + zdq
1540 zlglac=zdq*((1.-foealfa(ptu(jl,jk))) - &
1541 (1.-foealfa(ptu(jl,jk+1))))
1542 plu(jl,jk) = min(plu(jl,jk),5.e-3)
1543 dh(jl,jk) = pgeoh(jl,jk) + cpd*(ptu(jl,jk)+ralfdcp*zlglac)
1545 vptu(jl,jk) = ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk)-plu(jl,jk))+&
1547 vten(jl,jk) = ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))
1548 zbuo(jl,jk) = (vptu(jl,jk) - vten(jl,jk))/vten(jl,jk)
1549 abuoy(jl,jk)=(zbuo(jl,jk)+zbuo(jl,jk+1))*0.5*g
1550 atop1 = 1.0 - 2.*coef(jl)
1551 atop2 = 2.0*dz(jl)*abuoy(jl,jk)
1552 abot = 1.0 + 2.*coef(jl)
1553 kup(jl,jk) = (atop1*kup(jl,jk+1) + atop2) / abot
1556 if ( plu(jl,jk) > 0. .and. klab(jl,jk+1) == 1 )
then
1558 zqsu = foeewm(ptu(jl,ik))/paph(jl,ik)
1559 zqsu = min(0.5,zqsu)
1560 zcor = 1./(1.-vtmpc1*zqsu)
1562 zdq = min(0.,pqu(jl,ik)-zqsu)
1563 zalfaw = foealfa(ptu(jl,ik))
1564 zfacw = c5les/((ptu(jl,ik)-c4les)**2)
1565 zfaci = c5ies/((ptu(jl,ik)-c4ies)**2)
1566 zfac = zalfaw*zfacw + (1.-zalfaw)*zfaci
1567 zesdp = foeewm(ptu(jl,ik))/paph(jl,ik)
1568 zcor = 1./(1.-vtmpc1*zesdp)
1569 zdqsdt = zfac*zcor*zqsu
1570 zdtdp = rd*ptu(jl,ik)/(cpd*paph(jl,ik))
1571 zdp = zdq/(zdqsdt*zdtdp)
1572 zcbase(jl) = paph(jl,ik) + zdp
1574 zpdifftop = zcbase(jl) - paph(jl,jk)
1575 zpdiffbot = paph(jl,jk+1) - zcbase(jl)
1576 if ( zpdifftop > zpdiffbot .and. kup(jl,jk+1) > 0. )
then
1577 ikb = min(klev-1,jk+1)
1581 plu(jl,jk+1) = 1.0e-8
1582 else if ( zpdifftop <= zpdiffbot .and.kup(jl,jk) > 0. )
then
1588 if(kup(jl,jk) .lt. 0.)
then
1589 loflag(jl) = .false.
1590 if(plu(jl,jk+1) .gt. 0.)
then
1594 lldcum(jl) = .false.
1597 if(plu(jl,jk) .gt. 0.)
then
1611 if(paph(jl,ikb) - paph(jl,ikt) > zdnoprc) lldcum(jl) = .false.
1615 wbase(jl) = sqrt(max(2.*kup(jl,ikb),0.))
1632 culab(jl,jk) = klab(jl,jk)
1633 cutu(jl,jk) = ptu(jl,jk)
1634 cuqu(jl,jk) = pqu(jl,jk)
1635 culu(jl,jk) = plu(jl,jk)
1650 deepflag(jl) = .false.
1655 if((paph(jl,klev+1)-paph(jl,jk)) .lt. 350.e2) itoppacel(jl) = jk
1659 do levels=klevm1-1,klevm1-20,-1
1680 lldcum(jl) = .false.
1681 resetflag(jl)= .false.
1682 loflag(jl) = (.not. deepflag(jl)) .and. (levels.ge.itoppacel(jl))
1696 if(jk .eq. levels)
then
1699 if((paph(jl,klev+1)-paph(jl,jk)) < 60.e2)
then
1705 if(pmix < 50.e2)
then
1706 dp = paph(jl,nk) - paph(jl,nk-1)
1707 tmix=tmix+dp*ptenh(jl,nk)
1708 qmix=qmix+dp*pqenh(jl,nk)
1709 zmix=zmix+dp*pgeoh(jl,nk)
1722 pqu(jl,jk+1) = qmix + deltq
1723 dhen(jl,jk+1)= zmix + tmix*cpd
1724 dh(jl,jk+1) = dhen(jl,jk+1) + deltt*cpd
1725 ptu(jl,jk+1) = (dh(jl,jk+1)-pgeoh(jl,jk+1))*rcpd
1728 vptu(jl,jk+1)=ptu(jl,jk+1)*(1.+vtmpc1*pqu(jl,jk+1))
1729 vten(jl,jk+1)=ptenh(jl,jk+1)*(1.+vtmpc1*pqenh(jl,jk+1))
1730 zbuo(jl,jk+1)=(vptu(jl,jk+1)-vten(jl,jk+1))/vten(jl,jk+1)
1739 fscale = min(1.,(pqsen(jl,jk)/pqsen(jl,levels))**3)
1740 eta(jl) = 1.75e-3*fscale
1741 dz(jl) = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
1742 coef(jl)= 0.5*eta(jl)*dz(jl)
1743 dhen(jl,jk) = pgeoh(jl,jk) + cpd*ptenh(jl,jk)
1744 dh(jl,jk) = (coef(jl)*(dhen(jl,jk+1)+dhen(jl,jk))&
1745 & +(1.-coef(jl))*dh(jl,jk+1))/(1.+coef(jl))
1746 pqu(jl,jk) =(coef(jl)*(pqenh(jl,jk+1)+pqenh(jl,jk))&
1747 & +(1.-coef(jl))*pqu(jl,jk+1))/(1.+coef(jl))
1748 ptu(jl,jk) = (dh(jl,jk)-pgeoh(jl,jk))*rcpd
1749 zqold(jl) = pqu(jl,jk)
1756 call cuadjtqn(klon,klev,ik,zph,ptu,pqu,loflag,icall)
1759 if( loflag(jl) )
then
1760 zdq = max((zqold(jl) - pqu(jl,jk)),0.)
1761 plu(jl,jk) = plu(jl,jk+1) + zdq
1762 zlglac=zdq*((1.-foealfa(ptu(jl,jk))) - &
1763 (1.-foealfa(ptu(jl,jk+1))))
1764 plu(jl,jk) = 0.5*plu(jl,jk)
1765 dh(jl,jk) = pgeoh(jl,jk) + cpd*(ptu(jl,jk)+ralfdcp*zlglac)
1767 vptu(jl,jk) = ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk)-plu(jl,jk))+&
1769 vten(jl,jk) = ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))
1770 zbuo(jl,jk) = (vptu(jl,jk) - vten(jl,jk))/vten(jl,jk)
1771 abuoy(jl,jk)=(zbuo(jl,jk)+zbuo(jl,jk+1))*0.5*g
1772 atop1 = 1.0 - 2.*coef(jl)
1773 atop2 = 2.0*dz(jl)*abuoy(jl,jk)
1774 abot = 1.0 + 2.*coef(jl)
1775 kup(jl,jk) = (atop1*kup(jl,jk+1) + atop2) / abot
1777 if ( plu(jl,jk) > 0. .and. klab(jl,jk+1) == 1 )
then
1779 zqsu = foeewm(ptu(jl,ik))/paph(jl,ik)
1780 zqsu = min(0.5,zqsu)
1781 zcor = 1./(1.-vtmpc1*zqsu)
1783 zdq = min(0.,pqu(jl,ik)-zqsu)
1784 zalfaw = foealfa(ptu(jl,ik))
1785 zfacw = c5les/((ptu(jl,ik)-c4les)**2)
1786 zfaci = c5ies/((ptu(jl,ik)-c4ies)**2)
1787 zfac = zalfaw*zfacw + (1.-zalfaw)*zfaci
1788 zesdp = foeewm(ptu(jl,ik))/paph(jl,ik)
1789 zcor = 1./(1.-vtmpc1*zesdp)
1790 zdqsdt = zfac*zcor*zqsu
1791 zdtdp = rd*ptu(jl,ik)/(cpd*paph(jl,ik))
1792 zdp = zdq/(zdqsdt*zdtdp)
1793 zcbase(jl) = paph(jl,ik) + zdp
1795 zpdifftop = zcbase(jl) - paph(jl,jk)
1796 zpdiffbot = paph(jl,jk+1) - zcbase(jl)
1797 if ( zpdifftop > zpdiffbot .and. kup(jl,jk+1) > 0. )
then
1798 ikb = min(klev-1,jk+1)
1802 plu(jl,jk+1) = 1.0e-8
1803 else if ( zpdifftop <= zpdiffbot .and.kup(jl,jk) > 0. )
then
1809 if(kup(jl,jk) .lt. 0.)
then
1810 loflag(jl) = .false.
1811 if(plu(jl,jk+1) .gt. 0.)
then
1815 lldcum(jl) = .false.
1818 if(plu(jl,jk) .gt. 0.)
then
1833 if(paph(jl,ikb) - paph(jl,ikt) < zdnoprc) lldcum(jl) = .false.
1837 deepflag(jl) = .true.
1838 wbase(jl) = sqrt(max(2.*kup(jl,ikb),0.))
1843 resetflag(jl)= .true.
1850 if(resetflag(jl))
then
1853 if(jk .le. ikb .and. jk .ge. ikt )
then
1854 culab(jl,jk) = klab(jl,jk)
1855 cutu(jl,jk) = ptu(jl,jk)
1856 cuqu(jl,jk) = pqu(jl,jk)
1857 culu(jl,jk) = plu(jl,jk)
1860 cutu(jl,jk) = ptenh(jl,jk)
1861 cuqu(jl,jk) = pqenh(jl,jk)
1864 if ( jk .lt. ikt ) culab(jl,jk) = 0
1879 & (klon, klev, klevp1, klevm1, ptenh,&
1880 & pqenh, puen, pven, pten, pqen,&
1881 & pqsen, pgeo, pgeoh, pap, paph,&
1882 & pqte, pverv, klwmin, ldcum, phcbase,&
1883 & ktype, klab, ptu, pqu, plu,&
1884 & puu, pvu, pmfu, pmfub, &
1885 & pmfus, pmfuq, pmful, plude, pdmfup,&
1886 & kcbot, kctop, kctop0, kcum, ztmst,&
1887 & pqsenh, plglac, lndj, wup, wbase, kdpl, pmfude_rate)
1959 integer klev,klon,klevp1,klevm1
1960 real(kind=kind_phys) ptenh(klon,klev), pqenh(klon,klev), &
1961 & puen(klon,klev), pven(klon,klev),&
1962 & pten(klon,klev), pqen(klon,klev),&
1963 & pgeo(klon,klev), pgeoh(klon,klevp1),&
1964 & pap(klon,klev), paph(klon,klevp1),&
1965 & pqsen(klon,klev), pqte(klon,klev),&
1966 & pverv(klon,klev), pqsenh(klon,klev)
1967 real(kind=kind_phys) ptu(klon,klev), pqu(klon,klev),&
1968 & puu(klon,klev), pvu(klon,klev),&
1969 & pmfu(klon,klev), zph(klon),&
1971 & pmfus(klon,klev), pmfuq(klon,klev),&
1972 & plu(klon,klev), plude(klon,klev),&
1973 & pmful(klon,klev), pdmfup(klon,klev)
1974 real(kind=kind_phys) zdmfen(klon), zdmfde(klon),&
1975 & zmfuu(klon), zmfuv(klon),&
1976 & zpbase(klon), zqold(klon)
1977 real(kind=kind_phys) phcbase(klon), zluold(klon)
1978 real(kind=kind_phys) zprecip(klon), zlrain(klon,klev)
1979 real(kind=kind_phys) zbuo(klon,klev), kup(klon,klev)
1980 real(kind=kind_phys) wup(klon)
1981 real(kind=kind_phys) wbase(klon), zodetr(klon,klev)
1982 real(kind=kind_phys) plglac(klon,klev)
1984 real(kind=kind_phys) eta(klon),dz(klon)
1986 integer klwmin(klon), ktype(klon),&
1987 & klab(klon,klev), kcbot(klon),&
1988 & kctop(klon), kctop0(klon)
1990 logical ldcum(klon), loflag(klon)
1991 logical llo2,llo3, llo1(klon)
1994 real(kind=kind_phys) zoentr(klon), zdpmean(klon)
1995 real(kind=kind_phys) pdmfen(klon,klev), pmfude_rate(klon,klev)
1998 integer ikb,icum,itopm2,ik,icall,is,kcum,jlm,jll
2000 real(kind=kind_phys) ztmst,zcons2,zfacbuo,zprcdgw,z_cwdrag,z_cldmax,z_cwifrac,z_cprc2
2001 real(kind=kind_phys) zmftest,zmfmax,zqeen,zseen,zscde,zqude
2002 real(kind=kind_phys) zmfusk,zmfuqk,zmfulk
2003 real(kind=kind_phys) zbc,zbe,zkedke,zmfun,zwu,zprcon,zdt,zcbf,zzco
2004 real(kind=kind_phys) zlcrit,zdfi,zc,zd,zint,zlnew,zvw,zvi,zalfaw,zrold
2005 real(kind=kind_phys) zrnew,zz,zdmfeu,zdmfdu,dp
2006 real(kind=kind_phys) zfac,zbuoc,zdkbuo,zdken,zvv,zarg,zchange,zxe,zxs,zdshrd
2007 real(kind=kind_phys) atop1,atop2,abot
2012 zfacbuo = 0.5/(1.+0.5)
2013 zprcdgw = cprcon*zrg
2017 z_cwdrag = (3.0/8.0)*0.506/0.2
2027 if(.not.ldcum(jl))
then
2038 if(jk.ne.kcbot(jl)) plu(jl,jk)=0.
2050 pmfude_rate(jl,jk) = 0.
2051 if(.not.ldcum(jl).or.ktype(jl).eq.3) klab(jl,jk)=0
2052 if(.not.ldcum(jl).and.paph(jl,jk).lt.4.e4) kctop0(jl)=jk
2057 if ( ktype(jl) == 3 ) ldcum(jl) = .false.
2066 kup(jl,ikb) = 0.5*wbase(jl)**2
2067 pmfu(jl,ikb) = pmfub(jl)
2068 pmfus(jl,ikb) = pmfub(jl)*(cpd*ptu(jl,ikb)+pgeoh(jl,ikb))
2069 pmfuq(jl,ikb) = pmfub(jl)*pqu(jl,ikb)
2070 pmful(jl,ikb) = pmfub(jl)*plu(jl,ikb)
2087 & (klon, klev, klevm1, ik, pten,&
2088 & pqen, pqsen, puen, pven, pverv,&
2089 & pgeo, pgeoh, ldcum, ktype, klab, zlrain,&
2090 & pmfu, pmfub, kcbot, ptu,&
2091 & pqu, plu, puu, pvu, pmfus,&
2092 & pmfuq, pmful, pdmfup)
2096 loflag(jl) = .false.
2099 is = is + klab(jl,jk+1)
2100 if ( klab(jl,jk+1) == 0 ) klab(jl,jk) = 0
2101 if ( (ldcum(jl) .and. klab(jl,jk+1) == 2) .or. &
2102 (ktype(jl) == 3 .and. klab(jl,jk+1) == 1) )
then
2107 zph(jl) = paph(jl,jk)
2108 if ( ktype(jl) == 3 .and. jk == kcbot(jl) )
then
2109 zmfmax = (paph(jl,jk)-paph(jl,jk-1))*zcons2
2110 if ( pmfub(jl) > zmfmax )
then
2111 zfac = zmfmax/pmfub(jl)
2112 pmfu(jl,jk+1) = pmfu(jl,jk+1)*zfac
2113 pmfus(jl,jk+1) = pmfus(jl,jk+1)*zfac
2114 pmfuq(jl,jk+1) = pmfuq(jl,jk+1)*zfac
2117 pmfub(jl)=min(pmfub(jl),zmfmax)
2121 if(is.gt.0) llo3 = .true.
2126 call cuentrn(klon,klev,ik,kcbot,ldcum,llo3, &
2127 pgeoh,pmfu,zdmfen,zdmfde)
2138 zdmfde(jl) = min(zdmfde(jl),0.75*pmfu(jl,jk+1))
2139 if ( jk == kcbot(jl) )
then
2140 zoentr(jl) = -1.75e-3*(min(1.,pqen(jl,jk)/pqsen(jl,jk)) - &
2141 1.)*(pgeoh(jl,jk)-pgeoh(jl,jk+1))*zrg
2142 zoentr(jl) = min(0.4,zoentr(jl))*pmfu(jl,jk+1)
2144 if ( jk < kcbot(jl) )
then
2145 zmfmax = (paph(jl,jk)-paph(jl,jk-1))*zcons2
2146 zxs = max(pmfu(jl,jk+1)-zmfmax,0.)
2147 wup(jl) = wup(jl) + kup(jl,jk+1)*(pap(jl,jk+1)-pap(jl,jk))
2148 zdpmean(jl) = zdpmean(jl) + pap(jl,jk+1) - pap(jl,jk)
2149 zdmfen(jl) = zoentr(jl)
2150 if ( ktype(jl) >= 2 )
then
2151 zdmfen(jl) = 2.0*zdmfen(jl)
2152 zdmfde(jl) = zdmfen(jl)
2154 zdmfde(jl) = zdmfde(jl) * &
2155 (1.6-min(1.,pqen(jl,jk)/pqsen(jl,jk)))
2156 zmftest = pmfu(jl,jk+1) + zdmfen(jl) - zdmfde(jl)
2157 zchange = max(zmftest-zmfmax,0.)
2158 zxe = max(zchange-zxs,0.)
2159 zdmfen(jl) = zdmfen(jl) - zxe
2160 zchange = zchange - zxe
2161 zdmfde(jl) = zdmfde(jl) + zchange
2163 pdmfen(jl,jk) = zdmfen(jl) - zdmfde(jl)
2164 pmfu(jl,jk) = pmfu(jl,jk+1) + zdmfen(jl) - zdmfde(jl)
2165 zqeen = pqenh(jl,jk+1)*zdmfen(jl)
2166 zseen = (cpd*ptenh(jl,jk+1)+pgeoh(jl,jk+1))*zdmfen(jl)
2167 zscde = (cpd*ptu(jl,jk+1)+pgeoh(jl,jk+1))*zdmfde(jl)
2168 zqude = pqu(jl,jk+1)*zdmfde(jl)
2169 plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl)
2170 zmfusk = pmfus(jl,jk+1) + zseen - zscde
2171 zmfuqk = pmfuq(jl,jk+1) + zqeen - zqude
2172 zmfulk = pmful(jl,jk+1) - plude(jl,jk)
2173 plu(jl,jk) = zmfulk*(1./max(cmfcmin,pmfu(jl,jk)))
2174 pqu(jl,jk) = zmfuqk*(1./max(cmfcmin,pmfu(jl,jk)))
2175 ptu(jl,jk) = (zmfusk * &
2176 (1./max(cmfcmin,pmfu(jl,jk)))-pgeoh(jl,jk))*rcpd
2177 ptu(jl,jk) = max(100.,ptu(jl,jk))
2178 ptu(jl,jk) = min(400.,ptu(jl,jk))
2179 zqold(jl) = pqu(jl,jk)
2180 zlrain(jl,jk) = zlrain(jl,jk+1)*(pmfu(jl,jk+1)-zdmfde(jl)) * &
2181 (1./max(cmfcmin,pmfu(jl,jk)))
2182 zluold(jl) = plu(jl,jk)
2186 if ( jk > kdpl(jl) )
then
2187 ptu(jl,jk) = ptenh(jl,jk)
2188 pqu(jl,jk) = pqenh(jl,jk)
2190 zluold(jl) = plu(jl,jk)
2200 call cuadjtqn(klon,klev,ik,zph,ptu,pqu,loflag,icall)
2205 if ( pqu(jl,jk) /= zqold(jl) )
then
2206 plglac(jl,jk) = plu(jl,jk) * &
2207 ((1.-foealfa(ptu(jl,jk)))- &
2208 (1.-foealfa(ptu(jl,jk+1))))
2209 ptu(jl,jk) = ptu(jl,jk) + ralfdcp*plglac(jl,jk)
2214 if ( pqu(jl,jk) /= zqold(jl) )
then
2216 plu(jl,jk) = plu(jl,jk) + zqold(jl) - pqu(jl,jk)
2217 zbc = ptu(jl,jk)*(1.+vtmpc1*pqu(jl,jk)-plu(jl,jk+1) - &
2219 zbe = ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))
2220 zbuo(jl,jk) = zbc - zbe
2222 if ( ktype(jl) == 3 .and. klab(jl,jk+1) == 1 )
then
2223 if ( zbuo(jl,jk) > -0.5 )
then
2234 if ( klab(jl,jk+1) == 2 )
then
2235 if ( zbuo(jl,jk) < 0. )
then
2236 ptenh(jl,jk) = 0.5*(pten(jl,jk)+pten(jl,jk-1))
2237 pqenh(jl,jk) = 0.5*(pqen(jl,jk)+pqen(jl,jk-1))
2238 zbuo(jl,jk) = zbc - ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk))
2240 zbuoc = (zbuo(jl,jk) / &
2241 (ptenh(jl,jk)*(1.+vtmpc1*pqenh(jl,jk)))+zbuo(jl,jk+1) / &
2242 (ptenh(jl,jk+1)*(1.+vtmpc1*pqenh(jl,jk+1))))*0.5
2243 zdkbuo = (pgeoh(jl,jk)-pgeoh(jl,jk+1))*zfacbuo*zbuoc
2245 if ( zdmfen(jl) > 0. )
then
2246 zdken = min(1.,(1.+z_cwdrag)*zdmfen(jl) / &
2247 max(cmfcmin,pmfu(jl,jk+1)))
2249 zdken = min(1.,(1.+z_cwdrag)*zdmfde(jl) / &
2250 max(cmfcmin,pmfu(jl,jk+1)))
2252 kup(jl,jk) = (kup(jl,jk+1)*(1.-zdken)+zdkbuo) / &
2254 if ( zbuo(jl,jk) < 0. )
then
2255 zkedke = kup(jl,jk)/max(1.e-10,kup(jl,jk+1))
2256 zkedke = max(0.,min(1.,zkedke))
2257 zmfun = sqrt(zkedke)*pmfu(jl,jk+1)
2259 zdmfde(jl) = max(zdmfde(jl),pmfu(jl,jk+1)-zmfun)
2260 plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl)
2261 pmfu(jl,jk) = pmfu(jl,jk+1) + zdmfen(jl) - zdmfde(jl)
2263 if ( zbuo(jl,jk) > -0.2 )
then
2265 zoentr(jl) = 1.75e-3*(0.3-(min(1.,pqen(jl,jk-1) / &
2266 pqsen(jl,jk-1))-1.))*(pgeoh(jl,jk-1)-pgeoh(jl,jk)) * &
2267 zrg*min(1.,pqsen(jl,jk)/pqsen(jl,ikb))**3
2268 zoentr(jl) = min(0.4,zoentr(jl))*pmfu(jl,jk)
2273 if ( jk > kdpl(jl) )
then
2274 pmfu(jl,jk) = pmfu(jl,jk+1)
2277 if ( kup(jl,jk) > 0. .and. pmfu(jl,jk) > 0. )
then
2284 zdmfde(jl) = pmfu(jl,jk+1)
2285 plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl)
2288 if ( pmfu(jl,jk+1) > 0. ) pmfude_rate(jl,jk) = zdmfde(jl)
2290 else if ( ktype(jl) == 2 .and. pqu(jl,jk) == zqold(jl) )
then
2294 zdmfde(jl) = pmfu(jl,jk+1)
2295 plude(jl,jk) = plu(jl,jk+1)*zdmfde(jl)
2296 pmfude_rate(jl,jk) = zdmfde(jl)
2301 if ( llo1(jl) )
then
2305 if ( lndj(jl).eq.1 )
then
2311 if ( plu(jl,jk) > zdshrd )
then
2313 zwu = min(15.0,sqrt(2.*max(0.1,kup(jl,jk+1))))
2314 zprcon = zprcdgw/(0.75*zwu)
2316 zdt = min(rtber-rtice,max(rtber-ptu(jl,jk),0.))
2317 zcbf = 1. + z_cprc2*sqrt(zdt)
2319 zlcrit = zdshrd/zcbf
2320 zdfi = pgeoh(jl,jk) - pgeoh(jl,jk+1)
2321 zc = (plu(jl,jk)-zluold(jl))
2322 zarg = (plu(jl,jk)/zlcrit)**2
2323 if ( zarg < 25.0 )
then
2324 zd = zzco*(1.-exp(-zarg))*zdfi
2329 zlnew = zluold(jl)*zint + zc/zd*(1.-zint)
2330 zlnew = max(0.,min(plu(jl,jk),zlnew))
2331 zlnew = min(z_cldmax,zlnew)
2332 zprecip(jl) = max(0.,zluold(jl)+zc-zlnew)
2333 pdmfup(jl,jk) = zprecip(jl)*pmfu(jl,jk)
2334 zlrain(jl,jk) = zlrain(jl,jk) + zprecip(jl)
2340 if ( llo1(jl) )
then
2341 if ( zlrain(jl,jk) > 0. )
then
2342 zvw = 21.18*zlrain(jl,jk)**0.2
2344 zalfaw = foealfa(ptu(jl,jk))
2345 zvv = zalfaw*zvw + (1.-zalfaw)*zvi
2346 zrold = zlrain(jl,jk) - zprecip(jl)
2348 zwu = min(15.0,sqrt(2.*max(0.1,kup(jl,jk))))
2351 zrnew = zrold*zint + zc/zd*(1.-zint)
2352 zrnew = max(0.,min(zlrain(jl,jk),zrnew))
2353 zlrain(jl,jk) = zrnew
2359 pmful(jl,jk) = plu(jl,jk)*pmfu(jl,jk)
2360 pmfus(jl,jk) = (cpd*ptu(jl,jk)+pgeoh(jl,jk))*pmfu(jl,jk)
2361 pmfuq(jl,jk) = pqu(jl,jk)*pmfu(jl,jk)
2369 if ( kctop(jl) == -1 ) ldcum(jl) = .false.
2370 kcbot(jl) = max(kcbot(jl),kctop(jl))
2371 if ( ldcum(jl) )
then
2372 wup(jl) = max(1.e-2,wup(jl)/max(1.,zdpmean(jl)))
2373 wup(jl) = sqrt(2.*wup(jl))
2839 & ( klon, klev, ztmst &
2840 & , pten, pqen, pqsen, ptenh, pqenh &
2841 & , paph, pap, pgeoh, lndj, ldcum &
2842 & , kcbot, kctop, kdtop, ktopm2 &
2844 & , pmfu, pmfd, pmfus, pmfds &
2845 & , pmfuq, pmfdq, pmful, plude &
2846 & , pdmfup, pdmfdp, pdpmel, plglac &
2847 & , prain, pmfdde_rate, pmflxr, pmflxs )
2925 integer klev,klon,ktopm2
2926 real(kind=kind_phys) pten(klon,klev), ptenh(klon,klev), &
2927 & pqen(klon,klev), pqsen(klon,klev), &
2928 & pqenh(klon,klev), pap(klon,klev), &
2929 & paph(klon,klev+1), pgeoh(klon,klev+1), &
2932 real(kind=kind_phys) pmfu(klon,klev), pmfd(klon,klev), &
2933 & pmfus(klon,klev), pmfds(klon,klev), &
2934 & pmfuq(klon,klev), pmfdq(klon,klev), &
2935 & pdmfup(klon,klev), pdmfdp(klon,klev), &
2936 & pdpmel(klon,klev), prain(klon), &
2937 & pmful(klon,klev), plude(klon,klev), &
2938 & pmflxr(klon,klev+1), pmflxs(klon,klev+1)
2939 real(kind=kind_phys) pmfdde_rate(klon,klev)
2940 integer kcbot(klon), kctop(klon), &
2941 & kdtop(klon), ktype(klon)
2942 logical lddraf(klon), &
2947 integer is,ik,icall,ike,ikb
2948 real(kind=kind_phys) ztmst,ztaumel,zcons1a,zcons1,zcons2,zcucov,zcpecons
2949 real(kind=kind_phys) zalfaw,zrfl,zdrfl1,zrnew,zrmin,zrfln,zdrfl,zdenom
2950 real(kind=kind_phys) zpdr,zpds,zzp,zfac,zsnmlt
2951 real(kind=kind_phys) rhevap(klon)
2958 zcons1a=cpd/(alf*g*ztaumel)
2967 if(.not.ldcum(jl).or.kdtop(jl).lt.kctop(jl)) lddraf(jl)=.false.
2968 if(.not.ldcum(jl)) ktype(jl)=0
2970 if(lndj(jl) .eq. 1)
then
2978 ikb = min(jk+1,klev)
2983 if(ldcum(jl).and.jk.ge.kctop(jl))
then
2984 pmfus(jl,jk)=pmfus(jl,jk)-pmfu(jl,jk)* &
2985 & (cpd*ptenh(jl,jk)+pgeoh(jl,jk))
2986 pmfuq(jl,jk)=pmfuq(jl,jk)-pmfu(jl,jk)*pqenh(jl,jk)
2987 plglac(jl,jk)=pmfu(jl,jk)*plglac(jl,jk)
2988 llddraf = lddraf(jl) .and. jk >= kdtop(jl)
2989 if ( llddraf .and.jk.ge.kdtop(jl))
then
2990 pmfds(jl,jk) = pmfds(jl,jk)-pmfd(jl,jk) * &
2991 (cpd*ptenh(jl,jk)+pgeoh(jl,jk))
2992 pmfdq(jl,jk) = pmfdq(jl,jk)-pmfd(jl,jk)*pqenh(jl,jk)
2997 pdmfdp(jl,jk-1) = 0.
2999 if ( llddraf .and. pmfd(jl,jk) < 0. .and. &
3000 abs(pmfd(jl,ikb)) < 1.e-20 )
then
3020 pmflxr(jl,klev+1) = 0.
3021 pmflxs(jl,klev+1) = 0.
3027 zzp=((paph(jl,klev+1)-paph(jl,ik))/ &
3028 & (paph(jl,klev+1)-paph(jl,ikb)))
3029 if(ktype(jl).eq.3)
then
3032 pmfu(jl,ik)=pmfu(jl,ikb)*zzp
3033 pmfus(jl,ik)=(pmfus(jl,ikb)- &
3034 & foelhm(ptenh(jl,ikb))*pmful(jl,ikb))*zzp
3035 pmfuq(jl,ik)=(pmfuq(jl,ikb)+pmful(jl,ikb))*zzp
3042 if(ldcum(jl).and.jk.gt.kcbot(jl)+1)
then
3044 zzp=((paph(jl,klev+1)-paph(jl,jk))/ &
3045 & (paph(jl,klev+1)-paph(jl,ikb)))
3046 if(ktype(jl).eq.3)
then
3049 pmfu(jl,jk)=pmfu(jl,ikb)*zzp
3050 pmfus(jl,jk)=pmfus(jl,ikb)*zzp
3051 pmfuq(jl,jk)=pmfuq(jl,ikb)*zzp
3055 llddraf = lddraf(jl) .and. jk > ik .and. ik < klev
3056 if ( llddraf .and. ik == kcbot(jl)+1 )
then
3057 zzp = ((paph(jl,klev+1)-paph(jl,jk))/(paph(jl,klev+1)-paph(jl,ik)))
3058 if ( ktype(jl) == 3 ) zzp = zzp*zzp
3059 pmfd(jl,jk) = pmfd(jl,ik)*zzp
3060 pmfds(jl,jk) = pmfds(jl,ik)*zzp
3061 pmfdq(jl,jk) = pmfdq(jl,ik)*zzp
3062 pmfdde_rate(jl,jk) = -(pmfd(jl,jk-1)-pmfd(jl,jk))
3063 else if ( llddraf .and. ik /= kcbot(jl)+1 .and. jk == ik+1 )
then
3064 pmfdde_rate(jl,jk) = -(pmfd(jl,jk-1)-pmfd(jl,jk))
3075 if(ldcum(jl) .and. jk >=kctop(jl)-1 )
then
3076 prain(jl)=prain(jl)+pdmfup(jl,jk)
3077 if(pmflxs(jl,jk).gt.0..and.pten(jl,jk).gt.tmelt)
then
3078 zcons1=zcons1a*(1.+0.5*(pten(jl,jk)-tmelt))
3079 zfac=zcons1*(paph(jl,jk+1)-paph(jl,jk))
3080 zsnmlt=min(pmflxs(jl,jk),zfac*(pten(jl,jk)-tmelt))
3081 pdpmel(jl,jk)=zsnmlt
3082 pqsen(jl,jk)=foeewm(pten(jl,jk)-zsnmlt/zfac)/pap(jl,jk)
3084 zalfaw=foealfa(pten(jl,jk))
3088 if ( pten(jl,jk) < tmelt .and. zalfaw > 0. )
then
3089 plglac(jl,jk) = plglac(jl,jk)+zalfaw*(pdmfup(jl,jk)+pdmfdp(jl,jk))
3092 pmflxr(jl,jk+1)=pmflxr(jl,jk)+zalfaw* &
3093 & (pdmfup(jl,jk)+pdmfdp(jl,jk))+pdpmel(jl,jk)
3094 pmflxs(jl,jk+1)=pmflxs(jl,jk)+(1.-zalfaw)* &
3095 & (pdmfup(jl,jk)+pdmfdp(jl,jk))-pdpmel(jl,jk)
3096 if(pmflxr(jl,jk+1)+pmflxs(jl,jk+1).lt.0.0)
then
3097 pdmfdp(jl,jk)=-(pmflxr(jl,jk)+pmflxs(jl,jk)+pdmfup(jl,jk))
3101 else if ( pmflxr(jl,jk+1) < 0. )
then
3102 pmflxs(jl,jk+1) = pmflxs(jl,jk+1)+pmflxr(jl,jk+1)
3103 pmflxr(jl,jk+1) = 0.
3104 else if ( pmflxs(jl,jk+1) < 0. )
then
3105 pmflxr(jl,jk+1) = pmflxr(jl,jk+1)+pmflxs(jl,jk+1)
3106 pmflxs(jl,jk+1) = 0.
3113 if(ldcum(jl).and.jk.ge.kcbot(jl))
then
3114 zrfl=pmflxr(jl,jk)+pmflxs(jl,jk)
3115 if(zrfl.gt.1.e-20)
then
3116 zdrfl1=zcpecons*max(0.,pqsen(jl,jk)-pqen(jl,jk))*zcucov* &
3117 & (sqrt(paph(jl,jk)/paph(jl,klev+1))/5.09e-3* &
3118 & zrfl/zcucov)**0.5777* &
3119 & (paph(jl,jk+1)-paph(jl,jk))
3121 zrmin=zrfl-zcucov*max(0.,rhevap(jl)*pqsen(jl,jk) &
3122 & -pqen(jl,jk)) *zcons2*(paph(jl,jk+1)-paph(jl,jk))
3123 zrnew=max(zrnew,zrmin)
3125 zdrfl=min(0.,zrfln-zrfl)
3126 zdenom=1./max(1.e-20,pmflxr(jl,jk)+pmflxs(jl,jk))
3127 zalfaw=foealfa(pten(jl,jk))
3128 if ( pten(jl,jk) < tmelt ) zalfaw = 0.
3129 zpdr=zalfaw*pdmfdp(jl,jk)
3130 zpds=(1.-zalfaw)*pdmfdp(jl,jk)
3131 pmflxr(jl,jk+1)=pmflxr(jl,jk)+zpdr &
3132 & +pdpmel(jl,jk)+zdrfl*pmflxr(jl,jk)*zdenom
3133 pmflxs(jl,jk+1)=pmflxs(jl,jk)+zpds &
3134 & -pdpmel(jl,jk)+zdrfl*pmflxs(jl,jk)*zdenom
3135 pdmfup(jl,jk)=pdmfup(jl,jk)+zdrfl
3136 if ( pmflxr(jl,jk+1)+pmflxs(jl,jk+1) < 0. )
then
3137 pdmfup(jl,jk) = pdmfup(jl,jk)-(pmflxr(jl,jk+1)+pmflxs(jl,jk+1))
3138 pmflxr(jl,jk+1) = 0.
3139 pmflxs(jl,jk+1) = 0.
3141 else if ( pmflxr(jl,jk+1) < 0. )
then
3142 pmflxs(jl,jk+1) = pmflxs(jl,jk+1)+pmflxr(jl,jk+1)
3143 pmflxr(jl,jk+1) = 0.
3144 else if ( pmflxs(jl,jk+1) < 0. )
then
3145 pmflxr(jl,jk+1) = pmflxr(jl,jk+1)+pmflxs(jl,jk+1)
3146 pmflxs(jl,jk+1) = 0.