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

Update from master #322

Merged
merged 5 commits into from
Oct 27, 2024
Merged
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
2 changes: 1 addition & 1 deletion main/MOD_Forcing.F90
Original file line number Diff line number Diff line change
Expand Up @@ -487,7 +487,7 @@ SUBROUTINE read_forcing (idate, dir_forcing)
cosz = orb_coszen(calday, gforc%rlon(ilon), gforc%rlat(ilat))
cosz = max(0.001, cosz)
! 10/24/2024, yuan: deal with time log with backward or foreward
IF (timelog(ivar) == 'foreward') THEN
IF (trim(timelog(ivar)) == 'foreward') THEN
forcn(ivar)%blk(ib,jb)%val(i,j) = &
cosz / avgcos%blk(ib,jb)%val(i,j) * forcn_LB(ivar)%blk(ib,jb)%val(i,j)
ELSE
Expand Down
12 changes: 7 additions & 5 deletions main/MOD_GroundFluxes.F90
Original file line number Diff line number Diff line change
Expand Up @@ -32,18 +32,20 @@ SUBROUTINE GroundFluxes (zlnd, zsno, hu, ht, hq, hpbl, &
! Original author : Yongjiu Dai, 09/15/1999; 08/30/2002
!
! REVISIONS:
! Hua Yuan, 09/2019: removed sigf to be consistant with PFT runs, removed fsena,
! fevpa, renamed z0ma to z0m
! Shaofeng Liu, 05/2023: add option to call moninobuk_leddy, the LargeEddy
! surface turbulence scheme (LZD2022);
! make a proper update of um.
! 09/2019, Hua Yuan: removed sigf to be consistant with PFT runs, removed
! fsena, fevpa, renamed z0ma to z0m.
!
! 05/2023, Shaofeng Liu: add option to call moninobuk_leddy, the LargeEddy
! surface turbulence scheme (LZD2022); make a proper update of um.
!
!=======================================================================

USE MOD_Precision
USE MOD_Const_Physical, only: cpair,vonkar,grav
USE MOD_FrictionVelocity
USE mod_namelist, only: DEF_USE_CBL_HEIGHT,DEF_RSS_SCHEME
USE MOD_TurbulenceLEddy

IMPLICIT NONE

!----------------------- Dummy argument --------------------------------
Expand Down
32 changes: 26 additions & 6 deletions main/MOD_LeafTemperature.F90
Original file line number Diff line number Diff line change
Expand Up @@ -105,6 +105,7 @@ SUBROUTINE LeafTemperature ( &
USE MOD_CanopyLayerProfile
USE MOD_TurbulenceLEddy
USE MOD_AssimStomataConductance
USE MOD_UserSpecifiedForcing, only: HEIGHT_mode
USE MOD_Vars_TimeInvariants, only: patchclass
USE MOD_Const_LC, only: z0mr, displar
USE MOD_PlantHydraulic, only :PlantHydraulicStress_twoleaf, getvegwp_twoleaf
Expand Down Expand Up @@ -304,6 +305,9 @@ SUBROUTINE LeafTemperature ( &

real(r8) :: &
displa, &! displacement height [m]
hu_, &! adjusted observational height of wind [m]
ht_, &! adjusted observational height of temperature [m]
hq_, &! adjusted observational height of humidity [m]
zldis, &! reference height "minus" zero displacement heght [m]
zii, &! convective boundary layer height [m]
z0mv, &! roughness length, momentum [m]
Expand Down Expand Up @@ -518,7 +522,23 @@ SUBROUTINE LeafTemperature ( &
dth = thm - taf
dqh = qm - qaf
dthv = dth*(1.+0.61*qm) + 0.61*th*dqh
zldis = hu - displa

hu_ = hu; ht_ = ht; hq_ = hq;

IF (trim(HEIGHT_mode) == 'absolute') THEN
#ifndef SingPoint
! to ensure the obs height >= htop+10.
hu_ = max(htop+10., hu)
ht_ = max(htop+10., ht)
hq_ = max(htop+10., hq)
#endif
ELSE ! relative height
hu_ = htop + hu
ht_ = htop + ht
hq_ = htop + hq
ENDIF

zldis = hu_ - displa

IF(zldis <= 0.0) THEN
write(6,*) 'the obs height of u less than the zero displacement heght'
Expand Down Expand Up @@ -550,11 +570,11 @@ SUBROUTINE LeafTemperature ( &
! Evaluate stability-dependent variables using moz from prior iteration
IF (rd_opt == 3) THEN
IF (DEF_USE_CBL_HEIGHT) THEN
CALL moninobukm_leddy(hu,ht,hq,displa,z0mv,z0hv,z0qv,obu,um, &
CALL moninobukm_leddy(hu_,ht_,hq_,displa,z0mv,z0hv,z0qv,obu,um, &
displasink,z0mv,hpbl,ustar,fh2m,fq2m, &
htop,fmtop,fm,fh,fq,fht,fqt,phih)
ELSE
CALL moninobukm(hu,ht,hq,displa,z0mv,z0hv,z0qv,obu,um, &
CALL moninobukm(hu_,ht_,hq_,displa,z0mv,z0hv,z0qv,obu,um, &
displasink,z0mv,ustar,fh2m,fq2m, &
htop,fmtop,fm,fh,fq,fht,fqt,phih)
ENDIF
Expand All @@ -564,10 +584,10 @@ SUBROUTINE LeafTemperature ( &
raw = 1./(vonkar/(fq-fqt)*ustar)
ELSE
IF (DEF_USE_CBL_HEIGHT) THEN
CALL moninobuk_leddy(hu,ht,hq,displa,z0mv,z0hv,z0qv,obu,um,hpbl, &
CALL moninobuk_leddy(hu_,ht_,hq_,displa,z0mv,z0hv,z0qv,obu,um,hpbl, &
ustar,fh2m,fq2m,fm10m,fm,fh,fq)
ELSE
CALL moninobuk(hu,ht,hq,displa,z0mv,z0hv,z0qv,obu,um,&
CALL moninobuk(hu_,ht_,hq_,displa,z0mv,z0hv,z0qv,obu,um,&
ustar,fh2m,fq2m,fm10m,fm,fh,fq)
ENDIF
! Aerodynamic resistance
Expand Down Expand Up @@ -918,7 +938,7 @@ SUBROUTINE LeafTemperature ( &
um = max(ur,.1)
ELSE
IF (DEF_USE_CBL_HEIGHT) THEN !//TODO: Shaofeng, 2023.05.18
zii = max(5.*hu,hpbl)
zii = max(5.*hu_,hpbl)
ENDIF !//TODO: Shaofeng, 2023.05.18
wc = (-grav*ustar*thvstar*zii/thv)**(1./3.)
wc2 = beta*beta*(wc*wc)
Expand Down
35 changes: 28 additions & 7 deletions main/MOD_LeafTemperaturePC.F90
Original file line number Diff line number Diff line change
Expand Up @@ -91,14 +91,14 @@ SUBROUTINE LeafTemperaturePC ( &
!
! !REVISIONS:
!
! 01/2021, Xingjie Lu and Nan Wei: added plant hydraulic process interface
! 01/2021, Xingjie Lu and Nan Wei: added plant hydraulic process interface.
!
! 01/2021, Nan Wei: added interaction btw prec and canopy
! 01/2021, Nan Wei: added interaction btw prec and canopy.
!
! 05/2023, Shaofeng Liu: add option to call moninobuk_leddy, the LargeEddy
! surface turbulence scheme (LZD2022); make a proper update of um.
!
! 04/2024, Hua Yuan: add option to account for vegetation snow process
! 04/2024, Hua Yuan: add option to account for vegetation snow process.
!
!=======================================================================

Expand All @@ -114,6 +114,8 @@ SUBROUTINE LeafTemperaturePC ( &
USE MOD_AssimStomataConductance
USE MOD_PlantHydraulic, only: PlantHydraulicStress_twoleaf
USE MOD_Ozone, only: CalcOzoneStress
USE MOD_UserSpecifiedForcing, only: HEIGHT_mode

IMPLICIT NONE

!-----------------------Arguments---------------------------------------
Expand Down Expand Up @@ -325,6 +327,9 @@ SUBROUTINE LeafTemperaturePC ( &
rootfr(nl_soil,ps:pe) ! root fraction

real(r8) :: &
hu_, &! adjusted observational height of wind [m]
ht_, &! adjusted observational height of temperature [m]
hq_, &! adjusted observational height of humidity [m]
zldis, &! reference height "minus" zero displacement heght [m]
zii, &! convective boundary layer height [m]
z0mv, &! roughness length, momentum [m]
Expand Down Expand Up @@ -901,7 +906,23 @@ SUBROUTINE LeafTemperaturePC ( &
dth = thm - taf(toplay)
dqh = qm - qaf(toplay)
dthv = dth*(1.+0.61*qm) + 0.61*th*dqh
zldis = hu - displa_lays(3)

hu_ = hu; ht_ = ht; hq_ = hq;

IF (trim(HEIGHT_mode) == 'absolute') THEN
#ifndef SingPoint
! to ensure the obs height >= htop+10.
hu_ = max(htop_lay(toplay)+10., hu)
ht_ = max(htop_lay(toplay)+10., ht)
hq_ = max(htop_lay(toplay)+10., hq)
#endif
ELSE ! relative height
hu_ = htop_lay(toplay) + hu
ht_ = htop_lay(toplay) + ht
hq_ = htop_lay(toplay) + hq
ENDIF

zldis = hu_ - displa_lays(3)

IF(zldis <= 0.0) THEN
write(6,*) 'the obs height of u less than the zero displacement heght'
Expand Down Expand Up @@ -935,11 +956,11 @@ SUBROUTINE LeafTemperaturePC ( &
! Evaluate stability-dependent variables using moz from prior iteration

IF (DEF_USE_CBL_HEIGHT) THEN
CALL moninobukm_leddy(hu,ht,hq,displa_lays(toplay),z0mv,z0hv,z0qv,obu,um, &
CALL moninobukm_leddy(hu_,ht_,hq_,displa_lays(toplay),z0mv,z0hv,z0qv,obu,um, &
displa_lay(toplay),z0m_lay(toplay),hpbl,ustar,fh2m,fq2m, &
htop_lay(toplay),fmtop,fm,fh,fq,fht,fqt,phih)
ELSE
CALL moninobukm(hu,ht,hq,displa_lays(toplay),z0mv,z0hv,z0qv,obu,um, &
CALL moninobukm(hu_,ht_,hq_,displa_lays(toplay),z0mv,z0hv,z0qv,obu,um, &
displa_lay(toplay),z0m_lay(toplay),ustar,fh2m,fq2m, &
htop_lay(toplay),fmtop,fm,fh,fq,fht,fqt,phih)
ENDIF
Expand Down Expand Up @@ -1659,7 +1680,7 @@ SUBROUTINE LeafTemperaturePC ( &
um = max(ur,.1)
ELSE
IF (DEF_USE_CBL_HEIGHT) THEN !//TODO: Shaofeng, 2023.05.18
zii = max(5.*hu,hpbl)
zii = max(5.*hu_,hpbl)
ENDIF !//TODO: Shaofeng, 2023.05.18
wc = (-grav*ustar*thvstar*zii/thv)**(1./3.)
wc2 = beta*beta*(wc*wc)
Expand Down
3 changes: 3 additions & 0 deletions main/MOD_UserSpecifiedForcing.F90
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,8 @@ MODULE MOD_UserSpecifiedForcing
character(len=256) :: dataset

logical :: solarin_all_band ! whether solar radiation in all bands is available

character(len=256) :: HEIGHT_mode ! observation height mode
real(r8) :: HEIGHT_V ! observation height of wind speed
real(r8) :: HEIGHT_T ! observation height of air temperature
real(r8) :: HEIGHT_Q ! observation height of specific humidity
Expand Down Expand Up @@ -110,6 +112,7 @@ SUBROUTINE init_user_specified_forcing
allocate (tintalgo (NVAR))

solarin_all_band = DEF_forcing%solarin_all_band ! whether solar radiation in all bands is available
HEIGHT_mode = DEF_forcing%HEIGHT_mode ! observation height mode
HEIGHT_V = DEF_forcing%HEIGHT_V ! observation height of wind speed
HEIGHT_T = DEF_forcing%HEIGHT_T ! observation height of air temperature
HEIGHT_Q = DEF_forcing%HEIGHT_Q ! observation height of specific humidity
Expand Down
58 changes: 40 additions & 18 deletions main/URBAN/MOD_Urban_Flux.F90
Original file line number Diff line number Diff line change
Expand Up @@ -111,6 +111,7 @@ SUBROUTINE UrbanOnlyFlux ( &
USE MOD_Const_Physical, only: cpair,vonkar,grav,hvap
USE MOD_FrictionVelocity
USE MOD_CanopyLayerProfile
USE MOD_UserSpecifiedForcing, only: HEIGHT_mode
IMPLICIT NONE

!----------------------- Dummy argument --------------------------------
Expand Down Expand Up @@ -262,9 +263,9 @@ SUBROUTINE UrbanOnlyFlux ( &
numlay ! available layer number

real(r8) :: &
huu, &! observational height of wind [m]
htu, &! observational height of temperature [m]
hqu, &! observational height of humidity [m]
hu_, &! adjusted observational height of wind [m]
ht_, &! adjusted observational height of temperature [m]
hq_, &! adjusted observational height of humidity [m]
ktop, &! K value at a specific height
utop, &! u value at a specific height
fht, &! integral of profile function for heat at the top layer
Expand Down Expand Up @@ -509,12 +510,22 @@ SUBROUTINE UrbanOnlyFlux ( &
dqh = qm - qaf(2)
dthv = dth*(1.+0.61*qm) + 0.61*th*dqh

! to ensure the obs height >= hroof+10.
huu = max(hroof+10., hu)
htu = max(hroof+10., ht)
hqu = max(hroof+10., hq)
hu_ = hu; ht_ = ht; hq_ = hq;

zldis = huu - displa
IF (trim(HEIGHT_mode) == 'absolute') THEN
#ifndef SingPoint
! to ensure the obs height >= hroof+10.
hu_ = max(hroof+10., hu)
ht_ = max(hroof+10., ht)
hq_ = max(hroof+10., hq)
#endif
ELSE ! relative height
hu_ = hroof + hu
ht_ = hroof + ht
hq_ = hroof + hq
ENDIF

zldis = hu_ - displa

IF (zldis <= 0.0) THEN
write(6,*) 'the obs height of u less than the zero displacement heght'
Expand All @@ -538,7 +549,7 @@ SUBROUTINE UrbanOnlyFlux ( &

!NOTE: displat=hroof, z0mt=0, are set for roof
! fmtop is calculated at the same height of fht, fqt
CALL moninobukm(huu,htu,hqu,displa,z0m,z0h,z0q,obu,um, &
CALL moninobukm(hu_,ht_,hq_,displa,z0m,z0h,z0q,obu,um, &
hroof,0.,ustar,fh2m,fq2m,hroof,fmtop,fm,fh,fq,fht,fqt,phih)

! Aerodynamic resistance
Expand Down Expand Up @@ -885,6 +896,7 @@ SUBROUTINE UrbanVegFlux ( &
USE MOD_FrictionVelocity
USE MOD_CanopyLayerProfile
USE MOD_AssimStomataConductance
USE MOD_UserSpecifiedForcing, only: HEIGHT_mode
IMPLICIT NONE

!-----------------------Arguments---------------------------------------
Expand Down Expand Up @@ -1150,9 +1162,9 @@ SUBROUTINE UrbanVegFlux ( &
numlay ! available layer number

real(r8) :: &
huu, &! observational height of wind [m]
htu, &! observational height of temperature [m]
hqu, &! observational height of humidity [m]
hu_, &! adjusted observational height of wind [m]
ht_, &! adjusted observational height of temperature [m]
hq_, &! adjusted observational height of humidity [m]
ktop, &! K value at a specific height
utop, &! u value at a specific height
fht, &! integral of profile function for heat at the top layer
Expand Down Expand Up @@ -1487,12 +1499,22 @@ SUBROUTINE UrbanVegFlux ( &
dqh = qm - qaf(2)
dthv = dth*(1.+0.61*qm) + 0.61*th*dqh

! To ensure the obs height >= hroof+10.
huu = max(hroof+10., hu)
htu = max(hroof+10., ht)
hqu = max(hroof+10., hq)
hu_ = hu; ht_ = ht; hq_ = hq;

IF (trim(HEIGHT_mode) == 'absolute') THEN
#ifndef SingPoint
! to ensure the obs height >= hroof+10.
hu_ = max(hroof+10., hu)
ht_ = max(hroof+10., ht)
hq_ = max(hroof+10., hq)
#endif
ELSE ! relative height
hu_ = hroof + hu
ht_ = hroof + ht
hq_ = hroof + hq
ENDIF

zldis = huu - displa
zldis = hu_ - displa

IF (zldis <= 0.0) THEN
write(6,*) 'the obs height of u less than the zero displacement heght'
Expand All @@ -1517,7 +1539,7 @@ SUBROUTINE UrbanVegFlux ( &
!-----------------------------------------------------------------------
! Evaluate stability-dependent variables using moz from prior iteration

CALL moninobukm(huu,htu,hqu,displa,z0m,z0h,z0q,obu,um, &
CALL moninobukm(hu_,ht_,hq_,displa,z0m,z0h,z0q,obu,um, &
hroof,0.,ustar,fh2m,fq2m,hroof,fmtop,fm,fh,fq,fht,fqt,phih)

! Aerodynamic resistance
Expand Down
3 changes: 3 additions & 0 deletions run/forcing/ERA5.nml
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,9 @@

DEF_forcing%dataset = 'ERA5'
DEF_forcing%solarin_all_band = .true.
DEF_forcing%HEIGHT_mode = 'absolute' ! 'absolute' height to ground surface or
! 'relative' height to tree/building top
! if not speficied, default is absolute
DEF_forcing%HEIGHT_V = 50.0
DEF_forcing%HEIGHT_T = 40.0
DEF_forcing%HEIGHT_Q = 40.0
Expand Down
9 changes: 6 additions & 3 deletions run/forcing/ERA5LAND.nml
Original file line number Diff line number Diff line change
Expand Up @@ -5,9 +5,12 @@

DEF_forcing%dataset = 'ERA5LAND'
DEF_forcing%solarin_all_band = .true.
DEF_forcing%HEIGHT_V = 50.0
DEF_forcing%HEIGHT_T = 40.
DEF_forcing%HEIGHT_Q = 40.
DEF_forcing%HEIGHT_mode = 'relative' ! 'absolute' height to ground surface or
! 'relative' height to tree/building top
! if not speficied, default is absolute
DEF_forcing%HEIGHT_V = 10.0
DEF_forcing%HEIGHT_T = 10.0
DEF_forcing%HEIGHT_Q = 10.0

DEF_forcing%has_missing_value = .true.
DEF_forcing%missing_value_name = 'missing_value'
Expand Down
2 changes: 2 additions & 0 deletions share/MOD_Namelist.F90
Original file line number Diff line number Diff line change
Expand Up @@ -356,6 +356,7 @@ MODULE MOD_Namelist

character(len=256) :: dataset = 'CRUNCEP'
logical :: solarin_all_band = .true.
character(len=256) :: HEIGHT_mode = 'absolute'
real(r8) :: HEIGHT_V = 100.0
real(r8) :: HEIGHT_T = 50.
real(r8) :: HEIGHT_Q = 50.
Expand Down Expand Up @@ -1443,6 +1444,7 @@ SUBROUTINE read_namelist (nlfile)

CALL mpi_bcast (DEF_forcing%dataset ,256 ,mpi_character ,p_address_master ,p_comm_glb ,p_err)
CALL mpi_bcast (DEF_forcing%solarin_all_band ,1 ,mpi_logical ,p_address_master ,p_comm_glb ,p_err)
CALL mpi_bcast (DEF_forcing%HEIGHT_mode ,256 ,mpi_character ,p_address_master ,p_comm_glb ,p_err)
CALL mpi_bcast (DEF_forcing%HEIGHT_V ,1 ,mpi_real8 ,p_address_master ,p_comm_glb ,p_err)
CALL mpi_bcast (DEF_forcing%HEIGHT_T ,1 ,mpi_real8 ,p_address_master ,p_comm_glb ,p_err)
CALL mpi_bcast (DEF_forcing%HEIGHT_Q ,1 ,mpi_real8 ,p_address_master ,p_comm_glb ,p_err)
Expand Down