From 654b53dd55145bd374c0b3722485b23be4cff26f Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Mon, 7 Dec 2020 17:14:05 -0700 Subject: [PATCH 01/19] Grid cell-level error check for H2O: non-transient simulations only Works for non-transient simulations. In this attempt I convert most column-level variables to grid cell-level locally in subr. BalanceCheck. Using grid cell-level variables calculated in subr. lnd2atm did not work because subr. lnd2atm is called after subr. BalanceCheck. --- src/biogeophys/BalanceCheckMod.F90 | 312 +++++++++++++++++++++++++--- src/biogeophys/WaterBalanceType.F90 | 18 +- src/main/clm_driver.F90 | 18 +- src/main/lnd2atmMod.F90 | 4 +- 4 files changed, 315 insertions(+), 37 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index c2d58f711e..bc73b68265 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -11,7 +11,7 @@ module BalanceCheckMod use decompMod , only : bounds_type use abortutils , only : endrun use clm_varctl , only : iulog - use clm_varcon , only : namep, namec + use clm_varcon , only : namep, namec, nameg use clm_varpar , only : nlevsoi use GetGlobalValuesMod , only : GetGlobalIndex use atm2lndType , only : atm2lnd_type @@ -23,6 +23,7 @@ module BalanceCheckMod use WaterDiagnosticBulkType, only : waterdiagnosticbulk_type use WaterDiagnosticType, only : waterdiagnostic_type use Wateratm2lndType , only : wateratm2lnd_type +! use Waterlnd2atmType , only : waterlnd2atm_type ! slevis: place holder use WaterBalanceType , only : waterbalance_type use WaterFluxType , only : waterflux_type use WaterType , only : water_type @@ -43,8 +44,9 @@ module BalanceCheckMod ! !PUBLIC MEMBER FUNCTIONS: public :: BalanceCheckInit ! Initialization of Water and energy balance check - public :: BeginWaterBalance ! Initialize water balance check - public :: BalanceCheck ! Water and energy balance check + public :: BeginWaterGridcellBalance ! Initialize grid cell-level water balance check + public :: BeginWaterColumnBalance ! Initialize column-level water balance check + public :: BalanceCheck ! Water & energy balance checks public :: GetBalanceCheckSkipSteps ! Get the number of steps to skip for the balance check public :: BalanceCheckClean ! Clean up for BalanceCheck @@ -54,7 +56,8 @@ module BalanceCheckMod ! ! !PRIVATE MEMBER FUNCTIONS: - private :: BeginWaterBalanceSingle ! Initialize water balance check for bulk or a single tracer + private :: BeginWaterGridcellBalanceSingle ! Initialize grid cell-level water balance check for bulk or a single tracer + private :: BeginWaterColumnBalanceSingle ! Initialize column-level water balance check for bulk or a single tracer character(len=*), parameter, private :: sourcefile = & __FILE__ @@ -120,7 +123,46 @@ end function GetBalanceCheckSkipSteps !----------------------------------------------------------------------- !----------------------------------------------------------------------- - subroutine BeginWaterBalance(bounds, & + subroutine BeginWaterGridcellBalance(bounds, & + num_nolakec, filter_nolakec, num_lakec, filter_lakec, & + water_inst, soilhydrology_inst, & + use_aquifer_layer) + ! + ! !DESCRIPTION: + ! Initialize grid cell-level water balance at beginning of time step + ! for bulk water and each water tracer + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_nolakec ! number of column non-lake points in column filter + integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points + integer , intent(in) :: num_lakec ! number of column lake points in column filter + integer , intent(in) :: filter_lakec(:) ! column filter for lake points + type(water_type) , intent(inout) :: water_inst + type(soilhydrology_type), intent(in) :: soilhydrology_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run + ! + ! !LOCAL VARIABLES: + integer :: i + + character(len=*), parameter :: subname = 'BeginWaterGridcellBalance' + !----------------------------------------------------------------------- + + do i = water_inst%bulk_and_tracers_beg, water_inst%bulk_and_tracers_end + call BeginWaterGridcellBalanceSingle(bounds, & + num_nolakec, filter_nolakec, & + num_lakec, filter_lakec, & + soilhydrology_inst, & + water_inst%bulk_and_tracers(i)%waterstate_inst, & + water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & + water_inst%bulk_and_tracers(i)%waterbalance_inst, & + use_aquifer_layer = use_aquifer_layer) + end do + + end subroutine BeginWaterGridcellBalance + + !----------------------------------------------------------------------- + subroutine BeginWaterColumnBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & water_inst, soilhydrology_inst, & use_aquifer_layer) @@ -142,11 +184,11 @@ subroutine BeginWaterBalance(bounds, & ! !LOCAL VARIABLES: integer :: i - character(len=*), parameter :: subname = 'BeginWaterBalance' + character(len=*), parameter :: subname = 'BeginWaterColumnBalance' !----------------------------------------------------------------------- do i = water_inst%bulk_and_tracers_beg, water_inst%bulk_and_tracers_end - call BeginWaterBalanceSingle(bounds, & + call BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, & num_lakec, filter_lakec, & soilhydrology_inst, & @@ -156,10 +198,85 @@ subroutine BeginWaterBalance(bounds, & use_aquifer_layer = use_aquifer_layer) end do - end subroutine BeginWaterBalance + end subroutine BeginWaterColumnBalance !----------------------------------------------------------------------- - subroutine BeginWaterBalanceSingle(bounds, & + subroutine BeginWaterGridcellBalanceSingle(bounds, & + num_nolakec, filter_nolakec, num_lakec, filter_lakec, & + soilhydrology_inst, waterstate_inst, waterdiagnostic_inst, waterbalance_inst, & + use_aquifer_layer) + ! + ! !DESCRIPTION: + ! Initialize grid cell-level water balance at beginning of time step + ! for bulk or a single tracer + ! + ! !USES: + use subgridAveMod, only: c2g + ! + ! !ARGUMENTS: + type(bounds_type) , intent(in) :: bounds + integer , intent(in) :: num_nolakec ! number of column non-lake points in column filter + integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points + integer , intent(in) :: num_lakec ! number of column lake points in column filter + integer , intent(in) :: filter_lakec(:) ! column filter for lake points + type(soilhydrology_type) , intent(in) :: soilhydrology_inst + class(waterstate_type) , intent(inout) :: waterstate_inst + class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst + class(waterbalance_type) , intent(inout) :: waterbalance_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run + ! + ! !LOCAL VARIABLES: + integer :: c, j, fc ! indices + integer :: begc, endc, begg, endg ! bounds + !----------------------------------------------------------------------- + + associate( & + zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) + zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) + aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) + wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) + begwb_col => waterbalance_inst%begwb_col, & ! Output: [real(r8) (:) ] column-level water mass begining of the time step + begwb_grc => waterbalance_inst%begwb_grc & ! Output: [real(r8) (:) ] grid cell-level water mass begining of the time step + ) + + begc = bounds%begc + endc = bounds%endc + begg = bounds%begg + endg = bounds%endg + + ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake + if(use_aquifer_layer) then + do fc = 1, num_nolakec + c = filter_nolakec(fc) + if (col%hydrologically_active(c)) then + if(zwt(c) <= zi(c,nlevsoi)) then + wa(c) = aquifer_water_baseline + end if + end if + end do + endif + + ! NOTES subroutines Compute*Mass* are in TotalWaterAndHeatMod.F90 + ! endwb is calculated in HydrologyDrainageMod & LakeHydrologyMod + call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, & + waterstate_inst, waterdiagnostic_inst, & + subtract_dynbal_baselines = .false., & + water_mass = begwb_col(begc:endc)) + + call ComputeWaterMassLake(bounds, num_lakec, filter_lakec, & + waterstate_inst, & + subtract_dynbal_baselines = .false., & + water_mass = begwb_col(begc:endc)) + + call c2g(bounds, begwb_col(begc:endc), begwb_grc(begg:endg), & + c2l_scale_type='urbanf', l2g_scale_type='unity') + + end associate + + end subroutine BeginWaterGridcellBalanceSingle + + !----------------------------------------------------------------------- + subroutine BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & soilhydrology_inst, waterstate_inst, waterdiagnostic_inst, waterbalance_inst, & use_aquifer_layer) @@ -223,7 +340,7 @@ subroutine BeginWaterBalanceSingle(bounds, & end associate - end subroutine BeginWaterBalanceSingle + end subroutine BeginWaterColumnBalanceSingle !----------------------------------------------------------------------- subroutine BalanceCheck( bounds, & @@ -231,6 +348,7 @@ subroutine BalanceCheck( bounds, & atm2lnd_inst, solarabs_inst, waterflux_inst, waterstate_inst, & waterdiagnosticbulk_inst, waterbalance_inst, wateratm2lnd_inst, & surfalb_inst, energyflux_inst, canopystate_inst) +! waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst) ! ! !DESCRIPTION: ! This subroutine accumulates the numerical truncation errors of the water @@ -251,7 +369,7 @@ subroutine BalanceCheck( bounds, & use clm_time_manager , only : get_step_size_real, get_nstep use clm_time_manager , only : get_nstep_since_startup_or_lastDA_restart_or_pause use CanopyStateType , only : canopystate_type - use subgridAveMod + use subgridAveMod ! , only : c2g ? ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds @@ -263,6 +381,7 @@ subroutine BalanceCheck( bounds, & class(waterstate_type), intent(in) :: waterstate_inst type(waterdiagnosticbulk_type), intent(in) :: waterdiagnosticbulk_inst class(waterbalance_type), intent(inout) :: waterbalance_inst +! class(waterlnd2atm_type), intent(in) :: waterlnd2atm_inst class(wateratm2lnd_type) , intent(in) :: wateratm2lnd_inst type(surfalb_type) , intent(in) :: surfalb_inst type(energyflux_type) , intent(inout) :: energyflux_inst @@ -277,6 +396,18 @@ subroutine BalanceCheck( bounds, & real(r8) :: forc_rain_col(bounds%begc:bounds%endc) ! column level rain rate [mm/s] real(r8) :: forc_snow_col(bounds%begc:bounds%endc) ! column level snow rate [mm/s] real(r8) :: h2osno_total(bounds%begc:bounds%endc) ! total snow water [mm H2O] + real(r8) :: endwb_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_irrig_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_glcice_dyn_water_flux_locgrc(bounds%begg:bounds%endg) ! water flux needed for balance check due to glc_dyn_runoff_routing [mm H2O/s] (positive means addition of water to the system) + real(r8) :: qflx_evap_tot_locgrc(bounds%begg:bounds%endg) ! grid cell level total evapotranspiration [mm/s] + real(r8) :: qflx_surf_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_drain_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_drain_perched_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_qrgwl_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_ice_runoff_snwcp_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_ice_runoff_xs_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now + real(r8) :: qflx_snwcp_discarded_liq_locgrc(bounds%begg:bounds%endg) ! excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] + real(r8) :: qflx_snwcp_discarded_ice_locgrc(bounds%begg:bounds%endg) ! excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] real(r8) :: errh2o_max_val ! Maximum value of error in water conservation error over all columns [mm H2O] real(r8) :: errh2osno_max_val ! Maximum value of error in h2osno conservation error over all columns [kg m-2] @@ -294,6 +425,9 @@ subroutine BalanceCheck( bounds, & associate( & forc_solad => atm2lnd_inst%forc_solad_grc , & ! Input: [real(r8) (:,:) ] direct beam radiation (vis=forc_sols , nir=forc_soll ) forc_solai => atm2lnd_inst%forc_solai_grc , & ! Input: [real(r8) (:,:) ] diffuse radiation (vis=forc_solsd, nir=forc_solld) + forc_rain_grc => wateratm2lnd_inst%forc_rain_not_downscaled_grc , & ! Input: [real(r8) (:)] grid cell-level rain rate [mm/s] + forc_snow_grc => wateratm2lnd_inst%forc_snow_not_downscaled_grc , & ! Input: [real(r8) (:)] grid cell-level snow rate [mm/s] + qflx_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] total grid cell-level runoff from river model forc_rain => wateratm2lnd_inst%forc_rain_downscaled_col , & ! Input: [real(r8) (:) ] rain rate [mm/s] forc_snow => wateratm2lnd_inst%forc_snow_downscaled_col , & ! Input: [real(r8) (:) ] snow rate [mm/s] forc_lwrad => atm2lnd_inst%forc_lwrad_downscaled_col , & ! Input: [real(r8) (:) ] downward infrared (longwave) radiation (W/m**2) @@ -302,10 +436,13 @@ subroutine BalanceCheck( bounds, & frac_sno_eff => waterdiagnosticbulk_inst%frac_sno_eff_col , & ! Input: [real(r8) (:) ] effective snow fraction frac_sno => waterdiagnosticbulk_inst%frac_sno_col , & ! Input: [real(r8) (:) ] fraction of ground covered by snow (0 to 1) snow_depth => waterdiagnosticbulk_inst%snow_depth_col , & ! Input: [real(r8) (:) ] snow height (m) - begwb => waterbalance_inst%begwb_col , & ! Input: [real(r8) (:) ] water mass begining of the time step - errh2o => waterbalance_inst%errh2o_col , & ! Output: [real(r8) (:) ] water conservation error (mm H2O) + begwb_grc => waterbalance_inst%begwb_grc , & ! Input: [real(r8) (:) ] grid cell-level water mass begining of the time step +! endwb_grc => waterbalance_inst%endwb_grc , & ! Output: [real(r8) (:) ] grid cell-level water mass end of the time step + errh2o_grc => waterbalance_inst%errh2o_grc , & ! Output: [real(r8) (:) ] grid cell-level water conservation error (mm H2O) + begwb_col => waterbalance_inst%begwb_col , & ! Input: [real(r8) (:) ] column-level water mass begining of the time step + endwb_col => waterbalance_inst%endwb_col , & ! Output: [real(r8) (:) ] column-level water mass end of the time step + errh2o_col => waterbalance_inst%errh2o_col , & ! Output: [real(r8) (:) ] column-level water conservation error (mm H2O) errh2osno => waterbalance_inst%errh2osno_col , & ! Output: [real(r8) (:) ] error in h2osno (kg m-2) - endwb => waterbalance_inst%endwb_col , & ! Output: [real(r8) (:) ] water mass end of the time step snow_sources => waterbalance_inst%snow_sources_col , & ! Output: [real(r8) (:) ] snow sources (mm H2O /s) snow_sinks => waterbalance_inst%snow_sinks_col , & ! Output: [real(r8) (:) ] snow sinks (mm H2O /s) qflx_liq_grnd_col => waterflux_inst%qflx_liq_grnd_col , & ! Input: [real(r8) (:) ] liquid on ground after interception (mm H2O/s) [+] @@ -315,6 +452,7 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_liq => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:) ] excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` qflx_snwcp_discarded_ice => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:) ] excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` qflx_evap_tot => waterflux_inst%qflx_evap_tot_col , & ! Input: [real(r8) (:) ] qflx_evap_soi + qflx_evap_can + qflx_tran_veg +! qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc, & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_soliddew_to_top_layer => waterflux_inst%qflx_soliddew_to_top_layer_col , & ! Input: [real(r8) (:) ] rate of solid water deposited on top soil or snow layer (frost) (mm H2O /s) [+] qflx_solidevap_from_top_layer => waterflux_inst%qflx_solidevap_from_top_layer_col, & ! Input: [real(r8) (:) ] rate of ice evaporated from top soil or snow layer (sublimation) (mm H2O /s) [+] qflx_liqevap_from_top_layer => waterflux_inst%qflx_liqevap_from_top_layer_col , & ! Input: [real(r8) (:) ] rate of liquid water evaporated from top soil or snow layer (mm H2O/s) [+] @@ -323,16 +461,24 @@ subroutine BalanceCheck( bounds, & qflx_snow_h2osfc => waterflux_inst%qflx_snow_h2osfc_col , & ! Input: [real(r8) (:) ] snow falling on surface water (mm/s) qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice qflx_drain_perched => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) +! qflx_rofliq_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) qflx_floodc => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] total runoff due to flooding qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack +! qflx_liq_dynbal_grc => waterflux_inst%qflx_liq_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder +! qflx_ice_dynbal_grc => waterflux_inst%qflx_ice_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder +! qflx_runoff_col => waterflux_inst%qflx_runoff_col , & ! total runoff (mm H2O / s) slevis: place holder qflx_surf => waterflux_inst%qflx_surf_col , & ! Input: [real(r8) (:) ] surface runoff (mm H2O /s) +! qflx_rofliq_qsur_grc => waterlnd2atm_inst%qflx_rofliq_qsur_grc , & ! Input: [real(r8) (:) ] grid cell-level surface runoff (mm H20 /s) qflx_qrgwl => waterflux_inst%qflx_qrgwl_col , & ! Input: [real(r8) (:) ] qflx_surf at glaciers, wetlands, lakes +! qflx_rofliq_qgwl_grc => waterlnd2atm_inst%qflx_rofliq_qgwl_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_surf at glaciers, wetlands, lakes qflx_drain => waterflux_inst%qflx_drain_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) +! qflx_rofliq_qsub_grc => waterlnd2atm_inst%qflx_rofliq_qsub_grc , & ! Input: [real(r8) (:) ] grid cell-level drainage (mm H20 /s) qflx_ice_runoff_snwcp => waterflux_inst%qflx_ice_runoff_snwcp_col, & ! Input: [real(r8) (:) ] solid runoff from snow capping (mm H2O /s) qflx_ice_runoff_xs => waterflux_inst%qflx_ice_runoff_xs_col , & ! Input: [real(r8) (:) ] solid runoff from excess ice in soil (mm H2O /s) qflx_sl_top_soil => waterflux_inst%qflx_sl_top_soil_col , & ! Input: [real(r8) (:) ] liquid water + ice from layer above soil to top soil layer or sent to qflx_qrgwl (mm H2O/s) qflx_sfc_irrig => waterflux_inst%qflx_sfc_irrig_col , & ! Input: [real(r8) (:) ] irrigation flux (mm H2O /s) +! qirrig_grc => waterlnd2atm_inst%qirrig_grc , & ! Input: [real(r8) (:) ] grid cell-level irrigation flux (mm H20 /s) qflx_glcice_dyn_water_flux => waterflux_inst%qflx_glcice_dyn_water_flux_col, & ! Input: [real(r8) (:)] water flux needed for balance check due to glc_dyn_runoff_routing (mm H2O/s) (positive means addition of water to the system) eflx_lwrad_out => energyflux_inst%eflx_lwrad_out_patch , & ! Input: [real(r8) (:) ] emitted infrared (longwave) radiation (W/m**2) @@ -393,14 +539,14 @@ subroutine BalanceCheck( bounds, & end if end do - ! Water balance check + ! Water balance check at the column level do c = bounds%begc, bounds%endc ! add qflx_drain_perched and qflx_flood if (col%active(c)) then - errh2o(c) = endwb(c) - begwb(c) & + errh2o_col(c) = endwb_col(c) - begwb_col(c) & - (forc_rain_col(c) & + forc_snow_col(c) & + qflx_floodc(c) & @@ -418,31 +564,31 @@ subroutine BalanceCheck( bounds, & else - errh2o(c) = 0.0_r8 + errh2o_col(c) = 0.0_r8 end if end do - errh2o_max_val = maxval(abs(errh2o(bounds%begc:bounds%endc))) + errh2o_max_val = maxval(abs(errh2o_col(bounds%begc:bounds%endc))) if (errh2o_max_val > h2o_warning_thresh) then - indexc = maxloc( abs(errh2o(bounds%begc:bounds%endc)), 1 ) + bounds%begc -1 - write(iulog,*)'WARNING: water balance error ',& + indexc = maxloc( abs(errh2o_col(bounds%begc:bounds%endc)), 1 ) + bounds%begc -1 + write(iulog,*)'WARNING: column-level water balance error ',& ' nstep= ',nstep, & ' local indexc= ',indexc,& ! ' global indexc= ',GetGlobalIndex(decomp_index=indexc, clmlevel=namec), & - ' errh2o= ',errh2o(indexc) + ' errh2o= ',errh2o_col(indexc) if ((errh2o_max_val > error_thresh) .and. (DAnstep > skip_steps)) then write(iulog,*)'clm urban model is stopping - error is greater than 1e-5 (mm)' write(iulog,*)'nstep = ',nstep - write(iulog,*)'errh2o = ',errh2o(indexc) + write(iulog,*)'errh2o_col = ',errh2o_col(indexc) write(iulog,*)'forc_rain = ',forc_rain_col(indexc)*dtime write(iulog,*)'forc_snow = ',forc_snow_col(indexc)*dtime - write(iulog,*)'endwb = ',endwb(indexc) - write(iulog,*)'begwb = ',begwb(indexc) + write(iulog,*)'endwb_col = ',endwb_col(indexc) + write(iulog,*)'begwb_col = ',begwb_col(indexc) write(iulog,*)'qflx_evap_tot = ',qflx_evap_tot(indexc)*dtime write(iulog,*)'qflx_sfc_irrig = ',qflx_sfc_irrig(indexc)*dtime @@ -455,8 +601,8 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice(indexc)*dtime write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq(indexc)*dtime - write(iulog,*)'deltawb = ',endwb(indexc)-begwb(indexc) - write(iulog,*)'deltawb/dtime = ',(endwb(indexc)-begwb(indexc))/dtime + write(iulog,*)'deltawb = ',endwb_col(indexc)-begwb_col(indexc) + write(iulog,*)'deltawb/dtime = ',(endwb_col(indexc)-begwb_col(indexc))/dtime if (.not.(col%itype(indexc) == icol_roof .or. & col%itype(indexc) == icol_road_imperv .or. & @@ -473,7 +619,118 @@ subroutine BalanceCheck( bounds, & end if - ! Snow balance check + ! Water balance check at the grid cell level + + call c2g(bounds, & + endwb_col(bounds%begc:bounds%endc), & + endwb_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type='urbanf', l2g_scale_type='unity') + call c2g( bounds, & + qflx_sfc_irrig(bounds%begc:bounds%endc), & + qflx_irrig_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_glcice_dyn_water_flux(bounds%begc:bounds%endc), & + qflx_glcice_dyn_water_flux_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_evap_tot(bounds%begc:bounds%endc), & + qflx_evap_tot_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_surf(bounds%begc:bounds%endc), & + qflx_surf_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_qrgwl(bounds%begc:bounds%endc), & + qflx_qrgwl_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_drain(bounds%begc:bounds%endc), & + qflx_drain_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_drain_perched(bounds%begc:bounds%endc), & + qflx_drain_perched_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_ice_runoff_snwcp(bounds%begc:bounds%endc), & + qflx_ice_runoff_snwcp_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_ice_runoff_xs(bounds%begc:bounds%endc), & + qflx_ice_runoff_xs_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_snwcp_discarded_liq(bounds%begc:bounds%endc), & + qflx_snwcp_discarded_liq_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call c2g( bounds, & + qflx_snwcp_discarded_ice(bounds%begc:bounds%endc), & + qflx_snwcp_discarded_ice_locgrc(bounds%begg:bounds%endg), & + c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + + do g = bounds%begg, bounds%endg + errh2o_grc(g) = endwb_locgrc(g) - begwb_grc(g) & + - (forc_rain_grc(g) & + + forc_snow_grc(g) & + + qflx_flood_grc(g) & + + qflx_irrig_locgrc(g) & + + qflx_glcice_dyn_water_flux_locgrc(g) & + - qflx_evap_tot_locgrc(g) & + - qflx_surf_locgrc(g) & + - qflx_qrgwl_locgrc(g) & + - qflx_drain_locgrc(g) & + - qflx_drain_perched_locgrc(g) & + - qflx_ice_runoff_snwcp_locgrc(g) & + - qflx_ice_runoff_xs_locgrc(g) & + - qflx_snwcp_discarded_liq_locgrc(g) & + - qflx_snwcp_discarded_ice_locgrc(g)) * dtime + end do + + errh2o_max_val = maxval(abs(errh2o_grc(bounds%begg:bounds%endg))) + + if (errh2o_max_val > h2o_warning_thresh) then + + indexg = maxloc( abs(errh2o_grc(bounds%begg:bounds%endg)), 1 ) + bounds%begg -1 + write(iulog,*)'WARNING: grid cell-level water balance error ',& + ' nstep= ',nstep, & + ' local indexg= ',indexg,& + ' errh2o_grc= ',errh2o_grc(indexg) + if ((errh2o_max_val > error_thresh) .and. (DAnstep > skip_steps)) then + + write(iulog,*)'clm model is stopping - error is greater than 1e-5 (mm)' + write(iulog,*)'nstep = ',nstep + write(iulog,*)'errh2o_grc = ',errh2o_grc(indexg) + write(iulog,*)'forc_rain = ',forc_rain_grc(indexg)*dtime + write(iulog,*)'forc_snow = ',forc_snow_grc(indexg)*dtime + write(iulog,*)'endwb_loc = ',endwb_locgrc(indexg) + write(iulog,*)'begwb_grc = ',begwb_grc(indexg) + + write(iulog,*)'qflx_evap_tot_loc = ',qflx_evap_tot_locgrc(indexg)*dtime + write(iulog,*)'qflx_irrig_loc = ',qflx_irrig_locgrc(indexg)*dtime + write(iulog,*)'qflx_surf_loc = ',qflx_surf_locgrc(indexg)*dtime + write(iulog,*)'qflx_qrgwl_loc = ',qflx_qrgwl_locgrc(indexg)*dtime + write(iulog,*)'qflx_drain_loc = ',qflx_drain_locgrc(indexg)*dtime + write(iulog,*)'qflx_ice_runoff_snwcp_loc = ',qflx_ice_runoff_snwcp_locgrc(indexg)*dtime + write(iulog,*)'qflx_ice_runoff_xs_loc = ',qflx_ice_runoff_xs_locgrc(indexg)*dtime + write(iulog,*)'qflx_snwcp_discarded_ice_loc = ',qflx_snwcp_discarded_ice_locgrc(indexg)*dtime + write(iulog,*)'qflx_snwcp_discarded_liq_loc = ',qflx_snwcp_discarded_liq_locgrc(indexg)*dtime + write(iulog,*)'deltawb = ',endwb_locgrc(indexg)-begwb_grc(indexg) + write(iulog,*)'deltawb/dtime = ',(endwb_locgrc(indexg)-begwb_grc(indexg))/dtime + write(iulog,*)'qflx_drain_perched_loc = ',qflx_drain_perched_locgrc(indexg)*dtime + write(iulog,*)'qflx_flood = ',qflx_flood_grc(indexg)*dtime + write(iulog,*)'qflx_glcice_dyn_water_flux_grc_loc = ',qflx_glcice_dyn_water_flux_locgrc(indexg)*dtime + + write(iulog,*)'clm model is stopping' + call endrun(decomp_index=indexg, clmlevel=nameg, msg=errmsg(sourcefile, __LINE__)) + end if + + end if + + ! Snow balance check at the grid cell level. + ! Beginning snow balance variable h2osno_old is calculated once + ! for both the column-level and grid cell-level balance checks. call waterstate_inst%CalculateTotalH2osno(bounds, num_allc, filter_allc, & caller = 'BalanceCheck', & @@ -574,6 +831,7 @@ subroutine BalanceCheck( bounds, & end if + ! Energy balance checks do p = bounds%begp, bounds%endp diff --git a/src/biogeophys/WaterBalanceType.F90 b/src/biogeophys/WaterBalanceType.F90 index d8e6f3f8ba..4d747acebe 100644 --- a/src/biogeophys/WaterBalanceType.F90 +++ b/src/biogeophys/WaterBalanceType.F90 @@ -36,10 +36,13 @@ module WaterBalanceType ! Balance Checks - real(r8), pointer :: begwb_col (:) ! water mass begining of the time step - real(r8), pointer :: endwb_col (:) ! water mass end of the time step + real(r8), pointer :: begwb_grc (:) ! grid cell-level water mass begining of the time step + real(r8), pointer :: endwb_grc (:) ! grid cell-level water mass end of the time step + real(r8), pointer :: begwb_col (:) ! column-level water mass begining of the time step + real(r8), pointer :: endwb_col (:) ! column-level water mass end of the time step real(r8), pointer :: errh2o_patch (:) ! water conservation error (mm H2O) - real(r8), pointer :: errh2o_col (:) ! water conservation error (mm H2O) + real(r8), pointer :: errh2o_col (:) ! column-level water conservation error (mm H2O) + real(r8), pointer :: errh2o_grc (:) ! grid cell-level water conservation error (mm H2O) real(r8), pointer :: errh2osno_col (:) ! snow water conservation error(mm H2O) contains @@ -112,6 +115,12 @@ subroutine InitAllocate(this, bounds, tracer_vars) container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) + call AllocateVar1d(var = this%begwb_grc, name = 'begwb_grc', & + container = tracer_vars, & + bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL) + call AllocateVar1d(var = this%endwb_grc, name = 'endwb_grc', & + container = tracer_vars, & + bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL) call AllocateVar1d(var = this%begwb_col, name = 'begwb_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) @@ -124,6 +133,9 @@ subroutine InitAllocate(this, bounds, tracer_vars) call AllocateVar1d(var = this%errh2o_col, name = 'errh2o_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) + call AllocateVar1d(var = this%errh2o_grc, name = 'errh2o_grc', & + container = tracer_vars, & + bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL) call AllocateVar1d(var = this%errh2osno_col, name = 'errh2osno_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index f06c49bc12..c332f9a662 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -26,7 +26,7 @@ module clm_driver use abortutils , only : endrun ! use dynSubgridDriverMod , only : dynSubgrid_driver, dynSubgrid_wrapup_weight_changes - use BalanceCheckMod , only : BeginWaterBalance, BalanceCheck + use BalanceCheckMod , only : BeginWaterGridcellBalance, BeginWaterColumnBalance, BalanceCheck ! use BiogeophysPreFluxCalcsMod , only : BiogeophysPreFluxCalcs use SurfaceHumidityMod , only : CalculateSurfaceHumidity @@ -326,6 +326,13 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro end if call t_stopf('begcnbal_grc') + call t_startf('begwbal') + call BeginWaterGridcellBalance(bounds_clump, & + filter(nc)%num_nolakec, filter(nc)%nolakec, & + filter(nc)%num_lakec, filter(nc)%lakec, & + water_inst, soilhydrology_inst, & + use_aquifer_layer = use_aquifer_layer()) + call t_stopf('begwbal') end do !$OMP END PARALLEL DO @@ -361,9 +368,9 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro ! For water: Currently, I believe this needs to be done after weights are updated for ! prescribed transient patches or CNDV, because column-level water is not generally ! conserved when weights change (instead the difference is put in the grid cell-level - ! terms, qflx_liq_dynbal, etc.). In the future, we may want to change the balance - ! checks to ensure that the grid cell-level water is conserved, considering - ! qflx_liq_dynbal; in this case, the call to BeginWaterBalance should be moved to + ! terms, qflx_liq_dynbal, etc.). Grid cell-level balance + ! checks ensure that the grid cell-level water is conserved by considering + ! qflx_liq_dynbal and calling BeginWaterGridcellBalance ! before the weight updates. ! ! For carbon & nitrogen: This needs to be done after dynSubgrid_driver, because the @@ -381,7 +388,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call t_stopf('prescribed_sm') endif call t_startf('begwbal') - call BeginWaterBalance(bounds_clump, & + call BeginWaterColumnBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & water_inst, soilhydrology_inst, & @@ -1103,6 +1110,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, & water_inst%waterbalancebulk_inst, water_inst%wateratm2lndbulk_inst, & surfalb_inst, energyflux_inst, canopystate_inst) +! water_inst%waterlnd2atmbulk_inst, surfalb_inst, energyflux_inst, canopystate_inst) call t_stopf('balchk') ! ============================================================================ diff --git a/src/main/lnd2atmMod.F90 b/src/main/lnd2atmMod.F90 index 9a18428264..6c64043d8c 100644 --- a/src/main/lnd2atmMod.F90 +++ b/src/main/lnd2atmMod.F90 @@ -405,10 +405,10 @@ subroutine lnd2atm(bounds, & call c2g( bounds, & water_inst%waterbalancebulk_inst%endwb_col(bounds%begc:bounds%endc), & - water_inst%waterdiagnosticbulk_inst%tws_grc (bounds%begg:bounds%endg), & + water_inst%waterbalancebulk_inst%endwb_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg - water_inst%waterdiagnosticbulk_inst%tws_grc(g) = water_inst%waterdiagnosticbulk_inst%tws_grc(g) + water_inst%wateratm2lndbulk_inst%volr_grc(g) / grc%area(g) * 1.e-3_r8 + water_inst%waterdiagnosticbulk_inst%tws_grc(g) = water_inst%waterbalancebulk_inst%endwb_grc(g) + water_inst%wateratm2lndbulk_inst%volr_grc(g) / grc%area(g) * 1.e-3_r8 enddo end subroutine lnd2atm From 5a0c758d672fcae5ef383617663650d76f2348bd Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Tue, 8 Dec 2020 16:17:06 -0700 Subject: [PATCH 02/19] Deleted a blank line; committing this in prep for update to latest tag --- src/biogeophys/BalanceCheckMod.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index bc73b68265..206af92946 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -831,7 +831,6 @@ subroutine BalanceCheck( bounds, & end if - ! Energy balance checks do p = bounds%begp, bounds%endp From ff59a5fe62f57b9b7b795256254311bcda1aafe9 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Tue, 8 Dec 2020 17:51:35 -0700 Subject: [PATCH 03/19] Mods required in new code from conflicts not caught by git --- src/biogeophys/BalanceCheckMod.F90 | 12 ++++++++---- src/main/clm_driver.F90 | 8 ++++---- 2 files changed, 12 insertions(+), 8 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index ea3e4360ab..b9390c57b3 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -126,7 +126,7 @@ end function GetBalanceCheckSkipSteps !----------------------------------------------------------------------- subroutine BeginWaterGridcellBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, soilhydrology_inst, & + water_inst, soilhydrology_inst, lakestate_inst, & use_aquifer_layer) ! ! !DESCRIPTION: @@ -140,6 +140,7 @@ subroutine BeginWaterGridcellBalance(bounds, & integer , intent(in) :: num_lakec ! number of column lake points in column filter integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst + type(lakestate_type) , intent(in) :: lakestate_inst type(soilhydrology_type), intent(in) :: soilhydrology_inst logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! @@ -154,6 +155,7 @@ subroutine BeginWaterGridcellBalance(bounds, & num_nolakec, filter_nolakec, & num_lakec, filter_lakec, & soilhydrology_inst, & + lakestate_inst, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & water_inst%bulk_and_tracers(i)%waterbalance_inst, & @@ -206,7 +208,8 @@ end subroutine BeginWaterColumnBalance !----------------------------------------------------------------------- subroutine BeginWaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - soilhydrology_inst, waterstate_inst, waterdiagnostic_inst, waterbalance_inst, & + soilhydrology_inst, lakestate_inst, waterstate_inst, & + waterdiagnostic_inst, waterbalance_inst, & use_aquifer_layer) ! ! !DESCRIPTION: @@ -223,6 +226,7 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & integer , intent(in) :: num_lakec ! number of column lake points in column filter integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(soilhydrology_type) , intent(in) :: soilhydrology_inst + type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst @@ -267,8 +271,8 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & water_mass = begwb_col(begc:endc)) call ComputeWaterMassLake(bounds, num_lakec, filter_lakec, & - waterstate_inst, & - subtract_dynbal_baselines = .false., & + waterstate_inst, lakestate_inst, & + add_lake_water_and_subtract_dynbal_baselines = .false., & water_mass = begwb_col(begc:endc)) call c2g(bounds, begwb_col(begc:endc), begwb_grc(begg:endg), & diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index 3b619b1e45..943106bbd2 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -327,10 +327,10 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call t_stopf('begcnbal_grc') call t_startf('begwbal') - call BeginWaterGridcellBalance(bounds_clump, & - filter(nc)%num_nolakec, filter(nc)%nolakec, & - filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, soilhydrology_inst, & + call BeginWaterGridcellBalance(bounds_clump, & + filter(nc)%num_nolakec, filter(nc)%nolakec, & + filter(nc)%num_lakec, filter(nc)%lakec, & + water_inst, soilhydrology_inst, lakestate_inst, & use_aquifer_layer = use_aquifer_layer()) call t_stopf('begwbal') end do From 78a979fd989f1c82875e995b030ddb6ecf7f5465 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Thu, 17 Dec 2020 14:46:50 -0700 Subject: [PATCH 04/19] Move call BalanceCheck to after call lnd2atm; use existing _grc vars --- src/biogeophys/BalanceCheckMod.F90 | 139 +++++++++------------------- src/biogeophys/Waterlnd2atmType.F90 | 7 +- src/main/clm_driver.F90 | 32 ++++--- src/main/lnd2atmMod.F90 | 5 +- 4 files changed, 70 insertions(+), 113 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index b9390c57b3..fb73adfaaf 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -24,7 +24,7 @@ module BalanceCheckMod use WaterDiagnosticBulkType, only : waterdiagnosticbulk_type use WaterDiagnosticType, only : waterdiagnostic_type use Wateratm2lndType , only : wateratm2lnd_type -! use Waterlnd2atmType , only : waterlnd2atm_type ! slevis: place holder + use Waterlnd2atmType , only : waterlnd2atm_type use WaterBalanceType , only : waterbalance_type use WaterFluxType , only : waterflux_type use WaterType , only : water_type @@ -356,8 +356,7 @@ subroutine BalanceCheck( bounds, & num_allc, filter_allc, & atm2lnd_inst, solarabs_inst, waterflux_inst, waterstate_inst, & waterdiagnosticbulk_inst, waterbalance_inst, wateratm2lnd_inst, & - surfalb_inst, energyflux_inst, canopystate_inst) -! waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst) + waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst) ! ! !DESCRIPTION: ! This subroutine accumulates the numerical truncation errors of the water @@ -390,7 +389,7 @@ subroutine BalanceCheck( bounds, & class(waterstate_type), intent(in) :: waterstate_inst type(waterdiagnosticbulk_type), intent(in) :: waterdiagnosticbulk_inst class(waterbalance_type), intent(inout) :: waterbalance_inst -! class(waterlnd2atm_type), intent(in) :: waterlnd2atm_inst + class(waterlnd2atm_type), intent(in) :: waterlnd2atm_inst class(wateratm2lnd_type) , intent(in) :: wateratm2lnd_inst type(surfalb_type) , intent(in) :: surfalb_inst type(energyflux_type) , intent(inout) :: energyflux_inst @@ -405,18 +404,9 @@ subroutine BalanceCheck( bounds, & real(r8) :: forc_rain_col(bounds%begc:bounds%endc) ! column level rain rate [mm/s] real(r8) :: forc_snow_col(bounds%begc:bounds%endc) ! column level snow rate [mm/s] real(r8) :: h2osno_total(bounds%begc:bounds%endc) ! total snow water [mm H2O] - real(r8) :: endwb_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_irrig_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_glcice_dyn_water_flux_locgrc(bounds%begg:bounds%endg) ! water flux needed for balance check due to glc_dyn_runoff_routing [mm H2O/s] (positive means addition of water to the system) - real(r8) :: qflx_evap_tot_locgrc(bounds%begg:bounds%endg) ! grid cell level total evapotranspiration [mm/s] - real(r8) :: qflx_surf_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_drain_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_drain_perched_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_qrgwl_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_ice_runoff_snwcp_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_ice_runoff_xs_locgrc(bounds%begg:bounds%endg) ! slevis: using local for now - real(r8) :: qflx_snwcp_discarded_liq_locgrc(bounds%begg:bounds%endg) ! excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] - real(r8) :: qflx_snwcp_discarded_ice_locgrc(bounds%begg:bounds%endg) ! excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] + real(r8) :: qflx_glcice_dyn_water_flux_grc(bounds%begg:bounds%endg) ! grid cell-level water flux needed for balance check due to glc_dyn_runoff_routing [mm H2O/s] (positive means addition of water to the system) + real(r8) :: qflx_snwcp_discarded_liq_grc(bounds%begg:bounds%endg) ! grid cell-level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] + real(r8) :: qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg) ! grid cell-level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] real(r8) :: errh2o_max_val ! Maximum value of error in water conservation error over all columns [mm H2O] real(r8) :: errh2osno_max_val ! Maximum value of error in h2osno conservation error over all columns [kg m-2] @@ -446,7 +436,7 @@ subroutine BalanceCheck( bounds, & frac_sno => waterdiagnosticbulk_inst%frac_sno_col , & ! Input: [real(r8) (:) ] fraction of ground covered by snow (0 to 1) snow_depth => waterdiagnosticbulk_inst%snow_depth_col , & ! Input: [real(r8) (:) ] snow height (m) begwb_grc => waterbalance_inst%begwb_grc , & ! Input: [real(r8) (:) ] grid cell-level water mass begining of the time step -! endwb_grc => waterbalance_inst%endwb_grc , & ! Output: [real(r8) (:) ] grid cell-level water mass end of the time step + endwb_grc => waterbalance_inst%endwb_grc , & ! Output: [real(r8) (:) ] grid cell-level water mass end of the time step errh2o_grc => waterbalance_inst%errh2o_grc , & ! Output: [real(r8) (:) ] grid cell-level water conservation error (mm H2O) begwb_col => waterbalance_inst%begwb_col , & ! Input: [real(r8) (:) ] column-level water mass begining of the time step endwb_col => waterbalance_inst%endwb_col , & ! Output: [real(r8) (:) ] column-level water mass end of the time step @@ -461,7 +451,7 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_liq => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:) ] excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` qflx_snwcp_discarded_ice => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:) ] excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` qflx_evap_tot => waterflux_inst%qflx_evap_tot_col , & ! Input: [real(r8) (:) ] qflx_evap_soi + qflx_evap_can + qflx_tran_veg -! qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc, & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg + qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc, & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_soliddew_to_top_layer => waterflux_inst%qflx_soliddew_to_top_layer_col , & ! Input: [real(r8) (:) ] rate of solid water deposited on top soil or snow layer (frost) (mm H2O /s) [+] qflx_solidevap_from_top_layer => waterflux_inst%qflx_solidevap_from_top_layer_col, & ! Input: [real(r8) (:) ] rate of ice evaporated from top soil or snow layer (sublimation) (mm H2O /s) [+] qflx_liqevap_from_top_layer => waterflux_inst%qflx_liqevap_from_top_layer_col , & ! Input: [real(r8) (:) ] rate of liquid water evaporated from top soil or snow layer (mm H2O/s) [+] @@ -470,24 +460,23 @@ subroutine BalanceCheck( bounds, & qflx_snow_h2osfc => waterflux_inst%qflx_snow_h2osfc_col , & ! Input: [real(r8) (:) ] snow falling on surface water (mm/s) qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice qflx_drain_perched => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) -! qflx_rofliq_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) + qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) qflx_floodc => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] total runoff due to flooding qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack ! qflx_liq_dynbal_grc => waterflux_inst%qflx_liq_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder ! qflx_ice_dynbal_grc => waterflux_inst%qflx_ice_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder -! qflx_runoff_col => waterflux_inst%qflx_runoff_col , & ! total runoff (mm H2O / s) slevis: place holder qflx_surf => waterflux_inst%qflx_surf_col , & ! Input: [real(r8) (:) ] surface runoff (mm H2O /s) -! qflx_rofliq_qsur_grc => waterlnd2atm_inst%qflx_rofliq_qsur_grc , & ! Input: [real(r8) (:) ] grid cell-level surface runoff (mm H20 /s) + qflx_surf_grc => waterlnd2atm_inst%qflx_rofliq_qsur_grc , & ! Input: [real(r8) (:) ] grid cell-level surface runoff (mm H20 /s) qflx_qrgwl => waterflux_inst%qflx_qrgwl_col , & ! Input: [real(r8) (:) ] qflx_surf at glaciers, wetlands, lakes -! qflx_rofliq_qgwl_grc => waterlnd2atm_inst%qflx_rofliq_qgwl_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_surf at glaciers, wetlands, lakes + qflx_qrgwl_grc => waterlnd2atm_inst%qflx_rofliq_qgwl_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_surf at glaciers, wetlands, lakes qflx_drain => waterflux_inst%qflx_drain_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) -! qflx_rofliq_qsub_grc => waterlnd2atm_inst%qflx_rofliq_qsub_grc , & ! Input: [real(r8) (:) ] grid cell-level drainage (mm H20 /s) - qflx_ice_runoff_snwcp => waterflux_inst%qflx_ice_runoff_snwcp_col, & ! Input: [real(r8) (:) ] solid runoff from snow capping (mm H2O /s) - qflx_ice_runoff_xs => waterflux_inst%qflx_ice_runoff_xs_col , & ! Input: [real(r8) (:) ] solid runoff from excess ice in soil (mm H2O /s) + qflx_drain_grc => waterlnd2atm_inst%qflx_rofliq_qsub_grc , & ! Input: [real(r8) (:) ] grid cell-level drainage (mm H20 /s) + qflx_ice_runoff => waterlnd2atm_inst%qflx_ice_runoff_col , & ! Input: [real(r8) (:) ] column level solid runoff from snow capping and from excess ice in soil (mm H2O /s) + qflx_ice_runoff_grc => waterlnd2atm_inst%qflx_rofice_grc , & ! Input: [real(r8) (:) ] grid cell-level solid runoff from snow capping and from excess ice in soil (mm H2O /s) qflx_sl_top_soil => waterflux_inst%qflx_sl_top_soil_col , & ! Input: [real(r8) (:) ] liquid water + ice from layer above soil to top soil layer or sent to qflx_qrgwl (mm H2O/s) qflx_sfc_irrig => waterflux_inst%qflx_sfc_irrig_col , & ! Input: [real(r8) (:) ] irrigation flux (mm H2O /s) -! qirrig_grc => waterlnd2atm_inst%qirrig_grc , & ! Input: [real(r8) (:) ] grid cell-level irrigation flux (mm H20 /s) + qirrig_grc => waterlnd2atm_inst%qirrig_grc , & ! Input: [real(r8) (:) ] grid cell-level irrigation flux (mm H20 /s) qflx_glcice_dyn_water_flux => waterflux_inst%qflx_glcice_dyn_water_flux_col, & ! Input: [real(r8) (:)] water flux needed for balance check due to glc_dyn_runoff_routing (mm H2O/s) (positive means addition of water to the system) eflx_lwrad_out => energyflux_inst%eflx_lwrad_out_patch , & ! Input: [real(r8) (:) ] emitted infrared (longwave) radiation (W/m**2) @@ -566,8 +555,7 @@ subroutine BalanceCheck( bounds, & - qflx_qrgwl(c) & - qflx_drain(c) & - qflx_drain_perched(c) & - - qflx_ice_runoff_snwcp(c) & - - qflx_ice_runoff_xs(c) & + - qflx_ice_runoff(c) & - qflx_snwcp_discarded_liq(c) & - qflx_snwcp_discarded_ice(c)) * dtime @@ -605,8 +593,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_qrgwl = ',qflx_qrgwl(indexc)*dtime write(iulog,*)'qflx_drain = ',qflx_drain(indexc)*dtime - write(iulog,*)'qflx_ice_runoff_snwcp = ',qflx_ice_runoff_snwcp(indexc)*dtime - write(iulog,*)'qflx_ice_runoff_xs = ',qflx_ice_runoff_xs(indexc)*dtime + write(iulog,*)'qflx_ice_runoff = ',qflx_ice_runoff(indexc)*dtime write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice(indexc)*dtime write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq(indexc)*dtime @@ -630,71 +617,34 @@ subroutine BalanceCheck( bounds, & ! Water balance check at the grid cell level - call c2g(bounds, & - endwb_col(bounds%begc:bounds%endc), & - endwb_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type='urbanf', l2g_scale_type='unity') - call c2g( bounds, & - qflx_sfc_irrig(bounds%begc:bounds%endc), & - qflx_irrig_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & qflx_glcice_dyn_water_flux(bounds%begc:bounds%endc), & - qflx_glcice_dyn_water_flux_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_evap_tot(bounds%begc:bounds%endc), & - qflx_evap_tot_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_surf(bounds%begc:bounds%endc), & - qflx_surf_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_qrgwl(bounds%begc:bounds%endc), & - qflx_qrgwl_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_drain(bounds%begc:bounds%endc), & - qflx_drain_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_drain_perched(bounds%begc:bounds%endc), & - qflx_drain_perched_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_ice_runoff_snwcp(bounds%begc:bounds%endc), & - qflx_ice_runoff_snwcp_locgrc(bounds%begg:bounds%endg), & - c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call c2g( bounds, & - qflx_ice_runoff_xs(bounds%begc:bounds%endc), & - qflx_ice_runoff_xs_locgrc(bounds%begg:bounds%endg), & + qflx_glcice_dyn_water_flux_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & qflx_snwcp_discarded_liq(bounds%begc:bounds%endc), & - qflx_snwcp_discarded_liq_locgrc(bounds%begg:bounds%endg), & + qflx_snwcp_discarded_liq_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & qflx_snwcp_discarded_ice(bounds%begc:bounds%endc), & - qflx_snwcp_discarded_ice_locgrc(bounds%begg:bounds%endg), & + qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg - errh2o_grc(g) = endwb_locgrc(g) - begwb_grc(g) & + errh2o_grc(g) = endwb_grc(g) - begwb_grc(g) & - (forc_rain_grc(g) & + forc_snow_grc(g) & + qflx_flood_grc(g) & - + qflx_irrig_locgrc(g) & - + qflx_glcice_dyn_water_flux_locgrc(g) & - - qflx_evap_tot_locgrc(g) & - - qflx_surf_locgrc(g) & - - qflx_qrgwl_locgrc(g) & - - qflx_drain_locgrc(g) & - - qflx_drain_perched_locgrc(g) & - - qflx_ice_runoff_snwcp_locgrc(g) & - - qflx_ice_runoff_xs_locgrc(g) & - - qflx_snwcp_discarded_liq_locgrc(g) & - - qflx_snwcp_discarded_ice_locgrc(g)) * dtime + + qirrig_grc(g) & + + qflx_glcice_dyn_water_flux_grc(g) & + - qflx_evap_tot_grc(g) & + - qflx_surf_grc(g) & + - qflx_qrgwl_grc(g) & + - qflx_drain_grc(g) & + - qflx_drain_perched_grc(g) & + - qflx_ice_runoff_grc(g) & + - qflx_snwcp_discarded_liq_grc(g) & + - qflx_snwcp_discarded_ice_grc(g)) * dtime end do errh2o_max_val = maxval(abs(errh2o_grc(bounds%begg:bounds%endg))) @@ -713,23 +663,22 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errh2o_grc = ',errh2o_grc(indexg) write(iulog,*)'forc_rain = ',forc_rain_grc(indexg)*dtime write(iulog,*)'forc_snow = ',forc_snow_grc(indexg)*dtime - write(iulog,*)'endwb_loc = ',endwb_locgrc(indexg) + write(iulog,*)'endwb_grc = ',endwb_grc(indexg) write(iulog,*)'begwb_grc = ',begwb_grc(indexg) - write(iulog,*)'qflx_evap_tot_loc = ',qflx_evap_tot_locgrc(indexg)*dtime - write(iulog,*)'qflx_irrig_loc = ',qflx_irrig_locgrc(indexg)*dtime - write(iulog,*)'qflx_surf_loc = ',qflx_surf_locgrc(indexg)*dtime - write(iulog,*)'qflx_qrgwl_loc = ',qflx_qrgwl_locgrc(indexg)*dtime - write(iulog,*)'qflx_drain_loc = ',qflx_drain_locgrc(indexg)*dtime - write(iulog,*)'qflx_ice_runoff_snwcp_loc = ',qflx_ice_runoff_snwcp_locgrc(indexg)*dtime - write(iulog,*)'qflx_ice_runoff_xs_loc = ',qflx_ice_runoff_xs_locgrc(indexg)*dtime - write(iulog,*)'qflx_snwcp_discarded_ice_loc = ',qflx_snwcp_discarded_ice_locgrc(indexg)*dtime - write(iulog,*)'qflx_snwcp_discarded_liq_loc = ',qflx_snwcp_discarded_liq_locgrc(indexg)*dtime - write(iulog,*)'deltawb = ',endwb_locgrc(indexg)-begwb_grc(indexg) - write(iulog,*)'deltawb/dtime = ',(endwb_locgrc(indexg)-begwb_grc(indexg))/dtime - write(iulog,*)'qflx_drain_perched_loc = ',qflx_drain_perched_locgrc(indexg)*dtime + write(iulog,*)'qflx_evap_tot = ',qflx_evap_tot_grc(indexg)*dtime + write(iulog,*)'qirrig = ',qirrig_grc(indexg)*dtime + write(iulog,*)'qflx_surf = ',qflx_surf_grc(indexg)*dtime + write(iulog,*)'qflx_qrgwl = ',qflx_qrgwl_grc(indexg)*dtime + write(iulog,*)'qflx_drain = ',qflx_drain_grc(indexg)*dtime + write(iulog,*)'qflx_ice_runoff = ',qflx_ice_runoff_grc(indexg)*dtime + write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice_grc(indexg)*dtime + write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq_grc(indexg)*dtime + write(iulog,*)'deltawb = ',endwb_grc(indexg)-begwb_grc(indexg) + write(iulog,*)'deltawb/dtime = ',(endwb_grc(indexg)-begwb_grc(indexg))/dtime + write(iulog,*)'qflx_drain_perched = ',qflx_drain_perched_grc(indexg)*dtime write(iulog,*)'qflx_flood = ',qflx_flood_grc(indexg)*dtime - write(iulog,*)'qflx_glcice_dyn_water_flux_grc_loc = ',qflx_glcice_dyn_water_flux_locgrc(indexg)*dtime + write(iulog,*)'qflx_glcice_dyn_water_flux_grc = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime write(iulog,*)'clm model is stopping' call endrun(decomp_index=indexg, clmlevel=nameg, msg=errmsg(sourcefile, __LINE__)) diff --git a/src/biogeophys/Waterlnd2atmType.F90 b/src/biogeophys/Waterlnd2atmType.F90 index ed6e9ca0dd..fb59d1c83c 100644 --- a/src/biogeophys/Waterlnd2atmType.F90 +++ b/src/biogeophys/Waterlnd2atmType.F90 @@ -32,7 +32,8 @@ module Waterlnd2atmType real(r8), pointer :: qflx_rofliq_qsub_grc (:) ! rof liq -- subsurface runoff component real(r8), pointer :: qflx_rofliq_qgwl_grc (:) ! rof liq -- glacier, wetland and lakes water balance residual component real(r8), pointer :: qflx_rofliq_drain_perched_grc (:) ! rof liq -- perched water table runoff component - real(r8), pointer :: qflx_rofice_grc (:) ! rof ice forcing + real(r8), pointer :: qflx_ice_runoff_col(:) ! rof ice forcing, col level + real(r8), pointer :: qflx_rofice_grc (:) ! rof ice forcing, grc level real(r8), pointer :: qflx_liq_from_ice_col(:) ! liquid runoff from converted ice runoff real(r8), pointer :: qirrig_grc (:) ! irrigation flux @@ -119,6 +120,10 @@ subroutine InitAllocate(this, bounds, tracer_vars) container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL, & ival=ival) + call AllocateVar1d(var = this%qflx_ice_runoff_col, name = 'qflx_ice_runoff_col', & + container = tracer_vars, & + bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN, & + ival=ival) call AllocateVar1d(var = this%qflx_rofice_grc, name = 'qflx_rofice_grc', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL, & diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index 943106bbd2..eb94d5de83 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -1105,20 +1105,6 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro filter(nc)%num_soilp, filter(nc)%soilp, & filter(nc)%num_allc, filter(nc)%allc) - ! ============================================================================ - ! Check the energy and water balance - ! ============================================================================ - - call t_startf('balchk') - call BalanceCheck(bounds_clump, & - filter(nc)%num_allc, filter(nc)%allc, & - atm2lnd_inst, solarabs_inst, water_inst%waterfluxbulk_inst, & - water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, & - water_inst%waterbalancebulk_inst, water_inst%wateratm2lndbulk_inst, & - surfalb_inst, energyflux_inst, canopystate_inst) -! water_inst%waterlnd2atmbulk_inst, surfalb_inst, energyflux_inst, canopystate_inst) - call t_stopf('balchk') - ! ============================================================================ ! Check the carbon and nitrogen balance ! ============================================================================ @@ -1272,6 +1258,24 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro !$OMP END PARALLEL DO call t_stopf('lnd2glc') + ! ========================================================================== + ! Check the energy and water balance + ! ========================================================================== + + call t_startf('balchk') + !$OMP PARALLEL DO PRIVATE (nc, bounds_clump) + do nc = 1,nclumps + call get_clump_bounds(nc, bounds_clump) + call BalanceCheck(bounds_clump, & + filter(nc)%num_allc, filter(nc)%allc, & + atm2lnd_inst, solarabs_inst, water_inst%waterfluxbulk_inst, & + water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, & + water_inst%waterbalancebulk_inst, water_inst%wateratm2lndbulk_inst, & + water_inst%waterlnd2atmbulk_inst, surfalb_inst, energyflux_inst, canopystate_inst) + end do + !$OMP END PARALLEL DO + call t_stopf('balchk') + ! ============================================================================ ! Write global average diagnostics to standard output ! ============================================================================ diff --git a/src/main/lnd2atmMod.F90 b/src/main/lnd2atmMod.F90 index e7052b6020..d3eb22d610 100644 --- a/src/main/lnd2atmMod.F90 +++ b/src/main/lnd2atmMod.F90 @@ -180,7 +180,6 @@ subroutine lnd2atm(bounds, & ! ! !LOCAL VARIABLES: integer :: c, g ! indices - real(r8) :: qflx_ice_runoff_col(bounds%begc:bounds%endc) ! total column-level ice runoff real(r8) :: eflx_sh_ice_to_liq_grc(bounds%begg:bounds%endg) ! sensible heat flux generated from the ice to liquid conversion, averaged to gridcell real(r8), parameter :: amC = 12.0_r8 ! Atomic mass number for Carbon real(r8), parameter :: amO = 16.0_r8 ! Atomic mass number for Oxygen @@ -193,7 +192,7 @@ subroutine lnd2atm(bounds, & call handle_ice_runoff(bounds, water_inst%waterfluxbulk_inst, glc_behavior, & melt_non_icesheet_ice_runoff = lnd2atm_inst%params%melt_non_icesheet_ice_runoff, & - qflx_ice_runoff_col = qflx_ice_runoff_col(bounds%begc:bounds%endc), & + qflx_ice_runoff_col = water_inst%waterlnd2atmbulk_inst%qflx_ice_runoff_col(bounds%begc:bounds%endc), & qflx_liq_from_ice_col = water_inst%waterlnd2atmbulk_inst%qflx_liq_from_ice_col(bounds%begc:bounds%endc), & eflx_sh_ice_to_liq_col = lnd2atm_inst%eflx_sh_ice_to_liq_col(bounds%begc:bounds%endc)) @@ -391,7 +390,7 @@ subroutine lnd2atm(bounds, & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & - qflx_ice_runoff_col(bounds%begc:bounds%endc), & + water_inst%waterlnd2atmbulk_inst%qflx_ice_runoff_col(bounds%begc:bounds%endc), & water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg From 4ab9f3d54e7151b80d3b5b08f4e5ed5d857ed746 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Thu, 17 Dec 2020 19:21:53 -0700 Subject: [PATCH 05/19] Some clean-up before moving on to solving the transient h2o balance --- src/biogeophys/BalanceCheckMod.F90 | 108 ++++++++++++++--------------- 1 file changed, 54 insertions(+), 54 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index fb73adfaaf..91c6b9c0e9 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -377,7 +377,7 @@ subroutine BalanceCheck( bounds, & use clm_time_manager , only : get_step_size_real, get_nstep use clm_time_manager , only : get_nstep_since_startup_or_lastDA_restart_or_pause use CanopyStateType , only : canopystate_type - use subgridAveMod ! , only : c2g ? + use subgridAveMod , only : c2g ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds @@ -424,11 +424,10 @@ subroutine BalanceCheck( bounds, & associate( & forc_solad => atm2lnd_inst%forc_solad_grc , & ! Input: [real(r8) (:,:) ] direct beam radiation (vis=forc_sols , nir=forc_soll ) forc_solai => atm2lnd_inst%forc_solai_grc , & ! Input: [real(r8) (:,:) ] diffuse radiation (vis=forc_solsd, nir=forc_solld) - forc_rain_grc => wateratm2lnd_inst%forc_rain_not_downscaled_grc , & ! Input: [real(r8) (:)] grid cell-level rain rate [mm/s] - forc_snow_grc => wateratm2lnd_inst%forc_snow_not_downscaled_grc , & ! Input: [real(r8) (:)] grid cell-level snow rate [mm/s] - qflx_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] total grid cell-level runoff from river model - forc_rain => wateratm2lnd_inst%forc_rain_downscaled_col , & ! Input: [real(r8) (:) ] rain rate [mm/s] - forc_snow => wateratm2lnd_inst%forc_snow_downscaled_col , & ! Input: [real(r8) (:) ] snow rate [mm/s] + forc_rain => wateratm2lnd_inst%forc_rain_downscaled_col , & ! Input: [real(r8) (:) ] column level rain rate [mm/s] + forc_rain_grc => wateratm2lnd_inst%forc_rain_not_downscaled_grc, & ! Input: [real(r8) (:)] grid cell-level rain rate [mm/s] + forc_snow => wateratm2lnd_inst%forc_snow_downscaled_col , & ! Input: [real(r8) (:) ] column level snow rate [mm/s] + forc_snow_grc => wateratm2lnd_inst%forc_snow_not_downscaled_grc, & ! Input: [real(r8) (:)] grid cell-level snow rate [mm/s] forc_lwrad => atm2lnd_inst%forc_lwrad_downscaled_col , & ! Input: [real(r8) (:) ] downward infrared (longwave) radiation (W/m**2) h2osno_old => waterbalance_inst%h2osno_old_col , & ! Input: [real(r8) (:) ] snow water (mm H2O) at previous time step @@ -448,10 +447,10 @@ subroutine BalanceCheck( bounds, & qflx_snow_grnd_col => waterflux_inst%qflx_snow_grnd_col , & ! Input: [real(r8) (:) ] snow on ground after interception (mm H2O/s) [+] qflx_snwcp_liq => waterflux_inst%qflx_snwcp_liq_col , & ! Input: [real(r8) (:) ] excess liquid h2o due to snow capping (outgoing) (mm H2O /s) [+]` qflx_snwcp_ice => waterflux_inst%qflx_snwcp_ice_col , & ! Input: [real(r8) (:) ] excess solid h2o due to snow capping (outgoing) (mm H2O /s) [+]` - qflx_snwcp_discarded_liq => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:) ] excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` - qflx_snwcp_discarded_ice => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:) ] excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+]` - qflx_evap_tot => waterflux_inst%qflx_evap_tot_col , & ! Input: [real(r8) (:) ] qflx_evap_soi + qflx_evap_can + qflx_tran_veg - qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc, & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg + qflx_snwcp_discarded_liq_col => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:) ] column level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] + qflx_snwcp_discarded_ice_col => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:) ] column level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] + qflx_evap_tot_col => waterflux_inst%qflx_evap_tot_col , & ! Input: [real(r8) (:) ] column level qflx_evap_soi + qflx_evap_can + qflx_tran_veg + qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_soliddew_to_top_layer => waterflux_inst%qflx_soliddew_to_top_layer_col , & ! Input: [real(r8) (:) ] rate of solid water deposited on top soil or snow layer (frost) (mm H2O /s) [+] qflx_solidevap_from_top_layer => waterflux_inst%qflx_solidevap_from_top_layer_col, & ! Input: [real(r8) (:) ] rate of ice evaporated from top soil or snow layer (sublimation) (mm H2O /s) [+] qflx_liqevap_from_top_layer => waterflux_inst%qflx_liqevap_from_top_layer_col , & ! Input: [real(r8) (:) ] rate of liquid water evaporated from top soil or snow layer (mm H2O/s) [+] @@ -459,25 +458,26 @@ subroutine BalanceCheck( bounds, & qflx_prec_grnd => waterdiagnosticbulk_inst%qflx_prec_grnd_col, & ! Input: [real(r8) (:) ] water onto ground including canopy runoff [kg/(m2 s)] qflx_snow_h2osfc => waterflux_inst%qflx_snow_h2osfc_col , & ! Input: [real(r8) (:) ] snow falling on surface water (mm/s) qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice - qflx_drain_perched => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) - qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) - qflx_floodc => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] total runoff due to flooding + qflx_drain_perched_col => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) + qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) + qflx_flood_col => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] column level total runoff due to flooding + forc_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] grid cell-level total grid cell-level runoff from river model qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack ! qflx_liq_dynbal_grc => waterflux_inst%qflx_liq_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder ! qflx_ice_dynbal_grc => waterflux_inst%qflx_ice_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder - qflx_surf => waterflux_inst%qflx_surf_col , & ! Input: [real(r8) (:) ] surface runoff (mm H2O /s) + qflx_surf_col => waterflux_inst%qflx_surf_col , & ! Input: [real(r8) (:) ] column level surface runoff (mm H2O /s) qflx_surf_grc => waterlnd2atm_inst%qflx_rofliq_qsur_grc , & ! Input: [real(r8) (:) ] grid cell-level surface runoff (mm H20 /s) - qflx_qrgwl => waterflux_inst%qflx_qrgwl_col , & ! Input: [real(r8) (:) ] qflx_surf at glaciers, wetlands, lakes + qflx_qrgwl_col => waterflux_inst%qflx_qrgwl_col , & ! Input: [real(r8) (:) ] column level qflx_surf at glaciers, wetlands, lakes qflx_qrgwl_grc => waterlnd2atm_inst%qflx_rofliq_qgwl_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_surf at glaciers, wetlands, lakes - qflx_drain => waterflux_inst%qflx_drain_col , & ! Input: [real(r8) (:) ] sub-surface runoff (mm H2O /s) + qflx_drain_col => waterflux_inst%qflx_drain_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) qflx_drain_grc => waterlnd2atm_inst%qflx_rofliq_qsub_grc , & ! Input: [real(r8) (:) ] grid cell-level drainage (mm H20 /s) - qflx_ice_runoff => waterlnd2atm_inst%qflx_ice_runoff_col , & ! Input: [real(r8) (:) ] column level solid runoff from snow capping and from excess ice in soil (mm H2O /s) + qflx_ice_runoff_col => waterlnd2atm_inst%qflx_ice_runoff_col , & ! Input: [real(r8) (:) ] column level solid runoff from snow capping and from excess ice in soil (mm H2O /s) qflx_ice_runoff_grc => waterlnd2atm_inst%qflx_rofice_grc , & ! Input: [real(r8) (:) ] grid cell-level solid runoff from snow capping and from excess ice in soil (mm H2O /s) qflx_sl_top_soil => waterflux_inst%qflx_sl_top_soil_col , & ! Input: [real(r8) (:) ] liquid water + ice from layer above soil to top soil layer or sent to qflx_qrgwl (mm H2O/s) - qflx_sfc_irrig => waterflux_inst%qflx_sfc_irrig_col , & ! Input: [real(r8) (:) ] irrigation flux (mm H2O /s) - qirrig_grc => waterlnd2atm_inst%qirrig_grc , & ! Input: [real(r8) (:) ] grid cell-level irrigation flux (mm H20 /s) - qflx_glcice_dyn_water_flux => waterflux_inst%qflx_glcice_dyn_water_flux_col, & ! Input: [real(r8) (:)] water flux needed for balance check due to glc_dyn_runoff_routing (mm H2O/s) (positive means addition of water to the system) + qflx_sfc_irrig_col => waterflux_inst%qflx_sfc_irrig_col , & ! Input: [real(r8) (:) ] column level irrigation flux (mm H2O /s) + qflx_sfc_irrig_grc => waterlnd2atm_inst%qirrig_grc , & ! Input: [real(r8) (:) ] grid cell-level irrigation flux (mm H20 /s) + qflx_glcice_dyn_water_flux_col => waterflux_inst%qflx_glcice_dyn_water_flux_col, & ! Input: [real(r8) (:)] column level water flux needed for balance check due to glc_dyn_runoff_routing (mm H2O/s) (positive means addition of water to the system) eflx_lwrad_out => energyflux_inst%eflx_lwrad_out_patch , & ! Input: [real(r8) (:) ] emitted infrared (longwave) radiation (W/m**2) eflx_lwrad_net => energyflux_inst%eflx_lwrad_net_patch , & ! Input: [real(r8) (:) ] net infrared (longwave) rad (W/m**2) [+ = to atm] @@ -547,17 +547,17 @@ subroutine BalanceCheck( bounds, & errh2o_col(c) = endwb_col(c) - begwb_col(c) & - (forc_rain_col(c) & + forc_snow_col(c) & - + qflx_floodc(c) & - + qflx_sfc_irrig(c) & - + qflx_glcice_dyn_water_flux(c) & - - qflx_evap_tot(c) & - - qflx_surf(c) & - - qflx_qrgwl(c) & - - qflx_drain(c) & - - qflx_drain_perched(c) & - - qflx_ice_runoff(c) & - - qflx_snwcp_discarded_liq(c) & - - qflx_snwcp_discarded_ice(c)) * dtime + + qflx_flood_col(c) & + + qflx_sfc_irrig_col(c) & + + qflx_glcice_dyn_water_flux_col(c) & + - qflx_evap_tot_col(c) & + - qflx_surf_col(c) & + - qflx_qrgwl_col(c) & + - qflx_drain_col(c) & + - qflx_drain_perched_col(c) & + - qflx_ice_runoff_col(c) & + - qflx_snwcp_discarded_liq_col(c) & + - qflx_snwcp_discarded_ice_col(c)) * dtime else @@ -587,25 +587,25 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'endwb_col = ',endwb_col(indexc) write(iulog,*)'begwb_col = ',begwb_col(indexc) - write(iulog,*)'qflx_evap_tot = ',qflx_evap_tot(indexc)*dtime - write(iulog,*)'qflx_sfc_irrig = ',qflx_sfc_irrig(indexc)*dtime - write(iulog,*)'qflx_surf = ',qflx_surf(indexc)*dtime - write(iulog,*)'qflx_qrgwl = ',qflx_qrgwl(indexc)*dtime - write(iulog,*)'qflx_drain = ',qflx_drain(indexc)*dtime + write(iulog,*)'qflx_evap_tot = ',qflx_evap_tot_col(indexc)*dtime + write(iulog,*)'qflx_sfc_irrig = ',qflx_sfc_irrig_col(indexc)*dtime + write(iulog,*)'qflx_surf = ',qflx_surf_col(indexc)*dtime + write(iulog,*)'qflx_qrgwl = ',qflx_qrgwl_col(indexc)*dtime + write(iulog,*)'qflx_drain = ',qflx_drain_col(indexc)*dtime - write(iulog,*)'qflx_ice_runoff = ',qflx_ice_runoff(indexc)*dtime + write(iulog,*)'qflx_ice_runoff = ',qflx_ice_runoff_col(indexc)*dtime - write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice(indexc)*dtime - write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq(indexc)*dtime + write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice_col(indexc)*dtime + write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq_col(indexc)*dtime write(iulog,*)'deltawb = ',endwb_col(indexc)-begwb_col(indexc) write(iulog,*)'deltawb/dtime = ',(endwb_col(indexc)-begwb_col(indexc))/dtime if (.not.(col%itype(indexc) == icol_roof .or. & col%itype(indexc) == icol_road_imperv .or. & col%itype(indexc) == icol_road_perv)) then - write(iulog,*)'qflx_drain_perched = ',qflx_drain_perched(indexc)*dtime - write(iulog,*)'qflx_flood = ',qflx_floodc(indexc)*dtime - write(iulog,*)'qflx_glcice_dyn_water_flux = ', qflx_glcice_dyn_water_flux(indexc)*dtime + write(iulog,*)'qflx_drain_perched = ',qflx_drain_perched_col(indexc)*dtime + write(iulog,*)'qflx_flood = ',qflx_flood_col(indexc)*dtime + write(iulog,*)'qflx_glcice_dyn_water_flux = ', qflx_glcice_dyn_water_flux_col(indexc)*dtime end if write(iulog,*)'clm model is stopping' @@ -618,15 +618,15 @@ subroutine BalanceCheck( bounds, & ! Water balance check at the grid cell level call c2g( bounds, & - qflx_glcice_dyn_water_flux(bounds%begc:bounds%endc), & + qflx_glcice_dyn_water_flux_col(bounds%begc:bounds%endc), & qflx_glcice_dyn_water_flux_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & - qflx_snwcp_discarded_liq(bounds%begc:bounds%endc), & + qflx_snwcp_discarded_liq_col(bounds%begc:bounds%endc), & qflx_snwcp_discarded_liq_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) call c2g( bounds, & - qflx_snwcp_discarded_ice(bounds%begc:bounds%endc), & + qflx_snwcp_discarded_ice_col(bounds%begc:bounds%endc), & qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) @@ -634,8 +634,8 @@ subroutine BalanceCheck( bounds, & errh2o_grc(g) = endwb_grc(g) - begwb_grc(g) & - (forc_rain_grc(g) & + forc_snow_grc(g) & - + qflx_flood_grc(g) & - + qirrig_grc(g) & + + forc_flood_grc(g) & + + qflx_sfc_irrig_grc(g) & + qflx_glcice_dyn_water_flux_grc(g) & - qflx_evap_tot_grc(g) & - qflx_surf_grc(g) & @@ -667,7 +667,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'begwb_grc = ',begwb_grc(indexg) write(iulog,*)'qflx_evap_tot = ',qflx_evap_tot_grc(indexg)*dtime - write(iulog,*)'qirrig = ',qirrig_grc(indexg)*dtime + write(iulog,*)'qflx_sfc_irrig = ',qflx_sfc_irrig_grc(indexg)*dtime write(iulog,*)'qflx_surf = ',qflx_surf_grc(indexg)*dtime write(iulog,*)'qflx_qrgwl = ',qflx_qrgwl_grc(indexg)*dtime write(iulog,*)'qflx_drain = ',qflx_drain_grc(indexg)*dtime @@ -677,7 +677,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'deltawb = ',endwb_grc(indexg)-begwb_grc(indexg) write(iulog,*)'deltawb/dtime = ',(endwb_grc(indexg)-begwb_grc(indexg))/dtime write(iulog,*)'qflx_drain_perched = ',qflx_drain_perched_grc(indexg)*dtime - write(iulog,*)'qflx_flood = ',qflx_flood_grc(indexg)*dtime + write(iulog,*)'forc_flood = ',forc_flood_grc(indexg)*dtime write(iulog,*)'qflx_glcice_dyn_water_flux_grc = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime write(iulog,*)'clm model is stopping' @@ -709,7 +709,7 @@ subroutine BalanceCheck( bounds, & + qflx_liqdew_to_top_layer(c) snow_sinks(c) = qflx_solidevap_from_top_layer(c) + qflx_liqevap_from_top_layer(c) & + qflx_snow_drain(c) + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snwcp_discarded_ice(c) + qflx_snwcp_discarded_liq(c) & + + qflx_snwcp_discarded_ice_col(c) + qflx_snwcp_discarded_liq_col(c) & + qflx_sl_top_soil(c) if (lun%itype(l) == istdlak) then @@ -718,7 +718,7 @@ subroutine BalanceCheck( bounds, & + qflx_soliddew_to_top_layer(c) + qflx_liqdew_to_top_layer(c) ) snow_sinks(c) = frac_sno_eff(c) * (qflx_solidevap_from_top_layer(c) & + qflx_liqevap_from_top_layer(c) ) + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snwcp_discarded_ice(c) + qflx_snwcp_discarded_liq(c) & + + qflx_snwcp_discarded_ice_col(c) + qflx_snwcp_discarded_liq_col(c) & + qflx_snow_drain(c) + qflx_sl_top_soil(c) endif @@ -731,7 +731,7 @@ subroutine BalanceCheck( bounds, & + qflx_h2osfc_to_ice(c) snow_sinks(c) = frac_sno_eff(c) * (qflx_solidevap_from_top_layer(c) & + qflx_liqevap_from_top_layer(c)) + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snwcp_discarded_ice(c) + qflx_snwcp_discarded_liq(c) & + + qflx_snwcp_discarded_ice_col(c) + qflx_snwcp_discarded_liq_col(c) & + qflx_snow_drain(c) + qflx_sl_top_soil(c) endif @@ -780,8 +780,8 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_liqdew_to_top_layer = ',qflx_liqdew_to_top_layer(indexc)*dtime write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc)*dtime write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc)*dtime - write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice(indexc)*dtime - write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq(indexc)*dtime + write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice_col(indexc)*dtime + write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq_col(indexc)*dtime write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc)*dtime write(iulog,*)'clm model is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) From d267cebc28311fa6f19e85ab704fcc35b502d780 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Wed, 23 Dec 2020 10:13:08 -0700 Subject: [PATCH 06/19] Mods for successful grid cell-level H2O error checks in transient cases These tests PASS: ERP_D.f10_f10_musgs.IHistClm50Bgc.cheyenne_gnu.clm-decStart ERS_Ly3_Mmpi-serial.1x1_smallvilleIA.IHistClm50BgcCropQianRs.cheyenne_gnu.clm-cropMonthOutput --- src/biogeophys/BalanceCheckMod.F90 | 24 ++++++++++++++++++++++-- 1 file changed, 22 insertions(+), 2 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 91c6b9c0e9..6bad841b8f 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -159,6 +159,7 @@ subroutine BeginWaterGridcellBalance(bounds, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & water_inst%bulk_and_tracers(i)%waterbalance_inst, & + water_inst%bulk_and_tracers(i)%waterflux_inst, & use_aquifer_layer = use_aquifer_layer) end do @@ -209,7 +210,7 @@ end subroutine BeginWaterColumnBalance subroutine BeginWaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & soilhydrology_inst, lakestate_inst, waterstate_inst, & - waterdiagnostic_inst, waterbalance_inst, & + waterdiagnostic_inst, waterbalance_inst, waterflux_inst, & use_aquifer_layer) ! ! !DESCRIPTION: @@ -228,13 +229,16 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & type(soilhydrology_type) , intent(in) :: soilhydrology_inst type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst + class(waterflux_type) , intent(inout) :: waterflux_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: - integer :: c, j, fc ! indices + integer :: c, g, j, fc ! indices integer :: begc, endc, begg, endg ! bounds + real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at beginning of time step + real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc ice dynamic land cover change conversion runoff flux at beginning of time step !----------------------------------------------------------------------- associate( & @@ -278,6 +282,14 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & call c2g(bounds, begwb_col(begc:endc), begwb_grc(begg:endg), & c2l_scale_type='urbanf', l2g_scale_type='unity') + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_beg(bounds, qflx_liq_dynbal_left_to_dribble(begg:endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_beg(bounds, qflx_ice_dynbal_left_to_dribble(begg:endg)) + + do g = begg, endg + begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & + - qflx_ice_dynbal_left_to_dribble(g) + end do + end associate end subroutine BeginWaterGridcellBalanceSingle @@ -407,6 +419,8 @@ subroutine BalanceCheck( bounds, & real(r8) :: qflx_glcice_dyn_water_flux_grc(bounds%begg:bounds%endg) ! grid cell-level water flux needed for balance check due to glc_dyn_runoff_routing [mm H2O/s] (positive means addition of water to the system) real(r8) :: qflx_snwcp_discarded_liq_grc(bounds%begg:bounds%endg) ! grid cell-level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] real(r8) :: qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg) ! grid cell-level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] + real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step + real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step real(r8) :: errh2o_max_val ! Maximum value of error in water conservation error over all columns [mm H2O] real(r8) :: errh2osno_max_val ! Maximum value of error in h2osno conservation error over all columns [kg m-2] @@ -630,7 +644,13 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end(bounds, qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_end(bounds, qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg)) + do g = bounds%begg, bounds%endg + endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & + - qflx_ice_dynbal_left_to_dribble(g) + errh2o_grc(g) = endwb_grc(g) - begwb_grc(g) & - (forc_rain_grc(g) & + forc_snow_grc(g) & From 32a793082ee494232457bc35cf5207fa32953e03 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Wed, 23 Dec 2020 11:23:14 -0700 Subject: [PATCH 07/19] Minor clean-up --- src/biogeophys/BalanceCheckMod.F90 | 22 ++++++++++------------ 1 file changed, 10 insertions(+), 12 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 6bad841b8f..ad2cbd31dd 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -439,9 +439,9 @@ subroutine BalanceCheck( bounds, & forc_solad => atm2lnd_inst%forc_solad_grc , & ! Input: [real(r8) (:,:) ] direct beam radiation (vis=forc_sols , nir=forc_soll ) forc_solai => atm2lnd_inst%forc_solai_grc , & ! Input: [real(r8) (:,:) ] diffuse radiation (vis=forc_solsd, nir=forc_solld) forc_rain => wateratm2lnd_inst%forc_rain_downscaled_col , & ! Input: [real(r8) (:) ] column level rain rate [mm/s] - forc_rain_grc => wateratm2lnd_inst%forc_rain_not_downscaled_grc, & ! Input: [real(r8) (:)] grid cell-level rain rate [mm/s] + forc_rain_grc => wateratm2lnd_inst%forc_rain_not_downscaled_grc, & ! Input: [real(r8) (:) ] grid cell-level rain rate [mm/s] forc_snow => wateratm2lnd_inst%forc_snow_downscaled_col , & ! Input: [real(r8) (:) ] column level snow rate [mm/s] - forc_snow_grc => wateratm2lnd_inst%forc_snow_not_downscaled_grc, & ! Input: [real(r8) (:)] grid cell-level snow rate [mm/s] + forc_snow_grc => wateratm2lnd_inst%forc_snow_not_downscaled_grc, & ! Input: [real(r8) (:) ] grid cell-level snow rate [mm/s] forc_lwrad => atm2lnd_inst%forc_lwrad_downscaled_col , & ! Input: [real(r8) (:) ] downward infrared (longwave) radiation (W/m**2) h2osno_old => waterbalance_inst%h2osno_old_col , & ! Input: [real(r8) (:) ] snow water (mm H2O) at previous time step @@ -461,8 +461,8 @@ subroutine BalanceCheck( bounds, & qflx_snow_grnd_col => waterflux_inst%qflx_snow_grnd_col , & ! Input: [real(r8) (:) ] snow on ground after interception (mm H2O/s) [+] qflx_snwcp_liq => waterflux_inst%qflx_snwcp_liq_col , & ! Input: [real(r8) (:) ] excess liquid h2o due to snow capping (outgoing) (mm H2O /s) [+]` qflx_snwcp_ice => waterflux_inst%qflx_snwcp_ice_col , & ! Input: [real(r8) (:) ] excess solid h2o due to snow capping (outgoing) (mm H2O /s) [+]` - qflx_snwcp_discarded_liq_col => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:) ] column level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] - qflx_snwcp_discarded_ice_col => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:) ] column level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] + qflx_snwcp_discarded_liq_col => waterflux_inst%qflx_snwcp_discarded_liq_col, & ! Input: [real(r8) (:)] column level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] + qflx_snwcp_discarded_ice_col => waterflux_inst%qflx_snwcp_discarded_ice_col, & ! Input: [real(r8) (:)] column level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack (mm H2O /s) [+] qflx_evap_tot_col => waterflux_inst%qflx_evap_tot_col , & ! Input: [real(r8) (:) ] column level qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_evap_tot_grc => waterlnd2atm_inst%qflx_evap_tot_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_soliddew_to_top_layer => waterflux_inst%qflx_soliddew_to_top_layer_col , & ! Input: [real(r8) (:) ] rate of solid water deposited on top soil or snow layer (frost) (mm H2O /s) [+] @@ -473,20 +473,18 @@ subroutine BalanceCheck( bounds, & qflx_snow_h2osfc => waterflux_inst%qflx_snow_h2osfc_col , & ! Input: [real(r8) (:) ] snow falling on surface water (mm/s) qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice qflx_drain_perched_col => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) - qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:) ] grid cell-level sub-surface runoff (mm H2O /s) + qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:)] grid cell-level sub-surface runoff (mm H2O /s) qflx_flood_col => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] column level total runoff due to flooding forc_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] grid cell-level total grid cell-level runoff from river model qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack -! qflx_liq_dynbal_grc => waterflux_inst%qflx_liq_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder -! qflx_ice_dynbal_grc => waterflux_inst%qflx_ice_dynbal_grc , & ! Input: [real(r8) (:) ] slevis: place holder qflx_surf_col => waterflux_inst%qflx_surf_col , & ! Input: [real(r8) (:) ] column level surface runoff (mm H2O /s) qflx_surf_grc => waterlnd2atm_inst%qflx_rofliq_qsur_grc , & ! Input: [real(r8) (:) ] grid cell-level surface runoff (mm H20 /s) qflx_qrgwl_col => waterflux_inst%qflx_qrgwl_col , & ! Input: [real(r8) (:) ] column level qflx_surf at glaciers, wetlands, lakes qflx_qrgwl_grc => waterlnd2atm_inst%qflx_rofliq_qgwl_grc , & ! Input: [real(r8) (:) ] grid cell-level qflx_surf at glaciers, wetlands, lakes qflx_drain_col => waterflux_inst%qflx_drain_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) qflx_drain_grc => waterlnd2atm_inst%qflx_rofliq_qsub_grc , & ! Input: [real(r8) (:) ] grid cell-level drainage (mm H20 /s) - qflx_ice_runoff_col => waterlnd2atm_inst%qflx_ice_runoff_col , & ! Input: [real(r8) (:) ] column level solid runoff from snow capping and from excess ice in soil (mm H2O /s) - qflx_ice_runoff_grc => waterlnd2atm_inst%qflx_rofice_grc , & ! Input: [real(r8) (:) ] grid cell-level solid runoff from snow capping and from excess ice in soil (mm H2O /s) + qflx_ice_runoff_col => waterlnd2atm_inst%qflx_ice_runoff_col , & ! Input: [real(r8) (:) ] column level solid runoff from snow capping and from excess ice in soil (mm H2O /s) + qflx_ice_runoff_grc => waterlnd2atm_inst%qflx_rofice_grc , & ! Input: [real(r8) (:) ] grid cell-level solid runoff from snow capping and from excess ice in soil (mm H2O /s) qflx_sl_top_soil => waterflux_inst%qflx_sl_top_soil_col , & ! Input: [real(r8) (:) ] liquid water + ice from layer above soil to top soil layer or sent to qflx_qrgwl (mm H2O/s) qflx_sfc_irrig_col => waterflux_inst%qflx_sfc_irrig_col , & ! Input: [real(r8) (:) ] column level irrigation flux (mm H2O /s) @@ -585,7 +583,7 @@ subroutine BalanceCheck( bounds, & if (errh2o_max_val > h2o_warning_thresh) then - indexc = maxloc( abs(errh2o_col(bounds%begc:bounds%endc)), 1 ) + bounds%begc -1 + indexc = maxloc( abs(errh2o_col(bounds%begc:bounds%endc)), 1 ) + bounds%begc - 1 write(iulog,*)'WARNING: column-level water balance error ',& ' nstep= ',nstep, & ' local indexc= ',indexc,& @@ -671,7 +669,7 @@ subroutine BalanceCheck( bounds, & if (errh2o_max_val > h2o_warning_thresh) then - indexg = maxloc( abs(errh2o_grc(bounds%begg:bounds%endg)), 1 ) + bounds%begg -1 + indexg = maxloc( abs(errh2o_grc(bounds%begg:bounds%endg)), 1 ) + bounds%begg - 1 write(iulog,*)'WARNING: grid cell-level water balance error ',& ' nstep= ',nstep, & ' local indexg= ',indexg,& @@ -698,7 +696,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'deltawb/dtime = ',(endwb_grc(indexg)-begwb_grc(indexg))/dtime write(iulog,*)'qflx_drain_perched = ',qflx_drain_perched_grc(indexg)*dtime write(iulog,*)'forc_flood = ',forc_flood_grc(indexg)*dtime - write(iulog,*)'qflx_glcice_dyn_water_flux_grc = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime + write(iulog,*)'qflx_glcice_dyn_water_flux = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime write(iulog,*)'clm model is stopping' call endrun(decomp_index=indexg, clmlevel=nameg, msg=errmsg(sourcefile, __LINE__)) From 097ff9805c5d10e55c9be42e6ec19d066a3b963f Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Wed, 23 Dec 2020 12:08:39 -0700 Subject: [PATCH 08/19] First draft of ChangeLog --- doc/ChangeLog | 120 ++++++++++++++++++++++++++++++++++++++++++++++++++ doc/ChangeSum | 1 + 2 files changed, 121 insertions(+) diff --git a/doc/ChangeLog b/doc/ChangeLog index 5642709fc1..a9cbe04d1e 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,4 +1,124 @@ =============================================================== +Tag name: ctsm5.1.dev020 +Originator(s): slevis (Samuel Levis,303-665-1310) +Date: Wed Dec 23 11:29:33 MST 2020 +One-line Summary: Grid cell-level error check for H2O + +Purpose of changes +------------------ + + For more robust mass balance error checking, introduced + grid cell-level error check for H2O following the approach + of pull requests #984 and #1022 + + +Bugs fixed or introduced +------------------------ + +Issues fixed (include CTSM Issue #): #201 + + +Significant changes to scientifically-supported configurations +-------------------------------------------------------------- + +Does this tag change answers significantly for any of the following physics configurations? +(Details of any changes will be given in the "Answer changes" section below.) + + [Put an [X] in the box for any configuration with significant answer changes.] + +[X] clm5_1 + +[X] clm5_0 + +[X] ctsm5_0-nwp + +[X] clm4_5 + +Notes of particular relevance for users +--------------------------------------- + +Caveats for users (e.g., need to interpolate initial conditions): + None + +Changes to CTSM's user interface (e.g., new/renamed XML or namelist variables): + None + +Changes made to namelist defaults (e.g., changed parameter values): + None + +Changes to the datasets (e.g., parameter, surface or initial files): + None + +Substantial timing or memory changes: [For timing changes, can check PFS test(s) in the test suite] + None + +Notes of particular relevance for developers: (including Code reviews and testing) +--------------------------------------------- +NOTE: Be sure to review the steps in README.CHECKLIST.master_tags as well as the coding style in the Developers Guide + +Caveats for developers (e.g., code that is duplicated that requires double maintenance): + None + +Changes to tests or testing: + None + +CTSM testing: + + [PASS means all tests PASS and OK means tests PASS other than expected fails.] + + build-namelist tests: + + cheyenne - + + tools-tests (test/tools): + + cheyenne - + + PTCLM testing (tools/shared/PTCLM/test): + + cheyenne - + + python testing (see instructions in python/README.md; document testing done): + + (any machine) - + + regular tests (aux_clm): + + cheyenne ---- PEND (expect OK) + izumi ------- PEND (expect OK) + +If the tag used for baseline comparisons was NOT the previous tag, note that here: + + +Answer changes +-------------- + +Changes answers relative to baseline: + + Summarize any changes to answers, i.e., + - what code configurations: ALL + - what platforms/compilers: ALL + - nature of change (roundoff; larger than roundoff/same climate; new climate): + Specific example from running the single point test + ERI_D_Ld9.1x1_camdenNJ.I2000Clm50BgcCruRs.cheyenne_intel.clm-default: + RMS ERRH2O 6.0280E-21 NORMALIZED 7.6050E-06 + + Explanation: Moving call BalanceCheck to after the call lnd2glc in + subroutine clm_drv causes a change in order of operations that leads to + the above change in ERRH2O. + + +Detailed list of changes +------------------------ + +List any externals directories updated (cime, rtm, mosart, cism, fates, etc.): + None + +Pull Requests that document the changes (include PR ids): + https://github.com/ESCOMP/ctsm/pull/1228 + +=============================================================== +=============================================================== Tag name: ctsm5.1.dev019 Originator(s): sacks (Bill Sacks) Date: Sat Dec 19 06:55:46 MST 2020 diff --git a/doc/ChangeSum b/doc/ChangeSum index 3789dbee16..6cfae1392f 100644 --- a/doc/ChangeSum +++ b/doc/ChangeSum @@ -1,5 +1,6 @@ Tag Who Date Summary ============================================================================================================================ + ctsm5.1.dev020 slevis 12/23/2020 Grid cell-level error check for H2O ctsm5.1.dev019 sacks 12/19/2020 Fix ndep from coupler ctsm5.1.dev018 slevis 12/08/2020 Add ACTIVE (T/F) column to master hist fields table and alphabetize ctsm5.1.dev017 slevis 11/17/2020 Write history fields master list to separate optional file From cb6a84042a47e785e365ac3629715bd22a5810d0 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Sat, 23 Jan 2021 16:48:50 -0700 Subject: [PATCH 09/19] Revisions part 1 in response to code review --- src/biogeophys/BalanceCheckMod.F90 | 52 ++++++++++------------------- src/biogeophys/WaterBalanceType.F90 | 4 --- src/main/clm_driver.F90 | 3 +- src/main/lnd2atmMod.F90 | 21 ++++++++---- 4 files changed, 34 insertions(+), 46 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index ad2cbd31dd..733f24403d 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -168,8 +168,7 @@ end subroutine BeginWaterGridcellBalance !----------------------------------------------------------------------- subroutine BeginWaterColumnBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, soilhydrology_inst, lakestate_inst, & - use_aquifer_layer) + water_inst, lakestate_inst) ! ! !DESCRIPTION: ! Initialize column-level water balance at beginning of time step, for bulk water and @@ -183,8 +182,6 @@ subroutine BeginWaterColumnBalance(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst type(lakestate_type) , intent(in) :: lakestate_inst - type(soilhydrology_type) , intent(in) :: soilhydrology_inst - logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: integer :: i @@ -196,12 +193,10 @@ subroutine BeginWaterColumnBalance(bounds, & call BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, & num_lakec, filter_lakec, & - soilhydrology_inst, & lakestate_inst, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & - water_inst%bulk_and_tracers(i)%waterbalance_inst, & - use_aquifer_layer = use_aquifer_layer) + water_inst%bulk_and_tracers(i)%waterbalance_inst) end do end subroutine BeginWaterColumnBalance @@ -282,8 +277,12 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & call c2g(bounds, begwb_col(begc:endc), begwb_grc(begg:endg), & c2l_scale_type='urbanf', l2g_scale_type='unity') - call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_beg(bounds, qflx_liq_dynbal_left_to_dribble(begg:endg)) - call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_beg(bounds, qflx_ice_dynbal_left_to_dribble(begg:endg)) + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_beg( & + bounds, & + qflx_liq_dynbal_left_to_dribble(begg:endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_beg( & + bounds, & + qflx_ice_dynbal_left_to_dribble(begg:endg)) do g = begg, endg begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & @@ -297,9 +296,8 @@ end subroutine BeginWaterGridcellBalanceSingle !----------------------------------------------------------------------- subroutine BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - soilhydrology_inst, lakestate_inst, waterstate_inst, & - waterdiagnostic_inst, waterbalance_inst, & - use_aquifer_layer) + lakestate_inst, waterstate_inst, & + waterdiagnostic_inst, waterbalance_inst) ! ! !DESCRIPTION: ! Initialize column-level water balance at beginning of time step, for bulk or a @@ -311,37 +309,19 @@ subroutine BeginWaterColumnBalanceSingle(bounds, & integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points integer , intent(in) :: num_lakec ! number of column lake points in column filter integer , intent(in) :: filter_lakec(:) ! column filter for lake points - type(soilhydrology_type) , intent(in) :: soilhydrology_inst type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst - logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: - integer :: c, j, fc ! indices !----------------------------------------------------------------------- associate( & - zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) - zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) - aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) - wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) begwb => waterbalance_inst%begwb_col , & ! Output: [real(r8) (:) ] water mass begining of the time step h2osno_old => waterbalance_inst%h2osno_old_col & ! Output: [real(r8) (:) ] snow water (mm H2O) at previous time step ) - if(use_aquifer_layer) then - do fc = 1, num_nolakec - c = filter_nolakec(fc) - if (col%hydrologically_active(c)) then - if(zwt(c) <= zi(c,nlevsoi)) then - wa(c) = aquifer_water_baseline - end if - end if - end do - endif - call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, & waterstate_inst, waterdiagnostic_inst, & subtract_dynbal_baselines = .false., & @@ -413,6 +393,7 @@ subroutine BalanceCheck( bounds, & integer :: nstep ! time step number integer :: DAnstep ! time step number since last Data Assimilation (DA) integer :: indexp,indexc,indexl,indexg ! index of first found in search loop + real(r8) :: errh2o_grc(bounds%begg:bounds%endg) ! grid cell level water conservation error [mm H2O] real(r8) :: forc_rain_col(bounds%begc:bounds%endc) ! column level rain rate [mm/s] real(r8) :: forc_snow_col(bounds%begc:bounds%endc) ! column level snow rate [mm/s] real(r8) :: h2osno_total(bounds%begc:bounds%endc) ! total snow water [mm H2O] @@ -450,7 +431,6 @@ subroutine BalanceCheck( bounds, & snow_depth => waterdiagnosticbulk_inst%snow_depth_col , & ! Input: [real(r8) (:) ] snow height (m) begwb_grc => waterbalance_inst%begwb_grc , & ! Input: [real(r8) (:) ] grid cell-level water mass begining of the time step endwb_grc => waterbalance_inst%endwb_grc , & ! Output: [real(r8) (:) ] grid cell-level water mass end of the time step - errh2o_grc => waterbalance_inst%errh2o_grc , & ! Output: [real(r8) (:) ] grid cell-level water conservation error (mm H2O) begwb_col => waterbalance_inst%begwb_col , & ! Input: [real(r8) (:) ] column-level water mass begining of the time step endwb_col => waterbalance_inst%endwb_col , & ! Output: [real(r8) (:) ] column-level water mass end of the time step errh2o_col => waterbalance_inst%errh2o_col , & ! Output: [real(r8) (:) ] column-level water conservation error (mm H2O) @@ -642,8 +622,12 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end(bounds, qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg)) - call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_end(bounds, qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg)) + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end( & + bounds, & + qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_end( & + bounds, & + qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg)) do g = bounds%begg, bounds%endg endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & @@ -704,7 +688,7 @@ subroutine BalanceCheck( bounds, & end if - ! Snow balance check at the grid cell level. + ! Snow balance check at the column level. ! Beginning snow balance variable h2osno_old is calculated once ! for both the column-level and grid cell-level balance checks. diff --git a/src/biogeophys/WaterBalanceType.F90 b/src/biogeophys/WaterBalanceType.F90 index 4d747acebe..6b10895549 100644 --- a/src/biogeophys/WaterBalanceType.F90 +++ b/src/biogeophys/WaterBalanceType.F90 @@ -42,7 +42,6 @@ module WaterBalanceType real(r8), pointer :: endwb_col (:) ! column-level water mass end of the time step real(r8), pointer :: errh2o_patch (:) ! water conservation error (mm H2O) real(r8), pointer :: errh2o_col (:) ! column-level water conservation error (mm H2O) - real(r8), pointer :: errh2o_grc (:) ! grid cell-level water conservation error (mm H2O) real(r8), pointer :: errh2osno_col (:) ! snow water conservation error(mm H2O) contains @@ -133,9 +132,6 @@ subroutine InitAllocate(this, bounds, tracer_vars) call AllocateVar1d(var = this%errh2o_col, name = 'errh2o_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) - call AllocateVar1d(var = this%errh2o_grc, name = 'errh2o_grc', & - container = tracer_vars, & - bounds = bounds, subgrid_level = BOUNDS_SUBGRID_GRIDCELL) call AllocateVar1d(var = this%errh2osno_col, name = 'errh2osno_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index eb94d5de83..2340b5d139 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -391,8 +391,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call BeginWaterColumnBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, soilhydrology_inst, lakestate_inst, & - use_aquifer_layer = use_aquifer_layer()) + water_inst, lakestate_inst) call t_stopf('begwbal') diff --git a/src/main/lnd2atmMod.F90 b/src/main/lnd2atmMod.F90 index d3eb22d610..ed9d44f36f 100644 --- a/src/main/lnd2atmMod.F90 +++ b/src/main/lnd2atmMod.F90 @@ -358,8 +358,9 @@ subroutine lnd2atm(bounds, & ! qflx_runoff is the sum of a number of terms, including qflx_qrgwl. Since we ! are adjusting qflx_qrgwl above, we need to adjust qflx_runoff analogously. - water_inst%waterfluxbulk_inst%qflx_runoff_col(c) = water_inst%waterfluxbulk_inst%qflx_runoff_col(c) + & - water_inst%waterlnd2atmbulk_inst%qflx_liq_from_ice_col(c) + water_inst%waterfluxbulk_inst%qflx_runoff_col(c) = & + water_inst%waterfluxbulk_inst%qflx_runoff_col(c) + & + water_inst%waterlnd2atmbulk_inst%qflx_liq_from_ice_col(c) end if end do @@ -374,8 +375,12 @@ subroutine lnd2atm(bounds, & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg - water_inst%waterlnd2atmbulk_inst%qflx_rofliq_qgwl_grc(g) = water_inst%waterlnd2atmbulk_inst%qflx_rofliq_qgwl_grc(g) - water_inst%waterfluxbulk_inst%qflx_liq_dynbal_grc(g) - water_inst%waterlnd2atmbulk_inst%qflx_rofliq_grc(g) = water_inst%waterlnd2atmbulk_inst%qflx_rofliq_grc(g) - water_inst%waterfluxbulk_inst%qflx_liq_dynbal_grc(g) + water_inst%waterlnd2atmbulk_inst%qflx_rofliq_qgwl_grc(g) = & + water_inst%waterlnd2atmbulk_inst%qflx_rofliq_qgwl_grc(g) - & + water_inst%waterfluxbulk_inst%qflx_liq_dynbal_grc(g) + water_inst%waterlnd2atmbulk_inst%qflx_rofliq_grc(g) = & + water_inst%waterlnd2atmbulk_inst%qflx_rofliq_grc(g) - & + water_inst%waterfluxbulk_inst%qflx_liq_dynbal_grc(g) enddo call c2g( bounds, & @@ -394,7 +399,9 @@ subroutine lnd2atm(bounds, & water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg - water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(g) = water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(g) - water_inst%waterfluxbulk_inst%qflx_ice_dynbal_grc(g) + water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(g) = & + water_inst%waterlnd2atmbulk_inst%qflx_rofice_grc(g) - & + water_inst%waterfluxbulk_inst%qflx_ice_dynbal_grc(g) enddo ! calculate total water storage for history files @@ -407,7 +414,9 @@ subroutine lnd2atm(bounds, & water_inst%waterbalancebulk_inst%endwb_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg - water_inst%waterdiagnosticbulk_inst%tws_grc(g) = water_inst%waterbalancebulk_inst%endwb_grc(g) + water_inst%wateratm2lndbulk_inst%volr_grc(g) / grc%area(g) * 1.e-3_r8 + water_inst%waterdiagnosticbulk_inst%tws_grc(g) = & + water_inst%waterbalancebulk_inst%endwb_grc(g) + & + water_inst%wateratm2lndbulk_inst%volr_grc(g) / grc%area(g) * 1.e-3_r8 enddo end subroutine lnd2atm From 799e011fe71a02f4ddb8c97ed267cb6b5de6ebaf Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Sat, 23 Jan 2021 18:43:36 -0700 Subject: [PATCH 10/19] Revisions part 2: bypass `for_testing_zero_dynbal_fluxes` cases --- src/biogeophys/BalanceCheckMod.F90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 733f24403d..fae724f58d 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -370,6 +370,7 @@ subroutine BalanceCheck( bounds, & use clm_time_manager , only : get_nstep_since_startup_or_lastDA_restart_or_pause use CanopyStateType , only : canopystate_type use subgridAveMod , only : c2g + use dynSubgridControlMod, only : get_for_testing_zero_dynbal_fluxes ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds @@ -658,7 +659,8 @@ subroutine BalanceCheck( bounds, & ' nstep= ',nstep, & ' local indexg= ',indexg,& ' errh2o_grc= ',errh2o_grc(indexg) - if ((errh2o_max_val > error_thresh) .and. (DAnstep > skip_steps)) then + if (errh2o_max_val > error_thresh .and. DAnstep > skip_steps .and. & + .not. get_for_testing_zero_dynbal_fluxes()) then write(iulog,*)'clm model is stopping - error is greater than 1e-5 (mm)' write(iulog,*)'nstep = ',nstep From 23bf38baa6142b478ccff8c0cb7326d402edefb5 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Mon, 25 Jan 2021 14:48:21 -0700 Subject: [PATCH 11/19] Revisions part 3: Insert comment about negative dribble terms --- src/biogeophys/BalanceCheckMod.F90 | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index fae724f58d..f1df93fa4a 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -284,6 +284,10 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & bounds, & qflx_ice_dynbal_left_to_dribble(begg:endg)) + ! These dynbal dribblers store the delta state, (end - beg). Thus, the + ! amount dribbled out is the negative of the amount stored in the + ! dribblers. Therefore, conservation requires us to subtract the amount + ! remaining to dribble. do g = begg, endg begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - qflx_ice_dynbal_left_to_dribble(g) @@ -630,6 +634,10 @@ subroutine BalanceCheck( bounds, & bounds, & qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg)) + ! These dynbal dribblers store the delta state, (end - beg). Thus, the + ! amount dribbled out is the negative of the amount stored in the + ! dribblers. Therefore, conservation requires us to subtract the amount + ! remaining to dribble. do g = bounds%begg, bounds%endg endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - qflx_ice_dynbal_left_to_dribble(g) From c5b09617be65bab70f00820e1162713f50883058 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Mon, 25 Jan 2021 18:05:33 -0700 Subject: [PATCH 12/19] Addition of text to the comment that I added in the last commit --- src/biogeophys/BalanceCheckMod.F90 | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index f1df93fa4a..56dabca84c 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -288,6 +288,11 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & ! amount dribbled out is the negative of the amount stored in the ! dribblers. Therefore, conservation requires us to subtract the amount ! remaining to dribble. + ! This sign convention is opposite to the convention chosen for the + ! respective dribble terms used in the carbon balance. At some point + ! it may be worth making the two conventions consistent. + ! Bill Sacks states: I think the convention used for the water and + ! energy dribblers is counter-intuitive. do g = begg, endg begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - qflx_ice_dynbal_left_to_dribble(g) @@ -638,6 +643,11 @@ subroutine BalanceCheck( bounds, & ! amount dribbled out is the negative of the amount stored in the ! dribblers. Therefore, conservation requires us to subtract the amount ! remaining to dribble. + ! This sign convention is opposite to the convention chosen for the + ! respective dribble terms used in the carbon balance. At some point + ! it may be worth making the two conventions consistent. + ! Bill Sacks states: I think the convention used for the water and + ! energy dribblers is counter-intuitive. do g = bounds%begg, bounds%endg endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - qflx_ice_dynbal_left_to_dribble(g) From 6e38fc28ebab92d41981f4097c2f9ef039d14321 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Wed, 27 Jan 2021 15:46:07 -0700 Subject: [PATCH 13/19] Revisions part 4: Tracking wa_reset_nonconservation_gain --- src/biogeophys/BalanceCheckMod.F90 | 81 ++++++++++++++++++++--------- src/biogeophys/WaterBalanceType.F90 | 23 +++++++- src/main/clm_driver.F90 | 5 +- 3 files changed, 80 insertions(+), 29 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 56dabca84c..50575eb5c9 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -126,7 +126,7 @@ end function GetBalanceCheckSkipSteps !----------------------------------------------------------------------- subroutine BeginWaterGridcellBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, soilhydrology_inst, lakestate_inst, & + water_inst, lakestate_inst, & use_aquifer_layer) ! ! !DESCRIPTION: @@ -141,7 +141,6 @@ subroutine BeginWaterGridcellBalance(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst type(lakestate_type) , intent(in) :: lakestate_inst - type(soilhydrology_type), intent(in) :: soilhydrology_inst logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: @@ -154,7 +153,6 @@ subroutine BeginWaterGridcellBalance(bounds, & call BeginWaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, & num_lakec, filter_lakec, & - soilhydrology_inst, & lakestate_inst, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & @@ -168,7 +166,8 @@ end subroutine BeginWaterGridcellBalance !----------------------------------------------------------------------- subroutine BeginWaterColumnBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, lakestate_inst) + water_inst, soilhydrology_inst, lakestate_inst, & + use_aquifer_layer) ! ! !DESCRIPTION: ! Initialize column-level water balance at beginning of time step, for bulk water and @@ -182,6 +181,8 @@ subroutine BeginWaterColumnBalance(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst type(lakestate_type) , intent(in) :: lakestate_inst + type(soilhydrology_type) , intent(in) :: soilhydrology_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: integer :: i @@ -193,10 +194,12 @@ subroutine BeginWaterColumnBalance(bounds, & call BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, & num_lakec, filter_lakec, & + soilhydrology_inst, & lakestate_inst, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & - water_inst%bulk_and_tracers(i)%waterbalance_inst) + water_inst%bulk_and_tracers(i)%waterbalance_inst, & + use_aquifer_layer = use_aquifer_layer) end do end subroutine BeginWaterColumnBalance @@ -204,7 +207,7 @@ end subroutine BeginWaterColumnBalance !----------------------------------------------------------------------- subroutine BeginWaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - soilhydrology_inst, lakestate_inst, waterstate_inst, & + lakestate_inst, waterstate_inst, & waterdiagnostic_inst, waterbalance_inst, waterflux_inst, & use_aquifer_layer) ! @@ -221,7 +224,6 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points integer , intent(in) :: num_lakec ! number of column lake points in column filter integer , intent(in) :: filter_lakec(:) ! column filter for lake points - type(soilhydrology_type) , intent(in) :: soilhydrology_inst type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst class(waterflux_type) , intent(inout) :: waterflux_inst @@ -230,17 +232,15 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: - integer :: c, g, j, fc ! indices + integer :: g ! indices integer :: begc, endc, begg, endg ! bounds real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at beginning of time step real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc ice dynamic land cover change conversion runoff flux at beginning of time step + real(r8) :: wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) ! grc mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) !----------------------------------------------------------------------- associate( & - zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) - zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) - aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) - wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) + wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb_col => waterbalance_inst%begwb_col, & ! Output: [real(r8) (:) ] column-level water mass begining of the time step begwb_grc => waterbalance_inst%begwb_grc & ! Output: [real(r8) (:) ] grid cell-level water mass begining of the time step ) @@ -250,17 +250,20 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & begg = bounds%begg endg = bounds%endg - ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake - if(use_aquifer_layer) then - do fc = 1, num_nolakec - c = filter_nolakec(fc) - if (col%hydrologically_active(c)) then - if(zwt(c) <= zi(c,nlevsoi)) then - wa(c) = aquifer_water_baseline - end if - end if - end do - endif + if (use_aquifer_layer) then + ! wa_reset_nonconservation_gain_col should be non-zero only when + ! use_aquifer_layer is true. We do this c2g call only when needed + ! to avoid unnecessary calculations; by adding this term only when + ! use_aquifer_layer is true, we effectively let the balance checks + ! ensure that this term is zero when use_aquifer_layer is false, + ! as it should be. + call c2g( bounds, & + wa_reset_nonconservation_gain_col(begc:endc), & + wa_reset_nonconservation_gain_grc(begg:endg), & + c2l_scale_type='urbanf', l2g_scale_type='unity' ) + else + wa_reset_nonconservation_gain_grc(begg:endg) = 0._r8 + end if ! NOTES subroutines Compute*Mass* are in TotalWaterAndHeatMod.F90 ! endwb is calculated in HydrologyDrainageMod & LakeHydrologyMod @@ -295,7 +298,8 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & ! energy dribblers is counter-intuitive. do g = begg, endg begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - - qflx_ice_dynbal_left_to_dribble(g) + - qflx_ice_dynbal_left_to_dribble(g) & + - wa_reset_nonconservation_gain_grc(g) end do end associate @@ -305,8 +309,9 @@ end subroutine BeginWaterGridcellBalanceSingle !----------------------------------------------------------------------- subroutine BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - lakestate_inst, waterstate_inst, & - waterdiagnostic_inst, waterbalance_inst) + soilhydrology_inst, lakestate_inst, waterstate_inst, & + waterdiagnostic_inst, waterbalance_inst, & + use_aquifer_layer) ! ! !DESCRIPTION: ! Initialize column-level water balance at beginning of time step, for bulk or a @@ -318,19 +323,43 @@ subroutine BeginWaterColumnBalanceSingle(bounds, & integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points integer , intent(in) :: num_lakec ! number of column lake points in column filter integer , intent(in) :: filter_lakec(:) ! column filter for lake points + type(soilhydrology_type) , intent(in) :: soilhydrology_inst type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: + integer :: c, fc ! indices !----------------------------------------------------------------------- associate( & + zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) + zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) + aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) + wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) + wa_reset_nonconservation_gain => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb => waterbalance_inst%begwb_col , & ! Output: [real(r8) (:) ] water mass begining of the time step h2osno_old => waterbalance_inst%h2osno_old_col & ! Output: [real(r8) (:) ] snow water (mm H2O) at previous time step ) + ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake + if(use_aquifer_layer) then + do fc = 1, num_nolakec + c = filter_nolakec(fc) + if (col%hydrologically_active(c)) then + if(zwt(c) <= zi(c,nlevsoi)) then + wa_reset_nonconservation_gain(c) = aquifer_water_baseline - & + wa(c) + wa(c) = aquifer_water_baseline + else + wa_reset_nonconservation_gain(c) = 0._r8 + end if + end if + end do + endif + call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, & waterstate_inst, waterdiagnostic_inst, & subtract_dynbal_baselines = .false., & diff --git a/src/biogeophys/WaterBalanceType.F90 b/src/biogeophys/WaterBalanceType.F90 index 6b10895549..0bf0573913 100644 --- a/src/biogeophys/WaterBalanceType.F90 +++ b/src/biogeophys/WaterBalanceType.F90 @@ -33,6 +33,7 @@ module WaterBalanceType real(r8), pointer :: snow_sources_col (:) ! col snow sources (mm H2O/s) real(r8), pointer :: snow_sinks_col (:) ! col snow sinks (mm H2O/s) + real(r8), pointer :: wa_reset_nonconservation_gain_col(:) ! col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) ! Balance Checks @@ -49,6 +50,7 @@ module WaterBalanceType procedure :: Init procedure, private :: InitAllocate procedure, private :: InitHistory + procedure, private :: InitCold end type waterbalance_type @@ -70,8 +72,8 @@ subroutine Init(this, bounds, info, tracer_vars) this%info => info call this%InitAllocate(bounds, tracer_vars) - call this%InitHistory(bounds) + call this%InitCold(bounds) end subroutine Init @@ -113,6 +115,9 @@ subroutine InitAllocate(this, bounds, tracer_vars) call AllocateVar1d(var = this%snow_sinks_col, name = 'snow_sinks_col', & container = tracer_vars, & bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) + call AllocateVar1d(var = this%wa_reset_nonconservation_gain_col, name = 'wa_reset_nonconservation_gain_col', & + container = tracer_vars, & + bounds = bounds, subgrid_level = BOUNDS_SUBGRID_COLUMN) call AllocateVar1d(var = this%begwb_grc, name = 'begwb_grc', & container = tracer_vars, & @@ -233,4 +238,20 @@ subroutine InitHistory(this, bounds) ptr_col=this%errh2osno_col, c2l_scale_type='urbanf') end subroutine InitHistory + !----------------------------------------------------------------------- + subroutine InitCold(this, bounds) + ! + ! !USES: + ! + ! !ARGUMENTS: + class(waterbalance_type), intent(in) :: this + type(bounds_type) , intent(in) :: bounds + ! + ! !LOCAL VARIABLES: + !----------------------------------------------------------------------- + + this%wa_reset_nonconservation_gain_col(bounds%begc:bounds%endc) = 0.0_r8 + + end subroutine InitCold + end module WaterBalanceType diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index 2340b5d139..7504ae1c28 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -330,7 +330,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call BeginWaterGridcellBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, soilhydrology_inst, lakestate_inst, & + water_inst, lakestate_inst, & use_aquifer_layer = use_aquifer_layer()) call t_stopf('begwbal') end do @@ -391,7 +391,8 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call BeginWaterColumnBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, lakestate_inst) + water_inst, soilhydrology_inst, lakestate_inst, & + use_aquifer_layer = use_aquifer_layer()) call t_stopf('begwbal') From 6f0bf37db5e22422b7d58474989ee331957ee4c2 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Thu, 28 Jan 2021 15:12:05 -0700 Subject: [PATCH 14/19] Minor cleanup of comments and write statements --- src/biogeophys/BalanceCheckMod.F90 | 42 ++++++++++++++++-------------- 1 file changed, 23 insertions(+), 19 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 50575eb5c9..84f4bc9484 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -240,7 +240,7 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & !----------------------------------------------------------------------- associate( & - wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) + wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb_col => waterbalance_inst%begwb_col, & ! Output: [real(r8) (:) ] column-level water mass begining of the time step begwb_grc => waterbalance_inst%begwb_grc & ! Output: [real(r8) (:) ] grid cell-level water mass begining of the time step ) @@ -251,7 +251,7 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & endg = bounds%endg if (use_aquifer_layer) then - ! wa_reset_nonconservation_gain_col should be non-zero only when + ! wa_reset_nonconservation_gain may be non-zero only when ! use_aquifer_layer is true. We do this c2g call only when needed ! to avoid unnecessary calculations; by adding this term only when ! use_aquifer_layer is true, we effectively let the balance checks @@ -309,7 +309,7 @@ end subroutine BeginWaterGridcellBalanceSingle !----------------------------------------------------------------------- subroutine BeginWaterColumnBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - soilhydrology_inst, lakestate_inst, waterstate_inst, & + soilhydrology_inst, lakestate_inst, waterstate_inst, & waterdiagnostic_inst, waterbalance_inst, & use_aquifer_layer) ! @@ -335,16 +335,20 @@ subroutine BeginWaterColumnBalanceSingle(bounds, & !----------------------------------------------------------------------- associate( & - zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) - zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) + zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) + zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) - wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) - wa_reset_nonconservation_gain => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) + wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) + wa_reset_nonconservation_gain => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb => waterbalance_inst%begwb_col , & ! Output: [real(r8) (:) ] water mass begining of the time step h2osno_old => waterbalance_inst%h2osno_old_col & ! Output: [real(r8) (:) ] snow water (mm H2O) at previous time step ) ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake + ! wa_reset_nonconservation_gain is calculated for the grid cell-level + ! water balance check and may be non-zero only when + ! use_aquifer_layer is true. The grid cell-level balance check ensures + ! that this term is zero when use_aquifer_layer is false, as it should be. if(use_aquifer_layer) then do fc = 1, num_nolakec c = filter_nolakec(fc) @@ -610,7 +614,7 @@ subroutine BalanceCheck( bounds, & ' errh2o= ',errh2o_col(indexc) if ((errh2o_max_val > error_thresh) .and. (DAnstep > skip_steps)) then - write(iulog,*)'clm urban model is stopping - error is greater than 1e-5 (mm)' + write(iulog,*)'CTSM model stopping because errh2o > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2o_col = ',errh2o_col(indexc) write(iulog,*)'forc_rain = ',forc_rain_col(indexc)*dtime @@ -639,7 +643,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_glcice_dyn_water_flux = ', qflx_glcice_dyn_water_flux_col(indexc)*dtime end if - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if @@ -709,7 +713,7 @@ subroutine BalanceCheck( bounds, & if (errh2o_max_val > error_thresh .and. DAnstep > skip_steps .and. & .not. get_for_testing_zero_dynbal_fluxes()) then - write(iulog,*)'clm model is stopping - error is greater than 1e-5 (mm)' + write(iulog,*)'CTSM model stopping because errh2o > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2o_grc = ',errh2o_grc(indexg) write(iulog,*)'forc_rain = ',forc_rain_grc(indexg)*dtime @@ -731,7 +735,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'forc_flood = ',forc_flood_grc(indexg)*dtime write(iulog,*)'qflx_glcice_dyn_water_flux = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexg, clmlevel=nameg, msg=errmsg(sourcefile, __LINE__)) end if @@ -811,7 +815,7 @@ subroutine BalanceCheck( bounds, & ' errh2osno= ',errh2osno(indexc) if ((errh2osno_max_val > error_thresh) .and. (DAnstep > skip_steps) ) then - write(iulog,*)'clm model is stopping - error is greater than 1e-5 (mm)' + write(iulog,*)'CTSM model stopping because errh2osno > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2osno = ',errh2osno(indexc) write(iulog,*)'snl = ',col%snl(indexc) @@ -834,7 +838,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice_col(indexc)*dtime write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq_col(indexc)*dtime write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc)*dtime - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if @@ -906,7 +910,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errsol = ',errsol(indexp) if (errsol_max_val > error_thresh) then - write(iulog,*)'clm model is stopping - error is greater than 1e-5 (W/m2)' + write(iulog,*)'CTSM model stopping because errsol > ', error_thresh, ' W/m2' write(iulog,*)'fsa = ',fsa(indexp) write(iulog,*)'fsr = ',fsr(indexp) write(iulog,*)'forc_solad(1) = ',forc_solad(indexg,1) @@ -915,7 +919,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'forc_solai(2) = ',forc_solai(indexg,2) write(iulog,*)'forc_tot = ',forc_solad(indexg,1)+forc_solad(indexg,2) & +forc_solai(indexg,1)+forc_solai(indexg,2) - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if @@ -932,7 +936,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'nstep = ',nstep write(iulog,*)'errlon = ',errlon(indexp) if (errlon_max_val > error_thresh ) then - write(iulog,*)'clm model is stopping - error is greater than 1e-5 (W/m2)' + write(iulog,*)'CTSM model stopping because errlon > ', error_thresh, ' W/m2' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if end if @@ -952,7 +956,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errseb = ' ,errseb(indexp) if ( errseb_max_val > error_thresh ) then - write(iulog,*)'clm model is stopping - error is greater than 1e-5 (W/m2)' + write(iulog,*)'CTSM model stopping because errseb > ', error_thresh, ' W/m2' write(iulog,*)'sabv = ' ,sabv(indexp) write(iulog,*)'sabg = ' ,sabg(indexp), ((1._r8- frac_sno(indexc))*sabg_soil(indexp) + & frac_sno(indexc)*sabg_snow(indexp)),sabg_chk(indexp) @@ -968,7 +972,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'albd albi = ' ,albd(indexp,:), albi(indexp,:) write(iulog,*)'ftii ftdd ftid = ' ,ftii(indexp,:), ftdd(indexp,:),ftid(indexp,:) write(iulog,*)'elai esai = ' ,elai(indexp), esai(indexp) - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if @@ -985,7 +989,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errsoi_col = ',errsoi_col(indexc) if ((errsoi_col_max_val > 1.e-4_r8) .and. (DAnstep > skip_steps)) then - write(iulog,*)'clm model is stopping' + write(iulog,*)'CTSM model is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if end if From bf52b8ae0843dd7f3f3a6262847ebb103cc32ec8 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Fri, 29 Jan 2021 12:49:05 -0700 Subject: [PATCH 15/19] Cleanup of comments, write-statmts, and formatting --- doc/ChangeLog | 14 ++++----- src/biogeophys/BalanceCheckMod.F90 | 50 +++++++++++++++++++----------- 2 files changed, 39 insertions(+), 25 deletions(-) diff --git a/doc/ChangeLog b/doc/ChangeLog index a9cbe04d1e..dfa751387f 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,5 +1,5 @@ =============================================================== -Tag name: ctsm5.1.dev020 +Tag name: ctsm5.1.dev0?? Originator(s): slevis (Samuel Levis,303-665-1310) Date: Wed Dec 23 11:29:33 MST 2020 One-line Summary: Grid cell-level error check for H2O @@ -26,13 +26,13 @@ Does this tag change answers significantly for any of the following physics conf [Put an [X] in the box for any configuration with significant answer changes.] -[X] clm5_1 +[ ] clm5_1 -[X] clm5_0 +[ ] clm5_0 -[X] ctsm5_0-nwp +[ ] ctsm5_0-nwp -[X] clm4_5 +[ ] clm4_5 Notes of particular relevance for users --------------------------------------- @@ -93,12 +93,12 @@ If the tag used for baseline comparisons was NOT the previous tag, note that her Answer changes -------------- -Changes answers relative to baseline: +Changes answers relative to baseline: YES Summarize any changes to answers, i.e., - what code configurations: ALL - what platforms/compilers: ALL - - nature of change (roundoff; larger than roundoff/same climate; new climate): + - nature of change: ROUNDOFF Specific example from running the single point test ERI_D_Ld9.1x1_camdenNJ.I2000Clm50BgcCruRs.cheyenne_intel.clm-default: RMS ERRH2O 6.0280E-21 NORMALIZED 7.6050E-06 diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 84f4bc9484..b36f23893a 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -257,6 +257,9 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & ! use_aquifer_layer is true, we effectively let the balance checks ! ensure that this term is zero when use_aquifer_layer is false, ! as it should be. + ! The _col term converted to _grc here gets determined in + ! BeginWaterColumnBalanceSingle in the previous time step after any + ! dynamic landuse adjustments. call c2g( bounds, & wa_reset_nonconservation_gain_col(begc:endc), & wa_reset_nonconservation_gain_grc(begg:endg), & @@ -335,20 +338,33 @@ subroutine BeginWaterColumnBalanceSingle(bounds, & !----------------------------------------------------------------------- associate( & - zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) - zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) + zi => col%zi , & ! Input: [real(r8) (:,:) ] interface level below a "z" level (m) + zwt => soilhydrology_inst%zwt_col , & ! Input: [real(r8) (:) ] water table depth (m) aquifer_water_baseline => waterstate_inst%aquifer_water_baseline, & ! Input: [real(r8)] baseline value for water in the unconfined aquifer (wa_col) for this bulk / tracer (mm) - wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) + wa => waterstate_inst%wa_col , & ! Output: [real(r8) (:) ] water in the unconfined aquifer (mm) wa_reset_nonconservation_gain => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb => waterbalance_inst%begwb_col , & ! Output: [real(r8) (:) ] water mass begining of the time step h2osno_old => waterbalance_inst%h2osno_old_col & ! Output: [real(r8) (:) ] snow water (mm H2O) at previous time step ) - ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake + ! wa(c) gets added to liquid_mass in ComputeLiqIceMassNonLake called here. ! wa_reset_nonconservation_gain is calculated for the grid cell-level ! water balance check and may be non-zero only when ! use_aquifer_layer is true. The grid cell-level balance check ensures ! that this term is zero when use_aquifer_layer is false, as it should be. + ! In particular, we adjust wa back to the baseline under certain + ! conditions. The right way to do this might be to use explicit fluxes from + ! some other state, but in this case we don't have a source to pull from, + ! so we adjust wa without explicit fluxes. Because we do this before + ! initializing the column-level balance check, the column-level check is + ! unaware of the adjustment. However, since this adjustment happens after + ! initializing the gridcell-level balance check, we have to account for + ! it in the gridcell-level balance check. The normal way to account for an + ! adjustment like this would be to include the flux in the balance check. + ! Here we don't have an explicit flux, so instead we track the + ! non-conservation state. In principle, we could calculate an explicit flux + ! and use that, but we don't gain anything from using an explicit flux in + ! this case. if(use_aquifer_layer) then do fc = 1, num_nolakec c = filter_nolakec(fc) @@ -614,7 +630,7 @@ subroutine BalanceCheck( bounds, & ' errh2o= ',errh2o_col(indexc) if ((errh2o_max_val > error_thresh) .and. (DAnstep > skip_steps)) then - write(iulog,*)'CTSM model stopping because errh2o > ', error_thresh, ' mm' + write(iulog,*)'CTSM is stopping because errh2o > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2o_col = ',errh2o_col(indexc) write(iulog,*)'forc_rain = ',forc_rain_col(indexc)*dtime @@ -643,7 +659,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_glcice_dyn_water_flux = ', qflx_glcice_dyn_water_flux_col(indexc)*dtime end if - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if @@ -713,7 +729,7 @@ subroutine BalanceCheck( bounds, & if (errh2o_max_val > error_thresh .and. DAnstep > skip_steps .and. & .not. get_for_testing_zero_dynbal_fluxes()) then - write(iulog,*)'CTSM model stopping because errh2o > ', error_thresh, ' mm' + write(iulog,*)'CTSM is stopping because errh2o > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2o_grc = ',errh2o_grc(indexg) write(iulog,*)'forc_rain = ',forc_rain_grc(indexg)*dtime @@ -735,15 +751,13 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'forc_flood = ',forc_flood_grc(indexg)*dtime write(iulog,*)'qflx_glcice_dyn_water_flux = ',qflx_glcice_dyn_water_flux_grc(indexg)*dtime - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexg, clmlevel=nameg, msg=errmsg(sourcefile, __LINE__)) end if end if ! Snow balance check at the column level. - ! Beginning snow balance variable h2osno_old is calculated once - ! for both the column-level and grid cell-level balance checks. call waterstate_inst%CalculateTotalH2osno(bounds, num_allc, filter_allc, & caller = 'BalanceCheck', & @@ -815,7 +829,7 @@ subroutine BalanceCheck( bounds, & ' errh2osno= ',errh2osno(indexc) if ((errh2osno_max_val > error_thresh) .and. (DAnstep > skip_steps) ) then - write(iulog,*)'CTSM model stopping because errh2osno > ', error_thresh, ' mm' + write(iulog,*)'CTSM is stopping because errh2osno > ', error_thresh, ' mm' write(iulog,*)'nstep = ',nstep write(iulog,*)'errh2osno = ',errh2osno(indexc) write(iulog,*)'snl = ',col%snl(indexc) @@ -838,7 +852,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'qflx_snwcp_discarded_ice = ',qflx_snwcp_discarded_ice_col(indexc)*dtime write(iulog,*)'qflx_snwcp_discarded_liq = ',qflx_snwcp_discarded_liq_col(indexc)*dtime write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc)*dtime - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if @@ -910,7 +924,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errsol = ',errsol(indexp) if (errsol_max_val > error_thresh) then - write(iulog,*)'CTSM model stopping because errsol > ', error_thresh, ' W/m2' + write(iulog,*)'CTSM is stopping because errsol > ', error_thresh, ' W/m2' write(iulog,*)'fsa = ',fsa(indexp) write(iulog,*)'fsr = ',fsr(indexp) write(iulog,*)'forc_solad(1) = ',forc_solad(indexg,1) @@ -919,7 +933,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'forc_solai(2) = ',forc_solai(indexg,2) write(iulog,*)'forc_tot = ',forc_solad(indexg,1)+forc_solad(indexg,2) & +forc_solai(indexg,1)+forc_solai(indexg,2) - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if @@ -936,7 +950,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'nstep = ',nstep write(iulog,*)'errlon = ',errlon(indexp) if (errlon_max_val > error_thresh ) then - write(iulog,*)'CTSM model stopping because errlon > ', error_thresh, ' W/m2' + write(iulog,*)'CTSM is stopping because errlon > ', error_thresh, ' W/m2' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if end if @@ -956,7 +970,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errseb = ' ,errseb(indexp) if ( errseb_max_val > error_thresh ) then - write(iulog,*)'CTSM model stopping because errseb > ', error_thresh, ' W/m2' + write(iulog,*)'CTSM is stopping because errseb > ', error_thresh, ' W/m2' write(iulog,*)'sabv = ' ,sabv(indexp) write(iulog,*)'sabg = ' ,sabg(indexp), ((1._r8- frac_sno(indexc))*sabg_soil(indexp) + & frac_sno(indexc)*sabg_snow(indexp)),sabg_chk(indexp) @@ -972,7 +986,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'albd albi = ' ,albd(indexp,:), albi(indexp,:) write(iulog,*)'ftii ftdd ftid = ' ,ftii(indexp,:), ftdd(indexp,:),ftid(indexp,:) write(iulog,*)'elai esai = ' ,elai(indexp), esai(indexp) - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexp, clmlevel=namep, msg=errmsg(sourcefile, __LINE__)) end if @@ -989,7 +1003,7 @@ subroutine BalanceCheck( bounds, & write(iulog,*)'errsoi_col = ',errsoi_col(indexc) if ((errsoi_col_max_val > 1.e-4_r8) .and. (DAnstep > skip_steps)) then - write(iulog,*)'CTSM model is stopping' + write(iulog,*)'CTSM is stopping' call endrun(decomp_index=indexc, clmlevel=namec, msg=errmsg(sourcefile, __LINE__)) end if end if From 5dbbc81bbf4ab78dd671ad865266085a06ac30bf Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Fri, 29 Jan 2021 15:02:28 -0700 Subject: [PATCH 16/19] Revisions 4b: Move code from the beg. to the ending of h2o balance check --- src/biogeophys/BalanceCheckMod.F90 | 59 ++++++++++++++---------------- src/main/clm_driver.F90 | 6 +-- 2 files changed, 30 insertions(+), 35 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index b36f23893a..c1b672975a 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -126,8 +126,7 @@ end function GetBalanceCheckSkipSteps !----------------------------------------------------------------------- subroutine BeginWaterGridcellBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, lakestate_inst, & - use_aquifer_layer) + water_inst, lakestate_inst) ! ! !DESCRIPTION: ! Initialize grid cell-level water balance at beginning of time step @@ -141,7 +140,6 @@ subroutine BeginWaterGridcellBalance(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst type(lakestate_type) , intent(in) :: lakestate_inst - logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: integer :: i @@ -157,8 +155,7 @@ subroutine BeginWaterGridcellBalance(bounds, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & water_inst%bulk_and_tracers(i)%waterbalance_inst, & - water_inst%bulk_and_tracers(i)%waterflux_inst, & - use_aquifer_layer = use_aquifer_layer) + water_inst%bulk_and_tracers(i)%waterflux_inst) end do end subroutine BeginWaterGridcellBalance @@ -208,8 +205,7 @@ end subroutine BeginWaterColumnBalance subroutine BeginWaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & lakestate_inst, waterstate_inst, & - waterdiagnostic_inst, waterbalance_inst, waterflux_inst, & - use_aquifer_layer) + waterdiagnostic_inst, waterbalance_inst, waterflux_inst) ! ! !DESCRIPTION: ! Initialize grid cell-level water balance at beginning of time step @@ -229,18 +225,15 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & class(waterflux_type) , intent(inout) :: waterflux_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst - logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: integer :: g ! indices integer :: begc, endc, begg, endg ! bounds real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at beginning of time step real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc ice dynamic land cover change conversion runoff flux at beginning of time step - real(r8) :: wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) ! grc mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) !----------------------------------------------------------------------- associate( & - wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col , & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) begwb_col => waterbalance_inst%begwb_col, & ! Output: [real(r8) (:) ] column-level water mass begining of the time step begwb_grc => waterbalance_inst%begwb_grc & ! Output: [real(r8) (:) ] grid cell-level water mass begining of the time step ) @@ -250,24 +243,6 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & begg = bounds%begg endg = bounds%endg - if (use_aquifer_layer) then - ! wa_reset_nonconservation_gain may be non-zero only when - ! use_aquifer_layer is true. We do this c2g call only when needed - ! to avoid unnecessary calculations; by adding this term only when - ! use_aquifer_layer is true, we effectively let the balance checks - ! ensure that this term is zero when use_aquifer_layer is false, - ! as it should be. - ! The _col term converted to _grc here gets determined in - ! BeginWaterColumnBalanceSingle in the previous time step after any - ! dynamic landuse adjustments. - call c2g( bounds, & - wa_reset_nonconservation_gain_col(begc:endc), & - wa_reset_nonconservation_gain_grc(begg:endg), & - c2l_scale_type='urbanf', l2g_scale_type='unity' ) - else - wa_reset_nonconservation_gain_grc(begg:endg) = 0._r8 - end if - ! NOTES subroutines Compute*Mass* are in TotalWaterAndHeatMod.F90 ! endwb is calculated in HydrologyDrainageMod & LakeHydrologyMod call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, & @@ -301,8 +276,7 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & ! energy dribblers is counter-intuitive. do g = begg, endg begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - - qflx_ice_dynbal_left_to_dribble(g) & - - wa_reset_nonconservation_gain_grc(g) + - qflx_ice_dynbal_left_to_dribble(g) end do end associate @@ -406,7 +380,8 @@ subroutine BalanceCheck( bounds, & num_allc, filter_allc, & atm2lnd_inst, solarabs_inst, waterflux_inst, waterstate_inst, & waterdiagnosticbulk_inst, waterbalance_inst, wateratm2lnd_inst, & - waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst) + waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst, & + use_aquifer_layer) ! ! !DESCRIPTION: ! This subroutine accumulates the numerical truncation errors of the water @@ -445,6 +420,7 @@ subroutine BalanceCheck( bounds, & type(surfalb_type) , intent(in) :: surfalb_inst type(energyflux_type) , intent(inout) :: energyflux_inst type(canopystate_type), intent(inout) :: canopystate_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run ! ! !LOCAL VARIABLES: integer :: p,c,l,g,fc ! indices @@ -461,6 +437,7 @@ subroutine BalanceCheck( bounds, & real(r8) :: qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg) ! grid cell-level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step + real(r8) :: wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) ! grc mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) real(r8) :: errh2o_max_val ! Maximum value of error in water conservation error over all columns [mm H2O] real(r8) :: errh2osno_max_val ! Maximum value of error in h2osno conservation error over all columns [kg m-2] @@ -513,6 +490,7 @@ subroutine BalanceCheck( bounds, & qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice qflx_drain_perched_col => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:)] grid cell-level sub-surface runoff (mm H2O /s) + wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col, & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) qflx_flood_col => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] column level total runoff due to flooding forc_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] grid cell-level total grid cell-level runoff from river model qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack @@ -680,6 +658,22 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_ice_col(bounds%begc:bounds%endc), & qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) + if (use_aquifer_layer) then + ! wa_reset_nonconservation_gain may be non-zero only when + ! use_aquifer_layer is true. We do this c2g call only when needed + ! to avoid unnecessary calculations; by adding this term only when + ! use_aquifer_layer is true, we effectively let the balance checks + ! ensure that this term is zero when use_aquifer_layer is false, + ! as it should be. + ! The _col term was determined in BeginWaterColumnBalanceSingle + ! after any dynamic landuse adjustments. + call c2g( bounds, & + wa_reset_nonconservation_gain_col(bounds%begc:bounds%endc), & + wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg), & + c2l_scale_type='urbanf', l2g_scale_type='unity' ) + else + wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) = 0._r8 + end if call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end( & bounds, & @@ -699,7 +693,8 @@ subroutine BalanceCheck( bounds, & ! energy dribblers is counter-intuitive. do g = bounds%begg, bounds%endg endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - - qflx_ice_dynbal_left_to_dribble(g) + - qflx_ice_dynbal_left_to_dribble(g) & + - wa_reset_nonconservation_gain_grc(g) errh2o_grc(g) = endwb_grc(g) - begwb_grc(g) & - (forc_rain_grc(g) & diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index 7504ae1c28..4f88a5ef90 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -330,8 +330,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call BeginWaterGridcellBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, lakestate_inst, & - use_aquifer_layer = use_aquifer_layer()) + water_inst, lakestate_inst) call t_stopf('begwbal') end do !$OMP END PARALLEL DO @@ -1271,7 +1270,8 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro atm2lnd_inst, solarabs_inst, water_inst%waterfluxbulk_inst, & water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, & water_inst%waterbalancebulk_inst, water_inst%wateratm2lndbulk_inst, & - water_inst%waterlnd2atmbulk_inst, surfalb_inst, energyflux_inst, canopystate_inst) + water_inst%waterlnd2atmbulk_inst, surfalb_inst, energyflux_inst, & + canopystate_inst, use_aquifer_layer = use_aquifer_layer()) end do !$OMP END PARALLEL DO call t_stopf('balchk') From 6768534c40e9a185a45aff256af819160bd731ae Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Sat, 13 Feb 2021 17:28:40 -0700 Subject: [PATCH 17/19] Bypass grc-level h2o check when use_soil_moisture_streams = .true. --- src/biogeophys/BalanceCheckMod.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index c0aafe60e9..79eab50610 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -18,7 +18,6 @@ module BalanceCheckMod use EnergyFluxType , only : energyflux_type use SolarAbsorbedType , only : solarabs_type use SoilHydrologyType , only : soilhydrology_type - use SurfaceAlbedoType , only : surfalb_type use WaterStateType , only : waterstate_type use LakestateType , only : lakestate_type use WaterDiagnosticBulkType, only : waterdiagnosticbulk_type @@ -399,6 +398,7 @@ subroutine BalanceCheck( bounds, & ! ! !USES: use clm_varcon , only : spval + use clm_varctl , only : use_soil_moisture_streams use clm_time_manager , only : get_step_size_real, get_nstep use clm_time_manager , only : get_nstep_since_startup_or_lastDA_restart_or_pause use CanopyStateType , only : canopystate_type @@ -724,6 +724,7 @@ subroutine BalanceCheck( bounds, & ' local indexg= ',indexg,& ' errh2o_grc= ',errh2o_grc(indexg) if (errh2o_max_val > error_thresh .and. DAnstep > skip_steps .and. & + .not. use_soil_moisture_streams .and. & .not. get_for_testing_zero_dynbal_fluxes()) then write(iulog,*)'CTSM is stopping because errh2o > ', error_thresh, ' mm' From 09e35531f7422229a332fe7e972801981339e9f4 Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Thu, 18 Feb 2021 18:02:52 -0700 Subject: [PATCH 18/19] Revisions required for I*G tests and *clm-prescribed tests to PASS 1) Addressed the former with some code reorg. At the crux of it though was the need to pass subtract_dynbal_baselines = .true. and add_lake_water_and_subtract_dynbal_baselines = .true. instead of .false. when calling ComputeWaterMassNonLake and ComputeWaterMassLake at the beginning and end of the time step. 2) Addressed the *clm-prescribed test failure by bypassing the grc-level water balance check when use_soil_moisture_streams = .false. --- src/biogeophys/BalanceCheckMod.F90 | 183 ++++++++++++++++------------- src/main/clm_driver.F90 | 12 +- src/main/lnd2atmMod.F90 | 4 +- 3 files changed, 113 insertions(+), 86 deletions(-) diff --git a/src/biogeophys/BalanceCheckMod.F90 b/src/biogeophys/BalanceCheckMod.F90 index 79eab50610..e1d7497198 100644 --- a/src/biogeophys/BalanceCheckMod.F90 +++ b/src/biogeophys/BalanceCheckMod.F90 @@ -44,7 +44,7 @@ module BalanceCheckMod ! !PUBLIC MEMBER FUNCTIONS: public :: BalanceCheckInit ! Initialization of Water and energy balance check - public :: BeginWaterGridcellBalance ! Initialize grid cell-level water balance check + public :: WaterGridcellBalance ! Grid cell-level water balance check public :: BeginWaterColumnBalance ! Initialize column-level water balance check public :: BalanceCheck ! Water & energy balance checks public :: GetBalanceCheckSkipSteps ! Get the number of steps to skip for the balance check @@ -56,7 +56,7 @@ module BalanceCheckMod ! ! !PRIVATE MEMBER FUNCTIONS: - private :: BeginWaterGridcellBalanceSingle ! Initialize grid cell-level water balance check for bulk or a single tracer + private :: WaterGridcellBalanceSingle ! Grid cell-level water balance check for bulk or a single tracer private :: BeginWaterColumnBalanceSingle ! Initialize column-level water balance check for bulk or a single tracer character(len=*), parameter, private :: sourcefile = & @@ -123,13 +123,12 @@ end function GetBalanceCheckSkipSteps !----------------------------------------------------------------------- !----------------------------------------------------------------------- - subroutine BeginWaterGridcellBalance(bounds, & + subroutine WaterGridcellBalance(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - water_inst, lakestate_inst) + water_inst, lakestate_inst, use_aquifer_layer, flag) ! ! !DESCRIPTION: - ! Initialize grid cell-level water balance at beginning of time step - ! for bulk water and each water tracer + ! Grid cell-level water balance for bulk water and each water tracer ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds @@ -139,25 +138,28 @@ subroutine BeginWaterGridcellBalance(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(water_type) , intent(inout) :: water_inst type(lakestate_type) , intent(in) :: lakestate_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run + character(len=5) , intent(in) :: flag ! specifies begwb or endwb ! ! !LOCAL VARIABLES: integer :: i - character(len=*), parameter :: subname = 'BeginWaterGridcellBalance' + character(len=*), parameter :: subname = 'WaterGridcellBalance' !----------------------------------------------------------------------- do i = water_inst%bulk_and_tracers_beg, water_inst%bulk_and_tracers_end - call BeginWaterGridcellBalanceSingle(bounds, & - num_nolakec, filter_nolakec, & - num_lakec, filter_lakec, & + ! Obtain begwb_grc + call WaterGridcellBalanceSingle(bounds, & + num_nolakec, filter_nolakec, num_lakec, filter_lakec, & lakestate_inst, & water_inst%bulk_and_tracers(i)%waterstate_inst, & water_inst%bulk_and_tracers(i)%waterdiagnostic_inst, & water_inst%bulk_and_tracers(i)%waterbalance_inst, & - water_inst%bulk_and_tracers(i)%waterflux_inst) + water_inst%bulk_and_tracers(i)%waterflux_inst, & + use_aquifer_layer = use_aquifer_layer, flag = flag) end do - end subroutine BeginWaterGridcellBalance + end subroutine WaterGridcellBalance !----------------------------------------------------------------------- subroutine BeginWaterColumnBalance(bounds, & @@ -201,14 +203,14 @@ subroutine BeginWaterColumnBalance(bounds, & end subroutine BeginWaterColumnBalance !----------------------------------------------------------------------- - subroutine BeginWaterGridcellBalanceSingle(bounds, & + subroutine WaterGridcellBalanceSingle(bounds, & num_nolakec, filter_nolakec, num_lakec, filter_lakec, & - lakestate_inst, waterstate_inst, & - waterdiagnostic_inst, waterbalance_inst, waterflux_inst) + lakestate_inst, waterstate_inst, waterdiagnostic_inst, & + waterbalance_inst, waterflux_inst, use_aquifer_layer, flag) ! ! !DESCRIPTION: - ! Initialize grid cell-level water balance at beginning of time step - ! for bulk or a single tracer + ! Grid cell-level water balance for bulk or a single tracer + ! at beginning or end of time step as specified by "flag" ! ! !USES: use subgridAveMod, only: c2g @@ -221,20 +223,28 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & integer , intent(in) :: filter_lakec(:) ! column filter for lake points type(lakestate_type) , intent(in) :: lakestate_inst class(waterstate_type) , intent(inout) :: waterstate_inst - class(waterflux_type) , intent(inout) :: waterflux_inst class(waterdiagnostic_type), intent(in) :: waterdiagnostic_inst class(waterbalance_type) , intent(inout) :: waterbalance_inst + class(waterflux_type) , intent(inout) :: waterflux_inst + logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run + character(len=5) , intent(in) :: flag ! specifies begwb or endwb ! ! !LOCAL VARIABLES: integer :: g ! indices integer :: begc, endc, begg, endg ! bounds - real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at beginning of time step - real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc ice dynamic land cover change conversion runoff flux at beginning of time step + real(r8) :: wb_col(bounds%begc:bounds%endc) ! temporary column-level water mass + real(r8) :: wb_grc(bounds%begg:bounds%endg) ! temporary grid cell-level water mass + real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux + real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc ice dynamic land cover change conversion runoff flux + real(r8) :: wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) ! grc mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) + + character(len=*), parameter :: subname = 'WaterGridcellBalanceSingle' !----------------------------------------------------------------------- associate( & - begwb_col => waterbalance_inst%begwb_col, & ! Output: [real(r8) (:) ] column-level water mass begining of the time step - begwb_grc => waterbalance_inst%begwb_grc & ! Output: [real(r8) (:) ] grid cell-level water mass begining of the time step + begwb_grc => waterbalance_inst%begwb_grc, & ! Output: [real(r8) (:)] grid cell-level water mass begining of the time step + endwb_grc => waterbalance_inst%endwb_grc, & ! Output: [real(r8) (:)] grid cell-level water mass end of the time step + wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col & ! Input: [real(r8) (:)] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) ) begc = bounds%begc @@ -242,27 +252,40 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & begg = bounds%begg endg = bounds%endg - ! NOTES subroutines Compute*Mass* are in TotalWaterAndHeatMod.F90 - ! endwb is calculated in HydrologyDrainageMod & LakeHydrologyMod call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, & waterstate_inst, waterdiagnostic_inst, & - subtract_dynbal_baselines = .false., & - water_mass = begwb_col(begc:endc)) + subtract_dynbal_baselines = .true., & + water_mass = wb_col(begc:endc)) call ComputeWaterMassLake(bounds, num_lakec, filter_lakec, & waterstate_inst, lakestate_inst, & - add_lake_water_and_subtract_dynbal_baselines = .false., & - water_mass = begwb_col(begc:endc)) + add_lake_water_and_subtract_dynbal_baselines = .true., & + water_mass = wb_col(begc:endc)) - call c2g(bounds, begwb_col(begc:endc), begwb_grc(begg:endg), & + call c2g(bounds, wb_col(begc:endc), wb_grc(begg:endg), & c2l_scale_type='urbanf', l2g_scale_type='unity') - call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_beg( & - bounds, & - qflx_liq_dynbal_left_to_dribble(begg:endg)) - call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_beg( & - bounds, & - qflx_ice_dynbal_left_to_dribble(begg:endg)) + ! Call the beginning or ending version of the subroutine according + ! to flag value + if (flag == 'begwb') then + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_beg( & + bounds, & + qflx_liq_dynbal_left_to_dribble(begg:endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_beg( & + bounds, & + qflx_ice_dynbal_left_to_dribble(begg:endg)) + else if (flag == 'endwb') then + call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end( & + bounds, & + qflx_liq_dynbal_left_to_dribble(begg:endg)) + call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_end( & + bounds, & + qflx_ice_dynbal_left_to_dribble(begg:endg)) + else + write(iulog,*) 'Unknown flag passed into this subroutine.' + write(iulog,*) 'Expecting either begwb or endwb.' + call endrun(msg=errmsg(sourcefile, __LINE__)) + end if ! These dynbal dribblers store the delta state, (end - beg). Thus, the ! amount dribbled out is the negative of the amount stored in the @@ -271,16 +294,42 @@ subroutine BeginWaterGridcellBalanceSingle(bounds, & ! This sign convention is opposite to the convention chosen for the ! respective dribble terms used in the carbon balance. At some point ! it may be worth making the two conventions consistent. - ! Bill Sacks states: I think the convention used for the water and - ! energy dribblers is counter-intuitive. do g = begg, endg - begwb_grc(g) = begwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - - qflx_ice_dynbal_left_to_dribble(g) + wb_grc(g) = wb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & + - qflx_ice_dynbal_left_to_dribble(g) end do + ! Map wb_grc to beginning/ending water balance according to flag + if (flag == 'begwb') then + do g = begg, endg + begwb_grc(g) = wb_grc(g) + end do + else if (flag == 'endwb') then + ! endwb_grc requires one more step first + if (use_aquifer_layer) then + ! wa_reset_nonconservation_gain may be non-zero only when + ! use_aquifer_layer is true. We do this c2g call only when needed + ! to avoid unnecessary calculations; by adding this term only when + ! use_aquifer_layer is true, we effectively let the balance checks + ! ensure that this term is zero when use_aquifer_layer is false, + ! as it should be. + ! The _col term was determined in BeginWaterColumnBalanceSingle + ! after any dynamic landuse adjustments. + call c2g( bounds, & + wa_reset_nonconservation_gain_col(begc:endc), & + wa_reset_nonconservation_gain_grc(begg:endg), & + c2l_scale_type='urbanf', l2g_scale_type='unity' ) + else + wa_reset_nonconservation_gain_grc(begg:endg) = 0._r8 + end if + do g = begg, endg + endwb_grc(g) = wb_grc(g) - wa_reset_nonconservation_gain_grc(g) + end do + end if + end associate - end subroutine BeginWaterGridcellBalanceSingle + end subroutine WaterGridcellBalanceSingle !----------------------------------------------------------------------- subroutine BeginWaterColumnBalanceSingle(bounds, & @@ -377,6 +426,8 @@ end subroutine BeginWaterColumnBalanceSingle !----------------------------------------------------------------------- subroutine BalanceCheck( bounds, & num_allc, filter_allc, & + num_nolakec, filter_nolakec, num_lakec, filter_lakec, & + water_inst, lakestate_inst, & atm2lnd_inst, solarabs_inst, waterflux_inst, waterstate_inst, & waterdiagnosticbulk_inst, waterbalance_inst, wateratm2lnd_inst, & waterlnd2atm_inst, surfalb_inst, energyflux_inst, canopystate_inst, & @@ -410,6 +461,12 @@ subroutine BalanceCheck( bounds, & type(bounds_type) , intent(in) :: bounds integer , intent(in) :: num_allc ! number of columns in allc filter integer , intent(in) :: filter_allc(:) ! filter for all columns + integer , intent(in) :: num_nolakec ! number of column non-lake points in column filter + integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points + integer , intent(in) :: num_lakec ! number of column lake points in column filter + integer , intent(in) :: filter_lakec(:) ! column filter for lake points + type(water_type) , intent(inout) :: water_inst + type(lakestate_type) , intent(in) :: lakestate_inst type(atm2lnd_type) , intent(in) :: atm2lnd_inst type(solarabs_type) , intent(in) :: solarabs_inst class(waterflux_type) , intent(in) :: waterflux_inst @@ -436,9 +493,6 @@ subroutine BalanceCheck( bounds, & real(r8) :: qflx_glcice_dyn_water_flux_grc(bounds%begg:bounds%endg) ! grid cell-level water flux needed for balance check due to glc_dyn_runoff_routing [mm H2O/s] (positive means addition of water to the system) real(r8) :: qflx_snwcp_discarded_liq_grc(bounds%begg:bounds%endg) ! grid cell-level excess liquid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] real(r8) :: qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg) ! grid cell-level excess solid h2o due to snow capping, which we simply discard in order to reset the snow pack [mm H2O /s] - real(r8) :: qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step - real(r8) :: qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg) ! grc liq dynamic land cover change conversion runoff flux at end of time step - real(r8) :: wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) ! grc mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) real(r8) :: errh2o_max_val ! Maximum value of error in water conservation error over all columns [mm H2O] real(r8) :: errh2osno_max_val ! Maximum value of error in h2osno conservation error over all columns [kg m-2] @@ -491,7 +545,6 @@ subroutine BalanceCheck( bounds, & qflx_h2osfc_to_ice => waterflux_inst%qflx_h2osfc_to_ice_col , & ! Input: [real(r8) (:) ] conversion of h2osfc to ice qflx_drain_perched_col => waterflux_inst%qflx_drain_perched_col , & ! Input: [real(r8) (:) ] column level sub-surface runoff (mm H2O /s) qflx_drain_perched_grc => waterlnd2atm_inst%qflx_rofliq_drain_perched_grc, & ! Input: [real(r8) (:)] grid cell-level sub-surface runoff (mm H2O /s) - wa_reset_nonconservation_gain_col => waterbalance_inst%wa_reset_nonconservation_gain_col, & ! Output: [real(r8) (:) ] col mass gained from resetting water in the unconfined aquifer, wa_col (negative indicates mass lost) (mm) qflx_flood_col => waterflux_inst%qflx_floodc_col , & ! Input: [real(r8) (:) ] column level total runoff due to flooding forc_flood_grc => wateratm2lnd_inst%forc_flood_grc , & ! Input: [real(r8) (:) ] grid cell-level total grid cell-level runoff from river model qflx_snow_drain => waterflux_inst%qflx_snow_drain_col , & ! Input: [real(r8) (:) ] drainage from snow pack @@ -660,44 +713,14 @@ subroutine BalanceCheck( bounds, & qflx_snwcp_discarded_ice_col(bounds%begc:bounds%endc), & qflx_snwcp_discarded_ice_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) - if (use_aquifer_layer) then - ! wa_reset_nonconservation_gain may be non-zero only when - ! use_aquifer_layer is true. We do this c2g call only when needed - ! to avoid unnecessary calculations; by adding this term only when - ! use_aquifer_layer is true, we effectively let the balance checks - ! ensure that this term is zero when use_aquifer_layer is false, - ! as it should be. - ! The _col term was determined in BeginWaterColumnBalanceSingle - ! after any dynamic landuse adjustments. - call c2g( bounds, & - wa_reset_nonconservation_gain_col(bounds%begc:bounds%endc), & - wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg), & - c2l_scale_type='urbanf', l2g_scale_type='unity' ) - else - wa_reset_nonconservation_gain_grc(bounds%begg:bounds%endg) = 0._r8 - end if - call waterflux_inst%qflx_liq_dynbal_dribbler%get_amount_left_to_dribble_end( & - bounds, & - qflx_liq_dynbal_left_to_dribble(bounds%begg:bounds%endg)) - call waterflux_inst%qflx_ice_dynbal_dribbler%get_amount_left_to_dribble_end( & - bounds, & - qflx_ice_dynbal_left_to_dribble(bounds%begg:bounds%endg)) - - ! These dynbal dribblers store the delta state, (end - beg). Thus, the - ! amount dribbled out is the negative of the amount stored in the - ! dribblers. Therefore, conservation requires us to subtract the amount - ! remaining to dribble. - ! This sign convention is opposite to the convention chosen for the - ! respective dribble terms used in the carbon balance. At some point - ! it may be worth making the two conventions consistent. - ! Bill Sacks states: I think the convention used for the water and - ! energy dribblers is counter-intuitive. - do g = bounds%begg, bounds%endg - endwb_grc(g) = endwb_grc(g) - qflx_liq_dynbal_left_to_dribble(g) & - - qflx_ice_dynbal_left_to_dribble(g) & - - wa_reset_nonconservation_gain_grc(g) + ! Obtain endwb_grc + call WaterGridcellBalance(bounds, & + num_nolakec, filter_nolakec, num_lakec, filter_lakec, & + water_inst, lakestate_inst, & + use_aquifer_layer = use_aquifer_layer, flag = 'endwb') + do g = bounds%begg, bounds%endg errh2o_grc(g) = endwb_grc(g) - begwb_grc(g) & - (forc_rain_grc(g) & + forc_snow_grc(g) & diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90 index 1df3a3640c..1122455abf 100644 --- a/src/main/clm_driver.F90 +++ b/src/main/clm_driver.F90 @@ -26,7 +26,7 @@ module clm_driver use abortutils , only : endrun ! use dynSubgridDriverMod , only : dynSubgrid_driver, dynSubgrid_wrapup_weight_changes - use BalanceCheckMod , only : BeginWaterGridcellBalance, BeginWaterColumnBalance, BalanceCheck + use BalanceCheckMod , only : WaterGridcellBalance, BeginWaterColumnBalance, BalanceCheck ! use BiogeophysPreFluxCalcsMod , only : BiogeophysPreFluxCalcs use SurfaceHumidityMod , only : CalculateSurfaceHumidity @@ -328,10 +328,11 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call t_stopf('begcnbal_grc') call t_startf('begwbal') - call BeginWaterGridcellBalance(bounds_clump, & + call WaterGridcellBalance(bounds_clump, & filter(nc)%num_nolakec, filter(nc)%nolakec, & filter(nc)%num_lakec, filter(nc)%lakec, & - water_inst, lakestate_inst) + water_inst, lakestate_inst, & + use_aquifer_layer = use_aquifer_layer(), flag = 'begwb') call t_stopf('begwbal') end do !$OMP END PARALLEL DO @@ -370,7 +371,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro ! conserved when weights change (instead the difference is put in the grid cell-level ! terms, qflx_liq_dynbal, etc.). Grid cell-level balance ! checks ensure that the grid cell-level water is conserved by considering - ! qflx_liq_dynbal and calling BeginWaterGridcellBalance + ! qflx_liq_dynbal and calling WaterGridcellBalance ! before the weight updates. ! ! For carbon & nitrogen: This needs to be done after dynSubgrid_driver, because the @@ -1284,6 +1285,9 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro call get_clump_bounds(nc, bounds_clump) call BalanceCheck(bounds_clump, & filter(nc)%num_allc, filter(nc)%allc, & + filter(nc)%num_nolakec, filter(nc)%nolakec, & + filter(nc)%num_lakec, filter(nc)%lakec, & + water_inst, lakestate_inst, & atm2lnd_inst, solarabs_inst, water_inst%waterfluxbulk_inst, & water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, & water_inst%waterbalancebulk_inst, water_inst%wateratm2lndbulk_inst, & diff --git a/src/main/lnd2atmMod.F90 b/src/main/lnd2atmMod.F90 index ed9d44f36f..b05013ea09 100644 --- a/src/main/lnd2atmMod.F90 +++ b/src/main/lnd2atmMod.F90 @@ -411,11 +411,11 @@ subroutine lnd2atm(bounds, & call c2g( bounds, & water_inst%waterbalancebulk_inst%endwb_col(bounds%begc:bounds%endc), & - water_inst%waterbalancebulk_inst%endwb_grc(bounds%begg:bounds%endg), & + water_inst%waterdiagnosticbulk_inst%tws_grc(bounds%begg:bounds%endg), & c2l_scale_type= 'urbanf', l2g_scale_type='unity' ) do g = bounds%begg, bounds%endg water_inst%waterdiagnosticbulk_inst%tws_grc(g) = & - water_inst%waterbalancebulk_inst%endwb_grc(g) + & + water_inst%waterdiagnosticbulk_inst%tws_grc(g) + & water_inst%wateratm2lndbulk_inst%volr_grc(g) / grc%area(g) * 1.e-3_r8 enddo From 5a0ba100a094c7152fa0f7247d44f2e88a42823f Mon Sep 17 00:00:00 2001 From: Samuel Levis Date: Sat, 20 Feb 2021 14:51:46 -0700 Subject: [PATCH 19/19] Small updates to ChangeLog and ChangeSum --- doc/ChangeLog | 11 ++++------- doc/ChangeSum | 2 +- 2 files changed, 5 insertions(+), 8 deletions(-) diff --git a/doc/ChangeLog b/doc/ChangeLog index a5e4b0c493..a25df070f3 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,7 +1,7 @@ =============================================================== Tag name: ctsm5.1.dev024 Originator(s): slevis (Samuel Levis,303-665-1310) -Date: Wed Feb 11 10:35:33 MST 2021 +Date: Sat Feb 20 14:42:33 MST 2021 One-line Summary: Grid cell-level error check for H2O Purpose of changes @@ -84,8 +84,8 @@ CTSM testing: regular tests (aux_clm): - cheyenne ---- PEND (expect OK) - izumi ------- PEND (expect OK) + cheyenne ---- OK (comparisons to baseline fail as expected) + izumi ------- OK (comparisons to baseline fail as expected) If the tag used for baseline comparisons was NOT the previous tag, note that here: @@ -99,13 +99,10 @@ Changes answers relative to baseline: YES - what code configurations: ALL - what platforms/compilers: ALL - nature of change: ROUNDOFF - Specific example from running the single point test - ERI_D_Ld9.1x1_camdenNJ.I2000Clm50BgcCruRs.cheyenne_intel.clm-default: - RMS ERRH2O 6.0280E-21 NORMALIZED 7.6050E-06 Explanation: Moving call BalanceCheck to after the call lnd2glc in subroutine clm_drv causes a change in order of operations that leads to - the above change in ERRH2O. + roundoff change in ERRH2O. Confirmed by running ./summarize_cprnc_diffs Detailed list of changes diff --git a/doc/ChangeSum b/doc/ChangeSum index 3337e6e449..e179f3cd2d 100644 --- a/doc/ChangeSum +++ b/doc/ChangeSum @@ -1,6 +1,6 @@ Tag Who Date Summary ============================================================================================================================ - ctsm5.1.dev024 slevis 02/11/2021 Grid cell-level error check for H2O + ctsm5.1.dev024 slevis 02/20/2021 Grid cell-level error check for H2O ctsm5.1.dev023 erik 02/11/2021 Calculate leaf biomass for non-woody PFTS, and a few other small answer changes ctsm5.1.dev022 glemieux 02/05/2021 Merge fates_main_api into ctsm master ctsm5.1.dev021 erik 01/12/2021 Add option for biomass heat storage (BHS) to clm5_1 physics