Skip to content

Bugfix: Address an issue related to initialization of Noah-MP groundwater #463

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.

Sign up for GitHub

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

Merged
barlage merged 2 commits into from
May 4, 2018
Merged

Bugfix: Address an issue related to initialization of Noah-MP groundwater #463

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
330679e
Bugfix: Address an issue related to initialization of the groundwater…
Apr 26, 2018
bfca8ea
Wrap references and call to HALO_EM in subroutine groundwater_init wi…
May 3, 2018
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
8 changes: 4 additions & 4 deletions Registry/registry.noahmp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,14 +76,14 @@ state real chb2xy ij - 1 - i02rhd=(interp_mask_fiel
state real chstarxy ij - 1 - i02rhd=(interp_mask_field:lu_index,iswater)u=(copy_fcnm) "chstar" "dummy exchange coefficient" "m/s"
state real SMOISEQ ilj - 1 Z r "SMOISEQ" "EQ. SOIL MOISTURE" "m3 m-3"
state real smcwtdxy ij - 1 - rh "smcwtd" "deep soil moisture " "m3 m-3"
state real rechxy ij - 1 - h "rech" "water table recharge" "mm"
state real rechxy ij - 1 - rh "rech" "water table recharge" "mm"
state real deeprechxy ij - 1 - r "deeprech" "deep water table recharge" "mm"
state real areaxy ij - 1 - r "area" "area of grid boxes" "m2"
state real qrfxy ij - 1 - r "qrf" "baseflow " "m"
state real qrfsxy ij - 1 - h "qrfs" "sum baseflow " "mm"
state real qrfsxy ij - 1 - rh "qrfs" "sum baseflow " "mm"
state real qspringxy ij - 1 - r "qspring" "seeping water " "m"
state real qspringsxy ij - 1 - h "qsprings" "sum seeping water " "mm"
state real qslatxy ij - 1 - h "qslat" "sum lateral flow " "mm"
state real qspringsxy ij - 1 - rh "qsprings" "sum seeping water " "mm"
state real qslatxy ij - 1 - rh "qslat" "sum lateral flow " "mm"
state integer STEPWTD - misc 1 - r "STEPWTD" "NUMBER OF FUNDAMENTAL TIMESTEPS BETWEEN WTD CALLS" ""
state real pexpxy ij - 1 - r "pexp" "exponent for river conductance" "Kg m s-1"
state real rivercondxy ij - 1 - r "rivercond" "river conductance" "Kg m s-1"
Expand Down
20 changes: 20 additions & 0 deletions dyn_em/start_em.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read &
USE module_comm_dm
USE module_llxy, ONLY : proj_cassini
USE module_physics_init
USE module_sf_noahmpdrv, ONLY : groundwater_init
USE module_lightning_driver, ONLY : lightning_init
USE module_fr_fire_driver_wrf, ONLY : fire_driver_em_init
USE module_stoch, ONLY : setup_rand_perturb, rand_seed, update_stoch, initialize_stoch
Expand Down Expand Up @@ -1160,6 +1161,25 @@ SUBROUTINE start_domain_em ( grid, allowed_to_read &
,grid%aeropcu,grid%EFCG,grid%EFCS,grid%EFIG,grid%EFIS,grid%EFSG,grid%EFSS &
)
ENDDO ! loop of tiles for phy_init

if(.not.grid%restart) then

if(config_flags%sf_surface_physics .eq. NOAHMPSCHEME .and. config_flags%opt_run .eq. 5) then
call groundwater_init ( grid , &
config_flags%num_soil_layers, grid%dzs, grid%isltyp, grid%ivgtyp, config_flags%wtddt , &
grid%fdepthxy , grid%ht , grid%riverbedxy, grid%eqzwt , &
grid%rivercondxy, grid%pexpxy , grid%areaxy , grid%zwtxy , &
grid%smois , grid%sh2o , grid%smoiseq , grid%smcwtdxy , grid%deeprechxy, &
grid%rechxy , grid%qslatxy, grid%qrfsxy , grid%qspringsxy, &
grid%rechclim , &
ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
ips,ipe, jps,jpe, kps,kpe, &
its,ite, jts,jte, kts,kte )
endif

endif

ENDIF ! no phy_init for the backwards part of the DFI

CALL wrf_debug ( 100 , 'start_domain_em: After call to phy_init' )
Expand Down
67 changes: 44 additions & 23 deletions phys/module_sf_noahmpdrv.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1807,15 +1807,6 @@ SUBROUTINE NOAHMP_INIT ( MMINLU, SNOW , SNOWH , CANWAT , ISLTYP , IVGTYP, XLAT
STEPWTD = nint(WTDDT*60./DT)
STEPWTD = max(STEPWTD,1)

CALL groundwater_init ( &
& nsoil, zsoil , dzs ,isltyp, ivgtyp,wtddt , &
& fdepthxy, ht, riverbedxy, eqzwt, rivercondxy, pexpxy , areaxy, zwtxy, &
& smois,sh2o, smoiseq, smcwtdxy, deeprechxy, rechxy, qslatxy, qrfsxy, qspringsxy, &
& rechclim , &
& ids,ide, jds,jde, kds,kde, &
& ims,ime, jms,jme, kms,kme, &
& its,ite, jts,jte, kts,kte )

ELSE
CALL wrf_error_fatal ('Not enough fields to use groundwater option in Noah-MP')
END IF
Expand Down Expand Up @@ -1932,30 +1923,40 @@ END SUBROUTINE SNOW_INIT
! ==================================================================================================
! ----------------------------------------------------------------------
SUBROUTINE GROUNDWATER_INIT ( &
& NSOIL , ZSOIL , DZS, ISLTYP, IVGTYP, WTDDT , &
& GRID, NSOIL , DZS, ISLTYP, IVGTYP, WTDDT , &
& FDEPTH, TOPO, RIVERBED, EQWTD, RIVERCOND, PEXP , AREA ,WTD , &
& SMOIS,SH2O, SMOISEQ, SMCWTDXY, DEEPRECHXY, RECHXY , &
& QSLATXY, QRFSXY, QSPRINGSXY, &
& rechclim , &
& ids,ide, jds,jde, kds,kde, &
& ims,ime, jms,jme, kms,kme, &
& ips,ipe, jps,jpe, kps,kpe, &
& its,ite, jts,jte, kts,kte )


USE NOAHMP_TABLES, ONLY : BEXP_TABLE,SMCMAX_TABLE,PSISAT_TABLE,SMCWLT_TABLE,DWSAT_TABLE,DKSAT_TABLE, &
ISURBAN_TABLE, ISICE_TABLE ,ISWATER_TABLE
USE module_sf_noahmp_groundwater, ONLY : LATERALFLOW
USE module_domain, only: domain
#if (EM_CORE == 1)
#ifdef DM_PARALLEL
USE module_dm , ONLY : ntasks_x,ntasks_y,local_communicator,mytask,ntasks
USE module_comm_dm , ONLY : halo_em_hydro_noahmp_sub
#endif
#endif

! ----------------------------------------------------------------------
IMPLICIT NONE
! ----------------------------------------------------------------------

INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, &
& ims,ime, jms,jme, kms,kme, &
& ips,ipe, jps,jpe, kps,kpe, &
& its,ite, jts,jte, kts,kte
TYPE(domain) , TARGET :: grid ! state
INTEGER, INTENT(IN) :: NSOIL
REAL, INTENT(IN) :: WTDDT
REAL, INTENT(IN), DIMENSION(1:NSOIL) :: ZSOIL,DZS
REAL, INTENT(IN), DIMENSION(1:NSOIL) :: DZS
INTEGER, INTENT(IN), DIMENSION(ims:ime, jms:jme) :: ISLTYP, IVGTYP
REAL, INTENT(IN), DIMENSION(ims:ime, jms:jme) :: FDEPTH, TOPO , AREA
REAL, INTENT(IN), DIMENSION(ims:ime, jms:jme) :: rechclim
Expand All @@ -1973,15 +1974,22 @@ SUBROUTINE GROUNDWATER_INIT ( &
QRFSXY, &
QSPRINGSXY
! local
INTEGER :: I,J,K,ITER,itf,jtf, NITER, NCOUNT
INTEGER :: I,J,K,ITER,itf,jtf, NITER, NCOUNT,NS
REAL :: BEXP,SMCMAX,PSISAT,SMCWLT,DWSAT,DKSAT
REAL :: FRLIQ,SMCEQDEEP
REAL :: DELTAT,RCOND,TOTWATER
REAL :: AA,BBB,CC,DD,DX,FUNC,DFUNC,DDZ,EXPON,SMC,FLUX
REAL, DIMENSION(1:NSOIL) :: SMCEQ
REAL, DIMENSION(1:NSOIL) :: SMCEQ,ZSOIL
REAL, DIMENSION( ims:ime, jms:jme ) :: QLAT, QRF
INTEGER, DIMENSION( ims:ime, jms:jme ) :: LANDMASK !-1 for water (ice or no ice) and glacial areas, 1 for land where the LSM does its soil moisture calculations

! Given the soil layer thicknesses (in DZS), calculate the soil layer
! depths from the surface.
ZSOIL(1) = -DZS(1) ! negative
DO NS=2, NSOIL
ZSOIL(NS) = ZSOIL(NS-1) - DZS(NS)
END DO


itf=min0(ite,ide-1)
jtf=min0(jte,jde-1)
Expand All @@ -1999,33 +2007,51 @@ SUBROUTINE GROUNDWATER_INIT ( &

!readjust the raw aggregated water table from hires, so that it is better compatible with topography

DO NITER=1,500
!use WTD here, to use the lateral communication routine
WTD=EQWTD

NCOUNT=0

DO NITER=1,500

#if (EM_CORE == 1)
#ifdef DM_PARALLEL
# include "HALO_EM_HYDRO_NOAHMP.inc"
#endif
#endif

!Calculate lateral flow

IF(NCOUNT.GT.0.OR.NITER.eq.1)THEN
QLAT = 0.
CALL LATERALFLOW(ISLTYP,EQWTD,QLAT,FDEPTH,TOPO,LANDMASK,DELTAT,AREA &
CALL LATERALFLOW(ISLTYP,WTD,QLAT,FDEPTH,TOPO,LANDMASK,DELTAT,AREA &
,ids,ide,jds,jde,kds,kde &
,ims,ime,jms,jme,kms,kme &
,its,ite,jts,jte,kts,kte )

NCOUNT=0
DO J=jts,jtf
DO I=its,itf
IF(LANDMASK(I,J).GT.0)THEN
IF(QLAT(i,j).GT.1.e-2)THEN
NCOUNT=NCOUNT+1
EQWTD(i,j)=min(EQWTD(i,j)+0.25,0.)
WTD(i,j)=min(WTD(i,j)+0.25,0.)
ENDIF
ENDIF
ENDDO
ENDDO

IF(NCOUNT.EQ.0)EXIT
ENDIF

ENDDO

#if (EM_CORE == 1)
#ifdef DM_PARALLEL
# include "HALO_EM_HYDRO_NOAHMP.inc"
#endif
#endif

EQWTD=WTD

!after adjusting, where qlat > 1cm/year now wtd is at the surface.
!it may still happen that qlat + rech > 0 and eqwtd-rbed <0. There the wtd can
!rise to the surface (poor drainage) but the et will then increase.
Expand Down Expand Up @@ -2062,11 +2088,6 @@ SUBROUTINE GROUNDWATER_INIT ( &

RIVERBED = min( RIVERBED-TOPO, 0.)

!now inititalize wtd

WTD = EQWTD


!now recompute lateral flow and flow to rivers to initialize deep soil moisture

DELTAT = WTDDT * 60. !timestep in seconds for this calculation
Expand Down