Skip to content

Commit

Permalink
Add baseline approach for lake dynamical land units
Browse files Browse the repository at this point in the history
  • Loading branch information
@Ivanderkelen
Ivanderkelen committed Sep 19, 2019
1 parent f9318cc commit ceaa825
Show file tree
Hide file tree
Showing 3 changed files with 239 additions and 77 deletions.
208 changes: 156 additions & 52 deletions src/biogeophys/TotalWaterAndHeatMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,12 +37,14 @@ module TotalWaterAndHeatMod
! routines parallel with the water routines.
public :: ComputeWaterMassNonLake ! Compute total water mass of non-lake columns
public :: ComputeWaterMassLake ! Compute total water mass of lake columns
public :: AccumulateLiqIceMassLake ! Accumulate lake water mass of lake columns, separated into liquid and ice
public :: ComputeLiqIceMassNonLake ! Compute total water mass of non-lake columns, separated into liquid and ice
public :: AccumulateSoilLiqIceMassNonLake ! Accumulate soil water mass of non-lake columns, separated into liquid and ice
public :: ComputeLiqIceMassLake ! Compute total water mass of lake columns, separated into liquid and ice
public :: ComputeHeatNonLake ! Compute heat content of non-lake columns
public :: AccumulateSoilHeatNonLake ! Accumulate soil heat content of non-lake columns
public :: ComputeHeatLake ! Compute heat content of lake columns
public :: AccumulateHeatLake ! Accumulate heat content of lake water of lake columns
public :: AdjustDeltaHeatForDeltaLiq ! Adjusts the change in gridcell heat content due to land cover change to account for the implicit heat flux associated with delta_liq
public :: LiquidWaterHeat ! Get the total heat content of some mass of liquid water at a given temperature

