Sundials (IDA & KINSOL), Actors, BMI update

[ashleymedin]

This is the sundials merge. It uses cmake with compiler directives so that you can use backward compatible SUMMA, Sundials, BMI, or Actors and you don't have to download any of the extra libraries if you don't plan to use them.

[andywood]

Congrats on getting this massive PR in. Also kudos for making it backward compatible! It does seem large to point of being non-reviewable (183 changed files, ~800 commits). This is something we avoided in the past (I remember one of my PRs being split up and resubmitted by others). I wonder if we're still trying to have a protocol around PR submission with SUMMA.

[wknoben]

Would it make sense to first clean up existing PRs and build a new v3.x release on master from those, and next use this PR as the basis for SUMMA v4?

[seantrim]

Following a typical semantic versioning approach, I think it could make sense to advance to version 4, because this probably counts as a major update (e.g., actors being available seems to be a new major feature). I agree that the review process could be challenging given the length (although the end result would definitely be worthwhile). @ashleymedin, would it be possible to subdivide things into categories like IDA, KINSOL, Actors, BMI, etc (as suggested by @andywood )? That may help speed up the PR process. Also, I'm willing to assist if needed :)

[ashleymedin]

For splitting things up -- dividing out KINSOL, Actors, and BMI wouldn't do a ton. The biggest split would be to do the old summa (num_method = itertive or numrec) and then everything else. I suppose you could then add KINSOL, ACTORS, BMI ... you could do that. The spaces and bug fixes would be in the itertive pull. The major change in itertive Summa is that it has the ability to do the inner solver step more times so we can get BE32 for example. This is with parameter be_steps (=32.0). So, that took a lot of reworking on everything from coupled_em inside to the solver. The changes I made to the flux routines are mostly only for the purposes of fixing the Jacobian, which had some incorrect terms.

[ashleymedin]

These are the changes (I was going to put together a short presentation along with some N.America comparisons of the code results on BE1 before and after). But I should probably leave this here too:

1) eval8summa.f90 calling subroutine soilCmpres with mLayerMatricHead not mLayerMatricHeadLiq version, commit ce60d55d Aug 25 2022 2) eval8summa.f90 soil min water is now theta_res, commit c19473df Apr 12 2022 3) eval8summa.f90 updateCp (update heat capacity true or false) and needCm flag, with setup to use heat capacity based on enthalpy calculation in a flag that’s not working yet, commit b69592b0 Nov 30 2022 original Summa BE equivalent to updateCp=.false. as is now needCm to match sundials [coded to .false. as equivalent to original Summa BE] 4) summaSolve.f90 return not cycle, does not effect code results just efficiency, commit 6e4cec7c Oct 26 2022 5) paramCheck.f90 check k_macropore>k_soil, makes code fail right away if not true, commit 00a10bd7 Dec 7 2022 6) run_oneHRU.f90 and run_onGRU.f90 correct error readout id’s, does not effect code results, commit d3904b51 Dec 6 2022 7) Jacobian fixes, computed in various subroutines and then included in computJacob.f90: Fix bulk heat capacity depends on frac ice/liq if updateCp, changed so updates every step (did not used to, and gives a different solution) Fix thermal conductivity at snow soil layer interfaces depends on frac ice/liq (ssdNrgFlux) if updateCp Fix soil layer and aquifer transpiration depends on canopy nrg and wat (canopy transpiration) Fix scalarCanopyLiq derivatives were getting zeroed out in calculations Fix aquifer recharge depends on soil drainage from interface above Fix soil infiltration at surface depends on all layers below and above water and temp 8) throughout, made “indian bread” terminology for NaN say it’s not a number for advised clarity, and there might be other other spaces and comments changed (e.g. tabs deleted and comments deleted or clarified), does not effect code results 9) flxMapping.f90 flux mapping of soil resistance as an energy variable corrected (was missing), commit 315583df June 5, 2023 10) runOneGRU.f90 fixed basin aquifer recharge was summing incorrectly the HRU soil drainage instead of the HRU aquifer recharge, commit cd6f07f1 June 20, 2023 11) read_icond.f90 initial values, canopy water only initialized to be positive at the start of the simulation, commit c0f7fa26 Jan 30, 2023 12) summa_restart.f90 local column aquifer storage initialized at 0 instead of 1 so there is not a massive runoff peak at simulation start, commit cd6f07f1 June 20, 2023 13) SUNDIALS IDA and KINSOL, ACTORS, and BMI/NGEN -- adds new files and some compiler directives in code num_method [numrec or kinsol or ida] ! (07) choice of numerical method Choice 'itertive' is backwards compatible to numrec 14) Working on enthalpy conservation formulation, this work is not yet complete but is in there (and runs) The choice is in the howHeatCap, with backwards compatibility 'closedForm' which will always run if you choose 'itertive' howHeatCap [closedForm or enthalpyFD] ! (30) 15) build with cmake now 16) added possible parameters for BEXX and SUNDIALS, commit 0619a403 May 31, 2023 be_steps | 1.0000 | 1.0000 | 512.0000 relErrTol_ida | 1.0d-6 | 1.0d-10| 1.0d-3 absErrTol_ida | 1.0d-6 | 1.0d-10| 1.0d-3 This is backwards compatible to give default values if not put in. 17) ascii_util.f90 memory leak, commit 44933953 May 9, 2023 18) Took out all calculation of numerical derivatives in flux routines, commit 9e5b703 Jun 28, 202 We can do that better with Sundials and a lot of them were wrong/not completely calculated. It made some of the flux routines very long to include that. So now model decision choice fDerivMeth [analytic or numericl] ! (08) method used to calculate flux derivatives refers to whether or not you want Sundials to use the provided analytical Jacobian or a finite difference one that it calculates. (the numrec num_method choice will not have numerical derivatives as an option). commit 9e5b703, Jun 28, 2023

