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

use sfc_diag to diagnose 2m t/q inside NoahMP #43

Open
wants to merge 1 commit into
base: p8c
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
76 changes: 69 additions & 7 deletions physics/module_sf_noahmplsm.f90
Original file line number Diff line number Diff line change
Expand Up @@ -414,7 +414,7 @@ subroutine noahmp_sflx (parameters, &
pblhx , iz0tlnd , itime ,psi_opt ,&
prcpconv, prcpnonc, prcpshcv, prcpsnow, prcpgrpl, prcphail, & ! in : forcing
tbot , co2air , o2air , foln , ficeold , zlvl , & ! in : forcing
ep_1 , ep_2 , cp , & ! in : constants
ep_1 , ep_2 , epsm1 , cp , & ! in : constants
albold , sneqvo , & ! in/out :
stc , sh2o , smc , tah , eah , fwet , & ! in/out :
canliq , canice , tv , tg , qsfc, qsnow, qrain, & ! in/out :
Expand Down Expand Up @@ -462,6 +462,7 @@ subroutine noahmp_sflx (parameters, &
integer , intent(in) :: jloc !< grid index
real (kind=kind_phys) , intent(in) :: ep_1
real (kind=kind_phys) , intent(in) :: ep_2
real (kind=kind_phys) , intent(in) :: epsm1
real (kind=kind_phys) , intent(in) :: cp
real (kind=kind_phys) , intent(in) :: dt !< time step [sec]
real (kind=kind_phys), dimension( 1:nsoil), intent(in) :: zsoil !< layer-bottom depth from soil surf (m)
Expand Down Expand Up @@ -816,7 +817,7 @@ subroutine noahmp_sflx (parameters, &
fveg ,shdfac, pahv ,pahg ,pahb , & !in
qsnow ,dzsnso ,lat ,canliq ,canice ,iloc, jloc , & !in
thsfc_loc, prslkix,prsik1x,prslk1x,garea1, & !in
pblhx ,iz0tlnd, itime ,psi_opt, ep_1, ep_2, cp, &
pblhx ,iz0tlnd, itime ,psi_opt, ep_1, ep_2, epsm1, cp, &
z0wrf ,z0hwrf , & !out
imelt ,snicev ,snliqv ,epore ,t2m ,fsno , & !out
sav ,sag ,qmelt ,fsa ,fsr ,taux , & !out
Expand Down Expand Up @@ -1655,7 +1656,7 @@ subroutine energy (parameters,ice ,vegtyp ,ist ,nsnow ,nsoil , & !in
fveg ,shdfac, pahv ,pahg ,pahb , & !in
qsnow ,dzsnso ,lat ,canliq ,canice ,iloc , jloc, & !in
thsfc_loc, prslkix,prsik1x,prslk1x,garea1, & !in
pblhx , iz0tlnd, itime,psi_opt,ep_1, ep_2, cp, &
pblhx , iz0tlnd, itime,psi_opt,ep_1, ep_2, epsm1, cp, &
z0wrf ,z0hwrf , & !out
imelt ,snicev ,snliqv ,epore ,t2m ,fsno , & !out
sav ,sag ,qmelt ,fsa ,fsr ,taux , & !out
Expand Down Expand Up @@ -1739,6 +1740,7 @@ subroutine energy (parameters,ice ,vegtyp ,ist ,nsnow ,nsoil , & !in
real (kind=kind_phys) , intent(in) :: pblhx ! pbl height
real (kind=kind_phys) , intent(in) :: ep_1
real (kind=kind_phys) , intent(in) :: ep_2
real (kind=kind_phys) , intent(in) :: epsm1
real (kind=kind_phys) , intent(in) :: cp
integer , intent(in) :: iz0tlnd
integer , intent(in) :: itime
Expand Down Expand Up @@ -2202,7 +2204,7 @@ subroutine energy (parameters,ice ,vegtyp ,ist ,nsnow ,nsoil , & !in
foln ,co2air ,o2air ,btran ,sfcprs , & !in
rhsur ,iloc ,jloc ,q2 ,pahv ,pahg , & !in
thsfc_loc, prslkix,prsik1x,prslk1x, garea1, & !in
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, cp, &
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, epsm1,cp, &
eah ,tah ,tv ,tgv ,cmv, ustarx , & !inout
#ifdef CCPP
chv ,dx ,dz8w ,errmsg ,errflg , & !inout
Expand Down Expand Up @@ -2239,7 +2241,7 @@ subroutine energy (parameters,ice ,vegtyp ,ist ,nsnow ,nsoil , & !in
emg ,stc ,df ,rsurf ,latheag , & !in
gammag ,rhsur ,iloc ,jloc ,q2 ,pahb , & !in
thsfc_loc, prslkix,prsik1x,prslk1x,vegtyp,fveg,shdfac,garea1, & !in
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, cp, &
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, epsm1, cp, &
#ifdef CCPP
tgb ,cmb ,chb, ustarx,errmsg ,errflg , & !inout
#else
Expand Down Expand Up @@ -3644,7 +3646,7 @@ subroutine vege_flux(parameters,nsnow ,nsoil ,isnow ,vegtyp ,veg , &
foln ,co2air ,o2air ,btran ,sfcprs , & !in
rhsur ,iloc ,jloc ,q2 ,pahv ,pahg , & !in
thsfc_loc, prslkix,prsik1x,prslk1x, garea1, & !in
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, cp, &
pblhx ,iz0tlnd ,itime ,psi_opt ,ep_1, ep_2, epsm1, cp, &
eah ,tah ,tv ,tg ,cm,ustarx,& !inout
#ifdef CCPP
ch ,dx ,dz8w ,errmsg ,errflg , & !inout
Expand All @@ -3666,6 +3668,7 @@ subroutine vege_flux(parameters,nsnow ,nsoil ,isnow ,vegtyp ,veg , &
! -sav + irc[tv] + shc[tv] + evc[tv] + tr[tv] + canhs[tv] = 0
! -sag + irg[tg] + shg[tg] + evg[tg] + gh[tg] = 0
! --------------------------------------------------------------------------------------------------
use funcphys, only : fpvs
implicit none
! --------------------------------------------------------------------------------------------------
! input
Expand Down Expand Up @@ -3695,6 +3698,7 @@ subroutine vege_flux(parameters,nsnow ,nsoil ,isnow ,vegtyp ,veg , &
real (kind=kind_phys) , intent(in) :: pblhx ! pbl height
real (kind=kind_phys) , intent(in) :: ep_1
real (kind=kind_phys) , intent(in) :: ep_2
real (kind=kind_phys) , intent(in) :: epsm1
real (kind=kind_phys) , intent(in) :: cp
integer , intent(in) :: iz0tlnd
integer , intent(in) :: itime
Expand Down Expand Up @@ -3906,8 +3910,12 @@ subroutine vege_flux(parameters,nsnow ,nsoil ,isnow ,vegtyp ,veg , &

real (kind=kind_phys) :: t, tdc !kelvin to degree celsius with limit -50 to +50

real(kind=kind_phys) :: fhi, qss, wrk
real(kind=kind_phys), parameter :: qmin=1.0e-8

character(len=80) :: message


tdc(t) = min( 50., max(-50.,(t-tfrz)) )
! ---------------------------------------------------------------------------------------------

Expand Down Expand Up @@ -4306,6 +4314,30 @@ subroutine vege_flux(parameters,nsnow ,nsoil ,isnow ,vegtyp ,veg , &
q2v = qsfc - ((evc+tr)/fveg+evg)/(latheav*rhoair) * 1./cq2v
endif

! use sfc_diag to calculate t2mv and q2v for opt_sfc=1&3
if(opt_stc == 1 .or. opt_stc == 3) then

fhi = fh2/fh
wrk = 1.0 - fhi
if(thsfc_loc) then ! Use local potential temperature
t2mv = tv*wrk + sfctmp*prslkix*fhi - (grav+grav)/cp
else ! Use potential temperature referenced to 1000 hPa
t2mv = tv*wrk + sfctmp*fhi - (grav+grav)/cp
endif

if(evg >= 0.) then ! for evaporation>0, use inferred qsurf to deduce q2v
q2v = qsfc*wrk + max(qmin,qair)*fhi
else ! for dew formation, use saturated q at tskin
qss = fpvs(tv)
qss = ep_2 * qss / (psfc + epsm1 * qss)
q2v= qss*wrk + max(qmin,qair)*fhi
endif
qss = fpvs(t2mv)
qss = ep_2 * qss / (psfc + epsm1 * qss)
q2v = min(q2v,qss)

endif

! update ch for output
ch = cah
chleaf = cvh
Expand All @@ -4325,7 +4357,7 @@ subroutine bare_flux (parameters,nsnow ,nsoil ,isnow ,dt ,sag , &
emg ,stc ,df ,rsurf ,lathea , & !in
gamma ,rhsur ,iloc ,jloc ,q2 ,pahb , & !in
thsfc_loc, prslkix,prsik1x,prslk1x,vegtyp,fveg,shdfac,garea1, & !in
pblhx , iz0tlnd , itime ,psi_opt,ep_1,ep_2,cp ,&
pblhx , iz0tlnd , itime ,psi_opt,ep_1,ep_2,epsm1, cp ,&
#ifdef CCPP
tgb ,cm ,ch,ustarx,errmsg ,errflg , & !inout
#else
Expand All @@ -4344,6 +4376,7 @@ subroutine bare_flux (parameters,nsnow ,nsoil ,isnow ,dt ,sag , &
! bare soil:
! -sab + irb[tg] + shb[tg] + evb[tg] + ghb[tg] = 0
! ----------------------------------------------------------------------
use funcphys, only : fpvs
implicit none
! ----------------------------------------------------------------------
! input
Expand Down Expand Up @@ -4381,6 +4414,7 @@ subroutine bare_flux (parameters,nsnow ,nsoil ,isnow ,dt ,sag , &
real (kind=kind_phys), intent(in) :: pblhx !pbl height (m)
real (kind=kind_phys), intent(in) :: ep_1
real (kind=kind_phys), intent(in) :: ep_2
real (kind=kind_phys), intent(in) :: epsm1
real (kind=kind_phys), intent(in) :: cp
integer, intent(in) :: iz0tlnd
integer, intent(in) :: itime
Expand Down Expand Up @@ -4529,6 +4563,10 @@ subroutine bare_flux (parameters,nsnow ,nsoil ,isnow ,dt ,sag , &
real (kind=kind_phys) :: temptrs

real (kind=kind_phys) :: t, tdc !kelvin to degree celsius with limit -50 to +50

real(kind=kind_phys) :: fhi, qss, wrk
real(kind=kind_phys), parameter :: qmin=1.0e-8

tdc(t) = min( 50., max(-50.,(t-tfrz)) )

! -----------------------------------------------------------------
Expand Down Expand Up @@ -4765,8 +4803,32 @@ subroutine bare_flux (parameters,nsnow ,nsoil ,isnow ,dt ,sag , &
end if
endif ! 4

! use sfc_diag to calculate t2mv and q2v for opt_sfc=1&3
if(opt_stc == 1 .or. opt_stc == 3) then

fhi = fh2/fh
wrk = 1.0 - fhi
if(thsfc_loc) then ! Use local potential temperature
t2mb = tgb*wrk + sfctmp*prslkix*fhi - (grav+grav)/cp
else ! Use potential temperature referenced to 1000 hPa
t2mb = tgb*wrk + sfctmp*fhi - (grav+grav)/cp
endif

if(evb >= 0.) then ! for evaporation>0, use inferred qsurf to deduce q2v
q2b = qsfc*wrk + max(qmin,qair)*fhi
else ! for dew formation, use saturated q at tskin
qss = fpvs(tgb)
qss = ep_2 * qss / (psfc + epsm1 * qss)
q2b= qss*wrk + max(qmin,qair)*fhi
endif
qss = fpvs(t2mb)
qss = ep_2 * qss / (psfc + epsm1 * qss)
q2b = min(q2b,qss)

endif
if (parameters%urban_flag) q2b = qsfc


! update ch
ch = ehb

Expand Down
2 changes: 1 addition & 1 deletion physics/noahmpdrv.F90
Original file line number Diff line number Diff line change
Expand Up @@ -879,7 +879,7 @@ subroutine noahmpdrv_run &
precip_graupel ,precip_hail ,temperature_soil_bot , &
co2_air ,o2_air ,foliage_nitrogen , &
snow_ice_frac_old ,forcing_height , &
con_fvirt ,con_eps ,con_cp , &
con_fvirt ,con_eps,con_epsm1 ,con_cp , &
snow_albedo_old ,snow_water_equiv_old , &
temperature_snow_soil ,soil_liquid_vol ,soil_moisture_vol , &
temperature_canopy_air,vapor_pres_canopy_air ,canopy_wet_fraction , &
Expand Down