WIP: Add Lobatto IIIa-b and Radau II tableaus 2.0

[axla-io]

(New PR because, in the previous one, version history got messed up)
This PR adds Lobatto and Radau II RK methods to calculate the collocation loss for a BVP

Checklist:

• Add Lobatto IIIa tables
• Add Lobatto IIIb tables
• Add Radau II tables
• Implement RK collocation
• Adaptivity (Interpolated tables)
• Sparse Jacobian for RK collocation
• Add tests
• Add docstrings

Questions / comments

1. I renamed MIRKCache into an RKCache. I also added a new RKTableau struct for holding standard RK tableaus. Then I added a new collocation function that gets dispatched on the type of Tableau that is in the cache. Idk if this is the best way to set up dispatching, but it is pretty concise.
2. For a simple pendulum problem, the Lobatto IIIb collocation seems to match the MIRK results when it converges, but frequently runs into NaN solutions. I've tested the tableau values against RungeKutta.jl. The collocation implementation itself is based on the one for MIRKs with a small adaptation to make it work with a full Butcher's tableau. The NaN solution problem is not mitigated by decreasing dt, as far as I've checked (in the range 0.1 to 0.001).
3. Assuming that we can use the same machinery for adaptive stepping as for MIRK, the interpolation tables for the RK methods need to be derived. How can I do that?
4. Sparse Jacobian as previously implemented for MIRKs doesn't work. The method needs to be modified. Currently I'm using an autodiff Jacobian (the default in NonLinearSolve) so solves are slow.

[axla-io]

Added dispatch for nested solve, which currently is very slow. Also, the user should be able to pass arguments to the nested solve. I'd like some input on better performance for that one and also on the dispatch for nesting.

Started work on the sparse jacobian for the unnested version, but band width needs to be controlled.

[ChrisRackauckas]

For the nested version, are you using init and solve! so that the caches don't need to be reconstructed?

[axla-io]

No, will do so! (along with making the function inplace).

[avik-pal]

@axla-io do you have a cleaner version of this without the merge conflicts?

[axla-io]

@axla-io do you have a cleaner version of this without the merge conflicts?

@avik-pal working on it still. Sorry for the delays. I think that I'll be able to take more time this week and turn the PR around.

[avik-pal]

@axla-io do you have a cleaner version of this without the merge conflicts?

@avik-pal working on it still. Sorry for the delays. I think that I'll be able to take more time this week and turn the PR around.

Sounds good. I would still recommend opening a new PR, even if you are able to resolve the merge conflicts. Merging after resolving the conflicts, will mess up the commit history

[axla-io]

Thanks @avik-pal , speaking of the conflicts, it seems that a lot of them come from dispatching on the AbstractRK class. Should I remove that and to something like Union{MIRK, FIRK} instead?

[axla-io]

I'm a bit worried about the convergence tests.

For example LobattoIIIa I get this behavior:
For nested, the convergence is bad, below 1 for all orders
For non nested, the convergence is faster or equal than the method order (IIIa2 is correct, order 4 convergence for IIIa3, 6 for IIIa4, 11 for IIIa5)

The rest of the Lobatto and Radau methods behave similarily

The other tests look good.

[avik-pal]

can you check the states for each case, it is possible some reinit path is missing for the cached forwarddiff bindings which causes the problem?

One other way would be to use non-cached nested version and see if there is a difference between that and the cached version.

[axla-io]

can you check the states for each case, it is possible some reinit path is missing for the cached forwarddiff bindings which causes the problem?

Yes this is pretty much the problem as I remembered it.

One other way would be to use non-cached nested version and see if there is a difference between that and the cached version.

That's a good idea! I'll let you know how it goes

[axla-io]

I totally forgot that Lobatto and Radau methods are of order 2s-1, which fits my results, so aside from nested, we are good

[axla-io]

Okay so right now I need to fix a bug I noticed in the sparse differentiation for the expanded cache. Moving on to nested after that.
Other than that, the expanded cache passes all tests except for Problem 4 and 8