Skip to content

Commit

Permalink
Merge pull request ESCOMP#17 from NCAR/vic_no_infl_excess_runoff
Browse files Browse the repository at this point in the history
(This has been copied to pull request ESCOMP#190 from
ESCOMP/vic_no_inf_excess_runoff)

Remove infiltration excess runoff for VIC

Martyn Clark reviewed the VIC implementation, and felt that the current
implementation of infiltration excess runoff is inconsistent with the
standard VIC implementation. It appears that what was being called VIC's
infiltration excess runoff was actually just an attempt to give a better
numerical approximation to the solution for saturated surface excess
runoff. So deleting this leaves only a first-order approximation to
VIC's saturated surface excess runoff.

Eventually we may want to put in place a more accurate solution for
VIC's saturated surface excess runoff. But Martyn's feeling is that this
can come in with other changes we want to make regarding numerical
solutions in CTSM.
  • Loading branch information
billsacks committed Dec 28, 2017
2 parents dd8f549 + 5f9b7a2 commit fc06b68
Show file tree
Hide file tree
Showing 2 changed files with 23 additions and 78 deletions.
97 changes: 19 additions & 78 deletions src/biogeophys/InfiltrationExcessRunoffMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -53,13 +53,22 @@ module InfiltrationExcessRunoffMod
procedure, private :: InitCold

procedure, private, nopass :: ComputeQinmaxHksat
procedure, private, nopass :: ComputeQinmaxVic
end type infiltration_excess_runoff_type

! !PRIVATE DATA MEMBERS:

! For methods that don't generate any infiltration excess runoff, we get this end result
! by specifying a huge qinmax value - i.e., an effectively infinite max infiltration
! rate.
!
! 1e200 mm H2O /s seems large enough to be effectively infinite, while not being so
! large as to cause floating point overflows elsewhere.
real(r8), parameter :: qinmax_unlimited = 1.e200_r8 ! mm H2O /s

! The 'none' option avoids generating any infiltration excess runoff by setting qinmax
! to a huge value
integer, parameter :: QINMAX_METHOD_NONE = 0
integer, parameter :: QINMAX_METHOD_HKSAT = 1
integer, parameter :: QINMAX_METHOD_VIC = 2

character(len=*), parameter, private :: sourcefile = &
__FILE__
Expand Down Expand Up @@ -155,7 +164,7 @@ subroutine InitCold(this, bounds)

! TODO(wjs, 2017-08-14) We'll read qinmax_method from namelist.
if (use_vichydro) then
this%qinmax_method = QINMAX_METHOD_VIC
this%qinmax_method = QINMAX_METHOD_NONE
else
this%qinmax_method = QINMAX_METHOD_HKSAT
end if
Expand Down Expand Up @@ -203,16 +212,17 @@ subroutine InfiltrationExcessRunoff(this, bounds, num_hydrologyc, filter_hydrolo
)

