diff --git a/physics/progsigma_calc.f90 b/physics/progsigma_calc.f90 index c05af3003..fe74dc0c1 100644 --- a/physics/progsigma_calc.f90 +++ b/physics/progsigma_calc.f90 @@ -48,7 +48,7 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & real(kind=kind_phys) :: gcvalmx,epsilon,ZZ,cvg,mcon,buy2, & fdqb,dtdyn,dxlim,rmulacvg,tem, & - alpha,DEN,betascu,dp1 + alpha,DEN,betascu,dp1,invdelt !Parameters gcvalmx = 0.1 @@ -57,11 +57,11 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & km1=km-1 alpha=7000. betascu = 3.0 + invdelt = 1./delt !Initialization 2D do k = 1,km do i = 1,im - sigmaout(i,k)=0. inbu(i,k)=0. form(i,k)=0. dp(i,k)=0. @@ -70,7 +70,9 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & !Initialization 1D do i=1,im - sigmab(i)=0. + if(cnvflg(i))then + sigmab(i)=0. + endif sigmamax(i)=0.95 termA(i)=0. termB(i)=0. @@ -89,24 +91,16 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & enddo !Initial computations, place maximum sigmain in sigmab - if(flag_init .and. .not. flag_restart)then - do i=1,im - if(cnvflg(i))then - sigmab(i)=0.03 - endif - enddo - else - do i=1,im - if(cnvflg(i))then - do k=2,km - if(sigmain(i,k)>sigmab(i))then - sigmab(i)=sigmain(i,k) - endif - enddo - endif - enddo - endif - + do i=1,im + if(cnvflg(i))then + do k=2,km + if(sigmain(i,k)>sigmab(i))then + sigmab(i)=sigmain(i,k) + endif + enddo + endif + enddo + do i=1,im if(sigmab(i) < 1.E-5)then !after advection sigmab(i)=0. @@ -120,15 +114,20 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & enddo !Initial computations, dynamic q-tendency - do k = 1,km - do i = 1,im - if(flag_init .and. .not.flag_restart)then + if(flag_init .and. .not.flag_restart)then + do k = 1,km + do i = 1,im qadv(i,k)=0. - else - qadv(i,k)=(q(i,k) - prevsq(i,k))/delt - endif + enddo enddo - enddo + else + do k = 1,km + do i = 1,im + qadv(i,k)=(q(i,k) - prevsq(i,k))*invdelt + enddo + enddo + endif + !compute termD "The vertical integral of the latent heat convergence is limited to the !buoyant layers with positive moisture convergence (accumulated from the surface). @@ -173,7 +172,7 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & endif enddo enddo - + !termC do k = 2,km1 do i = 1,im @@ -188,33 +187,38 @@ subroutine progsigma_calc (im,km,flag_init,flag_restart, & enddo !sigmab - do i = 1,im - if(cnvflg(i))then - DEN=MIN(termC(i)+termB(i),1.E8) - cvg=termD(i)*delt - ZZ=MAX(0.0,SIGN(1.0,termA(i))) & - *MAX(0.0,SIGN(1.0,termB(i))) & - *MAX(0.0,SIGN(1.0,termC(i)-epsilon)) - cvg=MAX(0.0,cvg) - if(flag_init .and. .not. flag_restart)then - sigmab(i)=0.03 - else - sigmab(i)=(ZZ*(termA(i)+cvg))/(DEN+(1.0-ZZ)) - endif - if(sigmab(i)>0.)then - sigmab(i)=MIN(sigmab(i),sigmamax(i)) - sigmab(i)=MAX(sigmab(i),0.01) - endif - endif!cnvflg - enddo + if(flag_init .and. .not. flag_restart)then + do i = 1,im + if(cnvflg(i))then + sigmab(i)=0.03 + endif + enddo + else + do i = 1,im + if(cnvflg(i))then + DEN=MIN(termC(i)+termB(i),1.E8) + cvg=termD(i)*delt + ZZ=MAX(0.0,SIGN(1.0,termA(i))) & + *MAX(0.0,SIGN(1.0,termB(i))) & + *MAX(0.0,SIGN(1.0,termC(i)-epsilon)) + cvg=MAX(0.0,cvg) + sigmab(i)=(ZZ*(termA(i)+cvg))/(DEN+(1.0-ZZ)) + if(sigmab(i)>0.)then + sigmab(i)=MIN(sigmab(i),sigmamax(i)) + sigmab(i)=MAX(sigmab(i),0.01) + endif + endif!cnvflg + enddo + endif - do k=1,km - do i=1,im - if(cnvflg(i))then - sigmaout(i,k)=sigmab(i) - endif - enddo - enddo + do k=1,km + do i=1,im + if(cnvflg(i))then + sigmaout(i,k)=sigmab(i) + endif + enddo + enddo + !Since updraft velocity is much lower in shallow cu region, termC becomes small in shallow cu application, thus the area fraction !in this regime becomes too large compared with the deep cu region. To address this simply apply a scaling factor for shallow cu diff --git a/physics/samfdeepcnv.f b/physics/samfdeepcnv.f index 071bf0557..8398af769 100644 --- a/physics/samfdeepcnv.f +++ b/physics/samfdeepcnv.f @@ -85,8 +85,8 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & & QLCN, QICN, w_upi, cf_upi, CNV_MFD, & & CNV_DQLDT,CLCN,CNV_FICE,CNV_NDROP,CNV_NICE,mp_phys,mp_phys_mg,& & clam,c0s,c1,betal,betas,evef,pgcon,asolfac, & - & do_ca, ca_closure, ca_entr, ca_trigger, nthresh, ca_deep, & - & rainevap, sigmain, sigmaout, errmsg,errflg) + & do_ca, ca_closure, ca_entr, ca_trigger, nthresh,ca_deep, & + & ca_micro, rainevap, sigmain, sigmaout, errmsg,errflg) ! use machine , only : kind_phys use funcphys , only : fpvs @@ -107,7 +107,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & real(kind=kind_phys), intent(in) :: ca_deep(:) real(kind=kind_phys), intent(in) :: sigmain(:,:),qmicro(:,:), & & tmf(:,:),q(:,:), prevsq(:,:) - real(kind=kind_phys), intent(out) :: rainevap(:) + real(kind=kind_phys), intent(out) :: rainevap(:), ca_micro(:) real(kind=kind_phys), intent(out) :: sigmaout(:,:) logical, intent(in) :: do_ca,ca_closure,ca_entr,ca_trigger @@ -2894,6 +2894,16 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & !> - From Han et al.'s (2017) \cite han_et_al_2017 equation 6, calculate cloud base mass flux as a function of the mean updraft velcoity for the grid sizes where the quasi-equilibrium assumption of Arakawa-Schubert is not valid any longer. !! As discussed in Han et al. (2017) \cite han_et_al_2017 , when dtconv is larger than tauadv, the convective mixing is not fully conducted before the cumulus cloud is advected out of the grid cell. In this case, therefore, the cloud base mass flux is further reduced in proportion to the ratio of tauadv to dtconv. + do i=1,im + ca_micro(i)=0. + enddo + + do i=1,im + if(cnvflg(i))then + ca_micro(i)=sigmab(i) + endif + enddo + do i= 1, im if(cnvflg(i) .and. .not.asqecflg(i)) then k = kbcon(i) diff --git a/physics/samfdeepcnv.meta b/physics/samfdeepcnv.meta index 3f28035b6..1764a74fd 100644 --- a/physics/samfdeepcnv.meta +++ b/physics/samfdeepcnv.meta @@ -652,6 +652,14 @@ type = real kind = kind_phys intent = in +[ca_micro] + standard_name = output_prognostic_sigma_two + long_name = output of prognostic area fraction two + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = out [rainevap] standard_name = physics_field_for_coupling long_name = physics_field_for_coupling