Skip to content

Commit

Permalink
Merge pull request #487 from climbfuji/iacono_cloud_overlap_merged_wi…
Browse files Browse the repository at this point in the history
…th_rrtmgp_and_updated

RRTMG cloud overlap method update (contains #477)
  • Loading branch information
climbfuji committed Sep 25, 2020
2 parents 0808cc2 + a8d1b6b commit f91d1bf
Show file tree
Hide file tree
Showing 13 changed files with 707 additions and 162 deletions.
133 changes: 85 additions & 48 deletions physics/GFS_rrtmg_pre.F90
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
faerlw1, faerlw2, faerlw3, aerodp, &
clouds1, clouds2, clouds3, clouds4, clouds5, clouds6, &
clouds7, clouds8, clouds9, cldsa, cldfra, &
mtopa, mbota, de_lgth, alb1d, errmsg, errflg)
mtopa, mbota, de_lgth, alpha, alb1d, errmsg, errflg)

use machine, only: kind_phys
use GFS_typedefs, only: GFS_statein_type, &
Expand Down Expand Up @@ -128,6 +128,7 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
real(kind=kind_phys), dimension(size(Grid%xlon,1),5), intent(out) :: cldsa
integer, dimension(size(Grid%xlon,1),3), intent(out) :: mbota, mtopa
real(kind=kind_phys), dimension(size(Grid%xlon,1)), intent(out) :: de_lgth, alb1d
real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP), intent(out) :: alpha

character(len=*), intent(out) :: errmsg
integer, intent(out) :: errflg
Expand All @@ -137,22 +138,21 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input

integer :: i, j, k, k1, k2, lsk, lv, n, itop, ibtc, LP1, lla, llb, lya, lyb

real(kind=kind_phys) :: es, qs, delt, tem0d
real(kind=kind_phys) :: es, qs, delt, tem0d, pfac

real(kind=kind_phys), dimension(size(Grid%xlon,1)) :: cvt1, cvb1, tem1d, tskn

real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP) :: &
htswc, htlwc, gcice, grain, grime, htsw0, htlw0, &
rhly, tvly,qstl, vvel, clw, ciw, prslk1, tem2da, &
cldcov, deltaq, cnvc, cnvw, &
dzb, hzb, cldcov, deltaq, cnvc, cnvw, &
effrl, effri, effrr, effrs, rho, orho
! for Thompson MP
real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP) :: &
re_cloud, re_ice, re_snow, qv_mp, qc_mp, &
qi_mp, qs_mp, nc_mp, ni_mp, nwfa

real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP+1) :: tem2db
! real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP+1) :: hz
real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP+1) :: tem2db, hz

real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP,min(4,Model%ncnd)) :: ccnd
real(kind=kind_phys), dimension(size(Grid%xlon,1),lm+LTP,2:Model%ntrac) :: tracer1
Expand Down Expand Up @@ -424,17 +424,31 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
enddo

! --- ... level height and layer thickness (km)
! dz: Layer thickness between layer boundaries
! dzb: Layer thickness between layer centers (lowest is from surface to lowest layer center)
! hz: Height of each level (i.e. layer boundary)
! hzb: Height of each layer center

tem0d = 0.001 * rog
do i = 1, IM
do k = 1, LMK
dz(i,k) = tem0d * (tem2db(i,k+1) - tem2db(i,k)) * tvly(i,k)
enddo

! hz(i,LMP) = 0.0
! do k = LMK, 1, -1
! hz(i,k) = hz(i,k+1) + dz(i,k)
! enddo
hz(i,LMP) = 0.0
do k = LMK, 1, -1
hz(i,k) = hz(i,k+1) + dz(i,k)
enddo

do k = LMK, 1, -1
pfac = (tem2db(i,k+1) - tem2da(i,k)) / (tem2db(i,k+1) - tem2db(i,k))
hzb(i,k) = hz(i,k+1) + pfac * (hz(i,k) - hz(i,k+1))
enddo

do k = LMK-1, 1, -1
dzb(i,k) = hzb(i,k) - hzb(i,k+1)
enddo
dzb(i,LMK) = hzb(i,LMK) - hz(i,LMP)
enddo

else ! input data from sfc to toa
Expand Down Expand Up @@ -475,17 +489,31 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
enddo

! --- ... level height and layer thickness (km)
! dz: Layer thickness between layer boundaries
! dzb: Layer thickness between layer centers (lowest is from surface to lowest layer center)
! hz: Height of each level (i.e. layer boundary)
! hzb: Height of each layer center

tem0d = 0.001 * rog
do i = 1, IM
do k = LMK, 1, -1
dz(i,k) = tem0d * (tem2db(i,k) - tem2db(i,k+1)) * tvly(i,k)
enddo