select case (this%qinmax_method)
case (QINMAX_METHOD_NONE)
! NOTE(wjs, 2017-09-01) I'm putting this here for clarity, though it could be
! moved to initialization for efficiency
do fc = 1, num_hydrologyc
c = filter_hydrologyc(fc)
qinmax_on_unsaturated_area(c) = qinmax_unlimited
end do
case (QINMAX_METHOD_HKSAT)
call this%ComputeQinmaxHksat(bounds, num_hydrologyc, filter_hydrologyc, &
soilhydrology_inst, soilstate_inst, &
qinmax_on_unsaturated_area = qinmax_on_unsaturated_area(bounds%begc:bounds%endc))
case (QINMAX_METHOD_VIC)
call this%ComputeQinmaxVic(bounds, num_hydrologyc, filter_hydrologyc, &
soilhydrology_inst, &
fsat = fsat(bounds%begc:bounds%endc), &
qflx_in_soil = qflx_in_soil(bounds%begc:bounds%endc), &
qinmax_on_unsaturated_area = qinmax_on_unsaturated_area(bounds%begc:bounds%endc))
case default
write(iulog,*) subname//' ERROR: Unrecognized qinmax_method: ', this%qinmax_method
call endrun(subname//' ERROR: Unrecognized qinmax_method')
Expand Down Expand Up @@ -270,73 +280,4 @@ subroutine ComputeQinmaxHksat(bounds, num_hydrologyc, filter_hydrologyc, &

end subroutine ComputeQinmaxHksat

!-----------------------------------------------------------------------
subroutine ComputeQinmaxVic(bounds, num_hydrologyc, filter_hydrologyc, &
soilhydrology_inst, &
fsat, qflx_in_soil, qinmax_on_unsaturated_area)
!
! !DESCRIPTION:
! Compute qinmax using the VIC parameterization
!
! Citation: Wood et al. 1992, "A land-surface hydrology parameterization with subgrid
! variability for general circulation models", JGR 97(D3), 2717-2728.
!
! !ARGUMENTS:
type(bounds_type), intent(in) :: bounds
integer, intent(in) :: num_hydrologyc ! number of column soil points in column filter
integer, intent(in) :: filter_hydrologyc(:) ! column filter for soil points
type(soilhydrology_type) , intent(in) :: soilhydrology_inst
real(r8) , intent(in) :: fsat( bounds%begc: ) ! fractional area with water table at surface
real(r8) , intent(in) :: qflx_in_soil( bounds%begc: ) ! surface input to soil (mm/s)
real(r8) , intent(inout) :: qinmax_on_unsaturated_area( bounds%begc: ) ! maximum infiltration rate on the unsaturated fraction of the column (mm H2O /s)
!
! !LOCAL VARIABLES:
integer :: fc, c
real(r8) :: dtime ! land model time step (sec)
real(r8) :: top_icefrac ! ice fraction in top VIC layers
real(r8) :: max_infil ! max infiltration capacity using the VIC parameterization (mm)
real(r8) :: i_0 ! average soil moisture in top VIC layers (mm)
real(r8) :: rsurf_vic ! surface runoff based on the VIC parameterization
real(r8) :: basis ! variable soil moisture holding capacity in top VIC layers for runoff calculation

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

SHR_ASSERT_ALL((ubound(fsat) == (/bounds%endc/)), errMsg(sourcefile, __LINE__))
SHR_ASSERT_ALL((ubound(qflx_in_soil) == (/bounds%endc/)), errMsg(sourcefile, __LINE__))
SHR_ASSERT_ALL((ubound(qinmax_on_unsaturated_area) == (/bounds%endc/)), errMsg(sourcefile, __LINE__))

associate( &
top_max_moist => soilhydrology_inst%top_max_moist_col, & ! Input: [real(r8) (:) ] maximum soil moisture in top VIC layers
top_moist => soilhydrology_inst%top_moist_col , & ! Input: [real(r8) (:) ] soil moisture in top VIC layers
top_ice => soilhydrology_inst%top_ice_col , & ! Input: [real(r8) (:) ] ice len in top VIC layers
b_infil => soilhydrology_inst%b_infil_col & ! Input: [real(r8) (:) ] VIC b infiltration parameter
)

dtime = get_step_size()

do fc = 1, num_hydrologyc
c = filter_hydrologyc(fc)
top_icefrac = min(1._r8,top_ice(c)/top_max_moist(c))
max_infil = (1._r8+b_infil(c)) * top_max_moist(c)
i_0 = max_infil * (1._r8 - (1._r8 - fsat(c))**(1._r8/b_infil(c)))
if(qflx_in_soil(c) <= 0._r8) then
rsurf_vic = 0._r8
else if(max_infil <= 0._r8) then
rsurf_vic = qflx_in_soil(c)
else if((i_0 + qflx_in_soil(c)*dtime) > max_infil) then !(Eq.(3a) Wood et al. 1992)
rsurf_vic = (qflx_in_soil(c)*dtime - top_max_moist(c) + top_moist(c))/dtime
else !(Eq.(3b) Wood et al. 1992)
basis = 1._r8 - (i_0 + qflx_in_soil(c)*dtime)/max_infil
rsurf_vic = (qflx_in_soil(c)*dtime - top_max_moist(c) + top_moist(c) &
+ top_max_moist(c) * basis**(1._r8 + b_infil(c)))/dtime
end if
rsurf_vic = min(qflx_in_soil(c), rsurf_vic)
qinmax_on_unsaturated_area(c) = 10._r8**(-e_ice*top_icefrac)*(qflx_in_soil(c) - rsurf_vic)
end do

end associate

end subroutine ComputeQinmaxVic

end module InfiltrationExcessRunoffMod
4 changes: 4 additions & 0 deletions src/biogeophys/SaturatedExcessRunoffMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -334,6 +334,10 @@ subroutine ComputeFsatVic(bounds, num_hydrologyc, filter_hydrologyc, &
! Citation: Wood et al. 1992, "A land-surface hydrology parameterization with subgrid
! variability for general circulation models", JGR 97(D3), 2717-2728.
!
! This implementation gives a first-order approximation to saturated excess runoff.
! For now we're not including the more exact analytical solution, or even a better
! numerical approximation.
!
! !ARGUMENTS:
type(bounds_type), intent(in) :: bounds
integer, intent(in) :: num_hydrologyc ! number of column soil points in column filter
Expand Down

0 comments on commit fc06b68

Please sign in to comment.