Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

implementation of new orographic drag parameterization scheme to alleviate wind bias in E3SM #6667

Open
wants to merge 3 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 11 additions & 0 deletions components/eam/bld/build-namelist
Original file line number Diff line number Diff line change
Expand Up @@ -4089,13 +4089,24 @@ if ($waccm_phys or $cfg->get('nlev') >= 60) {
add_default($nl, 'use_gw_oro' , 'val'=>'.true.');
add_default($nl, 'use_gw_front' , 'val'=>'.true.');
add_default($nl, 'use_gw_convect', 'val'=>'.true.');
add_default($nl, 'use_od_ls', 'val'=>'.false.');
add_default($nl, 'use_od_bl', 'val'=>'.false.');
add_default($nl, 'use_od_ss', 'val'=>'.false.');
add_default($nl, 'use_od_fd', 'val'=>'.false.');
} else {
add_default($nl, 'use_gw_oro' , 'val'=>'.true.');
add_default($nl, 'use_gw_front' , 'val'=>'.false.');
add_default($nl, 'use_gw_convect', 'val'=>'.false.');
add_default($nl, 'use_od_ls', 'val'=>'.false.');
add_default($nl, 'use_od_bl', 'val'=>'.false.');
add_default($nl, 'use_od_ss', 'val'=>'.false.');
add_default($nl, 'use_od_fd', 'val'=>'.false.');
}
add_default($nl, 'pgwv', 'val'=>'32');
add_default($nl, 'gw_dc','val'=>'2.5D0');
add_default($nl, 'ncleff_ls', 'val'=>'3.D0');
add_default($nl, 'ncd_bl', 'val'=>'3.D0');
add_default($nl, 'sncleff_ss','val'=>'1.D0');

