Skip to content

Commit

Permalink
Change wa_col cold start initialization when use_aquifer_layer is false
Browse files Browse the repository at this point in the history
In the groundwater_irrigation branch
(ESCOMP#523), @swensosc has stopped
resetting wa_col each time step if use_aquifer_layer is false. However,
this leads to having a substantially different value of wa_col when
use_aquifer_layer is false: previously, it was reset to
aquifer_water_baseline each time step, but with the changes from
@swensosc, it stays close to its initial values, which have been 4000 in
most places. This commit changes the initial values to match the value
it was being reset to, so it simply starts at aquifer_water_baseline -
so the every-time-step-resetting to aquifer_water_baseline can be
removed without massively changing the value of wa_col.

In addition to changing the cold start initialization of wa_col, we are
also changing the cold start initialization of zwt in this case where
use_aquifer_layer is false: The old initialization of zwt wouldn't work
as intended now that we have changed the initial value of wa_col;
@swensosc suggested this new initialization method instead. This
initialization to zi(c,nbedrock(c)) is correct if there are no saturated
layers, and close enough for a decent cold start even if there are.
  • Loading branch information
billsacks committed Nov 27, 2018
1 parent 527f326 commit c600499
Show file tree
Hide file tree
Showing 5 changed files with 109 additions and 41 deletions.
41 changes: 34 additions & 7 deletions src/biogeophys/SoilHydrologyInitTimeConstMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,8 @@ module SoilHydrologyInitTimeConstMod
contains

!-----------------------------------------------------------------------
subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, waterstatebulk_inst)
subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, waterstatebulk_inst, &
use_aquifer_layer)
!
! !USES:
use shr_const_mod , only : shr_const_pi
Expand All @@ -53,6 +54,7 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, waterstatebulk
type(bounds_type) , intent(in) :: bounds
type(soilhydrology_type) , intent(inout) :: soilhydrology_inst
type(waterstatebulk_type) , intent(inout) :: waterstatebulk_inst
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run
!
! !LOCAL VARIABLES:
integer :: p,c,j,l,g,lev,nlevs
Expand Down Expand Up @@ -91,19 +93,44 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, waterstatebulk
if (.not. lun%lakpoi(l)) then !not lake
if (lun%urbpoi(l)) then
if (col%itype(c) == icol_road_perv) then
! Note that the following hard-coded constants (on the next line)
! seem implicitly related to aquifer_water_baseline
soilhydrology_inst%zwt_col(c) = (25._r8 + col%zi(c,nlevsoi)) - waterstatebulk_inst%wa_col(c)/0.2_r8 /1000._r8 ! One meter below soil column
if (use_aquifer_layer) then
! NOTE(wjs, 2018-11-27) There is no fundamental reason why zwt should
! be initialized differently based on use_aquifer_layer, but we (Bill
! Sacks and Sean Swenson) are changing the cold start initialization of
! wa_col when use_aquifer_layer is .false., and so need to come up with
! a different cold start initialization of zwt in that case, but we
! don't want to risk messing up the use_aquifer_layer = .true. case,
! so we're keeping that as it was before.

! Note that the following hard-coded constants (on the next line)
! seem implicitly related to the initial value of wa_col
soilhydrology_inst%zwt_col(c) = (25._r8 + col%zi(c,nlevsoi)) - waterstatebulk_inst%wa_col(c)/0.2_r8 /1000._r8 ! One meter below soil column
else
soilhydrology_inst%zwt_col(c) = col%zi(c,col%nbedrock(c))
end if
else
soilhydrology_inst%zwt_col(c) = spval
end if
! initialize frost_table, zwt_perched
soilhydrology_inst%zwt_perched_col(c) = spval
soilhydrology_inst%frost_table_col(c) = spval
else
! Note that the following hard-coded constants (on the next line) seem
! implicitly related to aquifer_water_baseline
soilhydrology_inst%zwt_col(c) = (25._r8 + col%zi(c,nlevsoi)) - waterstatebulk_inst%wa_col(c)/0.2_r8 /1000._r8 ! One meter below soil column
if (use_aquifer_layer) then
! NOTE(wjs, 2018-11-27) There is no fundamental reason why zwt should
! be initialized differently based on use_aquifer_layer, but we (Bill
! Sacks and Sean Swenson) are changing the cold start initialization of
! wa_col when use_aquifer_layer is .false., and so need to come up with
! a different cold start initialization of zwt in that case, but we
! don't want to risk messing up the use_aquifer_layer = .true. case,
! so we're keeping that as it was before.

