From f736477889046a905e2c0bb56eb843a1a7269771 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Sun, 14 Jun 2026 14:54:39 -0700 Subject: [PATCH 1/9] remove unused total canopy area variable --- biogeochem/EDCanopyStructureMod.F90 | 4 ---- 1 file changed, 4 deletions(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index d6ae9f52c8..85f803bd8e 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1352,14 +1352,12 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) type (fates_patch_type) , pointer :: currentPatch real(r8) :: bare_frac_area real(r8) :: total_patch_area - real(r8) :: total_canopy_area real(r8) :: total_patch_leaf_stem_area real(r8) :: weight ! Weighting for cohort variables in patch do s = 1,nsites total_patch_area = 0._r8 - total_canopy_area = 0._r8 bc_out(s)%canopy_fraction_pa(:) = 0._r8 bc_out(s)%dleaf_pa(:) = 0._r8 bc_out(s)%z0m_pa(:) = 0._r8 @@ -1457,8 +1455,6 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) total_patch_area = total_patch_area + bc_out(s)%canopy_fraction_pa(ifp) + bare_frac_area - total_canopy_area = total_canopy_area + bc_out(s)%canopy_fraction_pa(ifp) - bc_out(s)%nocomp_pft_label_pa(ifp) = currentPatch%nocomp_pft_label ! Calculate area indices for output boundary to HLM From 5c852bb730fb421590c9686e6052fd743609c4e2 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Sun, 14 Jun 2026 15:14:15 -0700 Subject: [PATCH 2/9] rename canopy_fraction_pa to patch_fraction This reflects that the output will be a fraction dependent on certain logic having to do with handling bareground --- main/FatesInterfaceMod.F90 | 4 ++-- main/FatesInterfaceTypesMod.F90 | 10 ++++++---- 2 files changed, 8 insertions(+), 6 deletions(-) diff --git a/main/FatesInterfaceMod.F90 b/main/FatesInterfaceMod.F90 index faf91b1da6..77e1419469 100644 --- a/main/FatesInterfaceMod.F90 +++ b/main/FatesInterfaceMod.F90 @@ -411,7 +411,7 @@ subroutine zero_bcs(fates,s) fates%bc_out(s)%dleaf_pa(:) = 0.0_r8 fates%bc_out(s)%nocomp_pft_label_pa(:) = 0 - fates%bc_out(s)%canopy_fraction_pa(:) = 0.0_r8 + fates%bc_out(s)%patch_fraction(:) = 0.0_r8 fates%bc_out(s)%frac_veg_nosno_alb_pa(:) = 0.0_r8 if (hlm_use_planthydro.eq.itrue) then @@ -744,7 +744,7 @@ subroutine allocate_bcout(bc_out, nlevsoil_in, nlevdecomp_in) allocate(bc_out%displa_pa(maxpatch_total)) allocate(bc_out%z0m_pa(maxpatch_total)) - allocate(bc_out%canopy_fraction_pa(maxpatch_total)) + allocate(bc_out%patch_fraction(maxpatch_total)) allocate(bc_out%frac_veg_nosno_alb_pa(maxpatch_total)) allocate(bc_out%nocomp_pft_label_pa(maxpatch_total)) diff --git a/main/FatesInterfaceTypesMod.F90 b/main/FatesInterfaceTypesMod.F90 index dcaf9cd972..a9b53ea970 100644 --- a/main/FatesInterfaceTypesMod.F90 +++ b/main/FatesInterfaceTypesMod.F90 @@ -766,10 +766,12 @@ module FatesInterfaceTypesMod real(r8), allocatable :: displa_pa(:) ! displacement height [m] real(r8), allocatable :: dleaf_pa(:) ! leaf characteristic dimension/width/diameter [m] - real(r8), allocatable :: canopy_fraction_pa(:) ! Area fraction of each patch in the site - ! Use most likely for weighting - ! This is currently the projected canopy - ! area of each patch [0-1] + real(r8), allocatable :: patch_fraction(:) ! Area fraction of each patch in the site + ! Use most likely for weighting + ! This is either the projected canopy + ! area of each patch [0-1] or the actual + ! are fraction that each patch contributes + ! to the site real(r8), allocatable :: frac_veg_nosno_alb_pa(:) ! This is not really a fraction ! this is actually binary based on if any From ac12b45942d3584110f1cd7737b5774785b2f910 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Sun, 14 Jun 2026 15:30:22 -0700 Subject: [PATCH 3/9] add logic to receive interstitial bareground flag from HLM --- main/FatesInterfaceMod.F90 | 7 +++++++ main/FatesInterfaceTypesMod.F90 | 3 +++ 2 files changed, 10 insertions(+) diff --git a/main/FatesInterfaceMod.F90 b/main/FatesInterfaceMod.F90 index 77e1419469..fe3ee9fa12 100644 --- a/main/FatesInterfaceMod.F90 +++ b/main/FatesInterfaceMod.F90 @@ -1574,6 +1574,7 @@ subroutine set_fates_ctrlparms(tag,ival,rval,cval) hlm_use_fixed_biogeog = unset_int hlm_use_nocomp = unset_int hlm_use_sp = unset_int + hlm_use_interstitial_bareground = unset_int hlm_use_inventory_init = unset_int hlm_use_dbh_init = unset_int hlm_inventory_ctrl_file = 'unset' @@ -2050,6 +2051,12 @@ subroutine set_fates_ctrlparms(tag,ival,rval,cval) write(fates_log(),*) 'Transfering hlm_use_sp= ',ival,' to FATES' end if + case('use_fates_interstitial_bareground') + hlm_use_interstitial_bareground = ival + if (fates_global_verbose()) then + write(fates_log(),*) 'Transfering hlm_use_interstitial_bareground= ',ival,' to FATES' + end if + case('use_planthydro') hlm_use_planthydro = ival if (fates_global_verbose()) then diff --git a/main/FatesInterfaceTypesMod.F90 b/main/FatesInterfaceTypesMod.F90 index a9b53ea970..0a3d92c1a2 100644 --- a/main/FatesInterfaceTypesMod.F90 +++ b/main/FatesInterfaceTypesMod.F90 @@ -225,6 +225,9 @@ module FatesInterfaceTypesMod integer, public :: hlm_use_sp ! Flag to use FATES satellite phenology (LAI) mode ! 1 = TRUE, 0 = FALSE + integer, public :: hlm_use_interstitial_bareground ! Flag to use count the FATES interstitial bareground area + ! with the HLM bareground patch + ! 1 = TRUE, 0 = FALSE ! Flag specifying what types of history fields to allocate and prepare ! The "_dynam" refers to history fields that can be updated on the dynamics (daily) step From 921589348ccaf9ce552f70ca63fe632888001d61 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Sun, 14 Jun 2026 15:35:55 -0700 Subject: [PATCH 4/9] update the boundary condition output for patch fraction --- biogeochem/EDCanopyStructureMod.F90 | 27 +++++++++++++++++++-------- 1 file changed, 19 insertions(+), 8 deletions(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 85f803bd8e..a4aad13f39 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -36,6 +36,7 @@ module EDCanopyStructureMod use FatesInterfaceTypesMod , only : hlm_use_planthydro use FatesInterfaceTypesMod , only : hlm_use_cohort_age_tracking use FatesInterfaceTypesMod , only : hlm_use_sp + use FatesInterfaceTypesMod , only : hlm_use_interstitial_bareground use FatesInterfaceTypesMod , only : numpft use FatesInterfaceTypesMod, only : bc_in_type use FatesPlantHydraulicsMod, only : UpdateH2OVeg,InitHydrCohort, RecruitWaterStorage @@ -1351,6 +1352,7 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) integer :: s, ifp, c, p type (fates_patch_type) , pointer :: currentPatch real(r8) :: bare_frac_area + real(r8) :: local_patch_fraction real(r8) :: total_patch_area real(r8) :: total_patch_leaf_stem_area real(r8) :: weight ! Weighting for cohort variables in patch @@ -1358,7 +1360,8 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) do s = 1,nsites total_patch_area = 0._r8 - bc_out(s)%canopy_fraction_pa(:) = 0._r8 + local_patch_fraction = 0._r8 + bc_out(s)%patch_fraction(:) = 0._r8 bc_out(s)%dleaf_pa(:) = 0._r8 bc_out(s)%z0m_pa(:) = 0._r8 bc_out(s)%displa_pa(:) = 0._r8 @@ -1443,17 +1446,25 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) ! currentPatch%total_canopy_area/currentPatch%area is fraction of this patch cover by plants ! currentPatch%area/AREA is the fraction of the soil covered by this patch. + ! If we are accounting for the interstitial bareground within the HLM bareground patch + ! we use the projected total canopy area fraction as the output boundary condition. + ! Otherwise we simply use the fraction of the area that the patch contributes to the given site. + if (hlm_use_interstitial_bareground == itrue) then + local_patch_fraction = min(1.0_r8,currentPatch%total_canopy_area/currentPatch%area) + else + local_patch_fraction = 1.0_r8 + end if + if(currentPatch%area.gt.0.0_r8)then - bc_out(s)%canopy_fraction_pa(ifp) = & - min(1.0_r8,currentPatch%total_canopy_area/currentPatch%area)*(currentPatch%area/AREA) + bc_out(s)%patch_fraction(ifp) = local_patch_fraction*(currentPatch%area/AREA) else - bc_out(s)%canopy_fraction_pa(ifp) = 0.0_r8 + bc_out(s)%patch_fraction(ifp) = 0.0_r8 endif - bare_frac_area = (1.0_r8 - min(1.0_r8,currentPatch%total_canopy_area/currentPatch%area)) * & + bare_frac_area = (1.0_r8 - local_patch_fraction) * & (currentPatch%area/AREA) - total_patch_area = total_patch_area + bc_out(s)%canopy_fraction_pa(ifp) + bare_frac_area + total_patch_area = total_patch_area + bc_out(s)%patch_fraction(ifp) + bare_frac_area bc_out(s)%nocomp_pft_label_pa(ifp) = currentPatch%nocomp_pft_label @@ -1480,7 +1491,7 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) bc_out(s)%frac_veg_nosno_alb_pa(ifp) = 0.0_r8 end if - else ! nocomp or SP, and currentPatch%nocomp_pft_label .eq. 0 + else ! nocomp or SP, (i.e.currentPatch%nocomp_pft_label .eq. 0 total_patch_area = total_patch_area + currentPatch%area/AREA @@ -1507,7 +1518,7 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) do while(associated(currentPatch)) ifp = currentPatch%patchno if(currentPatch%nocomp_pft_label.ne.nocomp_bareground)then ! for vegetated patches only - bc_out(s)%canopy_fraction_pa(ifp) = bc_out(s)%canopy_fraction_pa(ifp)/total_patch_area + bc_out(s)%patch_fraction(ifp) = bc_out(s)%patch_fraction(ifp)/total_patch_area endif ! veg patch currentPatch => currentPatch%younger end do From c20ceb38cd3f09732aad2248ff06266803451310 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Fri, 26 Jun 2026 16:13:38 -0700 Subject: [PATCH 5/9] add interstitial bareground logic to handle surface roughness and displacement height --- biogeochem/EDCanopyStructureMod.F90 | 32 +++++++++++++++++++++-------- 1 file changed, 23 insertions(+), 9 deletions(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index a4aad13f39..c597529a29 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1355,8 +1355,9 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) real(r8) :: local_patch_fraction real(r8) :: total_patch_area real(r8) :: total_patch_leaf_stem_area - real(r8) :: weight ! Weighting for cohort variables in patch - + real(r8) :: weight ! Weighting for cohort variables in patch + real(r8) :: weighting_area ! Area to normalize against depending on insterstitial bareground handling + do s = 1,nsites total_patch_area = 0._r8 @@ -1390,17 +1391,27 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) bc_out(s)%hbot_pa(ifp) = max(0._r8, min(0.2_r8, bc_out(s)%htop_pa(ifp)- 1.0_r8)) - ! Use canopy-only crown area weighting for all cohorts in the patch to define the characteristic - ! Roughness length and displacement height used by the HLM - ! use total LAI + SAI to weight the leaft characteristic dimension + ! Define the characteristic roughness length and displacement height used by the HLM + ! Nominally, FATES uses canopy-only crown area weighting for all cohorts in the patch + ! when the interstitial bareground is being collapsed into the HLM bareground patch. + ! When we account for the bareground as part of the patch, the full area is accounted + ! for in the calculation. + ! If there is no canopy and we are considering the interstitial bareground we do + ! something else. + ! Use total LAI + SAI to weight the leaft characteristic dimension ! Avoid this if running in satellite phenology mode ! ---------------------------------------------------------------------------- if (currentPatch%total_canopy_area > nearzero) then + if (hlm_use_interstitial_bareground .eq. itrue) then + weighting_area = currentPatch%total_canopy_area + else + weighting_area = currentPatch%area + end if currentCohort => currentPatch%shortest do while(associated(currentCohort)) if (currentCohort%canopy_layer .eq. 1) then - weight = min(1.0_r8,currentCohort%c_area/currentPatch%total_canopy_area) + weight = min(1.0_r8,currentCohort%c_area/weighting_area) bc_out(s)%z0m_pa(ifp) = bc_out(s)%z0m_pa(ifp) + & EDPftvarcon_inst%z0mr(currentCohort%pft) * currentCohort%height * weight bc_out(s)%displa_pa(ifp) = bc_out(s)%displa_pa(ifp) + & @@ -1434,9 +1445,12 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) endif else ! if no canopy, then use dummy values (first PFT) of aerodynamic properties - bc_out(s)%z0m_pa(ifp) = EDPftvarcon_inst%z0mr(1) * bc_out(s)%htop_pa(ifp) - bc_out(s)%displa_pa(ifp) = EDPftvarcon_inst%displar(1) * bc_out(s)%htop_pa(ifp) - bc_out(s)%dleaf_pa(ifp) = EDPftvarcon_inst%dleaf(1) + if (hlm_use_interstitial_bareground .eq. itrue) then + bc_out(s)%z0m_pa(ifp) = EDPftvarcon_inst%z0mr(1) * bc_out(s)%htop_pa(ifp) + bc_out(s)%displa_pa(ifp) = EDPftvarcon_inst%displar(1) * bc_out(s)%htop_pa(ifp) + bc_out(s)%dleaf_pa(ifp) = EDPftvarcon_inst%dleaf(1) + ! else do something else if interstitial bareground is not being used by the HLM + end if endif ! ----------------------------------------------------------------------------- From 08377484522ad8ce718233e1f80a86cc0883e9d5 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Mon, 10 Aug 2026 18:04:17 -0600 Subject: [PATCH 6/9] move weighting check outside of TCA check --- biogeochem/EDCanopyStructureMod.F90 | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index c597529a29..7a32254de1 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1402,12 +1402,15 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) ! Avoid this if running in satellite phenology mode ! ---------------------------------------------------------------------------- + ! Set the are to be used for weighting against the cohort area + ! If we are not using the interstitial bareground paradigm, simply + ! use the patch area. Otherwise use the total canopy area. + weighting_area = currentPatch%total_canopy_area + if (hlm_use_interstitial_bareground .eq. ifalse) then + weighting_area = currentPatch%area + end if + if (currentPatch%total_canopy_area > nearzero) then - if (hlm_use_interstitial_bareground .eq. itrue) then - weighting_area = currentPatch%total_canopy_area - else - weighting_area = currentPatch%area - end if currentCohort => currentPatch%shortest do while(associated(currentCohort)) if (currentCohort%canopy_layer .eq. 1) then From 8c2b11c3dc8fbd30a8a78dfad3667a54bcce562e Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Mon, 10 Aug 2026 18:05:29 -0600 Subject: [PATCH 7/9] remove interstitial check from zero TCA option These are effectively dummy values that are necessary to avoid NaN calcs in on the HLM side --- biogeochem/EDCanopyStructureMod.F90 | 3 --- 1 file changed, 3 deletions(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 7a32254de1..7799e6ad10 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1448,12 +1448,9 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) endif else ! if no canopy, then use dummy values (first PFT) of aerodynamic properties - if (hlm_use_interstitial_bareground .eq. itrue) then bc_out(s)%z0m_pa(ifp) = EDPftvarcon_inst%z0mr(1) * bc_out(s)%htop_pa(ifp) bc_out(s)%displa_pa(ifp) = EDPftvarcon_inst%displar(1) * bc_out(s)%htop_pa(ifp) bc_out(s)%dleaf_pa(ifp) = EDPftvarcon_inst%dleaf(1) - ! else do something else if interstitial bareground is not being used by the HLM - end if endif ! ----------------------------------------------------------------------------- From fdbca105d8039b3cfe6a2e60c377c0b33453d81e Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Mon, 24 Aug 2026 17:10:40 -0600 Subject: [PATCH 8/9] add input bc for bareground roughness from the HLM --- main/FatesInterfaceTypesMod.F90 | 3 +++ 1 file changed, 3 insertions(+) diff --git a/main/FatesInterfaceTypesMod.F90 b/main/FatesInterfaceTypesMod.F90 index 0a3d92c1a2..5559e03b24 100644 --- a/main/FatesInterfaceTypesMod.F90 +++ b/main/FatesInterfaceTypesMod.F90 @@ -506,6 +506,9 @@ module FatesInterfaceTypesMod ! soil temperature (Kelvin) real(r8), allocatable :: t_soisno_sl(:) + + ! surface roughness length, site bareground (m) + real(r8) :: z0mg ! Canopy Radiation Boundaries ! --------------------------------------------------------------------------------- From 4045b454ce4802555b7949ae1ed9c196ab5601d2 Mon Sep 17 00:00:00 2001 From: Gregory Lemieux Date: Mon, 24 Aug 2026 17:11:07 -0600 Subject: [PATCH 9/9] update roughness length calculation to include bareground --- biogeochem/EDCanopyStructureMod.F90 | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 7799e6ad10..bb67553699 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1329,7 +1329,7 @@ end subroutine leaf_area_profile ! ====================================================================================== - subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) + subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out,bc_in) ! ---------------------------------------------------------------------------------- ! The purpose of this routine is to package output boundary conditions related @@ -1346,6 +1346,7 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) type(ed_site_type), intent(inout), target :: sites(nsites) integer, intent(in) :: fcolumn(nsites) type(bc_out_type), intent(inout) :: bc_out(nsites) + type(bc_in_type), intent(in) :: bc_in(nsites) ! Locals type (fates_cohort_type) , pointer :: currentCohort @@ -1357,6 +1358,7 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) real(r8) :: total_patch_leaf_stem_area real(r8) :: weight ! Weighting for cohort variables in patch real(r8) :: weighting_area ! Area to normalize against depending on insterstitial bareground handling + real(r8) :: bareground_fraction ! fraction of the patch that is bareground do s = 1,nsites @@ -1422,6 +1424,15 @@ subroutine update_hlm_dynamics(nsites,sites,fcolumn,bc_out) endif currentCohort => currentCohort%taller end do + + ! If we are not using the interstitial bareground, we need to account for the bareground area + ! fraction in calculating the roughness length. Here we mimic ZengWang 2007. This should + ! eventually be replaced with a more robust approach that more accurately accounts for the + ! the bareground fraction of the patch. + if (hlm_use_interstitial_bareground .eq. ifalse) then + bareground_fraction = (currentPatch%area - currentPatch%total_canopy_area) / currentPatch%area + bc_out(s)%z0m_pa(ifp) = exp(log(bc_out(s)%z0m_pa(ifp)) + bareground_fraction * log(bc_in(s)%z0mg)) + end if ! for lai, scale to total LAI + SAI in patch. first add up all the LAI and SAI in the patch total_patch_leaf_stem_area = 0._r8