diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index c4ccb4f7e6..3a58ad19ad 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -186,6 +186,7 @@ subroutine trim_canopy( currentSite ) real(r8) :: vai ! leaf and stem area in this layer real(r8) :: kn ! nitrogen decay coefficient real(r8) :: sla_max ! Observational constraint on how large sla (m2/gC) can become + real(r8) :: tai_to_lai ! ratio of total area index (ie. sai + lai) to lai for individual tree !---------------------------------------------------------------------- @@ -226,19 +227,27 @@ subroutine trim_canopy( currentSite ) ! Observational constraint for maximum sla value (m2/gC): ! m2/gC = m2/gBiomass *kgC/kgBiomass - sla_max = sla_max_drymass * EDPftvarcon_inst%c2b(ipft) + sla_max = sla_max_drymass * EDPftvarcon_inst%c2b(ipft) + ! Ratio of total area index (ie. lai + sai) to lai for individual tree: + tai_to_lai = 1.0_r8 + EDPftvarcon_inst%allom_sai_scaler(ipft) !Leaf cost vs netuptake for each leaf layer. do z = 1,nlevleaf ! Vegetation area index - vai = (currentPatch%elai_profile(cl,ipft,z)+currentPatch%esai_profile(cl,ipft,z)) - if (z == 1) then - laican = laican + 0.5_r8 * vai - else - laican = laican + 0.5_r8 * (currentPatch%elai_profile(cl,ipft,z-1)+ & - currentPatch%esai_profile(cl,ipft,z-1)+vai) - end if + vai = (currentPatch%elai_profile(cl,ipft,z)+currentPatch%esai_profile(cl,ipft,z)) + if (vai > 0.0_r8) then + ! If there is vai in this leaf layer, add to laican + if (z == 1) then + ! If in first layer, add laican equal to laican value halfway through the layer + laican = laican + 0.5_r8 * tai_to_lai * dinc_ed + else + ! Since starting from halfway through first layer when z=1, + ! can add an entire dinc_ed to remain at halfway value for subsequent layers + laican = laican + tai_to_lai * dinc_ed + end if + else ! If vai = 0, do not add to laican + end if if (currentCohort%year_net_uptake(z) /= 999._r8)then !there was activity this year in this leaf layer.