Skip to content

Commit

Permalink
Added routine to calculate surface inundated fraction in polygonal tu…
Browse files Browse the repository at this point in the history
…ndra, but still needs work defining values for microtopographic parameters before it can be tested
  • Loading branch information
chuckaustin committed May 6, 2024
1 parent 735b5b2 commit fc8f584
Showing 1 changed file with 49 additions and 17 deletions.
66 changes: 49 additions & 17 deletions components/elm/src/biogeophys/CanopyHydrologyMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -783,7 +783,7 @@ subroutine FracH2OSfc(bounds, num_h2osfc, filter_h2osfc, &
! !USES:
use shr_const_mod , only : shr_const_pi
use shr_spfn_mod , only : erf => shr_spfn_erf
use landunit_varcon , only : istsoil, istcrop
use landunit_varcon , only : istsoil, istcrop, ispolygon
!
! !ARGUMENTS:
type(bounds_type) , intent(in) :: bounds
Expand All @@ -794,24 +794,28 @@ subroutine FracH2OSfc(bounds, num_h2osfc, filter_h2osfc, &
integer , intent(in), optional :: no_update ! flag to make calculation w/o updating variables
!
! !LOCAL VARIABLES:
integer :: c,f,l ! indices
integer :: c,f,l,k ! indices
real(r8):: d,fd,dfdd ! temporary variable for frac_h2oscs iteration
real(r8):: sigma ! microtopography pdf sigma in mm
real(r8):: swc ! surface water content in m
real(r8):: min_h2osfc
!-----------------------------------------------------------------------

associate( &
micro_sigma => col_pp%micro_sigma , & ! Input: [real(r8) (:) ] microtopography pdf sigma (m)

h2osno => col_ws%h2osno , & ! Input: [real(r8) (:) ] snow water (mm H2O)
excess_ice => col_ws%excess_ice , & ! Input: [real(r8) (:,:) ] excess ground ice [kg/m2]

h2osoi_liq => col_ws%h2osoi_liq , & ! Output: [real(r8) (:,:) ] liquid water (col,lyr) [kg/m2]
h2osfc => col_ws%h2osfc , & ! Output: [real(r8) (:) ] surface water (mm)
frac_sno => col_ws%frac_sno , & ! Output: [real(r8) (:) ] fraction of ground covered by snow (0 to 1)
frac_sno_eff => col_ws%frac_sno_eff , & ! Output: [real(r8) (:) ] eff. fraction of ground covered by snow (0 to 1)
frac_h2osfc => col_ws%frac_h2osfc , & ! Output: [real(r8) (:) ] col fractional area with surface water greater than zero
frac_h2osfc_act => col_ws%frac_h2osfc_act & ! Output: [real(r8) (:) ] col fractional area with surface water greater than zero
associate( &
micro_sigma => col_pp%micro_sigma , & ! Input: [real(r8) (:) ] microtopography pdf sigma (m)

iwp_microrel => col_pp%iwp_microrel , & ! Input: [real(r8) (:) ] ice wedge polygon microtopographic relief (m)
iwp_exclvol => col_pp%iwp_exclvol , & ! Input: [real(r8) (:) ] microtopography pdf sigma (m)
iwp_ddep => col_pp%iwp_ddep , & ! Input: [real(r8) (:) ] microtopography pdf sigma (m)

h2osno => col_ws%h2osno , & ! Input: [real(r8) (:) ] snow water (mm H2O)

h2osoi_liq => col_ws%h2osoi_liq , & ! Output: [real(r8) (:,:) ] liquid water (col,lyr) [kg/m2]
h2osfc => col_ws%h2osfc , & ! Output: [real(r8) (:) ] surface water (mm)
frac_sno => col_ws%frac_sno , & ! Output: [real(r8) (:) ] fraction of ground covered by snow (0 to 1)
frac_sno_eff => col_ws%frac_sno_eff , & ! Output: [real(r8) (:) ] eff. fraction of ground covered by snow (0 to 1)
frac_h2osfc => col_ws%frac_h2osfc , & ! Output: [real(r8) (:) ] col fractional area with surface water greater than zero
frac_h2osfc_act => col_ws%frac_h2osfc_act & ! Output: [real(r8) (:) ] col fractional area with surface water greater than zero
)

! arbitrary lower limit on h2osfc for safer numerics...
Expand All @@ -825,15 +829,42 @@ subroutine FracH2OSfc(bounds, num_h2osfc, filter_h2osfc, &
if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop) then

! Use newton-raphson method to iteratively determine frac_h20sfc
! based on amount of surface water storage (h2osfc) and
! microtopography variability (micro_sigma)
! based on amount of surface water storage (h2osfc) and
! microtopography variability (micro_sigma) if nonpolygonal, or
! amount of surface water (h2osfc) and microtopographic attributes
! (iwp_microrel, iwp_exclvol) if polygonal tundra

if (h2osfc(c) > min_h2osfc) then
! a cutoff is needed for numerical reasons...(nonconvergence after 5 iterations)

if lun_pp%ispolygon(l) then
! calculate water depth and inundation fraction if column is polygonal
swc = h2osfc(c)/1000 ! convert to m

if swc > iwp_microrel - iwp_exclvol then
d = swc + iwp_exclvol
else
d = 0.0
do k=1,10
fd = (2_r8*iwp_exclvol - iwp_microrel) * (d/iwp_microrel)**3_r8 &
+ (2_r8*iwp_microrel - 3_r8*iwp_exclvol) * (d/iwp_microrel)**2_r8 &
- swc
dfdd = (3_r8/iwp_microrel) * (2_r8*iwp_exclvol - iwp_microrel) * (d/iwp_microrel)**2_r8 &
+ (2_r8/iwp_microrel) * (2_r8*iwp_microrel - 3_r8*iwp_exclvol) * (d/iwp_microrel)
d = d - fd/dfdd
enddo
endif
!-- update the submerged areal fraction using the new d value
frac_h2osfc(c) = (3_r8/iwp_microrel) * (2_r8*iwp_exclvol - iwp_microrel) * (d/iwp_microrel)**2_r8 &
+ (2_r8/iwp_microrel) * (2_r8*iwp_microrel - 3_r8*iwp_exclvol) * (d/iwp_microrel)
endif

else
! calculate water depth and inudation fraction if column is non-polygonal
d=0.0

sigma=1.0e3 * micro_sigma(c) ! convert to mm
do l=1,10
do k=1,10
fd = 0.5*d*(1.0_r8+erf(d/(sigma*sqrt(2.0)))) &
+sigma/sqrt(2.0*shr_const_pi)*exp(-d**2/(2.0*sigma**2)) &
-h2osfc(c)
Expand All @@ -843,6 +874,7 @@ subroutine FracH2OSfc(bounds, num_h2osfc, filter_h2osfc, &
enddo
!-- update the submerged areal fraction using the new d value
frac_h2osfc(c) = 0.5*(1.0_r8+erf(d/(sigma*sqrt(2.0))))
endif

else
frac_h2osfc(c) = 0._r8
Expand Down

0 comments on commit fc8f584

Please sign in to comment.