hydro solver default

[rgknox]

The default solver method for plant-soil hydraulics is the Taylor 1D method. This is an explicit solve (ie first order taylor expansion, synonymous with an euler step) that breaks up the solution so that fluxes are are calculated independently for each root layer. For each solution, the plant compartments are partitioned as if fluxes through them are in parallel with other root-layer solves. So the final fluxes through the compartment is a sum of the indepenent solves including a weighting factor.

Yilin (@pnlfang) implemented the Picard solution. I know less about it though, but I do know that it solves the all root layers and plant compartments simultaneously.

There is also the Newton's method solve, which uses Newton's method to iteratively find an implicit solution that includes all root layers simultaneously.

In my recent tests from PR 1164, the Picard method was more stable. See thread here: #1164 (comment)

Can folks chime in with their preferred method?

[JunyanDing]

@rgknox That's great!
I am not an expertise on these numeric things, but it seems to me the Picard solution is the most advanced. Does it take longer time to find the solution? Let's keep all the options, but use Picard as the default.

[jennykowalcz]

I vote for Picard as well since I am also finding that PR 1164 with fates_hydro_solver=2 resolves the crashes I was getting in dry regions, but with that PR and fates_hydro_solver=1 I still get "Could not find a stable solution for hydro 1D solve" errors in dry gridcells.

[bchristo]

@rgknox agree this is great news.
Following Junyan's suggestion, I too vote for Picard as default, but keeping all options available -- with a some helpful comments/caveats about conditions in which each may encounter trouble. As the one who initially designed the 1D Taylor solve, it was intended as an initial efficient yet imperfect solve for hydro. My hope always was that having others would come along, such as Yilin, to implement more robust solvers and sounds like Picard now fits the bill.

[rgknox]

@bchristo , yup, we will be sure to keep Taylor 1D as an option.

[jennykowalcz]

I was reminded of this recently while running a test case -- running FATES-Hydro "out of the box" with API 36.0 and the default parameter file (which has fates_hydro_solver=1) crashes with a "Could not find a stable solution for hydro 1D solve" error. Setting fates_hydro_solver=2 allows FATES-Hydro to run successfully out of the box. So, I suggest that we set fates_hydro_solver=2 as the default in the FATES parameter file, but keep Taylor 1D as an option. If it would be helpful, I could incorporate this into the parameter file changes I'm making to address issue 1128

[rgknox]

I agree! That's a great idea @jennykowalcz . I think @pnlfang , who wrote both alternative solver's code, advocates for the Picard solver. Yilin?

[pnlfang]

Agree. :-)