! Note that the following hard-coded constants (on the next line) seem
! implicitly related to the initial value of wa_col
soilhydrology_inst%zwt_col(c) = (25._r8 + col%zi(c,nlevsoi)) - waterstatebulk_inst%wa_col(c)/0.2_r8 /1000._r8
else
soilhydrology_inst%zwt_col(c) = col%zi(c,col%nbedrock(c))
end if

! initialize frost_table, zwt_perched to bottom of soil column
soilhydrology_inst%zwt_perched_col(c) = col%zi(c,nlevsoi)
soilhydrology_inst%frost_table_col(c) = col%zi(c,nlevsoi)
Expand Down
12 changes: 9 additions & 3 deletions src/biogeophys/WaterStateBulkType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ module WaterStateBulkType

!------------------------------------------------------------------------
subroutine InitBulk(this, bounds, info, vars, &
h2osno_input_col, watsat_col, t_soisno_col)
h2osno_input_col, watsat_col, t_soisno_col, use_aquifer_layer)

class(waterstatebulk_type), intent(inout) :: this
type(bounds_type) , intent(in) :: bounds
Expand All @@ -59,10 +59,16 @@ subroutine InitBulk(this, bounds, info, vars, &
real(r8) , intent(in) :: h2osno_input_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run


call this%Init(bounds, info, vars, &
h2osno_input_col, watsat_col, t_soisno_col)
call this%Init(bounds = bounds, &
info = info, &
tracer_vars = vars, &
h2osno_input_col = h2osno_input_col, &
watsat_col = watsat_col, &
t_soisno_col = t_soisno_col, &
use_aquifer_layer = use_aquifer_layer)

call this%InitBulkAllocate(bounds)

Expand Down
53 changes: 33 additions & 20 deletions src/biogeophys/WaterStateType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -60,7 +60,7 @@ module WaterStateType

!------------------------------------------------------------------------
subroutine Init(this, bounds, info, tracer_vars, &
h2osno_input_col, watsat_col, t_soisno_col)
h2osno_input_col, watsat_col, t_soisno_col, use_aquifer_layer)

