Add USO stomatal model to FATES as alternative to ball-berry

[serbinsh]

I posted this question on Slack:

@charlie @Anthony Walker @rosie fisher @Ryan Knox et al, I forget what we decided with respect to Medlyn vs current B.B. slope?  Should be update to Medly to be consistent with CLM5?  Also how does it work if FATES had B.B. and CLM has medlyn?  Medlyn disabled with FATES?
Should we start an issue on git?

Charlies response:

charlie [11:01 AM]
I don’t think that’s something we ever decided on.  The main result of the work as I understand it in CLM5 is that the structural uncertainty (BB/Medlyn) is much less than the parametric uncertainty of what slope to use. Also the code resides in a completely different place so its not as simple as enabling/disabling the various options. Feel free to start a github issue on this though.

This also relates to the great leaf photosynthesis walkthrough led by @walkeranthonyp last fall. I am interested to hear from others about this, specifically if we want to replace Ball Berry in FATES with the Medlyn Optimal Stomatal Model approach.

My feeling is that, at least for tropical sites, models with VPD vs RH provide a better match to the data. Also, given that CLM5 has moved to Medlyn I wonder if it makes sense to be consistent but that enters in to my general ignorance on how much that would matter in terms of consistency.

OK, so any thoughts on this?

[rosiealice]

It's pretty easy (esp. now there is template code in CLM5) to replace BB with Medlyn. In answer to @serbinsh 's question, the entire photosynthesis routine is unique to FATES, on account of it needing to conduct multi-layer multi-PFT calculations.

In principle, this could be modularized such that basic physiology is common between FATES and the HLM models, but for better or worse, it's not... An upside of this is that within FATES we can play around with whatever we want. A downside is that there is some code duplication.

