Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add Atkin et al 2017 leaf respiration model #931

Merged
merged 34 commits into from
Mar 7, 2023
Merged

Add Atkin et al 2017 leaf respiration model #931

Show file tree
Hide file tree
Changes from 4 commits
Commits
Show all changes
34 commits
Select commit Hold shift + click to select a range
b907af7
starting to add hooks for multiple leaf maintenance respiration models
ckoven Nov 1, 2022
4e0d958
added Atkin et al 2017 respiration model
ckoven Nov 2, 2022
08d755c
bugfixes
ckoven Nov 3, 2022
bd091b0
fixed some Celsius vs Kelvin unit errors
ckoven Nov 3, 2022
5b372c9
changed vertical profile of leaf MR to be propotional to SLA and leaf…
ckoven Nov 15, 2022
a16e1b5
fixed bug in C4 respiration in prior version
ckoven Nov 15, 2022
1d90f6d
merged main and resolved conflicts
ckoven Dec 2, 2022
1f7d9af
merged main and resolved conflicts
ckoven Dec 2, 2022
8e582a6
merged two respiration branches
ckoven Dec 2, 2022
57572d3
added history diagnostic for running-mean Tgrowth variable
ckoven Jan 12, 2023
9fcdfa8
added history diagnostic for unreduced maintenance respiration
ckoven Jan 12, 2023
27d0c68
make Tgrowth and maintenance respiration reduction diagnostics inacti…
ckoven Jan 30, 2023
46fa764
reorganized C3/C4 logic to use separate subroutines
ckoven Jan 30, 2023
f5f1c9a
revert dark respiration model to the Ryan 1991 model
ckoven Jan 30, 2023
10db8b9
swapping hard-coded basal respiration rate in Ryan 1991 model with va…
ckoven Jan 30, 2023
27c5e78
removing PFT argument to Ryan 1991 respiration models as no longer ne…
ckoven Jan 31, 2023
22212b5
merged main and resolved conflicts
ckoven Feb 1, 2023
c2a1dd8
add xml patch file for update to api25.1.1
glemieux Feb 16, 2023
5c3fe60
fix missing third decimal place value in fates_base_mr_20
glemieux Feb 16, 2023
56b3266
Update xml file to old format to add new parameter value
glemieux Feb 22, 2023
10d686a
Update xml patch file for api25.1.1
glemieux Feb 23, 2023
8d680fa
Fix BatchPatchParams.py global history update call
glemieux Feb 24, 2023
c26f5c7
Use new Atkin parameter variable
glemieux Feb 24, 2023
6c00d7f
fix typo and comment out pft param registration
glemieux Feb 24, 2023
7e56708
Uncomment the atkin parameter registration and retrieval
glemieux Feb 27, 2023
8620f9e
Update maintenance respiration variable names
glemieux Feb 27, 2023
aa0375d
fix default param file name
glemieux Feb 28, 2023
59a865f
fix leaf model name assignment
glemieux Feb 28, 2023
5c97f6a
reverted C3/C4 logic to now use same model for all PFTs and same Atki…
ckoven Mar 1, 2023
fd9f5f6
Merge pull request #14 from glemieux/respiration-agnostic_vars
ckoven Mar 1, 2023
99ce45d
merged updates to parameter variable names and dimensions
ckoven Mar 1, 2023
ac8d58f
fixed bug from merge
ckoven Mar 1, 2023
a20b462
Rename the patch file api number
glemieux Mar 7, 2023
008f8f8
Minor instruction update in patch file
glemieux Mar 7, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
136 changes: 113 additions & 23 deletions biogeophys/FatesPlantRespPhotosynthMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,9 @@ module FATESPlantRespPhotosynthMod
use PRTGenericMod, only : max_nleafage
use EDTypesMod, only : do_fates_salinity
use EDParamsMod, only : q10_mr
use EDParamsMod, only : maintresp_model
use FatesConstantsMod, only : lmrmodel_ryan_1991
use FatesConstantsMod, only : lmrmodel_atkin_etal_2017
use PRTGenericMod, only : prt_carbon_allom_hyp
use PRTGenericMod, only : prt_cnp_flex_allom_hyp
use PRTGenericMod, only : all_carbon_elements
Expand Down Expand Up @@ -125,7 +128,6 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)
use FatesInterfaceTypesMod , only : bc_out_type
use EDCanopyStructureMod, only : calc_areaindex
use FatesConstantsMod, only : umolC_to_kgC
use FatesConstantsMod, only : g_per_kg
use FatesConstantsMod, only : umol_per_mmol
use FatesConstantsMod, only : rgas => rgas_J_K_kmol
use FatesConstantsMod, only : tfrz => t_water_freeze_k_1atm
Expand Down Expand Up @@ -269,12 +271,6 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)
integer :: iage ! loop counter for leaf age classes

! Parameters
! -----------------------------------------------------------------------
! Base maintenance respiration rate for plant tissues base_mr_20
! M. Ryan, 1991. Effects of climate change on plant respiration.
! Ecological Applications, 1(2), 157-167.
! Original expression is br = 0.0106 molC/(molN h)
! Conversion by molecular weights of C and N gives 2.525e-6 gC/(gN s)
!
! Base rate is at 20C. Adjust to 25C using the CN Q10 = 1.5
! (gC/gN/s)
Expand Down Expand Up @@ -520,17 +516,32 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)

end select

! MLO - Shouldn't these numbers be parameters too?
lmr25top = 2.525e-6_r8 * (1.5_r8 ** ((25._r8 - 20._r8)/10._r8))
lmr25top = lmr25top * lnc_top / (umolC_to_kgC * g_per_kg)
! Part VII: Calculate dark respiration (leaf maintenance) for this layer
select case (maintresp_model)

case (lmrmodel_ryan_1991)

! Part VII: Calculate dark respiration (leaf maintenance) for this layer
call LeafLayerMaintenanceRespiration( lmr25top, & ! in
nscaler, & ! in
ft, & ! in
bc_in(s)%t_veg_pa(ifp), & ! in
lmr_z(iv,ft,cl)) ! out
call LeafLayerMaintenanceRespiration_Ryan_1991( lnc_top, & ! in
nscaler, & ! in
ft, & ! in
bc_in(s)%t_veg_pa(ifp), & ! in
lmr_z(iv,ft,cl)) ! out

case (lmrmodel_atkin_etal_2017)

call LeafLayerMaintenanceRespiration_Atkin_etal_2017(lnc_top, & ! in
nscaler, & ! in
ft, & ! in
bc_in(s)%t_veg_pa(ifp), & ! in
currentPatch%tveg_lpa%GetMean(), & ! in
lmr_z(iv,ft,cl)) ! out

case default

write (fates_log(),*)'error, incorrect leaf respiration model specified'
call endrun(msg=errMsg(sourcefile, __LINE__))

end select

! Part VII: Calculate (1) maximum rate of carboxylation (vcmax),
! (2) maximum electron transport rate, (3) triose phosphate
Expand Down Expand Up @@ -1966,25 +1977,35 @@ end subroutine GetCanopyGasParameters

! ====================================================================================