! hz(i,1) = 0.0
! do k = 1, LMP
! hz(i,k+1) = hz(i,k) + dz(i,k)
! enddo
hz(i,1) = 0.0
do k = 1, LMK
hz(i,k+1) = hz(i,k) + dz(i,k)
enddo

do k = 1, LMK
pfac = (tem2db(i,k) - tem2da(i,k)) / (tem2db(i,k) - tem2db(i,k+1))
hzb(i,k) = hz(i,k) + pfac * (hz(i,k+1) - hz(i,k))
enddo

do k = 2, LMK
dzb(i,k) = hzb(i,k) - hzb(i,k-1)
enddo
dzb(i,1) = hzb(i,1) - hz(i,1)
enddo

endif ! end_if_ivflip
Expand Down Expand Up @@ -807,19 +835,21 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
! or unified cloud and/or with MG microphysics

if (Model%uni_cld .and. Model%ncld >= 2) then
call progclduni (plyr, plvl, tlyr, tvly, ccnd, ncndl, & ! --- inputs
Grid%xlat, Grid%xlon, Sfcprop%slmsk,dz,delp, &
IM, LMK, LMP, cldcov, &
effrl, effri, effrr, effrs, Model%effr_in, &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
call progclduni (plyr, plvl, tlyr, tvly, ccnd, ncndl, & ! --- inputs
Grid%xlat, Grid%xlon, Sfcprop%slmsk,dz,delp, &
IM, LMK, LMP, cldcov, &
effrl, effri, effrr, effrs, Model%effr_in, &
dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs
else
call progcld1 (plyr ,plvl, tlyr, tvly, qlyr, qstl, rhly, & ! --- inputs
ccnd(1:IM,1:LMK,1), Grid%xlat,Grid%xlon, &
Sfcprop%slmsk, dz, delp, IM, LMK, LMP, &
Model%uni_cld, Model%lmfshal, &
Model%lmfdeep2, cldcov, &
effrl, effri, effrr, effrs, Model%effr_in, &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
call progcld1 (plyr ,plvl, tlyr, tvly, qlyr, qstl, rhly, & ! --- inputs
ccnd(1:IM,1:LMK,1), Grid%xlat,Grid%xlon, &
Sfcprop%slmsk, dz, delp, IM, LMK, LMP, &
Model%uni_cld, Model%lmfshal, &
Model%lmfdeep2, cldcov, &
effrl, effri, effrr, effrs, Model%effr_in, &
dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs
endif

elseif(Model%imp_physics == 98) then ! zhao/moorthi's prognostic cloud+pdfcld
Expand All @@ -829,7 +859,8 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
cnvw, cnvc, Grid%xlat, Grid%xlon, &
Sfcprop%slmsk, dz, delp, im, lmk, lmp, deltaq, &
Model%sup, Model%kdt, me, &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs


elseif (Model%imp_physics == 11) then ! GFDL cloud scheme
Expand All @@ -839,21 +870,24 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
ccnd(1:IM,1:LMK,1), cnvw, cnvc, &
Grid%xlat, Grid%xlon, Sfcprop%slmsk, &
cldcov, dz, delp, im, lmk, lmp, &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs
else

call progclduni (plyr, plvl, tlyr, tvly, ccnd, ncndl, & ! --- inputs
Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz,delp, &
IM, LMK, LMP, cldcov, &
effrl, effri, effrr, effrs, Model%effr_in, &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
dzb, Grid%xlat_d, Model%julian, Model%yearlen,&
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs
! call progcld4o (plyr, plvl, tlyr, tvly, qlyr, qstl, rhly, & ! --- inputs
! tracer1, Grid%xlat, Grid%xlon, Sfcprop%slmsk, &
! dz, delp, &
! ntrac-1, Model%ntcw-1,Model%ntiw-1,Model%ntrw-1,&
! Model%ntsw-1,Model%ntgl-1,Model%ntclamt-1, &
! im, lmk, lmp, &
! clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
! dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
! clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs
endif

elseif(Model%imp_physics == 6 .or. Model%imp_physics == 15) then
Expand All @@ -863,15 +897,16 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
Tbd%phy_f3d(:,:,Model%nseffr) = 250.
endif

call progcld5 (plyr,plvl,tlyr,qlyr,qstl,rhly,tracer1, & ! --- inputs
Grid%xlat,Grid%xlon,Sfcprop%slmsk,dz,delp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lmk, lmp, Model%uni_cld, &
Model%lmfshal,Model%lmfdeep2, &
cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), &
Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), &
clouds,cldsa,mtopa,mbota, de_lgth) ! --- outputs
call progcld5 (plyr,plvl,tlyr,qlyr,qstl,rhly,tracer1, & ! --- inputs
Grid%xlat,Grid%xlon,Sfcprop%slmsk,dz,delp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lmk, lmp, Model%uni_cld, &
Model%lmfshal,Model%lmfdeep2, &
cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), &
Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), &
dzb, Grid%xlat_d, Model%julian, Model%yearlen,&
clouds,cldsa,mtopa,mbota, de_lgth, alpha) ! --- outputs


elseif(Model%imp_physics == Model%imp_physics_thompson) then ! Thompson MP
Expand All @@ -892,19 +927,21 @@ subroutine GFS_rrtmg_pre_run (Model, Grid, Sfcprop, Statein, & ! input
Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz,delp, &
IM, LMK, LMP, clouds(:,1:LMK,1), &
effrl, effri, effrr, effrs, Model%effr_in , &
clouds, cldsa, mtopa, mbota, de_lgth) ! --- outputs
dzb, Grid%xlat_d, Model%julian, Model%yearlen, &
clouds, cldsa, mtopa, mbota, de_lgth, alpha) ! --- outputs

else
! MYNN PBL or GF convective are not used
call progcld5 (plyr,plvl,tlyr,qlyr,qstl,rhly,tracer1, & ! --- inputs
Grid%xlat,Grid%xlon,Sfcprop%slmsk,dz,delp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lmk, lmp, Model%uni_cld, &
Model%lmfshal,Model%lmfdeep2, &
cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), &
Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), &
clouds,cldsa,mtopa,mbota, de_lgth) ! --- outputs
call progcld5 (plyr,plvl,tlyr,qlyr,qstl,rhly,tracer1, & ! --- inputs
Grid%xlat,Grid%xlon,Sfcprop%slmsk,dz,delp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lmk, lmp, Model%uni_cld, &
Model%lmfshal,Model%lmfdeep2, &
cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), &
Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), &
dzb, Grid%xlat_d, Model%julian, Model%yearlen,&
clouds,cldsa,mtopa,mbota, de_lgth, alpha) ! --- outputs
endif ! MYNN PBL or GF

