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.

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

parteh v1 #421

Merged
merged 46 commits into from
Oct 24, 2018
Merged

parteh v1 #421

merged 46 commits into from
Oct 24, 2018

Conversation

rgknox
Copy link
Contributor

@rgknox rgknox commented Oct 3, 2018

Description:

This changeset contains the first implementation of the Plant Allocation and Reactive Transport Extensible Hypotheses. The parteh system is designed to enable flexibility in which organs and species are enabled in a plant allocation hypothesis, and also how many discrete spatial positions for those organs are enabled.

Some notes:

  1. Only the carbon-only hypothesis, pre-existing, is implemented
  2. Generic/object looping for populating and reading restarts has been added, but nothing special for history output
  3. A functional unit testing framework has also been added, see the new folder functional_unit_testing/parteh
  4. New and changed parameters
  5. Re-translocation of C&N&P has been added
  6. Helper functions for event based leaf-drop, leaf-flush, burning and seed drop have been added

Synced with: ESCOMP/CTSM#537

Documentation PR: NGEET/fates-docs#13

Collaborators:

@ckoven, @bishtgautam , @rosiealice , @xuchongang , @jkshuman , @tompowell9 , @bchristo , @jenniferholm

Expectation of Answer Changes:

Answers are expected to change, but only slightly. So much code has changed, but it is intended to generate the same answers within reasonable precision.

Checklist:

  • My change requires a change to the documentation.
  • I have updated the in-code documentation .AND. (the technical note .OR. the wiki) accordingly.
  • I have read the CONTRIBUTING document.
  • FATES PASS/FAIL regression tests were run
  • If answers were expected to change, evaluation was performed and provided

Test Results:

30 Year regression test against master at BCI (no inventory): https://drive.google.com/file/d/1HwoBhMw1QBR_2qcCyzTmagouismT_C9b/view?usp=sharing

FATES test suite TBD

CTSM (or) E3SM (specify which) test hash-tag:

CTSM (or) E3SM (specify which) baseline hash-tag:

FATES baseline hash-tag:

Test Output:

rgknox and others added 23 commits September 18, 2018 15:19
@rgknox rgknox added enhancement PR status: Not Ready The author is signaling that this PR is a work in progress and not ready for integration. labels Oct 3, 2018
@rgknox
Copy link
Contributor Author

rgknox commented Oct 3, 2018

FYI: I plan to go through and do a few more passes, cleaning and commenting code. But I wanted to get this PR established so others can familiarize with changes at earliest convenience.

As you can see in the the acre tests, this seems to be generating fairly similar results compared to master, and seems to be stable. More testing is on its way, including a multi-decade f45 smoke test.

@rgknox
Copy link
Contributor Author

rgknox commented Oct 9, 2018

Conducted another test using the setup described in #424 (gridded run over South America, 80 year simulation, tropical evergreen and c3 stress deciduous grass).

TreeGrass-parteh-v2_plots.pdf

