diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 4c5985fe15..acf3d76577 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -1016,8 +1016,10 @@ subroutine canopy_summarization( nsites, sites, bc_in ) call carea_allom(currentCohort%dbh,currentCohort%n,sites(s)%spread,& currentCohort%pft,currentCohort%c_area) - currentCohort%treelai = tree_lai(currentCohort%bl, currentCohort%status_coh, & - currentCohort%pft, currentCohort%c_area, currentCohort%n ) + + currentCohort%treelai = tree_lai(currentCohort%bl, & + currentCohort%pft, currentCohort%c_area, currentCohort%n, & + currentCohort%canopy_layer, currentPatch%canopy_layer_tlai ) canopy_leaf_area = canopy_leaf_area + currentCohort%treelai *currentCohort%c_area @@ -1087,7 +1089,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) ! ! The following patch level diagnostics are updated here: ! - ! currentPatch%canopy_layer_tai(cl) ! TAI of each canopy layer + ! currentPatch%canopy_layer_tlai(cl) ! total leaf area index of canopy layer ! currentPatch%ncan(cl,ft) ! number of vegetation layers needed ! ! in this patch's pft/canopy-layer ! currentPatch%nrad(cl,ft) ! same as ncan, but does not include @@ -1157,7 +1159,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) ! calculate tree lai and sai. ! -------------------------------------------------------------------------------- - currentPatch%canopy_layer_tai(:) = 0._r8 + currentPatch%canopy_layer_tlai(:) = 0._r8 currentPatch%ncan(:,:) = 0 currentPatch%nrad(:,:) = 0 patch_lai = 0._r8 @@ -1176,38 +1178,38 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) if (currentPatch%total_canopy_area > nearzero ) then - currentCohort => currentPatch%shortest + + currentCohort => currentPatch%tallest do while(associated(currentCohort)) ft = currentCohort%pft cl = currentCohort%canopy_layer - currentCohort%treelai = tree_lai(currentCohort%bl, currentCohort%status_coh, currentCohort%pft, & - currentCohort%c_area, currentCohort%n ) - currentCohort%treesai = tree_sai(currentCohort%dbh, currentCohort%pft, currentCohort%canopy_trim, & - currentCohort%c_area, currentCohort%n) + ! Calculate LAI of layers above + ! Note that the canopy_layer_lai is also calculated in this loop + ! but since we go top down in terms of plant size, we should be okay + + currentCohort%treelai = tree_lai(currentCohort%bl, currentCohort%pft, currentCohort%c_area, & + currentCohort%n, currentCohort%canopy_layer, & + currentPatch%canopy_layer_tlai ) + + currentCohort%treesai = tree_sai(currentCohort%pft, currentCohort%dbh, currentCohort%canopy_trim, & + currentCohort%c_area, currentCohort%n, currentCohort%canopy_layer, & + currentPatch%canopy_layer_tlai, currentCohort%treelai ) currentCohort%lai = currentCohort%treelai *currentCohort%c_area/currentPatch%total_canopy_area currentCohort%sai = currentCohort%treesai *currentCohort%c_area/currentPatch%total_canopy_area ! Number of actual vegetation layers in this cohort's crown - currentCohort%NV = ceiling((currentCohort%treelai+currentCohort%treesai)/dinc_ed) + currentCohort%nv = ceiling((currentCohort%treelai+currentCohort%treesai)/dinc_ed) currentPatch%ncan(cl,ft) = max(currentPatch%ncan(cl,ft),currentCohort%NV) patch_lai = patch_lai + currentCohort%lai -! currentPatch%canopy_layer_tai(cl) = currentPatch%canopy_layer_tai(cl) + & -! currentCohort%lai + currentCohort%sai + currentPatch%canopy_layer_tlai(cl) = currentPatch%canopy_layer_tlai(cl) + currentCohort%lai - do cl = 1,nclmax-1 - if(currentCohort%canopy_layer == cl)then - currentPatch%canopy_layer_tai(cl) = currentPatch%canopy_layer_tai(cl) + & - currentCohort%lai + currentCohort%sai - endif - enddo - - currentCohort => currentCohort%taller + currentCohort => currentCohort%shorter enddo !currentCohort @@ -1384,7 +1386,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) if(iv==currentCohort%NV) then remainder = (currentCohort%treelai + currentCohort%treesai) - & - (dinc_ed*dble(currentCohort%NV-1.0_r8)) + (dinc_ed*dble(currentCohort%nv-1.0_r8)) if(remainder > dinc_ed )then write(fates_log(), *)'ED: issue with remainder', & currentCohort%treelai,currentCohort%treesai,dinc_ed, & @@ -1484,7 +1486,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) enddo !currentCohort call endrun(msg=errMsg(sourcefile, __LINE__)) end if - end do + end do do ft = 1,numpft do iv = 1,currentPatch%ncan(cl,ft) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index e102495a4e..d57000e18f 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -21,7 +21,6 @@ module EDCohortDynamicsMod use EDTypesMod , only : min_npm2, min_nppatch use EDTypesMod , only : min_n_safemath use EDTypesMod , only : nlevleaf - use EDTypesMod , only : dinc_ed use FatesInterfaceMod , only : hlm_use_planthydro use FatesPlantHydraulicsMod, only : FuseCohortHydraulics use FatesPlantHydraulicsMod, only : CopyCohortHydraulics @@ -152,11 +151,16 @@ subroutine create_cohort(currentSite, patchptr, pft, nn, hite, dbh, bleaf, bfine ! Assign canopy extent and depth call carea_allom(new_cohort%dbh,new_cohort%n,spread,new_cohort%pft,new_cohort%c_area) - new_cohort%treelai = tree_lai(new_cohort%bl, new_cohort%status_coh, new_cohort%pft, & - new_cohort%c_area, new_cohort%n) + new_cohort%treelai = tree_lai(new_cohort%bl, new_cohort%pft, new_cohort%c_area, & + new_cohort%n, new_cohort%canopy_layer, & + patchptr%canopy_layer_tlai ) + + new_cohort%treesai = tree_sai(new_cohort%pft, new_cohort%dbh, new_cohort%canopy_trim, & + new_cohort%c_area, new_cohort%n, new_cohort%canopy_layer, & + patchptr%canopy_layer_tlai, new_cohort%treelai ) + new_cohort%lai = new_cohort%treelai * new_cohort%c_area/patchptr%area - new_cohort%treesai = tree_sai(new_cohort%dbh, new_cohort%pft, new_cohort%canopy_trim, & - new_cohort%c_area, new_cohort%n) + ! Put cohort at the right place in the linked list storebigcohort => patchptr%tallest diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index a120835d4c..358f342297 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -1267,7 +1267,7 @@ subroutine zero_patch(cp_p) currentPatch%age = nan currentPatch%age_class = 1 currentPatch%area = nan - currentPatch%canopy_layer_tai(:) = nan + currentPatch%canopy_layer_tlai(:) = nan currentPatch%total_canopy_area = nan currentPatch%tlai_profile(:,:,:) = nan @@ -1347,7 +1347,7 @@ subroutine zero_patch(cp_p) currentPatch%burnt_frac_litter(:) = 0.0_r8 currentPatch%btran_ft(:) = 0.0_r8 - currentPatch%canopy_layer_tai(:) = 0.0_r8 + currentPatch%canopy_layer_tlai(:) = 0.0_r8 currentPatch%seeds_in(:) = 0.0_r8 currentPatch%seed_decay(:) = 0.0_r8 diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index 0065031905..52d7fb2b56 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -20,6 +20,7 @@ module EDPhysiologyMod use EDCohortDynamicsMod , only : create_cohort, sort_cohorts use FatesAllometryMod , only : tree_lai use FatesAllometryMod , only : tree_sai + use FatesAllometryMod , only : decay_coeff_kn use EDTypesMod , only : numWaterMem use EDTypesMod , only : dl_sf, dinc_ed @@ -162,8 +163,7 @@ subroutine trim_canopy( currentSite ) ! Canopy trimming / leaf optimisation. Removes leaves in negative annual carbon balance. ! ! !USES: - ! - ! + ! !ARGUMENTS type (ed_site_type),intent(inout), target :: currentSite ! @@ -171,31 +171,51 @@ subroutine trim_canopy( currentSite ) type (ed_cohort_type) , pointer :: currentCohort type (ed_patch_type) , pointer :: currentPatch - integer :: z ! leaf layer - integer :: ipft ! pft index - integer :: trimmed ! was this layer trimmed in this year? If not expand the canopy. - real(r8) :: tar_bl ! target leaf biomass (leaves flushed, trimmed) - real(r8) :: tar_bfr ! target fine-root biomass (leaves flushed, trimmed) - real(r8) :: bfr_per_bleaf ! ratio of fine root per leaf biomass + integer :: z ! leaf layer + integer :: ipft ! pft index + logical :: trimmed ! was this layer trimmed in this year? If not expand the canopy. + real(r8) :: tar_bl ! target leaf biomass (leaves flushed, trimmed) + real(r8) :: tar_bfr ! target fine-root biomass (leaves flushed, trimmed) + real(r8) :: bfr_per_bleaf ! ratio of fine root per leaf biomass + real(r8) :: sla_levleaf ! sla at leaf level z + real(r8) :: nscaler_levleaf ! nscaler value at leaf level z + integer :: cl ! canopy layer index + real(r8) :: kn ! nitrogen decay coefficient + real(r8) :: sla_max ! Observational constraint on how large sla (m2/gC) can become + + real(r8) :: leaf_inc ! LAI-only portion of the vegetation increment of dinc_ed + real(r8) :: lai_canopy_above ! the LAI in the canopy layers above the layer of interest + real(r8) :: lai_layers_above ! the LAI in the leaf layers, within the current canopy, + ! above the leaf layer of interest + real(r8) :: lai_current ! the LAI in the current leaf layer + real(r8) :: cumulative_lai ! the cumulative LAI, top down, to the leaf layer of interest !---------------------------------------------------------------------- currentPatch => currentSite%youngest_patch - do while(associated(currentPatch)) + currentCohort => currentPatch%tallest do while (associated(currentCohort)) - trimmed = 0 + + trimmed = .false. ipft = currentCohort%pft call carea_allom(currentCohort%dbh,currentCohort%n,currentSite%spread,currentCohort%pft,currentCohort%c_area) - currentCohort%treelai = tree_lai(currentCohort%bl, currentCohort%status_coh, currentCohort%pft, & - currentCohort%c_area, currentCohort%n ) - currentCohort%treesai = tree_sai(currentCohort%dbh, currentCohort%pft, currentCohort%canopy_trim, & - currentCohort%c_area, currentCohort%n) - currentCohort%nv = ceiling((currentCohort%treelai+currentCohort%treesai)/dinc_ed) + currentCohort%treelai = tree_lai(currentCohort%bl, currentCohort%pft, currentCohort%c_area, & + currentCohort%n, currentCohort%canopy_layer, & + currentPatch%canopy_layer_tlai ) + + currentCohort%treesai = tree_sai(currentCohort%pft, currentCohort%dbh, currentCohort%canopy_trim, & + currentCohort%c_area, currentCohort%n, currentCohort%canopy_layer, & + currentPatch%canopy_layer_tlai, currentCohort%treelai ) + + currentCohort%nv = ceiling((currentCohort%treelai+currentCohort%treesai)/dinc_ed) + if (currentCohort%nv > nlevleaf)then - write(fates_log(),*) 'nv > nlevleaf',currentCohort%nv,currentCohort%treelai,currentCohort%treesai, & - currentCohort%c_area,currentCohort%n,currentCohort%bl + write(fates_log(),*) 'nv > nlevleaf',currentCohort%nv, & + currentCohort%treelai,currentCohort%treesai, & + currentCohort%c_area,currentCohort%n,currentCohort%bl + call endrun(msg=errMsg(sourcefile, __LINE__)) endif call bleaf(currentcohort%dbh,ipft,currentcohort%canopy_trim,tar_bl) @@ -206,35 +226,68 @@ subroutine trim_canopy( currentSite ) bfr_per_bleaf = tar_bfr/tar_bl endif + ! Identify current canopy layer (cl) + cl = currentCohort%canopy_layer + + ! PFT-level maximum SLA value, even if under a thick canopy (same units as slatop) + sla_max = EDPftvarcon_inst%slamax(ipft) + !Leaf cost vs netuptake for each leaf layer. - do z = 1,nlevleaf - if (currentCohort%year_net_uptake(z) /= 999._r8)then !there was activity this year in this leaf layer. + do z = 1, currentCohort%nv + + ! Calculate the cumulative total vegetation area index (no snow occlusion, stems and leaves) + + leaf_inc = dinc_ed * & + currentCohort%treelai/(currentCohort%treelai+currentCohort%treesai) + + ! Now calculate the cumulative top-down lai of the current layer's midpoint + lai_canopy_above = sum(currentPatch%canopy_layer_tlai(1:cl-1)) + lai_layers_above = leaf_inc * (z-1) + lai_current = min(leaf_inc, currentCohort%treelai - lai_layers_above) + cumulative_lai = lai_canopy_above + lai_layers_above + 0.5*lai_current + + if (currentCohort%year_net_uptake(z) /= 999._r8)then !there was activity this year in this leaf layer. + + + ! Calculate sla_levleaf following the sla profile with overlying leaf area + ! Scale for leaf nitrogen profile + kn = decay_coeff_kn(ipft) + ! Nscaler value at leaf level z + nscaler_levleaf = exp(-kn * cumulative_lai) + ! Sla value at leaf level z after nitrogen profile scaling (m2/gC) + sla_levleaf = EDPftvarcon_inst%slatop(ipft)/nscaler_levleaf + + if(sla_levleaf > sla_max)then + sla_levleaf = sla_max + end if + !Leaf Cost kgC/m2/year-1 !decidous costs. if (EDPftvarcon_inst%season_decid(ipft) == 1.or. & EDPftvarcon_inst%stress_decid(ipft) == 1)then - - currentCohort%leaf_cost = 1._r8/(EDPftvarcon_inst%slatop(ipft)*1000.0_r8) + ! Leaf cost at leaf level z accounting for sla profile (kgC/m2) + currentCohort%leaf_cost = 1._r8/(sla_levleaf*1000.0_r8) if ( int(EDPftvarcon_inst%allom_fmode(ipft)) .eq. 1 ) then ! if using trimmed leaf for fine root biomass allometry, add the cost of the root increment ! to the leaf increment; otherwise do not. currentCohort%leaf_cost = currentCohort%leaf_cost + & - 1.0_r8/(EDPftvarcon_inst%slatop(ipft)*1000.0_r8) * & + 1.0_r8/(sla_levleaf*1000.0_r8) * & bfr_per_bleaf / EDPftvarcon_inst%root_long(ipft) endif currentCohort%leaf_cost = currentCohort%leaf_cost * & (EDPftvarcon_inst%grperc(ipft) + 1._r8) else !evergreen costs - currentCohort%leaf_cost = 1.0_r8/(EDPftvarcon_inst%slatop(ipft)* & + ! Leaf cost at leaf level z accounting for sla profile + currentCohort%leaf_cost = 1.0_r8/(sla_levleaf* & EDPftvarcon_inst%leaf_long(ipft)*1000.0_r8) !convert from sla in m2g-1 to m2kg-1 if ( int(EDPftvarcon_inst%allom_fmode(ipft)) .eq. 1 ) then ! if using trimmed leaf for fine root biomass allometry, add the cost of the root increment ! to the leaf increment; otherwise do not. currentCohort%leaf_cost = currentCohort%leaf_cost + & - 1.0_r8/(EDPftvarcon_inst%slatop(ipft)*1000.0_r8) * & + 1.0_r8/(sla_levleaf*1000.0_r8) * & bfr_per_bleaf / EDPftvarcon_inst%root_long(ipft) endif currentCohort%leaf_cost = currentCohort%leaf_cost * & @@ -256,7 +309,7 @@ subroutine trim_canopy( currentSite ) currentCohort%laimemory = currentCohort%laimemory * & (1.0_r8 - EDPftvarcon_inst%trim_inc(ipft)) endif - trimmed = 1 + trimmed = .true. endif endif endif @@ -264,7 +317,7 @@ subroutine trim_canopy( currentSite ) enddo !z currentCohort%year_net_uptake(:) = 999.0_r8 - if (trimmed == 0.and.currentCohort%canopy_trim < 1.0_r8)then + if ( (.not.trimmed) .and.currentCohort%canopy_trim < 1.0_r8)then currentCohort%canopy_trim = currentCohort%canopy_trim + EDPftvarcon_inst%trim_inc(ipft) endif @@ -2416,4 +2469,9 @@ subroutine flux_into_litter_pools(nsites, sites, bc_in, bc_out) end subroutine flux_into_litter_pools + ! =================================================================================== + + + + end module EDPhysiologyMod diff --git a/biogeochem/FatesAllometryMod.F90 b/biogeochem/FatesAllometryMod.F90 index 432e1efb1b..34236ca93a 100644 --- a/biogeochem/FatesAllometryMod.F90 +++ b/biogeochem/FatesAllometryMod.F90 @@ -93,6 +93,8 @@ module FatesAllometryMod use FatesGlobals , only : fates_log use FatesGlobals , only : endrun => fates_endrun use EDTypesMod , only : nlevleaf, dinc_ed + use EDTypesMod , only : nclmax + implicit none @@ -104,13 +106,14 @@ module FatesAllometryMod public :: bleaf ! Generic actual leaf biomass wrapper public :: storage_fraction_of_target ! storage as fraction of leaf biomass public :: tree_lai ! Calculate tree-level LAI from actual leaf biomass - public :: tree_sai ! Calculate tree-level SAI from target leaf biomass + public :: tree_sai ! Calculate tree-level SAI from tree-level LAI public :: bsap_allom ! Generic sapwood wrapper public :: bbgw_allom ! Generic coarse root wrapper public :: bfineroot ! Generic actual fine root biomass wrapper public :: bdead_allom ! Generic bdead wrapper public :: carea_allom ! Generic crown area wrapper public :: bstore_allom ! Generic maximum storage carbon wrapper + public :: decay_coeff_kn public :: StructureResetOfDH ! Method to set DBH to sync with structure biomass public :: CheckIntegratedAllometries public :: set_root_fraction ! Generic wrapper to calculate normalized @@ -540,85 +543,169 @@ subroutine storage_fraction_of_target(b_leaf, bstore, frac) end subroutine storage_fraction_of_target ! ===================================================================================== - - real(r8) function tree_lai( bl, status_coh, pft, c_area, n ) - ! ============================================================================ - ! LAI of individual trees is a function of the total leaf area and the total canopy area. - ! ============================================================================ + real(r8) function tree_lai( bl, pft, c_area, nplant, cl, canopy_lai) + + ! ----------------------------------------------------------------------------------- + ! LAI of individual trees is a function of the total leaf area and the total + ! canopy area. + ! ---------------------------------------------------------------------------------- - real(r8), intent(in) :: bl ! plant leaf biomass [kg] - integer, intent(in) :: status_coh ! growth status of plant (2 = leaves on , 1 = leaves off) + ! !ARGUMENTS + real(r8), intent(in) :: bl ! plant leaf biomass [kg] integer, intent(in) :: pft - real(r8), intent(in) :: c_area ! areal extent of canopy (m2) - real(r8), intent(in) :: n ! number of individuals in cohort per 'area' (10000m2 default) + real(r8), intent(in) :: c_area ! areal extent of canopy (m2) + real(r8), intent(in) :: nplant ! number of individuals in cohort per ha + integer, intent(in) :: cl ! canopy layer index + real(r8), intent(in) :: canopy_lai(nclmax) ! total leaf area index of + ! each canopy layer + ! !LOCAL VARIABLES: real(r8) :: leafc_per_unitarea ! KgC of leaf per m2 area of ground. real(r8) :: slat ! the sla of the top leaf layer. m2/kgC + real(r8) :: canopy_lai_above ! total LAI of canopy layer overlying this tree + real(r8) :: vai_per_lai ! ratio of vegetation area index (ie. sai+lai) + ! to lai for individual tree + real(r8) :: kn ! coefficient for exponential decay of 1/sla and + ! vcmax with canopy depth + real(r8) :: sla_max ! Observational constraint on how large sla + ! (m2/gC) can become + real(r8) :: leafc_slamax ! Leafc_per_unitarea at which sla_max is reached + real(r8) :: clim ! Upper limit for leafc_per_unitarea in exponential + ! tree_lai function + !---------------------------------------------------------------------- if( bl < 0._r8 .or. pft == 0 ) then write(fates_log(),*) 'problem in treelai',bl,pft endif slat = g_per_kg * EDPftvarcon_inst%slatop(pft) ! m2/g to m2/kg - leafc_per_unitarea = bl/(c_area/n) !KgC/m2 + leafc_per_unitarea = bl/(c_area/nplant) !KgC/m2 + if(leafc_per_unitarea > 0.0_r8)then - tree_lai = leafc_per_unitarea * slat !kg/m2 * m2/kg = unitless LAI - else - tree_lai = 0.0_r8 - endif - ! here, if the LAI exceeeds the maximum size of the possible array, then we have no way of accomodating it - ! at the moments nlevleaf default is 40, which is very large, so exceeding this would clearly illustrate a - ! huge error - if(tree_lai > nlevleaf*dinc_ed)then - write(fates_log(),*) 'too much lai' , tree_lai , pft , nlevleaf * dinc_ed - write(fates_log(),*) 'Aborting' - call endrun(msg=errMsg(sourcefile, __LINE__)) - endif + if (cl==1) then ! if in we are in the canopy (top) layer) + canopy_lai_above = 0._r8 + else + canopy_lai_above = sum(canopy_lai(1:cl-1)) + end if - return + ! Coefficient for exponential decay of 1/sla with canopy depth: + kn = decay_coeff_kn(pft) + + ! take PFT-level maximum SLA value, even if under a thick canopy (which has units of m2/gC), + ! and put into units of m2/kgC + sla_max = g_per_kg *EDPftvarcon_inst%slamax(pft) + ! Leafc_per_unitarea at which sla_max is reached due to exponential sla profile in canopy: + leafc_slamax = (slat - sla_max * exp(-1.0_r8 * kn * canopy_lai_above)) / & + (-1.0_r8 * kn * slat * sla_max) + if(leafc_slamax < 0.0_r8)then + leafc_slamax = 0.0_r8 + endif + + ! Calculate tree_lai (m2 leaf area /m2 ground) = unitless LAI + !---------------------------------------------------------------------- + ! If leafc_per_unitarea is less than leafc_slamax, + ! sla with depth in the canopy will not exceed sla_max. + ! In this case, we can use an exponential profile for sla throughout the entire canopy. + ! The exponential profile for sla is given by: + ! sla(at a given canopy depth) = slat / exp(-kn (canopy_lai_above + tree_lai) + ! + ! We can solve for tree_lai using the above function for the sla profile and first setting + ! leafc_per_unitarea = integral of e^(-kn(x + canopy_lai_above)) / slatop + ! over x = 0 to tree_lai + ! Then, rearranging the equation to solve for tree_lai. + + if (leafc_per_unitarea <= leafc_slamax)then + tree_lai = (log(exp(-1.0_r8 * kn * canopy_lai_above) - & + kn * slat * leafc_per_unitarea) + & + (kn * canopy_lai_above)) / (-1.0_r8 * kn) + + ! If leafc_per_unitarea becomes too large, tree_lai becomes an imaginary number + ! (because the tree_lai equation requires us to take the natural log of something >0) + ! Thus, we include the following error message in case leafc_per_unitarea becomes too large. + clim = (exp(-1.0_r8 * kn * canopy_lai_above)) / (kn * slat) + if (leafc_per_unitarea >= clim) then + write(fates_log(),*) 'too much leafc_per_unitarea' , leafc_per_unitarea, clim, pft, canopy_lai_above + write(fates_log(),*) 'Aborting' + call endrun(msg=errMsg(sourcefile, __LINE__)) + endif + + ! When leafc_per_unitarea is greater than leafc_slamax, + ! tree_lai could become so great that the sla profile surpasses sla_max at depth. + ! In this case, we use the exponential profile to calculate tree_lai until + ! we reach the maximum allowed value for sla (sla_max). + ! Then, calculate the remaining tree_lai using a linear function of sla_max and the remaining leafc. + + else if(leafc_per_unitarea > leafc_slamax)then + + ! Add exponential and linear portions of tree_lai + ! Exponential term for leafc = leafc_slamax; + ! Linear term (static sla = sla_max) for portion of leafc > leafc_slamax + tree_lai = ((log(exp(-1.0_r8 * kn * canopy_lai_above) - & + kn * slat * leafc_slamax) + & + (kn * canopy_lai_above)) / (-1.0_r8 * kn)) + & + (leafc_per_unitarea - leafc_slamax) * sla_max + + ! if leafc_slamax becomes too large, tree_lai_exp becomes an imaginary number + ! (because the tree_lai equation requires us to take the natural log of something >0) + ! Thus, we include the following error message in case leafc_slamax becomes too large. + clim = (exp(-1.0_r8 * kn * canopy_lai_above)) / (kn * slat) + if(leafc_slamax >= clim)then + write(fates_log(),*) 'too much leafc_slamax' , & + leafc_per_unitarea, leafc_slamax, clim, pft, canopy_lai_above + write(fates_log(),*) 'Aborting' + call endrun(msg=errMsg(sourcefile, __LINE__)) + endif + end if ! (leafc_per_unitarea > leafc_slamax) + else + tree_lai = 0.0_r8 + endif ! (leafc_per_unitarea > 0.0_r8) + return end function tree_lai ! ============================================================================ - real(r8) function tree_sai( dbh, pft, canopy_trim, c_area, n ) + real(r8) function tree_sai( pft, dbh, canopy_trim, c_area, nplant, cl, canopy_lai, treelai ) ! ============================================================================ - ! SAI of individual trees is a function of the target leaf biomass + ! SAI of individual trees is a function of the LAI of individual trees ! ============================================================================ - real(r8),intent(in) :: dbh - integer, intent(in) :: pft - real(r8),intent(in) :: canopy_trim - real(r8), intent(in) :: c_area ! areal extent of canopy (m2) - real(r8), intent(in) :: n ! number of individuals in cohort per 'area' (10000m2 default) + integer, intent(in) :: pft + real(r8), intent(in) :: dbh + real(r8), intent(in) :: canopy_trim ! trimming function (0-1) + real(r8), intent(in) :: c_area ! crown area (m2) + real(r8), intent(in) :: nplant ! number of plants + integer, intent(in) :: cl ! canopy layer index + real(r8), intent(in) :: canopy_lai(nclmax) ! total leaf area index of + ! each canopy layer + real(r8), intent(in) :: treelai ! tree LAI for checking purposes only - real(r8) :: leafc_per_unitarea ! KgC of target leaf per m2 area of ground. - real(r8) :: sai_scaler - real(r8) :: b_leaf + real(r8) :: target_bleaf + real(r8) :: target_lai - sai_scaler = g_per_kg * EDPftvarcon_inst%allom_sai_scaler(pft) ! m2/g to m2/kg + call bleaf(dbh,pft,canopy_trim,target_bleaf) - call bleaf(dbh,pft,canopy_trim,b_leaf) + target_lai = tree_lai( target_bleaf, pft, c_area, nplant, cl, canopy_lai) - leafc_per_unitarea = b_leaf/(c_area/n) !KgC/m2 + tree_sai = EDPftvarcon_inst%allom_sai_scaler(pft) * target_lai - tree_sai = leafc_per_unitarea * sai_scaler !kg/m2 * m2/kg = unitless SAI - ! here, if the LAI exceeeds the maximum size of the possible array, then we have no way of accomodating it - ! at the moments nlevleaf default is 40, which is very large, so exceeding this would clearly illustrate a - ! huge error - if(tree_sai > nlevleaf*dinc_ed)then - write(fates_log(),*) 'too much sai' , tree_sai , pft , nlevleaf * dinc_ed - write(fates_log(),*) 'Aborting' + if( (treelai + tree_sai) > (nlevleaf*dinc_ed) )then + write(fates_log(),*) 'The leaf and stem are predicted for a cohort, maxed out the array size' + write(fates_log(),*) 'lai: ',treelai + write(fates_log(),*) 'sai: ',tree_sai + write(fates_log(),*) 'lai+sai: ',treelai+tree_sai + write(fates_log(),*) 'nlevleaf,dinc_ed,nlevleaf*dinc_ed :',nlevleaf,dinc_ed,nlevleaf*dinc_ed call endrun(msg=errMsg(sourcefile, __LINE__)) - endif + end if + - return + return end function tree_sai ! ============================================================================ @@ -1975,7 +2062,36 @@ subroutine jackson_beta_root_profile(root_fraction, ft, zi) return end subroutine jackson_beta_root_profile + ! ===================================================================================== + + real(r8) function decay_coeff_kn(pft) + + ! --------------------------------------------------------------------------------- + ! This function estimates the decay coefficient used to estimate vertical + ! attenuation of properties in the canopy. + ! + ! Decay coefficient (kn) is a function of vcmax25top for each pft. + ! + ! Currently, this decay is applied to vcmax attenuation, and SLA (optionally) + ! + ! --------------------------------------------------------------------------------- + + !ARGUMENTS + integer, intent(in) :: pft + + !LOCAL VARIABLES + ! ----------------------------------------------------------------------------------- + + ! Bonan et al (2011) JGR, 116, doi:10.1029/2010JG001593 used + ! kn = 0.11. Here, we derive kn from vcmax25 as in Lloyd et al + ! (2010) Biogeosciences, 7, 1833-1859 + + decay_coeff_kn = exp(0.00963_r8 * EDPftvarcon_inst%vcmax25top(pft) - 2.43_r8) + + return + end function decay_coeff_kn + ! ===================================================================================== subroutine StructureResetOfDH( bdead, ipft, canopy_trim, d, h ) diff --git a/biogeophys/EDSurfaceAlbedoMod.F90 b/biogeophys/EDSurfaceAlbedoMod.F90 index 95c37de148..ee525cd479 100644 --- a/biogeophys/EDSurfaceAlbedoMod.F90 +++ b/biogeophys/EDSurfaceAlbedoMod.F90 @@ -112,7 +112,7 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out ) integer, parameter :: diag_nlevleaf = min(nlevleaf,max_diag_nlevleaf) ! for diagnostics, write a small number of leaf layers real(r8) :: denom - real(r8) :: lai_reduction(2) + real(r8) :: lai_reduction(nclmax) integer :: fp,iv,s ! array indices integer :: ib ! waveband number diff --git a/biogeophys/FatesPlantRespPhotosynthMod.F90 b/biogeophys/FatesPlantRespPhotosynthMod.F90 index 7a0950689b..ebbbdc37e5 100644 --- a/biogeophys/FatesPlantRespPhotosynthMod.F90 +++ b/biogeophys/FatesPlantRespPhotosynthMod.F90 @@ -70,6 +70,7 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) use EDTypesMod , only : ed_cohort_type use EDTypesMod , only : ed_site_type use EDTypesMod , only : maxpft + use EDTypesMod , only : dinc_ed use FatesInterfaceMod , only : bc_in_type use FatesInterfaceMod , only : bc_out_type use EDCanopyStructureMod, only : calc_areaindex @@ -84,6 +85,7 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) use FatesAllometryMod, only : storage_fraction_of_target use FatesAllometryMod, only : set_root_fraction use FatesAllometryMod, only : i_hydro_rootprof_context + use FatesAllometryMod, only : decay_coeff_kn ! ARGUMENTS: ! ----------------------------------------------------------------------------------- @@ -148,8 +150,6 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) real(r8) :: ceair ! vapor pressure of air, constrained (Pa) real(r8) :: nscaler ! leaf nitrogen scaling coefficient real(r8) :: leaf_frac ! ratio of to leaf biomass to total alive biomass - real(r8) :: laican ! canopy sum of lai_z - real(r8) :: vai ! leaf and steam area in ths layer. real(r8) :: tcsoi ! Temperature response function for root respiration. real(r8) :: tcwood ! Temperature response function for wood @@ -169,12 +169,21 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) real(r8) :: r_stomata ! Mean stomatal resistance across all leaves in the patch [s/m] - real(r8) :: maintresp_reduction_factor ! factor by which to reduce maintenance respiration when storage pools are low + real(r8) :: maintresp_reduction_factor ! factor by which to reduce maintenance + ! respiration when storage pools are low real(r8) :: b_leaf ! leaf biomass kgC real(r8) :: frac ! storage pool as a fraction of target leaf biomass real(r8) :: check_elai ! This is a check on the effective LAI that is calculated ! over each cohort x layer. real(r8) :: cohort_eleaf_area ! This is the effective leaf area [m2] reported by each cohort + + real(r8) :: leaf_inc ! LAI-only portion of the vegetation increment of dinc_ed + real(r8) :: lai_canopy_above ! the LAI in the canopy layers above the layer of interest + real(r8) :: lai_layers_above ! the LAI in the leaf layers, within the current canopy, + ! above the leaf layer of interest + real(r8) :: lai_current ! the LAI in the current leaf layer + real(r8) :: cumulative_lai ! the cumulative LAI, top down, to the leaf layer of interest + ! ----------------------------------------------------------------------------------- ! Keeping these two definitions in case they need to be added later @@ -215,7 +224,8 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) ! projected area basis [m^2/gC] woody => EDPftvarcon_inst%woody , & ! Is vegetation woody or not? leafcn => EDPftvarcon_inst%leafcn , & ! leaf C:N (gC/gN) - frootcn => EDPftvarcon_inst%frootcn, & ! froot C:N (gc/gN) ! slope of BB relationship + frootcn => EDPftvarcon_inst%frootcn, & ! froot C:N (gc/gN) + woodcn => EDPftvarcon_inst%woodcn, & ! wood C:N (gc/gN) q10 => FatesSynchronizedParamsInst%Q10 ) bbbopt(0) = ED_val_bbopt_c4 @@ -242,9 +252,6 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) g_sb_leaves = 0._r8 check_elai = 0._r8 - !psncanopy_pa = 0._r8 - !lmrcanopy_pa = 0._r8 - ! Part II. Filter out patches ! Patch level filter flag for photosynthesis calculations ! has a short memory, flags: @@ -293,13 +300,9 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) ! Bonan et al (2011) JGR, 116, doi:10.1029/2010JG001593 used ! kn = 0.11. Here, derive kn from vcmax25 as in Lloyd et al ! (2010) Biogeosciences, 7, 1833-1859 - ! Remove daylength factor from vcmax25 so that kn is based on maximum vcmax25 - if (bc_in(s)%dayl_factor_pa(ifp) == 0._r8) then - kn(ft) = 0._r8 - else - kn(ft) = exp(0.00963_r8 * EDPftvarcon_inst%vcmax25top(ft) - 2.43_r8) - end if + kn(ft) = decay_coeff_kn(ft) + ! This is probably unnecessary and already calculated ! ALSO, THIS ROOTING PROFILE IS USED TO CALCULATE RESPIRATION @@ -353,14 +356,6 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) ! are there any leaves of this pft in this layer? if(currentPatch%canopy_mask(cl,ft) == 1)then - if(cl==1)then !are we in the top canopy layer or a shaded layer? - laican = 0._r8 - else - - laican = sum(currentPatch%canopy_layer_tai(1:cl-1)) - - end if - ! Loop over leaf-layers do iv = 1,currentCohort%nv @@ -377,31 +372,41 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) if ( .not.rate_mask_z(iv,ft,cl) .or. (hlm_use_planthydro.eq.itrue) ) then if (hlm_use_planthydro.eq.itrue) then -! write(fates_log(),*) 'hlm_use_planthydro' -! write(fates_log(),*) 'has been set to true. You have inadvertently' -! write(fates_log(),*) 'turned on a future feature that is not in the' -! write(fates_log(),*) 'FATES codeset yet. Please set this to' -! write(fates_log(),*) 'false and re-compile.' -! call endrun(msg=errMsg(sourcefile, __LINE__)) - - bbb = max (bbbopt(nint(c3psn(ft)))*currentCohort%co_hydr%btran(1), 1._r8) + + + bbb = max (bbbopt(nint(c3psn(ft)))*currentCohort%co_hydr%btran(1), 1._r8) btran_eff = currentCohort%co_hydr%btran(1) + + ! dinc_ed is the total vegetation area index of each "leaf" layer + ! we convert to the leaf only portion of the increment + ! ------------------------------------------------------ + leaf_inc = dinc_ed * & + currentCohort%treelai/(currentCohort%treelai+currentCohort%treesai) + + ! Now calculate the cumulative top-down lai of the current layer's midpoint + lai_canopy_above = sum(currentPatch%canopy_layer_tlai(1:cl-1)) + lai_layers_above = leaf_inc * (iv-1) + lai_current = min(leaf_inc, currentCohort%treelai - lai_layers_above) + cumulative_lai = lai_canopy_above + lai_layers_above + 0.5*lai_current + else + bbb = max (bbbopt(nint(c3psn(ft)))*currentPatch%btran_ft(ft), 1._r8) btran_eff = currentPatch%btran_ft(ft) + ! For consistency sake, we use total LAI here, and not exposed + ! if the plant is under-snow, it will be effectively dormant for + ! the purposes of nscaler + + cumulative_lai = sum(currentPatch%canopy_layer_tlai(1:cl-1)) + & + sum(currentPatch%tlai_profile(cl,ft,1:iv-1)) + & + 0.5*currentPatch%tlai_profile(cl,ft,iv) + + end if - - ! Vegetation area index - vai = (currentPatch%elai_profile(cl,ft,iv)+currentPatch%esai_profile(cl,ft,iv)) - if (iv == 1) then - laican = laican + 0.5_r8 * vai - else - laican = laican + 0.5_r8 * (currentPatch%elai_profile(cl,ft,iv-1)+ & - currentPatch%esai_profile(cl,ft,iv-1)+vai) - end if + ! Scale for leaf nitrogen profile - nscaler = exp(-kn(ft) * laican) + nscaler = exp(-kn(ft) * cumulative_lai) ! Part VII: Calculate dark respiration (leaf maintenance) for this layer call LeafLayerMaintenanceRespiration( param_derived%lmr25top(ft),& ! in @@ -529,9 +534,9 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) ! Units are in (kgN/plant) ! ------------------------------------------------------------------ live_stem_n = EDPftvarcon_inst%allom_agb_frac(currentCohort%pft) * currentCohort%bsw / & - frootcn(currentCohort%pft) + woodcn(currentCohort%pft) live_croot_n = (1.0_r8-EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)) * currentCohort%bsw / & - frootcn(currentCohort%pft) + woodcn(currentCohort%pft) froot_n = currentCohort%br / frootcn(currentCohort%pft) @@ -1105,7 +1110,6 @@ subroutine ScaleLeafLayerFluxToCohort(nv, & ! in currentCohort%nv ! ------------------------------------------------------------------------------------ use FatesConstantsMod, only : umolC_to_kgC - use EDTypesMod, only : dinc_ed ! Arguments integer, intent(in) :: nv ! number of active leaf layers diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index c1f98bf061..a5c70ca2a4 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -51,6 +51,7 @@ module EDPftvarcon real(r8), allocatable :: stress_decid(:) real(r8), allocatable :: season_decid(:) real(r8), allocatable :: evergreen(:) + real(r8), allocatable :: slamax(:) real(r8), allocatable :: slatop(:) real(r8), allocatable :: leaf_long(:) real(r8), allocatable :: roota_par(:) @@ -68,6 +69,7 @@ module EDPftvarcon real(r8), allocatable :: vcmax25top(:) real(r8), allocatable :: leafcn(:) real(r8), allocatable :: frootcn(:) + real(r8), allocatable :: woodcn(:) real(r8), allocatable :: smpso(:) real(r8), allocatable :: smpsc(:) real(r8), allocatable :: grperc(:) @@ -357,6 +359,10 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_leaf_slamax' + call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & + dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_leaf_slatop' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -422,6 +428,10 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_wood_cn_ratio' + call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & + dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_smpso' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -782,6 +792,10 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%evergreen) + name = 'fates_leaf_slamax' + call fates_params%RetreiveParameterAllocate(name=name, & + data=this%slamax) + name = 'fates_leaf_slatop' call fates_params%RetreiveParameterAllocate(name=name, & data=this%slatop) @@ -846,6 +860,10 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%frootcn) + name = 'fates_wood_cn_ratio' + call fates_params%RetreiveParameterAllocate(name=name, & + data=this%woodcn) + name = 'fates_smpso' call fates_params%RetreiveParameterAllocate(name=name, & data=this%smpso) @@ -1487,7 +1505,8 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'stress_decid = ',EDPftvarcon_inst%stress_decid write(fates_log(),fmt0) 'season_decid = ',EDPftvarcon_inst%season_decid write(fates_log(),fmt0) 'evergreen = ',EDPftvarcon_inst%evergreen - write(fates_log(),fmt0) 'slatop = ',EDPftvarcon_inst%slatop + write(fates_log(),fmt0) 'slamax = ',EDPftvarcon_inst%slamax + write(fates_log(),fmt0) 'slatop = ',EDPftvarcon_inst%slatop write(fates_log(),fmt0) 'leaf_long = ',EDPftvarcon_inst%leaf_long write(fates_log(),fmt0) 'roota_par = ',EDPftvarcon_inst%roota_par write(fates_log(),fmt0) 'rootb_par = ',EDPftvarcon_inst%rootb_par @@ -1503,6 +1522,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'vcmax25top = ',EDPftvarcon_inst%vcmax25top write(fates_log(),fmt0) 'leafcn = ',EDPftvarcon_inst%leafcn write(fates_log(),fmt0) 'frootcn = ',EDPftvarcon_inst%frootcn + write(fates_log(),fmt0) 'woodcn = ',EDPftvarcon_inst%woodcn write(fates_log(),fmt0) 'smpso = ',EDPftvarcon_inst%smpso write(fates_log(),fmt0) 'smpsc = ',EDPftvarcon_inst%smpsc write(fates_log(),fmt0) 'grperc = ',EDPftvarcon_inst%grperc diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index d9a2083482..7d4c2e2841 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -11,7 +11,7 @@ module EDTypesMod integer, parameter :: maxPatchesPerSite = 10 ! maximum number of patches to live on a site integer, parameter :: maxCohortsPerPatch = 160 ! maximum number of cohorts per patch - + integer, parameter :: nclmax = 2 ! Maximum number of canopy layers integer, parameter :: ican_upper = 1 ! Nominal index for the upper canopy integer, parameter :: ican_ustory = 2 ! Nominal index for understory in two-canopy system @@ -22,7 +22,6 @@ module EDTypesMod ! are used, but this helps allocate scratch ! space and output arrays. - ! ------------------------------------------------------------------------------------- ! Radiation parameters ! These should be part of the radiation module, but since we only have one option @@ -81,7 +80,8 @@ module EDTypesMod ! BIOLOGY/BIOGEOCHEMISTRY integer , parameter :: external_recruitment = 0 ! external recruitment flag 1=yes integer , parameter :: SENES = 10 ! Window of time over which we track temp for cold sensecence (days) - real(r8), parameter :: DINC_ED = 1.0_r8 ! size of LAI bins. + real(r8), parameter :: dinc_ed = 1.0_r8 ! size of VAI bins (LAI+SAI) [CHANGE THIS NAME WITH NEXT INTERFACE + ! UPDATE] integer , parameter :: N_DIST_TYPES = 3 ! Disturbance Modes 1) tree-fall, 2) fire, 3) logging integer , parameter :: dtype_ifall = 1 ! index for naturally occuring tree-fall generated event integer , parameter :: dtype_ifire = 2 ! index for fire generated disturbance event @@ -312,7 +312,7 @@ module EDTypesMod ! LEAF ORGANIZATION real(r8) :: pft_agb_profile(maxpft,n_dbh_bins) ! binned above ground biomass, for patch fusion: KgC/m2 - real(r8) :: canopy_layer_tai(nclmax) ! total area index of each canopy layer + real(r8) :: canopy_layer_tlai(nclmax) ! total leaf area index of each canopy layer ! used to determine attenuation of parameters during ! photosynthesis m2 veg / m2 of canopy area (patch without bare ground) real(r8) :: total_canopy_area ! area that is covered by vegetation : m2 diff --git a/parameter_files/fates_params_14pfts.cdl b/parameter_files/fates_params_14pfts.cdl index 6ac2539a26..c513356358 100644 --- a/parameter_files/fates_params_14pfts.cdl +++ b/parameter_files/fates_params_14pfts.cdl @@ -192,8 +192,8 @@ variables: fates_allom_lmode:units = "index" ; fates_allom_lmode:long_name = "leaf biomass allometry function index" ; float fates_allom_sai_scaler(fates_pft) ; - fates_allom_sai_scaler:units = "m2/gC" ; - fates_allom_sai_scaler:long_name = "allometric ratio of SAI to target bleaf" ; + fates_allom_sai_scaler:units = "m2/m2" ; + fates_allom_sai_scaler:long_name = "allometric ratio of SAI per LAI" ; float fates_allom_smode(fates_pft) ; fates_allom_smode:units = "index" ; fates_allom_smode:long_name = "sapwood allometry function index" ; @@ -308,6 +308,9 @@ variables: float fates_leaf_long(fates_pft) ; fates_leaf_long:units = "yr" ; fates_leaf_long:long_name = "Leaf longevity (ie turnover timescale)" ; + float fates_leaf_slamax(fates_pft) ; + fates_leaf_slamax:units = "m^2/gC" ; + fates_leaf_slamax:long_name = "Maximum Specific Leaf Area (SLA), even if under a dense canopy" ; float fates_leaf_slatop(fates_pft) ; fates_leaf_slatop:units = "m^2/gC" ; fates_leaf_slatop:long_name = "Specific Leaf Area (SLA) at top of canopy, projected area basis" ; @@ -470,6 +473,9 @@ variables: float fates_trim_limit(fates_pft) ; fates_trim_limit:units = "m2/m2" ; fates_trim_limit:long_name = "Arbitrary limit to reductions in leaf area with stress" ; + float fates_wood_cn_ratio(fates_pft) ; + fates_wood_cn_ratio:units = "gC/gN" ; + fates_wood_cn_ratio:long_name = "Wood C:N" ; float fates_wood_density(fates_pft) ; fates_wood_density:units = "g/cm3" ; fates_wood_density:long_name = "mean density of woody tissue in plant" ; @@ -708,8 +714,8 @@ data: fates_allom_lmode = 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ; - fates_allom_sai_scaler = 0.0012, 0.0012, 0.0012, 0.0012, 0.0012, 0.0012, - 0.0012, 0.0012, 0.0012, 0.0012, 0.0012, 0.0012, 0.0012, 0.0012 ; + fates_allom_sai_scaler = 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, + 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1 ; fates_allom_smode = 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ; @@ -867,6 +873,9 @@ data: fates_leaf_long = 1.5, 4, 6, 1, 1.5, 1, 1, 1, 1.5, 1, 1, 1, 1, 1 ; + fates_leaf_slamax = 0.0954, 0.0954, 0.0954, 0.0954, 0.0954, 0.0954, 0.0954, 0.0954, + 0.012, 0.03, 0.03, 0.03, 0.03, 0.03 ; + fates_leaf_slatop = 0.012, 0.01, 0.008, 0.024, 0.012, 0.03, 0.03, 0.03, 0.012, 0.03, 0.03, 0.03, 0.03, 0.03 ; @@ -1022,6 +1031,9 @@ data: fates_trim_limit = 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3, 0.3 ; + fates_wood_cn_ratio = 210, 210, 210, 210, 210, 210, 210, 210, 210, 210, 210, + 210, 210, 210 ; + fates_wood_density = 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7 ; diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index 3d30509506..8f58d4cefa 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -159,8 +159,8 @@ variables: fates_allom_lmode:units = "index" ; fates_allom_lmode:long_name = "leaf biomass allometry function index" ; float fates_allom_sai_scaler(fates_pft) ; - fates_allom_sai_scaler:units = "m2/gC" ; - fates_allom_sai_scaler:long_name = "allometric ratio of SAI to target bleaf" ; + fates_allom_sai_scaler:units = "m2/m2" ; + fates_allom_sai_scaler:long_name = "allometric ratio of SAI per LAI" ; float fates_allom_smode(fates_pft) ; fates_allom_smode:units = "index" ; fates_allom_smode:long_name = "sapwood allometry function index" ; @@ -248,6 +248,9 @@ variables: float fates_leaf_long(fates_pft) ; fates_leaf_long:units = "yr" ; fates_leaf_long:long_name = "Leaf longevity (ie turnover timescale)" ; + float fates_leaf_slamax(fates_pft) ; + fates_leaf_slamax:units = "m^2/gC" ; + fates_leaf_slamax:long_name = "Maximum Specific Leaf Area (SLA), even if under a dense canopy" ; float fates_leaf_slatop(fates_pft) ; fates_leaf_slatop:units = "m^2/gC" ; fates_leaf_slatop:long_name = "Specific Leaf Area (SLA) at top of canopy, projected area basis" ; @@ -407,6 +410,9 @@ variables: float fates_trim_limit(fates_pft) ; fates_trim_limit:units = "m2/m2" ; fates_trim_limit:long_name = "Arbitrary limit to reductions in leaf area with stress" ; + float fates_wood_cn_ratio(fates_pft) ; + fates_wood_cn_ratio:units = "gC/gN" ; + fates_wood_cn_ratio:long_name = "Wood C:N" ; float fates_wood_density(fates_pft) ; fates_wood_density:units = "g/cm3" ; fates_wood_density:long_name = "mean density of woody tissue in plant" ; @@ -695,7 +701,7 @@ data: fates_allom_lmode = 1, 1 ; - fates_allom_sai_scaler = 0.0012, 0.0012 ; + fates_allom_sai_scaler = 0.10, 0.10 ; fates_allom_smode = 1, 1 ; @@ -755,6 +761,8 @@ data: fates_leaf_long = 1.5, 1.5 ; + fates_leaf_slamax = 0.0954, 0.0954 ; + fates_leaf_slatop = 0.012, 0.012 ; fates_leaf_stor_priority = 0.8, 0.8 ; @@ -861,6 +869,8 @@ data: fates_trim_limit = 0.3, 0.3 ; + fates_wood_cn_ratio = 210, 210 ; + fates_wood_density = 0.7, 0.7 ; fates_woody = 1, 1 ;