From 44183112cc55b30e5940f78b2d56a59654c528c5 Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 10:07:02 -0600 Subject: [PATCH 01/11] Merge effectively identical code for calculation of qflx_in_soil Previously, hydrologically-active urban columns were doing their own calculation of what is essentially qflx_in_soil. But it turns out that they were doing essentially the same thing as other landunits (once you note that frac_h2osfc is 0 for urban). So I'm removing this duplication. --- src/biogeophys/HydrologyNoDrainageMod.F90 | 2 +- src/biogeophys/SoilHydrologyMod.F90 | 59 +++++++---------------- 2 files changed, 19 insertions(+), 42 deletions(-) diff --git a/src/biogeophys/HydrologyNoDrainageMod.F90 b/src/biogeophys/HydrologyNoDrainageMod.F90 index 5cfca79a6a..c5bf0d6320 100644 --- a/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -196,7 +196,7 @@ subroutine HydrologyNoDrainage(bounds, & call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc,& infiltration_excess_runoff_inst, & - energyflux_inst, soilhydrology_inst, saturated_excess_runoff_inst, & + energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) call TotalSurfaceRunoff(bounds, num_hydrologyc, filter_hydrologyc, & diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index cf300064eb..832ffe8bca 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -223,37 +223,29 @@ subroutine SetQflxInputs(bounds, num_hydrologyc, filter_hydrologyc, & do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) + qflx_top_soil(c) = qflx_rain_plus_snomelt(c) + qflx_snow_h2osfc(c) + qflx_floodc(c) + ! ------------------------------------------------------------------------ ! Partition surface inputs between soil and h2osfc - ! - ! This is only done for soil & crop landunits, because other - ! hydrologically-active landunits (in particular, urban pervious road) do not use - ! h2osfc. - if (lun%itype(col%landunit(c)) == istsoil .or. lun%itype(col%landunit(c))==istcrop) then + ! ------------------------------------------------------------------------ - ! explicitly use frac_sno=0 if snl=0 - if (snl(c) >= 0) then - fsno=0._r8 - ! if no snow layers, sublimation is removed from h2osoi_ice in drainage - qflx_evap(c)=qflx_evap_grnd(c) - else - fsno=frac_sno(c) - qflx_evap(c)=qflx_ev_soil(c) - endif + if (snl(c) >= 0) then + fsno=0._r8 + ! if no snow layers, sublimation is removed from h2osoi_ice in drainage + qflx_evap(c)=qflx_evap_grnd(c) + else + fsno=frac_sno(c) + qflx_evap(c)=qflx_ev_soil(c) + endif - qflx_in_soil(c) = (1._r8 - frac_h2osfc(c)) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) - qflx_top_soil_to_h2osfc(c) = frac_h2osfc(c) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) + qflx_in_soil(c) = (1._r8 - frac_h2osfc(c)) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) + qflx_top_soil_to_h2osfc(c) = frac_h2osfc(c) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) - ! remove evaporation (snow treated in SnowHydrology) - qflx_in_soil(c) = qflx_in_soil(c) - (1.0_r8 - fsno - frac_h2osfc(c))*qflx_evap(c) - qflx_top_soil_to_h2osfc(c) = qflx_top_soil_to_h2osfc(c) - frac_h2osfc(c) * qflx_ev_h2osfc(c) + ! remove evaporation (snow treated in SnowHydrology) + qflx_in_soil(c) = qflx_in_soil(c) - (1.0_r8 - fsno - frac_h2osfc(c))*qflx_evap(c) + qflx_top_soil_to_h2osfc(c) = qflx_top_soil_to_h2osfc(c) - frac_h2osfc(c) * qflx_ev_h2osfc(c) - else - ! FIXME(wjs, 2017-08-15) Rather than this kludge, instead set qflx_in_soil and - ! qflx_top_soil_to_h2osfc to reasonable values for other landunits - qflx_in_soil(c) = 0._r8 - end if end do end associate @@ -264,7 +256,7 @@ end subroutine SetQflxInputs !----------------------------------------------------------------------- subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc, & infiltration_excess_runoff_inst, & - energyflux_inst, soilhydrology_inst, saturated_excess_runoff_inst, & + energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) ! ! !DESCRIPTION: @@ -284,7 +276,6 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f type(infiltration_excess_runoff_type), intent(in) :: infiltration_excess_runoff_inst type(energyflux_type) , intent(in) :: energyflux_inst type(soilhydrology_type) , intent(in) :: soilhydrology_inst - type(saturated_excess_runoff_type), intent(in) :: saturated_excess_runoff_inst type(waterstate_type) , intent(inout) :: waterstate_inst type(waterflux_type) , intent(inout) :: waterflux_inst ! @@ -298,27 +289,19 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f !----------------------------------------------------------------------- associate( & - snl => col%snl , & ! Input: [integer (:) ] minus number of snow layers - qflx_infl_excess => infiltration_excess_runoff_inst%qflx_infl_excess_col , & ! Input: [real(r8) (:)] infiltration excess runoff (mm H2O /s) qinmax => infiltration_excess_runoff_inst%qinmax_col , & ! Input: [real(r8) (:)] maximum infiltration capacity (mm H2O /s) frac_h2osfc => waterstate_inst%frac_h2osfc_col , & ! Input: [real(r8) (:) ] fraction of ground covered by surface water (0 to 1) frac_h2osfc_nosnow => waterstate_inst%frac_h2osfc_nosnow_col, & ! Input: [real(r8) (:) ] col fractional area with surface water greater than zero (if no snow present) - frac_sno => waterstate_inst%frac_sno_eff_col , & ! Input: [real(r8) (:) ] fraction of ground covered by snow (0 to 1) h2osfc => waterstate_inst%h2osfc_col , & ! Output: [real(r8) (:) ] surface water (mm) qflx_in_soil => waterflux_inst%qflx_in_soil_col , & ! Input: [real(r8) (:) ] surface input to soil (mm/s) qflx_top_soil_to_h2osfc => waterflux_inst%qflx_top_soil_to_h2osfc_col, & ! Input: [real(r8) (:)] portion of qflx_top_soil going to h2osfc, minus evaporation (mm/s) - qflx_ev_soil => waterflux_inst%qflx_ev_soil_col , & ! Input: [real(r8) (:) ] evaporation flux from soil (W/m**2) [+ to atm] - qflx_evap_grnd => waterflux_inst%qflx_evap_grnd_col , & ! Input: [real(r8) (:) ] ground surface evaporation rate (mm H2O/s) [+] - qflx_top_soil => waterflux_inst%qflx_top_soil_col , & ! Input: [real(r8) (:) ] net water input into soil from top (mm/s) qflx_infl_excess_surf => waterflux_inst%qflx_infl_excess_surf_col, & ! Output: [real(r8) (:) ] surface runoff due to infiltration excess (mm H2O /s) qflx_h2osfc_surf => waterflux_inst%qflx_h2osfc_surf_col , & ! Output: [real(r8) (:) ] surface water runoff (mm H2O /s) qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) - qflx_sat_excess_surf => saturated_excess_runoff_inst%qflx_sat_excess_surf_col, & ! Input: [real(r8) (:) ] surface runoff due to saturated surface (mm H2O /s) - h2osfc_thresh => soilhydrology_inst%h2osfc_thresh_col, & ! Input: [real(r8) (:) ] level at which h2osfc "percolates" h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: logical ) @@ -389,13 +372,7 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) else ! non-vegetated landunits (i.e. urban) use original CLM4 code - if (snl(c) >= 0) then - ! when no snow present, sublimation is removed in Drainage - qflx_infl(c) = qflx_top_soil(c) - qflx_sat_excess_surf(c) - qflx_evap_grnd(c) - else - qflx_infl(c) = qflx_top_soil(c) - qflx_sat_excess_surf(c) & - - (1.0_r8 - frac_sno(c)) * qflx_ev_soil(c) - end if + qflx_infl(c) = qflx_in_soil(c) qflx_h2osfc_surf(c) = 0._r8 qflx_infl_excess_surf(c) = 0._r8 endif From d69be4ec5d0d1a5b16b841fabf4db3bb21e17595 Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 11:42:13 -0600 Subject: [PATCH 02/11] Turn an array into a scalar Also add some fixme notes --- src/biogeophys/SoilHydrologyMod.F90 | 15 ++++++++++----- 1 file changed, 10 insertions(+), 5 deletions(-) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index 832ffe8bca..c316e61fad 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -196,8 +196,8 @@ subroutine SetQflxInputs(bounds, num_hydrologyc, filter_hydrologyc, & ! ! !LOCAL VARIABLES: integer :: fc, c - real(r8) :: qflx_evap(bounds%begc:bounds%endc) ! local evaporation array - real(r8) :: fsno ! copy of frac_sno + real(r8) :: qflx_evap ! evaporation for this column + real(r8) :: fsno ! copy of frac_sno character(len=*), parameter :: subname = 'SetQflxInputs' !----------------------------------------------------------------------- @@ -233,17 +233,17 @@ subroutine SetQflxInputs(bounds, num_hydrologyc, filter_hydrologyc, & if (snl(c) >= 0) then fsno=0._r8 ! if no snow layers, sublimation is removed from h2osoi_ice in drainage - qflx_evap(c)=qflx_evap_grnd(c) + qflx_evap=qflx_evap_grnd(c) else fsno=frac_sno(c) - qflx_evap(c)=qflx_ev_soil(c) + qflx_evap=qflx_ev_soil(c) endif qflx_in_soil(c) = (1._r8 - frac_h2osfc(c)) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) qflx_top_soil_to_h2osfc(c) = frac_h2osfc(c) * (qflx_top_soil(c) - qflx_sat_excess_surf(c)) ! remove evaporation (snow treated in SnowHydrology) - qflx_in_soil(c) = qflx_in_soil(c) - (1.0_r8 - fsno - frac_h2osfc(c))*qflx_evap(c) + qflx_in_soil(c) = qflx_in_soil(c) - (1.0_r8 - fsno - frac_h2osfc(c))*qflx_evap qflx_top_soil_to_h2osfc(c) = qflx_top_soil_to_h2osfc(c) - frac_h2osfc(c) * qflx_ev_h2osfc(c) end do @@ -336,6 +336,10 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f endif ! limit runoff to value of storage above S(pc) + + ! FIXME(wjs, 2017-08-22) Change >= to >: in the equals case, + ! qflx_h2osfc_surf will end up being 0 anyway. (This will avoid doing these + ! calculations for points with h2osfc = h2osfc_thresh = 0.) if(h2osfc(c) >= h2osfc_thresh(c) .and. h2osfcflag/=0) then ! spatially variable k_wet k_wet=1.0e-4_r8 * sin((rpi/180.) * col%topo_slope(c)) @@ -380,6 +384,7 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f enddo ! No infiltration for impervious urban surfaces + ! FIXME(wjs, 2017-08-22) Move the following to initialization? do fc = 1, num_urbanc c = filter_urbanc(fc) From feb92ef908f097258ba9f0fcc92729f7cc8430aa Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 14:47:10 -0600 Subject: [PATCH 03/11] Add routine, RouteInfiltrationExcess --- src/biogeophys/HydrologyNoDrainageMod.F90 | 5 +- src/biogeophys/SoilHydrologyMod.F90 | 83 +++++++++++++++++------ src/biogeophys/WaterfluxType.F90 | 2 + 3 files changed, 68 insertions(+), 22 deletions(-) diff --git a/src/biogeophys/HydrologyNoDrainageMod.F90 b/src/biogeophys/HydrologyNoDrainageMod.F90 index c5bf0d6320..9cb0e1d6ef 100644 --- a/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -62,7 +62,7 @@ subroutine HydrologyNoDrainage(bounds, & use SnowHydrologyMod , only : SnowCompaction, CombineSnowLayers, DivideSnowLayers, SnowCapping use SnowHydrologyMod , only : SnowWater, BuildSnowFilter use SoilHydrologyMod , only : CLMVICMap, SetSoilWaterFractions - use SoilHydrologyMod , only : SetQflxInputs, Infiltration, TotalSurfaceRunoff + use SoilHydrologyMod , only : SetQflxInputs, RouteInfiltrationExcess, Infiltration, TotalSurfaceRunoff use SoilHydrologyMod , only : UpdateUrbanPonding use SoilHydrologyMod , only : WaterTable, PerchedWaterTable use SoilHydrologyMod , only : ThetaBasedWaterTable, RenewCondensation @@ -194,6 +194,9 @@ subroutine HydrologyNoDrainage(bounds, & soilhydrology_inst, soilstate_inst, saturated_excess_runoff_inst, waterflux_inst, & waterstate_inst) + call RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & + waterflux_inst, infiltration_excess_runoff_inst, soilhydrology_inst) + call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc,& infiltration_excess_runoff_inst, & energyflux_inst, soilhydrology_inst, & diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index c316e61fad..e71c4a5979 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -252,6 +252,64 @@ subroutine SetQflxInputs(bounds, num_hydrologyc, filter_hydrologyc, & end subroutine SetQflxInputs + !----------------------------------------------------------------------- + subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & + waterflux_inst, infiltration_excess_runoff_inst, soilhydrology_inst) + ! + ! !DESCRIPTION: + ! Route the infiltration excess runoff flux + ! + ! This adjusts infiltration, input to h2osfc, and the associated surface runoff flux. + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_hydrologyc ! number of column soil points in column filter + integer , intent(in) :: filter_hydrologyc(:) ! column filter for soil points + type(waterflux_type) , intent(inout) :: waterflux_inst + type(infiltration_excess_runoff_type), intent(in) :: infiltration_excess_runoff_inst + type(soilhydrology_type) , intent(in) :: soilhydrology_inst + ! + ! !LOCAL VARIABLES: + integer :: fc, c + + character(len=*), parameter :: subname = 'RouteInfiltrationExcess' + !----------------------------------------------------------------------- + + associate( & + qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) + qflx_in_h2osfc => waterflux_inst%qflx_in_h2osfc_col , & ! Output: [real(r8) (:) ] total surface input to h2osfc + qflx_infl_excess_surf => waterflux_inst%qflx_infl_excess_surf_col , & ! Output: [real(r8) (:) ] surface runoff due to infiltration excess (mm H2O /s) + qflx_in_soil => waterflux_inst%qflx_in_soil_col , & ! Input: [real(r8) (:) ] surface input to soil (mm/s) + qflx_top_soil_to_h2osfc => waterflux_inst%qflx_top_soil_to_h2osfc_col , & ! Input: [real(r8) (:) ] portion of qflx_top_soil going to h2osfc, minus evaporation (mm/s) + + qflx_infl_excess => infiltration_excess_runoff_inst%qflx_infl_excess_col , & ! Input: [real(r8) (:) ] infiltration excess runoff (mm H2O /s) + + h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: logical + ) + + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) + if (lun%itype(col%landunit(c)) == istsoil .or. lun%itype(col%landunit(c))==istcrop) then + qflx_infl(c) = qflx_in_soil(c) - qflx_infl_excess(c) + if (h2osfcflag /= 0) then + qflx_in_h2osfc(c) = qflx_top_soil_to_h2osfc(c) + qflx_infl_excess(c) + qflx_infl_excess_surf(c) = 0._r8 + else + ! No h2osfc pool, so qflx_infl_excess goes directly to surface runoff + qflx_in_h2osfc(c) = qflx_top_soil_to_h2osfc(c) + qflx_infl_excess_surf(c) = qflx_infl_excess(c) + end if + else + ! non-vegetated landunits (i.e. urban) use original CLM4 code + qflx_infl(c) = qflx_in_soil(c) + qflx_in_h2osfc(c) = 0._r8 + qflx_infl_excess_surf(c) = 0._r8 + end if + end do + + end associate + + end subroutine RouteInfiltrationExcess !----------------------------------------------------------------------- subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc, & @@ -283,22 +341,18 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f integer :: c,l,fc ! indices real(r8) :: dtime ! land model time step (sec) real(r8) :: qflx_h2osfc_drain(bounds%begc:bounds%endc) ! bottom drainage from h2osfc - real(r8) :: qflx_in_h2osfc(bounds%begc:bounds%endc) ! net surface input to h2osfc real(r8) :: frac_infclust ! fraction of submerged area that is connected real(r8) :: k_wet ! linear reservoir coefficient for h2osfc !----------------------------------------------------------------------- associate( & - qflx_infl_excess => infiltration_excess_runoff_inst%qflx_infl_excess_col , & ! Input: [real(r8) (:)] infiltration excess runoff (mm H2O /s) qinmax => infiltration_excess_runoff_inst%qinmax_col , & ! Input: [real(r8) (:)] maximum infiltration capacity (mm H2O /s) frac_h2osfc => waterstate_inst%frac_h2osfc_col , & ! Input: [real(r8) (:) ] fraction of ground covered by surface water (0 to 1) frac_h2osfc_nosnow => waterstate_inst%frac_h2osfc_nosnow_col, & ! Input: [real(r8) (:) ] col fractional area with surface water greater than zero (if no snow present) h2osfc => waterstate_inst%h2osfc_col , & ! Output: [real(r8) (:) ] surface water (mm) - qflx_in_soil => waterflux_inst%qflx_in_soil_col , & ! Input: [real(r8) (:) ] surface input to soil (mm/s) - qflx_top_soil_to_h2osfc => waterflux_inst%qflx_top_soil_to_h2osfc_col, & ! Input: [real(r8) (:)] portion of qflx_top_soil going to h2osfc, minus evaporation (mm/s) - qflx_infl_excess_surf => waterflux_inst%qflx_infl_excess_surf_col, & ! Output: [real(r8) (:) ] surface runoff due to infiltration excess (mm H2O /s) + qflx_in_h2osfc => waterflux_inst%qflx_in_h2osfc_col , & ! Input: [real(r8) (:) ] total surface input to h2osfc qflx_h2osfc_surf => waterflux_inst%qflx_h2osfc_surf_col , & ! Output: [real(r8) (:) ] surface water runoff (mm H2O /s) qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) @@ -312,19 +366,9 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) + ! FIXME(wjs, 2017-08-22) Remove this conditional if (lun%itype(col%landunit(c)) == istsoil .or. lun%itype(col%landunit(c))==istcrop) then - !4. soil infiltration and h2osfc "run-on" - qflx_infl(c) = qflx_in_soil(c) - qflx_infl_excess(c) - if (h2osfcflag /= 0) then - qflx_in_h2osfc(c) = qflx_top_soil_to_h2osfc(c) + qflx_infl_excess(c) - qflx_infl_excess_surf(c) = 0._r8 - else - ! No h2osfc pool, so qflx_infl_excess goes directly to surface runoff - qflx_in_h2osfc(c) = qflx_top_soil_to_h2osfc(c) - qflx_infl_excess_surf(c) = qflx_infl_excess(c) - end if - !5. surface runoff from h2osfc if (h2osfcflag==1) then ! calculate runoff from h2osfc ------------------------------------- @@ -353,10 +397,8 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f ! cutoff lower limit if ( qflx_h2osfc_surf(c) < 1.0e-8) qflx_h2osfc_surf(c) = 0._r8 - qflx_in_h2osfc(c) = qflx_in_h2osfc(c) - qflx_h2osfc_surf(c) - !6. update h2osfc prior to calculating bottom drainage from h2osfc - h2osfc(c) = h2osfc(c) + qflx_in_h2osfc(c) * dtime + h2osfc(c) = h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime !-- if all water evaporates, there will be no bottom drainage if (h2osfc(c) < 0.0) then @@ -376,9 +418,8 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) else ! non-vegetated landunits (i.e. urban) use original CLM4 code - qflx_infl(c) = qflx_in_soil(c) + ! FIXME(wjs, 2017-08-22) Remove this qflx_h2osfc_surf(c) = 0._r8 - qflx_infl_excess_surf(c) = 0._r8 endif enddo diff --git a/src/biogeophys/WaterfluxType.F90 b/src/biogeophys/WaterfluxType.F90 index 33c84bf908..2398763a35 100644 --- a/src/biogeophys/WaterfluxType.F90 +++ b/src/biogeophys/WaterfluxType.F90 @@ -82,6 +82,7 @@ module WaterfluxType real(r8), pointer :: qflx_top_soil_col (:) ! col net water input into soil from top (mm/s) real(r8), pointer :: qflx_in_soil_col (:) ! col surface input to soil (mm/s) real(r8), pointer :: qflx_top_soil_to_h2osfc_col(:) ! col portion of qflx_top_soil going to h2osfc, minus evaporation (mm/s) + real(r8), pointer :: qflx_in_h2osfc_col(:) ! col total surface input to h2osfc real(r8), pointer :: qflx_h2osfc_to_ice_col (:) ! col conversion of h2osfc to ice real(r8), pointer :: qflx_snow_h2osfc_col (:) ! col snow falling on surface water real(r8), pointer :: qflx_drain_perched_col (:) ! col sub-surface runoff from perched wt (mm H2O /s) @@ -226,6 +227,7 @@ subroutine InitAllocate(this, bounds) allocate(this%qflx_top_soil_col (begc:endc)) ; this%qflx_top_soil_col (:) = nan allocate(this%qflx_in_soil_col (begc:endc)) ; this%qflx_in_soil_col (:) = nan allocate(this%qflx_top_soil_to_h2osfc_col(begc:endc)) ; this%qflx_top_soil_to_h2osfc_col(:) = nan + allocate(this%qflx_in_h2osfc_col (begc:endc)) ; this%qflx_in_h2osfc_col(:) = nan allocate(this%qflx_h2osfc_to_ice_col (begc:endc)) ; this%qflx_h2osfc_to_ice_col (:) = nan allocate(this%qflx_infl_excess_surf_col(begc:endc)) ; this%qflx_infl_excess_surf_col(:) = nan allocate(this%qflx_h2osfc_surf_col (begc:endc)) ; this%qflx_h2osfc_surf_col (:) = nan From f5abae2f8a5c8438594bd573833ea9ae223ece3b Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 14:57:31 -0600 Subject: [PATCH 04/11] Move setting of qflx_infl for urban Main point is so that subroutine Infiltration doesn't need to bother with an urban filter --- src/biogeophys/HydrologyDrainageMod.F90 | 1 + src/biogeophys/HydrologyNoDrainageMod.F90 | 2 +- src/biogeophys/SoilHydrologyMod.F90 | 14 +------------- 3 files changed, 3 insertions(+), 14 deletions(-) diff --git a/src/biogeophys/HydrologyDrainageMod.F90 b/src/biogeophys/HydrologyDrainageMod.F90 index 226d85dd40..e6b1b8ee27 100644 --- a/src/biogeophys/HydrologyDrainageMod.F90 +++ b/src/biogeophys/HydrologyDrainageMod.F90 @@ -185,6 +185,7 @@ subroutine HydrologyDrainage(bounds, & qflx_drain_perched(c) = 0._r8 qflx_rsub_sat(c) = spval + qflx_infl(c) = 0._r8 end if diff --git a/src/biogeophys/HydrologyNoDrainageMod.F90 b/src/biogeophys/HydrologyNoDrainageMod.F90 index 9cb0e1d6ef..3dce903d26 100644 --- a/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -197,7 +197,7 @@ subroutine HydrologyNoDrainage(bounds, & call RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & waterflux_inst, infiltration_excess_runoff_inst, soilhydrology_inst) - call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc,& + call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & infiltration_excess_runoff_inst, & energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index e71c4a5979..16723f9238 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -312,7 +312,7 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & end subroutine RouteInfiltrationExcess !----------------------------------------------------------------------- - subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filter_urbanc, & + subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & infiltration_excess_runoff_inst, & energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) @@ -329,8 +329,6 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f type(bounds_type) , intent(in) :: bounds integer , intent(in) :: num_hydrologyc ! number of column soil points in column filter integer , intent(in) :: filter_hydrologyc(:) ! column filter for soil points - integer , intent(in) :: num_urbanc ! number of column urban points in column filter - integer , intent(in) :: filter_urbanc(:) ! column filter for urban points type(infiltration_excess_runoff_type), intent(in) :: infiltration_excess_runoff_inst type(energyflux_type) , intent(in) :: energyflux_inst type(soilhydrology_type) , intent(in) :: soilhydrology_inst @@ -424,16 +422,6 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, f enddo - ! No infiltration for impervious urban surfaces - ! FIXME(wjs, 2017-08-22) Move the following to initialization? - - do fc = 1, num_urbanc - c = filter_urbanc(fc) - if (col%itype(c) /= icol_road_perv) then - qflx_infl(c) = 0._r8 - end if - end do - end associate end subroutine Infiltration From 5a8be1a4ef776a4703f3231d982b341ad59907e8 Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 15:14:52 -0600 Subject: [PATCH 05/11] Treat urban like other landunits in h2osfc calculations This should give the same answers as before simply because the relevant quantities (e.g., fluxes in) are already 0. This simplifies the code, and is a step towards reconciling urban pervious road with other hydrologically-active columns. Also, change >= to > in conditional for computing qflx_h2osfc_surf: in the equals case, qflx_h2osfc_surf will end up being 0 anyway. (This will avoid doing these calculations for points with h2osfc = h2osfc_thresh = 0.) --- src/biogeophys/SoilHydrologyMod.F90 | 81 ++++++++++++----------------- 1 file changed, 34 insertions(+), 47 deletions(-) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index 16723f9238..de7c69527e 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -360,66 +360,53 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & dtime = get_step_size() - ! Infiltration into surface soil layer (minus the evaporation) - do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) - ! FIXME(wjs, 2017-08-22) Remove this conditional - if (lun%itype(col%landunit(c)) == istsoil .or. lun%itype(col%landunit(c))==istcrop) then - - !5. surface runoff from h2osfc - if (h2osfcflag==1) then - ! calculate runoff from h2osfc ------------------------------------- - if (frac_h2osfc_nosnow(c) <= pc) then - frac_infclust=0.0_r8 - else - frac_infclust=(frac_h2osfc_nosnow(c)-pc)**mu - endif - endif - ! limit runoff to value of storage above S(pc) - - ! FIXME(wjs, 2017-08-22) Change >= to >: in the equals case, - ! qflx_h2osfc_surf will end up being 0 anyway. (This will avoid doing these - ! calculations for points with h2osfc = h2osfc_thresh = 0.) - if(h2osfc(c) >= h2osfc_thresh(c) .and. h2osfcflag/=0) then - ! spatially variable k_wet - k_wet=1.0e-4_r8 * sin((rpi/180.) * col%topo_slope(c)) - qflx_h2osfc_surf(c) = k_wet * frac_infclust * (h2osfc(c) - h2osfc_thresh(c)) - - qflx_h2osfc_surf(c)=min(qflx_h2osfc_surf(c),(h2osfc(c) - h2osfc_thresh(c))/dtime) + !5. surface runoff from h2osfc + if (h2osfcflag==1) then + ! calculate runoff from h2osfc ------------------------------------- + if (frac_h2osfc_nosnow(c) <= pc) then + frac_infclust=0.0_r8 else - qflx_h2osfc_surf(c)= 0._r8 + frac_infclust=(frac_h2osfc_nosnow(c)-pc)**mu endif + endif - ! cutoff lower limit - if ( qflx_h2osfc_surf(c) < 1.0e-8) qflx_h2osfc_surf(c) = 0._r8 + ! limit runoff to value of storage above S(pc) + if(h2osfc(c) > h2osfc_thresh(c) .and. h2osfcflag/=0) then + ! spatially variable k_wet + k_wet=1.0e-4_r8 * sin((rpi/180.) * col%topo_slope(c)) + qflx_h2osfc_surf(c) = k_wet * frac_infclust * (h2osfc(c) - h2osfc_thresh(c)) - !6. update h2osfc prior to calculating bottom drainage from h2osfc - h2osfc(c) = h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime + qflx_h2osfc_surf(c)=min(qflx_h2osfc_surf(c),(h2osfc(c) - h2osfc_thresh(c))/dtime) + else + qflx_h2osfc_surf(c)= 0._r8 + endif - !-- if all water evaporates, there will be no bottom drainage - if (h2osfc(c) < 0.0) then - qflx_infl(c) = qflx_infl(c) + h2osfc(c)/dtime - h2osfc(c) = 0.0 - qflx_h2osfc_drain(c)= 0._r8 - else - qflx_h2osfc_drain(c)=min(frac_h2osfc(c)*qinmax(c),h2osfc(c)/dtime) - endif + ! cutoff lower limit + if ( qflx_h2osfc_surf(c) < 1.0e-8) qflx_h2osfc_surf(c) = 0._r8 - if(h2osfcflag==0) then - qflx_h2osfc_drain(c)= max(0._r8,h2osfc(c)/dtime) !ensure no h2osfc - endif + !6. update h2osfc prior to calculating bottom drainage from h2osfc + h2osfc(c) = h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime - !7. remove drainage from h2osfc and add to qflx_infl - h2osfc(c) = h2osfc(c) - qflx_h2osfc_drain(c) * dtime - qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) + !-- if all water evaporates, there will be no bottom drainage + if (h2osfc(c) < 0.0) then + qflx_infl(c) = qflx_infl(c) + h2osfc(c)/dtime + h2osfc(c) = 0.0 + qflx_h2osfc_drain(c)= 0._r8 else - ! non-vegetated landunits (i.e. urban) use original CLM4 code - ! FIXME(wjs, 2017-08-22) Remove this - qflx_h2osfc_surf(c) = 0._r8 + qflx_h2osfc_drain(c)=min(frac_h2osfc(c)*qinmax(c),h2osfc(c)/dtime) endif + if(h2osfcflag==0) then + qflx_h2osfc_drain(c)= max(0._r8,h2osfc(c)/dtime) !ensure no h2osfc + endif + + !7. remove drainage from h2osfc and add to qflx_infl + h2osfc(c) = h2osfc(c) - qflx_h2osfc_drain(c) * dtime + qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) + enddo end associate From 955f4917dc0b109433cb62a378615dc2d33265ea Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 15:19:36 -0600 Subject: [PATCH 06/11] Rename Infiltration to UpdateH2osfc since that's what it does now --- src/biogeophys/HydrologyNoDrainageMod.F90 | 4 ++-- src/biogeophys/SoilHydrologyMod.F90 | 8 ++++---- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/biogeophys/HydrologyNoDrainageMod.F90 b/src/biogeophys/HydrologyNoDrainageMod.F90 index 3dce903d26..22def389ec 100644 --- a/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -62,7 +62,7 @@ subroutine HydrologyNoDrainage(bounds, & use SnowHydrologyMod , only : SnowCompaction, CombineSnowLayers, DivideSnowLayers, SnowCapping use SnowHydrologyMod , only : SnowWater, BuildSnowFilter use SoilHydrologyMod , only : CLMVICMap, SetSoilWaterFractions - use SoilHydrologyMod , only : SetQflxInputs, RouteInfiltrationExcess, Infiltration, TotalSurfaceRunoff + use SoilHydrologyMod , only : SetQflxInputs, RouteInfiltrationExcess, UpdateH2osfc, TotalSurfaceRunoff use SoilHydrologyMod , only : UpdateUrbanPonding use SoilHydrologyMod , only : WaterTable, PerchedWaterTable use SoilHydrologyMod , only : ThetaBasedWaterTable, RenewCondensation @@ -197,7 +197,7 @@ subroutine HydrologyNoDrainage(bounds, & call RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & waterflux_inst, infiltration_excess_runoff_inst, soilhydrology_inst) - call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & + call UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & infiltration_excess_runoff_inst, & energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index de7c69527e..15486ab25b 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -36,7 +36,7 @@ module SoilHydrologyMod public :: SoilHydReadNML ! Read in the Soil hydrology namelist public :: SetSoilWaterFractions ! Set diagnostic variables related to the fraction of water and ice in each layer public :: SetQflxInputs ! Set the flux of water into the soil from the top - public :: Infiltration ! Calculate infiltration into surface soil layer + public :: UpdateH2osfc ! Calculate fluxes out of h2osfc and update the h2osfc state public :: TotalSurfaceRunoff ! Calculate total surface runoff public :: UpdateUrbanPonding ! Update the state variable representing ponding on urban surfaces public :: WaterTable ! Calculate water table before imposing drainage @@ -312,13 +312,13 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & end subroutine RouteInfiltrationExcess !----------------------------------------------------------------------- - subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & + subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & infiltration_excess_runoff_inst, & energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) ! ! !DESCRIPTION: - ! Calculate infiltration into surface soil layer (minus the evaporation) + ! Calculate fluxes out of h2osfc and update the h2osfc state ! ! !USES: use shr_const_mod , only : shr_const_pi @@ -411,7 +411,7 @@ subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & end associate - end subroutine Infiltration + end subroutine UpdateH2osfc !----------------------------------------------------------------------- subroutine TotalSurfaceRunoff(bounds, num_hydrologyc, filter_hydrologyc, & From b586cebd88939540863319eec542e9049e62511d Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Tue, 22 Aug 2017 20:40:26 -0600 Subject: [PATCH 07/11] Separate flux calculations from state updates for h2osfc --- src/biogeophys/SoilHydrologyMod.F90 | 240 ++++++++++++++++++++++------ 1 file changed, 191 insertions(+), 49 deletions(-) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index 15486ab25b..44d35bf3f1 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -4,6 +4,7 @@ module SoilHydrologyMod ! !DESCRIPTION: ! Calculate soil hydrology ! +#include "shr_assert.h" use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use abortutils , only : endrun @@ -48,6 +49,10 @@ module SoilHydrologyMod public :: LateralFlowPowerLaw ! Calculate lateral flow based on power law drainage function public :: RenewCondensation ! Misc. corrections + ! !PRIVATE MEMBER FUNCTIONS: + private :: QflxH2osfcSurf ! Compute qflx_h2osfc_surf + private :: QflxH2osfcDrain ! Compute qflx_h2osfc_drain + !----------------------------------------------------------------------- real(r8), private :: baseflow_scalar = 1.e-2_r8 @@ -284,7 +289,7 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & qflx_infl_excess => infiltration_excess_runoff_inst%qflx_infl_excess_col , & ! Input: [real(r8) (:) ] infiltration excess runoff (mm H2O /s) - h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: logical + h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: integer ) do fc = 1, num_hydrologyc @@ -322,7 +327,7 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & ! ! !USES: use shr_const_mod , only : shr_const_pi - use clm_varcon , only : denh2o, denice, roverg, wimp, pc, mu, tfrz + use clm_varcon , only : denh2o, denice, roverg, wimp, tfrz use column_varcon , only : icol_roof, icol_road_imperv, icol_sunwall, icol_shadewall, icol_road_perv ! ! !ARGUMENTS: @@ -338,9 +343,11 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & ! !LOCAL VARIABLES: integer :: c,l,fc ! indices real(r8) :: dtime ! land model time step (sec) + real(r8) :: h2osfc_partial(bounds%begc:bounds%endc) ! partially-updated h2osfc + ! FIXME(wjs, 2017-08-22) put qflx_h2osfc_drain in waterflux_type (this will be + ! needed when we move the update of qflx_infl to a separate routine) real(r8) :: qflx_h2osfc_drain(bounds%begc:bounds%endc) ! bottom drainage from h2osfc - real(r8) :: frac_infclust ! fraction of submerged area that is connected - real(r8) :: k_wet ! linear reservoir coefficient for h2osfc + logical :: truncate_h2osfc_to_zero(bounds%begc:bounds%endc) !----------------------------------------------------------------------- associate( & @@ -355,63 +362,198 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) h2osfc_thresh => soilhydrology_inst%h2osfc_thresh_col, & ! Input: [real(r8) (:) ] level at which h2osfc "percolates" - h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: logical + h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: integer ) - dtime = get_step_size() + dtime = get_step_size() - do fc = 1, num_hydrologyc - c = filter_hydrologyc(fc) + call QflxH2osfcSurf(bounds, num_hydrologyc, filter_hydrologyc, & + h2osfcflag = h2osfcflag, & + h2osfc = h2osfc(bounds%begc:bounds%endc), & + h2osfc_thresh = h2osfc_thresh(bounds%begc:bounds%endc), & + frac_h2osfc_nosnow = frac_h2osfc_nosnow(bounds%begc:bounds%endc), & + topo_slope = col%topo_slope(bounds%begc:bounds%endc), & + qflx_h2osfc_surf = qflx_h2osfc_surf(bounds%begc:bounds%endc)) + + ! Update h2osfc prior to calculating bottom drainage from h2osfc. + ! + ! This could be removed if we wanted to do a straight forward Euler, and/or set + ! things up for a more sophisticated solution method. In the latter, rather than + ! having h2osfc_partial, we'd have some more sophisticated state estimate here + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) + h2osfc_partial(c) = h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime + end do - !5. surface runoff from h2osfc - if (h2osfcflag==1) then - ! calculate runoff from h2osfc ------------------------------------- - if (frac_h2osfc_nosnow(c) <= pc) then - frac_infclust=0.0_r8 - else - frac_infclust=(frac_h2osfc_nosnow(c)-pc)**mu - endif - endif + call QflxH2osfcDrain(bounds, num_hydrologyc, filter_hydrologyc, & + h2osfcflag = h2osfcflag, & + h2osfc = h2osfc_partial(bounds%begc:bounds%endc), & + frac_h2osfc = frac_h2osfc(bounds%begc:bounds%endc), & + qinmax = qinmax(bounds%begc:bounds%endc), & + qflx_h2osfc_drain = qflx_h2osfc_drain(bounds%begc:bounds%endc), & + truncate_h2osfc_to_zero = truncate_h2osfc_to_zero(bounds%begc:bounds%endc)) - ! limit runoff to value of storage above S(pc) - if(h2osfc(c) > h2osfc_thresh(c) .and. h2osfcflag/=0) then - ! spatially variable k_wet - k_wet=1.0e-4_r8 * sin((rpi/180.) * col%topo_slope(c)) - qflx_h2osfc_surf(c) = k_wet * frac_infclust * (h2osfc(c) - h2osfc_thresh(c)) + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) - qflx_h2osfc_surf(c)=min(qflx_h2osfc_surf(c),(h2osfc(c) - h2osfc_thresh(c))/dtime) - else - qflx_h2osfc_surf(c)= 0._r8 - endif + ! The parenthesization of this expression is just needed to maintain bfb answers + ! in the major refactor. Note that the first parenthesized expression is + ! h2osfc_partial(c), but I'm writing it out explicitly to facilitate a possible + ! future removal of h2osfc_partial. + h2osfc(c) = (h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime) & + - qflx_h2osfc_drain(c) * dtime + + ! TODO(wjs, 2017-08-22) This is here to maintain bit-for-bit answers with the old + ! code. If we're okay changing answers, we could remove it. Or maybe we want to + ! put a more general truncation in here, like: + ! if (abs(h2osfc(c)) < 1.e-14_r8 * abs(h2osfc_orig)) h2osfc(c) = 0._r8 + ! But I'm not sure what the general philosophy is regarding whether and when we + ! want to do truncations like this. + if (truncate_h2osfc_to_zero(c)) then + h2osfc(c) = 0._r8 + end if - ! cutoff lower limit - if ( qflx_h2osfc_surf(c) < 1.0e-8) qflx_h2osfc_surf(c) = 0._r8 + end do - !6. update h2osfc prior to calculating bottom drainage from h2osfc - h2osfc(c) = h2osfc(c) + (qflx_in_h2osfc(c) - qflx_h2osfc_surf(c)) * dtime + ! FIXME(wjs, 2017-08-22) Move this update to a separate routine that calculations + ! qflx_infl from its components. + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) + qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) + end do - !-- if all water evaporates, there will be no bottom drainage - if (h2osfc(c) < 0.0) then - qflx_infl(c) = qflx_infl(c) + h2osfc(c)/dtime - h2osfc(c) = 0.0 - qflx_h2osfc_drain(c)= 0._r8 - else - qflx_h2osfc_drain(c)=min(frac_h2osfc(c)*qinmax(c),h2osfc(c)/dtime) - endif + end associate - if(h2osfcflag==0) then - qflx_h2osfc_drain(c)= max(0._r8,h2osfc(c)/dtime) !ensure no h2osfc - endif + end subroutine UpdateH2osfc - !7. remove drainage from h2osfc and add to qflx_infl - h2osfc(c) = h2osfc(c) - qflx_h2osfc_drain(c) * dtime - qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) + !----------------------------------------------------------------------- + subroutine QflxH2osfcSurf(bounds, num_hydrologyc, filter_hydrologyc, & + h2osfcflag, h2osfc, h2osfc_thresh, frac_h2osfc_nosnow, topo_slope, & + qflx_h2osfc_surf) + ! + ! !DESCRIPTION: + ! Compute qflx_h2osfc_surf + ! + ! !USES: + use clm_varcon, only : pc, mu + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_hydrologyc ! number of column soil points in column filter + integer , intent(in) :: filter_hydrologyc(:) ! column filter for soil points + integer , intent(in) :: h2osfcflag ! true => surface water is active + real(r8) , intent(in) :: h2osfc( bounds%begc: ) ! surface water (mm) + real(r8) , intent(in) :: h2osfc_thresh( bounds%begc: ) ! level at which h2osfc "percolates" + real(r8) , intent(in) :: frac_h2osfc_nosnow( bounds%begc: ) ! fractional area with surface water greater than zero (if no snow present) + real(r8) , intent(in) :: topo_slope( bounds%begc: ) ! topographic slope + real(r8) , intent(inout) :: qflx_h2osfc_surf( bounds%begc: ) ! surface water runoff (mm H2O /s) + ! + ! !LOCAL VARIABLES: + integer :: fc, c + real(r8) :: dtime ! land model time step (sec) + real(r8) :: frac_infclust ! fraction of submerged area that is connected + real(r8) :: k_wet ! linear reservoir coefficient for h2osfc - enddo + character(len=*), parameter :: subname = 'QflxH2osfcSurf' + !----------------------------------------------------------------------- - end associate + SHR_ASSERT_ALL((ubound(h2osfc) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(h2osfc_thresh) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(frac_h2osfc_nosnow) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(topo_slope) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(qflx_h2osfc_surf) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + + dtime = get_step_size() + + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) + + if (h2osfcflag==1) then + if (frac_h2osfc_nosnow(c) <= pc) then + frac_infclust=0.0_r8 + else + frac_infclust=(frac_h2osfc_nosnow(c)-pc)**mu + endif + endif + + ! limit runoff to value of storage above S(pc) + if(h2osfc(c) > h2osfc_thresh(c) .and. h2osfcflag/=0) then + ! spatially variable k_wet + k_wet=1.0e-4_r8 * sin((rpi/180.) * topo_slope(c)) + qflx_h2osfc_surf(c) = k_wet * frac_infclust * (h2osfc(c) - h2osfc_thresh(c)) + + qflx_h2osfc_surf(c)=min(qflx_h2osfc_surf(c),(h2osfc(c) - h2osfc_thresh(c))/dtime) + else + qflx_h2osfc_surf(c)= 0._r8 + endif + + ! cutoff lower limit + if ( qflx_h2osfc_surf(c) < 1.0e-8) then + qflx_h2osfc_surf(c) = 0._r8 + end if + + end do + + end subroutine QflxH2osfcSurf + + !----------------------------------------------------------------------- + subroutine QflxH2osfcDrain(bounds, num_hydrologyc, filter_hydrologyc, & + h2osfcflag, h2osfc, frac_h2osfc, qinmax, & + qflx_h2osfc_drain, truncate_h2osfc_to_zero) + ! + ! !DESCRIPTION: + ! Compute qflx_h2osfc_drain + ! + ! Note that, if h2osfc is negative, then qflx_h2osfc_drain will be negative - acting + ! to exactly restore h2osfc to 0. In this case, truncate_h2osfc_to_zero will be set + ! to true; in all other cases, truncate_h2osfc_to_zero will be set to false. This + ! particular behavior of truncate_h2osfc_to_zero is just like this to maintain + ! bit-for-bit answers with the old code (see also the TODO note in UpdateH2osfc). + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_hydrologyc ! number of column soil points in column filter + integer , intent(in) :: filter_hydrologyc(:) ! column filter for soil points + integer , intent(in) :: h2osfcflag ! true => surface water is active + real(r8) , intent(in) :: h2osfc( bounds%begc: ) ! surface water (mm) + real(r8) , intent(in) :: frac_h2osfc( bounds%begc: ) ! fraction of ground covered by surface water (0 to 1) + real(r8) , intent(in) :: qinmax( bounds%begc: ) ! maximum infiltration capacity (mm H2O /s) + real(r8) , intent(inout) :: qflx_h2osfc_drain( bounds%begc: ) ! bottom drainage from h2osfc (mm H2O /s) + logical , intent(inout) :: truncate_h2osfc_to_zero( bounds%begc: ) ! whether h2osfc should be truncated to 0 to correct for roundoff errors, in order to maintain bit-for-bit the same answers as the old code + ! + ! !LOCAL VARIABLES: + integer :: fc, c + real(r8) :: dtime ! land model time step (sec) + + character(len=*), parameter :: subname = 'QflxH2osfcDrain' + !----------------------------------------------------------------------- + + SHR_ASSERT_ALL((ubound(h2osfc) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(frac_h2osfc) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(qinmax) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(qflx_h2osfc_drain) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + SHR_ASSERT_ALL((ubound(truncate_h2osfc_to_zero) == (/bounds%endc/)), errMsg(sourcefile, __LINE__)) + + dtime = get_step_size() + + do fc = 1, num_hydrologyc + c = filter_hydrologyc(fc) + + if (h2osfc(c) < 0.0) then + qflx_h2osfc_drain(c) = h2osfc(c)/dtime + truncate_h2osfc_to_zero(c) = .true. + else + qflx_h2osfc_drain(c)=min(frac_h2osfc(c)*qinmax(c),h2osfc(c)/dtime) + if(h2osfcflag==0) then + ! ensure no h2osfc + qflx_h2osfc_drain(c)= max(0._r8,h2osfc(c)/dtime) + end if + truncate_h2osfc_to_zero(c) = .false. + end if + end do + + end subroutine QflxH2osfcDrain - end subroutine UpdateH2osfc !----------------------------------------------------------------------- subroutine TotalSurfaceRunoff(bounds, num_hydrologyc, filter_hydrologyc, & @@ -986,7 +1128,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ice => soilhydrology_inst%ice_col , & ! Input: [real(r8) (:,:) ] soil layer moisture (mm) qcharge => soilhydrology_inst%qcharge_col , & ! Input: [real(r8) (:) ] aquifer recharge rate (mm/s) origflag => soilhydrology_inst%origflag , & ! Input: logical - h2osfcflag => soilhydrology_inst%h2osfcflag , & ! Input: logical + h2osfcflag => soilhydrology_inst%h2osfcflag , & ! Input: integer qflx_snwcp_liq => waterflux_inst%qflx_snwcp_liq_col , & ! Output: [real(r8) (:) ] excess liquid h2o due to snow capping (outgoing) (mm H2O /s) [+] qflx_ice_runoff_xs => waterflux_inst%qflx_ice_runoff_xs_col , & ! Output: [real(r8) (:) ] solid runoff from excess ice in soil (mm H2O /s) [+] @@ -2032,7 +2174,7 @@ subroutine LateralFlowPowerLaw(bounds, num_hydrologyc, filter_hydrologyc, & ice => soilhydrology_inst%ice_col , & ! Input: [real(r8) (:,:) ] soil layer moisture (mm) qcharge => soilhydrology_inst%qcharge_col , & ! Input: [real(r8) (:) ] aquifer recharge rate (mm/s) origflag => soilhydrology_inst%origflag , & ! Input: logical - h2osfcflag => soilhydrology_inst%h2osfcflag , & ! Input: logical + h2osfcflag => soilhydrology_inst%h2osfcflag , & ! Input: integer qflx_snwcp_liq => waterflux_inst%qflx_snwcp_liq_col , & ! Output: [real(r8) (:) ] excess rainfall due to snow capping (mm H2O /s) [+] qflx_ice_runoff_xs => waterflux_inst%qflx_ice_runoff_xs_col , & ! Output: [real(r8) (:) ] solid runoff from excess ice in soil (mm H2O /s) [+] From 00352314ed9e23dad2bc1759a1656ff45883193d Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Wed, 23 Aug 2017 14:25:53 -0600 Subject: [PATCH 08/11] Change comment for qinmax Martyn Clark points out that this is maximum infiltration rate, not maximum infiltration capacity --- src/biogeophys/InfiltrationExcessRunoffMod.F90 | 10 +++++----- src/biogeophys/SaturatedExcessRunoffMod.F90 | 2 +- src/biogeophys/SoilHydrologyMod.F90 | 4 ++-- 3 files changed, 8 insertions(+), 8 deletions(-) diff --git a/src/biogeophys/InfiltrationExcessRunoffMod.F90 b/src/biogeophys/InfiltrationExcessRunoffMod.F90 index 1d482f1d84..690ec7c10b 100644 --- a/src/biogeophys/InfiltrationExcessRunoffMod.F90 +++ b/src/biogeophys/InfiltrationExcessRunoffMod.F90 @@ -36,7 +36,7 @@ module InfiltrationExcessRunoffMod ! Both of these give averages over the entire column. However, qinmax is implicitly ! 0 over the fraction of the column given by fsat, and qflx_infl_excess is ! implicitly 0 over both fsat and frac_h2osfc. - real(r8), pointer, public :: qinmax_col(:) ! maximum infiltration capacity (mm H2O /s) + real(r8), pointer, public :: qinmax_col(:) ! maximum infiltration rate (mm H2O /s) real(r8), pointer, public :: qflx_infl_excess_col(:) ! infiltration excess runoff (mm H2O /s) ! Private data members @@ -186,13 +186,13 @@ subroutine InfiltrationExcessRunoff(this, bounds, num_hydrologyc, filter_hydrolo ! ! !LOCAL VARIABLES: integer :: fc, c - real(r8) :: qinmax_on_unsaturated_area(bounds%begc:bounds%endc) ! maximum infiltration capacity on the unsaturated fraction of the column (mm H2O /s) + real(r8) :: qinmax_on_unsaturated_area(bounds%begc:bounds%endc) ! maximum infiltration rate on the unsaturated fraction of the column (mm H2O /s) character(len=*), parameter :: subname = 'InfiltrationExcessRunoff' !----------------------------------------------------------------------- associate( & - qinmax => this%qinmax_col , & ! Output: [real(r8) (:) ] maximum infiltration capacity (mm H2O /s) + qinmax => this%qinmax_col , & ! Output: [real(r8) (:) ] maximum infiltration rate (mm H2O /s) qflx_infl_excess => this%qflx_infl_excess_col , & ! Output: [real(r8) (:) ] infiltration excess runoff (mm H2O /s) fsat => saturated_excess_runoff_inst%fsat_col, & ! Input: [real(r8) (:) ] fractional area with water table at surface @@ -245,7 +245,7 @@ subroutine ComputeQinmaxHksat(bounds, num_hydrologyc, filter_hydrologyc, & integer, intent(in) :: filter_hydrologyc(:) ! column filter for soil points type(soilhydrology_type) , intent(in) :: soilhydrology_inst type(soilstate_type), intent(in) :: soilstate_inst - real(r8), intent(inout) :: qinmax_on_unsaturated_area( bounds%begc: ) ! maximum infiltration capacity on the unsaturated fraction of the column (mm H2O /s) + real(r8), intent(inout) :: qinmax_on_unsaturated_area( bounds%begc: ) ! maximum infiltration rate on the unsaturated fraction of the column (mm H2O /s) ! ! !LOCAL VARIABLES: integer :: fc, c @@ -288,7 +288,7 @@ subroutine ComputeQinmaxVic(bounds, num_hydrologyc, filter_hydrologyc, & type(soilhydrology_type) , intent(in) :: soilhydrology_inst real(r8) , intent(in) :: fsat( bounds%begc: ) ! fractional area with water table at surface real(r8) , intent(in) :: qflx_in_soil( bounds%begc: ) ! surface input to soil (mm/s) - real(r8) , intent(inout) :: qinmax_on_unsaturated_area( bounds%begc: ) ! maximum infiltration capacity on the unsaturated fraction of the column (mm H2O /s) + real(r8) , intent(inout) :: qinmax_on_unsaturated_area( bounds%begc: ) ! maximum infiltration rate on the unsaturated fraction of the column (mm H2O /s) ! ! !LOCAL VARIABLES: integer :: fc, c diff --git a/src/biogeophys/SaturatedExcessRunoffMod.F90 b/src/biogeophys/SaturatedExcessRunoffMod.F90 index bb24aea1d6..21dcedd5f6 100644 --- a/src/biogeophys/SaturatedExcessRunoffMod.F90 +++ b/src/biogeophys/SaturatedExcessRunoffMod.F90 @@ -231,7 +231,7 @@ subroutine SaturatedExcessRunoff (this, bounds, num_hydrologyc, filter_hydrology ! ------------------------------------------------------------------------ ! Compute qflx_sat_excess_surf ! - ! assume qinmax (maximum infiltration capacity) is large relative to + ! assume qinmax (maximum infiltration rate) is large relative to ! qflx_rain_plus_snomelt in control ! ------------------------------------------------------------------------ diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index 44d35bf3f1..a733f97212 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -351,7 +351,7 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & !----------------------------------------------------------------------- associate( & - qinmax => infiltration_excess_runoff_inst%qinmax_col , & ! Input: [real(r8) (:)] maximum infiltration capacity (mm H2O /s) + qinmax => infiltration_excess_runoff_inst%qinmax_col , & ! Input: [real(r8) (:)] maximum infiltration rate (mm H2O /s) frac_h2osfc => waterstate_inst%frac_h2osfc_col , & ! Input: [real(r8) (:) ] fraction of ground covered by surface water (0 to 1) frac_h2osfc_nosnow => waterstate_inst%frac_h2osfc_nosnow_col, & ! Input: [real(r8) (:) ] col fractional area with surface water greater than zero (if no snow present) @@ -517,7 +517,7 @@ subroutine QflxH2osfcDrain(bounds, num_hydrologyc, filter_hydrologyc, & integer , intent(in) :: h2osfcflag ! true => surface water is active real(r8) , intent(in) :: h2osfc( bounds%begc: ) ! surface water (mm) real(r8) , intent(in) :: frac_h2osfc( bounds%begc: ) ! fraction of ground covered by surface water (0 to 1) - real(r8) , intent(in) :: qinmax( bounds%begc: ) ! maximum infiltration capacity (mm H2O /s) + real(r8) , intent(in) :: qinmax( bounds%begc: ) ! maximum infiltration rate (mm H2O /s) real(r8) , intent(inout) :: qflx_h2osfc_drain( bounds%begc: ) ! bottom drainage from h2osfc (mm H2O /s) logical , intent(inout) :: truncate_h2osfc_to_zero( bounds%begc: ) ! whether h2osfc should be truncated to 0 to correct for roundoff errors, in order to maintain bit-for-bit the same answers as the old code ! From 6c979ad9d57e27e21b113dfa06661a2ea9d1dc5e Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Wed, 23 Aug 2017 16:26:21 -0600 Subject: [PATCH 09/11] Tweak comments based on feedback from Martyn Clark --- src/biogeophys/InfiltrationExcessRunoffMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/biogeophys/InfiltrationExcessRunoffMod.F90 b/src/biogeophys/InfiltrationExcessRunoffMod.F90 index 690ec7c10b..ddd50136d5 100644 --- a/src/biogeophys/InfiltrationExcessRunoffMod.F90 +++ b/src/biogeophys/InfiltrationExcessRunoffMod.F90 @@ -294,9 +294,9 @@ subroutine ComputeQinmaxVic(bounds, num_hydrologyc, filter_hydrologyc, & integer :: fc, c real(r8) :: dtime ! land model time step (sec) real(r8) :: top_icefrac ! ice fraction in top VIC layers - real(r8) :: max_infil ! max infiltration capacity in VIC (mm) + real(r8) :: max_infil ! max infiltration capacity using the VIC parameterization (mm) real(r8) :: i_0 ! average soil moisture in top VIC layers (mm) - real(r8) :: rsurf_vic ! VIC surface runoff + real(r8) :: rsurf_vic ! surface runoff based on the VIC parameterization real(r8) :: basis ! variable soil moisture holding capacity in top VIC layers for runoff calculation character(len=*), parameter :: subname = 'ComputeQinmaxVic' From 5c75419d917feef8d2ef17c93d4e30ee6ff716c5 Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Wed, 23 Aug 2017 20:21:28 -0600 Subject: [PATCH 10/11] Move setting of qflx_infl to its own subroutine This is cleaner and safer than having it set initially in one routine and then updated in a later routine. --- src/biogeophys/HydrologyNoDrainageMod.F90 | 6 ++- src/biogeophys/SoilHydrologyMod.F90 | 49 +++++++++++++++++------ src/biogeophys/WaterfluxType.F90 | 4 ++ 3 files changed, 46 insertions(+), 13 deletions(-) diff --git a/src/biogeophys/HydrologyNoDrainageMod.F90 b/src/biogeophys/HydrologyNoDrainageMod.F90 index 22def389ec..299c43b004 100644 --- a/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -62,7 +62,8 @@ subroutine HydrologyNoDrainage(bounds, & use SnowHydrologyMod , only : SnowCompaction, CombineSnowLayers, DivideSnowLayers, SnowCapping use SnowHydrologyMod , only : SnowWater, BuildSnowFilter use SoilHydrologyMod , only : CLMVICMap, SetSoilWaterFractions - use SoilHydrologyMod , only : SetQflxInputs, RouteInfiltrationExcess, UpdateH2osfc, TotalSurfaceRunoff + use SoilHydrologyMod , only : SetQflxInputs, RouteInfiltrationExcess, UpdateH2osfc + use SoilHydrologyMod , only : Infiltration, TotalSurfaceRunoff use SoilHydrologyMod , only : UpdateUrbanPonding use SoilHydrologyMod , only : WaterTable, PerchedWaterTable use SoilHydrologyMod , only : ThetaBasedWaterTable, RenewCondensation @@ -202,6 +203,9 @@ subroutine HydrologyNoDrainage(bounds, & energyflux_inst, soilhydrology_inst, & waterflux_inst, waterstate_inst) + call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & + waterflux_inst) + call TotalSurfaceRunoff(bounds, num_hydrologyc, filter_hydrologyc, & num_urbanc, filter_urbanc, & waterflux_inst, soilhydrology_inst, saturated_excess_runoff_inst, waterstate_inst) diff --git a/src/biogeophys/SoilHydrologyMod.F90 b/src/biogeophys/SoilHydrologyMod.F90 index a733f97212..c16cef0e51 100644 --- a/src/biogeophys/SoilHydrologyMod.F90 +++ b/src/biogeophys/SoilHydrologyMod.F90 @@ -38,6 +38,7 @@ module SoilHydrologyMod public :: SetSoilWaterFractions ! Set diagnostic variables related to the fraction of water and ice in each layer public :: SetQflxInputs ! Set the flux of water into the soil from the top public :: UpdateH2osfc ! Calculate fluxes out of h2osfc and update the h2osfc state + public :: Infiltration ! Calculate total infiltration public :: TotalSurfaceRunoff ! Calculate total surface runoff public :: UpdateUrbanPonding ! Update the state variable representing ponding on urban surfaces public :: WaterTable ! Calculate water table before imposing drainage @@ -281,7 +282,7 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & !----------------------------------------------------------------------- associate( & - qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) + qflx_in_soil_limited => waterflux_inst%qflx_in_soil_limited_col , & ! Output: [real(r8) (:) ] surface input to soil, limited by max infiltration rate (mm H2O /s) qflx_in_h2osfc => waterflux_inst%qflx_in_h2osfc_col , & ! Output: [real(r8) (:) ] total surface input to h2osfc qflx_infl_excess_surf => waterflux_inst%qflx_infl_excess_surf_col , & ! Output: [real(r8) (:) ] surface runoff due to infiltration excess (mm H2O /s) qflx_in_soil => waterflux_inst%qflx_in_soil_col , & ! Input: [real(r8) (:) ] surface input to soil (mm/s) @@ -295,7 +296,7 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) if (lun%itype(col%landunit(c)) == istsoil .or. lun%itype(col%landunit(c))==istcrop) then - qflx_infl(c) = qflx_in_soil(c) - qflx_infl_excess(c) + qflx_in_soil_limited(c) = qflx_in_soil(c) - qflx_infl_excess(c) if (h2osfcflag /= 0) then qflx_in_h2osfc(c) = qflx_top_soil_to_h2osfc(c) + qflx_infl_excess(c) qflx_infl_excess_surf(c) = 0._r8 @@ -306,7 +307,7 @@ subroutine RouteInfiltrationExcess(bounds, num_hydrologyc, filter_hydrologyc, & end if else ! non-vegetated landunits (i.e. urban) use original CLM4 code - qflx_infl(c) = qflx_in_soil(c) + qflx_in_soil_limited(c) = qflx_in_soil(c) qflx_in_h2osfc(c) = 0._r8 qflx_infl_excess_surf(c) = 0._r8 end if @@ -344,9 +345,6 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & integer :: c,l,fc ! indices real(r8) :: dtime ! land model time step (sec) real(r8) :: h2osfc_partial(bounds%begc:bounds%endc) ! partially-updated h2osfc - ! FIXME(wjs, 2017-08-22) put qflx_h2osfc_drain in waterflux_type (this will be - ! needed when we move the update of qflx_infl to a separate routine) - real(r8) :: qflx_h2osfc_drain(bounds%begc:bounds%endc) ! bottom drainage from h2osfc logical :: truncate_h2osfc_to_zero(bounds%begc:bounds%endc) !----------------------------------------------------------------------- @@ -359,7 +357,7 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & qflx_in_h2osfc => waterflux_inst%qflx_in_h2osfc_col , & ! Input: [real(r8) (:) ] total surface input to h2osfc qflx_h2osfc_surf => waterflux_inst%qflx_h2osfc_surf_col , & ! Output: [real(r8) (:) ] surface water runoff (mm H2O /s) - qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) + qflx_h2osfc_drain => waterflux_inst%qflx_h2osfc_drain_col , & ! Output: [real(r8) (:) ] bottom drainage from h2osfc (mm H2O /s) h2osfc_thresh => soilhydrology_inst%h2osfc_thresh_col, & ! Input: [real(r8) (:) ] level at which h2osfc "percolates" h2osfcflag => soilhydrology_inst%h2osfcflag & ! Input: integer @@ -415,16 +413,43 @@ subroutine UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, & end do - ! FIXME(wjs, 2017-08-22) Move this update to a separate routine that calculations - ! qflx_infl from its components. + end associate + + end subroutine UpdateH2osfc + + !----------------------------------------------------------------------- + subroutine Infiltration(bounds, num_hydrologyc, filter_hydrologyc, & + waterflux_inst) + ! + ! !DESCRIPTION: + ! Calculate total infiltration + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_hydrologyc ! number of column soil points in column filter + integer , intent(in) :: filter_hydrologyc(:) ! column filter for soil points + type(waterflux_type) , intent(inout) :: waterflux_inst + ! + ! !LOCAL VARIABLES: + integer :: fc, c + + character(len=*), parameter :: subname = 'Infiltration' + !----------------------------------------------------------------------- + + associate( & + qflx_infl => waterflux_inst%qflx_infl_col , & ! Output: [real(r8) (:) ] infiltration (mm H2O /s) + qflx_in_soil_limited => waterflux_inst%qflx_in_soil_limited_col , & ! Input: [real(r8) (:) ] surface input to soil, limited by max infiltration rate (mm H2O /s) + qflx_h2osfc_drain => waterflux_inst%qflx_h2osfc_drain_col & ! Input: [real(r8) (:) ] bottom drainage from h2osfc (mm H2O /s) + ) + do fc = 1, num_hydrologyc c = filter_hydrologyc(fc) - qflx_infl(c) = qflx_infl(c) + qflx_h2osfc_drain(c) + qflx_infl(c) = qflx_in_soil_limited(c) + qflx_h2osfc_drain(c) end do - end associate + end associate - end subroutine UpdateH2osfc + end subroutine Infiltration !----------------------------------------------------------------------- subroutine QflxH2osfcSurf(bounds, num_hydrologyc, filter_hydrologyc, & diff --git a/src/biogeophys/WaterfluxType.F90 b/src/biogeophys/WaterfluxType.F90 index 2398763a35..ff55279596 100644 --- a/src/biogeophys/WaterfluxType.F90 +++ b/src/biogeophys/WaterfluxType.F90 @@ -81,6 +81,8 @@ module WaterfluxType real(r8), pointer :: qflx_rain_plus_snomelt_col(:) ! col rain plus snow melt falling on the soil (mm/s) real(r8), pointer :: qflx_top_soil_col (:) ! col net water input into soil from top (mm/s) real(r8), pointer :: qflx_in_soil_col (:) ! col surface input to soil (mm/s) + real(r8), pointer :: qflx_in_soil_limited_col (:) ! col surface input to soil, limited by max infiltration rate (mm/s) + real(r8), pointer :: qflx_h2osfc_drain_col (:) ! col bottom drainage from h2osfc (mm/s) real(r8), pointer :: qflx_top_soil_to_h2osfc_col(:) ! col portion of qflx_top_soil going to h2osfc, minus evaporation (mm/s) real(r8), pointer :: qflx_in_h2osfc_col(:) ! col total surface input to h2osfc real(r8), pointer :: qflx_h2osfc_to_ice_col (:) ! col conversion of h2osfc to ice @@ -226,6 +228,8 @@ subroutine InitAllocate(this, bounds) allocate(this%qflx_rain_plus_snomelt_col(begc:endc)) ; this%qflx_rain_plus_snomelt_col(:) = nan allocate(this%qflx_top_soil_col (begc:endc)) ; this%qflx_top_soil_col (:) = nan allocate(this%qflx_in_soil_col (begc:endc)) ; this%qflx_in_soil_col (:) = nan + allocate(this%qflx_in_soil_limited_col (begc:endc)) ; this%qflx_in_soil_limited_col (:) = nan + allocate(this%qflx_h2osfc_drain_col (begc:endc)) ; this%qflx_h2osfc_drain_col (:) = nan allocate(this%qflx_top_soil_to_h2osfc_col(begc:endc)) ; this%qflx_top_soil_to_h2osfc_col(:) = nan allocate(this%qflx_in_h2osfc_col (begc:endc)) ; this%qflx_in_h2osfc_col(:) = nan allocate(this%qflx_h2osfc_to_ice_col (begc:endc)) ; this%qflx_h2osfc_to_ice_col (:) = nan From 406d86bd33c6d023209670d9fe3204b65ed74032 Mon Sep 17 00:00:00 2001 From: Bill Sacks Date: Thu, 24 Aug 2017 12:55:05 -0600 Subject: [PATCH 11/11] Update changelog --- ChangeLog_branch | 80 ++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/ChangeLog_branch b/ChangeLog_branch index 0a2617eef8..8f65aed1f1 100644 --- a/ChangeLog_branch +++ b/ChangeLog_branch @@ -1,4 +1,84 @@ =============================================================== +Tag name: ctsm_n08_clm4_5_16_r249 +Originator(s): sacks +Date: Aug 24, 2017 +One-line Summary: Finish Infiltration modularization - mainly h2osfc updates + +Purpose of changes +------------------ + +Finish the modularization of subroutine Infiltration. The main focus here is +modularizing the calculation of fluxes out of h2osfc and the updates of the +h2osfc state variable. Flux calculations and state updates are now done +separately. Where a flux calculation depends on a partially-updated version of +h2osfc, this is now made explicit in a way that could facilitate removing this +partial update later. + + +Notes of particular relevance for developers: (including Code reviews and testing) +--------------------------------------------- + +CLM testing: + + build-namelist tests: + + yellowstone - not run + + unit-tests (components/clm/src): + + yellowstone - pass + + tools-tests (components/clm/test/tools): + + yellowstone - not run + + PTCLM testing (components/clm/tools/shared/PTCLM/test): + + yellowstone - not run + + regular tests (aux_clm): + + yellowstone_intel - pass + yellowstone_pgi - pass + yellowstone_gnu - pass + cheyenne_intel - pass + hobart_nag - pass + hobart_pgi - pass + hobart_intel - pass + +CLM tag used for the baseline comparisons: ctsm_n07_clm4_5_16_r249 + + +Answer changes +-------------- + +Changes answers relative to baseline: NO + + +Detailed list of changes +------------------------ + +List any svn externals directories updated (cime, rtm, mosart, cism, etc.): none + +List all files eliminated: none + +List all files added and what they do: none + +List all existing files that have been modified, and describe the changes: + +========= Main changes, as described above +M components/clm/src/biogeophys/SoilHydrologyMod.F90 +M components/clm/src/biogeophys/HydrologyNoDrainageMod.F90 +M components/clm/src/biogeophys/HydrologyDrainageMod.F90 +M components/clm/src/biogeophys/WaterfluxType.F90 + +========= Tweak some comments based on feedback on + https://github.com/NCAR/clm-ctsm/pull/2 +M components/clm/src/biogeophys/InfiltrationExcessRunoffMod.F90 +M components/clm/src/biogeophys/SaturatedExcessRunoffMod.F90 + +=============================================================== +=============================================================== Tag name: ctsm_n07_clm4_5_16_r249 Originator(s): sacks Date: Aug 22, 2017