diff --git a/drivers/hrldas/ConfigVarInTransferMod.F90 b/drivers/hrldas/ConfigVarInTransferMod.F90 index 4fa054b2..8ea389e1 100644 --- a/drivers/hrldas/ConfigVarInTransferMod.F90 +++ b/drivers/hrldas/ConfigVarInTransferMod.F90 @@ -60,6 +60,7 @@ subroutine ConfigVarInTransfer(noahmp, NoahmpIO) noahmp%config%nmlist%OptRunoffSubsurface = NoahmpIO%IOPT_RUNSUB noahmp%config%nmlist%OptGlacierTreatment = NoahmpIO%IOPT_GLA noahmp%config%nmlist%OptSnowCompaction = NoahmpIO%IOPT_COMPACT + noahmp%config%nmlist%OptWetlandModel = NoahmpIO%IOPT_WETLAND ! config domain variable noahmp%config%domain%SurfaceType = 1 diff --git a/drivers/hrldas/NoahmpDriverMainMod.F90 b/drivers/hrldas/NoahmpDriverMainMod.F90 index fd23b55f..dd593a8b 100644 --- a/drivers/hrldas/NoahmpDriverMainMod.F90 +++ b/drivers/hrldas/NoahmpDriverMainMod.F90 @@ -70,6 +70,7 @@ subroutine NoahmpDriverMain(NoahmpIO) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 endif if ( NoahmpIO%SOIL_UPDATE_STEPS > 1 ) then @@ -83,6 +84,7 @@ subroutine NoahmpDriverMain(NoahmpIO) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 end if endif diff --git a/drivers/hrldas/NoahmpIOVarInitMod.F90 b/drivers/hrldas/NoahmpIOVarInitMod.F90 index b09a0cd7..640a599f 100644 --- a/drivers/hrldas/NoahmpIOVarInitMod.F90 +++ b/drivers/hrldas/NoahmpIOVarInitMod.F90 @@ -297,6 +297,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) if ( .not. allocated (NoahmpIO%ACC_QINSURXY)) allocate ( NoahmpIO%ACC_QINSURXY (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%ACC_QSEVAXY) ) allocate ( NoahmpIO%ACC_QSEVAXY (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%ACC_ETRANIXY)) allocate ( NoahmpIO%ACC_ETRANIXY (XSTART:XEND,1:NSOIL,YSTART:YEND) ) + if ( .not. allocated (NoahmpIO%ACC_GLAFLWXY)) allocate ( NoahmpIO%ACC_GLAFLWXY (XSTART:XEND, YSTART:YEND) ) ! Needed for MMF_RUNOFF (IOPT_RUN = 5); not part of MP driver in WRF if ( .not. allocated (NoahmpIO%MSFTX) ) allocate ( NoahmpIO%MSFTX (XSTART:XEND,YSTART:YEND) ) @@ -325,6 +326,16 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) if ( .not. allocated (NoahmpIO%SEASON_GDD)) allocate ( NoahmpIO%SEASON_GDD (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%CROPTYPE) ) allocate ( NoahmpIO%CROPTYPE (XSTART:XEND,5,YSTART:YEND) ) + ! Needed for Zhang et al. 2022 wetland model (OPT_WETLAND=1 or 2) + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + if ( .not. allocated (NoahmpIO%FSATXY) ) allocate ( NoahmpIO%FSATXY (XSTART:XEND, YSTART:YEND) ) ! saturated fraction of the grid (-) + if ( .not. allocated (NoahmpIO%WSURFXY)) allocate ( NoahmpIO%WSURFXY (XSTART:XEND, YSTART:YEND) ) ! wetland water storage [mm] + endif + if ( NoahmpIO%IOPT_WETLAND == 2 ) then + if ( .not. allocated (NoahmpIO%FSATMX) ) allocate ( NoahmpIO%FSATMX (XSTART:XEND, YSTART:YEND) ) ! maximum saturated fraction + if ( .not. allocated (NoahmpIO%WCAP) ) allocate ( NoahmpIO%WCAP (XSTART:XEND, YSTART:YEND) ) ! maximum wetland capacity [m] + endif + ! Single- and Multi-layer Urban Models if ( NoahmpIO%SF_URBAN_PHYSICS > 0 ) then if ( .not. allocated (NoahmpIO%sh_urb2d) ) allocate ( NoahmpIO%sh_urb2d (XSTART:XEND,YSTART:YEND) ) @@ -339,7 +350,6 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) if ( .not. allocated (NoahmpIO%lb_urb2d) ) allocate ( NoahmpIO%lb_urb2d (XSTART:XEND,YSTART:YEND) ) if ( .not. allocated (NoahmpIO%hgt_urb2d) ) allocate ( NoahmpIO%hgt_urb2d (XSTART:XEND,YSTART:YEND) ) if ( .not. allocated (NoahmpIO%ust) ) allocate ( NoahmpIO%ust (XSTART:XEND,YSTART:YEND) ) - !ENDIF !IF(NoahmpIO%SF_URBAN_PHYSICS == 1 ) THEN ! single layer urban model if ( .not. allocated (NoahmpIO%cmr_sfcdif) ) allocate ( NoahmpIO%cmr_sfcdif (XSTART:XEND, YSTART:YEND) ) @@ -651,6 +661,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 ! MMF Groundwater NoahmpIO%TERRAIN = undefined_real @@ -701,6 +712,16 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) NoahmpIO%IRRSPLH = 0.0 NoahmpIO%LOCTIM = undefined_real + ! wetland model (Zhang et al. 2022) + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + NoahmpIO%FSATXY = undefined_real + NoahmpIO%WSURFXY = undefined_real + endif + if ( NoahmpIO%IOPT_WETLAND == 2 ) then + NoahmpIO%FSATMX = undefined_real + NoahmpIO%WCAP = undefined_real + endif + ! spatial varying soil texture if ( NoahmpIO%IOPT_SOIL > 1 ) then NoahmpIO%SOILCL1 = undefined_real diff --git a/drivers/hrldas/NoahmpIOVarType.F90 b/drivers/hrldas/NoahmpIOVarType.F90 index e0be132c..12581dbf 100644 --- a/drivers/hrldas/NoahmpIOVarType.F90 +++ b/drivers/hrldas/NoahmpIOVarType.F90 @@ -64,6 +64,7 @@ module NoahmpIOVarType integer :: IOPT_INFDV ! infiltration options for dynamic VIC (1->Philip; 2-> Green-Ampt;3->Smith-Parlange) integer :: IOPT_TDRN ! drainage option (0->off; 1->simple scheme; 2->Hooghoudt's scheme) integer :: IOPT_COMPACT ! snowpack compaction (1->Anderson1976; 2->Abolafia-Rosenzweig2024) + integer :: IOPT_WETLAND ! wetland model option (0->off; 1->Zhang2022 fixed parameter; 2->Zhang2022 read in 2D parameter) real(kind=kind_noahmp) :: XICE_THRESHOLD ! fraction of grid determining seaice real(kind=kind_noahmp) :: JULIAN ! Julian day real(kind=kind_noahmp) :: DTBL ! timestep [s] @@ -331,6 +332,7 @@ module NoahmpIOVarType real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_ETRANXY ! accumulated transpiration per soil timestep [mm] real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_EDIRXY ! accumulated net ground (soil/snow) evaporation per soil timestep [mm] real(kind=kind_noahmp), allocatable, dimension(:,:,:) :: ACC_ETRANIXY ! accumualted transpiration rate within soil timestep [m/s * dt_soil/dt_main] + real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_GLAFLWXY ! accumulated glacier excessive flow [mm] per soil timestep !------------------------------------------------------------------------ ! Needed for MMF_RUNOFF (IOPT_RUN = 5); not part of MP driver in WRF @@ -372,6 +374,14 @@ module NoahmpIOVarType real(kind=kind_noahmp), allocatable, dimension(:,:) :: SEASON_GDD ! seasonal GDD real(kind=kind_noahmp), allocatable, dimension(:,:,:) :: CROPTYPE ! crop type +!------------------------------------------------------------------------ +! Needed for wetland model (OPT_WETLAND=1 or 2) +!------------------------------------------------------------------------ + real(kind=kind_noahmp), allocatable, dimension(:,:) :: FSATXY ! saturated fraction of the grid (-) + real(kind=kind_noahmp), allocatable, dimension(:,:) :: WSURFXY ! wetland water storage [mm] + real(kind=kind_noahmp), allocatable, dimension(:,:) :: FSATMX ! maximum saturated fraction + real(kind=kind_noahmp), allocatable, dimension(:,:) :: WCAP ! maximum wetland capacity [m] + !------------------------------------------------------------------------ ! Single- and Multi-layer Urban Models !------------------------------------------------------------------------ @@ -878,6 +888,9 @@ module NoahmpIOVarType real(kind=kind_noahmp), allocatable, dimension(:,:) :: RTCT_TABLE ! fraction of carbohydrate translocation from root to grain real(kind=kind_noahmp), allocatable, dimension(:) :: BIO2LAI_TABLE ! leaf area per living leaf biomass [m2/kg] + ! wetland parameter (OPT_WETLAND=1) + real(kind=kind_noahmp) :: WCAP_TABLE ! maximum surface wetland capacity + ! soil parameters integer :: SLCATS_TABLE ! number of soil categories real(kind=kind_noahmp), allocatable, dimension(:) :: BEXP_TABLE ! soil B parameter diff --git a/drivers/hrldas/NoahmpInitMainMod.F90 b/drivers/hrldas/NoahmpInitMainMod.F90 index 7a590791..a1ab2fe1 100644 --- a/drivers/hrldas/NoahmpInitMainMod.F90 +++ b/drivers/hrldas/NoahmpInitMainMod.F90 @@ -151,6 +151,10 @@ subroutine NoahmpInitMain(NoahmpIO) NoahmpIO%QSNOWXY(I,J) = 0.0 NoahmpIO%QRAINXY(I,J) = 0.0 NoahmpIO%WSLAKEXY(I,J) = 0.0 + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + NoahmpIO%FSATXY(I,J) = 0.0 + NoahmpIO%WSURFXY(I,J) = 0.0 + endif if ( NoahmpIO%IOPT_RUNSUB /= 5 ) then NoahmpIO%WAXY(I,J) = 4900.0 NoahmpIO%WTXY(I,J) = NoahmpIO%WAXY(i,j) diff --git a/drivers/hrldas/NoahmpReadNamelistMod.F90 b/drivers/hrldas/NoahmpReadNamelistMod.F90 index 1c0fe477..34701629 100644 --- a/drivers/hrldas/NoahmpReadNamelistMod.F90 +++ b/drivers/hrldas/NoahmpReadNamelistMod.F90 @@ -93,6 +93,7 @@ subroutine NoahmpReadNamelist(NoahmpIO) integer :: irrigation_method = 0 integer :: dvic_infiltration_option = 1 integer :: tile_drainage_option = 0 + integer :: wetland_option = 0 integer :: split_output_count = 1 logical :: skip_first_output = .false. integer :: khour = -9999 @@ -126,7 +127,7 @@ subroutine NoahmpReadNamelist(NoahmpIO) btr_option, surface_drag_option, supercooled_water_option, & frozen_soil_option, radiative_transfer_option, snow_albedo_option, & snow_thermal_conductivity, surface_runoff_option, subsurface_runoff_option, & - pcp_partition_option, tbot_option, temp_time_scheme_option, & + pcp_partition_option, tbot_option, temp_time_scheme_option, wetland_option, & glacier_option, surface_resistance_option, snow_compaction_option, & irrigation_option, irrigation_method, dvic_infiltration_option, & tile_drainage_option,soil_data_option, pedotransfer_option, crop_option, & @@ -345,6 +346,7 @@ subroutine NoahmpReadNamelist(NoahmpIO) NoahmpIO%IOPT_INFDV = dvic_infiltration_option NoahmpIO%IOPT_TDRN = tile_drainage_option NoahmpIO%IOPT_COMPACT = snow_compaction_option + NoahmpIO%IOPT_WETLAND = wetland_option ! basic model setup variables NoahmpIO%indir = indir NoahmpIO%forcing_timestep = forcing_timestep diff --git a/drivers/hrldas/NoahmpReadTableMod.F90 b/drivers/hrldas/NoahmpReadTableMod.F90 index eda242fc..b9f9a694 100644 --- a/drivers/hrldas/NoahmpReadTableMod.F90 +++ b/drivers/hrldas/NoahmpReadTableMod.F90 @@ -114,7 +114,7 @@ subroutine NoahmpReadTable(NoahmpIO) SNOWCOMPACTm_AR24,SNOWCOMPACTb_AR24,SNOWCOMPACT_P1_AR24, SNOWCOMPACT_P2_AR24,& SNOWCOMPACT_P3_AR24, SNOWCOMPACT_Up_AR24, & WSLMAX, ROUS, CMIC, SNOWDEN_MAX, CLASS_ALB_REF, CLASS_SNO_AGE, CLASS_ALB_NEW,& - PSIWLT, Z0SOIL, Z0LAKE + PSIWLT, Z0SOIL, Z0LAKE, WCAP namelist / noahmp_global_parameters / CO2, O2, TIMEAN, FSATMX, Z0SNO, SSI, SNOW_RET_FAC ,SNOW_EMIS, SWEMX, TAU0, & GRAIN_GROWTH, EXTRA_GROWTH, DIRT_SOOT, BATS_COSZ, BATS_VIS_NEW, & BATS_NIR_NEW, BATS_VIS_AGE, BATS_NIR_AGE, BATS_VIS_DIR, BATS_NIR_DIR, & @@ -123,7 +123,7 @@ subroutine NoahmpReadTable(NoahmpIO) SNOWCOMPACTm_AR24,SNOWCOMPACTb_AR24,SNOWCOMPACT_P1_AR24, SNOWCOMPACT_P2_AR24,& SNOWCOMPACT_P3_AR24, SNOWCOMPACT_Up_AR24, & WSLMAX, ROUS, CMIC, SNOWDEN_MAX, CLASS_ALB_REF, CLASS_SNO_AGE, CLASS_ALB_NEW,& - PSIWLT, Z0SOIL, Z0LAKE + PSIWLT, Z0SOIL, Z0LAKE, WCAP ! irrigation parameters integer :: IRR_HAR @@ -531,6 +531,7 @@ subroutine NoahmpReadTable(NoahmpIO) NoahmpIO%PSIWLT_TABLE = undefined_real NoahmpIO%Z0SOIL_TABLE = undefined_real NoahmpIO%Z0LAKE_TABLE = undefined_real + NoahmpIO%WCAP_TABLE = undefined_real ! irrigation parameters NoahmpIO%IRR_HAR_TABLE = undefined_int @@ -933,6 +934,7 @@ subroutine NoahmpReadTable(NoahmpIO) NoahmpIO%PSIWLT_TABLE = PSIWLT NoahmpIO%Z0SOIL_TABLE = Z0SOIL NoahmpIO%Z0LAKE_TABLE = Z0LAKE + NoahmpIO%WCAP_TABLE = WCAP !---------------- NoahmpTable.TBL irrigation parameters inquire( file='NoahmpTable.TBL', exist=file_named ) diff --git a/drivers/hrldas/WaterVarInTransferMod.F90 b/drivers/hrldas/WaterVarInTransferMod.F90 index 98324ea4..f604da63 100644 --- a/drivers/hrldas/WaterVarInTransferMod.F90 +++ b/drivers/hrldas/WaterVarInTransferMod.F90 @@ -77,6 +77,10 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO) noahmp%water%state%SoilMoistureEqui (1:NumSoilLayer) = NoahmpIO%SMOISEQ (I,1:NumSoilLayer,J) noahmp%water%state%RechargeGwDeepWT = 0.0 noahmp%water%state%RechargeGwShallowWT = 0.0 + if ( noahmp%config%nmlist%OptWetlandModel > 0 ) then + noahmp%water%state%SoilSaturateFrac = NoahmpIO%FSATXY (I,J) + noahmp%water%state%WaterStorageWetland = NoahmpIO%WSURFXY (I,J) + endif #ifdef WRF_HYDRO noahmp%water%state%WaterTableHydro = NoahmpIO%ZWATBLE2D (I,J) noahmp%water%state%WaterHeadSfc = NoahmpIO%sfcheadrt (I,J) @@ -91,6 +95,7 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO) noahmp%water%flux%TranspirationAcc = NoahmpIO%ACC_ETRANXY (I,J) noahmp%water%flux%EvapGroundNetAcc = NoahmpIO%ACC_EDIRXY (I,J) noahmp%water%flux%TranspWatLossSoilAcc(1:NumSoilLayer)= NoahmpIO%ACC_ETRANIXY(I,1:NumSoilLayer,J) + noahmp%water%flux%GlacierExcessFlowAcc = NoahmpIO%ACC_GLAFLWXY(I,J) ! water parameter variables noahmp%water%param%DrainSoilLayerInd = NoahmpIO%DRAIN_LAYER_OPT_TABLE @@ -152,6 +157,7 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO) noahmp%water%param%DrainWatDepToImperv = NoahmpIO%TD_D_TABLE(SoilType(1)) noahmp%water%param%NumSoilLayerRoot = NoahmpIO%NROOT_TABLE(VegType) noahmp%water%param%SoilDrainSlope = NoahmpIO%SLOPE_TABLE(RunoffSlopeType) + noahmp%water%param%WetlandCapMax = NoahmpIO%WCAP_TABLE do IndexSoilLayer = 1, size(SoilType) noahmp%water%param%SoilMoistureSat (IndexSoilLayer) = NoahmpIO%SMCMAX_TABLE(SoilType(IndexSoilLayer)) @@ -205,6 +211,12 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO) noahmp%water%param%DrainTubeDist = NoahmpIO%TD_SPAC (I,J) ! tile spacing endif + ! spatial varying wetland parameters from input + if ( noahmp%config%nmlist%OptWetlandModel == 2 ) then + noahmp%water%param%SoilSfcSatFracMax = NoahmpIO%FSATMX(I,J) + noahmp%water%param%WetlandCapMax = NoahmpIO%WCAP(I,J) + endif + ! derived water parameters noahmp%water%param%SoilInfilMaxCoeff = noahmp%water%param%SoilInfilFacRef * & noahmp%water%param%SoilWatConductivitySat(1) / & diff --git a/drivers/hrldas/WaterVarOutTransferMod.F90 b/drivers/hrldas/WaterVarOutTransferMod.F90 index 7a30df02..94c75c22 100644 --- a/drivers/hrldas/WaterVarOutTransferMod.F90 +++ b/drivers/hrldas/WaterVarOutTransferMod.F90 @@ -120,6 +120,7 @@ subroutine WaterVarOutTransfer(noahmp, NoahmpIO) NoahmpIO%ACC_ECANXY (I,J) = noahmp%water%flux%EvapCanopyNetAcc NoahmpIO%ACC_ETRANXY (I,J) = noahmp%water%flux%TranspirationAcc NoahmpIO%ACC_EDIRXY (I,J) = noahmp%water%flux%EvapGroundNetAcc + NoahmpIO%ACC_GLAFLWXY(I,J) = noahmp%water%flux%GlacierExcessFlowAcc NoahmpIO%RECHXY (I,J) = NoahmpIO%RECHXY(I,J) + (noahmp%water%state%RechargeGwShallowWT*1.0e3) NoahmpIO%DEEPRECHXY (I,J) = NoahmpIO%DEEPRECHXY(I,J) + noahmp%water%state%RechargeGwDeepWT NoahmpIO%SMCWTDXY (I,J) = noahmp%water%state%SoilMoistureToWT @@ -141,6 +142,12 @@ subroutine WaterVarOutTransfer(noahmp, NoahmpIO) NoahmpIO%IRFIVOL (I,J) = NoahmpIO%IRFIVOL(I,J)+(noahmp%water%flux%IrrigationRateFlood*1000.0) NoahmpIO%IRELOSS (I,J) = NoahmpIO%IRELOSS(I,J)+(noahmp%water%flux%EvapIrriSprinkler*NoahmpIO%DTBL) + ! wetland (Zhang2022) + if ( noahmp%config%nmlist%OptWetlandModel > 0 ) then + NoahmpIO%WSURFXY(I,J) = noahmp%water%state%WaterStorageWetland + NoahmpIO%FSATXY (I,J) = noahmp%water%state%SoilSaturateFrac + endif + #ifdef WRF_HYDRO NoahmpIO%infxsrt (I,J) = max(noahmp%water%flux%RunoffSurface, 0.0) ! mm, surface runoff NoahmpIO%soldrain (I,J) = max(noahmp%water%flux%RunoffSubsurface, 0.0) ! mm, underground runoff diff --git a/drivers/lis/ConfigVarInTransferMod.F90 b/drivers/lis/ConfigVarInTransferMod.F90 index 30a592ff..b0cde6a5 100644 --- a/drivers/lis/ConfigVarInTransferMod.F90 +++ b/drivers/lis/ConfigVarInTransferMod.F90 @@ -60,6 +60,7 @@ subroutine ConfigVarInTransfer(noahmp, NoahmpIO) noahmp%config%nmlist%OptRunoffSubsurface = NoahmpIO%IOPT_RUNSUB noahmp%config%nmlist%OptGlacierTreatment = NoahmpIO%IOPT_GLA noahmp%config%nmlist%OptSnowCompaction = NoahmpIO%IOPT_COMPACT + noahmp%config%nmlist%OptWetlandModel = NoahmpIO%IOPT_WETLAND ! config domain variable noahmp%config%domain%SurfaceType = 1 diff --git a/drivers/lis/LisNoahmpParamType.F90 b/drivers/lis/LisNoahmpParamType.F90 index a0678d4e..9f834401 100644 --- a/drivers/lis/LisNoahmpParamType.F90 +++ b/drivers/lis/LisNoahmpParamType.F90 @@ -210,6 +210,9 @@ module LisNoahmpParamType real(kind=kind_noahmp) :: RTCT(NSTAGE) ! fraction of carbohydrate translocation from root to grain real(kind=kind_noahmp) :: BIO2LAI ! leaf area per living leaf biomass [m2/kg] + ! wetland parameters + real(kind=kind_noahmp) :: WCAP ! maximum wetland water holding capacity [m] (tunable) for opt_wetland=1 + ! soil parameters real(kind=kind_noahmp) :: BEXP(NSOIL) ! soil B parameter real(kind=kind_noahmp) :: SMCDRY(NSOIL) ! dry soil moisture threshold diff --git a/drivers/lis/NoahmpDriverMainMod.F90 b/drivers/lis/NoahmpDriverMainMod.F90 index a5dbe856..7ede99fc 100644 --- a/drivers/lis/NoahmpDriverMainMod.F90 +++ b/drivers/lis/NoahmpDriverMainMod.F90 @@ -68,6 +68,7 @@ subroutine NoahmpDriverMain(NoahmpIO, LISparam) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 endif if ( NoahmpIO%SOIL_UPDATE_STEPS > 1 ) then @@ -81,6 +82,7 @@ subroutine NoahmpDriverMain(NoahmpIO, LISparam) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 end if endif diff --git a/drivers/lis/NoahmpIOVarInitMod.F90 b/drivers/lis/NoahmpIOVarInitMod.F90 index e2439223..1bec48c2 100644 --- a/drivers/lis/NoahmpIOVarInitMod.F90 +++ b/drivers/lis/NoahmpIOVarInitMod.F90 @@ -300,6 +300,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) if ( .not. allocated (NoahmpIO%ACC_QINSURXY)) allocate ( NoahmpIO%ACC_QINSURXY (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%ACC_QSEVAXY) ) allocate ( NoahmpIO%ACC_QSEVAXY (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%ACC_ETRANIXY)) allocate ( NoahmpIO%ACC_ETRANIXY (XSTART:XEND,1:NSOIL,YSTART:YEND) ) + if ( .not. allocated (NoahmpIO%ACC_GLAFLWXY)) allocate ( NoahmpIO%ACC_GLAFLWXY (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%FGEV_PET) ) allocate ( NoahmpIO%FGEV_PET (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%FCEV_PET) ) allocate ( NoahmpIO%FCEV_PET (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%FCTR_PET) ) allocate ( NoahmpIO%FCTR_PET (XSTART:XEND, YSTART:YEND) ) @@ -331,6 +332,16 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) if ( .not. allocated (NoahmpIO%SEASON_GDD)) allocate ( NoahmpIO%SEASON_GDD (XSTART:XEND, YSTART:YEND) ) if ( .not. allocated (NoahmpIO%CROPTYPE) ) allocate ( NoahmpIO%CROPTYPE (XSTART:XEND,5,YSTART:YEND) ) + ! Needed for Zhang et al. 2022 wetland model (OPT_WETLAND=1 or 2) + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + if ( .not. allocated (NoahmpIO%FSATXY) ) allocate ( NoahmpIO%FSATXY (XSTART:XEND, YSTART:YEND) ) ! saturated fraction of the grid (-) + if ( .not. allocated (NoahmpIO%WSURFXY)) allocate ( NoahmpIO%WSURFXY (XSTART:XEND, YSTART:YEND) ) ! wetland water storage [mm] + endif + if ( NoahmpIO%IOPT_WETLAND == 2 ) then + if ( .not. allocated (NoahmpIO%FSATMX) ) allocate ( NoahmpIO%FSATMX (XSTART:XEND, YSTART:YEND) ) ! maximum saturated fraction + if ( .not. allocated (NoahmpIO%WCAP) ) allocate ( NoahmpIO%WCAP (XSTART:XEND, YSTART:YEND) ) ! maximum wetland capacity + endif + ! Single- and Multi-layer Urban Models if ( NoahmpIO%SF_URBAN_PHYSICS > 0 ) then if ( .not. allocated (NoahmpIO%sh_urb2d) ) allocate ( NoahmpIO%sh_urb2d (XSTART:XEND,YSTART:YEND) ) @@ -656,6 +667,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) NoahmpIO%ACC_ECANXY = 0.0 NoahmpIO%ACC_ETRANXY = 0.0 NoahmpIO%ACC_EDIRXY = 0.0 + NoahmpIO%ACC_GLAFLWXY = 0.0 NoahmpIO%FGEV_PET = undefined_real NoahmpIO%FCEV_PET = undefined_real NoahmpIO%FCTR_PET = undefined_real @@ -709,6 +721,16 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO) NoahmpIO%IRRSPLH = 0.0 NoahmpIO%LOCTIM = undefined_real + ! wetland model (Zhang et al. 2022) + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + NoahmpIO%FSATXY = undefined_real + NoahmpIO%WSURFXY = undefined_real + endif + if ( NoahmpIO%IOPT_WETLAND == 2 ) then + NoahmpIO%FSATMX = undefined_real + NoahmpIO%WCAP = undefined_real + endif + ! spatial varying soil texture if ( NoahmpIO%IOPT_SOIL > 1 ) then NoahmpIO%SOILCL1 = undefined_real diff --git a/drivers/lis/NoahmpIOVarType.F90 b/drivers/lis/NoahmpIOVarType.F90 index bf7c5eab..81ee9ab4 100644 --- a/drivers/lis/NoahmpIOVarType.F90 +++ b/drivers/lis/NoahmpIOVarType.F90 @@ -58,6 +58,7 @@ module NoahmpIOVarType integer :: IOPT_INFDV ! infiltration options for dynamic VIC (1->Philip; 2-> Green-Ampt;3->Smith-Parlange) integer :: IOPT_TDRN ! drainage option (0->off; 1->simple scheme; 2->Hooghoudt's scheme) integer :: IOPT_COMPACT ! snowpack compaction (1->Anderson1976; 2->Abolafia-Rosenzweig2024) + integer :: IOPT_WETLAND ! wetland model option (0->off; 1->Zhang2022 scheme fixed parameter; 2->Zhang2022 read in 2D parameter) integer :: sf_urban_physics ! urban physics option ! basic model info @@ -316,6 +317,7 @@ module NoahmpIOVarType real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_ETRANXY ! accumulated transpiration per soil timestep [mm] real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_EDIRXY ! accumulated net ground (soil/snow) evaporation per soil timestep [mm] real(kind=kind_noahmp), allocatable, dimension(:,:,:) :: ACC_ETRANIXY ! accumualted transpiration rate within soil timestep [m/s * dt_soil/dt_main] + real(kind=kind_noahmp), allocatable, dimension(:,:) :: ACC_GLAFLWXY ! accumulated glacier excessive flow [mm] per soil timestep real(kind=kind_noahmp), allocatable, dimension(:,:) :: FGEV_PET ! ground evapo flux due to PET [W/m2] real(kind=kind_noahmp), allocatable, dimension(:,:) :: FCEV_PET ! canopy evapo flux due to PET [W/m2] real(kind=kind_noahmp), allocatable, dimension(:,:) :: FCTR_PET ! transpiration flux due to PET [W/m2] @@ -345,6 +347,13 @@ module NoahmpIOVarType real(kind=kind_noahmp) :: WTDDT = 30.0 ! frequency of groundwater call [minutes] integer :: STEPWTD ! step of groundwater call +!------------------------------------------------------------------------ +! Needed for wetland model (OPT_WETLAND = 1 or 2) +!------------------------------------------------------------------------ + real(kind=kind_noahmp), allocatable, dimension(:,:) :: FSATXY ! saturated fraction of the grid (-) + real(kind=kind_noahmp), allocatable, dimension(:,:) :: WSURFXY ! wetland water storage + real(kind=kind_noahmp), allocatable, dimension(:,:) :: FSATMX ! maximum saturated fraction + real(kind=kind_noahmp), allocatable, dimension(:,:) :: WCAP ! maximum wetland capacity !------------------------------------------------------------------------ ! Single- and Multi-layer Urban Models @@ -853,6 +862,9 @@ module NoahmpIOVarType real(kind=kind_noahmp), allocatable, dimension(:,:) :: RTCT_TABLE ! fraction of carbohydrate translocation from root to grain real(kind=kind_noahmp), allocatable, dimension(:) :: BIO2LAI_TABLE ! leaf area per living leaf biomass [m2/kg] + ! wetland parameter (OPT_WETLAND=1) + real(kind=kind_noahmp) :: WCAP_TABLE ! maximum surface wetland capacity + ! soil parameters integer :: SLCATS_TABLE ! number of soil categories real(kind=kind_noahmp), allocatable, dimension(:) :: BEXP_TABLE ! soil B parameter diff --git a/drivers/lis/NoahmpInitMainMod.F90 b/drivers/lis/NoahmpInitMainMod.F90 index 7a590791..a1ab2fe1 100644 --- a/drivers/lis/NoahmpInitMainMod.F90 +++ b/drivers/lis/NoahmpInitMainMod.F90 @@ -151,6 +151,10 @@ subroutine NoahmpInitMain(NoahmpIO) NoahmpIO%QSNOWXY(I,J) = 0.0 NoahmpIO%QRAINXY(I,J) = 0.0 NoahmpIO%WSLAKEXY(I,J) = 0.0 + if ( NoahmpIO%IOPT_WETLAND > 0 ) then + NoahmpIO%FSATXY(I,J) = 0.0 + NoahmpIO%WSURFXY(I,J) = 0.0 + endif if ( NoahmpIO%IOPT_RUNSUB /= 5 ) then NoahmpIO%WAXY(I,J) = 4900.0 NoahmpIO%WTXY(I,J) = NoahmpIO%WAXY(i,j) diff --git a/drivers/lis/NoahmpReadTableMod.F90 b/drivers/lis/NoahmpReadTableMod.F90 index f92f9647..aff96124 100644 --- a/drivers/lis/NoahmpReadTableMod.F90 +++ b/drivers/lis/NoahmpReadTableMod.F90 @@ -114,7 +114,7 @@ subroutine NoahmpReadTable(LANDDATA_name,MPTABLE_file) SNOWCOMPACTm_AR24,SNOWCOMPACTb_AR24,SNOWCOMPACT_P1_AR24, SNOWCOMPACT_P2_AR24,& SNOWCOMPACT_P3_AR24, SNOWCOMPACT_Up_AR24, & WSLMAX, ROUS, CMIC, SNOWDEN_MAX, CLASS_ALB_REF, CLASS_SNO_AGE, CLASS_ALB_NEW,& - PSIWLT, Z0SOIL, Z0LAKE + PSIWLT, Z0SOIL, Z0LAKE, WCAP namelist / noahmp_global_parameters / CO2, O2, TIMEAN, FSATMX, Z0SNO, SSI, SNOW_RET_FAC ,SNOW_EMIS, SWEMX, TAU0, & GRAIN_GROWTH, EXTRA_GROWTH, DIRT_SOOT, BATS_COSZ, BATS_VIS_NEW, & BATS_NIR_NEW, BATS_VIS_AGE, BATS_NIR_AGE, BATS_VIS_DIR, BATS_NIR_DIR, & @@ -123,7 +123,7 @@ subroutine NoahmpReadTable(LANDDATA_name,MPTABLE_file) SNOWCOMPACTm_AR24,SNOWCOMPACTb_AR24,SNOWCOMPACT_P1_AR24, SNOWCOMPACT_P2_AR24,& SNOWCOMPACT_P3_AR24, SNOWCOMPACT_Up_AR24, & WSLMAX, ROUS, CMIC, SNOWDEN_MAX, CLASS_ALB_REF, CLASS_SNO_AGE, CLASS_ALB_NEW,& - PSIWLT, Z0SOIL, Z0LAKE + PSIWLT, Z0SOIL, Z0LAKE, WCAP ! irrigation parameters integer :: IRR_HAR @@ -531,6 +531,7 @@ subroutine NoahmpReadTable(LANDDATA_name,MPTABLE_file) NoahmpIO%PSIWLT_TABLE = undefined_real NoahmpIO%Z0SOIL_TABLE = undefined_real NoahmpIO%Z0LAKE_TABLE = undefined_real + NoahmpIO%WCAP_TABLE = undefined_real ! irrigation parameters NoahmpIO%IRR_HAR_TABLE = undefined_int @@ -939,6 +940,7 @@ subroutine NoahmpReadTable(LANDDATA_name,MPTABLE_file) NoahmpIO%PSIWLT_TABLE = PSIWLT NoahmpIO%Z0SOIL_TABLE = Z0SOIL NoahmpIO%Z0LAKE_TABLE = Z0LAKE + NoahmpIO%WCAP_TABLE = WCAP !---------------- NoahmpTable.TBL irrigation parameters inquire( file=trim(MPTABLE_file), exist=file_named ) diff --git a/drivers/lis/WaterVarInTransferMod.F90 b/drivers/lis/WaterVarInTransferMod.F90 index fb8671f4..96415601 100644 --- a/drivers/lis/WaterVarInTransferMod.F90 +++ b/drivers/lis/WaterVarInTransferMod.F90 @@ -79,6 +79,10 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO, LISparam) noahmp%water%state%SoilMoistureEqui (1:NumSoilLayer) = NoahmpIO%SMOISEQ (I,1:NumSoilLayer,J) noahmp%water%state%RechargeGwDeepWT = 0.0 noahmp%water%state%RechargeGwShallowWT = 0.0 + if ( noahmp%config%nmlist%OptWetlandModel > 0 ) then + noahmp%water%state%SoilSaturateFrac = NoahmpIO%FSATXY (I,J) + noahmp%water%state%WaterStorageWetland = NoahmpIO%WSURFXY (I,J) + endif #ifdef WRF_HYDRO noahmp%water%state%WaterTableHydro = NoahmpIO%ZWATBLE2D (I,J) noahmp%water%state%WaterHeadSfc = NoahmpIO%sfcheadrt (I,J) @@ -93,6 +97,7 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO, LISparam) noahmp%water%flux%TranspirationAcc = NoahmpIO%ACC_ETRANXY (I,J) noahmp%water%flux%EvapGroundNetAcc = NoahmpIO%ACC_EDIRXY (I,J) noahmp%water%flux%TranspWatLossSoilAcc(1:NumSoilLayer)= NoahmpIO%ACC_ETRANIXY(I,1:NumSoilLayer,J) + noahmp%water%flux%GlacierExcessFlowAcc = NoahmpIO%ACC_GLAFLWXY(I,J) ! water parameter variables noahmp%water%param%DrainSoilLayerInd = LISparam%DRAIN_LAYER_OPT @@ -154,6 +159,7 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO, LISparam) noahmp%water%param%DrainWatDepToImperv = LISparam%TD_D noahmp%water%param%NumSoilLayerRoot = LISparam%NROOT noahmp%water%param%SoilDrainSlope = LISparam%SLOPE + noahmp%water%param%WetlandCapMax = LISparam%WCAP do IndexSoilLayer = 1, size(SoilType) noahmp%water%param%SoilMoistureSat (IndexSoilLayer) = LISparam%SMCMAX(IndexSoilLayer) @@ -207,6 +213,12 @@ subroutine WaterVarInTransfer(noahmp, NoahmpIO, LISparam) noahmp%water%param%DrainTubeDist = NoahmpIO%TD_SPAC (I,J) ! tile spacing endif + ! spatial varying wetland parameters from input + if ( noahmp%config%nmlist%OptWetlandModel == 2 ) then + noahmp%water%param%SoilSfcSatFracMax = NoahmpIO%FSATMX(I,J) + noahmp%water%param%WetlandCapMax = NoahmpIO%WCAP(I,J) + endif + ! derived water parameters noahmp%water%param%SoilInfilMaxCoeff = noahmp%water%param%SoilInfilFacRef * & noahmp%water%param%SoilWatConductivitySat(1) / & diff --git a/drivers/lis/WaterVarOutTransferMod.F90 b/drivers/lis/WaterVarOutTransferMod.F90 index 3682e8b5..e14851b6 100644 --- a/drivers/lis/WaterVarOutTransferMod.F90 +++ b/drivers/lis/WaterVarOutTransferMod.F90 @@ -118,6 +118,7 @@ subroutine WaterVarOutTransfer(noahmp, NoahmpIO) NoahmpIO%ACC_ECANXY (I,J) = noahmp%water%flux%EvapCanopyNetAcc NoahmpIO%ACC_ETRANXY (I,J) = noahmp%water%flux%TranspirationAcc NoahmpIO%ACC_EDIRXY (I,J) = noahmp%water%flux%EvapGroundNetAcc + NoahmpIO%ACC_GLAFLWXY(I,J) = noahmp%water%flux%GlacierExcessFlowAcc NoahmpIO%RECHXY (I,J) = NoahmpIO%RECHXY(I,J) + (noahmp%water%state%RechargeGwShallowWT*1.0e3) NoahmpIO%DEEPRECHXY (I,J) = NoahmpIO%DEEPRECHXY(I,J) + noahmp%water%state%RechargeGwDeepWT NoahmpIO%SMCWTDXY (I,J) = noahmp%water%state%SoilMoistureToWT @@ -142,6 +143,12 @@ subroutine WaterVarOutTransfer(noahmp, NoahmpIO) NoahmpIO%IRFIVOL (I,J) = NoahmpIO%IRFIVOL(I,J)+(noahmp%water%flux%IrrigationRateFlood*1000.0) NoahmpIO%IRELOSS (I,J) = NoahmpIO%IRELOSS(I,J)+(noahmp%water%flux%EvapIrriSprinkler*NoahmpIO%DTBL) + ! wetland (Zhang2022) + if ( noahmp%config%nmlist%OptWetlandModel > 0 ) then + NoahmpIO%WSURFXY(I,J) = noahmp%water%state%WaterStorageWetland + NoahmpIO%FSATXY (I,J) = noahmp%water%state%SoilSaturateFrac + endif + #ifdef WRF_HYDRO NoahmpIO%infxsrt (I,J) = max(noahmp%water%flux%RunoffSurface, 0.0) ! mm, surface runoff NoahmpIO%soldrain (I,J) = max(noahmp%water%flux%RunoffSubsurface, 0.0) ! mm, underground runoff diff --git a/parameters/NoahmpTable.TBL b/parameters/NoahmpTable.TBL index 59a3fb44..cad98ba4 100644 --- a/parameters/NoahmpTable.TBL +++ b/parameters/NoahmpTable.TBL @@ -476,6 +476,8 @@ PSIWLT = -150.0 ! metric potential for wilting point (m) Z0SOIL = 0.002 ! Bare-soil roughness length (m) (i.e., under the canopy) Z0LAKE = 0.01 ! Lake surface roughness length (m) + ! wetland parameter + WCAP = 0.10 ! maximum wetland water holding capacity [m] (tunable) from Zhang et al. 2022 / &noahmp_irrigation_parameters diff --git a/src/BalanceErrorCheckMod.F90 b/src/BalanceErrorCheckMod.F90 index f076e2a5..f698821e 100644 --- a/src/BalanceErrorCheckMod.F90 +++ b/src/BalanceErrorCheckMod.F90 @@ -36,13 +36,14 @@ subroutine BalanceWaterInit(noahmp) SnowWaterEquiv => noahmp%water%state%SnowWaterEquiv ,& ! in, snow water equivalent [mm] SoilMoisture => noahmp%water%state%SoilMoisture ,& ! in, total soil moisture [m3/m3] WaterStorageAquifer => noahmp%water%state%WaterStorageAquifer ,& ! in, water storage in aquifer [mm] + WaterStorageWetland => noahmp%water%state%WaterStorageWetland ,& ! in, water storage in wetland [mm] WaterStorageTotBeg => noahmp%water%state%WaterStorageTotBeg & ! out, total water storage [mm] at the beginning ) ! ---------------------------------------------------------------------- ! compute total water storage before NoahMP processes if ( SurfaceType == 1 ) then ! soil - WaterStorageTotBeg = CanopyLiqWater + CanopyIce + SnowWaterEquiv + WaterStorageAquifer + WaterStorageTotBeg = CanopyLiqWater + CanopyIce + SnowWaterEquiv + WaterStorageAquifer + WaterStorageWetland do LoopInd = 1, NumSoilLayer WaterStorageTotBeg = WaterStorageTotBeg + SoilMoisture(LoopInd) * ThicknessSnowSoilLayer(LoopInd) * 1000.0 enddo @@ -87,6 +88,7 @@ subroutine BalanceWaterCheck(noahmp) SnowWaterEquiv => noahmp%water%state%SnowWaterEquiv ,& ! in, snow water equivalent [mm] SoilMoisture => noahmp%water%state%SoilMoisture ,& ! in, total soil moisture [m3/m3] WaterStorageAquifer => noahmp%water%state%WaterStorageAquifer ,& ! in, water storage in aquifer [mm] + WaterStorageWetland => noahmp%water%state%WaterStorageWetland ,& ! in, water storage in wetland [mm] WaterStorageTotBeg => noahmp%water%state%WaterStorageTotBeg ,& ! in, total water storage [mm] at the beginning PrecipTotRefHeight => noahmp%water%flux%PrecipTotRefHeight ,& ! in, total precipitation [mm/s] at reference height EvapCanopyNet => noahmp%water%flux%EvapCanopyNet ,& ! in, evaporation of intercepted water [mm/s] @@ -116,7 +118,7 @@ subroutine BalanceWaterCheck(noahmp) ! only water balance check for every soil timestep ! Error in water balance should be < 0.1 mm if ( SurfaceType == 1 ) then ! soil - WaterStorageTotEnd = CanopyLiqWater + CanopyIce + SnowWaterEquiv + WaterStorageAquifer + WaterStorageTotEnd = CanopyLiqWater + CanopyIce + SnowWaterEquiv + WaterStorageAquifer + WaterStorageWetland do LoopInd = 1, NumSoilLayer WaterStorageTotEnd = WaterStorageTotEnd + SoilMoisture(LoopInd) * ThicknessSnowSoilLayer(LoopInd) * 1000.0 enddo @@ -143,11 +145,11 @@ subroutine BalanceWaterCheck(noahmp) write(*, & '(" GridIndexI GridIndexJ SfcWaterTotChgAcc PrecipTotRefHeightAcc IrrigationRateMicro & IrrigationRateFlood EvapCanopyNetAcc EvapGroundNetAcc TranspirationAcc RunoffSurface & - RunoffSubsurface WaterTableDepth TileDrain")') - write(*,'(i6,i6,f10.3,10f10.5)') GridIndexI, GridIndexJ, SfcWaterTotChgAcc, PrecipTotAcc, & + RunoffSubsurface WaterTableDepth TileDrain WaterStorageWetland ")') + write(*,'(i6,i6,f10.3,11f10.5)') GridIndexI, GridIndexJ, SfcWaterTotChgAcc, PrecipTotAcc, & IrrigationRateMicro*1000.0, IrrigationRateFlood*1000.0, & EvapCanopyNetAcc, EvapGroundNetAcc, TranspirationAcc, RunoffSurface, & - RunoffSubsurface, WaterTableDepth, TileDrain + RunoffSubsurface, WaterTableDepth, TileDrain, WaterStorageWetland stop "Error: Water budget problem in NoahMP LSM" endif #endif diff --git a/src/ConfigVarInitMod.F90 b/src/ConfigVarInitMod.F90 index 99a96e35..af9ec5b3 100644 --- a/src/ConfigVarInitMod.F90 +++ b/src/ConfigVarInitMod.F90 @@ -47,11 +47,13 @@ subroutine ConfigVarInitDefault(noahmp) noahmp%config%nmlist%OptPedotransfer = undefined_int noahmp%config%nmlist%OptGlacierTreatment = undefined_int noahmp%config%nmlist%OptSnowCompaction = undefined_int + noahmp%config%nmlist%OptWetlandModel = undefined_int ! config domain variable noahmp%config%domain%LandUseDataName = "MODIFIED_IGBP_MODIS_NOAH" noahmp%config%domain%FlagUrban = .false. noahmp%config%domain%FlagCropland = .false. + noahmp%config%domain%FlagWetland = .false. noahmp%config%domain%FlagDynamicCrop = .false. noahmp%config%domain%FlagDynamicVeg = .false. noahmp%config%domain%FlagSoilProcess = .false. diff --git a/src/ConfigVarType.F90 b/src/ConfigVarType.F90 index 25bfb0d6..31dbc861 100644 --- a/src/ConfigVarType.F90 +++ b/src/ConfigVarType.F90 @@ -123,7 +123,10 @@ module ConfigVarType integer :: OptSnowCompaction ! options for ground snow compaction ! 1 -> original scheme from Anderson (1976) ! 2 -> new scheme from Abolafia-Rosenzweig et al. (2024) - + integer :: OptWetlandModel ! option for wetland model + ! 0 -> No Wetland model (default) + ! 1 -> Single-point/uniform parameter (Zhang, et al. 2022 WRR) + ! 2 -> 2-D regional parameter input (Zhang, et al. 2022 WRR) end type namelist_type @@ -133,6 +136,7 @@ module ConfigVarType character(len=256) :: LandUseDataName ! landuse dataset name (USGS or MODIFIED_IGBP_MODIS_NOAH) logical :: FlagUrban ! flag for urban grid logical :: FlagCropland ! flag to identify croplands + logical :: FlagWetland ! flag to identify wetlands logical :: FlagDynamicCrop ! flag to activate dynamic crop model logical :: FlagDynamicVeg ! flag to activate dynamic vegetation scheme logical :: FlagSoilProcess ! flag to determine if calculating soil processes diff --git a/src/GeneralInitMod.F90 b/src/GeneralInitMod.F90 index 551e0176..f26fe71f 100644 --- a/src/GeneralInitMod.F90 +++ b/src/GeneralInitMod.F90 @@ -27,12 +27,16 @@ subroutine GeneralInit(noahmp) ! -------------------------------------------------------------------- associate( & + LandUseDataName => noahmp%config%domain%LandUseDataName ,& ! in, landuse data name (USGS or MODIS_IGBP) + VegType => noahmp%config%domain%VegType ,& ! in, vegetation type NumSoilLayer => noahmp%config%domain%NumSoilLayer ,& ! in, number of soil layers DepthSoilLayer => noahmp%config%domain%DepthSoilLayer ,& ! in, depth [m] of layer-bottom from soil surface NumSoilLayerRoot => noahmp%water%param%NumSoilLayerRoot ,& ! in, number of soil layers with root present NumSnowLayerNeg => noahmp%config%domain%NumSnowLayerNeg ,& ! in, actual number of snow layers (negative) DepthSnowSoilLayer => noahmp%config%domain%DepthSnowSoilLayer ,& ! in, depth of snow/soil layer-bottom [m] TemperatureSoilSnow => noahmp%energy%state%TemperatureSoilSnow ,& ! in, snow and soil layer temperature [K] + FlagCropland => noahmp%config%domain%FlagCropland ,& ! out, flag to identify croplands + FlagWetland => noahmp%config%domain%FlagWetland ,& ! out, flag to identify wetlands ThicknessSnowSoilLayer => noahmp%config%domain%ThicknessSnowSoilLayer ,& ! out, thickness of snow/soil layers [m] TemperatureRootZone => noahmp%energy%state%TemperatureRootZone & ! out, root-zone averaged temperature [K] ) @@ -54,6 +58,17 @@ subroutine GeneralInit(noahmp) TemperatureSoilSnow(LoopInd) * ThicknessSnowSoilLayer(LoopInd) / (-DepthSoilLayer(NumSoilLayerRoot)) enddo + ! initialize special land type flags + FlagCropland = .false. + FlagWetland = .false. + if ( trim(LandUseDataName) == "USGS" ) then + if ( (VegType >= 3 ) .and. (VegType <= 6 ) ) FlagCropland = .true. + if ( (VegType >= 17) .and. (VegType <= 18) ) FlagWetland = .true. + elseif ( trim(LandUseDataName) == "MODIFIED_IGBP_MODIS_NOAH") then + if ( (VegType == 12) .or. (VegType == 14) ) FlagCropland = .true. + if ( (VegType == 11) ) FlagWetland = .true. + endif + end associate end subroutine GeneralInit diff --git a/src/IrrigationPrepareMod.F90 b/src/IrrigationPrepareMod.F90 index 108bbe68..2996f107 100644 --- a/src/IrrigationPrepareMod.F90 +++ b/src/IrrigationPrepareMod.F90 @@ -27,8 +27,7 @@ subroutine IrrigationPrepare(noahmp) ! ---------------------------------------------------------------------- associate( & - LandUseDataName => noahmp%config%domain%LandUseDataName ,& ! in, landuse data name (USGS or MODIS_IGBP) - VegType => noahmp%config%domain%VegType ,& ! in, vegetation type + FlagCropland => noahmp%config%domain%FlagCropland ,& ! in, flag to identify croplands FlagSoilProcess => noahmp%config%domain%FlagSoilProcess ,& ! in, flag to calculate soil processes OptIrrigationMethod => noahmp%config%nmlist%OptIrrigationMethod ,& ! in, irrigation method option IrriFracThreshold => noahmp%water%param%IrriFracThreshold ,& ! in, irrigation fraction threshold @@ -38,23 +37,12 @@ subroutine IrrigationPrepare(noahmp) IrrigationAmtFlood => noahmp%water%state%IrrigationAmtFlood ,& ! inout, flood irrigation water amount [m] IrrigationAmtMicro => noahmp%water%state%IrrigationAmtMicro ,& ! inout, micro irrigation water amount [m] RainfallRefHeight => noahmp%water%flux%RainfallRefHeight ,& ! inout, rainfall [mm/s] at reference height - FlagCropland => noahmp%config%domain%FlagCropland ,& ! out, flag to identify croplands IrrigationFracSprinkler => noahmp%water%state%IrrigationFracSprinkler ,& ! out, sprinkler irrigation fraction (0 to 1) IrrigationFracMicro => noahmp%water%state%IrrigationFracMicro ,& ! out, fraction of grid under micro irrigation (0 to 1) IrrigationFracFlood => noahmp%water%state%IrrigationFracFlood & ! out, fraction of grid under flood irrigation (0 to 1) ) ! ---------------------------------------------------------------------- - ! initialize - FlagCropland = .false. - - ! determine cropland - if ( trim(LandUseDataName) == "USGS" ) then - if ( (VegType >= 3) .and. (VegType <= 6) ) FlagCropland = .true. - elseif ( trim(LandUseDataName) == "MODIFIED_IGBP_MODIS_NOAH") then - if ( (VegType == 12) .or. (VegType == 14) ) FlagCropland = .true. - endif - ! if OptIrrigationMethod = 0 and if methods are unknown for certain area, then use sprinkler irrigation method if ( (OptIrrigationMethod == 0) .and. (IrrigationFracSprinkler == 0.0) .and. (IrrigationFracMicro == 0.0) & .and. (IrrigationFracFlood == 0.0) .and. (IrrigationFracGrid >= IrriFracThreshold) ) then diff --git a/src/Makefile b/src/Makefile index f9fb49d7..da2d8569 100644 --- a/src/Makefile +++ b/src/Makefile @@ -36,6 +36,8 @@ OBJS = ConstantDefineMod.o \ RunoffSurfaceTopModelMmfMod.o \ RunoffSurfaceVicMod.o \ RunoffSurfaceXinAnJiangMod.o \ + RunoffSurfaceWetlandMod.o \ + WetlandWaterZhang22Mod.o \ ShallowWaterTableMmfMod.o \ SnowfallBelowCanopyMod.o \ SnowLayerCombineMod.o \ @@ -198,6 +200,8 @@ RunoffSurfaceTopModelGrdMod.o: ../utility/Machine.o NoahmpVarType.o Const RunoffSurfaceTopModelMmfMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o RunoffSurfaceVicMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o RunoffSurfaceXinAnJiangMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o +RunoffSurfaceWetlandMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o +WetlandWaterZhang22Mod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o ShallowWaterTableMmfMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o SnowfallBelowCanopyMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o SnowLayerCombineMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o \ @@ -231,7 +235,8 @@ SoilWaterMainMod.o: ../utility/Machine.o NoahmpVarType.o Const RunoffSubSurfaceEquiWaterTableMod.o RunoffSubSurfaceGroundWaterMod.o \ RunoffSubSurfaceDrainageMod.o RunoffSubSurfaceShallowMmfMod.o \ SoilWaterDiffusionRichardsMod.o SoilMoistureSolverMod.o \ - TileDrainageSimpleMod.o TileDrainageHooghoudtMod.o + TileDrainageSimpleMod.o TileDrainageHooghoudtMod.o \ + RunoffSurfaceWetlandMod.o TileDrainageEquiDepthMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o TileDrainageHooghoudtMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o \ TileDrainageEquiDepthMod.o WaterTableDepthSearchMod.o \ @@ -239,7 +244,7 @@ TileDrainageHooghoudtMod.o: ../utility/Machine.o NoahmpVarType.o Const TileDrainageSimpleMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o WaterMainMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o \ CanopyHydrologyMod.o SnowWaterMainMod.o IrrigationFloodMod.o \ - IrrigationMicroMod.o SoilWaterMainMod.o + IrrigationMicroMod.o SoilWaterMainMod.o WetlandWaterZhang22Mod.o WaterTableDepthSearchMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o WaterTableEquilibriumMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o IrrigationTriggerMod.o: ../utility/Machine.o NoahmpVarType.o ConstantDefineMod.o diff --git a/src/RunoffSurfaceWetlandMod.F90 b/src/RunoffSurfaceWetlandMod.F90 new file mode 100644 index 00000000..f66e5dcd --- /dev/null +++ b/src/RunoffSurfaceWetlandMod.F90 @@ -0,0 +1,50 @@ +module RunoffSurfaceWetlandMod + +!!! Calculate surface runoff based on first-layer soil moisture with wetland scheme (Zhang et al., 2022) + + use Machine + use NoahmpVarType + use ConstantDefineMod + + implicit none + +contains + + subroutine RunoffSurfaceWetland(noahmp) + +! ------------------------ Code history -------------------------------------------------- +! Refactered code: C. He, P. Valayamkunnath, & refactor team (He et al. 2023) +! Implemented by Z. Zhang (Zhang et al. 2022) +! ---------------------------------------------------------------------------------------- + + implicit none + + type(noahmp_type), intent(inout) :: noahmp + +! -------------------------------------------------------------------- + associate( & + SoilSfcInflowMean => noahmp%water%flux%SoilSfcInflowMean ,& ! in, mean water input on soil surface [m/s] + SoilSfcSatFracMax => noahmp%water%param%SoilSfcSatFracMax ,& ! in, maximum surface saturated fraction (global mean) + SoilImpervFrac => noahmp%water%state%SoilImpervFrac ,& ! in, impervious fraction due to frozen soil + SoilLiqWater => noahmp%water%state%SoilLiqWater ,& ! in, soil water content [m3/m3] + SoilMoistureSat => noahmp%water%param%SoilMoistureSat ,& ! in, saturated value of soil moisture [m3/m3] + SoilSaturateFrac => noahmp%water%state%SoilSaturateFrac ,& ! out, fractional saturated area for soil moisture + RunoffSurface => noahmp%water%flux%RunoffSurface ,& ! out, surface runoff [m/s] + InfilRateSfc => noahmp%water%flux%InfilRateSfc & ! out, infiltration rate at surface [m/s] + ) +! ---------------------------------------------------------------------- + + ! compute saturated area fraction, based on top-layer soil moisture + SoilSaturateFrac = SoilSfcSatFracMax * (SoilLiqWater(1)/SoilMoistureSat(1)) + + ! compute surface runoff and infiltration m/s + if ( SoilSfcInflowMean > 0.0 ) then + RunoffSurface = SoilSfcInflowMean * ((1.0-SoilImpervFrac(1)) * SoilSaturateFrac + SoilImpervFrac(1)) + InfilRateSfc = SoilSfcInflowMean - RunoffSurface + endif + + end associate + + end subroutine RunoffSurfaceWetland + +end module RunoffSurfaceWetlandMod diff --git a/src/SoilWaterMainMod.F90 b/src/SoilWaterMainMod.F90 index a03a983b..cae67801 100644 --- a/src/SoilWaterMainMod.F90 +++ b/src/SoilWaterMainMod.F90 @@ -14,6 +14,7 @@ module SoilWaterMainMod use RunoffSurfaceVicMod, only : RunoffSurfaceVIC use RunoffSurfaceXinAnJiangMod, only : RunoffSurfaceXinAnJiang use RunoffSurfaceDynamicVicMod, only : RunoffSurfaceDynamicVic + use RunoffSurfaceWetlandMod, only : RunoffSurfaceWetland use RunoffSubSurfaceEquiWaterTableMod, only : RunoffSubSurfaceEquiWaterTable use RunoffSubSurfaceGroundWaterMod, only : RunoffSubSurfaceGroundWater use RunoffSubSurfaceDrainageMod, only : RunoffSubSurfaceDrainage @@ -65,9 +66,11 @@ subroutine SoilWaterMain(noahmp) SoilTimeStep => noahmp%config%domain%SoilTimeStep ,& ! in, noahmp soil time step [s] ThicknessSnowSoilLayer => noahmp%config%domain%ThicknessSnowSoilLayer ,& ! in, thickness of snow/soil layers [m] FlagUrban => noahmp%config%domain%FlagUrban ,& ! in, logical flag for urban grid + FlagWetland => noahmp%config%domain%FlagWetland ,& ! in, logical flag for wetland grid OptRunoffSurface => noahmp%config%nmlist%OptRunoffSurface ,& ! in, options for surface runoff OptRunoffSubsurface => noahmp%config%nmlist%OptRunoffSubsurface ,& ! in, options for subsurface runoff OptTileDrainage => noahmp%config%nmlist%OptTileDrainage ,& ! in, options for tile drainage + OptWetlandModel => noahmp%config%nmlist%OptWetlandModel ,& ! in, options for wetland model SoilIce => noahmp%water%state%SoilIce ,& ! in, soil ice content [m3/m3] TileDrainFrac => noahmp%water%state%TileDrainFrac ,& ! in, tile drainage map (fraction) SoilSfcInflowMean => noahmp%water%flux%SoilSfcInflowMean ,& ! in, mean water input on soil surface [m/s] @@ -149,6 +152,10 @@ subroutine SoilWaterMain(noahmp) if ( OptRunoffSurface == 7 ) call RunoffSurfaceXinAnJiang(noahmp,SoilTimeStep) if ( OptRunoffSurface == 8 ) call RunoffSurfaceDynamicVic(noahmp,SoilTimeStep,InfilSfcAcc) + ! special treatment for wetland points (due to subgrid wetland treatment, currently no flag control) + !if ( (FlagWetland .eqv. .true.) .and. (OptWetlandModel > 0) ) call RunoffSurfaceWetland(noahmp) + if ( OptWetlandModel > 0 ) call RunoffSurfaceWetland(noahmp) + ! determine iteration times to solve soil water diffusion and moisture NumIterSoilWat = 3 if ( (InfilRateSfc*SoilTimeStep) > (ThicknessSnowSoilLayer(1)*SoilMoistureSat(1)) ) then diff --git a/src/WaterMainMod.F90 b/src/WaterMainMod.F90 index d737e81e..988961f5 100644 --- a/src/WaterMainMod.F90 +++ b/src/WaterMainMod.F90 @@ -6,11 +6,12 @@ module WaterMainMod use Machine use NoahmpVarType use ConstantDefineMod - use CanopyHydrologyMod, only : CanopyHydrology - use SnowWaterMainMod, only : SnowWaterMain - use IrrigationFloodMod, only : IrrigationFlood - use IrrigationMicroMod, only : IrrigationMicro - use SoilWaterMainMod, only : SoilWaterMain + use CanopyHydrologyMod, only : CanopyHydrology + use SnowWaterMainMod, only : SnowWaterMain + use IrrigationFloodMod, only : IrrigationFlood + use IrrigationMicroMod, only : IrrigationMicro + use SoilWaterMainMod, only : SoilWaterMain + use WetlandWaterZhang22Mod, only : WetlandWaterZhang22 implicit none @@ -37,9 +38,11 @@ subroutine WaterMain(noahmp) SoilTimeStep => noahmp%config%domain%SoilTimeStep ,& ! in, soil process timestep [s] SurfaceType => noahmp%config%domain%SurfaceType ,& ! in, surface type 1-soil; 2-lake FlagCropland => noahmp%config%domain%FlagCropland ,& ! in, flag to identify croplands + FlagWetland => noahmp%config%domain%FlagWetland ,& ! in, flag to identify wetlands FlagUrban => noahmp%config%domain%FlagUrban ,& ! in, urban point flag FlagSoilProcess => noahmp%config%domain%FlagSoilProcess ,& ! in, flag to calculate soil processes NumSoilTimeStep => noahmp%config%domain%NumSoilTimeStep ,& ! in, number of timesteps for soil process calculation + OptWetlandModel => noahmp%config%nmlist%OptWetlandModel ,& ! in, options for wetland model VaporizeGrd => noahmp%water%flux%VaporizeGrd ,& ! in, ground vaporize rate total (evap+sublim) [mm/s] CondenseVapGrd => noahmp%water%flux%CondenseVapGrd ,& ! in, ground vapor condense rate total (dew+frost) [mm/s] RainfallGround => noahmp%water%flux%RainfallGround ,& ! in, ground surface rain rate [mm/s] @@ -71,6 +74,7 @@ subroutine WaterMain(noahmp) SublimSnowSfcIce => noahmp%water%flux%SublimSnowSfcIce ,& ! inout, snow surface sublimation rate[mm/s] TranspWatLossSoil => noahmp%water%flux%TranspWatLossSoil ,& ! inout, transpiration water loss from soil layers [m/s] GlacierExcessFlow => noahmp%water%flux%GlacierExcessFlow ,& ! inout, glacier excess flow [mm/s] + GlacierExcessFlowAcc => noahmp%water%flux%GlacierExcessFlowAcc ,& ! inout, accumulated glacier excess flow [mm] SoilSfcInflowAcc => noahmp%water%flux%SoilSfcInflowAcc ,& ! inout, accumulated water flux into soil during soil timestep [m/s * dt_soil/dt_main] EvapSoilSfcLiqAcc => noahmp%water%flux%EvapSoilSfcLiqAcc ,& ! inout, accumulated soil surface evaporation during soil timestep [m/s * dt_soil/dt_main] TranspWatLossSoilAcc => noahmp%water%flux%TranspWatLossSoilAcc ,& ! inout, accumualted transpiration water loss during soil timestep [m/s * dt_soil/dt_main] @@ -94,10 +98,9 @@ subroutine WaterMain(noahmp) ! initialize TranspWatLossSoil = 0.0 - GlacierExcessFlow = 0.0 - RunoffSubsurface = 0.0 - RunoffSurface = 0.0 SoilSfcInflow = 0.0 + RunoffSurface = 0.0 + RunoffSubsurface = 0.0 TileDrain = 0.0 ! prepare for water process @@ -128,6 +131,9 @@ subroutine WaterMain(noahmp) ! snowpack water processs call SnowWaterMain(noahmp) + ! accumulate glacier excessive flow [mm] + GlacierExcessFlowAcc = GlacierExcessFlowAcc + GlacierExcessFlow * MainTimeStep + ! treat frozen ground/soil if ( FlagFrozenGround .eqv. .true. ) then SoilIce(1) = SoilIce(1) + (DewSoilSfcLiq-EvapSoilSfcLiq) * MainTimeStep / & @@ -190,10 +196,10 @@ subroutine WaterMain(noahmp) call SoilWaterMain(noahmp) endif - endif ! FlagSoilProcess soil timestep + ! merge excess glacier snow flow to subsurface runoff + RunoffSubsurface = RunoffSubsurface + GlacierExcessFlowAcc ! mm per soil timestep - ! merge excess glacier snow flow to subsurface runoff - RunoffSubsurface = RunoffSubsurface + GlacierExcessFlow * MainTimeStep ! mm per soil timestep + endif ! FlagSoilProcess soil timestep ! update surface water vapor flux ! urban - jref WaterToAtmosTotal = Transpiration + EvapCanopyNet + EvapGroundNet @@ -202,6 +208,12 @@ subroutine WaterMain(noahmp) SpecHumidity2mBare = SpecHumiditySfc endif + ! call surface wetland scheme (due to subgrid wetland treatment, currently no flag control) + !if ( (FlagWetland .eqv. .true.) .and. (OptWetlandModel > 0) ) then + if ( OptWetlandModel > 0 ) then + call WetlandWaterZhang22(noahmp,MainTimeStep) + endif + end associate end subroutine WaterMain diff --git a/src/WaterVarInitMod.F90 b/src/WaterVarInitMod.F90 index f60af81b..c07122c9 100644 --- a/src/WaterVarInitMod.F90 +++ b/src/WaterVarInitMod.F90 @@ -75,6 +75,7 @@ subroutine WaterVarInitDefault(noahmp) noahmp%water%state%PrecipAreaFrac = undefined_real noahmp%water%state%TileDrainFrac = undefined_real noahmp%water%state%FrozenPrecipFrac = undefined_real + noahmp%water%state%WaterStorageWetland = undefined_real if ( .not. allocated(noahmp%water%state%IndexPhaseChange) ) & allocate( noahmp%water%state%IndexPhaseChange(-NumSnowLayerMax+1:NumSoilLayer) ) @@ -181,6 +182,7 @@ subroutine WaterVarInitDefault(noahmp) noahmp%water%flux%WaterToAtmosTotal = undefined_real noahmp%water%flux%EvapSoilSfcLiqAcc = undefined_real noahmp%water%flux%SoilSfcInflowAcc = undefined_real + noahmp%water%flux%GlacierExcessFlowAcc = undefined_real noahmp%water%flux%SfcWaterTotChgAcc = undefined_real noahmp%water%flux%PrecipTotAcc = undefined_real noahmp%water%flux%EvapCanopyNetAcc = undefined_real @@ -282,6 +284,7 @@ subroutine WaterVarInitDefault(noahmp) noahmp%water%param%SoilMatPotentialWilt = undefined_real noahmp%water%param%SnowMeltFac = undefined_real noahmp%water%param%SnowCoverFac = undefined_real + noahmp%water%param%WetlandCapMax = undefined_real if ( .not. allocated(noahmp%water%param%SoilMoistureSat) ) & allocate( noahmp%water%param%SoilMoistureSat(1:NumSoilLayer) ) diff --git a/src/WaterVarType.F90 b/src/WaterVarType.F90 index ed788c18..b23d7966 100644 --- a/src/WaterVarType.F90 +++ b/src/WaterVarType.F90 @@ -70,6 +70,7 @@ module WaterVarType real(kind=kind_noahmp) :: EvapCanopyNetAcc ! accumulated net evaporation of canopy intercepted water per soil timestep [mm] real(kind=kind_noahmp) :: TranspirationAcc ! accumulated transpiration per soil timestep [mm] real(kind=kind_noahmp) :: EvapGroundNetAcc ! accumulated net ground (soil/snow) evaporation per soil timestep [mm] + real(kind=kind_noahmp) :: GlacierExcessFlowAcc ! accumulated glacier excessive flow [mm] per soil timestep real(kind=kind_noahmp) :: EvapSoilSfcLiqMean ! mean soil surface water evaporation during soil timestep [m/s] real(kind=kind_noahmp) :: SoilSfcInflowMean ! mean water input on soil surface during soil timestep [m/s] @@ -123,6 +124,7 @@ module WaterVarType real(kind=kind_noahmp) :: WaterStorageAquifer ! water storage in aquifer [mm] real(kind=kind_noahmp) :: WaterStorageSoilAqf ! water storage in aquifer + saturated soil [mm] real(kind=kind_noahmp) :: WaterStorageLake ! water storage in lake (can be negative) [mm] + real(kind=kind_noahmp) :: WaterStorageWetland ! water storage in wetland [mm] real(kind=kind_noahmp) :: WaterHeadSfc ! surface water head [mm] real(kind=kind_noahmp) :: IrrigationFracGrid ! total irrigation fraction from input for a grid real(kind=kind_noahmp) :: PrecipAreaFrac ! fraction of the gridcell that receives precipitation @@ -225,6 +227,7 @@ module WaterVarType real(kind=kind_noahmp) :: SoilMatPotentialWilt ! soil metric potential for wilting point [m] real(kind=kind_noahmp) :: SnowMeltFac ! snowmelt m parameter in snow cover fraction calculation real(kind=kind_noahmp) :: SnowCoverFac ! snow cover factor [m] (originally hard-coded 2.5*z0 in SCF formulation) + real(kind=kind_noahmp) :: WetlandCapMax ! maximum wetland capacity [m] real(kind=kind_noahmp), allocatable, dimension(:) :: SoilMoistureSat ! saturated value of soil moisture [m3/m3] real(kind=kind_noahmp), allocatable, dimension(:) :: SoilMoistureWilt ! wilting point soil moisture [m3/m3] diff --git a/src/WetlandWaterZhang22Mod.F90 b/src/WetlandWaterZhang22Mod.F90 new file mode 100644 index 00000000..d223de31 --- /dev/null +++ b/src/WetlandWaterZhang22Mod.F90 @@ -0,0 +1,102 @@ +module WetlandWaterZhang22Mod + +!!! Calculate wetland water processes (Zhang et al., 2022), evaporation is estimated based on Priestley Taylor (P-T) Method + + use Machine + use NoahmpVarType + use ConstantDefineMod + + implicit none + +contains + + subroutine WetlandWaterZhang22(noahmp,TimeStep) + +! ------------------------ Code history -------------------------------------------------- +! Refactered code: C. He, P. Valayamkunnath, & refactor team (He et al. 2023) +! Implementation by Z. Zhang (Zhang et al. 2022) +! ---------------------------------------------------------------------------------------- + + implicit none + + type(noahmp_type), intent(inout) :: noahmp + real(kind=kind_noahmp), intent(in) :: TimeStep ! timestep (may not be the same as main model timestep) + +! local variables + real(kind=kind_noahmp) :: LatHeatSpec ! latent heat vap./sublimation [J/kg] + real(kind=kind_noahmp) :: VapPresSlope ! slope of saturation vapor pressure curve + real(kind=kind_noahmp) :: LatHeatPot ! evaporation from P-T method, energy flux [W/m2] + real(kind=kind_noahmp) :: PsychroConst ! gamma psychrometic constant + real(kind=kind_noahmp) :: EvapLatentHeat ! evaporation heat from surface [W/m2] + real(kind=kind_noahmp) :: EvapWaterFlux ! evaporation water from surface [mm/s] + +! -------------------------------------------------------------------- + associate( & + TemperatureSfc => noahmp%energy%state%TemperatureSfc ,& ! in, surface air temperature [K] + RadSwAbsSfc => noahmp%energy%flux%RadSwAbsSfc ,& ! in, total absorbed solar radiation [W/m2] + RadSwReflSfc => noahmp%energy%flux%RadSwReflSfc ,& ! in, total reflected solar radiation [W/m2] + RadLwNetSfc => noahmp%energy%flux%RadLwNetSfc ,& ! in, total net longwave rad [W/m2] (+ to atm) + HeatGroundTot => noahmp%energy%flux%HeatGroundTot ,& ! in, total ground heat flux [W/m2] (+ to soil/snow) + HeatSensibleSfc => noahmp%energy%flux%HeatSensibleSfc ,& ! in, total sensible heat [W/m2] (+ to atm) + SoilSaturateFrac => noahmp%water%state%SoilSaturateFrac ,& ! in, fractional saturated area for soil moisture + WetlandCapMax => noahmp%water%param%WetlandCapMax ,& ! in, maximum wetland capacity [m] + WaterStorageWetland => noahmp%water%state%WaterStorageWetland ,& ! inout, wetland water storage [mm] + RunoffSurface => noahmp%water%flux%RunoffSurface ,& ! inout, surface runoff [mm] per soil timestep + EvapGroundNet => noahmp%water%flux%EvapGroundNet ,& ! inout, accumulated net ground evaporation per soil timestep [mm] + HeatLatentGrd => noahmp%energy%flux%HeatLatentGrd & ! inout, ground evaporation heat flux [W/m2] (+ to atm) + ) +! ---------------------------------------------------------------------- + + ! set initial value + EvapWaterFlux = 0.0 + EvapLatentHeat = 0.0 + + ! set psychrometric constant and VapPresSlope + if ( TemperatureSfc > ConstFreezePoint ) then + LatHeatSpec = ConstLatHeatEvap + else + LatHeatSpec = ConstLatHeatSublim + endif + + ! determine psychrometic constant + if ( TemperatureSfc < (273.15+26.85) ) then + PsychroConst = 0.00040 + else + PsychroConst = 0.00041 + endif + + ! calculate slope VapPresSlope based on surface temperature + if ( TemperatureSfc < (273.15+6.85) ) then + VapPresSlope = 0.00022 + elseif ( TemperatureSfc < (273.15+16.85) ) then + VapPresSlope = 0.00042 + elseif ( TemperatureSfc < (273.15+26.85) ) then + VapPresSlope = 0.00078 + else + VapPresSlope = 0.00132 + endif + + ! compute wetland water balance + LatHeatPot = 1.26 * VapPresSlope * (RadSwAbsSfc - RadLwNetSfc - HeatGroundTot) / & + (VapPresSlope + PsychroConst) ! LatHeatPot Potential latent heat W/m2 P-T method + WaterStorageWetland = WaterStorageWetland + RunoffSurface + if ( WaterStorageWetland > (LatHeatPot*TimeStep/LatHeatSpec*SoilSaturateFrac) ) then ! if current wetland storage is larger than PET rate + EvapWaterFlux = max(LatHeatPot / LatHeatSpec * SoilSaturateFrac, 0.0) ! evaporation rate mm/s + else ! if current wetland storage is less than PET rate + EvapWaterFlux = max(WaterStorageWetland/TimeStep, 0.0) ! use all of it + endif + WaterStorageWetland = max(WaterStorageWetland-EvapWaterFlux*TimeStep,0.0) ! adjust surface wetland storage + + ! adjuest energy and water balance + EvapLatentHeat = EvapWaterFlux * LatHeatSpec ! convert evaporation to latent heat flux + HeatSensibleSfc = HeatSensibleSfc - EvapLatentHeat ! reduce sensible heat flux + HeatLatentGrd = HeatLatentGrd + EvapLatentHeat ! increase direct evaporation + RunoffSurface = max(WaterStorageWetland - WetlandCapMax*1000.0, 0.0) ! excessive storage becomes runoff + EvapGroundNet = EvapGroundNet + EvapWaterFlux ! increase direct evaporation + WaterStorageWetland = min(WetlandCapMax*1000.0, WaterStorageWetland) + + end associate + + end subroutine WetlandWaterZhang22 + +end module WetlandWaterZhang22Mod