Skip to content

Commit

Permalink
excess ice restart from old files
Browse files Browse the repository at this point in the history
  • Loading branch information
mvdebolskiy committed Apr 29, 2022
1 parent f99106e commit 54ccc3b
Show file tree
Hide file tree
Showing 4 changed files with 123 additions and 22 deletions.
18 changes: 13 additions & 5 deletions src/biogeophys/WaterStateBulkType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -188,32 +188,40 @@ end subroutine InitBulkCold

!------------------------------------------------------------------------
subroutine RestartBulk(this, bounds, ncid, flag, &
watsat_col)
watsat_col, NLFilename, t_soisno_col)
!
! !DESCRIPTION:
! Read/Write module information to/from restart file.
!
! !USES:
use ncdio_pio , only : file_desc_t, ncd_double
use restUtilMod
use controlMod , only : use_excess_ice
!
! !ARGUMENTS:
class(waterstatebulk_type), intent(in) :: this
type(bounds_type), intent(in) :: bounds
type(file_desc_t), intent(inout) :: ncid ! netcdf id
character(len=*) , intent(in) :: flag ! 'read' or 'write'
real(r8) , intent(in) :: watsat_col (bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
character(len=*) , intent(in), optional :: NLFilename ! Namelist filename
real(r8) , intent(in), optional :: t_soisno_col(bounds%begc:bounds%endc,:) ! Soil column temperature (Kelvin)
!
! !LOCAL VARIABLES:
integer :: c,l,j
logical :: readvar
!------------------------------------------------------------------------

SHR_ASSERT_ALL_FL((ubound(watsat_col) == (/bounds%endc,nlevmaxurbgrnd/)) , sourcefile, __LINE__)

call this%restart (bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:))

if(use_excess_ice .and. flag == 'read') then
call this%restart (bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:), &
NLFilename= NLFilename, &
t_soisno_col=t_soisno_col(bounds%begc:bounds%endc,:))
else
call this%restart (bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:))
endif

call restartvar(ncid=ncid, flag=flag, &
varname=this%info%fname('INT_SNOW'), &
Expand Down
86 changes: 80 additions & 6 deletions src/biogeophys/WaterStateType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -353,7 +353,7 @@ subroutine InitCold(this, bounds, &
integer :: c,j,l,nlevs,g
integer :: nbedrock
real(r8) :: ratio
real(r8) :: exice_bulk_init(bounds%begc:bounds%endc)
real(r8) :: exice_bulk_init(bounds%begc:bounds%endc) ! initial excess ice amount (m/m)
!-----------------------------------------------------------------------

SHR_ASSERT_ALL_FL((ubound(h2osno_input_col) == (/bounds%endc/)) , sourcefile, __LINE__)
Expand Down Expand Up @@ -603,15 +603,16 @@ end subroutine InitCold

!------------------------------------------------------------------------
subroutine Restart(this, bounds, ncid, flag, &
watsat_col)
watsat_col, NLFilename, t_soisno_col)
!
! !DESCRIPTION:
! Read/Write module information to/from restart file.
!
! !USES:
use clm_varcon , only : denice, denh2o, pondmx, watmin
use landunit_varcon , only : istcrop, istdlak, istsoil
use clm_varcon , only : denice, denh2o, pondmx, watmin, tfrz
use landunit_varcon , only : istcrop, istdlak, istsoil, istice
use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall
use column_varcon , only : icol_road_imperv, icol_road_perv
use clm_time_manager , only : is_first_step, is_restart
use clm_varctl , only : bound_h2osoi
use ncdio_pio , only : file_desc_t, ncd_double
Expand All @@ -623,16 +624,28 @@ subroutine Restart(this, bounds, ncid, flag, &
type(file_desc_t), intent(inout) :: ncid ! netcdf id
character(len=*) , intent(in) :: flag ! 'read' or 'write'
real(r8) , intent(in) :: watsat_col (bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
character(len=*) , intent(in), optional :: NLFilename ! Namelist filename
real(r8) , intent(in), optional :: t_soisno_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd) ! Soil column temperature (Kelvin)
!
! !LOCAL VARIABLES:
integer :: p,c,l,j,nlevs
integer :: p,c,l,j,nlevs,g
integer :: nbedrock
logical :: readvar
real(r8) :: maxwatsat ! maximum porosity
real(r8) :: excess ! excess volumetric soil water
real(r8) :: totwat ! total soil water (mm)
character(len=256) :: l_NLFilename !local namelist filename var
real(r8) :: l_t_soisno_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd) !local soil column temperature (Kelvin)
real(r8) :: exice_bulk_init(bounds%begc:bounds%endc) !Initial excess ice ammout (m/m)
!------------------------------------------------------------------------

l_NLFilename = ''
l_t_soisno_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd) = tfrz + 1.0_r8
if (present(NLFilename)) then
l_NLFilename = NLFilename
l_t_soisno_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd) = t_soisno_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)
end if
SHR_ASSERT_ALL_FL((ubound(watsat_col) == (/bounds%endc,nlevmaxurbgrnd/)) , sourcefile, __LINE__)
SHR_ASSERT_ALL_FL((ubound(t_soisno_col) == (/bounds%endc,nlevmaxurbgrnd/)) , sourcefile, __LINE__)

call restartvar(ncid=ncid, flag=flag, &
varname=this%info%fname('H2OSFC'), &
Expand Down Expand Up @@ -730,7 +743,68 @@ subroutine Restart(this, bounds, ncid, flag, &
this%excess_ice_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%exice_melt_lev(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%exice_melt(bounds%begc:bounds%endc)=0.0_r8
else if (use_excess_ice .and. flag == 'read') then
! try to read the first variable
call restartvar(ncid=ncid, flag=flag, varname=this%info%fname('EXCESS_ICE'), xtype=ncd_double, &
dim1name='column', dim2name='levtot', switchdim=.true., &
long_name=this%info%lname('excess soil ice (vegetated landunits only)'), units='kg/m2', &
scale_by_thickness=.true., &
interpinic_flag='interp', readvar=readvar, data=this%excess_ice_col)
if (.not. readvar) then
! did not read the first variable
this%init_exice(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%excess_ice_col(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%exice_melt_lev(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%exice_melt(bounds%begc:bounds%endc)=0.0_r8

call this%exicestream%Init(bounds, l_NLFilename) ! get initial fraction of excess ice per column
call this%exicestream%CalcExcessIce(bounds, exice_bulk_init)
do c = bounds%begc,bounds%endc
g = col%gridcell(c)
l = col%landunit(c)
if (.not. lun%lakpoi(l)) then !not lake
if (lun%itype(l) == istsoil .or. lun%itype(l) == istcrop) then
if (use_bedrock) then
nbedrock = col%nbedrock(c)
else
nbedrock = nlevsoi
endif
do j = 2, nlevmaxurbgrnd ! ignore first layer
if (j<nbedrock .and. l_t_soisno_col(c,j) <= tfrz ) then
this%excess_ice_col(c,j) = col%dz(c,j)*denice*(exice_bulk_init(c))
this%init_exice(c,j) = col%dz(c,j)*denice*exice_bulk_init(c)
else
this%excess_ice_col(c,j) = 0.0_r8
this%init_exice(c,j) = 0.0_r8
endif
end do
endif
else ! just in case zeros for lakes and other columns
this%excess_ice_col(c,-nlevsno+1:nlevmaxurbgrnd) = 0.0_r8
this%init_exice(c,-nlevsno+1:nlevmaxurbgrnd) = 0.0_r8
endif
enddo
this%exice_melt_lev(bounds%begc:bounds%endc,-nlevsno+1:nlevmaxurbgrnd)=0.0_r8
this%exice_melt(bounds%begc:bounds%endc)=0.0_r8
else
! first var read succesfully, read the rest
call restartvar(ncid=ncid, flag=flag, varname=this%info%fname('EXICE_INIT'), xtype=ncd_double, &
dim1name='column', dim2name='levtot', switchdim=.true., &
long_name=this%info%lname('inital excess soil ice (vegetated landunits only)'), units='m/m', &
scale_by_thickness=.false., &
interpinic_flag='interp', readvar=readvar, data=this%init_exice)

call restartvar(ncid=ncid, flag=flag, varname=this%info%fname('EXCESS_MELT_LEV'), xtype=ncd_double, &
dim1name='column', dim2name='levtot', switchdim=.true., &
long_name=this%info%lname('melt from excess ice per layer (vegetated landunits only)'), units='kg/m2', &
scale_by_thickness=.true., &
interpinic_flag='interp', readvar=readvar, data=this%exice_melt_lev)

call restartvar(ncid=ncid, flag=flag, varname=this%info%fname('EXICE_MELT'), xtype=ncd_double, &
dim1name='column', &
long_name=this%info%lname('melt from excess ice (vegetated landunits only)'), units='m', &
interpinic_flag='interp', readvar=readvar, data=this%exice_melt)
endif
else
! have to at least define them
call restartvar(ncid=ncid, flag=flag, varname=this%info%fname('EXCESS_ICE'), xtype=ncd_double, &
Expand Down
26 changes: 20 additions & 6 deletions src/biogeophys/WaterType.F90
Original file line number Diff line number Diff line change
Expand Up @@ -727,18 +727,22 @@ end subroutine UpdateAccVars

!-----------------------------------------------------------------------
subroutine Restart(this, bounds, ncid, flag, writing_finidat_interp_dest_file, &
watsat_col)
watsat_col, NLFilename, t_soisno_col)
!
! !DESCRIPTION:
! Read/write information to/from restart file for all water variables
!
! USES:
use controlMod, only : use_excess_ice
! !ARGUMENTS:
class(water_type), intent(inout) :: this
type(bounds_type), intent(in) :: bounds
type(file_desc_t), intent(inout) :: ncid ! netcdf id
character(len=*) , intent(in) :: flag ! 'read', 'write' or 'define'
logical , intent(in) :: writing_finidat_interp_dest_file ! true if we are writing a finidat_interp_dest file (ignored for flag=='read')
real(r8) , intent(in) :: watsat_col (bounds%begc:, 1:) ! volumetric soil water at saturation (porosity)
character(len=*) , intent(in), optional :: NLFilename ! Namelist filename
real(r8) , intent(in), optional :: t_soisno_col(bounds%begc:bounds%endc,:) ! Soil column temperature (Kelvin)
!
! !LOCAL VARIABLES:
integer :: i
Expand All @@ -749,9 +753,14 @@ subroutine Restart(this, bounds, ncid, flag, writing_finidat_interp_dest_file, &
SHR_ASSERT_ALL_FL((ubound(watsat_col, 1) == bounds%endc), sourcefile, __LINE__)

call this%waterfluxbulk_inst%restartBulk (bounds, ncid, flag=flag)

call this%waterstatebulk_inst%restartBulk (bounds, ncid, flag=flag, &
if( use_excess_ice .and. flag == 'read') then
call this%waterstatebulk_inst%restartBulk (bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:), &
NLFilename=NLFilename, t_soisno_col=t_soisno_col(bounds%begc:bounds%endc,:))
else
call this%waterstatebulk_inst%restartBulk (bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:))
endif

call this%waterdiagnosticbulk_inst%restartBulk (bounds, ncid, flag=flag, &
writing_finidat_interp_dest_file=writing_finidat_interp_dest_file, &
Expand All @@ -760,9 +769,14 @@ subroutine Restart(this, bounds, ncid, flag, writing_finidat_interp_dest_file, &
do i = this%tracers_beg, this%tracers_end

call this%bulk_and_tracers(i)%waterflux_inst%Restart(bounds, ncid, flag=flag)

call this%bulk_and_tracers(i)%waterstate_inst%Restart(bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:))
if (use_excess_ice .and. flag == 'read') then
call this%bulk_and_tracers(i)%waterstate_inst%Restart(bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:), &
NLFilename=NLFilename, t_soisno_col=t_soisno_col(bounds%begc:bounds%endc,:))
else
call this%bulk_and_tracers(i)%waterstate_inst%Restart(bounds, ncid, flag=flag, &
watsat_col=watsat_col(bounds%begc:bounds%endc,:))
endif

call this%bulk_and_tracers(i)%waterdiagnostic_inst%Restart(bounds, ncid, flag=flag)

Expand Down
15 changes: 10 additions & 5 deletions src/main/clm_instMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -528,11 +528,16 @@ subroutine clm_instRest(bounds, ncid, flag, writing_finidat_interp_dest_file)
is_simple_buildtemp=IsSimpleBuildTemp(), is_prog_buildtemp=IsProgBuildTemp())

call soilstate_inst%restart (bounds, ncid, flag=flag)

call water_inst%restart(bounds, ncid, flag=flag, &
writing_finidat_interp_dest_file = writing_finidat_interp_dest_file, &
watsat_col = soilstate_inst%watsat_col(bounds%begc:bounds%endc,:))

if(use_excess_ice .and. flag == 'read') then
call water_inst%restart(bounds, ncid, flag=flag, &
writing_finidat_interp_dest_file = writing_finidat_interp_dest_file, &
watsat_col = soilstate_inst%watsat_col(bounds%begc:bounds%endc,:), &
NLFilename = NLFilename, t_soisno_col = temperature_inst%t_soisno_col(begc:endc, -nlevsno+1:))
else
call water_inst%restart(bounds, ncid, flag=flag, &
writing_finidat_interp_dest_file = writing_finidat_interp_dest_file, &
watsat_col = soilstate_inst%watsat_col(bounds%begc:bounds%endc,:))
endif
call irrigation_inst%restart (bounds, ncid, flag=flag)

call aerosol_inst%restart (bounds, ncid, flag=flag, &
Expand Down

0 comments on commit 54ccc3b

Please sign in to comment.