subroutine LeafLayerMaintenanceRespiration(lmr25top_ft, &
subroutine LeafLayerMaintenanceRespiration_Ryan_1991(lnc_top, &
nscaler, &
ft, &
veg_tempk, &
lmr)

use FatesConstantsMod, only : tfrz => t_water_freeze_k_1atm
use FatesConstantsMod, only : umolC_to_kgC
use FatesConstantsMod, only : g_per_kg
use EDPftvarcon , only : EDPftvarcon_inst

! -----------------------------------------------------------------------
! Base maintenance respiration rate for plant tissues base_mr_20
! M. Ryan, 1991. Effects of climate change on plant respiration.
! Ecological Applications, 1(2), 157-167.
! Original expression is br = 0.0106 molC/(molN h)
! Conversion by molecular weights of C and N gives 2.525e-6 gC/(gN s)

! Arguments
real(r8), intent(in) :: lmr25top_ft ! canopy top leaf maint resp rate at 25C
! for this pft (umol CO2/m**2/s)
real(r8), intent(in) :: lnc_top ! Leaf nitrogen content per unit area at canopy top [gN/m2]

integer, intent(in) :: ft ! (plant) Functional Type Index
real(r8), intent(in) :: nscaler ! Scale for leaf nitrogen profile
real(r8), intent(in) :: veg_tempk ! vegetation temperature
real(r8), intent(out) :: lmr ! Leaf Maintenance Respiration (umol CO2/m**2/s)

! Locals
real(r8) :: lmr25 ! leaf layer: leaf maintenance respiration rate at 25C (umol CO2/m**2/s)
real(r8) :: lmr25top ! canopy top leaf maint resp rate at 25C for this pft (umol CO2/m**2/s)

! Parameter
real(r8), parameter :: lmrha = 46390._r8 ! activation energy for lmr (J/mol)
Expand All @@ -1994,14 +2015,16 @@ subroutine LeafLayerMaintenanceRespiration(lmr25top_ft, &
! temperature inhibition (25 C = 1.0)



! MLO - Shouldn't these numbers be parameters too?
lmr25top = 2.525e-6_r8 * (1.5_r8 ** ((25._r8 - 20._r8)/10._r8))
lmr25top = lmr25top * lnc_top / (umolC_to_kgC * g_per_kg)


! Part I: Leaf Maintenance respiration: umol CO2 / m**2 [leaf] / s
! ----------------------------------------------------------------------------------
lmr25 = lmr25top_ft * nscaler
lmr25 = lmr25top * nscaler

if ( nint(EDpftvarcon_inst%c3psn(ft)) == 1)then
if ( nint(EDPftvarcon_inst%c3psn(ft)) == 1)then
lmr = lmr25 * ft1_f(veg_tempk, lmrha) * &
fth_f(veg_tempk, lmrhd, lmrse, lmrc)
else
Expand All @@ -2012,7 +2035,74 @@ subroutine LeafLayerMaintenanceRespiration(lmr25top_ft, &
! Any hydrodynamic limitations could go here, currently none
! lmr = lmr * (nothing)

end subroutine LeafLayerMaintenanceRespiration
end subroutine LeafLayerMaintenanceRespiration_Ryan_1991

! ====================================================================================

subroutine LeafLayerMaintenanceRespiration_Atkin_etal_2017(lnc_top, &
nscaler, &
ft, &
veg_tempk, &
tgrowth, &
lmr)

use FatesConstantsMod, only : tfrz => t_water_freeze_k_1atm
use FatesConstantsMod, only : umolC_to_kgC
use FatesConstantsMod, only : g_per_kg
use EDPftvarcon , only : EDPftvarcon_inst

! Arguments
real(r8), intent(in) :: lnc_top ! Leaf nitrogen content per unit area at canopy top [gN/m2]
integer, intent(in) :: ft ! (plant) Functional Type Index
real(r8), intent(in) :: nscaler ! Scale for leaf nitrogen profile
real(r8), intent(in) :: veg_tempk ! vegetation temperature (degrees K)
real(r8), intent(in) :: tgrowth ! lagged vegetation temperature averaged over acclimation timescale (degrees K)
real(r8), intent(out) :: lmr ! Leaf Maintenance Respiration (umol CO2/m**2/s)

! Locals
real(r8) :: lmr25 ! leaf layer: leaf maintenance respiration rate at 25C (umol CO2/m**2/s)
real(r8) :: r_0 ! base respiration rate, PFT-dependent (umol CO2/m**2/s)
real(r8) :: r_t_ref ! acclimated ref respiration rate (umol CO2/m**2/s)
real(r8) :: lnc ! Leaf nitrogen content per unit area at this level [gN/m2]
glemieux marked this conversation as resolved.
Show resolved Hide resolved
real(r8) :: lmr25top ! canopy top leaf maint resp rate at 25C for this pft (umol CO2/m**2/s)

! Parameters
! values from Atkin et al., 2017 https://doi.org/10.1007/978-3-319-68703-2_6
! and Heskel et al., 2016 https://doi.org/10.1073/pnas.1520282113
real(r8), parameter :: b = 0.1012_r8 ! (degrees C**-1)
real(r8), parameter :: c = -0.0005_r8 ! (degrees C**-2)
real(r8), parameter :: TrefC = 25._r8 ! (degrees C)
real(r8), parameter :: r_1 = 0.2061_r8 ! (umol CO2/m**2/s / (gN/(m2 leaf)))
real(r8), parameter :: r_2 = -0.0402_r8 ! (umol CO2/m**2/s/degree C)

! parameter values of r_0 as listed in Atkin et al 2017: (umol CO2/m**2/s)
! Broad-leaved trees 1.7560
! Needle-leaf trees 1.4995
! Shrubs 2.0749
! C3 herbs/grasses 2.1956

lnc = lnc_top * nscaler

if ( nint(EDPftvarcon_inst%c3psn(ft)) == 1)then

! r_0 currently put into the EDPftvarcon_inst%dev_arbitrary_pft
! all figs in Atkin et al 2017 stop at zero Celsius so we will assume acclimation is fixed below that
r_0 = EDPftvarcon_inst%dev_arbitrary_pft(ft)
r_t_ref = r_0 + r_1 * lnc + r_2 * max(0._r8, (tgrowth - tfrz) )

lmr = r_t_ref * exp(b * (veg_tempk - tfrz - TrefC) + c * ((veg_tempk-tfrz)**2 - TrefC**2))

else
! revert to Q10 model for C4 plants, parameter values as described above in Ryan 1991 method

lmr25top = 2.525e-6_r8 * (1.5_r8 ** ((25._r8 - 20._r8)/10._r8))
lmr25top = lmr25top * lnc / (umolC_to_kgC * g_per_kg)

lmr = lmr25 * 2._r8**((veg_tempk-(tfrz+25._r8))/10._r8)
lmr = lmr / (1._r8 + exp( 1.3_r8*(veg_tempk-(tfrz+55._r8)) ))
end if

end subroutine LeafLayerMaintenanceRespiration_Atkin_etal_2017

! ====================================================================================

Expand Down
2 changes: 1 addition & 1 deletion main/EDParamsMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ module EDParamsMod
! temperature acclimation (NOT YET IMPLEMENTED)

integer,protected, public :: maintresp_model ! switch for choosing between leaf maintenance
! respiration model. 1=Ryan (1991) (NOT YET IMPLEMENTED)
! respiration model. 1=Ryan (1991), 2=Atkin et al (2017)
integer,protected, public :: photo_tempsens_model ! switch for choosing the model that defines the temperature
! sensitivity of photosynthetic parameters (vcmax, jmax).
! 1=non-acclimating (NOT YET IMPLEMENTED)
Expand Down
3 changes: 3 additions & 0 deletions main/FatesConstantsMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -74,6 +74,9 @@ module FatesConstantsMod
integer, parameter, public :: hlm_harvest_area_fraction = 1 ! Code for harvesting by area
integer, parameter, public :: hlm_harvest_carbon = 2 ! Code for harvesting based on carbon extracted.

! integer labels for specifying leaf maintenance respiration models
integer, parameter, public :: lmrmodel_ryan_1991 = 1
integer, parameter, public :: lmrmodel_atkin_etal_2017 = 2

! Error Tolerances

Expand Down
5 changes: 3 additions & 2 deletions parameter_files/fates_params_default.cdl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1020,7 +1020,8 @@ data:

fates_damage_recovery_scalar = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ;

fates_dev_arbitrary_pft = _, _, _, _, _, _, _, _, _, _, _, _ ;
fates_dev_arbitrary_pft = 1.7560, 1.4995, 1.4995, 1.7560, 1.7560, 1.7560,
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ckoven what would we call this parameter? I see it matches r_0: "! base respiration rate, PFT-dependent (umol CO2/m**2/s)"

Also, is there a similar parameter in the Ryan model that could also leverage the parameter?

I could add this parameter to the nutrient coupling PR (which has other parameter changes, and assuming it gets through queue more quickly).

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good question, and good idea. Maybe fates_maintresp_atkinetal2017model_baserate? There isn't a direct analogue in the Ryan 1991 model.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

how about: fates_maintresp_atkin2017_baserate

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@rgknox that works for me, thanks!

2.0749, 2.0749, 2.0749, 2.1956, 2.1956, _ ;
glemieux marked this conversation as resolved.
Show resolved Hide resolved

fates_fire_alpha_SH = 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2,
0.2 ;
Expand Down Expand Up @@ -1490,7 +1491,7 @@ data:

fates_leaf_theta_cj_c4 = 0.999 ;

fates_maintresp_model = 1 ;
fates_maintresp_model = 2 ;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also, what is our staging for making the change in the default file? My sense from discussion with others is that the Atkin model is the path forward. However, all of the other parameters are calibrated around the Ryan model, and this will undoubtedly have an impact on output.

Was the plan to switch over to Atkins first, and accept the differences good or bad to output, and then continue with other calibrations? Or was the plan to switch over to Atkins during a more comprehensive calibration exercise perhaps with the satphen simuluation or other biogeographically constrained testbed?

In either situations, we should probably run this through the ilamb mill to get a sense of the changes.

Copy link
Contributor Author

@ckoven ckoven Nov 9, 2022 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think for staging purposes, we should bring in the code first with this switch reverted to the current value, and if we subsequently decide to make it the default model, that would likely need to be accompanied by other default parameter changes as well.


fates_maxcohort = 100 ;

Expand Down