[ashleymedin]

I'll add to this a major change is I took out all calculation of numerical derivatives since we can do that better with Sundials and a lot of them were wrong/not completely calculated. It made some of the flux routines very long to have all that stuff in --- so now model decision choice fDerivMeth [analytic or numericl] ! (08) method used to calculate flux derivatives refers to whether or not you want Sundials to use the provided analytical Jacobian or a finite difference one that it calculates (the itertive num_method choice will not have numerical derivatives as an option). I can't imagine this will affect users since the old numerical derivatives probably would have failed.

Note, this should be the same as what Sean just put in as a pull request.

[andywood]

Does this restrict options for users who chose not to use Sundials/Prime? If Sundials slows the model then it may not be a preferred choice for some applications. Sorry, the description wasn't clear to me.

On Tue, Jun 27, 2023 at 10:05 PM ashleymedin @.***> wrote:

I'll add to this a major change is I took out all calculation of numerical derivatives since we can do that better with Sundials and a lot of them were wrong/not completely calculated. It made some of the flux routines very long to have all that stuff in --- so now model decision choice fDerivMeth [analytic or numericl] ! (08) method used to calculate flux derivatives refers to whether or not you want Sundials to use the provided analytical Jacobian or a finite difference one that it calculates (the itertive num_method choice will not have numerical derivatives as an option). I can't imagine this will affect users since the old numerical derivatives probably would have failed.

— Reply to this email directly, view it on GitHub https://github.com/CH-Earth/summa/pull/526#issuecomment-1610665828, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABIKARNPW7TPJIHV6EQ2SLTXNOURDANCNFSM6AAAAAAZP44XJ4 . You are receiving this because you were mentioned.Message ID: @.***>

[ashleymedin]

No it doesn't slow the non-sundials stuff down. You cannot use the fDerivMeth=numericl choice in non-sundials summa --- it was always there but only half implemented (and would throw a lot of warnings if you tried to use it, probably failing). Now it is just removed and will say you can't use it. Sean submitted the same pull request in the non-sundials Summa a few days ago. The only difference here is that I have made it so you can use that choice if you are using Sundials, since Sundials will calculate numerical derivatives. (But we would recommend the analytic choice in all cases as it should be faster; this is just helpful for debugging. )

[andywood]

Are there cases where the analytical choice is not possible for the non-Sundials solver? Now or in the future, eg, someone wants to add a parameterization for which an analytical derivative is not possible and a numerical one is needed? Just trying to guess use cases from the notes and it's been a while since I was working on summa code.

If it's just that a little used stub method was present in the code, even if buggy in current form, it's good to present pros/cons on whether to strip out out or leave it for future development or use.

On Wed, Jun 28, 2023 at 6:20 PM ashleymedin @.***> wrote:

No it doesn't slow the non-sundials stuff down. You cannot use the fDerivMeth=numericl choice in non-sundials summa --- it was always there but only half implemented (and would throw a lot of warnings if you tried to use it, probably failing). Now it is just removed and will say you can't use it. Sean submitted the same pull request in the non-sundials Summa a few days ago. The only difference here is that I have made it so you can use that choice if you are using Sundials, since Sundials will calculate numerical derivatives. (But we would recommend the analytic choice in all cases as it should be faster; this is just helpful for debugging. )

— Reply to this email directly, view it on GitHub https://github.com/CH-Earth/summa/pull/526#issuecomment-1612266344, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABIKAROPQIPRWM3CSXDTFLLXNTC57ANCNFSM6AAAAAAZP44XJ4 . You are receiving this because you were mentioned.Message ID: @.***>

[ashleymedin]

I would say no, or maybe the more exact thing to say is it could be possible, but the old calculations of numerical derivatives in the flux routines wouldn't have solved the scenario you describe. If we want numerical derivatives in the old Summa code, I believe we would need to write a separate routine/module that calculated them with a small perturbation to the state variables, instead of the way we were doing it in each routine with only some derivatives. For example, a lot of the cross derivatives were missed in the old way (that has been stripped out). It would be better to do it outside the flux routines in a more systematic way. I don't think the past way of doing it would help debug anything.

Are there cases where the analytical choice is not possible for the non-Sundials solver? Now or in the future, eg, someone wants to add a parameterization for which an analytical derivative is not possible and a numerical one is needed? Just trying to guess use cases from the notes and it's been a while since I was working on summa code. If it's just that a little used stub method was present in the code, even if buggy in current form, it's good to present pros/cons on whether to strip out out or leave it for future development or use. … On Wed, Jun 28, 2023 at 6:20 PM ashleymedin @.> wrote: No it doesn't slow the non-sundials stuff down. You cannot use the fDerivMeth=numericl choice in non-sundials summa --- it was always there but only half implemented (and would throw a lot of warnings if you tried to use it, probably failing). Now it is just removed and will say you can't use it. Sean submitted the same pull request in the non-sundials Summa a few days ago. The only difference here is that I have made it so you can use that choice if you are using Sundials, since Sundials will calculate numerical derivatives. (But we would recommend the analytic choice in all cases as it should be faster; this is just helpful for debugging. ) — Reply to this email directly, view it on GitHub <#526 (comment)>, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABIKAROPQIPRWM3CSXDTFLLXNTC57ANCNFSM6AAAAAAZP44XJ4 . You are receiving this because you were mentioned.Message ID: @.>

[ashleymedin]

This pull request has been moved to merge into branch develop_sundials, so I'm closing it. @wknoben