-
Notifications
You must be signed in to change notification settings - Fork 92
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.
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
Drought deciduous phenology updates #958
Conversation
… representation. A summary of the changes: 1. The soil moisture threshold is now defined for the rooting layer instead of at an arbitrary soil depth. 2. Because the rooting depth is PFT-dependent, most variables driving phenology became PFT variables too. 3. The code now accepts thresholds both in terms of liquid volume (m3/m3) or matric potential (mm) 4. Drought phenology can now be defined either as a "hard" deciduous or semi-deciduous. The difference is that in semi-deciduous, plants may reduce foliage when soil is moderately dry, and only shed all leaves under very dry conditions. 5. The code now allows PFTs to actively abscise fine roots as a strategy to reduce carbon costs. There will be additional implementations in upcoming commits.
@mpaiao , with the elongation factors, can we get rid of the leaf status variable? |
@rgknox I don't think so, because of the semi-deciduous phenology. For these PFTs, the elongation factor tells the maximum fraction of leaves that the PFT can have, relative to on-allometry under well-watered conditions, whereas the leaf status informs whether conditions are improving (elongation factor is increasing or at maximum) or deteriorating (elongation factor is decreasing or at minimum). The leaf status is still used to decide how to allocate carbon: when elongation factor is deteriorating, plants stop allocating carbon to live tissues, growth and reproduction, even if their NPP is above zero, which is now possible for semi-deciduous PFTs. |
For now, semi-deciduous trees should not attempt to grow or reproduce if the conditions are deteriorating (in reality many drought deciduous trees do bloom during the dry season, but this would require seed phenology).
Some elongation factor parameters were missing when initialising intgr_params, which likely caused the elongation factors to be unrealistic. I also had forgotten to update one deciduous test to check for both hard and semi-deciduous. Some of
biogeochem/EDPhysiologyMod.F90
Outdated
! MLO: Need to double check that this will not cause unintended problems | ||
pft_leaf_lifespan = min( canopy_leaf_lifespan, & | ||
nint(365*maxval(prt_params%leaf_long(ipft,:))) ) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Flagging this part for revision...
Parameter canopy_leaf_lifespan
(365 days) is defined locally, and this is used to force leaf shedding when leaves have been flushed for longer than their leaf longevity. However, we also have leaf longevity in the parameter file (fates_turnover_leaf
), and I wonder if this wouldn't make the leaf dynamics inconsistent with the parameter set. I proposed this change so the time for shedding is now the shorter between fates_turnover_leaf
and canopy_leaf_lifespan
. Does this make sense?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this makes sense. I would only suggest that we explain the dual meaning in the parameter file and the parameter definition in-code.
fates_turnover_leaf:long_name = "Leaf longevity (ie turnover timescale)" ; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also, FYI: the leaf lifespan also impacts leaf aging. My sense is that many simulations out there do not use variable leaf age impacts on vcmax, but it is an option. I'm not saying what you are proposing is in conflict with it, just wanted to point it out.
https://github.com/NGEET/fates/blob/main/parteh/PRTGenericMod.F90#L1312
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this makes sense. I would only suggest that we explain the dual meaning in the parameter file and the parameter definition in-code.
fates_turnover_leaf:long_name = "Leaf longevity (ie turnover timescale)" ;
Sure, I will add a brief explanation in both places. By the way, I find the parameter name fates_turnover_leaf
slightly confusing, because it reminds me of leaf turnover rate (yr-1)...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also, FYI: the leaf lifespan also impacts leaf aging. My sense is that many simulations out there do not use variable leaf age impacts on vcmax, but it is an option. I'm not saying what you are proposing is in conflict with it, just wanted to point it out.
https://github.com/NGEET/fates/blob/main/parteh/PRTGenericMod.F90#L1312
Thanks for the heads up, but my impression is that this change would not interfere with the leaf age impacts on Vcmax. The code already forced leaf shedding once the cohort's leaf age reached 1 year. The change here was simply to shorten the leaf "expiration date" to be the maximum leaf longevity across all leaf age classes when this value was < 1 year.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
shouldn't that be sum() (or something similar) instead of maxval then?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Indeed, thanks! I will submit an update shortly.
@@ -434,12 +434,12 @@ variables: | |||
double fates_phen_evergreen(fates_pft) ; | |||
fates_phen_evergreen:units = "logical flag" ; | |||
fates_phen_evergreen:long_name = "Binary flag for evergreen leaf habit" ; | |||
double fates_phen_fnrt_drop_fraction(fates_pft) ; | |||
fates_phen_fnrt_drop_fraction:units = "fraction" ; | |||
fates_phen_fnrt_drop_fraction:long_name = "fraction of fine roots to drop during drought/cold" ; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@rgknox Quick follow-up to our Zoom discussion yesterday on parameter names.
If you think this is fine, I would rather keep the new name fates_phen_fnrt_drop_fraction
instead of fates_phen_fnrt_drop_frac
just so the parameter names are consistent with similarly named variables (fates_phen_stem_drop_fraction
and fates_phen_flush_fraction
).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I see that we did not maintain this name change request, sorry about that. We can check-in about this.
…rameters from EDPftvarcon.F90 to PRTParametersMod following discussion with Ryan.
…reshold is fixed.
…ught deciduous. Added the description that the leaf longevity parameter is also the maximum time that leaves can be "on" when running drought-deciduous PFT, and fixed the check in sub-routine "phenology" to compare the length with the growing season with the sum of leaf longevities across ages (as opposed to the maximum value).
biogeochem/EDPhysiologyMod.F90
Outdated
drought_smoist_ifelse: if (model_day_int <= numWaterMem) then | ||
! Too early in the simulation. Do not change phenology status as we need to | ||
! populate the soil moisture memory. | ||
continue |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm not too familiar with the continue statement in modern fortran, and documentation is sparse. In some documentation its an alternative to an end do statement. It seems more of a legacy for fortran 77. TLDR, what does this continue do, and do we really need it?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Continue does nothing in this case. It may be my incorrect understanding, but I think each if block needs at least one command in Fortran, otherwise the compiler doesn't work. I could rewrite this block by putting all the other cases in a nested block like this:
drought_initial_ifelse: if (model_day_int > numWaterMem) then
drought_smoist_ifelse: if ( prolonged_off_period .and. ( .not. smoist_below_threshold ) ) then
...
end if drought_smoist_ifelse
end if drought_initial_ifelse
The only reason why I hadn't done so was because I find nested if
blocks a bit more difficult to read.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I found this on stack overflow. The first method on the upvoted solution seems reasonable:
https://stackoverflow.com/questions/25676793/break-from-select-case
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, the first method looks more elegant than my solution, I will update this in my next commit.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I changed this in my branch and I will push the change once I incorporate yours and Rosie's suggestions.
I will just add the caveat that exit
may not work for older Fortran compilers (before Fortran 2008) because exit used to work exclusively with do
loops. Do you think this would be a problem?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could you replace the exit command with some sort of 'do while' thing?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think I can replace it with a nested if statement, this way we don't use continue
and also make sure the solution works with any fortran compiler.
@mpaiao I curious why you didn't put the elongation factors in the allometry routines. I suppose it could be a argument about definition, and what we mean by a "target". On the other hand though, I'm worried that we will be applying the elongation factors in some places, but not in others, and if they are part of the allometry routine, then it should protect us from not applying them ubiquitously. |
Oh, I guess I didn't include it because it never occurred to me... I think the only potential danger of including the elongation factor in the allometry is that anything that depends on the fully flushed leaves, fine roots etc. could get the total incorrect. But I may be overthinking this, what do you think? |
IIRC, the sapwood allometry is proportional to the leaf allometry. So for that one we would assume fully elongated leaves. My sense is that it would be safe to apply these in the allometry formulations. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Hi @mpaiao, I finally had a look through this. I am glad I did because now I know how it works! Hard to predict what will happen a priori from the code, so I guess I am infested to see how it performs in the pesky edge cases of the drought deciduous world (and indeed if any of the previous biannual oscillating behaviour of the old algorithm is a feature of this new code). Most of my comments are pretty minor. I am not sure I am a good reviewer of code these days! Very happy to see this happening though :)
@@ -551,6 +573,9 @@ subroutine nan_cohort(cc_p) | |||
currentCohort%canopy_layer_yesterday = nan ! recent canopy status of cohort (1 = canopy, 2 = understorey, etc.) | |||
currentCohort%NV = fates_unset_int ! Number of leaf layers: - | |||
currentCohort%status_coh = fates_unset_int ! growth status of plant (2 = leaves on , 1 = leaves off) | |||
currentCohort%efleaf_coh = nan ! leaf elongation factor (fraction from 0 (fully abscissed) to 1 (fully flushed) | |||
currentCohort%effnrt_coh = nan ! fine-root "elongation factor" (fraction from 0 (fully abscissed) to 1 (fully flushed) | |||
currentCohort%efstem_coh = nan ! stem "elongation factor" (fraction from 0 (fully abscissed) to 1 (fully flushed) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
What type of stem tissue?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Both sapwood and heartwood. This is something that was already implemented in the code, my understanding is that this is needed for non-woody plants when running FATES-Hydro. Both fine-root and woody phenology are controlled by parameters (fates_phen_fnrt_drop_fraction and fates_phen_stem_drop_fraction, respectively). If these parameters are set to zero, then the drought deciduous code will not actively abscise fine roots and woody tissues due to phenology.
! This is the shortest between the PFT leaf lifespan | ||
! and the maximum lifespan of drought deciduous | ||
! (see canopy_leaf_lifespan below) | ||
real(r8) :: phen_drought_threshold ! For drought hard-deciduous, this is the threshold |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Are these the same for all PFTs?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
They used to be, but not anymore ;).
This is just a local variable that takes the value of prt_params%phen_drought_threshold(ipft)
, to make the if
statements a bit shorter.
|
||
! PFT leaf lifespan in days. This is the shortest between the leaf longevity | ||
! (defined as a PFT parameter) and the maximum canopy leaf life span allowed | ||
! for drought deciduous (local parameter). The sum term accounts for the |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Local as in site specific?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Actually, local as in defined locally in this subroutine (as opposed to a FATES parameter). I didn't define this parameter, it was already there when I started working on the sub-routine, but I do wonder whether this should be made a more "formal" FATES parameter or not.
|
||
|
||
! Find elongation factor by comparing the moisture with the thresholds. | ||
case_drought_phen: select case (prt_params%stress_decid(ipft)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can you maybe add a comment on what this case statement is doing? I for one still find them a bit unintuitive...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Sure, I will include more comments in the updated pull request. This case is just deciding whether the PFT is a "hard" deciduous or a semi-deciduous PFT. If the PFT is a hard deciduous, elongation factors are going to be either 0 (completely leafless) or 1 (fully flushed). If the plant is semi-deciduous, then elongation factors can be any value between 0 and 1, depending on how dry the rooting zone soils are.
biogeochem/EDPhysiologyMod.F90
Outdated
!------ | ||
|
||
|
||
! First guess elongation factor |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Maybe describe why this needs a first guess?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I added the comment. The first guess is entirely based on soil water availability. But we may force elongation factors to be zero or one depending on the time since last flushing or abscising event.
integer :: cleafon ! DOY for cold-decid leaf-on, initial guess | ||
integer :: cleafoff ! DOY for cold-decid leaf-off, initial guess | ||
integer :: dleafoff ! DOY for drought-decid leaf-off, initial guess | ||
integer :: dleafon ! DOY for drought-decid leaf-on, initial guess | ||
integer :: cndleafon ! days since leaf on (cold), initial guess | ||
integer :: cndleafoff ! days since leaf off (cold), initial guess | ||
integer :: dndleafon ! days since leaf on (drought), initial guess |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why is it dn at the beginning ?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The "d" at the beginning is for "drought deciduous", and "c" is for cold deciduous, the "nd" is short for number of days.
main/EDTypesMod.F90
Outdated
@@ -83,12 +83,32 @@ module EDTypesMod | |||
! can be approximated to be equal to the visible band | |||
|
|||
|
|||
integer, parameter, public :: leaves_pshed = 3 ! Flag specifying that a deciduous plant has leaves |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Could we call this leaves_shedding ?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, I will update it.
@@ -138,6 +158,7 @@ module EDTypesMod | |||
integer, parameter, public :: phen_dstat_moistoff = 1 ! Leaves off due to moisture avail (drought phenology) | |||
integer, parameter, public :: phen_dstat_moiston = 2 ! Leaves on due to moisture avail (drought phenology) | |||
integer, parameter, public :: phen_dstat_timeon = 3 ! Leaves on due to time exceedance (drought phenology) | |||
integer, parameter, public :: phen_dstat_pshed = 4 ! Leaves partially abscissing (drought phenology) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Is there a status for the leaves partly flushing? Do we need one?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I originally had the partial flushing status, but it turned out that we didn't need one. As long as the plants are growing, the elongation factors should control the maximum amount of biomass allowed. The partial shedding status was needed for something, but to be honest, now I can't remember why...
UPDATE: I checked the code for the difference. This is needed when making decisions about flushing or shedding leaves. Plants that only partially abscised leaves can adjust their elongation factors up and down, but they shouldn't do it too frequently when they are completely leafless.
@@ -1246,12 +1243,20 @@ data: | |||
|
|||
fates_phen_cold_size_threshold = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; | |||
|
|||
fates_phen_drought_threshold = 0.15, 0.15, 0.15, 0.15, 0.15, 0.15, 0.15, 0.15, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I feel like it makes more sense for this to default to being in pressure units (even if those are mm!) as that is more physiologically relevant...?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I didn't change the defaults just because they already existed, but I agree with you that defining them based on pressure is more intuitive.
parteh/PRTParamsFATESMod.F90
Outdated
@@ -1054,23 +1105,108 @@ subroutine PRTCheckParams(is_master) | |||
pftloop: do ipft = 1,npft | |||
|
|||
! Check to see if evergreen, deciduous flags are mutually exclusive | |||
! MLO. Changed the check because season_decid can be 1 or 2. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
! MLO. Changed the check because season_decid can be 1 or 2. |
Thank you very much for reviewing the code @rosiealice, really appreciate your inputs as you were the one who started implementing this feature! I will be working on two proposals for the next 10 days, but I will get back to this once I submit them. |
This also updates the comment header identifying which changes are associated with which pull request
these will be utilitized in a future pull request
…paiao-pr-drgt-decid
Revised the lower threshold for soil potential (and turned into a local constant).
smp and parameter updates
This also updates the output default parameter file
Regression testing for the fates suite against
I will look more into this shortly. Folder location: |
Retesting the above ERP test individually resulted in it passing RUN. I neglected to compare it against the baseline. Folder location: UPDATE: Resubmitting the original failed case passes the run as well. As with most of the other tests, it is not B4B. |
Regression testing on izumi is complete. All expected tests pass against Folder location: |
This is corrected via ESCOMP/CTSM@74f9bd1. The NEON testmods were rerun and pass through all test phases successfully (except the baseline comparison as there is nothing to compare against as the baseline didn't have this fix). This is ready to integrate. New izumi test results: |
Multiple updates and revision of how drought deciduousness is handled in FATES. This pull request supersedes PR #801.
Description:
fates_phen_drought_threshold
: when the parameter is positive, it continues to use the soil volumetric water content (m3/m3), and when it is negative, it uses soil matric potential (in mm).fates_phen_stress_decid
to 2 (it will be revisited in a future pull request). With this option, parameterfates_phen_drought_threshold
becomes the lower soil availability threshold, below which PFTs are completely leafless. A similar parameterfates_phen_moist_threshold
defines the limit below which semi-deciduous PFTs start abscising leaves and is defined similarly tofates_phen_drought_threshold
. The maximum fraction of remaining leaves (elongation factor) is linearly defined between 0 (when soil moisture/potential isfates_phen_drought_threshold
) and 1 (when soil moisture/potential isfates_phen_moist_threshold
).fates_phen_fnrt_drop_fraction
. This parameter allows drought-deciduous PFTs to actively abscise fine roots during drought conditions, allowing trees to further reduce carbon losses when leaves are off.This builds on now defunct pull request #801.
Fixes #1046
This pull request should be coordinated with ESCOMP/CTSM#2009.
Collaborators:
Discussions with @rgknox @glemieux @ckoven @rosiealice @jkshuman @XiulinGao
Expectation of Answer Changes:
These changes should affect the code when running with drought deciduous phenology. The code now uses a different approach to find moisture threshold (rooting zone and optionally, matric potential), and the carbon allocation strategies changed during the leaf-off periods.
Checklist:
Test Results:
CTSM (or) E3SM (specify which) test hash-tag:
CTSM (or) E3SM (specify which) baseline hash-tag:
FATES baseline hash-tag:
Test Output: