From ffc5e40311b7446030da6406c8dd1aeb502f0889 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:06:14 -0600 Subject: [PATCH 01/21] Rename internal area-related history variables for clarity. FATES_PATCHAREA_LU and FATES_PATCHAREA_AP are not in units of m2, but rather m2/m2. It is thus misleading for them to be named "AREA". In the interest of not messing up user workflows, I'm not changing those---at least for now. However, this commit renames related internal variables in the history code for clarity. Specifically, "area" is replaced with "fracarea" (fractional area, relative to site area) in the following: - ih_area_si_landuse - ih_area_si_age - hio_area_si_landuse - hio_fracarea_si_age --- main/FatesHistoryInterfaceMod.F90 | 26 +++++++++++++------------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index c5b48735f6..7207e256e0 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -358,7 +358,7 @@ module FatesHistoryInterfaceMod integer :: ih_primaryland_fusion_error_si - integer :: ih_area_si_landuse + integer :: ih_fracarea_si_landuse integer :: ih_disturbance_rate_si_lulu integer :: ih_transition_matrix_si_lulu @@ -648,7 +648,7 @@ module FatesHistoryInterfaceMod integer :: ih_seeds_in_gc_si_pft ! indices to (site x patch-age) variables - integer :: ih_area_si_age + integer :: ih_fracarea_si_age integer :: ih_lai_si_age integer :: ih_canopy_area_si_age integer :: ih_gpp_si_age @@ -3258,7 +3258,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_scorch_height_si_agepft => this%hvars(ih_scorch_height_si_agepft)%r82d, & hio_yesterdaycanopylevel_canopy_si_scls => this%hvars(ih_yesterdaycanopylevel_canopy_si_scls)%r82d, & hio_yesterdaycanopylevel_understory_si_scls => this%hvars(ih_yesterdaycanopylevel_understory_si_scls)%r82d, & - hio_area_si_age => this%hvars(ih_area_si_age)%r82d, & + hio_fracarea_si_age => this%hvars(ih_fracarea_si_age)%r82d, & hio_lai_si_age => this%hvars(ih_lai_si_age)%r82d, & hio_canopy_area_si_age => this%hvars(ih_canopy_area_si_age)%r82d, & hio_ncl_si_age => this%hvars(ih_ncl_si_age)%r82d, & @@ -3267,7 +3267,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_biomass_si_age => this%hvars(ih_biomass_si_age)%r82d, & hio_agesince_anthrodist_si_age => this%hvars(ih_agesince_anthrodist_si_age)%r82d, & hio_secondarylands_area_si_age => this%hvars(ih_secondarylands_area_si_age)%r82d, & - hio_area_si_landuse => this%hvars(ih_area_si_landuse)%r82d, & + hio_fracarea_si_landuse => this%hvars(ih_fracarea_si_landuse)%r82d, & hio_area_burnt_si_age => this%hvars(ih_area_burnt_si_age)%r82d, & ! hio_fire_rate_of_spread_front_si_age => this%hvars(ih_fire_rate_of_spread_front_si_age)%r82d, & hio_fire_intensity_si_age => this%hvars(ih_fire_intensity_si_age)%r82d, & @@ -3395,13 +3395,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) cpatch%age_class = get_age_class_index(cpatch%age) ! Increment the fractional area in each age class bin - hio_area_si_age(io_si,cpatch%age_class) = hio_area_si_age(io_si,cpatch%age_class) & + hio_fracarea_si_age(io_si,cpatch%age_class) = hio_fracarea_si_age(io_si,cpatch%age_class) & + cpatch%area * AREA_INV ! ignore land use info on nocomp bareground (where landuse label = 0) if (cpatch%land_use_label .gt. nocomp_bareground_land) then - hio_area_si_landuse(io_si, cpatch%land_use_label) = & - hio_area_si_landuse(io_si, cpatch%land_use_label) & + hio_fracarea_si_landuse(io_si, cpatch%land_use_label) = & + hio_fracarea_si_landuse(io_si, cpatch%land_use_label) & + cpatch%area * AREA_INV end if @@ -4305,13 +4305,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! divide so-far-just-summed but to-be-averaged patch-age-class ! variables by patch-age-class area to get mean values do ipa2 = 1, nlevage - if (hio_area_si_age(io_si, ipa2) .gt. nearzero) then - hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / (hio_area_si_age(io_si, ipa2)*AREA) - hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / (hio_area_si_age(io_si, ipa2)*AREA) + if (hio_fracarea_si_age(io_si, ipa2) .gt. nearzero) then + hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / (hio_fracarea_si_age(io_si, ipa2)*AREA) + hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / (hio_fracarea_si_age(io_si, ipa2)*AREA) do ft = 1, numpft iagepft = ipa2 + (ft-1) * nlevage hio_scorch_height_si_agepft(io_si, iagepft) = & - hio_scorch_height_si_agepft(io_si, iagepft) / (hio_area_si_age(io_si, ipa2)*AREA) + hio_scorch_height_si_agepft(io_si, iagepft) / (hio_fracarea_si_age(io_si, ipa2)*AREA) enddo else hio_lai_si_age(io_si, ipa2) = 0._r8 @@ -6794,7 +6794,7 @@ subroutine define_history_vars(this, initialize_variables) call this%set_history_var(vname='FATES_PATCHAREA_LU', units='m2 m-2', & long='patch area by land use type', use_default='active', & avgflag='A', vtype=site_landuse_r8, hlms='CLM:ALM', upfreq=group_dyna_complx, & - ivar=ivar, initialize=initialize_variables, index=ih_area_si_landuse) + ivar=ivar, initialize=initialize_variables, index=ih_fracarea_si_landuse) call this%set_history_var(vname='FATES_TRANSITION_MATRIX_LULU', units='m2 m-2 yr-1', & long='land use transition matrix', use_default='active', & @@ -6977,7 +6977,7 @@ subroutine define_history_vars(this, initialize_variables) call this%set_history_var(vname='FATES_PATCHAREA_AP', units='m2 m-2', & long='patch area by age bin per m2 land area', use_default='active', & avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, & - initialize=initialize_variables, index=ih_area_si_age) + initialize=initialize_variables, index=ih_fracarea_si_age) call this%set_history_var(vname='FATES_LAI_AP', units='m2 m-2', & long='total leaf area index by age bin per m2 land area', & From fcf2b54cb5c688b0a3cef8d215b24b4e65ad1179 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:12:40 -0600 Subject: [PATCH 02/21] FatesHistoryInterfaceMod: Rename frac_area to patch_fracarea. 1) Makes it clearer that this is a patch-level variable. 2) Uses "fracarea" instead of "area" for clarity and consistency with previous commit. --- main/FatesHistoryInterfaceMod.F90 | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 7207e256e0..ff4f5ed7dd 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -2389,7 +2389,7 @@ subroutine update_history_dyn1(this,nc,nsites,sites,bc_in) real(r8) :: struct_m_net_alloc ! mass allocated to structure [kg/yr] real(r8) :: repro_m_net_alloc ! mass allocated to reproduction [kg/yr] real(r8) :: n_perm2 ! abundance per m2 - real(r8) :: area_frac ! Fraction of area for this patch + real(r8) :: patch_fracarea ! Fraction of area for this patch associate( hio_npatches_si => this%hvars(ih_npatches_si)%r81d, & hio_npatches_sec_si => this%hvars(ih_npatches_sec_si)%r81d, & @@ -2685,7 +2685,7 @@ subroutine update_history_dyn1(this,nc,nsites,sites,bc_in) litt => cpatch%litter(element_pos(carbon12_element)) - area_frac = cpatch%area * AREA_INV + patch_fracarea = cpatch%area * AREA_INV ! Sum up all output fluxes (fragmentation) kgC/m2/day -> kgC/m2/s hio_litter_out_si(io_si) = hio_litter_out_si(io_si) + & @@ -2695,29 +2695,29 @@ subroutine update_history_dyn1(this,nc,nsites,sites,bc_in) sum(litt%bg_cwd_frag(:,:)) + & sum(litt%seed_decay(:)) + & sum(litt%seed_germ_decay(:))) * & - area_frac * days_per_sec + patch_fracarea * days_per_sec ! Sum up total seed bank (germinated and ungerminated) hio_seed_bank_si(io_si) = hio_seed_bank_si(io_si) + & (sum(litt%seed(:))+sum(litt%seed_germ(:))) * & - area_frac + patch_fracarea ! Sum up total seed bank (just ungerminated) hio_ungerm_seed_bank_si(io_si) = hio_ungerm_seed_bank_si(io_si) + & - sum(litt%seed(:)) * area_frac + sum(litt%seed(:)) * patch_fracarea ! Sum up total seedling pool hio_seedling_pool_si(io_si) = hio_seedling_pool_si(io_si) + & - sum(litt%seed_germ(:)) * area_frac + sum(litt%seed_germ(:)) * patch_fracarea ! Sum up the input flux into the seed bank (local and external) hio_seeds_in_si(io_si) = hio_seeds_in_si(io_si) + & (sum(litt%seed_in_local(:)) + sum(litt%seed_in_extern(:))) * & - area_frac * days_per_sec + patch_fracarea * days_per_sec hio_seeds_in_local_si(io_si) = hio_seeds_in_local_si(io_si) + & sum(litt%seed_in_local(:)) * & - area_frac * days_per_sec + patch_fracarea * days_per_sec ! loop through cohorts on patch ccohort => cpatch%shortest @@ -3058,7 +3058,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: n_perm2 ! abundance per m2 integer :: ageclass_since_anthrodist ! what is the equivalent age class for ! time-since-anthropogenic-disturbance of secondary forest - real(r8) :: area_frac ! Fraction of area for this patch + real(r8) :: patch_fracarea ! Fraction of area for this patch real(r8) :: frac_canopy_in_bin ! fraction of a leaf's canopy that is within a given height bin real(r8) :: binbottom,bintop ! edges of height bins integer :: height_bin_max, height_bin_min ! which height bin a given cohort's canopy is in @@ -4633,7 +4633,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) litt => cpatch%litter(el) - area_frac = cpatch%area * AREA_INV + patch_fracarea = cpatch%area * AREA_INV ! Sum up all output fluxes (fragmentation) hio_litter_out_elem(io_si,el) = hio_litter_out_elem(io_si,el) + & From 4adb6655f573f4f8c5da142187795d5b150e4f28 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:23:46 -0600 Subject: [PATCH 03/21] Replace use of hio_fracarea_si_age with Site class area_by_age. --- main/FatesHistoryInterfaceMod.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index ff4f5ed7dd..7494084235 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -4305,13 +4305,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! divide so-far-just-summed but to-be-averaged patch-age-class ! variables by patch-age-class area to get mean values do ipa2 = 1, nlevage - if (hio_fracarea_si_age(io_si, ipa2) .gt. nearzero) then - hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / (hio_fracarea_si_age(io_si, ipa2)*AREA) - hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / (hio_fracarea_si_age(io_si, ipa2)*AREA) + if (sites(s)%area_by_age(ipa2) .gt. nearzero) then + hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / sites(s)%area_by_age(ipa2) + hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / sites(s)%area_by_age(ipa2) do ft = 1, numpft iagepft = ipa2 + (ft-1) * nlevage hio_scorch_height_si_agepft(io_si, iagepft) = & - hio_scorch_height_si_agepft(io_si, iagepft) / (hio_fracarea_si_age(io_si, ipa2)*AREA) + hio_scorch_height_si_agepft(io_si, iagepft) / sites(s)%area_by_age(ipa2) enddo else hio_lai_si_age(io_si, ipa2) = 0._r8 From d7a8d60f9e27e96f93358ca56e746ac49cab36bb Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:37:12 -0600 Subject: [PATCH 04/21] Refactor ageclass-area weighting for FATES_LAI_AP, FATES_NCL_AP. --- main/FatesHistoryInterfaceMod.F90 | 27 +++++++++++++++------------ 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 7494084235..d4cbf74144 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3083,6 +3083,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: storep_understory_scpf(numpft*nlevsclass) real(r8) :: storec_canopy_scpf(numpft*nlevsclass) real(r8) :: storec_understory_scpf(numpft*nlevsclass) + real(r8) :: weight integer :: i_dist, j_dist @@ -3406,12 +3407,20 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) end if ! Increment some patch-age-resolved diagnostics - hio_lai_si_age(io_si,cpatch%age_class) = hio_lai_si_age(io_si,cpatch%age_class) & - + sum(cpatch%tlai_profile(:,:,:) * cpatch%canopy_area_profile(:,:,:) ) * cpatch%total_canopy_area - - hio_ncl_si_age(io_si,cpatch%age_class) = hio_ncl_si_age(io_si,cpatch%age_class) & - + cpatch%ncl_p * cpatch%area - + if (sites(s)%area_by_age(cpatch%age_class) .gt. nearzero) then + weight = cpatch%total_canopy_area / sites(s)%area_by_age(cpatch%age_class) + hio_lai_si_age(io_si,cpatch%age_class) = hio_lai_si_age(io_si,cpatch%age_class) & + + sum(cpatch%tlai_profile(:,:,:) * cpatch%canopy_area_profile(:,:,:) ) & + * weight + + weight = cpatch%area / sites(s)%area_by_age(cpatch%age_class) + hio_ncl_si_age(io_si,cpatch%age_class) = hio_ncl_si_age(io_si,cpatch%age_class) & + + cpatch%ncl_p * weight + else + hio_lai_si_age(io_si,cpatch%age_class) = 0._r8 + hio_ncl_si_age(io_si,cpatch%age_class) = 0._r8 + end if + hio_npatches_si_age(io_si,cpatch%age_class) = hio_npatches_si_age(io_si,cpatch%age_class) + 1._r8 @@ -4306,17 +4315,11 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! variables by patch-age-class area to get mean values do ipa2 = 1, nlevage if (sites(s)%area_by_age(ipa2) .gt. nearzero) then - hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / sites(s)%area_by_age(ipa2) - hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / sites(s)%area_by_age(ipa2) do ft = 1, numpft iagepft = ipa2 + (ft-1) * nlevage hio_scorch_height_si_agepft(io_si, iagepft) = & hio_scorch_height_si_agepft(io_si, iagepft) / sites(s)%area_by_age(ipa2) enddo - else - hio_lai_si_age(io_si, ipa2) = 0._r8 - hio_ncl_si_age(io_si, ipa2) = 0._r8 - endif end do From 6950195bb86ae36b9f187795156655a1f28e7e05 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:42:49 -0600 Subject: [PATCH 05/21] Refactor ageclass-area weighting for FATES_SCORCH_HEIGHT_APPF. --- main/FatesHistoryInterfaceMod.F90 | 27 ++++++++------------------- 1 file changed, 8 insertions(+), 19 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index d4cbf74144..7d85382da5 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3451,9 +3451,14 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! for scorch height, weight the value by patch area within any ! given age class - in the event that there is more than one ! patch per age class. - iagepft = cpatch%age_class + (ft-1) * nlevage - hio_scorch_height_si_agepft(io_si,iagepft) = hio_scorch_height_si_agepft(io_si,iagepft) + & - cpatch%Scorch_ht(ft) * cpatch%area + if (sites(s)%area_by_age(cpatch%age_class) .gt. nearzero) then + iagepft = cpatch%age_class + (ft-1) * nlevage + weight = cpatch%area / sites(s)%area_by_age(cpatch%age_class) + hio_scorch_height_si_agepft(io_si,iagepft) = hio_scorch_height_si_agepft(io_si,iagepft) + & + cpatch%Scorch_ht(ft) * weight + else + hio_scorch_height_si_agepft(io_si,iagepft) = 0._r8 + end if ! and also pft-labeled patch areas in the event that we are in nocomp mode if ( hlm_use_nocomp .eq. itrue .and. cpatch%nocomp_pft_label .eq. ft) then @@ -4308,22 +4313,6 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) cpatch => cpatch%younger end do patchloop !patch loop - - - - ! divide so-far-just-summed but to-be-averaged patch-age-class - ! variables by patch-age-class area to get mean values - do ipa2 = 1, nlevage - if (sites(s)%area_by_age(ipa2) .gt. nearzero) then - do ft = 1, numpft - iagepft = ipa2 + (ft-1) * nlevage - hio_scorch_height_si_agepft(io_si, iagepft) = & - hio_scorch_height_si_agepft(io_si, iagepft) / sites(s)%area_by_age(ipa2) - enddo - end do - - - ! pass the cohort termination mortality as a flux to the history, and then reset the termination mortality buffer ! note there are various ways of reporting the total mortality, so pass to these as well do i_pft = 1, numpft From 76a70624420c7adb49600f78393afdbed56f4374 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 12:49:17 -0600 Subject: [PATCH 06/21] Use age_class_area variable in patch loop. --- main/FatesHistoryInterfaceMod.F90 | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 7d85382da5..ea8e4538e8 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3084,6 +3084,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: storec_canopy_scpf(numpft*nlevsclass) real(r8) :: storec_understory_scpf(numpft*nlevsclass) real(r8) :: weight + real(r8) :: age_class_area ! [m2] integer :: i_dist, j_dist @@ -3394,6 +3395,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) cpatch%age_class = get_age_class_index(cpatch%age) + age_class_area = sites(s)%area_by_age(cpatch%age_class) ! Increment the fractional area in each age class bin hio_fracarea_si_age(io_si,cpatch%age_class) = hio_fracarea_si_age(io_si,cpatch%age_class) & @@ -3407,13 +3409,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) end if ! Increment some patch-age-resolved diagnostics - if (sites(s)%area_by_age(cpatch%age_class) .gt. nearzero) then - weight = cpatch%total_canopy_area / sites(s)%area_by_age(cpatch%age_class) + if (age_class_area .gt. nearzero) then + weight = cpatch%total_canopy_area / age_class_area hio_lai_si_age(io_si,cpatch%age_class) = hio_lai_si_age(io_si,cpatch%age_class) & + sum(cpatch%tlai_profile(:,:,:) * cpatch%canopy_area_profile(:,:,:) ) & * weight - weight = cpatch%area / sites(s)%area_by_age(cpatch%age_class) + weight = cpatch%area / age_class_area hio_ncl_si_age(io_si,cpatch%age_class) = hio_ncl_si_age(io_si,cpatch%age_class) & + cpatch%ncl_p * weight else @@ -3451,9 +3453,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! for scorch height, weight the value by patch area within any ! given age class - in the event that there is more than one ! patch per age class. - if (sites(s)%area_by_age(cpatch%age_class) .gt. nearzero) then + if (age_class_area .gt. nearzero) then iagepft = cpatch%age_class + (ft-1) * nlevage - weight = cpatch%area / sites(s)%area_by_age(cpatch%age_class) + weight = cpatch%area / age_class_area hio_scorch_height_si_agepft(io_si,iagepft) = hio_scorch_height_si_agepft(io_si,iagepft) + & cpatch%Scorch_ht(ft) * weight else From 57723c71f8180ef445dc1d0fd8ddac81edd42361 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 13:05:16 -0600 Subject: [PATCH 07/21] Refactor (but don't fix) bad age-class weighting of several outputs. - FATES_ZSTAR_AP - FATES_SECONDAREA_ANTHRODIST_AP - FATES_SECONDAREA_DIST_AP - FATES_BURNFRAC_AP - FATES_FIRE_INTENSITY_BURNFRAC_AP - FATES_VEGC_AP - FATES_FUEL_AMOUNT_AP - FATES_CANOPYAREA_AP --- main/FatesHistoryInterfaceMod.F90 | 28 +++++++++++++++++++--------- 1 file changed, 19 insertions(+), 9 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index ea8e4538e8..9bb3172cc1 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3059,6 +3059,8 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) integer :: ageclass_since_anthrodist ! what is the equivalent age class for ! time-since-anthropogenic-disturbance of secondary forest real(r8) :: patch_fracarea ! Fraction of area for this patch + real(r8) :: ageclass_fracarea ! Fraction of area for this age class + real(r8) :: ageclass_canopy_fracarea ! Fraction of canopy area for this age class real(r8) :: frac_canopy_in_bin ! fraction of a leaf's canopy that is within a given height bin real(r8) :: binbottom,bintop ! edges of height bins integer :: height_bin_max, height_bin_min ! which height bin a given cohort's canopy is in @@ -3084,6 +3086,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: storec_canopy_scpf(numpft*nlevsclass) real(r8) :: storec_understory_scpf(numpft*nlevsclass) real(r8) :: weight + real(r8) :: denominator real(r8) :: age_class_area ! [m2] integer :: i_dist, j_dist @@ -3428,8 +3431,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) if ( ED_val_comp_excln .lt. 0._r8 ) then ! only valid when "strict ppa" enabled + weight = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area hio_zstar_si_age(io_si,cpatch%age_class) = hio_zstar_si_age(io_si,cpatch%age_class) & - + cpatch%zstar * cpatch%area * AREA_INV + + cpatch%zstar * weight endif ! some diagnostics on secondary forest area and its age distribution @@ -3437,13 +3441,15 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ageclass_since_anthrodist = get_age_class_index(cpatch%age_since_anthro_disturbance) + ageclass_fracarea = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / sites(s)%age_class_area(ageclass_since_anthrodist) hio_agesince_anthrodist_si_age(io_si,ageclass_since_anthrodist) = & hio_agesince_anthrodist_si_age(io_si,ageclass_since_anthrodist) & - + cpatch%area * AREA_INV + + ageclass_fracarea + ageclass_fracarea = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area hio_secondarylands_area_si_age(io_si,cpatch%age_class) = & hio_secondarylands_area_si_age(io_si,cpatch%age_class) & - + cpatch%area * AREA_INV + + ageclass_fracarea endif @@ -3477,20 +3483,23 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) end do + ! Weight for the following per-age-class fire variables + weight = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area + ! fractional area burnt [frac/day] -> [frac/sec] hio_area_burnt_si_age(io_si,cpatch%age_class) = hio_area_burnt_si_age(io_si,cpatch%age_class) + & - cpatch%frac_burnt * cpatch%area * AREA_INV / sec_per_day + cpatch%frac_burnt * weight / sec_per_day ! hio_fire_rate_of_spread_front_si_age(io_si, cpatch%age_class) = hio_fire_rate_of_spread_si_age(io_si, cpatch%age_class) + & - ! cpatch%ros_front * cpatch*frac_burnt * cpatch%area * AREA_INV + ! cpatch%ros_front * cpatch*frac_burnt * weight ! Fire intensity weighted by burned fraction [kJ/m/s] -> [J/m/s] hio_fire_intensity_si_age(io_si, cpatch%age_class) = hio_fire_intensity_si_age(io_si, cpatch%age_class) + & - cpatch%FI * cpatch%frac_burnt * cpatch%area * AREA_INV * J_per_kJ + cpatch%FI * cpatch%frac_burnt * weight * J_per_kJ ! Fuel sum [kg/m2] hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) = hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) + & - cpatch%sum_fuel * cpatch%area * AREA_INV + cpatch%sum_fuel * weight @@ -3510,8 +3519,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) n_perm2 = ccohort%n * AREA_INV + ageclass_canopy_fracarea = ccohort%c_area * AREA_INV ! WRONG; *AREA_INV should be /age_class_area hio_canopy_area_si_age(io_si,cpatch%age_class) = hio_canopy_area_si_age(io_si,cpatch%age_class) & - + ccohort%c_area * AREA_INV + + ageclass_canopy_fracarea ! calculate leaf height distribution, assuming leaf area is evenly distributed thru crown depth call CrownDepth(ccohort%height,ft,crown_depth) @@ -3605,7 +3615,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update total biomass per age bin hio_biomass_si_age(io_si,cpatch%age_class) = hio_biomass_si_age(io_si,cpatch%age_class) & - + total_m * ccohort%n * AREA_INV + + total_m * ccohort%n * AREA_INV ! Wrong; *AREA_INV should be /age_class_area if (ccohort%canopy_layer .eq. 1) then storec_canopy_scpf(i_scpf) = & From 395732d7d25fef5421bdee95875ecaa3f2892bd6 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 13:40:35 -0600 Subject: [PATCH 08/21] Change "area" to "fracarea" in *secondarylands_fracarea_si_age vars. --- main/FatesHistoryInterfaceMod.F90 | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 9bb3172cc1..3e24023c5b 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -660,7 +660,7 @@ module FatesHistoryInterfaceMod integer :: ih_c_stomata_si_age integer :: ih_c_lblayer_si_age integer :: ih_agesince_anthrodist_si_age - integer :: ih_secondarylands_area_si_age + integer :: ih_secondarylands_fracarea_si_age integer :: ih_area_burnt_si_age ! integer :: ih_fire_rate_of_spread_front_si_age integer :: ih_fire_intensity_si_age @@ -3271,7 +3271,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_zstar_si_age => this%hvars(ih_zstar_si_age)%r82d, & hio_biomass_si_age => this%hvars(ih_biomass_si_age)%r82d, & hio_agesince_anthrodist_si_age => this%hvars(ih_agesince_anthrodist_si_age)%r82d, & - hio_secondarylands_area_si_age => this%hvars(ih_secondarylands_area_si_age)%r82d, & + hio_secondarylands_fracarea_si_age => this%hvars(ih_secondarylands_fracarea_si_age)%r82d, & hio_fracarea_si_landuse => this%hvars(ih_fracarea_si_landuse)%r82d, & hio_area_burnt_si_age => this%hvars(ih_area_burnt_si_age)%r82d, & ! hio_fire_rate_of_spread_front_si_age => this%hvars(ih_fire_rate_of_spread_front_si_age)%r82d, & @@ -3447,8 +3447,8 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) + ageclass_fracarea ageclass_fracarea = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area - hio_secondarylands_area_si_age(io_si,cpatch%age_class) = & - hio_secondarylands_area_si_age(io_si,cpatch%age_class) & + hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) = & + hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) & + ageclass_fracarea endif @@ -7050,7 +7050,7 @@ subroutine define_history_vars(this, initialize_variables) long='secondary forest patch area age distribution since any kind of disturbance', & use_default='inactive', avgflag='A', vtype=site_age_r8, & hlms='CLM:ALM', upfreq=group_dyna_complx, ivar=ivar, initialize=initialize_variables, & - index=ih_secondarylands_area_si_age) + index=ih_secondarylands_fracarea_si_age) call this%set_history_var(vname='FATES_FRAGMENTATION_SCALER_SL', units='', & long='factor (0-1) by which litter/cwd fragmentation proceeds relative to max rate by soil layer', & From cdb560055d30f5d480f65089069a5fc2d80d9482 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 13:50:12 -0600 Subject: [PATCH 09/21] Align a comment. --- main/FatesHistoryInterfaceMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 3e24023c5b..e99c8ce3ee 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -5382,7 +5382,7 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) end associate endif -!!! canopy leaf carbon balance + !!! canopy leaf carbon balance ican = ccohort%canopy_layer do ileaf=1,ccohort%nv cnlf_indx = ileaf + (ican-1) * nlevleaf From 79044a7231a2148e8dd2507d8d5b931d0c2bc213 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 13:59:32 -0600 Subject: [PATCH 10/21] Refactor calculation of FATES_GPP_AP, FATES_NPP_AP. Using fates_patch_type%area_by_age allows consolidation into one block, rather than splitting the calculation across two parts of update_history_hifrq2. --- main/FatesHistoryInterfaceMod.F90 | 30 +++++++++++------------------- 1 file changed, 11 insertions(+), 19 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index e99c8ce3ee..f1509b4036 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -5189,7 +5189,7 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) real(r8) :: npp ! npp for this time-step (adjusted for g resp) [kgC/indiv/step] real(r8) :: aresp ! autotrophic respiration (adjusted for g resp) [kgC/indiv/step] real(r8) :: n_perm2 ! individuals per m2 for the whole column - real(r8) :: patch_area_by_age(nlevage) ! patch area in each bin for normalizing purposes + real(r8) :: ageclass_area ! [m2] real(r8) :: canopy_area_by_age(nlevage) ! canopy area in each bin for normalizing purposes real(r8) :: site_area_veg_inv ! 1/area of the site that is not bare-ground integer :: ipa2 ! patch incrementer @@ -5271,7 +5271,6 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ipa = 0 - patch_area_by_age(1:nlevage) = 0._r8 canopy_area_by_age(1:nlevage) = 0._r8 cpatch => sites(s)%oldest_patch @@ -5279,9 +5278,6 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ipa = ipa + 1 - patch_area_by_age(cpatch%age_class) = & - patch_area_by_age(cpatch%age_class) + cpatch%area - canopy_area_by_age(cpatch%age_class) = & canopy_area_by_age(cpatch%age_class) + cpatch%total_canopy_area @@ -5341,11 +5337,17 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ccohort%froot_mr * n_perm2 ! accumulate fluxes per patch age bin - hio_gpp_si_age(io_si,cpatch%age_class) = hio_gpp_si_age(io_si,cpatch%age_class) & - + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv + ageclass_area = sites(s)%area_by_age(cpatch%age_class) + if (ageclass_area .gt. nearzero) then + hio_gpp_si_age(io_si,cpatch%age_class) = hio_gpp_si_age(io_si,cpatch%age_class) & + + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv / ageclass_area - hio_npp_si_age(io_si,cpatch%age_class) = hio_npp_si_age(io_si,cpatch%age_class) & - + npp * ccohort%n * dt_tstep_inv + hio_npp_si_age(io_si,cpatch%age_class) = hio_npp_si_age(io_si,cpatch%age_class) & + + npp * ccohort%n * dt_tstep_inv / ageclass_area + else + hio_gpp_si_age(io_si,cpatch%age_class) = 0._r8 + hio_npp_si_age(io_si,cpatch%age_class) = 0._r8 + end if ! accumulate fluxes on canopy- and understory- separated fluxes if (ccohort%canopy_layer .eq. 1) then @@ -5578,16 +5580,6 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ! Normalize age stratified diagnostics ! ---------------------------------------------------------------- do ipa2 = 1, nlevage - if (patch_area_by_age(ipa2) .gt. nearzero) then - hio_gpp_si_age(io_si, ipa2) = & - hio_gpp_si_age(io_si, ipa2) / (patch_area_by_age(ipa2)) - hio_npp_si_age(io_si, ipa2) = & - hio_npp_si_age(io_si, ipa2) / (patch_area_by_age(ipa2)) - else - hio_gpp_si_age(io_si, ipa2) = 0._r8 - hio_npp_si_age(io_si, ipa2) = 0._r8 - endif - ! Normalize resistance diagnostics if (canopy_area_by_age(ipa2) .gt. nearzero) then hio_c_stomata_si_age(io_si,ipa2) = & From ff87ce15012a285513710d995292052c4a3fa5b2 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 17 Jul 2024 14:08:25 -0600 Subject: [PATCH 11/21] Refactor calculation of FATES_LBLAYER_COND_AP, FATES_STOMATAL_COND_AP. Moving calculation of canopy_area_by_age to an earlier patch loop allows consolidation into one block, rather than splitting the calculation across two parts of update_history_hifrq2. --- main/FatesHistoryInterfaceMod.F90 | 57 ++++++++++++------------------- 1 file changed, 22 insertions(+), 35 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index f1509b4036..fc524b181a 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -5190,6 +5190,7 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) real(r8) :: aresp ! autotrophic respiration (adjusted for g resp) [kgC/indiv/step] real(r8) :: n_perm2 ! individuals per m2 for the whole column real(r8) :: ageclass_area ! [m2] + real(r8) :: ageclass_canopy_area ! [m2] real(r8) :: canopy_area_by_age(nlevage) ! canopy area in each bin for normalizing purposes real(r8) :: site_area_veg_inv ! 1/area of the site that is not bare-ground integer :: ipa2 ! patch incrementer @@ -5256,9 +5257,12 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) call this%zero_site_hvars(sites(s), upfreq_in=group_hifr_complx) site_area_veg_inv = 0._r8 + canopy_area_by_age(1:nlevage) = 0._r8 cpatch => sites(s)%oldest_patch do while(associated(cpatch)) site_area_veg_inv = site_area_veg_inv + cpatch%total_canopy_area + canopy_area_by_age(cpatch%age_class) = & + canopy_area_by_age(cpatch%age_class) + cpatch%total_canopy_area cpatch => cpatch%younger end do !patch loop @@ -5271,26 +5275,28 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ipa = 0 - canopy_area_by_age(1:nlevage) = 0._r8 - cpatch => sites(s)%oldest_patch do while(associated(cpatch)) ipa = ipa + 1 - - canopy_area_by_age(cpatch%age_class) = & - canopy_area_by_age(cpatch%age_class) + cpatch%total_canopy_area - - - - ! Canopy resitance terms - hio_c_stomata_si_age(io_si,cpatch%age_class) = & - hio_c_stomata_si_age(io_si,cpatch%age_class) + & - cpatch%c_stomata * cpatch%total_canopy_area * mol_per_umol - - hio_c_lblayer_si_age(io_si,cpatch%age_class) = & - hio_c_lblayer_si_age(io_si,cpatch%age_class) + & - cpatch%c_lblayer * cpatch%total_canopy_area * mol_per_umol + ageclass_area = sites(s)%area_by_age(cpatch%age_class) + ageclass_canopy_area = canopy_area_by_age(cpatch%age_class) + + ! Canopy resistance terms + if (ageclass_canopy_area .gt. nearzero) then + hio_c_stomata_si_age(io_si,cpatch%age_class) = & + hio_c_stomata_si_age(io_si,cpatch%age_class) + & + cpatch%c_stomata * cpatch%total_canopy_area * mol_per_umol / & + ageclass_canopy_area + + hio_c_lblayer_si_age(io_si,cpatch%age_class) = & + hio_c_lblayer_si_age(io_si,cpatch%age_class) + & + cpatch%c_lblayer * cpatch%total_canopy_area * mol_per_umol / & + ageclass_canopy_area + else + hio_c_stomata_si_age(io_si,cpatch%age_class) = 0._r8 + hio_c_lblayer_si_age(io_si,cpatch%age_class) = 0._r8 + end if ccohort => cpatch%shortest do while(associated(ccohort)) @@ -5337,7 +5343,6 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) ccohort%froot_mr * n_perm2 ! accumulate fluxes per patch age bin - ageclass_area = sites(s)%area_by_age(cpatch%age_class) if (ageclass_area .gt. nearzero) then hio_gpp_si_age(io_si,cpatch%age_class) = hio_gpp_si_age(io_si,cpatch%age_class) & + ccohort%gpp_tstep * ccohort%n * dt_tstep_inv / ageclass_area @@ -5576,24 +5581,6 @@ subroutine update_history_hifrq2(this,nc,nsites,sites,bc_in,bc_out,dt_tstep) end do do_ican2 - - ! Normalize age stratified diagnostics - ! ---------------------------------------------------------------- - do ipa2 = 1, nlevage - ! Normalize resistance diagnostics - if (canopy_area_by_age(ipa2) .gt. nearzero) then - hio_c_stomata_si_age(io_si,ipa2) = & - hio_c_stomata_si_age(io_si,ipa2) / canopy_area_by_age(ipa2) - - hio_c_lblayer_si_age(io_si,ipa2) = & - hio_c_lblayer_si_age(io_si,ipa2) / canopy_area_by_age(ipa2) - else - hio_c_stomata_si_age(io_si,ipa2) = 0._r8 - hio_c_lblayer_si_age(io_si,ipa2) = 0._r8 - end if - - end do - enddo do_sites ! site loop end associate From 5942a0df173a7bca5b367bbe005cc6d7cd7fd433 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 12:46:56 -0600 Subject: [PATCH 12/21] Fix weighting of most _AP variables. --- main/FatesHistoryInterfaceMod.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index fc524b181a..13f347b6b0 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3431,7 +3431,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) if ( ED_val_comp_excln .lt. 0._r8 ) then ! only valid when "strict ppa" enabled - weight = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area + weight = cpatch%area / age_class_area hio_zstar_si_age(io_si,cpatch%age_class) = hio_zstar_si_age(io_si,cpatch%age_class) & + cpatch%zstar * weight endif @@ -3484,7 +3484,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) end do ! Weight for the following per-age-class fire variables - weight = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area + weight = cpatch%area / age_class_area ! fractional area burnt [frac/day] -> [frac/sec] hio_area_burnt_si_age(io_si,cpatch%age_class) = hio_area_burnt_si_age(io_si,cpatch%age_class) + & @@ -3519,7 +3519,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) n_perm2 = ccohort%n * AREA_INV - ageclass_canopy_fracarea = ccohort%c_area * AREA_INV ! WRONG; *AREA_INV should be /age_class_area + ageclass_canopy_fracarea = ccohort%c_area / age_class_area hio_canopy_area_si_age(io_si,cpatch%age_class) = hio_canopy_area_si_age(io_si,cpatch%age_class) & + ageclass_canopy_fracarea @@ -3615,7 +3615,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update total biomass per age bin hio_biomass_si_age(io_si,cpatch%age_class) = hio_biomass_si_age(io_si,cpatch%age_class) & - + total_m * ccohort%n * AREA_INV ! Wrong; *AREA_INV should be /age_class_area + + total_m * ccohort%n / age_class_area if (ccohort%canopy_layer .eq. 1) then storec_canopy_scpf(i_scpf) = & From 66cc4f81b5ec20f54a6d635aee0143138af36158 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 12:54:26 -0600 Subject: [PATCH 13/21] FATES_SECONDAREA_(ANTHRO)DIST_AP are actually fine. --- main/FatesHistoryInterfaceMod.F90 | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 13f347b6b0..41c5d0d9ee 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3441,12 +3441,18 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ageclass_since_anthrodist = get_age_class_index(cpatch%age_since_anthro_disturbance) - ageclass_fracarea = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / sites(s)%age_class_area(ageclass_since_anthrodist) + ! The value for each age class is the fraction of the site that is secondary land + ! of that age class. So summing across age classes in the output file should give + ! you the fraction of each site that's secondary land. If instead we want this to + ! sum to 1 (so that each age class gives the fraction of secondary land in the + ! age class), then instead of multiplying by AREA_INV, divide by total secondary + ! land area. + ageclass_fracarea = cpatch%area * AREA_INV + hio_agesince_anthrodist_si_age(io_si,ageclass_since_anthrodist) = & hio_agesince_anthrodist_si_age(io_si,ageclass_since_anthrodist) & + ageclass_fracarea - ageclass_fracarea = cpatch%area * AREA_INV ! WRONG; should be: cpatch%area / age_class_area hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) = & hio_secondarylands_fracarea_si_age(io_si,cpatch%age_class) & + ageclass_fracarea From 99adff994af7deb5ce610a97dc4a5ee02e19134d Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 15:10:36 -0600 Subject: [PATCH 14/21] Fix denominator of FATES_NPP_APPF and FATES_VEGC_APPF. Instead of dividing by total site area, divide by age class area. --- main/FatesHistoryInterfaceMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 41c5d0d9ee..5bbe0fd646 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3971,10 +3971,10 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) iagepft = get_agepft_class_index(cpatch%age,ccohort%pft) hio_npp_si_agepft(io_si,iagepft) = hio_npp_si_agepft(io_si,iagepft) + & - ccohort%n * ccohort%npp_acc_hold * AREA_INV / days_per_year / sec_per_day + ccohort%n * ccohort%npp_acc_hold / age_class_area / days_per_year / sec_per_day hio_biomass_si_agepft(io_si,iagepft) = hio_biomass_si_agepft(io_si,iagepft) + & - total_m * ccohort%n * AREA_INV + total_m * ccohort%n / age_class_area ! update SCPF/SCLS- and canopy/subcanopy- partitioned quantities canlayer: if (ccohort%canopy_layer .eq. 1) then From 2ee7502ee3f63badb4e5fb5d1ec61bbc6b9f4b3f Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 15:28:13 -0600 Subject: [PATCH 15/21] Refactor size/age output variables to share scag_denominator_area. --- main/FatesHistoryInterfaceMod.F90 | 18 +++++++++++------- 1 file changed, 11 insertions(+), 7 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 5bbe0fd646..376155525e 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3088,6 +3088,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: weight real(r8) :: denominator real(r8) :: age_class_area ! [m2] + real(r8) :: scag_denominator_area ! Denominator for variables partitioned by size and age class [m2] integer :: i_dist, j_dist @@ -3400,6 +3401,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) cpatch%age_class = get_age_class_index(cpatch%age) age_class_area = sites(s)%area_by_age(cpatch%age_class) + ! Get denominators for non-site-level variables + scag_denominator_area = m2_per_ha + ! Increment the fractional area in each age class bin hio_fracarea_si_age(io_si,cpatch%age_class) = hio_fracarea_si_age(io_si,cpatch%age_class) & + cpatch%area * AREA_INV @@ -3957,7 +3961,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) iscag = get_sizeage_class_index(ccohort%dbh,cpatch%age) - hio_nplant_si_scag(io_si,iscag) = hio_nplant_si_scag(io_si,iscag) + ccohort%n / m2_per_ha + hio_nplant_si_scag(io_si,iscag) = hio_nplant_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_nplant_si_scls(io_si,scls) = hio_nplant_si_scls(io_si,scls) + ccohort%n / m2_per_ha @@ -3978,12 +3982,12 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update SCPF/SCLS- and canopy/subcanopy- partitioned quantities canlayer: if (ccohort%canopy_layer .eq. 1) then - hio_nplant_canopy_si_scag(io_si,iscag) = hio_nplant_canopy_si_scag(io_si,iscag) + ccohort%n / m2_per_ha + hio_nplant_canopy_si_scag(io_si,iscag) = hio_nplant_canopy_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & (ccohort%bmort + ccohort%hmort + ccohort%cmort + & ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha hio_ddbh_canopy_si_scag(io_si,iscag) = hio_ddbh_canopy_si_scag(io_si,iscag) + & - ccohort%ddbhdt*ccohort%n * m_per_cm / m2_per_ha + ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area hio_bstor_canopy_si_scpf(io_si,scpf) = hio_bstor_canopy_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_canopy_si_scpf(io_si,scpf) = hio_bleaf_canopy_si_scpf(io_si,scpf) + & @@ -4104,12 +4108,12 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) else canlayer - hio_nplant_understory_si_scag(io_si,iscag) = hio_nplant_understory_si_scag(io_si,iscag) + ccohort%n / m2_per_ha + hio_nplant_understory_si_scag(io_si,iscag) = hio_nplant_understory_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & (ccohort%bmort + ccohort%hmort + ccohort%cmort + & ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha hio_ddbh_understory_si_scag(io_si,iscag) = hio_ddbh_understory_si_scag(io_si,iscag) + & - ccohort%ddbhdt*ccohort%n * m_per_cm / m2_per_ha + ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area hio_bstor_understory_si_scpf(io_si,scpf) = hio_bstor_understory_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_understory_si_scpf(io_si,scpf) = hio_bleaf_understory_si_scpf(io_si,scpf) + & @@ -4456,9 +4460,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! for scag variables, also treat as happening in the newly-disurbed patch hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_canopy(i_scls, ft) / m2_per_ha + sites(s)%fmort_rate_canopy(i_scls, ft) / scag_denominator_area hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_ustory(i_scls, ft) / m2_per_ha + sites(s)%fmort_rate_ustory(i_scls, ft) / scag_denominator_area ! while in this loop, pass the fusion-induced growth rate flux to history hio_growthflux_fusion_si_scpf(io_si,i_scpf) = hio_growthflux_fusion_si_scpf(io_si,i_scpf) + & From 87d3bfe8abca044a427fd2ba3f9a311dbb6f262b Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 15:32:30 -0600 Subject: [PATCH 16/21] Tiny refactor of FATES_FUEL_AMOUNT_APFC. --- main/FatesHistoryInterfaceMod.F90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 376155525e..57429dfa4d 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3089,6 +3089,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: denominator real(r8) :: age_class_area ! [m2] real(r8) :: scag_denominator_area ! Denominator for variables partitioned by size and age class [m2] + real(r8) :: apfc_denominator_area ! Denominator for variables partitioned by age and fuel class [m2] integer :: i_dist, j_dist @@ -3403,6 +3404,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! Get denominators for non-site-level variables scag_denominator_area = m2_per_ha + apfc_denominator_area = AREA ! Increment the fractional area in each age class bin hio_fracarea_si_age(io_si,cpatch%age_class) = hio_fracarea_si_age(io_si,cpatch%age_class) & @@ -4294,7 +4296,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) i_agefuel = get_agefuel_class_index(cpatch%age,i_fuel) hio_fuel_amount_age_fuel(io_si,i_agefuel) = hio_fuel_amount_age_fuel(io_si,i_agefuel) + & - cpatch%fuel_frac(i_fuel) * cpatch%sum_fuel * cpatch%area * AREA_INV + cpatch%fuel_frac(i_fuel) * cpatch%sum_fuel * cpatch%area / apfc_denominator_area hio_litter_moisture_si_fuel(io_si, i_fuel) = hio_litter_moisture_si_fuel(io_si, i_fuel) + & cpatch%litter_moisture(i_fuel) * cpatch%area * AREA_INV From 0bb23b9d86b27b30aec42a279515baad44ecffe4 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 15:34:40 -0600 Subject: [PATCH 17/21] Tiny refactor of FATES_NPLANT_SZAPPF. --- main/FatesHistoryInterfaceMod.F90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 57429dfa4d..bc3dcd0ccf 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3089,6 +3089,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) real(r8) :: denominator real(r8) :: age_class_area ! [m2] real(r8) :: scag_denominator_area ! Denominator for variables partitioned by size and age class [m2] + real(r8) :: scagpft_denominator_area ! Denominator for variables partitioned by size class, age class, and PFT [m2] real(r8) :: apfc_denominator_area ! Denominator for variables partitioned by age and fuel class [m2] integer :: i_dist, j_dist @@ -3404,6 +3405,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! Get denominators for non-site-level variables scag_denominator_area = m2_per_ha + scagpft_denominator_area = m2_per_ha apfc_denominator_area = AREA ! Increment the fractional area in each age class bin @@ -3971,7 +3973,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update size, age, and PFT - indexed quantities iscagpft = get_sizeagepft_class_index(ccohort%dbh,cpatch%age,ccohort%pft) - hio_nplant_si_scagpft(io_si,iscagpft) = hio_nplant_si_scagpft(io_si,iscagpft) + ccohort%n / m2_per_ha + hio_nplant_si_scagpft(io_si,iscagpft) = hio_nplant_si_scagpft(io_si,iscagpft) + ccohort%n / scagpft_denominator_area ! update age and PFT - indexed quantities iagepft = get_agepft_class_index(cpatch%age,ccohort%pft) From 1ec6d6eb9bc1e73ed9af6b810edb0aa7e7172a31 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Thu, 18 Jul 2024 15:37:42 -0600 Subject: [PATCH 18/21] All history vars with age class multiplexed should divide by age class area, not site area. --- main/FatesHistoryInterfaceMod.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index bc3dcd0ccf..d99302a60e 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3404,9 +3404,9 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) age_class_area = sites(s)%area_by_age(cpatch%age_class) ! Get denominators for non-site-level variables - scag_denominator_area = m2_per_ha - scagpft_denominator_area = m2_per_ha - apfc_denominator_area = AREA + scag_denominator_area = age_class_area + scagpft_denominator_area = age_class_area + apfc_denominator_area = age_class_area ! Increment the fractional area in each age class bin hio_fracarea_si_age(io_si,cpatch%age_class) = hio_fracarea_si_age(io_si,cpatch%age_class) & From 103fdc96397b340f2c14a1f5b47bd8dbe479019a Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Tue, 23 Jul 2024 13:59:37 -0600 Subject: [PATCH 19/21] Ensure that no area denominator is zero. --- main/FatesHistoryInterfaceMod.F90 | 123 +++++++++++++++++++++--------- 1 file changed, 86 insertions(+), 37 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index d99302a60e..d0dddc45d1 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -3439,9 +3439,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) if ( ED_val_comp_excln .lt. 0._r8 ) then ! only valid when "strict ppa" enabled - weight = cpatch%area / age_class_area - hio_zstar_si_age(io_si,cpatch%age_class) = hio_zstar_si_age(io_si,cpatch%age_class) & - + cpatch%zstar * weight + if (age_class_area .gt. nearzero) then + weight = cpatch%area / age_class_area + hio_zstar_si_age(io_si,cpatch%age_class) = hio_zstar_si_age(io_si,cpatch%age_class) & + + cpatch%zstar * weight + else + hio_zstar_si_age(io_si,cpatch%age_class) = 0._r8 + end if endif ! some diagnostics on secondary forest area and its age distribution @@ -3498,22 +3502,29 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) end do ! Weight for the following per-age-class fire variables - weight = cpatch%area / age_class_area + if (age_class_area .gt. nearzero) then + weight = cpatch%area / age_class_area - ! fractional area burnt [frac/day] -> [frac/sec] - hio_area_burnt_si_age(io_si,cpatch%age_class) = hio_area_burnt_si_age(io_si,cpatch%age_class) + & - cpatch%frac_burnt * weight / sec_per_day + ! fractional area burnt [frac/day] -> [frac/sec] + hio_area_burnt_si_age(io_si,cpatch%age_class) = hio_area_burnt_si_age(io_si,cpatch%age_class) + & + cpatch%frac_burnt * weight / sec_per_day - ! hio_fire_rate_of_spread_front_si_age(io_si, cpatch%age_class) = hio_fire_rate_of_spread_si_age(io_si, cpatch%age_class) + & - ! cpatch%ros_front * cpatch*frac_burnt * weight + ! hio_fire_rate_of_spread_front_si_age(io_si, cpatch%age_class) = hio_fire_rate_of_spread_si_age(io_si, cpatch%age_class) + & + ! cpatch%ros_front * cpatch*frac_burnt * weight - ! Fire intensity weighted by burned fraction [kJ/m/s] -> [J/m/s] - hio_fire_intensity_si_age(io_si, cpatch%age_class) = hio_fire_intensity_si_age(io_si, cpatch%age_class) + & - cpatch%FI * cpatch%frac_burnt * weight * J_per_kJ + ! Fire intensity weighted by burned fraction [kJ/m/s] -> [J/m/s] + hio_fire_intensity_si_age(io_si, cpatch%age_class) = hio_fire_intensity_si_age(io_si, cpatch%age_class) + & + cpatch%FI * cpatch%frac_burnt * weight * J_per_kJ - ! Fuel sum [kg/m2] - hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) = hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) + & - cpatch%sum_fuel * weight + ! Fuel sum [kg/m2] + hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) = hio_fire_sum_fuel_si_age(io_si, cpatch%age_class) + & + cpatch%sum_fuel * weight + else + hio_area_burnt_si_age(io_si,cpatch%age_class) = 0._r8 + ! hio_fire_rate_of_spread_front_si_age(io_si,cpatch%age_class) = 0._r8 + hio_fire_intensity_si_age(io_si,cpatch%age_class) = 0._r8 + hio_fire_sum_fuel_si_age(io_si,cpatch%age_class) = 0._r8 + end if @@ -3533,7 +3544,11 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) n_perm2 = ccohort%n * AREA_INV - ageclass_canopy_fracarea = ccohort%c_area / age_class_area + if (age_class_area .gt. nearzero) then + ageclass_canopy_fracarea = ccohort%c_area / age_class_area + else + ageclass_canopy_fracarea = 0._r8 + end if hio_canopy_area_si_age(io_si,cpatch%age_class) = hio_canopy_area_si_age(io_si,cpatch%age_class) & + ageclass_canopy_fracarea @@ -3628,8 +3643,12 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) (ccohort%n * AREA_INV) * total_m ! update total biomass per age bin - hio_biomass_si_age(io_si,cpatch%age_class) = hio_biomass_si_age(io_si,cpatch%age_class) & - + total_m * ccohort%n / age_class_area + if (age_class_area .gt. nearzero) then + hio_biomass_si_age(io_si,cpatch%age_class) = hio_biomass_si_age(io_si,cpatch%age_class) & + + total_m * ccohort%n / age_class_area + else + hio_biomass_si_age(io_si,cpatch%age_class) = 0._r8 + end if if (ccohort%canopy_layer .eq. 1) then storec_canopy_scpf(i_scpf) = & @@ -3965,7 +3984,11 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) iscag = get_sizeage_class_index(ccohort%dbh,cpatch%age) - hio_nplant_si_scag(io_si,iscag) = hio_nplant_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area + if (scag_denominator_area .gt. nearzero) then + hio_nplant_si_scag(io_si,iscag) = hio_nplant_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area + else + hio_nplant_si_scag(io_si,iscag) = 0._r8 + end if hio_nplant_si_scls(io_si,scls) = hio_nplant_si_scls(io_si,scls) + ccohort%n / m2_per_ha @@ -3973,32 +3996,45 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update size, age, and PFT - indexed quantities iscagpft = get_sizeagepft_class_index(ccohort%dbh,cpatch%age,ccohort%pft) - hio_nplant_si_scagpft(io_si,iscagpft) = hio_nplant_si_scagpft(io_si,iscagpft) + ccohort%n / scagpft_denominator_area + if (scagpft_denominator_area .gt. nearzero) then + hio_nplant_si_scagpft(io_si,iscagpft) = hio_nplant_si_scagpft(io_si,iscagpft) + ccohort%n / scagpft_denominator_area + else + hio_nplant_si_scagpft(io_si,iscagpft) = 0._r8 + end if ! update age and PFT - indexed quantities iagepft = get_agepft_class_index(cpatch%age,ccohort%pft) - hio_npp_si_agepft(io_si,iagepft) = hio_npp_si_agepft(io_si,iagepft) + & - ccohort%n * ccohort%npp_acc_hold / age_class_area / days_per_year / sec_per_day + if (age_class_area .gt. nearzero) then + hio_npp_si_agepft(io_si,iagepft) = hio_npp_si_agepft(io_si,iagepft) + & + ccohort%n * ccohort%npp_acc_hold / age_class_area / days_per_year / sec_per_day - hio_biomass_si_agepft(io_si,iagepft) = hio_biomass_si_agepft(io_si,iagepft) + & - total_m * ccohort%n / age_class_area + hio_biomass_si_agepft(io_si,iagepft) = hio_biomass_si_agepft(io_si,iagepft) + & + total_m * ccohort%n / age_class_area + else + hio_npp_si_agepft(io_si,iagepft) = 0._r8 + hio_biomass_si_agepft(io_si,iagepft) = 0._r8 + end if ! update SCPF/SCLS- and canopy/subcanopy- partitioned quantities canlayer: if (ccohort%canopy_layer .eq. 1) then - hio_nplant_canopy_si_scag(io_si,iscag) = hio_nplant_canopy_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & (ccohort%bmort + ccohort%hmort + ccohort%cmort + & ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha - hio_ddbh_canopy_si_scag(io_si,iscag) = hio_ddbh_canopy_si_scag(io_si,iscag) + & - ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area hio_bstor_canopy_si_scpf(io_si,scpf) = hio_bstor_canopy_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_canopy_si_scpf(io_si,scpf) = hio_bleaf_canopy_si_scpf(io_si,scpf) + & leaf_m * ccohort%n / m2_per_ha hio_lai_canopy_si_scpf(io_si,scpf) = hio_lai_canopy_si_scpf(io_si,scpf) + & ccohort%treelai*ccohort%c_area * AREA_INV - + if (scag_denominator_area .gt. nearzero) then + hio_nplant_canopy_si_scag(io_si,iscag) = hio_nplant_canopy_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area + hio_ddbh_canopy_si_scag(io_si,iscag) = hio_ddbh_canopy_si_scag(io_si,iscag) + & + ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area + else + hio_nplant_canopy_si_scag(io_si,iscag) = 0._r8 + hio_ddbh_canopy_si_scag(io_si,iscag) = 0._r8 + end if hio_crownarea_canopy_si_scpf(io_si,scpf) = hio_crownarea_canopy_si_scpf(io_si,scpf) + & ccohort%c_area * AREA_INV !hio_mortality_canopy_si_scpf(io_si,scpf) = hio_mortality_canopy_si_scpf(io_si,scpf)+ & @@ -4112,16 +4148,19 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) else canlayer - hio_nplant_understory_si_scag(io_si,iscag) = hio_nplant_understory_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & (ccohort%bmort + ccohort%hmort + ccohort%cmort + & ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha - hio_ddbh_understory_si_scag(io_si,iscag) = hio_ddbh_understory_si_scag(io_si,iscag) + & - ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area hio_bstor_understory_si_scpf(io_si,scpf) = hio_bstor_understory_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_understory_si_scpf(io_si,scpf) = hio_bleaf_understory_si_scpf(io_si,scpf) + & leaf_m * ccohort%n / m2_per_ha + if (scag_denominator_area .gt. nearzero) then + hio_nplant_understory_si_scag(io_si,iscag) = hio_nplant_understory_si_scag(io_si,iscag) + & + ccohort%n / scag_denominator_area + hio_ddbh_understory_si_scag(io_si,iscag) = hio_ddbh_understory_si_scag(io_si,iscag) + & + ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area + end if hio_lai_understory_si_scpf(io_si,scpf) = hio_lai_understory_si_scpf(io_si,scpf) + & ccohort%treelai*ccohort%c_area * AREA_INV @@ -4297,8 +4336,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) do i_fuel = 1,nfsc i_agefuel = get_agefuel_class_index(cpatch%age,i_fuel) - hio_fuel_amount_age_fuel(io_si,i_agefuel) = hio_fuel_amount_age_fuel(io_si,i_agefuel) + & - cpatch%fuel_frac(i_fuel) * cpatch%sum_fuel * cpatch%area / apfc_denominator_area + + if (apfc_denominator_area .gt. nearzero) then + hio_fuel_amount_age_fuel(io_si,i_agefuel) = hio_fuel_amount_age_fuel(io_si,i_agefuel) + & + cpatch%fuel_frac(i_fuel) * cpatch%sum_fuel * cpatch%area / apfc_denominator_area + else + hio_fuel_amount_age_fuel(io_si,i_agefuel) = 0._r8 + end if hio_litter_moisture_si_fuel(io_si, i_fuel) = hio_litter_moisture_si_fuel(io_si, i_fuel) + & cpatch%litter_moisture(i_fuel) * cpatch%area * AREA_INV @@ -4463,10 +4507,15 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! ! for scag variables, also treat as happening in the newly-disurbed patch - hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_canopy(i_scls, ft) / scag_denominator_area - hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_ustory(i_scls, ft) / scag_denominator_area + if (scag_denominator_area .gt. nearzero) then + hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & + sites(s)%fmort_rate_canopy(i_scls, ft) / scag_denominator_area + hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & + sites(s)%fmort_rate_ustory(i_scls, ft) / scag_denominator_area + else + hio_mortality_canopy_si_scag(io_si,iscag) = 0._r8 + hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 + end if ! while in this loop, pass the fusion-induced growth rate flux to history hio_growthflux_fusion_si_scpf(io_si,i_scpf) = hio_growthflux_fusion_si_scpf(io_si,i_scpf) + & From a0881c5363eb4d52d296c6962677074ae82f9a66 Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 24 Jul 2024 10:11:13 -0600 Subject: [PATCH 20/21] Fix FATES_MORTALITY_CANOPY_SZAP and FATES_MORTALITY_USTORY_SZAP. --- main/FatesHistoryInterfaceMod.F90 | 25 +++++++++++++++++-------- 1 file changed, 17 insertions(+), 8 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index d0dddc45d1..ed4db1c73f 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -4018,9 +4018,6 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update SCPF/SCLS- and canopy/subcanopy- partitioned quantities canlayer: if (ccohort%canopy_layer .eq. 1) then - hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & - (ccohort%bmort + ccohort%hmort + ccohort%cmort + & - ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha hio_bstor_canopy_si_scpf(io_si,scpf) = hio_bstor_canopy_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_canopy_si_scpf(io_si,scpf) = hio_bleaf_canopy_si_scpf(io_si,scpf) + & @@ -4028,12 +4025,16 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_lai_canopy_si_scpf(io_si,scpf) = hio_lai_canopy_si_scpf(io_si,scpf) + & ccohort%treelai*ccohort%c_area * AREA_INV if (scag_denominator_area .gt. nearzero) then + hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & + (ccohort%bmort + ccohort%hmort + ccohort%cmort + & + ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / scag_denominator_area hio_nplant_canopy_si_scag(io_si,iscag) = hio_nplant_canopy_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area hio_ddbh_canopy_si_scag(io_si,iscag) = hio_ddbh_canopy_si_scag(io_si,iscag) + & ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area else hio_nplant_canopy_si_scag(io_si,iscag) = 0._r8 hio_ddbh_canopy_si_scag(io_si,iscag) = 0._r8 + hio_mortality_canopy_si_scag(io_si,iscag) = 0._r8 end if hio_crownarea_canopy_si_scpf(io_si,scpf) = hio_crownarea_canopy_si_scpf(io_si,scpf) + & ccohort%c_area * AREA_INV @@ -4148,9 +4149,6 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) else canlayer - hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - (ccohort%bmort + ccohort%hmort + ccohort%cmort + & - ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / m2_per_ha hio_bstor_understory_si_scpf(io_si,scpf) = hio_bstor_understory_si_scpf(io_si,scpf) + & store_m * ccohort%n / m2_per_ha hio_bleaf_understory_si_scpf(io_si,scpf) = hio_bleaf_understory_si_scpf(io_si,scpf) + & @@ -4160,6 +4158,13 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ccohort%n / scag_denominator_area hio_ddbh_understory_si_scag(io_si,iscag) = hio_ddbh_understory_si_scag(io_si,iscag) + & ccohort%ddbhdt*ccohort%n * m_per_cm / scag_denominator_area + hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & + (ccohort%bmort + ccohort%hmort + ccohort%cmort + & + ccohort%frmort + ccohort%smort + ccohort%asmort + ccohort%dgmort) * ccohort%n / scag_denominator_area + else + hio_nplant_understory_si_scag(io_si,iscag) = 0._r8 + hio_ddbh_understory_si_scag(io_si,iscag) = 0._r8 + hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 end if hio_lai_understory_si_scpf(io_si,scpf) = hio_lai_understory_si_scpf(io_si,scpf) + & @@ -4471,8 +4476,12 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) sites(s)%imort_rate(i_scls, ft) / m2_per_ha ! iscag = i_scls ! since imort is by definition something that only happens in newly disturbed patches, treat as such - hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - sites(s)%imort_rate(i_scls, ft) / m2_per_ha + if (scag_denominator_area .gt. nearzero) then + hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & + sites(s)%imort_rate(i_scls, ft) / scag_denominator_area + else + hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 + end if ! fire mortality from the site-level diagnostic rates hio_m5_si_scpf(io_si,i_scpf) = (sites(s)%fmort_rate_canopy(i_scls, ft) + & From f21fa95bbc4d76a48aceaa7738e6dcba68b7be1c Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 24 Jul 2024 12:48:51 -0600 Subject: [PATCH 21/21] scag_denominator_area needs to be in patchloop. Should fix FATES_MORTALITY_CANOPY_SZAP and FATES_MORTALITY_USTORY_SZAP. --- main/FatesHistoryInterfaceMod.F90 | 47 +++++++++++++++-------------- main/FatesSizeAgeTypeIndicesMod.F90 | 23 +++++++++++--- 2 files changed, 44 insertions(+), 26 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index ed4db1c73f..ef67717d01 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -125,7 +125,8 @@ module FatesHistoryInterfaceMod use FatesLitterMod , only : ncwd use FatesConstantsMod , only : ican_upper use FatesConstantsMod , only : ican_ustory - use FatesSizeAgeTypeIndicesMod, only : get_sizeage_class_index + use FatesSizeAgeTypeIndicesMod, only : get_sizeage_class_index_from_classes + use FatesSizeAgeTypeIndicesMod, only : get_sizeage_class_index_from_dbh_age use FatesSizeAgeTypeIndicesMod, only : get_sizeagepft_class_index use FatesSizeAgeTypeIndicesMod, only : get_agepft_class_index use FatesSizeAgeTypeIndicesMod, only : get_agefuel_class_index @@ -3471,6 +3472,28 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) endif + ! sizeclass-patchage-resolved outputs representing mortality variables specific to size class and PFT + do i_scls = 1,nlevsclass + iscag = get_sizeage_class_index_from_classes(i_scls, cpatch%age_class) + if (scag_denominator_area .gt. nearzero) then + do ft = 1,numpft + ! Impact mortality + hio_mortality_understory_si_scag(io_si,iscag) = & + hio_mortality_understory_si_scag(io_si,iscag) + & + sites(s)%imort_rate(i_scls, ft) / scag_denominator_area + ! Fire mortality + hio_mortality_canopy_si_scag(io_si,iscag) = & + hio_mortality_canopy_si_scag(io_si,iscag) + & + sites(s)%fmort_rate_canopy(i_scls, ft) / scag_denominator_area + hio_mortality_understory_si_scag(io_si,iscag) = & + hio_mortality_understory_si_scag(io_si,iscag) + & + sites(s)%fmort_rate_ustory(i_scls, ft) / scag_denominator_area + end do + else + hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 + hio_mortality_canopy_si_scag(io_si,iscag) = 0._r8 + end if + end do ! patch-age-resolved fire variables do ft = 1,numpft @@ -3982,7 +4005,7 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) ! update size-class x patch-age related quantities - iscag = get_sizeage_class_index(ccohort%dbh,cpatch%age) + iscag = get_sizeage_class_index_from_dbh_age(ccohort%dbh,cpatch%age) if (scag_denominator_area .gt. nearzero) then hio_nplant_si_scag(io_si,iscag) = hio_nplant_si_scag(io_si,iscag) + ccohort%n / scag_denominator_area @@ -4475,13 +4498,6 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_mortality_understory_si_scls(io_si,i_scls) = hio_mortality_understory_si_scls(io_si,i_scls) + & sites(s)%imort_rate(i_scls, ft) / m2_per_ha ! - iscag = i_scls ! since imort is by definition something that only happens in newly disturbed patches, treat as such - if (scag_denominator_area .gt. nearzero) then - hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - sites(s)%imort_rate(i_scls, ft) / scag_denominator_area - else - hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 - end if ! fire mortality from the site-level diagnostic rates hio_m5_si_scpf(io_si,i_scpf) = (sites(s)%fmort_rate_canopy(i_scls, ft) + & @@ -4513,19 +4529,6 @@ subroutine update_history_dyn2(this,nc,nsites,sites,bc_in) hio_mortality_understory_si_scls(io_si,i_scls) = hio_mortality_understory_si_scls(io_si,i_scls) + & sites(s)%fmort_rate_ustory(i_scls, ft) / m2_per_ha - ! - ! for scag variables, also treat as happening in the newly-disurbed patch - - if (scag_denominator_area .gt. nearzero) then - hio_mortality_canopy_si_scag(io_si,iscag) = hio_mortality_canopy_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_canopy(i_scls, ft) / scag_denominator_area - hio_mortality_understory_si_scag(io_si,iscag) = hio_mortality_understory_si_scag(io_si,iscag) + & - sites(s)%fmort_rate_ustory(i_scls, ft) / scag_denominator_area - else - hio_mortality_canopy_si_scag(io_si,iscag) = 0._r8 - hio_mortality_understory_si_scag(io_si,iscag) = 0._r8 - end if - ! while in this loop, pass the fusion-induced growth rate flux to history hio_growthflux_fusion_si_scpf(io_si,i_scpf) = hio_growthflux_fusion_si_scpf(io_si,i_scpf) + & sites(s)%growthflux_fusion(i_scls, ft) * days_per_year / m2_per_ha diff --git a/main/FatesSizeAgeTypeIndicesMod.F90 b/main/FatesSizeAgeTypeIndicesMod.F90 index 7945bef035..4a14f46eda 100644 --- a/main/FatesSizeAgeTypeIndicesMod.F90 +++ b/main/FatesSizeAgeTypeIndicesMod.F90 @@ -19,7 +19,8 @@ module FatesSizeAgeTypeIndicesMod ! Make public necessary subroutines and functions public :: get_age_class_index - public :: get_sizeage_class_index + public :: get_sizeage_class_index_from_classes + public :: get_sizeage_class_index_from_dbh_age public :: sizetype_class_index public :: get_size_class_index public :: get_height_index @@ -75,7 +76,21 @@ end function get_layersizetype_class_index ! ===================================================================================== - function get_sizeage_class_index(dbh,age) result(size_by_age_class) + function get_sizeage_class_index_from_classes(size_class, age_class) result(size_by_age_class) + + ! Arguments + integer, intent(in) :: size_class + integer, intent(in) :: age_class + + integer :: size_by_age_class + + size_by_age_class = (age_class-1)*nlevsclass + size_class + + end function get_sizeage_class_index_from_classes + + ! ===================================================================================== + + function get_sizeage_class_index_from_dbh_age(dbh,age) result(size_by_age_class) ! Arguments real(r8),intent(in) :: dbh @@ -89,9 +104,9 @@ function get_sizeage_class_index(dbh,age) result(size_by_age_class) age_class = get_age_class_index(age) - size_by_age_class = (age_class-1)*nlevsclass + size_class + size_by_age_class = get_sizeage_class_index_from_classes(size_class, age_class) - end function get_sizeage_class_index + end function get_sizeage_class_index_from_dbh_age !======================================================================================