And yes, most of the uncertainty is in the slope, with the slight exception that Medlyn-operating plants fare better (don't close their stomata as much) in very low humidity compared to Ball-Berry operating plants... (not that that is a real thing, but I find it vaguely amusing to imagine them duking it out nonetheless)

[serbinsh]

Thanks @rosiealice I actually agree that having FATES use its own internal schemes means that, at least for now, we can play around with whatever we want and dont need to impact HLM at all. Perhaps later on we standardize but for now I dont mind the duplication.

The question would be allow for both approaches (with a flag) or just one?

Also, I didnt see any parameter in the FATES list that would be equivalent to a g0 or a cuticular_conductance, is that correct?

[rosiealice]

It's fates_BB_slope, and fates_bbopt_c3/4 for the intercept.

We have a flag in CLM5. The differences are quite small, so that's probably what we might ultimately want to do.

[serbinsh]

Thanks @rosiealice

One thing to note, a slight mistake in the param file metadata

[sserbin@modex fates_parameter_files]$ ncdump fates_params_2troppftclones.c171018_sps.nc | grep fates_bbopt_c*
    double fates_bbopt_c3(fates_scalar) ;
        fates_bbopt_c3:long_name = "Ball-Berry minimum unstressed leaf conductance for C4" ;
        fates_bbopt_c3:units = "umol H2O/m**2/s" ;
    double fates_bbopt_c4(fates_scalar) ;
        fates_bbopt_c4:long_name = "Ball-Berry minimum unstressed leaf conductance for C4" ;
        fates_bbopt_c4:units = "umol H2O/m**2/s" ;
 fates_bbopt_c3 = 10000 ;
 fates_bbopt_c4 = 40000 ;

both say for "C4" probably a good time to fix that since @rgknox is working on the updated param file

[rgknox]

thanks @serbinsh , noted

[serbinsh]

@Qianyuxuan Here is the current issue on adding the USO model. Feel free to expand on this with notes, updates, etc as you make progress

[ckoven]

Hi All on this thread: it's great to see @Qianyuxuan starting to work on this via #604.

One question I had is how shall we handle the new parameter arguments to the code? Given our current parameter infrastructure, I think we basically have a choice to make—we can either (1) generalize our existing stomatal arguments (bbslope and bbopt) to be something like stomatal_conductance_param1, stomatal_conductance_param2, etc. or (2) we can carry around both sets of the specific named BB and USO parameters, and either turn them on or off depending on the stomatal model to use.

This problem is pretty analogous to what we've done for allometry, where we've gone with the first option, and so I guess I'd advocate doing that here too; because as as we add more switches to the model, doing so inevitably complicates running the model and so I think we need to be consistent in how we do so. The downside of that approach is the meaning of the parameters becomes less clear, so we also need to be very consistent in how we then document these options both in the parameter file metadata and in the tech note.

So, anyway: do people broadly agree with this or not? If so, do we want to move the bbopt values to be PFT-indexed rather than C3/C4 at this point as well? Or do we just rename from bbopt to something like cuticular_conductance if it applies to both stomatal models, and keep the values fixed for all C3 and C4 PFTs?

[jkshuman]

@ckoven I agree with you on doing this in a way that is similar to allometry. maintaining consistent practice makes things easier to explain and understand for new and expanding users. Though your suggestion about renaming bbopt to a more verbose explanatory name is also nice. I am a fan of the verbose names.

[alistairrogers]

@ckoven I agree that it's important to keep the nomenclature obvious, either in the parameter file etc. or through option 2 - but don't fully understand the consequences of either approach so happy to defer. Naming the g0 term cuticular conductance - whilst assumed - adds an extra level of confusion. I would advocate for something obvious e.g. BallBerrySlope and BallBerryIntercept. I also think we want to move the values to PFT-indexed. Just to be clear the intercept or g0 values will differ between formulations and not just the slopes (m or g1).

[serbinsh]

Hey thanks for all the great comments, suggestions! I sat with @Qianyuxuan today and looked at what she had coded up thus far and some prelim results. @alistairrogers and I will be sitting with her more tomorrow for a more detailed review of current comparisons with BB and model behavior at the diurnal to monthly scales

From what I saw, Cherry has put much of what is needed already in her version of the FATES code and is successfully reading the new Medlyn param values from the parameter file. I did however suggest that she moves the scalar cuticular conductance to the ParamsMods like the bboptc3 and bbopt4 variables; for now, though at some point we may want it in the EDPftvarcon(?) to allow for PFT specific g0, if we think we can parameterize this. So for the most part close to an operational prototype.

One issue I found was that Cherry was trying to utilize namelist functions from the HLM and I pointed out this all had to occur on the FATES side. As such, one last piece of the puzzle, I think, is that we need to set a switch somewhere on the FATES side to toggle gs type. I know we were hesitant for this to be a param option, so where would be the place to set this switch as a namelist (?) or similar option somewhere in the model setup/configuration step or prior to compilation such that the gs model type is set in the compiled code. Is that correct? Or perhaps that does have to be configured on the HLM side, though my understanding is that we should only have to modify the FATES side.

So, it seems we potentially have one more step before we can share a version that allows you to select/toggle and set the g1/slope (PFT specific) and g0 parameter (scalar)in the fates param file. Make sense?

[rosiealice]

I guess I'd advocate for having all the stomatal conductance parameters as PFT specific. In @danicalombardozzi 's minimum conductance paper, for example, we modified g0 as a PFT specific term and it would be annoying for someone to have to reimplement that capacity again.. .https://www.researchonline.mq.edu.au/vital/access/services/Download/mq:69411/DS01

[walkeranthonyp]

@ckoven on the whole I support generalisation of parameter names, or at least a common naming convention. Though also trying to maintain some level of identification. I guess this is a not very helpful comment. g is the common identifier in the code for conductance related terms, perhaps we could go with a prefix something like: gs_paramX... . The rest somewhat depends on how unified the parameters can be. If they're all different (i.e. separate slope and intercept for each gs model), maybe gs_param1_g0_ballberry and gs_param1_g0_medlyn or gs_param2_g0_medlyn . g0 could be replaced with intercept if preferred or if these parameters can be unified then the final suffix can be dropped.

@alistairrogers Why do you say g0 / bbopt will be different between medlyn and BB? Because when fitting the model empirically they come out different? I hadn't considered that but that makes sense. I think I would argue that conceptually they're the same, but maybe I'm mistaken. Also that the fitting means that g0 and g1 are correlated and g0 is estimated outside of the range of the data so we're less confident in g0, plus are they really different statistically? My guess would be that the differences would be within the range of uncertainty.

Also, based on those arguments isn't the g1 estimate for the medlyn model generally fit now assuming a zero intercept?

This brings up an issue that strictly speaking the intercept term shouldn't feature in the quadratic soln for gs with Medlyn ... and perhaps not for ball berry either. But maybe that's a discussion for another day ... it's unlikely to make much of a difference.

[danicalombardozzi]

Having PFT-level parameter makes the code more flexible, but I suggest being careful with changing the g0 parameters. Our paper shows that changing g0 can have a large impact because it's a value that is added to the equation at all times, and changes will require also changing the g1 parameter. We worked on this in the context of nighttime conductance, which is not accounted for in the BB or Medlyn conductance equations.

The Medlyn g1 parameter is fit assuming a zero intercept, but that caused CLM to crash because conductance then equaled zero at some times. Our solution was to include a very small g0 to keep the model from crashing.

On Tue, Jan 28, 2020 at 6:12 PM walkeranthonyp notifications@github.com wrote:

@ckoven https://github.com/ckoven on the whole I support generalisation of parameter names, or at least a common naming convention. Though also trying to maintain some level of identification. I guess this is a not very helpful comment. g is the common identifier in the code for conductance related terms, perhaps we could go with a prefix something like: gs_paramX... . The rest somewhat depends on how unified the parameters can be. If they're all different (i.e. separate slope and intercept for each gs model), maybe gs_param1_g0_ballberry and gs_param1_g0_medlyn or gs_param2_g0_medlyn . g0 could be replaced with intercept if preferred or if these parameters can be unified then the final suffix can be dropped.

@alistairrogers https://github.com/alistairrogers Why do you say g0 / bbopt will be different between medlyn and BB? Because when fitting the model empirically they come out different? I hadn't considered that but that makes sense. I think I would argue that conceptually they're the same, but maybe I'm mistaken. Also that the fitting means that g0 and g1 are correlated and g0 is estimated outside of the range of the data so we're less confident in g0, plus are they really different statistically? My guess would be that the differences would be within the range of uncertainty.

Also, based on those arguments isn't the g1 estimate for the medlyn model generally fit now assuming a zero intercept?

This brings up an issue that strictly speaking the intercept term shouldn't feature in the quadratic soln for gs with Medlyn ... and perhaps not for ball berry either. But maybe that's a discussion for another day ... it's unlikely to make much of a difference.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/NGEET/fates/issues/350?email_source=notifications&email_token=AGHW2QI7N4SMHMH52MSJE5TRADJXPA5CNFSM4EUSFCTKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOEKFTGBQ#issuecomment-579547910, or unsubscribe https://github.com/notifications/unsubscribe-auth/AGHW2QOHME6ICZSHL6N4CE3RADJXPANCNFSM4EUSFCTA .

-- Dr. Danica Lombardozzi she/her/hers Terrestrial Sciences Section Climate and Global Dynamics National Center for Atmospheric Research Boulder, CO 80305

email: dll@ucar.edu office: (303) 497-1777

[walkeranthonyp]

@danicalombardozzi just to be clear - I'm not advocating to make g0 = 0. I'm advocating, or at least pointing out, if g1 is fit assuming g0 is zero, then g0 should not feature in the quadratic solve for gs. If it was decided that was the way forward then I would still advocate for g0 to be applied whenever gs < g0 or when anet < 0 .

My guess is that would also be a more appropriate way to account for night-time gs.

[danicalombardozzi]

@walkeranthonyp, thanks for the clarification. The method you suggest -- g0 can be applied whenever gs < g0 -- is one method we tested for use in nighttime conductance. One thing that we found when implementing measured values of nighttime conductance this way was that gs < g0 during midday in several instances, leading us to question whether measured nighttime conductance was equivalent to daytime minimum conductance. To the best of my knowledge, we still do not have measurements to determine whether nighttime conductance is higher (under the same moisture conditions) than minimum daytime conductance. You can find the full discussion here: https://www.geosci-model-dev.net/10/321/2017/

On Thu, Jan 30, 2020 at 8:24 AM walkeranthonyp notifications@github.com wrote:

@danicalombardozzi https://github.com/danicalombardozzi just to be clear - I'm not advocating to make g0 = 0. I'm advocating, or at least pointing out, if g1 is fit assuming g0 is zero, then g0 should not feature in the quadratic solve for gs. If it was decided that was the way forward then I would still advocate for g0 to be applied whenever gs < g0 or when anet < 0 .

My guess is that would also be a more appropriate way to account for night-time gs.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/NGEET/fates/issues/350?email_source=notifications&email_token=AGHW2QMYTYAXIFZ4IUH2PM3RALWJHA5CNFSM4EUSFCTKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOEKLMBBQ#issuecomment-580305030, or unsubscribe https://github.com/notifications/unsubscribe-auth/AGHW2QMCTW5PIOHUUS6OULTRALWJHANCNFSM4EUSFCTA .

[serbinsh]

This is addressed by PR https://github.com/NGEET/fates/pull/608