class(waterstate_type), intent(inout) :: this
type(bounds_type) , intent(in) :: bounds
Expand All @@ -69,15 +69,19 @@ subroutine Init(this, bounds, info, tracer_vars, &
real(r8) , intent(in) :: h2osno_input_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run

this%info => info

call this%InitAllocate(bounds, tracer_vars)

call this%InitHistory(bounds)

call this%InitCold(bounds, &
h2osno_input_col, watsat_col, t_soisno_col)
call this%InitCold(bounds = bounds, &
h2osno_input_col = h2osno_input_col, &
watsat_col = watsat_col, &
t_soisno_col = t_soisno_col, &
use_aquifer_layer = use_aquifer_layer)

end subroutine Init

Expand Down Expand Up @@ -279,7 +283,7 @@ end subroutine InitHistory

!-----------------------------------------------------------------------
subroutine InitCold(this, bounds, &
h2osno_input_col, watsat_col, t_soisno_col)
h2osno_input_col, watsat_col, t_soisno_col, use_aquifer_layer)
!
! !DESCRIPTION:
! Initialize time constant variables and cold start conditions
Expand All @@ -298,6 +302,7 @@ subroutine InitCold(this, bounds, &
real(r8) , intent(in) :: h2osno_input_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run
!
! !LOCAL VARIABLES:
integer :: p,c,j,l,g,lev,nlevs
Expand Down Expand Up @@ -455,25 +460,33 @@ subroutine InitCold(this, bounds, &
end do


this%wa_col(bounds%begc:bounds%endc) = aquifer_water_baseline
do c = bounds%begc,bounds%endc
l = col%landunit(c)
if (.not. lun%lakpoi(l)) then !not lake
if (lun%urbpoi(l)) then
if (col%itype(c) == icol_road_perv) then
! Note that the following hard-coded constant (on the next line)
! seems implicitly related to aquifer_water_baseline
this%wa_col(c) = 4800._r8
this%wa_col(bounds%begc:bounds%endc) = aquifer_water_baseline
if (use_aquifer_layer) then
! NOTE(wjs, 2018-11-27) There is no fundamental reason why wa_col should be
! initialized differently based on use_aquifer_layer, but we (Bill Sacks and Sean
! Swenson) want to change the cold start initialization of wa_col to be
! aquifer_water_baseline everywhere for use_aquifer_layer .false., and we aren't
! sure of the implications of this change for use_aquifer_layer .true., so are
! maintaining the old cold start initialization in the latter case.
do c = bounds%begc,bounds%endc
l = col%landunit(c)
if (.not. lun%lakpoi(l)) then !not lake
if (lun%urbpoi(l)) then
if (col%itype(c) == icol_road_perv) then
! Note that the following hard-coded constant (on the next line)
! seems implicitly related to aquifer_water_baseline
this%wa_col(c) = 4800._r8
else
this%wa_col(c) = spval
end if
else
this%wa_col(c) = spval
! Note that the following hard-coded constant (on the next line) seems
! implicitly related to aquifer_water_baseline
this%wa_col(c) = 4000._r8
end if
else
! Note that the following hard-coded constant (on the next line) seems
! implicitly related to aquifer_water_baseline
this%wa_col(c) = 4000._r8
end if
end if
end do
end do
end if

end associate

Expand Down
37 changes: 28 additions & 9 deletions src/biogeophys/WaterType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -167,7 +167,7 @@ end function water_params_constructor

!-----------------------------------------------------------------------
subroutine Init(this, bounds, NLFilename, &
h2osno_col, snow_depth_col, watsat_col, t_soisno_col)
h2osno_col, snow_depth_col, watsat_col, t_soisno_col, use_aquifer_layer)
!
! !DESCRIPTION:
! Initialize all water variables
Expand All @@ -180,21 +180,26 @@ subroutine Init(this, bounds, NLFilename, &
real(r8) , intent(in) :: snow_depth_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run
!
! !LOCAL VARIABLES:

character(len=*), parameter :: subname = 'Init'
!-----------------------------------------------------------------------

call this%ReadNamelist(NLFilename)
call this%DoInit(bounds, &
h2osno_col, snow_depth_col, watsat_col, t_soisno_col)
call this%DoInit(bounds = bounds, &
h2osno_col = h2osno_col, &
snow_depth_col = snow_depth_col, &
watsat_col = watsat_col, &
t_soisno_col = t_soisno_col, &
use_aquifer_layer = use_aquifer_layer)

end subroutine Init

!-----------------------------------------------------------------------
subroutine InitForTesting(this, bounds, params, &
h2osno_col, snow_depth_col, watsat_col, t_soisno_col)
h2osno_col, snow_depth_col, watsat_col, t_soisno_col, use_aquifer_layer)
!
! !DESCRIPTION:
! Version of Init routine just for unit tests
Expand All @@ -209,21 +214,32 @@ subroutine InitForTesting(this, bounds, params, &
real(r8) , intent(in) :: snow_depth_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in), optional :: use_aquifer_layer ! whether an aquifer layer is used in this run (false by default)
!
! !LOCAL VARIABLES:
logical :: l_use_aquifer_layer

character(len=*), parameter :: subname = 'InitForTesting'
!-----------------------------------------------------------------------

l_use_aquifer_layer = .false.
if (present(use_aquifer_layer)) then
l_use_aquifer_layer = use_aquifer_layer
end if

this%params = params
call this%DoInit(bounds, &
h2osno_col, snow_depth_col, watsat_col, t_soisno_col)
call this%DoInit(bounds = bounds, &
h2osno_col = h2osno_col, &
snow_depth_col = snow_depth_col, &
watsat_col = watsat_col, &
t_soisno_col = t_soisno_col, &
use_aquifer_layer = l_use_aquifer_layer)

end subroutine InitForTesting

!-----------------------------------------------------------------------
subroutine DoInit(this, bounds, &
h2osno_col, snow_depth_col, watsat_col, t_soisno_col)
h2osno_col, snow_depth_col, watsat_col, t_soisno_col, use_aquifer_layer)
!
! !DESCRIPTION:
! Actually do the initialization (shared between main Init routine and InitForTesting)
Expand All @@ -237,6 +253,7 @@ subroutine DoInit(this, bounds, &
real(r8) , intent(in) :: snow_depth_col(bounds%begc:)
real(r8) , intent(in) :: watsat_col(bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
real(r8) , intent(in) :: t_soisno_col(bounds%begc:, -nlevsno+1:) ! col soil temperature (Kelvin)
logical , intent(in) :: use_aquifer_layer ! whether an aquifer layer is used in this run
!
! !LOCAL VARIABLES:
integer :: begc, endc
Expand All @@ -261,7 +278,8 @@ subroutine DoInit(this, bounds, &
this%bulk_vars, &
h2osno_input_col = h2osno_col(begc:endc), &
watsat_col = watsat_col(begc:endc, 1:), &
t_soisno_col = t_soisno_col(begc:endc, -nlevsno+1:) )
t_soisno_col = t_soisno_col(begc:endc, -nlevsno+1:), &
use_aquifer_layer = use_aquifer_layer)

call this%waterdiagnosticbulk_inst%InitBulk(bounds, &
this%bulk_info, &
Expand Down Expand Up @@ -308,7 +326,8 @@ subroutine DoInit(this, bounds, &
this%tracer_vars(i), &
h2osno_input_col = h2osno_col(begc:endc), &
watsat_col = watsat_col(begc:endc, 1:), &
t_soisno_col = t_soisno_col(begc:endc, -nlevsno+1:) )
t_soisno_col = t_soisno_col(begc:endc, -nlevsno+1:), &
use_aquifer_layer = use_aquifer_layer)

call this%waterdiagnostic_tracer_inst(i)%Init(bounds, &
this%tracer_info(i)%info, &
Expand Down
7 changes: 5 additions & 2 deletions src/main/clm_instMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -82,6 +82,7 @@ module clm_instMod
use SurfaceAlbedoMod , only : SurfaceAlbedoInitTimeConst
use LakeCon , only : LakeConInit
use SoilBiogeochemPrecisionControlMod, only: SoilBiogeochemPrecisionControlInit
use SoilWaterMovementMod , only : use_aquifer_layer
!
implicit none
public ! By default everything is public
Expand Down Expand Up @@ -280,7 +281,8 @@ subroutine clm_instInit(bounds)
h2osno_col = h2osno_col(begc:endc), &
snow_depth_col = snow_depth_col(begc:endc), &
watsat_col = soilstate_inst%watsat_col(begc:endc, 1:), &
t_soisno_col = temperature_inst%t_soisno_col(begc:endc, -nlevsno+1:))
t_soisno_col = temperature_inst%t_soisno_col(begc:endc, -nlevsno+1:), &
use_aquifer_layer = use_aquifer_layer())

call glacier_smb_inst%Init(bounds)

Expand All @@ -303,7 +305,8 @@ subroutine clm_instInit(bounds)
call photosyns_inst%Init(bounds)

call soilhydrology_inst%Init(bounds, nlfilename)
call SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, water_inst%waterstatebulk_inst) ! sets time constant properties
call SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, water_inst%waterstatebulk_inst, &
use_aquifer_layer = use_aquifer_layer())

call saturated_excess_runoff_inst%Init(bounds)
call infiltration_excess_runoff_inst%Init(bounds)
Expand Down

0 comments on commit c600499

Please sign in to comment.