Expand Down Expand Up @@ -147,7 +149,7 @@ subroutine ComputeWaterMassLake(bounds, num_lakec, filter_lakec, &
integer , intent(in) :: num_lakec ! number of column lake points in column filter
integer , intent(in) :: filter_lakec(:) ! column filter for lake points
type(waterstate_type) , intent(in) :: waterstate_inst
type(lakestate_type) , intent(in) :: lakestate_inst
type(lakestate_type) , intent(in) :: lakestate_inst

! BUG(wjs, 2019-03-12, ESCOMP/ctsm#659) When we can accept answer changes to methane,
! remove this argument, always assuming it's true.
Expand Down Expand Up @@ -413,9 +415,7 @@ subroutine ComputeLiqIceMassLake(bounds, num_lakec, filter_lakec, &
!
! !LOCAL VARIABLES:
integer :: c, j, fc ! indices
real(r8) :: h2olak_liq(bounds%begc:bounds%endc,1:nlevlak) ! liquid water content of lake layer [kg/m²]
real(r8) :: h2olak_ice(bounds%begc:bounds%endc,1:nlevlak) ! ice water content of lake layer [kg/m²]


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

Expand All @@ -424,11 +424,9 @@ subroutine ComputeLiqIceMassLake(bounds, num_lakec, filter_lakec, &

associate( &
snl => col%snl , & ! Input: [integer (:) ] negative number of snow layers
dz_lake => col%dz_lake , & ! Input: [real(r8) (:,:) ] lake depth (m)
h2osno => waterstate_inst%h2osno_col , & ! Input: [real(r8) (:) ] snow water (mm H2O)
h2osoi_ice => waterstate_inst%h2osoi_ice_col , & ! Input: [real(r8) (:,:) ] ice lens (kg/m2)
h2osoi_liq => waterstate_inst%h2osoi_liq_col , & ! Input: [real(r8) (:,:) ] liquid water (kg/m2)
lake_icefrac => lakestate_inst%lake_icefrac_col, & ! Input: [real(r8) (:,:) ] lake ice fraction
dynbal_baseline_liq => waterstate_inst%dynbal_baseline_liq_col, & ! Input: [real(r8) (:) ] baseline liquid water content subtracted from each column's total liquid water calculation (mm H2O)
dynbal_baseline_ice => waterstate_inst%dynbal_baseline_ice_col & ! Input: [real(r8) (:) ] baseline ice content subtracted from each column's total ice calculation (mm H2O)
)
Expand All @@ -454,19 +452,10 @@ subroutine ComputeLiqIceMassLake(bounds, num_lakec, filter_lakec, &
ice_mass(c) = ice_mass(c) + h2osno(c)
end if
end do

! Lake water content
do j = 1, nlevlak
do fc = 1, num_lakec
c = filter_lakec(fc)
! calculate lake liq and ice content per lake layer first
h2olak_liq(c,j) = dz_lake(c,j) * denh2o * (1 - lake_icefrac(c,j))
h2olak_ice(c,j) = dz_lake(c,j) * denice * lake_icefrac(c,j)

liquid_mass(c) = liquid_mass(c) + h2olak_liq(c,j)
ice_mass(c) = ice_mass(c) + h2olak_ice(c,j)
end do
end do

! Lake water content
call AccumulateLiqIceMassLake(bounds, num_lakec, filter_lakec, &
lakestate_inst, liquid_mass, ice_mass)

! Soil water content of the soil under the lake
do j = 1, nlevgrnd
Expand All @@ -490,6 +479,60 @@ subroutine ComputeLiqIceMassLake(bounds, num_lakec, filter_lakec, &

end subroutine ComputeLiqIceMassLake

!-----------------------------------------------------------------------
subroutine AccumulateLiqIceMassLake(bounds, num_c, filter_c, &
lakestate_inst, liquid_mass, ice_mass)
!
! !DESCRIPTION:
! Accumulate lake water mass of lake columns, separated into liquid and ice.
!
! Adds to any existing values in liquid_mass and ice_mass.
!
! Note: Changes to this routine should generally be accompanied by similar changes to
! AccumulateSoilHeatLake
!
! !ARGUMENTS:
type(bounds_type) , intent(in) :: bounds
integer , intent(in) :: num_c ! number of column points in column filter (should not include lake points; can be a subset of the no-lake filter)
integer , intent(in) :: filter_c(:) ! column filter (should not include lake points; can be a subset of the no-lake filter)
type(lakestate_type) , intent(in) :: lakestate_inst
real(r8) , intent(inout) :: liquid_mass( bounds%begc: ) ! accumulated liquid mass (kg m-2)
real(r8) , intent(inout) :: ice_mass( bounds%begc: ) ! accumulated ice mass (kg m-2)
!
! !LOCAL VARIABLES:
integer :: c, j, fc ! indices
real(r8) :: h2olak_liq(bounds%begc:bounds%endc,1:nlevlak) ! liquid water content of lake layer [kg/m²]
real(r8) :: h2olak_ice(bounds%begc:bounds%endc,1:nlevlak) ! ice water content of lake layer [kg/m²]

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

SHR_ASSERT_ALL((ubound(liquid_mass) == [bounds%endc]), errMsg(sourcefile, __LINE__))
SHR_ASSERT_ALL((ubound(ice_mass) == [bounds%endc]), errMsg(sourcefile, __LINE__))

associate( &
lake_icefrac => lakestate_inst%lake_icefrac_col, & ! Input: [real(r8) (:,:) ] lake ice fraction
dz_lake => col%dz_lake & ! Input: [real(r8) (:,:) ] lake depth (m)
)

! Lake water content
do j = 1, nlevlak
do fc = 1, num_c
c = filter_c(fc)
! calculate lake liq and ice content per lake layer first
h2olak_liq(c,j) = dz_lake(c,j) * denh2o * (1 - lake_icefrac(c,j))
h2olak_ice(c,j) = dz_lake(c,j) * denice * lake_icefrac(c,j)

liquid_mass(c) = liquid_mass(c) + h2olak_liq(c,j)
ice_mass(c) = ice_mass(c) + h2olak_ice(c,j)
end do
end do

end associate

end subroutine AccumulateLiqIceMassLake


!-----------------------------------------------------------------------
subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, &
urbanparams_inst, soilstate_inst, &
Expand Down Expand Up @@ -869,9 +912,7 @@ subroutine ComputeHeatLake(bounds, num_lakec, filter_lakec, &
real(r8) :: heat_dry_mass(bounds%begc:bounds%endc) ! sum of heat content: dry mass [J/m^2]
real(r8) :: heat_ice(bounds%begc:bounds%endc) ! sum of heat content: ice [J/m^2]
real(r8) :: latent_heat_liquid(bounds%begc:bounds%endc) ! sum of latent heat content of liquid water [J/m^2]
real(r8) :: h2olak_liq(bounds%begc:bounds%endc,1:nlevlak) ! liquid water content of lake layer [kg/m²]
real(r8) :: h2olak_ice(bounds%begc:bounds%endc,1:nlevlak) ! ice water content of lake layer [kg/m²]


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

Expand All @@ -882,16 +923,13 @@ subroutine ComputeHeatLake(bounds, num_lakec, filter_lakec, &
associate( &
snl => col%snl, & ! number of snow layers
dz => col%dz, & ! layer depth (m)
dz_lake => col%dz_lake, & ! lake layer depth (m)
watsat => soilstate_inst%watsat_col, & ! volumetric soil water at saturation (porosity)
csol => soilstate_inst%csol_col, & ! heat capacity, soil solids (J/m**3/Kelvin)
t_lake => temperature_inst%t_lake_col, & ! lake temperature (K)
t_soisno => temperature_inst%t_soisno_col, & ! soil temperature (Kelvin)
dynbal_baseline_heat => temperature_inst%dynbal_baseline_heat_col, & ! Input: [real(r8) (:) ] baseline heat content subtracted from each column's total heat calculation (J/m2)
h2osoi_liq => waterstate_inst%h2osoi_liq_col, & ! liquid water (kg/m2)
h2osoi_ice => waterstate_inst%h2osoi_ice_col, & ! frozen water (kg/m2)
h2osno => waterstate_inst%h2osno_col, & ! snow water (mm H2O)
lake_icefrac => lakestate_inst%lake_icefrac_col & ! Input: [real(r8) (:,:) ] lake ice fraction
h2osno => waterstate_inst%h2osno_col & ! snow water (mm H2O)
)

do fc = 1, num_lakec
Expand Down Expand Up @@ -927,7 +965,7 @@ subroutine ComputeHeatLake(bounds, num_lakec, filter_lakec, &
end if
end do

! Soil water content of the soil under the lake
! Soil heat content of the soil under the lake
do j = 1,nlevgrnd
do fc = 1, num_lakec
c = filter_lakec(fc)
Expand All @@ -944,35 +982,15 @@ subroutine ComputeHeatLake(bounds, num_lakec, filter_lakec, &
TempToHeat(temp = t_soisno(c,j), cv = (h2osoi_ice(c,j)*cpice))
end do
end do

! TODO(wjs, 2017-03-11) Include heat content of water in lakes, once we include
! lake water as an explicit water state (https://github.com/ESCOMP/ctsm/issues/200)

! calculate heat content of lake itself
do j = 1, nlevlak

do fc = 1, num_lakec
c = filter_lakec(fc)
! liquid heat
h2olak_liq(c,j) = dz_lake(c,j) * denh2o * (1 - lake_icefrac(c,j))
call AccumulateLiquidWaterHeat( &
temp = t_lake(c,j), &
h2o = h2olak_liq(c,j), &
cv_liquid = cv_liquid(c), &
heat_liquid = heat_liquid(c), &
latent_heat_liquid = latent_heat_liquid(c))
! ice heat
h2olak_ice(c,j) = dz_lake(c,j) * denice * lake_icefrac(c,j)
heat_ice(c) = heat_ice(c) + &
TempToHeat(temp=t_lake(c,j), cv = (h2olak_ice(c,j) * cpice))
end do
end do

! write(iulog,*) 'lake heat (J/m^2)', heat_lake(c)+latent_heat_liquid(c)

do fc = 1, num_lakec
c = filter_lakec(fc)
heat(c) = heat_dry_mass(c) + heat_ice(c) + heat_liquid(c) + latent_heat_liquid(c)
end do

! Lake water heat content
call AccumulateHeatLake(bounds, num_lakec, filter_lakec, temperature_inst, lakestate_inst, &
heat, heat_liquid, cv_liquid)

! Subtract baselines set in initialization
!
Expand All @@ -996,6 +1014,92 @@ subroutine ComputeHeatLake(bounds, num_lakec, filter_lakec, &

end subroutine ComputeHeatLake

!-----------------------------------------------------------------------
subroutine AccumulateHeatLake(bounds, num_c, filter_c, &
temperature_inst, lakestate_inst, &
heat, heat_liquid, cv_liquid)
!
! !DESCRIPTION:
! Accumulate heat of lake water in lake columns.
!
! Adds to any existing values in heat, heat_liquid and cv_liquid.
!
! Note: Changes to this routine should generally be accompanied by similar changes to
! AccumulateLiqIceMassLake
!
! !ARGUMENTS:
type(bounds_type) , intent(in) :: bounds
integer , intent(in) :: num_c ! number of column points in column filter (should not include lake points; can be a subset of the no-lake filter)
integer , intent(in) :: filter_c(:) ! column filter (should not include lake points; can be a subset of the no-lake filter)
type(temperature_type) , intent(in) :: temperature_inst
type(lakestate_type) , intent(in) :: lakestate_inst
real(r8) , intent(inout) :: heat( bounds%begc: ) ! accumulated heat content [J/m^2]
real(r8) , intent(inout) :: heat_liquid( bounds%begc: ) ! accumulated heat content: liquid water, excluding latent heat [J/m^2]
real(r8) , intent(inout) :: cv_liquid( bounds%begc: ) ! accumulated liquid water heat capacity [J/(m^2 K)]
!
! !LOCAL VARIABLES:
integer :: fc
integer :: l, c, j
real(r8) :: h2olak_liq(bounds%begc:bounds%endc,1:nlevlak) ! liquid water content of lake layer [kg/m²]
real(r8) :: h2olak_ice(bounds%begc:bounds%endc,1:nlevlak) ! ice water content of lake layer [kg/m²]
real(r8) :: lake_heat_liquid(bounds%begc:bounds%endc) ! sum of heat content: liquid water in lake, excluding latent heat [J/m^2]
real(r8) :: lake_heat_ice(bounds%begc:bounds%endc) ! sum of heat content: ice in lake [J/m^2]
real(r8) :: lake_latent_heat_liquid(bounds%begc:bounds%endc) ! sum of heat content: latent heat of liquid water in lake [J/m^2]


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

SHR_ASSERT_ALL((ubound(heat) == [bounds%endc]), errMsg(sourcefile, __LINE__))
SHR_ASSERT_ALL((ubound(heat_liquid) == [bounds%endc]), errMsg(sourcefile, __LINE__))
SHR_ASSERT_ALL((ubound(cv_liquid) == [bounds%endc]), errMsg(sourcefile, __LINE__))

associate( &
dz_lake => col%dz_lake, & ! lake layer depth (m)
t_lake => temperature_inst%t_lake_col, & ! lake temperature (K)
lake_icefrac => lakestate_inst%lake_icefrac_col & ! Input: [real(r8) (:,:) ] lake ice fraction
)

do fc = 1, num_c
c = filter_c(fc)

lake_heat_liquid(c) = 0._r8
lake_heat_ice(c) = 0._r8
lake_latent_heat_liquid(c) = 0._r8
end do


! calculate heat content of lake itself
do j = 1, nlevlak
do fc = 1, num_c
c = filter_c(fc)
! liquid heat
h2olak_liq(c,j) = dz_lake(c,j) * denh2o * (1 - lake_icefrac(c,j))
call AccumulateLiquidWaterHeat( &
temp = t_lake(c,j), &
h2o = h2olak_liq(c,j), &
cv_liquid = cv_liquid(c), &
heat_liquid = lake_heat_liquid(c), &
latent_heat_liquid = lake_latent_heat_liquid(c))
! ice heat
h2olak_ice(c,j) = dz_lake(c,j) * denice * lake_icefrac(c,j)
lake_heat_ice(c) = lake_heat_ice(c) + &
TempToHeat(temp=t_lake(c,j), cv = (h2olak_ice(c,j) * cpice))
end do
end do

! add ice heat and liquid heat together (look at this part)
do fc = 1, num_c
c = filter_c(fc)
heat_liquid(c) = heat_liquid(c) + lake_heat_liquid(c)
heat(c) = heat(c) + lake_heat_ice(c) + &
lake_heat_liquid(c) + lake_latent_heat_liquid(c)
end do

end associate

end subroutine AccumulateHeatLake

!-----------------------------------------------------------------------
subroutine AdjustDeltaHeatForDeltaLiq(bounds, delta_liq, &
liquid_water_temp1, liquid_water_temp2, &
Expand Down
Loading

0 comments on commit ceaa825

Please sign in to comment.