diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90
index 1f391bd67b1f..d75856411da6 100755
--- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90
+++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90
@@ -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
@@ -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...
@@ -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)
@@ -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