diff --git a/physics/CONV/SAMF/samfdeepcnv.f b/physics/CONV/SAMF/samfdeepcnv.f index 3ad067657..05b801c7e 100644 --- a/physics/CONV/SAMF/samfdeepcnv.f +++ b/physics/CONV/SAMF/samfdeepcnv.f @@ -76,7 +76,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & & tmf,qmicro,itc,ntc,cliq,cp,cvap, & & eps,epsm1,fv,grav,hvap,rd,rv, & & t0c,delt,ntk,ntr,delp, & - & prslp,psp,phil,qtr,prevsq,q,q1,t1,u1,v1,fscav, & + & prslp,psp,phil,tkeh,qtr,prevsq,q,q1,t1,u1,v1,fscav, & & hwrf_samfdeep,progsigma,cldwrk,rn,kbot,ktop,kcnv, & & islimsk,garea,dot,ncloud,hpbl,ud_mf,dd_mf,dt_mf,cnvw,cnvc, & & QLCN, QICN, w_upi, cf_upi, CNV_MFD, & @@ -97,7 +97,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & & fv, grav, hvap, rd, rv, t0c real(kind=kind_phys), intent(in) :: delt real(kind=kind_phys), intent(in) :: psp(:), delp(:,:), & - & prslp(:,:), garea(:), hpbl(:), dot(:,:), phil(:,:) + & prslp(:,:), garea(:), hpbl(:), dot(:,:), phil(:,:), tkeh(:,:) real(kind=kind_phys), dimension(:), intent(in) :: fscav logical, intent(in) :: first_time_step,restart,hwrf_samfdeep, & & progsigma,do_mynnedmf @@ -953,8 +953,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & if(cnvflg(i)) then if(k >= kb(i) .and. k < kbcon(i)) then dz = zo(i,k+1) - zo(i,k) - tem = 0.5 * (qtr(i,k,ntk)+qtr(i,k+1,ntk)) - tkemean(i) = tkemean(i) + tem * dz + tkemean(i) = tkemean(i) + tkeh(i,k) * dz sumx(i) = sumx(i) + dz endif endif @@ -1286,6 +1285,22 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & enddo enddo enddo + kk = ntk -2 + do k = 2, km1 + do i = 1, im + if (cnvflg(i)) then + if(k > kb(i) .and. k < kmax(i)) then + dz = zi(i,k) - zi(i,k-1) + tem = 0.25 * (xlamue(i,k)+xlamue(i,k-1)) * dz + tem = cq * tem + factor = 1. + tem + ecko(i,k,kk) = ((1. - tem) * ecko(i,k-1,kk) + tem * + & (ctro(i,k,kk) + ctro(i,k-1,kk))) / factor + ercko(i,k,kk) = ecko(i,k,kk) + endif + endif + enddo + enddo endif endif c @@ -2941,7 +2956,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, & flag_mid = .false. call progsigma_calc(im,km,first_time_step,restart,flag_shallow, & flag_mid,del,tmfq,qmicro,dbyo1,zdqca,omega_u,zeta,hvap, - & delt,qadv,kbcon1,ktcon,cnvflg,betascu,betamcu,betadcu, + & delt,qadv,kb,kbcon1,ktcon,cnvflg,betascu,betamcu,betadcu, & sigmind,sigminm,sigmins,sigmain,sigmaout,sigmab) endif diff --git a/physics/CONV/SAMF/samfdeepcnv.meta b/physics/CONV/SAMF/samfdeepcnv.meta index 3652a0d27..756c95ba5 100644 --- a/physics/CONV/SAMF/samfdeepcnv.meta +++ b/physics/CONV/SAMF/samfdeepcnv.meta @@ -242,6 +242,14 @@ type = real kind = kind_phys intent = in +[tkeh] + standard_name = vertical_turbulent_kinetic_energy_at_interface + long_name = vertical turbulent kinetic energy at model layer interfaces + units = m2 s-2 + dimensions = (horizontal_loop_extent,vertical_layer_dimension) + type = real + kind = kind_phys + intent = in [qtr] standard_name = convective_transportable_tracers long_name = array to contain cloud water and other convective trans. tracers diff --git a/physics/CONV/SAMF/samfshalcnv.f b/physics/CONV/SAMF/samfshalcnv.f index ce783ea15..1d3b3dbd2 100644 --- a/physics/CONV/SAMF/samfshalcnv.f +++ b/physics/CONV/SAMF/samfshalcnv.f @@ -53,7 +53,7 @@ subroutine samfshalcnv_run(im,km,itc,ntc,cliq,cp,cvap, & & eps,epsm1,fv,grav,hvap,rd,rv, & & t0c,delt,ntk,ntr,delp,first_time_step,restart, & & tmf,qmicro,progsigma, & - & prslp,psp,phil,qtr,prevsq,q,q1,t1,u1,v1,fscav, & + & prslp,psp,phil,tkeh,qtr,prevsq,q,q1,t1,u1,v1,fscav, & & rn,kbot,ktop,kcnv,islimsk,garea, & & dot,ncloud,hpbl,ud_mf,dt_mf,cnvw,cnvc, & & clam,c0s,c1,evef,pgcon,asolfac,hwrf_samfshal, & @@ -71,7 +71,7 @@ subroutine samfshalcnv_run(im,km,itc,ntc,cliq,cp,cvap, & & betamcu real(kind=kind_phys), intent(in) :: delt real(kind=kind_phys), intent(in) :: psp(:), delp(:,:), & - & prslp(:,:), garea(:), hpbl(:), dot(:,:), phil(:,:), & + & prslp(:,:), garea(:), hpbl(:), dot(:,:), phil(:,:), tkeh(:,:), & & tmf(:,:,:), q(:,:) real(kind=kind_phys), intent(in), optional :: qmicro(:,:), & & prevsq(:,:) @@ -872,14 +872,13 @@ subroutine samfshalcnv_run(im,km,itc,ntc,cliq,cp,cvap, & tkemean(i) = 0. endif enddo - +! do k = 1, km1 do i = 1, im if(cnvflg(i)) then if(k >= kb(i) .and. k < kbcon(i)) then dz = zo(i,k+1) - zo(i,k) - tem = 0.5 * (qtr(i,k,ntk)+qtr(i,k+1,ntk)) - tkemean(i) = tkemean(i) + tem * dz + tkemean(i) = tkemean(i) + tkeh(i,k) * dz sumx(i) = sumx(i) + dz endif endif @@ -1093,6 +1092,22 @@ subroutine samfshalcnv_run(im,km,itc,ntc,cliq,cp,cvap, & enddo enddo enddo + kk = ntk -2 + do k = 2, km1 + do i = 1, im + if (cnvflg(i)) then + if(k > kb(i) .and. k < kmax(i)) then + dz = zi(i,k) - zi(i,k-1) + tem = 0.25 * (xlamue(i,k)+xlamue(i,k-1)) * dz + tem = cq * tem + factor = 1. + tem + ecko(i,k,kk) = ((1. - tem) * ecko(i,k-1,kk) + tem * + & (ctro(i,k,kk) + ctro(i,k-1,kk))) / factor + ercko(i,k,kk) = ecko(i,k,kk) + endif + endif + enddo + enddo endif endif c @@ -1982,7 +1997,7 @@ subroutine samfshalcnv_run(im,km,itc,ntc,cliq,cp,cvap, & flag_mid = .false. call progsigma_calc(im,km,first_time_step,restart,flag_shallow, & flag_mid,del,tmfq,qmicro,dbyo1,zdqca,omega_u,zeta,hvap, - & delt,qadv,kbcon1,ktcon,cnvflg,betascu,betamcu,betadcu, + & delt,qadv,kb,kbcon1,ktcon,cnvflg,betascu,betamcu,betadcu, & sigmind,sigminm,sigmins,sigmain,sigmaout,sigmab) endif diff --git a/physics/CONV/SAMF/samfshalcnv.meta b/physics/CONV/SAMF/samfshalcnv.meta index 4dfa8ac20..fcf964e2b 100644 --- a/physics/CONV/SAMF/samfshalcnv.meta +++ b/physics/CONV/SAMF/samfshalcnv.meta @@ -242,6 +242,14 @@ type = real kind = kind_phys intent = in +[tkeh] + standard_name = vertical_turbulent_kinetic_energy_at_interface + long_name = vertical turbulent kinetic energy at model layer interfaces + units = m2 s-2 + dimensions = (horizontal_loop_extent,vertical_layer_dimension) + type = real + kind = kind_phys + intent = in [qtr] standard_name = convective_transportable_tracers long_name = array to contain cloud water and other convective trans. tracers