Skip to content

Commit

Permalink
addressing some review comments
Browse files Browse the repository at this point in the history
  • Loading branch information
lisa-bengtsson committed May 5, 2022
1 parent 8b815e0 commit be534e7
Show file tree
Hide file tree
Showing 3 changed files with 80 additions and 58 deletions.
114 changes: 59 additions & 55 deletions physics/progsigma_calc.f90
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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.
Expand All @@ -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.
Expand All @@ -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.
Expand All @@ -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).
Expand Down Expand Up @@ -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
Expand All @@ -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
Expand Down
16 changes: 13 additions & 3 deletions physics/samfdeepcnv.f
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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

Expand Down Expand Up @@ -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)
Expand Down
8 changes: 8 additions & 0 deletions physics/samfdeepcnv.meta
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down

0 comments on commit be534e7

Please sign in to comment.