From baadc328fa0e05021a955645aac4f7933c570f02 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Fri, 9 May 2025 14:41:48 -0700 Subject: [PATCH 01/10] added biophys vars by land use: LH, SH, LW, SW, TV, TSA --- main/FatesHistoryInterfaceMod.F90 | 183 ++++++++++++++++++++++++++++++ main/FatesInterfaceMod.F90 | 12 ++ main/FatesInterfaceTypesMod.F90 | 7 ++ 3 files changed, 202 insertions(+) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index f7785cc3c7..5c2f215d2e 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -9,6 +9,7 @@ module FatesHistoryInterfaceMod use FatesConstantsMod , only : mg_per_kg use FatesConstantsMod , only : pi_const use FatesConstantsMod , only : nearzero + use FatesConstantsMod , only : rsnbl_math_prec use FatesConstantsMod , only : t_water_freeze_k_1atm use FatesConstantsMod , only : n_term_mort_types use FatesConstantsMod , only : i_term_mort_type_cstarv @@ -138,6 +139,7 @@ module FatesHistoryInterfaceMod use FatesSizeAgeTypeIndicesMod, only : get_cdamagesize_class_index use FatesSizeAgeTypeIndicesMod, only : get_cdamagesizepft_class_index use FatesSizeAgeTypeIndicesMod, only : coagetype_class_index + use FatesInterfaceTypesMod , only : hlm_use_luh implicit none private ! By default everything is private @@ -362,6 +364,12 @@ module FatesHistoryInterfaceMod integer :: ih_burnedarea_si_landuse integer :: ih_gpp_si_landuse integer :: ih_npp_si_landuse + integer :: ih_tveg_si_landuse + integer :: ih_tsa_si_landuse + integer :: ih_sw_abs_si_landuse + integer :: ih_lw_net_si_landuse + integer :: ih_shflux_si_landuse + integer :: ih_lhflux_si_landuse ! land use by land use variables integer :: ih_disturbance_rate_si_lulu @@ -824,6 +832,7 @@ module FatesHistoryInterfaceMod procedure :: update_history_dyn2 procedure :: update_history_hifrq procedure :: update_history_hifrq1 + procedure :: update_history_hifrq_landuse procedure :: update_history_hifrq2 procedure :: update_history_hydraulics procedure :: update_history_nutrflux @@ -4897,6 +4906,9 @@ subroutine update_history_hifrq(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) if(hlm_hist_level_hifrq>0) then call update_history_hifrq1(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) + if (hlm_use_luh .eq. itrue) then + call update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) + end if if(hlm_hist_level_hifrq>1) then call update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) end if @@ -5124,6 +5136,138 @@ end subroutine update_history_hifrq1 ! =============================================================================================== + subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) + + ! + ! Arguments + class(fates_history_interface_type) :: this + integer , intent(in) :: nc ! clump index + integer , intent(in) :: nsites + type(ed_site_type) , intent(inout), target :: sites(nsites) + type(bc_in_type) , intent(in) :: bc_in(nsites) + + ! Locals + integer :: s ! The local site index + integer :: io_si ! The site index of the IO array + + real(r8) :: landuse_statevector(n_landuse_cats) + real(r8) :: canopy_area_bylanduse(n_landuse_cats) + integer :: i_lu + logical :: foundbaregroundpatch + + type(fates_patch_type),pointer :: cpatch + + associate( hio_tveg_si_landuse => this%hvars(ih_tveg_si_landuse)%r82d,& + hio_tsa_si_landuse => this%hvars(ih_tsa_si_landuse)%r82d,& + hio_sw_abs_si_landuse => this%hvars(ih_sw_abs_si_landuse)%r82d,& + hio_lw_net_si_landuse => this%hvars(ih_lw_net_si_landuse)%r82d,& + hio_shflux_si_landuse => this%hvars(ih_shflux_si_landuse)%r82d,& + hio_lhflux_si_landuse => this%hvars(ih_lhflux_si_landuse)%r82d) + + do_sites: do s = 1,nsites + + io_si = sites(s)%h_gid + + ! biophysical properties that are indexed by land use + landuse_statevector(:) = sites(s)%get_current_landuse_statevector() * AREA + + ! get the total canopy area for each land use type + canopy_area_bylanduse(:) = 0._r8 + cpatch => sites(s)%oldest_patch + do while(associated(cpatch)) + if ( cpatch%land_use_label .ne. nocomp_bareground_land ) then + canopy_area_bylanduse(cpatch%land_use_label) = canopy_area_bylanduse(cpatch%land_use_label) + & + cpatch%total_canopy_area + endif + cpatch => cpatch%younger + end do + + cpatch => sites(s)%oldest_patch + do while(associated(cpatch)) + if (cpatch%total_canopy_area .gt. rsnbl_math_prec) then + ! for TVEG, since it is only defined on vegetated area of vegetated patches, normalize by the total vegetated area + hio_tveg_si_landuse(io_si,cpatch%land_use_label) = hio_tveg_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%t_veg_pa(cpatch%patchno) * cpatch%total_canopy_area/canopy_area_bylanduse(cpatch%land_use_label) + + ! for the rest of these, first weight by the vegetated area of each patch over the total patch area for each land use type + hio_tsa_si_landuse(io_si,cpatch%land_use_label) = hio_tsa_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%t2m_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + + hio_sw_abs_si_landuse(io_si,cpatch%land_use_label) = hio_sw_abs_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%swabs_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + + hio_lw_net_si_landuse(io_si,cpatch%land_use_label) = hio_lw_net_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%netlw_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + + hio_shflux_si_landuse(io_si,cpatch%land_use_label) = hio_shflux_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%shflux_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + + hio_lhflux_si_landuse(io_si,cpatch%land_use_label) = hio_lhflux_si_landuse(io_si,cpatch%land_use_label) + & + bc_in(s)%lhflux_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + endif + cpatch => cpatch%younger + end do + + ! for all the land-use indexed variables, excpet for TVEG, also add in the component for the unvegetated area of each land use + foundbaregroundpatch = .false. + cpatch => sites(s)%oldest_patch + do while(associated(cpatch)) + if (cpatch%land_use_label .eq. nocomp_bareground_land .and. .not. foundbaregroundpatch) then + foundbaregroundpatch = .true. + do i_lu = 1, n_landuse_cats + if ( landuse_statevector(i_lu) .gt. rsnbl_math_prec ) then + hio_tsa_si_landuse(io_si,i_lu) = hio_tsa_si_landuse(io_si,i_lu) + & + bc_in(s)%t2m_pa(cpatch%patchno) * & + (landuse_statevector(i_lu) - canopy_area_bylanduse(i_lu)) / landuse_statevector(i_lu) + + hio_sw_abs_si_landuse(io_si,i_lu) = hio_sw_abs_si_landuse(io_si,i_lu) + & + bc_in(s)%swabs_pa(cpatch%patchno) * & + (landuse_statevector(i_lu) - canopy_area_bylanduse(i_lu)) / landuse_statevector(i_lu) + + hio_lw_net_si_landuse(io_si,i_lu) = hio_lw_net_si_landuse(io_si,i_lu) + & + bc_in(s)%netlw_pa(cpatch%patchno) * & + (landuse_statevector(i_lu) - canopy_area_bylanduse(i_lu)) / landuse_statevector(i_lu) + + hio_shflux_si_landuse(io_si,i_lu) = hio_shflux_si_landuse(io_si,i_lu) + & + bc_in(s)%shflux_pa(cpatch%patchno) * & + (landuse_statevector(i_lu) - canopy_area_bylanduse(i_lu)) / landuse_statevector(i_lu) + + hio_lhflux_si_landuse(io_si,i_lu) = hio_lhflux_si_landuse(io_si,i_lu) + & + bc_in(s)%lhflux_pa(cpatch%patchno) * & + (landuse_statevector(i_lu) - canopy_area_bylanduse(i_lu)) / landuse_statevector(i_lu) + end if + end do + end if + cpatch => cpatch%younger + end do + + ! instead of leaving the values for unoccupied areas as zero, set as missing values + do i_lu = 1, n_landuse_cats + + ! if a given land use type is not present, set the value as missing + if ( landuse_statevector(i_lu) .le. rsnbl_math_prec ) then + hio_tsa_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + hio_sw_abs_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + hio_lw_net_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + hio_shflux_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + hio_lhflux_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + end if + + ! for tveg, ignore if there is no vegetation present on any patches of a given land use type + if ( canopy_area_bylanduse(i_lu) .le. rsnbl_math_prec ) then + hio_tveg_si_landuse(io_si,i_lu) = hlm_hio_ignore_val + end if + + end do + + end do do_sites + + end associate + return + end subroutine update_history_hifrq_landuse + + ! =============================================================================================== + subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ! --------------------------------------------------------------------------------- @@ -8507,6 +8651,45 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & ivar=ivar, initialize=initialize_variables, index = ih_tveg_si ) + if (hlm_use_luh .eq. itrue) then + ! biophysics variables that are indexed by land use type + call this%set_history_var(vname='FATES_TVEG_LU', units='degrees Kelvin', & + long='fates instantaneous mean vegetation temperature by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_tveg_si_landuse ) + + call this%set_history_var(vname='FATES_TSA_LU', units='degrees Kelvin', & + long='fates instantaneous mean near-surface (2m) air temperature by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_tsa_si_landuse ) + + call this%set_history_var(vname='FATES_SWABS_LU', units='W m-2', & + long='fates absorbed shortwave radiation by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_sw_abs_si_landuse ) + + call this%set_history_var(vname='FATES_NETLW_LU', units='W m-2', & + long='fates net longwave flux by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_lw_net_si_landuse ) + + call this%set_history_var(vname='FATES_SHFLUX_LU', units='W m-2', & + long='fates sensible heat flux by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_shflux_si_landuse ) + + call this%set_history_var(vname='FATES_LHFLUX_LU', units='W m-2', & + long='fates latent heat flux by land use type', & + use_default='active', & + avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & + ivar=ivar, initialize=initialize_variables, index = ih_lhflux_si_landuse ) + endif + call this%set_history_var(vname='FATES_VIS_RAD_ERROR', units='-', & long='mean two-stream solver error for VIS', use_default='active', & avgflag='A', vtype=site_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & diff --git a/main/FatesInterfaceMod.F90 b/main/FatesInterfaceMod.F90 index 555fa526c1..17557d6964 100644 --- a/main/FatesInterfaceMod.F90 +++ b/main/FatesInterfaceMod.F90 @@ -321,6 +321,12 @@ subroutine zero_bcs(fates,s) fates%bc_in(s)%hksat_sisl(:) = 0.0_r8 end if + ! Diagnostic quantities for outputting FATES patch-resolved + fates%bc_in(s)%lhflux_pa(:) = 0._r8 + fates%bc_in(s)%shflux_pa(:) = 0._r8 + fates%bc_in(s)%swabs_pa(:) = 0._r8 + fates%bc_in(s)%netlw_pa(:) = 0._r8 + fates%bc_in(s)%t2m_pa(:) = 0._r8 ! Output boundaries fates%bc_out(s)%active_suction_sl(:) = .false. @@ -552,6 +558,12 @@ subroutine allocate_bcin(bc_in, nlevsoil_in, nlevdecomp_in, num_lu_harvest_cats, end if ! Land use + ! Diagnostic quantities for outputting FATES patch-resolved + allocate(bc_in%lhflux_pa(0:maxpatch_total)) + allocate(bc_in%shflux_pa(0:maxpatch_total)) + allocate(bc_in%swabs_pa(0:maxpatch_total)) + allocate(bc_in%netlw_pa(0:maxpatch_total)) + allocate(bc_in%t2m_pa(0:maxpatch_total)) ! harvest flag denote data from hlm, ! while the logging flag signifies only that logging is occurring (which could just be FATES logging) diff --git a/main/FatesInterfaceTypesMod.F90 b/main/FatesInterfaceTypesMod.F90 index 416af2728a..a49cf824bd 100644 --- a/main/FatesInterfaceTypesMod.F90 +++ b/main/FatesInterfaceTypesMod.F90 @@ -593,6 +593,13 @@ module FatesInterfaceTypesMod real(r8),allocatable :: h2o_liq_sisl(:) ! Liquid water mass in each layer (kg/m2) real(r8) :: smpmin_si ! restriction for min of soil potential (mm) + ! Diagnostic quantities for outputting FATES patch-resolved + real(r8),allocatable :: lhflux_pa(:) ! latent heat flux + real(r8),allocatable :: shflux_pa(:) ! sensible heat flux + real(r8),allocatable :: swabs_pa(:) ! shortwave absorbed radiation + real(r8),allocatable :: netlw_pa(:) ! longwave net radiation + real(r8),allocatable :: t2m_pa(:) ! 2m air temeprature + ! Land use ! --------------------------------------------------------------------------------- real(r8),allocatable :: hlm_harvest_rates(:) ! annual harvest rate per cat from hlm for a site From 125c01b2932a1c198013c95b095cd36b06de158a Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Fri, 9 May 2025 15:52:56 -0700 Subject: [PATCH 02/10] moved FATES_GPP_LU hist var to the hifrq landuse call --- main/FatesHistoryInterfaceMod.F90 | 44 +++++++++++++++++++++---------- 1 file changed, 30 insertions(+), 14 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 701a68318a..65fa6c841a 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -5017,7 +5017,7 @@ subroutine update_history_hifrq(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) if(hlm_hist_level_hifrq>0) then call update_history_hifrq_sitelevel(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) if (hlm_use_luh .eq. itrue) then - call update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) + call update_history_hifrq_landuse(this,nc,nsites,sites,bc_in,dt_tstep) end if if(hlm_hist_level_hifrq>1) then call update_history_hifrq_subsite(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) @@ -5252,7 +5252,7 @@ end subroutine update_history_hifrq_sitelevel ! =============================================================================================== - subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) + subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in,dt_tstep) ! ! Arguments @@ -5261,6 +5261,7 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) integer , intent(in) :: nsites type(ed_site_type) , intent(inout), target :: sites(nsites) type(bc_in_type) , intent(in) :: bc_in(nsites) + real(r8) , intent(in) :: dt_tstep ! Locals integer :: s ! The local site index @@ -5272,8 +5273,11 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) logical :: foundbaregroundpatch type(fates_patch_type),pointer :: cpatch + type(fates_cohort_type),pointer :: ccohort + real(r8) :: dt_tstep_inv ! Time step in frequency units (/s) associate( hio_tveg_si_landuse => this%hvars(ih_tveg_si_landuse)%r82d,& + hio_gpp_si_landuse => this%hvars(ih_gpp_si_landuse)%r82d, & hio_tsa_si_landuse => this%hvars(ih_tsa_si_landuse)%r82d,& hio_sw_abs_si_landuse => this%hvars(ih_sw_abs_si_landuse)%r82d,& hio_lw_net_si_landuse => this%hvars(ih_lw_net_si_landuse)%r82d,& @@ -5284,6 +5288,8 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) io_si = sites(s)%h_gid + dt_tstep_inv = 1.0_r8/dt_tstep + ! biophysical properties that are indexed by land use landuse_statevector(:) = sites(s)%get_current_landuse_statevector() * AREA @@ -5376,6 +5382,22 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in) end do + + cpatch => sites(s)%oldest_patch + do while(associated(cpatch)) + ccohort => cpatch%shortest + do while(associated(ccohort)) + if ( .not. ccohort%isnew ) then + if (cpatch%land_use_label .gt. nocomp_bareground_land) then + hio_gpp_si_landuse(io_si,cpatch%land_use_label) = hio_gpp_si_landuse(io_si,cpatch%land_use_label) & + + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv + end if + end if + ccohort => ccohort%taller + end do + cpatch => cpatch%younger + end do + end do do_sites end associate @@ -5442,7 +5464,6 @@ subroutine update_history_hifrq_subsite(this,nc,nsites,sites,bc_in,bc_out,dt_tst hio_froot_mr_understory_si_scls => this%hvars(ih_froot_mr_understory_si_scls)%r82d, & hio_resp_g_understory_si_scls => this%hvars(ih_resp_g_understory_si_scls)%r82d, & hio_resp_m_understory_si_scls => this%hvars(ih_resp_m_understory_si_scls)%r82d, & - hio_gpp_si_landuse => this%hvars(ih_gpp_si_landuse)%r82d, & hio_parsun_z_si_cnlf => this%hvars(ih_parsun_z_si_cnlf)%r82d, & hio_parsha_z_si_cnlf => this%hvars(ih_parsha_z_si_cnlf)%r82d, & hio_ts_net_uptake_si_cnlf => this%hvars(ih_ts_net_uptake_si_cnlf)%r82d, & @@ -5531,11 +5552,6 @@ subroutine update_history_hifrq_subsite(this,nc,nsites,sites,bc_in,bc_out,dt_tst hio_ar_frootm_si_scpf(io_si,scpf) = hio_ar_frootm_si_scpf(io_si,scpf) + & ccohort%froot_mr * n_perm2 - if (cpatch%land_use_label .gt. nocomp_bareground_land) then - hio_gpp_si_landuse(io_si,cpatch%land_use_label) = hio_gpp_si_landuse(io_si,cpatch%land_use_label) & - + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv - end if - ! accumulate fluxes on canopy- and understory- separated fluxes if (ccohort%canopy_layer .eq. 1) then @@ -8924,6 +8940,12 @@ subroutine define_history_vars(this, initialize_variables) use_default='active', & avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_hifr_simple, & ivar=ivar, initialize=initialize_variables, index = ih_lhflux_si_landuse ) + + call this%set_history_var(vname='FATES_GPP_LU', units='kg m-2 s-1', & + long='gross primary productivity by age bin in kg carbon per m2 per second', & + use_default='inactive', avgflag='A', vtype=site_landuse_r8, & + hlms='CLM:ALM', upfreq=group_hifr_simple, ivar=ivar, initialize=initialize_variables, & + index = ih_gpp_si_landuse) endif call this%set_history_var(vname='FATES_VIS_RAD_ERROR', units='-', & @@ -9085,12 +9107,6 @@ subroutine define_history_vars(this, initialize_variables) hlms='CLM:ALM', upfreq=group_hifr_complx, ivar=ivar, initialize=initialize_variables, & index = ih_gpp_si_age) - call this%set_history_var(vname='FATES_GPP_LU', units='kg m-2 s-1', & - long='gross primary productivity by land use type in kg carbon per m2 per second', & - use_default='inactive', avgflag='A', vtype=site_landuse_r8, & - hlms='CLM:ALM', upfreq=group_hifr_complx, ivar=ivar, initialize=initialize_variables, & - index = ih_gpp_si_landuse) - call this%set_history_var(vname='FATES_RDARK_USTORY_SZ', & units = 'kg m-2 s-1', & long='dark respiration for understory plants in kg carbon per m2 per second by size', & From 4fb1c648514eb08d89e92a74ce92fbeae06a73f1 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Mon, 14 Jul 2025 09:49:03 -0700 Subject: [PATCH 03/10] resolve ecosystem-level fire fluxes by dist type, and add hist var for LUC-associated fire --- biogeochem/EDPatchDynamicsMod.F90 | 26 ++++++++++++++------------ main/ChecksBalancesMod.F90 | 2 +- main/EDMainMod.F90 | 4 ++-- main/EDTypesMod.F90 | 4 ++-- main/FatesHistoryInterfaceMod.F90 | 17 +++++++++++++++-- 5 files changed, 34 insertions(+), 19 deletions(-) diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index a513f7717b..62a21ae5f0 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -1134,8 +1134,8 @@ subroutine spawn_patches( currentSite, bc_in ) endif - currentSite%mass_balance(el)%burn_flux_to_atm = & - currentSite%mass_balance(el)%burn_flux_to_atm + & + currentSite%mass_balance(el)%burn_flux_to_atm(i_disturbance_type) = & + currentSite%mass_balance(el)%burn_flux_to_atm(i_disturbance_type) + & leaf_burn_frac * leaf_m * nc%n ! This term increments the loss flux from surviving trees @@ -2042,7 +2042,7 @@ subroutine TransLitterNewPatch(currentSite, & if (debug) then - burn_flux0 = site_mass%burn_flux_to_atm + burn_flux0 = site_mass%burn_flux_to_atm(dist_type) litter_stock0 = curr_litt%GetTotalLitterMass()*currentPatch%area + & new_litt%GetTotalLitterMass()*newPatch%area end if @@ -2063,7 +2063,7 @@ subroutine TransLitterNewPatch(currentSite, & new_litt%ag_cwd(c) = new_litt%ag_cwd(c) + donatable_mass*donate_m2 curr_litt%ag_cwd(c) = curr_litt%ag_cwd(c) + donatable_mass*retain_m2 - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dist_type) = site_mass%burn_flux_to_atm(dist_type) + burned_mass ! Transfer below ground CWD (none burns) do sl = 1,currentSite%nlevsoil @@ -2091,7 +2091,7 @@ subroutine TransLitterNewPatch(currentSite, & new_litt%leaf_fines(dcmpy) = new_litt%leaf_fines(dcmpy) + donatable_mass*donate_m2 curr_litt%leaf_fines(dcmpy) = curr_litt%leaf_fines(dcmpy) + donatable_mass*retain_m2 - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dist_type) = site_mass%burn_flux_to_atm(dist_type) + burned_mass ! Transfer root fines (none burns) do sl = 1,currentSite%nlevsoil @@ -2122,7 +2122,7 @@ subroutine TransLitterNewPatch(currentSite, & ! EDMainMod start triggering. ! -------------------------------------------------------------------------- if (debug) then - burn_flux1 = site_mass%burn_flux_to_atm + burn_flux1 = site_mass%burn_flux_to_atm(dist_type) litter_stock1 = curr_litt%GetTotalLitterMass()*remainder_area + & new_litt%GetTotalLitterMass()*newPatch%area error = (litter_stock1 - litter_stock0) + (burn_flux1-burn_flux0) @@ -2301,7 +2301,7 @@ subroutine fire_litter_fluxes(currentSite, currentPatch, & donatable_mass*retain_m2*dcmpy_frac end do - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dtype_ifire) = site_mass%burn_flux_to_atm(dtype_ifire) + burned_mass call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & bc_in%max_rooting_depth_index_col) @@ -2363,7 +2363,7 @@ subroutine fire_litter_fluxes(currentSite, currentPatch, & donatable_mass = donatable_mass * (1.0_r8-currentCohort%fraction_crown_burned) burned_mass = num_dead_trees * SF_val_CWD_frac_adj(c) * bstem * & currentCohort%fraction_crown_burned - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dtype_ifire) = site_mass%burn_flux_to_atm(dtype_ifire) + burned_mass endif new_litt%ag_cwd(c) = new_litt%ag_cwd(c) + donatable_mass * donate_m2 curr_litt%ag_cwd(c) = curr_litt%ag_cwd(c) + donatable_mass * retain_m2 @@ -2768,7 +2768,8 @@ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & donatable_mass*retain_m2*dcmpy_frac end do - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dtype_ilandusechange) = & + site_mass%burn_flux_to_atm(dtype_ilandusechange) + burned_mass call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & bc_in%max_rooting_depth_index_col) @@ -2828,8 +2829,8 @@ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & burned_mass = num_dead_trees * SF_val_CWD_frac(c) * bstem * & EDPftvarcon_inst%landusechange_frac_burned(pft) - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass - + site_mass%burn_flux_to_atm(dtype_ilandusechange) = & + site_mass%burn_flux_to_atm(dtype_ilandusechange) + burned_mass else ! all other pools can end up as timber products or burn or go to litter donatable_mass = donatable_mass * (1.0_r8-EDPftvarcon_inst%landusechange_frac_exported(pft)) * & (1.0_r8-EDPftvarcon_inst%landusechange_frac_burned(pft)) @@ -2841,7 +2842,8 @@ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & woodproduct_mass = num_dead_trees * SF_val_CWD_frac(c) * bstem * & EDPftvarcon_inst%landusechange_frac_exported(pft) - site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass + site_mass%burn_flux_to_atm(dtype_ilandusechange) = & + site_mass%burn_flux_to_atm(dtype_ilandusechange) + burned_mass ! Amount of trunk mass exported off site [kg/m2] elflux_diags%exported_harvest = elflux_diags%exported_harvest + & diff --git a/main/ChecksBalancesMod.F90 b/main/ChecksBalancesMod.F90 index 325a1089e8..cb0df00c6b 100644 --- a/main/ChecksBalancesMod.F90 +++ b/main/ChecksBalancesMod.F90 @@ -331,7 +331,7 @@ subroutine CheckIntegratedMassPools(site) ibal%iflux_litter = ibal%iflux_litter + & tot_litter_input - & (site_mass%frag_out*area_inv - ediag%tot_seed_turnover) - & - (site_mass%burn_flux_to_atm*area_inv - ediag%burned_liveveg) + (sum(site_mass%burn_flux_to_atm(:))*area_inv - ediag%burned_liveveg) ediag%err_liveveg = ibal%iflux_liveveg - ibal%state_liveveg diff --git a/main/EDMainMod.F90 b/main/EDMainMod.F90 index 7a5bd21840..a2c527d417 100644 --- a/main/EDMainMod.F90 +++ b/main/EDMainMod.F90 @@ -1006,7 +1006,7 @@ subroutine TotalBalanceCheck (currentSite, call_index ) flux_out = sum(site_mass%wood_product_harvest(:)) + & sum(site_mass%wood_product_landusechange(:)) + & - site_mass%burn_flux_to_atm + & + sum(site_mass%burn_flux_to_atm(:)) + & site_mass%seed_out + & site_mass%flux_generic_out + & site_mass%frag_out + & @@ -1039,7 +1039,7 @@ subroutine TotalBalanceCheck (currentSite, call_index ) write(fates_log(),*) 'wood_product_harvest: ',site_mass%wood_product_harvest(:) write(fates_log(),*) 'wood_product_landusechange: ',site_mass%wood_product_landusechange(:) write(fates_log(),*) 'error from patch resizing: ',site_mass%patch_resize_err - write(fates_log(),*) 'burn_flux_to_atm: ',site_mass%burn_flux_to_atm + write(fates_log(),*) 'burn_flux_to_atm: ',site_mass%burn_flux_to_atm(:) write(fates_log(),*) 'seed_out: ',site_mass%seed_out write(fates_log(),*) 'flux_generic_out: ',site_mass%flux_generic_out write(fates_log(),*) 'frag_out: ',site_mass%frag_out diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 90be4df5ec..22e81073c3 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -301,7 +301,7 @@ module EDTypesMod real(r8) :: wood_product_landusechange(maxpft) ! Total mass exported as wood product from land use change [kg/site/day] - real(r8) :: burn_flux_to_atm ! Total mass burned and exported to the atmosphere [kg/site/day] + real(r8) :: burn_flux_to_atm(n_dist_types) ! Total mass burned and exported to the atmosphere [kg/site/day] real(r8) :: flux_generic_in ! Used for prescribed or artificial input fluxes ! and initialization [kg/site/day] @@ -726,7 +726,7 @@ subroutine ZeroMassBalFlux(this) this%frag_out = 0._r8 this%wood_product_harvest(:) = 0._r8 this%wood_product_landusechange(:) = 0._r8 - this%burn_flux_to_atm = 0._r8 + this%burn_flux_to_atm(:) = 0._r8 this%flux_generic_in = 0._r8 this%flux_generic_out = 0._r8 this%patch_resize_err = 0._r8 diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 0ac13af5a1..a9e633002e 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -41,6 +41,7 @@ module FatesHistoryInterfaceMod use FatesConstantsMod , only : dtype_ifall use FatesConstantsMod , only : dtype_ifire use FatesConstantsMod , only : dtype_ilog + use FatesConstantsMod , only : dtype_ilandusechange use FatesIODimensionsMod , only : fates_io_dimension_type use FatesIOVariableKindMod , only : fates_io_variable_kind_type use FatesIOVariableKindMod , only : site_int @@ -406,6 +407,7 @@ module FatesHistoryInterfaceMod integer :: ih_vis_rad_err_si integer :: ih_nir_rad_err_si integer :: ih_fire_c_to_atm_si + integer :: ih_fire_c_to_atm_landusechange_si integer :: ih_interr_liveveg_elem integer :: ih_interr_litter_elem integer :: ih_cbal_err_fates_si @@ -2496,6 +2498,7 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) hio_canopy_mortality_crownarea_si => this%hvars(ih_canopy_mortality_crownarea_si)%r81d, & hio_ustory_mortality_crownarea_si => this%hvars(ih_understory_mortality_crownarea_si)%r81d, & hio_fire_c_to_atm_si => this%hvars(ih_fire_c_to_atm_si)%r81d, & + hio_fire_c_to_atm_landusechange_si => this%hvars(ih_fire_c_to_atm_landusechange_si)%r81d, & hio_demotion_carbonflux_si => this%hvars(ih_demotion_carbonflux_si)%r81d, & hio_promotion_carbonflux_si => this%hvars(ih_promotion_carbonflux_si)%r81d, & hio_canopy_mortality_carbonflux_si => this%hvars(ih_canopy_mortality_carbonflux_si)%r81d, & @@ -2540,7 +2543,11 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) ! Total carbon lost to atmosphere from burning (kgC/site/day -> kgC/m2/s) hio_fire_c_to_atm_si(io_si) = & - sites(s)%mass_balance(element_pos(carbon12_element))%burn_flux_to_atm * & + sum(sites(s)%mass_balance(element_pos(carbon12_element))%burn_flux_to_atm(:)) * & + ha_per_m2 * days_per_sec + + hio_fire_c_to_atm_landusechange_si(io_si) = & + sites(s)%mass_balance(element_pos(carbon12_element))%burn_flux_to_atm(dtype_ilandusechange) * & ha_per_m2 * days_per_sec ! damage variables - site level - this needs to be OUT of the patch loop @@ -3395,7 +3402,7 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) ! Total element lost to atmosphere from burning (kg/site/day -> kg/m2/s) hio_burn_flux_elem(io_si,el) = & - sites(s)%mass_balance(el)%burn_flux_to_atm * ha_per_m2 * & + sum(sites(s)%mass_balance(el)%burn_flux_to_atm(:)) * ha_per_m2 * & days_per_sec end do @@ -7004,6 +7011,12 @@ subroutine define_history_vars(this, initialize_variables) upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & index = ih_fire_c_to_atm_si) + call this%set_history_var(vname='FATES_FIRE_CLOSS_LANDUSECHANGE', units='kg m-2 s-1', & + long='carbon loss to atmosphere from fire in kg carbon per m2 per second from land use change only', & + use_default='active', avgflag='A', vtype=site_r8, hlms='CLM:ALM', & + upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & + index = ih_fire_c_to_atm_landusechange_si) + call this%set_history_var(vname='FATES_CBALANCE_ERROR', & units='kg s-1', & long='total carbon error in kg carbon per second', & From 10d02e0ad24f24a7cab86b475d2e5217fd996f16 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Tue, 29 Jul 2025 10:20:09 -0700 Subject: [PATCH 04/10] fixed handoff to fire C flux BC out --- main/EDMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/EDMainMod.F90 b/main/EDMainMod.F90 index a2c527d417..3563ba0105 100644 --- a/main/EDMainMod.F90 +++ b/main/EDMainMod.F90 @@ -926,8 +926,8 @@ subroutine ed_update_site( currentSite, bc_in, bc_out, is_restarting ) ! Set boundary condition to HLM for carbon loss to atm from fires and grazing ! [kgC/ha/day]*[ha/m2]*[day/s] = [kg/m2/s] - - bc_out%fire_closs_to_atm_si = site_cmass%burn_flux_to_atm * ha_per_m2 * days_per_sec + + bc_out%fire_closs_to_atm_si = sum(site_cmass%burn_flux_to_atm(:)) * ha_per_m2 * days_per_sec bc_out%grazing_closs_to_atm_si = site_cmass%herbivory_flux_out * ha_per_m2 * days_per_sec From d0d2f1c702929b2415438b21384ed5bf84f7ab98 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Tue, 29 Jul 2025 14:43:22 -0700 Subject: [PATCH 05/10] added hist var for fire C loss from live fuels only --- main/FatesHistoryInterfaceMod.F90 | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index a9e633002e..d0c2f0c364 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -408,6 +408,7 @@ module FatesHistoryInterfaceMod integer :: ih_nir_rad_err_si integer :: ih_fire_c_to_atm_si integer :: ih_fire_c_to_atm_landusechange_si + integer :: ih_fire_livec_to_atm_si integer :: ih_interr_liveveg_elem integer :: ih_interr_litter_elem integer :: ih_cbal_err_fates_si @@ -2498,6 +2499,7 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) hio_canopy_mortality_crownarea_si => this%hvars(ih_canopy_mortality_crownarea_si)%r81d, & hio_ustory_mortality_crownarea_si => this%hvars(ih_understory_mortality_crownarea_si)%r81d, & hio_fire_c_to_atm_si => this%hvars(ih_fire_c_to_atm_si)%r81d, & + hio_fire_livec_to_atm_si => this%hvars(ih_fire_livec_to_atm_si)%r81d, & hio_fire_c_to_atm_landusechange_si => this%hvars(ih_fire_c_to_atm_landusechange_si)%r81d, & hio_demotion_carbonflux_si => this%hvars(ih_demotion_carbonflux_si)%r81d, & hio_promotion_carbonflux_si => this%hvars(ih_promotion_carbonflux_si)%r81d, & @@ -2684,6 +2686,9 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) sum(elflux_diags_c%root_litter_input(:))) * & AREA_INV * days_per_sec + hio_fire_livec_to_atm_si(io_si) = & + elflux_diags_c%burned_liveveg * ha_per_m2 * days_per_sec + ! Loop through patches to sum up diagonistics cpatch => sites(s)%oldest_patch patchloop: do while(associated(cpatch)) @@ -7011,6 +7016,12 @@ subroutine define_history_vars(this, initialize_variables) upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & index = ih_fire_c_to_atm_si) + call this%set_history_var(vname='FATES_FIRE_CLOSS_LIVEFUELS', units='kg m-2 s-1', & + long='carbon loss to atmosphere from live fuels only via fire in kg carbon per m2 per second', & + use_default='active', avgflag='A', vtype=site_r8, hlms='CLM:ALM', & + upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & + index = ih_fire_livec_to_atm_si) + call this%set_history_var(vname='FATES_FIRE_CLOSS_LANDUSECHANGE', units='kg m-2 s-1', & long='carbon loss to atmosphere from fire in kg carbon per m2 per second from land use change only', & use_default='active', avgflag='A', vtype=site_r8, hlms='CLM:ALM', & From 7ef58a8b472cf6b3a81aa3713100f0b54d54f4df Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Mon, 14 Jul 2025 15:01:08 -0700 Subject: [PATCH 06/10] added hist variable for 95th percentile height across patches --- main/FatesHistoryInterfaceMod.F90 | 29 ++++++++++++++++++++++++++++- 1 file changed, 28 insertions(+), 1 deletion(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index d0c2f0c364..0dbe4e3edf 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -307,6 +307,7 @@ module FatesHistoryInterfaceMod integer :: ih_seedling_pool_si ! carbon only integer :: ih_ba_weighted_height_si integer :: ih_ca_weighted_height_si + integer :: ih_patch_weighted_95thpctile_height_si integer :: ih_seeds_in_local_elem integer :: ih_seeds_in_extern_elem integer :: ih_seed_decay_elem @@ -2432,7 +2433,8 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) real(r8) :: leaf_herbivory ! mass of leaves eaten by herbivores [kg/yr] real(r8) :: n_perm2 ! abundance per m2 real(r8) :: patch_fracarea ! Fraction of area for this patch - + real(r8) :: crown_area_covered ! accumulator variable for patch crown area + associate( hio_npatches_si => this%hvars(ih_npatches_si)%r81d, & hio_ncohorts_si => this%hvars(ih_ncohorts_si)%r81d, & hio_ncl_si => this%hvars(ih_ncl_si)%r81d, & @@ -2443,6 +2445,7 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) hio_fates_fraction_si => this%hvars(ih_fates_fraction_si)%r81d, & hio_ba_weighted_height_si => this%hvars(ih_ba_weighted_height_si)%r81d, & hio_ca_weighted_height_si => this%hvars(ih_ca_weighted_height_si)%r81d, & + hio_patch_weighted_95thpctile_height_si => this%hvars(ih_patch_weighted_95thpctile_height_si)%r81d, & hio_canopy_spread_si => this%hvars(ih_canopy_spread_si)%r81d, & hio_nesterov_fire_danger_si => this%hvars(ih_nesterov_fire_danger_si)%r81d, & hio_rx_burn_window_si => this%hvars(ih_rx_burn_window_si)%r81d, & @@ -3018,6 +3021,24 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) ccohort => ccohort%taller enddo cohortloop ! cohort loop + ! mean 95th percentile height. loop through cohorts on patch again, this time from tallest to shortest + crown_area_covered = 0._r8 + ccohort => cpatch%tallest + do while(associated(ccohort)) + if (ccohort%canopy_layer .eq. 1) then ! ignore anything not in the canopy + if (((ccohort%c_area + crown_area_covered)/cpatch%area) .ge. 0.05_r8 ) then + hio_patch_weighted_95thpctile_height_si(io_si) = & + hio_patch_weighted_95thpctile_height_si(io_si) + ccohort%height * cpatch%area/AREA + ccohort => null() ! exit the cohort loop + else + crown_area_covered = crown_area_covered + ccohort%c_area + ccohort => ccohort%shorter + endif + else + ccohort => ccohort%shorter + endif + enddo ! cohort loop + cpatch => cpatch%younger end do patchloop !patch loop @@ -6414,6 +6435,12 @@ subroutine define_history_vars(this, initialize_variables) upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & index=ih_ca_weighted_height_si) + call this%set_history_var(vname='FATES_MEAN_95PCTILE_HEIGHT', units='m', & + long='The patch-area-weighted mean of 95th percentile height of canopy plants within a given patch', use_default='active', & + avgflag='A', vtype=site_r8, hlms='CLM:ALM', & + upfreq=group_dyna_simple, ivar=ivar, initialize=initialize_variables, & + index=ih_patch_weighted_95thpctile_height_si) + call this%set_history_var(vname='FATES_COLD_STATUS', units='', & long='site-level cold status, 0=not cold-dec, 1=too cold for leaves, 2=not too cold', & use_default='active', avgflag='A', vtype=site_r8, hlms='CLM:ALM', & From f96296de1089afef5c634167798e271b1b252d39 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Wed, 10 Sep 2025 15:22:30 -0700 Subject: [PATCH 07/10] added and subtracted some history vars --- main/FatesHistoryInterfaceMod.F90 | 67 +++++++++---------------------- 1 file changed, 20 insertions(+), 47 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 0dbe4e3edf..10d043466e 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -679,14 +679,9 @@ module FatesHistoryInterfaceMod integer :: ih_biomass_si_age integer :: ih_c_stomata_si_age integer :: ih_c_lblayer_si_age - integer :: ih_agesince_anthrodist_si integer :: ih_agesince_anthrodist_si_age - integer :: ih_secondarylands_area_si_age - integer :: ih_primarylands_area_si_age - integer :: ih_area_burnt_si_age - integer :: ih_primarylands_fracarea_si - integer :: ih_secondarylands_fracarea_si integer :: ih_secondarylands_fracarea_si_age + integer :: ih_secondary_agb_si_agesinceanthro integer :: ih_primarylands_fracarea_si_age integer :: ih_fracarea_burnt_si_age integer :: ih_rx_fracarea_burnt_si_age @@ -697,6 +692,7 @@ module FatesHistoryInterfaceMod integer :: ih_rx_intensity_si_age integer :: ih_nonrx_intensity_si_age + ! indices to (site x height) variables integer :: ih_canopy_height_dist_si_height integer :: ih_leaf_height_dist_si_height @@ -3333,9 +3329,6 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) hio_yesterdaycanopylevel_understory_si_scls => this%hvars(ih_yesterdaycanopylevel_understory_si_scls)%r82d, & hio_fracarea_si => this%hvars(ih_fracarea_si)%r81d, & hio_canopy_fracarea_si => this%hvars(ih_canopy_fracarea_si)%r81d, & - hio_agesince_anthrodist_si => this%hvars(ih_agesince_anthrodist_si)%r81d, & - hio_primarylands_fracarea_si => this%hvars(ih_primarylands_fracarea_si)%r81d, & - hio_secondarylands_fracarea_si => this%hvars(ih_secondarylands_fracarea_si)%r81d, & hio_fracarea_si_landuse => this%hvars(ih_fracarea_si_landuse)%r82d, & hio_npp_si_landuse => this%hvars(ih_npp_si_landuse)%r82d, & hio_biomass_si_landuse => this%hvars(ih_biomass_si_landuse)%r82d, & @@ -3370,6 +3363,7 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) hio_cstarvmortality_continuous_carbonflux_si_pft => this%hvars(ih_cstarvmortality_continuous_carbonflux_si_pft)%r82d, & hio_transition_matrix_si_lulu => this%hvars(ih_transition_matrix_si_lulu)%r82d, & hio_scorch_height_si_pft => this%hvars(ih_scorch_height_si_pft)%r82d, & + hio_secondary_agb_si_agesinceanthro => this%hvars(ih_secondary_agb_si_agesinceanthro)%r82d, & hio_sapwood_area_scpf => this%hvars(ih_sapwood_area_scpf)%r82d) model_day_int = nint(hlm_model_day) @@ -3476,25 +3470,6 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) cpatch%frac_burnt * cpatch%area * AREA_INV / sec_per_day end if - ! some diagnostics on secondary forest area and its age distribution - if ( cpatch%land_use_label .eq. secondaryland ) then - - hio_agesince_anthrodist_si(io_si) = & - hio_agesince_anthrodist_si(io_si) & - + cpatch%area * AREA_INV - - hio_secondarylands_fracarea_si(io_si) = & - hio_secondarylands_fracarea_si(io_si) & - + cpatch%area * AREA_INV - - else if ( cpatch%land_use_label .eq. primaryland ) then - - hio_primarylands_fracarea_si(io_si) = & - hio_primarylands_fracarea_si(io_si) & - + cpatch%area * AREA_INV - - endif - do ft = 1,numpft hio_scorch_height_si_pft(io_si,ft) = hio_scorch_height_si_pft(io_si,ft) + & cpatch%Scorch_ht(ft) * cpatch%area * AREA_INV @@ -3923,6 +3898,16 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) hio_biomass_si_scls(io_si,scls) = hio_biomass_si_scls(io_si,scls) + & total_m * ccohort%n * AREA_INV + ! cdkcdk agb by patch age since anthropogenic disturbance + if ( cpatch%land_use_label .eq. secondaryland ) then + + iscag_anthrodist = get_age_class_index(cpatch%age_since_anthro_disturbance) + + hio_secondary_agb_si_agesinceanthro(io_si,iscag_anthrodist) = & + hio_secondary_agb_si_agesinceanthro(io_si,iscag_anthrodist) & + + total_m * ccohort%n * prt_params%allom_agb_frac(ccohort%pft) * AREA_INV + endif + ! update size-class quantities hio_nplant_si_scls(io_si,scls) = hio_nplant_si_scls(io_si,scls) + ccohort%n / m2_per_ha @@ -4911,6 +4896,8 @@ subroutine update_history_dyn_subsite_ageclass(this,nc,nsites,sites) hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) = & hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) & + patch_area_div_site_area + + else if ( cpatch%land_use_label .eq. primaryland) then hio_primarylands_fracarea_si_age(io_si,cpatch%age_class) = & hio_primarylands_fracarea_si_age(io_si,cpatch%age_class) & @@ -7446,20 +7433,6 @@ subroutine define_history_vars(this, initialize_variables) hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & index=ih_agesince_anthrodist_si_age) - call this%set_history_var(vname='FATES_SECONDARY_AREA_ANTHRO', & - units='m2 m-2', & - long='secondary forest patch area since anthropgenic disturbance', & - use_default='inactive', avgflag='A', vtype=site_r8, & - hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & - index=ih_agesince_anthrodist_si) - - call this%set_history_var(vname='FATES_SECONDARY_AREA', & - units='m2 m-2', & - long='secondary forest patch area since any kind of disturbance', & - use_default='inactive', avgflag='A', vtype=site_r8, & - hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & - index=ih_secondarylands_fracarea_si) - call this%set_history_var(vname='FATES_SECONDARY_AREA_AP', & units='m2 m-2', & long='secondary forest patch area age distribution since any kind of disturbance', & @@ -7467,12 +7440,12 @@ subroutine define_history_vars(this, initialize_variables) hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & index=ih_secondarylands_fracarea_si_age) - call this%set_history_var(vname='FATES_PRIMARY_AREA', & - units='m2 m-2', & - long='primary forest patch area since any kind of disturbance', & - use_default='inactive', avgflag='A', vtype=site_r8, & + call this%set_history_var(vname='FATES_SECONDARY_AGB_ANTHROAGE_AP', & + units='kg m-2', & + long='secondary forest patch agb as resolved by age since anthropogenic disturbance', & + use_default='inactive', avgflag='A', vtype=site_age_r8, & hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & - index=ih_primarylands_fracarea_si) + index=ih_secondary_agb_si_agesinceanthro) call this%set_history_var(vname='FATES_PRIMARY_AREA_AP', & units='m2 m-2', & From bd40d256bb643ec41c9f369dc6df3db5c22f8ced Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Fri, 26 Sep 2025 15:24:08 -0700 Subject: [PATCH 08/10] added units to grazing term in coment --- main/EDTypesMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 22e81073c3..26d39bcf45 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -311,7 +311,7 @@ module EDTypesMod ! due to re-sizing patches when area math starts to lose ! precision - real(r8) :: herbivory_flux_out ! loss of element due to grazing (and/or browsing) by herbivores + real(r8) :: herbivory_flux_out ! loss of element due to grazing (and/or browsing) by herbivores [kg/site/day] contains From 641a2cf10cbfa7c1c4a44ddeb205afdd0fe758cb Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Fri, 10 Oct 2025 11:21:34 -0700 Subject: [PATCH 09/10] cleaning up some comments --- main/FatesHistoryInterfaceMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 10d043466e..2f2627fb13 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -2964,7 +2964,7 @@ subroutine update_history_dyn_sitelevel(this,nc,nsites,sites) hio_npp_stor_si(io_si) = hio_npp_stor_si(io_si) + & store_m_net_alloc * n_perm2 / days_per_year / sec_per_day - leaf_herbivory = ccohort%prt%GetHerbivory(leaf_organ, carbon12_element) * days_per_year !cdkcdk + leaf_herbivory = ccohort%prt%GetHerbivory(leaf_organ, carbon12_element) * days_per_year hio_grazing_si(io_si) = hio_grazing_si(io_si) + leaf_herbivory * n_perm2 / days_per_year / sec_per_day ! Woody State Variables (basal area growth increment) @@ -3898,7 +3898,7 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) hio_biomass_si_scls(io_si,scls) = hio_biomass_si_scls(io_si,scls) + & total_m * ccohort%n * AREA_INV - ! cdkcdk agb by patch age since anthropogenic disturbance + ! agb by patch age since anthropogenic disturbance if ( cpatch%land_use_label .eq. secondaryland ) then iscag_anthrodist = get_age_class_index(cpatch%age_since_anthro_disturbance) From deed8889f6437f4b6effe5714c9532ef67f89af9 Mon Sep 17 00:00:00 2001 From: Charles D Koven Date: Fri, 10 Oct 2025 14:09:00 -0700 Subject: [PATCH 10/10] updates per code review comments --- main/FatesHistoryInterfaceMod.F90 | 27 ++++++++++++++------------- 1 file changed, 14 insertions(+), 13 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 65fa6c841a..65752ddca1 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -5275,6 +5275,7 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in,dt_tstep) type(fates_patch_type),pointer :: cpatch type(fates_cohort_type),pointer :: ccohort real(r8) :: dt_tstep_inv ! Time step in frequency units (/s) + real(r8) :: vegarea_per_patcharea ! temporary weighting variable (unitless) associate( hio_tveg_si_landuse => this%hvars(ih_tveg_si_landuse)%r82d,& hio_gpp_si_landuse => this%hvars(ih_gpp_si_landuse)%r82d, & @@ -5312,25 +5313,27 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in,dt_tstep) bc_in(s)%t_veg_pa(cpatch%patchno) * cpatch%total_canopy_area/canopy_area_bylanduse(cpatch%land_use_label) ! for the rest of these, first weight by the vegetated area of each patch over the total patch area for each land use type + vegarea_per_patcharea = cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + hio_tsa_si_landuse(io_si,cpatch%land_use_label) = hio_tsa_si_landuse(io_si,cpatch%land_use_label) + & - bc_in(s)%t2m_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + bc_in(s)%t2m_pa(cpatch%patchno) * vegarea_per_patcharea hio_sw_abs_si_landuse(io_si,cpatch%land_use_label) = hio_sw_abs_si_landuse(io_si,cpatch%land_use_label) + & - bc_in(s)%swabs_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + bc_in(s)%swabs_pa(cpatch%patchno) * vegarea_per_patcharea hio_lw_net_si_landuse(io_si,cpatch%land_use_label) = hio_lw_net_si_landuse(io_si,cpatch%land_use_label) + & - bc_in(s)%netlw_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + bc_in(s)%netlw_pa(cpatch%patchno) * vegarea_per_patcharea hio_shflux_si_landuse(io_si,cpatch%land_use_label) = hio_shflux_si_landuse(io_si,cpatch%land_use_label) + & - bc_in(s)%shflux_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + bc_in(s)%shflux_pa(cpatch%patchno) * vegarea_per_patcharea hio_lhflux_si_landuse(io_si,cpatch%land_use_label) = hio_lhflux_si_landuse(io_si,cpatch%land_use_label) + & - bc_in(s)%lhflux_pa(cpatch%patchno) * cpatch%total_canopy_area/landuse_statevector(cpatch%land_use_label) + bc_in(s)%lhflux_pa(cpatch%patchno) * vegarea_per_patcharea endif cpatch => cpatch%younger end do - ! for all the land-use indexed variables, excpet for TVEG, also add in the component for the unvegetated area of each land use + ! for all the land-use indexed variables, except for TVEG, also add in the component for the unvegetated area of each land use foundbaregroundpatch = .false. cpatch => sites(s)%oldest_patch do while(associated(cpatch)) @@ -5382,16 +5385,14 @@ subroutine update_history_hifrq_landuse(this,nc,nsites,sites,bc_in,dt_tstep) end do - + ! for GPP by land use, we need to loop over both patches and cohorts cpatch => sites(s)%oldest_patch do while(associated(cpatch)) ccohort => cpatch%shortest do while(associated(ccohort)) - if ( .not. ccohort%isnew ) then - if (cpatch%land_use_label .gt. nocomp_bareground_land) then - hio_gpp_si_landuse(io_si,cpatch%land_use_label) = hio_gpp_si_landuse(io_si,cpatch%land_use_label) & - + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv - end if + if ( (.not. ccohort%isnew) .and. (cpatch%land_use_label .gt. nocomp_bareground_land) ) then + hio_gpp_si_landuse(io_si,cpatch%land_use_label) = hio_gpp_si_landuse(io_si,cpatch%land_use_label) & + + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv end if ccohort => ccohort%taller end do @@ -8942,7 +8943,7 @@ subroutine define_history_vars(this, initialize_variables) ivar=ivar, initialize=initialize_variables, index = ih_lhflux_si_landuse ) call this%set_history_var(vname='FATES_GPP_LU', units='kg m-2 s-1', & - long='gross primary productivity by age bin in kg carbon per m2 per second', & + long='gross primary productivity by land use type in kg carbon per m2 per second', & use_default='inactive', avgflag='A', vtype=site_landuse_r8, & hlms='CLM:ALM', upfreq=group_hifr_simple, ivar=ivar, initialize=initialize_variables, & index = ih_gpp_si_landuse)