Skip to content

Commit

Permalink
[Size- and age-dependent mortality functionality ]
Browse files Browse the repository at this point in the history
[These changes introduce two new mortality terms: M9 -
size dependent and M10 age dependent mortality. They are
logistic functions that gradually increase the rate of
mortality with size or age. Each is controlled by two
parameters a rate and an inflection point. To turn off
size dependent mortality the inflection point can be set
to a very large number (>10000). To turn off age mortality
set the fates_cohort_age_fusion_tol parameter to 0. Cohort
age will simply stay at 0. Otherwise, cohort age is tracked.
Cohorts cannot be fused if they are too different in age. ]

Fixes: [NGT-ED Github issue #]

User interface changes?: [Yes - five new parameters. User has to
specify age and size mortality as either on or off by setting
fates_mort_ip_size_senescence and fates_cohort_age_fusion_tol. If
they decide to have size or age on, then sensible parameter values
need to be selected for inflection point and rate.]

Code review: [Names]

Test suite: [suite name, machine, compilers]
Test baseline:
Test namelist changes:
Test answer changes: [bit for bit, roundoff, climate changing]
Test summary: No testing yet.
  • Loading branch information
Jessica Needham committed Dec 17, 2019
1 parent 5e81d7c commit ed84697
Show file tree
Hide file tree
Showing 21 changed files with 161 additions and 8,506 deletions.
13 changes: 9 additions & 4 deletions biogeochem/EDCanopyStructureMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,8 @@ module EDCanopyStructureMod
use FatesInterfaceMod , only : numpft
use FatesPlantHydraulicsMod, only : UpdateH2OVeg,InitHydrCohort, RecruitWaterStorage
use EDTypesMod , only : maxCohortsPerPatch

use EDParamsMod , only : ED_val_cohort_age_fusion_tol

use PRTGenericMod, only : leaf_organ
use PRTGenericMod, only : all_carbon_elements
use PRTGenericMod, only : leaf_organ
Expand Down Expand Up @@ -1261,7 +1262,8 @@ subroutine canopy_summarization( nsites, sites, bc_in )
use FatesSizeAgeTypeIndicesMod, only : coagetype_class_index
use EDtypesMod , only : area
use EDPftvarcon , only : EDPftvarcon_inst

use EDParamsMod , only : ED_val_cohort_age_fusion_tol

! !ARGUMENTS
integer , intent(in) :: nsites
type(ed_site_type) , intent(inout), target :: sites(nsites)
Expand All @@ -1280,8 +1282,9 @@ subroutine canopy_summarization( nsites, sites, bc_in )
real(r8) :: sapw_c ! sapwood carbon [kg]
real(r8) :: store_c ! storage carbon [kg]
real(r8) :: struct_c ! structure carbon [kg]
!----------------------------------------------------------------------

!----------------------------------------------------------------------

if ( debug ) then
write(fates_log(),*) 'in canopy_summarization'
endif
Expand Down Expand Up @@ -1322,9 +1325,11 @@ subroutine canopy_summarization( nsites, sites, bc_in )
call sizetype_class_index(currentCohort%dbh,currentCohort%pft, &
currentCohort%size_class,currentCohort%size_by_pft_class)

if (ED_val_cohort_age_fusion_tol > 0.0_r8) then
call coagetype_class_index(currentCohort%coage,currentCohort%pft, &
currentCohort%coage_class,currentCohort%coage_by_pft_class)

end if

call carea_allom(currentCohort%dbh,currentCohort%n,sites(s)%spread,&
currentCohort%pft,currentCohort%c_area)

Expand Down
39 changes: 24 additions & 15 deletions biogeochem/EDCohortDynamicsMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@ module EDCohortDynamicsMod
use EDTypesMod , only : ican_upper
use EDTypesMod , only : site_fluxdiags_type
use EDTypesMod , only : num_elements
use EDParamsMod , only : ED_val_cohort_age_fusion_tol
use FatesInterfaceMod , only : hlm_use_planthydro
use FatesInterfaceMod , only : hlm_parteh_mode
use FatesPlantHydraulicsMod, only : FuseCohortHydraulics
Expand Down Expand Up @@ -129,7 +130,7 @@ module EDCohortDynamicsMod
!-------------------------------------------------------------------------------------!


subroutine create_cohort(currentSite, patchptr, pft, nn, hite,coage, dbh, &
subroutine create_cohort(currentSite, patchptr, pft, nn, hite, coage, dbh, &
prt, laimemory, status, recruitstatus,ctrim, &
clayer, spread, bc_in)
!
Expand Down Expand Up @@ -183,7 +184,7 @@ subroutine create_cohort(currentSite, patchptr, pft, nn, hite,coage, dbh, &
integer :: nlevsoi_hyd ! number of hydraulically active soil layers

!----------------------------------------------------------------------

allocate(new_cohort)

call nan_cohort(new_cohort) ! Make everything in the cohort not-a-number
Expand Down Expand Up @@ -225,9 +226,13 @@ subroutine create_cohort(currentSite, patchptr, pft, nn, hite,coage, dbh, &
call sizetype_class_index(new_cohort%dbh, new_cohort%pft, &
new_cohort%size_class,new_cohort%size_by_pft_class)

! If cohort age trackign is off we call this here once
! just so everythin is in the first bin -
! this makes it easier to copy and terminate cohorts later
! we don't need to update this ever if cohort age tracking is off
call coagetype_class_index(new_cohort%coage, new_cohort%pft, &
new_cohort%coage_class,new_cohort%coage_by_pft_class)

new_cohort%coage_class,new_cohort%coage_by_pft_class)
! This routine may be called during restarts, and at this point in the call sequence
! the actual cohort data is unknown, as this is really only used for allocation
! In these cases, testing if things like biomass are reasonable is pre-mature
Expand Down Expand Up @@ -931,7 +936,7 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
! Join similar cohorts to reduce total number
!
! !USES:
use EDParamsMod , only : ED_val_cohort_fusion_tol
use EDParamsMod , only : ED_val_cohort_size_fusion_tol
use EDParamsMod , only : ED_val_cohort_age_fusion_tol
use FatesConstantsMod, only : days_per_year

Expand Down Expand Up @@ -960,7 +965,7 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
real(r8) :: leaf_c_next ! Leaf carbon * plant density of current (for weighting)
real(r8) :: leaf_c_curr ! Leaf carbon * plant density of next (for weighting)
real(r8) :: leaf_c_target
real(r8) :: dynamic_fusion_tolerance
real(r8) :: dynamic_size_fusion_tolerance
real(r8) :: dynamic_age_fusion_tolerance
real(r8) :: dbh
real(r8) :: leaf_c ! leaf carbon [kg]
Expand All @@ -976,7 +981,7 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
!----------------------------------------------------------------------

!set initial fusion tolerance
dynamic_fusion_tolerance = ED_val_cohort_fusion_tol
dynamic_size_fusion_tolerance = ED_val_cohort_size_fusion_tol
! set the cohort age fusion tolerance (this is in days - remains constant)
dynamic_age_fusion_tolerance = ED_val_cohort_age_fusion_tol

Expand Down Expand Up @@ -1013,18 +1018,20 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)

!Criteria used to divide up the height continuum into different cohorts.

if (diff < dynamic_fusion_tolerance) then
if (diff < dynamic_size_fusion_tolerance) then

! Only fuse if the cohorts are within x years of each other
! if they are the same age we make diff - to avoid errors divding by zero
! if they are the same age we make diff 0- to avoid errors divding by zero
!NB if cohort age tracking is off then the age of both should be 0
! and hence the age fusion criterion is met
if (currentCohort%coage .eq. nextc%coage) then
coage_diff = 0.0_r8
else
coage_diff = abs((currentCohort%coage - nextc%coage)/ &
(0.5*(currentCohort%coage + nextc%coage)))
end if

if (coage_diff < dynamic_age_fusion_tolerance ) then
if (coage_diff <= dynamic_age_fusion_tolerance ) then

! Don't fuse a cohort with itself!
if (.not.associated(currentCohort,nextc) ) then
Expand Down Expand Up @@ -1073,9 +1080,11 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
(nextc%coage * (nextc%n/(currentCohort%n + nextc%n)))

! update the cohort age again
call coagetype_class_index(currentCohort%coage, currentCohort%pft, &
currentCohort%coage_class, currentCohort%coage_by_pft_class)

if (ED_val_cohort_age_fusion_tol > 0.0_r8) then
call coagetype_class_index(currentCohort%coage, currentCohort%pft, &
currentCohort%coage_class, currentCohort%coage_by_pft_class)
end if

! Fuse all mass pools
call currentCohort%prt%WeightedFusePRTVartypes(nextc%prt, &
currentCohort%n/newn )
Expand Down Expand Up @@ -1395,7 +1404,7 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
!---------------------------------------------------------------------!
! Making profile tolerance larger means that more fusion will happen !
!---------------------------------------------------------------------!
dynamic_fusion_tolerance = dynamic_fusion_tolerance * 1.1_r8
dynamic_size_fusion_tolerance = dynamic_size_fusion_tolerance * 1.1_r8
dynamic_age_fusion_tolerance = dynamic_age_fusion_tolerance * 1.1_r8
!write(fates_log(),*) 'maxcohorts exceeded',dynamic_fusion_tolerance

Expand All @@ -1404,7 +1413,7 @@ subroutine fuse_cohorts(currentSite, currentPatch, bc_in)
iterate = 0
endif

if ( dynamic_fusion_tolerance .gt. 100._r8) then
if ( dynamic_size_fusion_tolerance .gt. 100._r8) then
! something has gone terribly wrong and we need to report what
write(fates_log(),*) 'exceeded reasonable expectation of cohort fusion.'
currentCohort => currentPatch%tallest
Expand Down
55 changes: 35 additions & 20 deletions biogeochem/EDMortalityFunctionsMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ module EDMortalityFunctionsMod
use FatesInterfaceMod , only : hlm_use_planthydro
use EDLoggingMortalityMod , only : LoggingMortality_frac
use EDParamsMod , only : fates_mortality_disturbance_fraction
use EDParamsMod , only : ED_val_cohort_age_fusion_tol
use FatesInterfaceMod , only : bc_in_type

use PRTGenericMod, only : all_carbon_elements
Expand Down Expand Up @@ -66,8 +67,8 @@ subroutine mortality_rates( cohort_in,bc_in,cmort,hmort,bmort,frmort,smort,asmor
real(r8) :: store_c
real(r8) :: hf_sm_threshold ! hydraulic failure soil moisture threshold
real(r8) :: hf_flc_threshold ! hydraulic failure fractional loss of conductivity threshold
real(r8) :: mort_ip_senescence ! inflection point for increase in mortality with dbh
real(r8) :: mort_r_senescence ! rate of mortality increase with dbh in senesence term
real(r8) :: mort_ip_size_senescence ! inflection point for increase in mortality with dbh
real(r8) :: mort_r_size_senescence ! rate of mortality increase with dbh in senesence term
real(r8) :: mort_ip_age_senescence ! inflection point for increase in mortality with age
real(r8) :: mort_r_age_senescence ! rate of mortality increase with age in senescence term
real(r8) :: temp_dep_fraction ! Temp. function (freezing mortality)
Expand All @@ -81,22 +82,36 @@ subroutine mortality_rates( cohort_in,bc_in,cmort,hmort,bmort,frmort,smort,asmor
logical, parameter :: test_zero_mortality = .false. ! Developer test which
! may help to debug carbon imbalances
! and the like

! Size Dependent Senescence
! rate (r) and inflection point (ip) define the increase in mortality rate with dbh
mort_r_senescence = EDPftvarcon_inst%mort_r_senescence(cohort_in%pft)
mort_ip_senescence = EDPftvarcon_inst%mort_ip_senescence(cohort_in%pft)
! smort = 1.0_r8 / ( 1.0_r8 + exp( -1.0_r8 * mort_r_senescence * &
! (cohort_in%dbh - mort_ip_senescence) ) )
smort = 0.0_r8

! Age Dependent Senescence
! rate and inflection point define the change in mortality with age
mort_r_age_senescence = EDPftvarcon_inst%mort_r_age_senescence(cohort_in%pft)
mort_ip_age_senescence = EDPftvarcon_inst%mort_ip_age_senescence(cohort_in%pft)
asmort = 1.0_r8 / (1.0_r8 + exp(-1.0_r8 * mort_r_age_senescence * &
(cohort_in%coage - mort_ip_age_senescence ) ) )
mort_r_size_senescence = EDPftvarcon_inst%mort_r_size_senescence(cohort_in%pft)
mort_ip_size_senescence = EDPftvarcon_inst%mort_ip_size_senescence(cohort_in%pft)

! if the user has set the inflection point to be a huge number then
! we zero smort (even though it would essentially be zero anyway)
if (mort_ip_size_senescence < 10000 ) then
smort = 1.0_r8 / ( 1.0_r8 + exp( -1.0_r8 * mort_r_size_senescence * &
(cohort_in%dbh - mort_ip_size_senescence) ) )
else
smort = 0.0_r8
end if

! if user has set cohort age fusion param > 0 we calculate age
! dependent mortality
if (ED_val_cohort_age_fusion_tol > 0.0_r8) then
! Age Dependent Senescence
! rate and inflection point define the change in mortality with age
mort_r_age_senescence = EDPftvarcon_inst%mort_r_age_senescence(cohort_in%pft)
mort_ip_age_senescence = EDPftvarcon_inst%mort_ip_age_senescence(cohort_in%pft)
asmort = 1.0_r8 / (1.0_r8 + exp(-1.0_r8 * mort_r_age_senescence * &
(cohort_in%coage - mort_ip_age_senescence ) ) )
else
asmort = 0.0_r8
end if



if (hlm_use_ed_prescribed_phys .eq. ifalse) then

! 'Background' mortality (can vary as a function of
Expand Down Expand Up @@ -237,10 +252,10 @@ subroutine Mortality_Derivative( currentSite, currentCohort, bc_in)
currentCohort%lmort_collateral + &
currentCohort%lmort_infra)/hlm_freq_day

! this check caps asmort so that daily mortality cannot exceed 0.995
! this caps bmort so that daily mortality cannot exceed 0.995
if ((cmort+hmort+bmort+frmort+smort+asmort + dndt_logging)*hlm_freq_day &
> 0.995_r8)then
asmort = (0.995_r8 - ((cmort+hmort+bmort+frmort+smort+dndt_logging)*hlm_freq_day)) &
bmort = (0.995_r8 - ((cmort+hmort+asmort+frmort+smort+dndt_logging)*hlm_freq_day)) &
/hlm_freq_day
endif

Expand All @@ -249,12 +264,12 @@ subroutine Mortality_Derivative( currentSite, currentCohort, bc_in)
* currentCohort%n
else

! cap on smort for canopy layers - smort is adjusted so that total mortality
! cap on bmort for canopy layers - bmort is adjusted so that total mortality
! including disturbance fraction does not exceed 0.995
if(((cmort+hmort+bmort+frmort+smort+asmort)*hlm_freq_day) > &
0.995_r8 - fates_mortality_disturbance_fraction)then
asmort = (0.995_r8 - fates_mortality_disturbance_fraction - &
((cmort+hmort+bmort+frmort+smort)*hlm_freq_day))/hlm_freq_day
bmort = (0.995_r8 - fates_mortality_disturbance_fraction - &
((cmort+hmort+asmort+frmort+smort)*hlm_freq_day))/hlm_freq_day
endif


Expand Down
Loading

0 comments on commit ed84697

Please sign in to comment.