if ($nl->get_value('use_gw_oro') =~ /$TRUE/io) {
add_default($nl, 'effgw_oro');
Expand Down
10 changes: 6 additions & 4 deletions components/eam/bld/namelist_files/namelist_defaults_eam.xml
Original file line number Diff line number Diff line change
Expand Up @@ -126,14 +126,13 @@

<!-- Topography -->
<bnd_topo hgrid="64x128" >atm/cam/topo/USGS-gtopo30_64x128_c050520.nc</bnd_topo>

<bnd_topo hgrid="ne4np4" npg="0">atm/cam/topo/USGS-gtopo30_ne4np4_16x.c20160612.nc</bnd_topo>
<bnd_topo hgrid="ne4np4" npg="2">atm/cam/topo/USGS-gtopo30_ne4np4pg2_16x_converted.c20200527.nc</bnd_topo>
<bnd_topo hgrid="ne4np4" npg="0">atm/cam/topo/USGS-gtopo30_ne4np4_16x_forOroDrag.c20241019.nc</bnd_topo>
<bnd_topo hgrid="ne4np4" npg="2">atm/cam/topo/USGS-gtopo30_ne4np4pg2_16x_converted_forOroDrag.c20241019.nc</bnd_topo>
<bnd_topo hgrid="ne11np4" npg="0">atm/cam/topo/USGS-gtopo30_ne11np4_16xconsistentSGH.c20160612.nc</bnd_topo>
<bnd_topo hgrid="ne16np4" npg="0">atm/cam/topo/USGS-gtopo30_ne16np4_16xconsistentSGH.c20160612.nc</bnd_topo>
<bnd_topo hgrid="ne16np4" npg="2">atm/cam/topo/USGS-gtopo30_ne16np4pg2_16xdel2_20200527.nc</bnd_topo>
<bnd_topo hgrid="ne30np4" npg="0">atm/cam/topo/USGS-gtopo30_ne30np4_16xdel2-PFC-consistentSGH.nc</bnd_topo>
<bnd_topo hgrid="ne30np4" npg="2">atm/cam/topo/USGS-gtopo30_ne30np4pg2_x6t-SGH.c20210614.nc</bnd_topo>
<bnd_topo hgrid="ne30np4" npg="2">atm/cam/topo/USGS-gtopo30_ne30np4pg2_x6t-SGH_forOroDrag.c20241001.nc</bnd_topo>
<bnd_topo hgrid="ne30np4" npg="3">atm/cam/topo/USGS-gtopo30_ne30np4pg3_16xdel2.c20200504.nc</bnd_topo>
<bnd_topo hgrid="ne30np4" npg="4">atm/cam/topo/USGS-gtopo30_ne30np4pg4_16xdel2.c20200504.nc</bnd_topo>
<bnd_topo hgrid="ne45np4" npg="2">atm/cam/topo/USGS-gtopo30_ne45np4pg2_16xdel2.c20200615.nc</bnd_topo>
Expand Down Expand Up @@ -1884,6 +1883,9 @@ with se_tstep, dt_remap_factor, dt_tracer_factor set to -1
<hdepth_scaling_factor phys="default" use_MMF="1" > 1.0 </hdepth_scaling_factor>
<effgw_oro phys="default"> 0.375 </effgw_oro>
<use_gw_energy_fix phys="default"> .true. </use_gw_energy_fix>
<ncleff_ls phys="default"> 3 </ncleff_ls>
<ncd_bl phys="default"> 3 </ncd_bl>
<sncleff_ss phys="default"> 1 </sncleff_ss>
<clubb_C14 phys="default"> 2.5D0 </clubb_C14>
<clubb_tk1 phys="default"> 268.15D0 </clubb_tk1>
<dust_emis_fact phys="default"> 13.8D0 </dust_emis_fact>
Expand Down
42 changes: 42 additions & 0 deletions components/eam/bld/namelist_files/namelist_definition.xml
Original file line number Diff line number Diff line change
Expand Up @@ -1078,6 +1078,48 @@ Whether or not to enable GWD brute-force energy fix.
Default: set by build-namelist.
</entry>

<entry id="use_od_ls" type="logical" category="gw_drag"
group="phys_ctl_nl" valid_values="" >
Whether or not to enable nonlinear orographic gravity wave drag (oGWD).
Default: set by build-namelist.
</entry>

<entry id="use_od_bl" type="logical" category="gw_drag"
group="phys_ctl_nl" valid_values="" >
Whether or not to enable flow-blocking drag (FBD).
Default: set by build-namelist.
</entry>

<entry id="use_od_ss" type="logical" category="gw_drag"
group="phys_ctl_nl" valid_values="" >
Whether or not to enable small-scale orographic GWD drag (sGWD).
Default: set by build-namelist.
</entry>

<entry id="use_od_fd" type="logical" category="pbl"
group="phys_ctl_nl" valid_values="" >
Whether or not to enable turbulent orographic form drag (TOFD).
Default: set by build-namelist.
</entry>

<entry id="ncleff_ls" type="real" category="gw_drag"
group="oro_drag_nl" valid_values="" >
Tuning parameter of orographic GWD (oGWD). See <varname>use_od_ls</varname>.
Default: set by build-namelist.
</entry>

<entry id="ncd_bl" type="real" category="gw_drag"
group="oro_drag_nl" valid_values="" >
Tuning parameter of flow-blocking drag (FBD). See <varname>use_od_bl</varname>.
Default: set by build-namelist.
</entry>

<entry id="sncleff_ss" type="real" category="gw_drag"
group="oro_drag_nl" valid_values="" >
Tuning parameter of small-scale GWD (sGWD). See <varname>use_od_ss</varname>.
Default: set by build-namelist.
</entry>

<entry id="pgwv" type="integer" category="gw_drag"
group="gw_drag_nl" valid_values="" >
Gravity wave spectrum dimension (wave numbers are from -pgwv to pgwv).
Expand Down
41 changes: 41 additions & 0 deletions components/eam/src/control/startup_initialconds.F90
Original file line number Diff line number Diff line change
Expand Up @@ -5,16 +5,28 @@ module startup_initialconds
!
!-----------------------------------------------------------------------

use pio, only: file_desc_t

implicit none
private
save

public :: initial_conds ! Read in initial conditions (dycore dependent)
!added for orographic drag
public topoGWD_file_get_id
public setup_initialGWD
public close_initial_fileGWD
type(file_desc_t), pointer :: ncid_topoGWD

!=======================================================================
contains
!=======================================================================

function topoGWD_file_get_id()
type(file_desc_t), pointer :: topoGWD_file_get_id
topoGWD_file_get_id => ncid_topoGWD
end function topoGWD_file_get_id

subroutine initial_conds(dyn_in)

! This routine does some initializing of buffers that should move to a
Expand Down Expand Up @@ -62,4 +74,33 @@ end subroutine initial_conds

!=======================================================================

subroutine setup_initialGWD()
use filenames, only: bnd_topo
use ioFileMod, only: getfil
use cam_pio_utils, only: cam_pio_openfile
use pio, only: pio_nowrite
!
! Input arguments
!
!-----------------------------------------------------------------------
include 'netcdf.inc'
!-----------------------------------------------------------------------
character(len=256) :: bnd_topo_loc ! filepath of topo file on local disk
allocate(ncid_topoGWD)
call getfil(bnd_topo, bnd_topo_loc)
call cam_pio_openfile(ncid_topoGWD, bnd_topo_loc, PIO_NOWRITE)
end subroutine setup_initialGWD

subroutine close_initial_fileGWD
use pio, only: pio_closefile
call pio_closefile(ncid_topoGWD)
deallocate(ncid_topoGWD)
nullify(ncid_topoGWD)
end subroutine close_initial_fileGWD
!=======================================================================





end module startup_initialconds
125 changes: 115 additions & 10 deletions components/eam/src/physics/cam/clubb_intr.F90
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ module clubb_intr
use shr_kind_mod, only: r8=>shr_kind_r8
use shr_log_mod , only: errMsg => shr_log_errMsg
use ppgrid, only: pver, pverp
use phys_control, only: phys_getopts
use phys_control, only: phys_getopts,use_od_ss,use_od_fd,ncleff_ls,ncd_bl,sncleff_ss
use physconst, only: rair, cpair, gravit, latvap, latice, zvir, rh2o, karman, &
tms_orocnst, tms_z0fac, pi
use cam_logfile, only: iulog
Expand Down Expand Up @@ -927,7 +927,15 @@ subroutine clubb_ini_cam(pbuf2d, dp1_in)
call addfld ('VMAGDP', horiz_only, 'A', '-', 'ZM gustiness enhancement')
call addfld ('VMAGCL', horiz_only, 'A', '-', 'CLUBB gustiness enhancement')
call addfld ('TPERTBLT', horiz_only, 'A', 'K', 'perturbation temperature at PBL top')

!!added for TOFD output
call addfld ('DTAUX3_FD',(/'lev'/),'A','m/s2','U tendency - fd orographic drag')
call addfld ('DTAUY3_FD',(/'lev'/),'A','m/s2','V tendency - fd orographic drag')
call addfld ('DUSFC_FD',horiz_only,'A','N/m2','fd zonal oro surface stress')
call addfld ('DVSFC_FD',horiz_only,'A','N/m2','fd merio oro surface stress')
call add_default('DTAUX3_FD', 1, ' ')
call add_default('DTAUY3_FD', 1, ' ')
call add_default('DUSFC_FD', 1, ' ')
call add_default('DVSFC_FD', 1, ' ')
! Initialize statistics, below are dummy variables
dum1 = 300._r8
dum2 = 1200._r8
Expand Down Expand Up @@ -1155,7 +1163,10 @@ subroutine clubb_tend_cam( &
use model_flags, only: ipdf_call_placement
use advance_clubb_core_module, only: ipdf_post_advance_fields
#endif

use gw_common, only: grid_size,gw_oro_interface
use hycoef, only: etamid
use physconst, only: rh2o,pi,rearth,r_universal
!!get the znu,znw,p_top set to 0
implicit none

! --------------- !
Expand Down Expand Up @@ -1518,8 +1529,30 @@ subroutine clubb_tend_cam( &

real(r8) :: sfc_v_diff_tau(pcols) ! Response to tau perturbation, m/s
real(r8), parameter :: pert_tau = 0.1_r8 ! tau perturbation, Pa


!add par for tofd
real(r8) :: dtaux3_fd(pcols,pver)
real(r8) :: dtauy3_fd(pcols,pver)
real(r8) :: dusfc_fd(pcols)
real(r8) :: dvsfc_fd(pcols)
logical :: gwd_ls,gwd_bl,gwd_ss,gwd_fd
real(r8) :: dummy_nm(pcols,pver)
real(r8) :: dummy_utgw(pcols,pver)
real(r8) :: dummy_vtgw(pcols,pver)
real(r8) :: dummy_ttgw(pcols,pver)
!
real(r8) :: dummx_ls(pcols,pver)
real(r8) :: dummx_bl(pcols,pver)
real(r8) :: dummx_ss(pcols,pver)
real(r8) :: dummy_ls(pcols,pver)
real(r8) :: dummy_bl(pcols,pver)
real(r8) :: dummy_ss(pcols,pver)
real(r8) :: dummx3_ls(pcols,pver)
real(r8) :: dummx3_bl(pcols,pver)
real(r8) :: dummx3_ss(pcols,pver)
real(r8) :: dummy3_ls(pcols,pver)
real(r8) :: dummy3_bl(pcols,pver)
real(r8) :: dummy3_ss(pcols,pver)
!
real(r8) :: inv_exner_clubb_surf


Expand Down Expand Up @@ -1946,7 +1979,36 @@ subroutine clubb_tend_cam( &
tautmsx, tautmsy, cam_in%landfrac )
call t_stopf('compute_tms')
endif

!
if (use_od_fd) then
gwd_ls=.false.
gwd_bl=.false.
gwd_ss=.false.
gwd_fd=use_od_fd
dummy_nm=0.0_r8
dummy_utgw=0.0_r8
dummy_vtgw=0.0_r8
dummy_ttgw=0.0_r8
!sgh30 as the input for TOFD instead of sgh
call gw_oro_interface(state,cam_in,sgh30,pbuf,hdtime,dummy_nm,&
gwd_ls,gwd_bl,gwd_ss,gwd_fd,&
ncleff_ls,ncd_bl,sncleff_ss,&
dummy_utgw,dummy_vtgw,dummy_ttgw,&
dtaux3_ls=dummx3_ls,dtauy3_ls=dummy3_ls,&
dtaux3_bl=dummx3_bl,dtauy3_bl=dummy3_bl,&
dtaux3_ss=dummx3_ss,dtauy3_ss=dummy3_ss,&
dtaux3_fd=dtaux3_fd,dtauy3_fd=dtauy3_fd,&
dusfc_ls=dummx_ls,dvsfc_ls=dummy_ls,&
dusfc_bl=dummx_bl,dvsfc_bl=dummy_bl,&
dusfc_ss=dummx_ss,dvsfc_ss=dummy_ss,&
dusfc_fd=dusfc_fd,dvsfc_fd=dvsfc_fd)
!
call outfld ('DTAUX3_FD', dtaux3_fd, pcols, lchnk)
call outfld ('DTAUY3_FD', dtauy3_fd, pcols, lchnk)
call outfld ('DUSFC_FD', dusfc_fd, pcols, lchnk)
call outfld ('DVSFC_FD', dvsfc_fd, pcols, lchnk)
endif
!
if (micro_do_icesupersat) then
call physics_ptend_init(ptend_loc,state%psetcols, 'clubb_ice3', ls=.true., lu=.true., lv=.true., lq=lq)
endif
Expand Down Expand Up @@ -2067,7 +2129,14 @@ subroutine clubb_tend_cam( &
dum_core_rknd = real((ksrftms(i)*state1%v(i,pver)), kind = core_rknd)
vpwp_sfc = vpwp_sfc-(dum_core_rknd/rho_ds_zm(1))
endif

!----------------------------------------------------!
!Apply TOFD
!----------------------------------------------------!
!tendency is flipped already
if (use_od_fd) then
um_forcing(2:pverp)=dtaux3_fd(i,pver:1:-1)
vm_forcing(2:pverp)=dtauy3_fd(i,pver:1:-1)
endif
! Need to flip arrays around for CLUBB core
do k=1,pverp
um_in(k) = real(um(i,pverp-k+1), kind = core_rknd)
Expand Down Expand Up @@ -3108,10 +3177,12 @@ subroutine clubb_surface (state, cam_in, ustar, obklen)
!-------------------------------------------------------------------------------

use physics_types, only: physics_state
use physconst, only: zvir
use physconst, only: zvir,gravit
use ppgrid, only: pver, pcols
use constituents, only: cnst_get_ind
use camsrfexch, only: cam_in_t
use hb_diff, only: pblintd_ri


implicit none

Expand All @@ -3136,13 +3207,17 @@ subroutine clubb_surface (state, cam_in, ustar, obklen)
! --------------- !

integer :: i ! indicees
integer :: k
integer :: ncol ! # of atmospheric columns

real(r8) :: th(pcols) ! surface potential temperature
real(r8) :: thv(pcols) ! surface virtual potential temperature
real(r8) :: th_lv(pcols,pver) ! level potential temperature
real(r8) :: thv_lv(pcols,pver) ! level virtual potential temperature
real(r8) :: kinheat ! kinematic surface heat flux
real(r8) :: kinwat ! kinematic surface vapor flux
real(r8) :: kbfs ! kinematic surface buoyancy flux
real(r8) :: kbfs_pcol(pcols)
integer :: ixq,ixcldliq !PMA fix for thv
real(r8) :: rrho ! Inverse air density

Expand Down Expand Up @@ -3173,14 +3248,44 @@ subroutine clubb_surface (state, cam_in, ustar, obklen)
thv(i) = th(i)*(1._r8+zvir*state%q(i,pver,ixq)) ! diagnose virtual potential temperature
end if
enddo

!
do i = 1, ncol
call calc_ustar( state%t(i,pver), state%pmid(i,pver), cam_in%wsx(i), cam_in%wsy(i), &
rrho, ustar(i) )
call calc_obklen( th(i), thv(i), cam_in%cflx(i,1), cam_in%shf(i), rrho, ustar(i), &
kinheat, kinwat, kbfs, obklen(i) )
enddo

!
if (use_od_ss) then
!add calculation of bulk richardson number here
!
!compute the whole level th and thv for diagnose of bulk richardson number
thv_lv=0.0_r8
th_lv=0.0_r8
!
do i=1,ncol
do k=1,pver
th_lv(i,k) = state%t(i,k)*state%exner(i,k)
if (use_sgv) then
thv_lv(i,k) = th_lv(i,k)*(1.0_r8+zvir*state%q(i,k,ixq) &
- state%q(i,k,ixcldliq)) !PMA corrects thv formula
else
thv_lv(i,k) = th_lv(i,k)*(1.0_r8+zvir*state%q(i,k,ixq))
end if
enddo
enddo
!
kbfs_pcol=0.0_r8
do i=1,ncol
call calc_obklen( th(i), thv(i), cam_in%cflx(i,1), cam_in%shf(i), rrho, ustar(i), &
kinheat, kinwat, kbfs, obklen(i) )
kbfs_pcol(i)=kbfs
enddo
!
call pblintd_ri(ncol, gravit, thv_lv, state%zm, state%u, state%v, &
ustar, obklen, kbfs_pcol, state%ribulk)
endif
!
return

#endif
Expand Down
Loading