Skip to content

Commit

Permalink
Second step in allowing for nlevgrnd to be less than nlevurb.
Browse files Browse the repository at this point in the history 
Code passes as bfb with ctsm1.0.dev035 for this test:
ERP_Ld3.f09_g17.I1850Clm50BgcCropCru.cheyenne_intel.clm-ciso
  • Loading branch information
@olyson
olyson committed Sep 18, 2019
1 parent 5ef0675 commit 998473e
Show file tree
Hide file tree
Showing 8 changed files with 225 additions and 91 deletions.
41 changes: 37 additions & 4 deletions src/biogeophys/CanopyTemperatureMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,9 @@ module CanopyTemperatureMod
!------------------------------------------------------------------------------
subroutine CanopyTemperature(bounds, &
num_nolakec, filter_nolakec, num_nolakep, filter_nolakep, &
!KO
num_urbanc, filter_urbanc, &
!KO
clm_fates, &
atm2lnd_inst, canopystate_inst, soilstate_inst, frictionvel_inst, &
waterstatebulk_inst, waterdiagnosticbulk_inst, waterfluxbulk_inst, &
Expand Down Expand Up @@ -90,6 +93,10 @@ subroutine CanopyTemperature(bounds, &
integer , intent(in) :: filter_nolakec(:) ! column filter for non-lake points
integer , intent(in) :: num_nolakep ! number of column non-lake points in patch filter
integer , intent(in) :: filter_nolakep(:) ! patch filter for non-lake points
!KO
integer , intent(in) :: num_urbanc ! number of urban columns in clump
integer , intent(in) :: filter_urbanc(:) ! urban column filter
!KO
type(hlm_fates_interface_type), intent(inout) :: clm_fates
type(atm2lnd_type) , intent(in) :: atm2lnd_inst
type(canopystate_type) , intent(inout) :: canopystate_inst
Expand All @@ -107,6 +114,9 @@ subroutine CanopyTemperature(bounds, &
integer :: j ! soil/snow level index
integer :: fp ! lake filter patch index
integer :: fc ! lake filter column index
!KO
integer :: nlev ! greater of nlevgrnd and nlevurb
!KO
real(r8) :: qred ! soil surface relative humidity
real(r8) :: avmuir ! ir inverse optical depth per unit leaf area
real(r8) :: eg ! water vapor pressure at temperature T [pa]
Expand Down Expand Up @@ -218,17 +228,40 @@ subroutine CanopyTemperature(bounds, &
do j = -nlevsno+1, nlevgrnd
do fc = 1,num_nolakec
c = filter_nolakec(fc)
if ((col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall &
.or. col%itype(c) == icol_roof) .and. j > nlevurb) then
tssbef(c,j) = spval
else
!KO if ((col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall &
!KO .or. col%itype(c) == icol_roof) .and. j > nlevurb) then
!KO tssbef(c,j) = spval
!KO else
!KO
if (col%itype(c) /= icol_sunwall .and. col%itype(c) /= icol_shadewall &
.and. col%itype(c) /= icol_roof) then
!KO
tssbef(c,j) = t_soisno(c,j)
end if
! record t_h2osfc prior to updating
t_h2osfc_bef(c) = t_h2osfc(c)
end do
end do

!KO
nlev = max0(nlevgrnd,nlevurb)
do j = -nlevsno+1, nlev
do fc = 1,num_urbanc
c = filter_urbanc(fc)
if (col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall &
.or. col%itype(c) == icol_roof) then
if (j > nlevurb) then
tssbef(c,j) = spval
else
tssbef(c,j) = t_soisno(c,j)
end if
end if
! record t_h2osfc prior to updating
t_h2osfc_bef(c) = t_h2osfc(c)
end do
end do
!KO

! calculate moisture stress/resistance for soil evaporation
call calc_soilevap_resis(bounds, num_nolakec, filter_nolakec, &
soilstate_inst, waterstatebulk_inst, waterdiagnosticbulk_inst, temperature_inst)
Expand Down
22 changes: 19 additions & 3 deletions src/biogeophys/HydrologyDrainageMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -151,14 +151,30 @@ subroutine HydrologyDrainage(bounds, &
do j = 1, nlevgrnd
do fc = 1, num_nolakec
c = filter_nolakec(fc)
if ((ctype(c) == icol_sunwall .or. ctype(c) == icol_shadewall &
.or. ctype(c) == icol_roof) .and. j > nlevurb) then
else
!KO if ((ctype(c) == icol_sunwall .or. ctype(c) == icol_shadewall &
!KO .or. ctype(c) == icol_roof) .and. j > nlevurb) then
!KO
if (ctype(c) /= icol_sunwall .and. ctype(c) /= icol_shadewall &
.and. ctype(c) /= icol_roof) then
!KO
!KO else
h2osoi_vol(c,j) = h2osoi_liq(c,j)/(dz(c,j)*denh2o) + h2osoi_ice(c,j)/(dz(c,j)*denice)
end if
end do
end do

!KO
do j = 1, nlevurb
do fc = 1, num_urbanc
c = filter_urbanc(fc)
if (col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall &
.or. col%itype(c) == icol_roof) then
h2osoi_vol(c,j) = h2osoi_liq(c,j)/(dz(c,j)*denh2o) + h2osoi_ice(c,j)/(dz(c,j)*denice)
end if
end do
end do
!KO

call ComputeWaterMassNonLake(bounds, num_nolakec, filter_nolakec, &
waterstatebulk_inst, waterdiagnosticbulk_inst, &
subtract_dynbal_baselines = .false., &
Expand Down
22 changes: 19 additions & 3 deletions src/biogeophys/HydrologyNoDrainageMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -514,13 +514,29 @@ subroutine HydrologyNoDrainage(bounds, &
do j = 1, nlevgrnd
do fc = 1, num_nolakec
c = filter_nolakec(fc)
if ((ctype(c) == icol_sunwall .or. ctype(c) == icol_shadewall &
.or. ctype(c) == icol_roof) .and. j > nlevurb) then
else
!KO if ((ctype(c) == icol_sunwall .or. ctype(c) == icol_shadewall &
!KO .or. ctype(c) == icol_roof) .and. j > nlevurb) then
!KO
if (ctype(c) /= icol_sunwall .and. ctype(c) /= icol_shadewall &
.and. ctype(c) /= icol_roof) then
!KO
!KO else
h2osoi_vol(c,j) = h2osoi_liq(c,j)/(dz(c,j)*denh2o) + h2osoi_ice(c,j)/(dz(c,j)*denice)
end if
end do
end do

!KO
do j = 1, nlevurb
do fc = 1, num_urbanc
c = filter_urbanc(fc)
if (col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall &
.or. col%itype(c) == icol_roof) then
h2osoi_vol(c,j) = h2osoi_liq(c,j)/(dz(c,j)*denh2o) + h2osoi_ice(c,j)/(dz(c,j)*denice)
end if
end do
end do
!KO

! if (use_cn) then
! Update soilpsi.
Expand Down
78 changes: 47 additions & 31 deletions src/biogeophys/SoilWaterMovementMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2188,43 +2188,59 @@ subroutine TridiagonalCol (ci, lbj, ubj, jtop, a, b, c, r, u)
bet = b(jtop)

do j = lbj, ubj
if ((col%itype(ci) == icol_sunwall .or. col%itype(ci) == icol_shadewall &
.or. col%itype(ci) == icol_roof) .and. j <= nlevurb) then
if (j >= jtop) then
if (j == jtop) then
u(j) = r(j) / bet
else
gam(j) = c(j-1) / bet
bet = b(j) - a(j) * gam(j)
u(j) = (r(j) - a(j)*u(j-1)) / bet
end if
end if
else if (col%itype(ci) /= icol_sunwall .and. col%itype(ci) /= icol_shadewall &
.and. col%itype(ci) /= icol_roof) then
if (j >= jtop) then
if (j == jtop) then
u(j) = r(j) / bet
else
gam(j) = c(j-1) / bet
bet = b(j) - a(j) * gam(j)
u(j) = (r(j) - a(j)*u(j-1)) / bet
end if
!KO if ((col%itype(ci) == icol_sunwall .or. col%itype(ci) == icol_shadewall &
!KO .or. col%itype(ci) == icol_roof) .and. j <= nlevurb) then
!KO if (j >= jtop) then
!KO if (j == jtop) then
!KO u(j) = r(j) / bet
!KO else
!KO gam(j) = c(j-1) / bet
!KO bet = b(j) - a(j) * gam(j)
!KO u(j) = (r(j) - a(j)*u(j-1)) / bet
!KO end if
!KO end if
!KO else if (col%itype(ci) /= icol_sunwall .and. col%itype(ci) /= icol_shadewall &
!KO .and. col%itype(ci) /= icol_roof) then
!KO if (j >= jtop) then
!KO if (j == jtop) then
!KO u(j) = r(j) / bet
!KO else
!KO gam(j) = c(j-1) / bet
!KO bet = b(j) - a(j) * gam(j)
!KO u(j) = (r(j) - a(j)*u(j-1)) / bet
!KO end if
!KO end if
!KO end if
!KO
if (j >= jtop) then
if (j == jtop) then
u(j) = r(j) / bet
else
gam(j) = c(j-1) / bet
bet = b(j) - a(j) * gam(j)
u(j) = (r(j) - a(j)*u(j-1)) / bet
end if
end if
!KO
end do

do j = ubj-1,lbj,-1
if ((col%itype(ci) == icol_sunwall .or. col%itype(ci) == icol_shadewall &
.or. col%itype(ci) == icol_roof) .and. j <= nlevurb-1) then
if (j >= jtop) then
u(j) = u(j) - gam(j+1) * u(j+1)
end if
else if (col%itype(ci) /= icol_sunwall .and. col%itype(ci) /= icol_shadewall &
.and. col%itype(ci) /= icol_roof) then
if (j >= jtop) then
u(j) = u(j) - gam(j+1) * u(j+1)
end if
!KO if ((col%itype(ci) == icol_sunwall .or. col%itype(ci) == icol_shadewall &
!KO .or. col%itype(ci) == icol_roof) .and. j <= nlevurb-1) then
!KO if (j >= jtop) then
!KO u(j) = u(j) - gam(j+1) * u(j+1)
!KO end if
!KO else if (col%itype(ci) /= icol_sunwall .and. col%itype(ci) /= icol_shadewall &
!KO .and. col%itype(ci) /= icol_roof) then
!KO if (j >= jtop) then
!KO u(j) = u(j) - gam(j+1) * u(j+1)
!KO end if
!KO end if
!KO
if (j >= jtop) then
u(j) = u(j) - gam(j+1) * u(j+1)
end if
!KO
end do

end subroutine TridiagonalCol
Expand Down
51 changes: 39 additions & 12 deletions src/biogeophys/TotalWaterAndHeatMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -348,6 +348,9 @@ subroutine AccumulateSoilLiqIceMassNonLake(bounds, num_c, filter_c, &
! !LOCAL VARIABLES:
integer :: c, j, fc ! indices
logical :: has_h2o ! whether this point potentially has water to add
!KO
integer :: nlev ! greater of nlevgrnd and nlevurb
!KO

character(len=*), parameter :: subname = 'AccumulateSoilLiqIceMassNonLake'
!-----------------------------------------------------------------------
Expand All @@ -360,7 +363,13 @@ subroutine AccumulateSoilLiqIceMassNonLake(bounds, num_c, filter_c, &
h2osoi_liq => waterstate_inst%h2osoi_liq_col & ! Input: [real(r8) (:,:) ] liquid water (kg/m2)
)

do j = 1, nlevgrnd
!KO
nlev = max0(nlevgrnd,nlevurb)
!KO
!KO do j = 1, nlevgrnd
!KO
do j = 1, nlev
!KO
do fc = 1, num_c
c = filter_c(fc)
if (col%itype(c) == icol_sunwall .or. col%itype(c) == icol_shadewall) then
Expand All @@ -372,7 +381,12 @@ subroutine AccumulateSoilLiqIceMassNonLake(bounds, num_c, filter_c, &
has_h2o = .false.
end if
else
has_h2o = .true.
!KO has_h2o = .true.
!KO
if (j <= nlevgrnd) then
has_h2o = .true.
end if
!KO
end if

if (has_h2o) then
Expand Down Expand Up @@ -719,6 +733,9 @@ subroutine AccumulateSoilHeatNonLake(bounds, num_c, filter_c, &
! !LOCAL VARIABLES:
integer :: fc
integer :: l, c, j
!KO
integer :: nlev ! greater of nlevgrnd and nlevurb
!KO
logical :: has_h2o ! whether this point potentially has water to add

real(r8) :: soil_heat_liquid(bounds%begc:bounds%endc) ! sum of heat content: liquid water in soil, excluding latent heat [J/m^2]
Expand Down Expand Up @@ -755,7 +772,11 @@ subroutine AccumulateSoilHeatNonLake(bounds, num_c, filter_c, &
soil_latent_heat_liquid(c) = 0._r8
end do

do j = 1, nlevgrnd
!KO
nlev = max0(nlevgrnd,nlevurb)
do j = 1, nlev
!KO
!KO do j = 1, nlevgrnd
do fc = 1, num_c
c = filter_c(fc)
l = col%landunit(c)
Expand All @@ -777,17 +798,23 @@ subroutine AccumulateSoilHeatNonLake(bounds, num_c, filter_c, &
end if

else
has_h2o = .true.
!KO
if (j <= nlevgrnd) then
!KO
has_h2o = .true.

if (col%itype(c) == icol_road_imperv .and. j <= nlev_improad(l)) then
soil_heat_dry_mass(c) = soil_heat_dry_mass(c) + &
TempToHeat(temp = t_soisno(c,j), cv = (cv_improad(l,j) * dz(c,j)))
else if (lun%itype(l) /= istwet .and. lun%itype(l) /= istice_mec) then
! Note that this also includes impervious roads below nlev_improad (where
! we have soil)
soil_heat_dry_mass(c) = soil_heat_dry_mass(c) + &
TempToHeat(temp = t_soisno(c,j), cv = (csol(c,j)*(1-watsat(c,j))*dz(c,j)))
if (col%itype(c) == icol_road_imperv .and. j <= nlev_improad(l)) then
soil_heat_dry_mass(c) = soil_heat_dry_mass(c) + &
TempToHeat(temp = t_soisno(c,j), cv = (cv_improad(l,j) * dz(c,j)))
else if (lun%itype(l) /= istwet .and. lun%itype(l) /= istice_mec) then
! Note that this also includes impervious roads below nlev_improad (where
! we have soil)
soil_heat_dry_mass(c) = soil_heat_dry_mass(c) + &
TempToHeat(temp = t_soisno(c,j), cv = (csol(c,j)*(1-watsat(c,j))*dz(c,j)))
end if
!KO
end if
!KO
end if

if (has_h2o) then
Expand Down
Loading

0 comments on commit 998473e

Please sign in to comment.