endif ! end if_imp_physics
Expand Down
9 changes: 9 additions & 0 deletions physics/GFS_rrtmg_pre.meta
Original file line number Diff line number Diff line change
Expand Up @@ -564,6 +564,15 @@
kind = kind_phys
intent = out
optional = F
[alpha]
standard_name = cloud_overlap_decorrelation_parameter
long_name = cloud overlap decorrelation parameter
units = frac
dimensions = (horizontal_dimension,adjusted_vertical_layer_dimension_for_radiation)
type = real
kind = kind_phys
intent = out
optional = F
[alb1d]
standard_name = surface_albedo_perturbation
long_name = surface albedo perturbation
Expand Down
4 changes: 2 additions & 2 deletions physics/GFS_rrtmg_setup.meta
Original file line number Diff line number Diff line change
Expand Up @@ -113,15 +113,15 @@
intent = in
optional = F
[iovr_sw]
standard_name = flag_for_max_random_overlap_clouds_for_shortwave_radiation
standard_name = flag_for_cloud_overlap_method_for_shortwave_radiation
long_name = sw: max-random overlap clouds
units = flag
dimensions = ()
type = integer
intent = in
optional = F
[iovr_lw]
standard_name = flag_for_max_random_overlap_clouds_for_longwave_radiation
standard_name = flag_for_cloud_overlap_method_for_longwave_radiation
long_name = lw: max-random overlap clouds
units = flag
dimensions = ()
Expand Down
4 changes: 2 additions & 2 deletions physics/GFS_rrtmgp_setup.meta
Original file line number Diff line number Diff line change
Expand Up @@ -161,15 +161,15 @@
intent = in
optional = F
[iovr_sw]
standard_name = flag_for_max_random_overlap_clouds_for_shortwave_radiation
standard_name = flag_for_cloud_overlap_method_for_shortwave_radiation
long_name = sw: max-random overlap clouds
units = flag
dimensions = ()
type = integer
intent = in
optional = F
[iovr_lw]
standard_name = flag_for_max_random_overlap_clouds_for_longwave_radiation
standard_name = flag_for_cloud_overlap_method_for_longwave_radiation
long_name = lw: max-random overlap clouds
units = flag
dimensions = ()
Expand Down
2 changes: 1 addition & 1 deletion physics/GFS_rrtmgp_zhaocarr_pre.F90
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ module GFS_rrtmgp_zhaocarr_pre
use machine, only: kind_phys
use rrtmgp_aux, only: check_error_msg
use funcphys, only: fpvs
use module_radiation_clouds, only: get_alpha_dcorr
use module_radiation_clouds, only: get_alpha_dcorr

! Zhao-Carr MP parameters.
real(kind_phys), parameter :: &
Expand Down
34 changes: 31 additions & 3 deletions physics/module_SGSCloud_RadPre.F90
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,7 @@ subroutine sgscloud_radpre_run( &
nlay, plyr, xlat, dz,de_lgth, &
cldsa,mtopa,mbota, &
imp_physics, imp_physics_gfdl,&
iovr, &
errmsg, errflg )

! should be moved to inside the mynn:
Expand Down Expand Up @@ -81,6 +82,7 @@ subroutine sgscloud_radpre_run( &
real(kind=kind_phys), dimension(im,nlay), intent(in) :: plyr, dz
real(kind=kind_phys), dimension(im,5), intent(inout) :: cldsa
integer, dimension(im,3), intent(inout) :: mbota, mtopa
integer, intent(in) :: iovr
character(len=*), intent(out) :: errmsg
integer, intent(out) :: errflg
! Local variables
Expand All @@ -93,6 +95,9 @@ subroutine sgscloud_radpre_run( &
real(kind=kind_phys), dimension(im) :: rxlat
real (kind=kind_phys):: Tc, iwc
integer :: i, k, id
! DH* 20200723 - see comment at the end of this routine around 'gethml'
real(kind=kind_phys), dimension(im,nlay) :: alpha_dummy
! *DH

! PARAMETERS FOR RANDALL AND XU (1996) CLOUD FRACTION
REAL, PARAMETER :: coef_p = 0.25, coef_gamm = 0.49, coef_alph = 100.
Expand Down Expand Up @@ -123,7 +128,7 @@ subroutine sgscloud_radpre_run( &

if (h2oliq > clwt) then
onemrh= max( 1.e-10, 1.0-rhgrid )
tem1 = min(max((onemrh*qsat)**0.49,0.0001),1.0) !jhan
tem1 = min(max((onemrh*qsat)**0.49,0.0001),1.0) !jhan
tem1 = 100.0 / tem1
value = max( min( tem1*(h2oliq-clwt), 50.0 ), 0.0 )
tem2 = sqrt( sqrt(rhgrid) )
Expand Down Expand Up @@ -304,12 +309,35 @@ subroutine sgscloud_radpre_run( &

cldcnv = 0.

! DH* 20200723
! iovr == 4 or 5 requires alpha, which is computed in GFS_rrmtg_pre,
! which comes after SGSCloud_RadPre. Computing alpha here requires
! a lot more input variables and computations (dzlay etc.), and
! recomputing it in GFS_rrmtg_pre is a waste of time. Workaround:
! pass a dummy array initialized to zero to gethml for other values of iovr.
if ( iovr == 4 .or. iovr == 5 ) then
errmsg = 'Logic error in sgscloud_radpre: iovr==4 or 5 not implemented'
errflg = 1
return
end if
!! Call subroutine get_alpha_exp to define alpha parameter for EXP and ER cloud overlap options
! if ( iovr == 4 .or. iovr == 5 ) then
! call get_alpha_exp &
!! --- inputs:
! (im, nlay, dzlay, iovr, latdeg, julian, yearlen, clouds1, &
!! --- outputs:
! alpha &
! )
! endif
alpha_dummy = 0.0
! *DH 2020723

!> - Recompute the diagnostic high, mid, low, total and bl cloud fraction
call gethml &
! --- inputs:
( plyr, ptop1, clouds1, cldcnv, dz, de_lgth, im, nlay, &
( plyr, ptop1, clouds1, cldcnv, dz, de_lgth, alpha_dummy, &
! --- outputs:
cldsa, mtopa, mbota)
im, nlay, cldsa, mtopa, mbota)

!print*,"===Finished adding subgrid clouds to the resolved-scale clouds"
!print*,"qc_save:",qc_save(1,1)," qi_save:",qi_save(1,1)
Expand Down
8 changes: 8 additions & 0 deletions physics/module_SGSCloud_RadPre.meta
Original file line number Diff line number Diff line change
Expand Up @@ -337,6 +337,14 @@
type = integer
intent = in
optional = F
[iovr]
standard_name = flag_for_cloud_overlap_method_for_radiation
long_name = max-random overlap clouds
units = flag
dimensions = ()
type = integer
intent = in
optional = F
[errmsg]
standard_name = ccpp_error_message
long_name = error message for error handling in CCPP
Expand Down
Loading

0 comments on commit f91d1bf

Please sign in to comment.