@@ -1017,7 +1080,7 @@ 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%treelai = tree_lai(leaf_c, &
currentCohort%pft, currentCohort%c_area, currentCohort%n, &

Choose a reason for hiding this comment

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

Apologies for this additional naive question but I dont see SLA in this subroutine/function to calculate LAI for each cohort, or is that the c_area?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

@kovenock updated the LAI calculations to use SLA at the top of the plant, and optionally a higher SLA that increases as one decreases in elevation in the canopy. But at its heart, the tree_lai() routine should be making use of these parameters. It is a parameter, not a dynamic quantity, so it does not need to be passed in as an argument.

Choose a reason for hiding this comment

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

Got it, thanks!

Choose a reason for hiding this comment

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

@rgknox as far as my comment about bleaf vs leaf_c, here is one of the examples using leaf_c. Just wanted to make sure whether you had planned for everything to be leaf_c or if in the other case bleaf is appropriate?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

My take is that we should drop the b_ convention altogether, except we want to describe things in terms of total biomass (including not just carbon, but lignans, water, nutrients, etc). We do have a subroutine called bleaf(), and I want to change that too. That subroutine is actually calculating the plants target leaf carbon for a fully flushed, on-allometry, plant of a given size and trimming. So bleaf(), is not a very concise subroutine name.

@@ -150,7 +204,7 @@ 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%pft, new_cohort%c_area, &
new_cohort%treelai = tree_lai(bleaf, new_cohort%pft, new_cohort%c_area, &

Choose a reason for hiding this comment

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

Just curious why the calculation of LAI here is slightly different than how it is done in biogeochem/EDCanopyStructureMod.F90 where here we have bleaf (leaf biomass) and above leaf_c which is also leaf biomass? Just trying to understand the steps....

Copy link
Contributor Author

Choose a reason for hiding this comment

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

@serbinsh , the stuff that mostly happens in EDCanopyStructureMod, is the tallying of LAI from each cohort, and assembling those LAI's into canopy, leaf and PFT bins within each patch. But that part of the code shouldn't be calculating the individual plant's LAI any differently. Could you point to the line in question so that I can double check?

call SetState(currentCohort%prt,sapw_organ,carbon12_species,0.0_r8)
call SetState(currentCohort%prt,struct_organ,carbon12_species,0.0_r8)
call SetState(currentCohort%prt,repro_organ,carbon12_species,0.0_r8)
call SetState(currentCohort%prt,store_organ,carbon12_species,0.0_r8)

Choose a reason for hiding this comment

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

So here we we set these pools we are using carbon12 (specific isotopic version) and above we are just using "all_carbon_species" Here we would allow for isotopic differentiation?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Right, good question. When we are using SetState(), we are initializing a specific variable, we can't set a set of them. But in most other cases, when we are querying, we just want to know the total mass, so we use "all_carbon_species" to signal that we want all carbon isotopes.

If you have a scenario, where on the FATES side of things you want to ask PARTEH for the leaf C13 for instance, you could call prt%GetState(leaf_organ,carbon13_species). ... Assuming that you have that isotope enabled in the specify PARTEH hypothesis you have turned on.

@@ -57,7 +59,8 @@ subroutine mortality_rates( cohort_in,bc_in,cmort,hmort,bmort,frmort )
real(r8),intent(out) :: frmort ! freezing stress mortality

real(r8) :: frac ! relativised stored carbohydrate
real(r8) :: b_leaf ! target leaf biomass kgC
real(r8) :: leaf_c_target ! target leaf biomass kgC

Choose a reason for hiding this comment

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

so again to clear up some possible confusion, we seem to still be interchanging b_leaf with leaf_c in some places, is that correct? Or should we replace all previous b_leaf with leaf_c to represent leaf biomass?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

these naming conventions should definitely be cleared up. I'm pretty sure that anywhere in the code you see "b_leaf", it actually is referring to the target leaf biomass, ie the value of the pool when it is on allometry and fully flushed.

Choose a reason for hiding this comment

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

OK got it

call bleaf(cohort_in%dbh,cohort_in%pft,cohort_in%canopy_trim,b_leaf)
call storage_fraction_of_target(b_leaf, cohort_in%bstore, frac)

call bleaf(cohort_in%dbh,cohort_in%pft,cohort_in%canopy_trim,leaf_c_target)

Choose a reason for hiding this comment

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

Or maybe I am accidentally conflating things since I also see this sub is called bleaf.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

that subroutine should be renamed to convey that it is returning the allometric target leaf carbon

Choose a reason for hiding this comment

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

aha! that makes sense then. I agree renaming would reduce confusion

if(EDPftvarcon_inst%woody(currentCohort%pft) == 1)then
burned_leaves = currentCohort%bl * currentCohort%fraction_crown_burned
burned_leaves = leaf_c * currentCohort%fraction_crown_burned

Choose a reason for hiding this comment

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

So here leaf_c is a mass and this X area fraction gives remaining mass of burned leaves?

Copy link
Contributor

Choose a reason for hiding this comment

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

It's not an area fraction, it's more like a fraction of the crown burned, calculated from flame height.

Choose a reason for hiding this comment

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

so a 1D fraction so say, 5kg of leaf biomass * 40% burned crown = 2 kg of living biomass / 3 kg of burned? where burned_frac = f(flame height)?

Copy link
Contributor

Choose a reason for hiding this comment

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

Yes, that sort of thing.

canopy_mortality_root_litter(p) = canopy_mortality_root_litter(p) + &
canopy_dead*(currentCohort%br+currentCohort%bstore)
canopy_dead*(fnrt_c + store_c)

Choose a reason for hiding this comment

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

@rgknox FYI use of "dead" here - justing pointing you to it

Copy link
Contributor Author

Choose a reason for hiding this comment

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

thanks @serbinsh , in this case, dead really refers to dead trees, specifically the number density of dead trees


leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_species)

currentCohort%treelai = tree_lai(leaf_c, currentCohort%pft, currentCohort%c_area, &

Choose a reason for hiding this comment

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

sorry just continuing this thread, here leaf_c in the calc of LAI for each cohort.....

@@ -366,6 +410,12 @@ subroutine phenology( currentSite, bc_in )
integer :: day ! day of month (1, ..., 31)
integer :: sec ! seconds of the day

real(r8) :: leaf_c ! leaf carbon [kg]

Choose a reason for hiding this comment

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

minor, comment alignment is off (my OCD is showing)

Copy link
Contributor

Choose a reason for hiding this comment

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

You know you can suggest edits with the new button for these things?

Choose a reason for hiding this comment

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

I couldnt figure out how to "align" it with the other comments in that block otherwise I would have done that ;)

Copy link
Contributor Author

Choose a reason for hiding this comment

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

ok, will change

currentCohort%bl = 0.0_r8

! Remember what the lai (leaf mass actually) was for next year
currentCohort%laimemory = leaf_c

Choose a reason for hiding this comment

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

Huh? "next" year? Also why do we consider leaf mass also and area? Personally I think that is confusing

Copy link
Contributor

Choose a reason for hiding this comment

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

I guess technically it should be leafc_memory... Do we still need this, @rgknox, in the presence of the new flushing parameter and the target-centric allocation scheme. It might be better, although not bfb, to have leaf flushing as a fraction of the target and forget about lai_memory?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I am fine with that @rosiealice . It seems a little cleaner to me, but I didn't want to mess with existing stuff. If it has your blessing I will remove laimemory and just use the target.

@@ -580,7 +580,10 @@ real(r8) function tree_lai( bl, pft, c_area, nplant, cl, canopy_lai)
!----------------------------------------------------------------------

if( bl < 0._r8 .or. pft == 0 ) then
write(fates_log(),*) 'problem in treelai',bl,pft

Choose a reason for hiding this comment

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

Is that a sufficient warning/error message? Idea is to then go to that portion of the output?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I am not a fan of this error message either @serbinsh . In fact, it triggers quite a bit, because machine math precision will generate small negatives that are effectively zero, like 1e-19. But it will complain. I think this comparison should be against, not zero, but a some number close to zero but higher than machine precision, like our carbon allocation error tolerance. If it triggers, it should fail, and be more verbose.

Actually, I just realized I had deleted the warning.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

i will remove the if clause (although compilers will remove it anyway I think)

@rgknox
Copy link
Contributor Author

rgknox commented Oct 19, 2018

FYI, github is not updating the code changes I make to rectify many of your comments. Or at least the updates are not appearing "inline" on this thread. The git side of things (the branch itself) seems to updating correctly though.

@rgknox rgknox added PR status: Ready This PR is ready for integration after having passed through peer review and final testing. and removed PR status: Not Ready The author is signaling that this PR is a work in progress and not ready for integration. labels Oct 22, 2018
@rgknox
Copy link
Contributor Author

rgknox commented Oct 24, 2018

Performed a final comparison with master. 100 year simulation, BCI, spin-up from bare ground, default 2 tropical clones.

acre report here: https://drive.google.com/file/d/1B2Asdo9obQZ0WvOWDbkZ7dlQ-dp-858z/view?usp=sharing

I compared their timing information, the parteh version was about 3% slower. I was expecting it to be slower by a small margin, I suppose that is about what 3% is. The flexibility that the code offers in connecting arbitrary hypotheses to the rest of the code comes at a marginal performance cost.

Incidentally, @jkshuman and I have both not experienced noticable performance hits on global gridded runs.

@rgknox
Copy link
Contributor Author

rgknox commented Oct 24, 2018

We have a few other PR's submitted right now that are fairly straight forward.
I will integrate those that can be done without deliberations or checks. Then I will integrate PARTEH. For those PRs that do not get integrated before parteh, I will help with the conflict resolutions.

@rgknox
Copy link
Contributor Author

rgknox commented Oct 24, 2018

I went through the list of PRs. Only one other PR is ready for testing, and the changes are quite small, so I will integrate parteh now.

@rgknox rgknox merged commit 5a9b1d6 into NGEET:master Oct 24, 2018
@rgknox rgknox deleted the rgknox-parteh branch October 31, 2023 18:57
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement PR status: Ready This PR is ready for integration after having passed through peer review and final testing.
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants