Why both tableau and explicit solver in DifferentialEquations.jl?

Question

I am looking through the DifferentialEquations.jl package. In DiffEqDevTools/src/ode_tableaus.jl I can see the tableaux for Midpoint and RK4.

But I can also see explicit code for these schemes in OrdinaryDiffEq/src/integrators/fixed_timestep_integrators.jl.

I sort of expected some code to use the tableaux rather than there being an explicit solver.

I am not sure how to check whether the tableau is being used. I tried removing OrdinaryDiffEq.jl but then my example would not run.

This rather suggests that the explicit code is being used. But in that case why does the tableau exist at all?

Answer 1

Yes, in most cases the tableaus are not used. In fact, the tableaus are only used if you use the ExplicitRK(tableau=...) method. You can see the convergence tests on each of these in DiffEqDevTools' tests, but other than that they generally aren't used.

The reason for this is because tableau-based implementations don't tend to be as efficient. Getting rid of the indirection caused by using pointers to the constants has a measurable effect on the runtime. Sure, in the asymptotic limit where the user's f is extremely costly the implementation details won't matter, but most of the real world cases aren't in that limit as evidenced by real benchmarks (in that limit, you should be using multistep methods aynways). So there are hard-coded versions of the most efficient Runge-Kutta methods with hard-coded higher order interpolation schemes as these are the workhorse methods and should get maximum efficiency.

So then why do these tableaus exist? I think it's important to note that DifferentialEquations.jl is not just a software package for using differential equations, but it's also a software package for the development and testing of new methods. This is evidenced by the testing features in the devdocs. For algorithms which have a more efficient implementation, the tableaus still have a development use because the tableaus all have the same implementation, and thus this gives an easy scientific way to determine the true efficiency between methods. Not only does having the tableaus allow comparing efficiency, but then you can also compare stability regions using the plot recipes:

plot(constructRK4())

This large tableau library was used to comb through all of the RK methods and create modernized choices. I posted some haphazard notes written about it, and documented some parts more extensively in a CompSci SO post. These were all experiments done using the development tools.

In the end, DifferentialEquations.jl is really unique because it's not just a re-implementation of what you've seen before. One clear change from before is that the order 4/5 Runge-Kutta method of choice is not the Dormand-Prince pair DP5 like in MATLAB or SciPy (which is just a dopri5 wrapper), but a modern algorithm: the Tsit5() method. This is because this more recent method theoretically achieves lower error when computational cost is even, and the devtools confirmed this. Other things in DifferentialEquations.jl which are unique are its adaptivity for stochastic differential equations, high-order methods for delay differential equations, etc., and there's lots more research coming (things in private repositories until publication). Most of this was made possible (or at least easy to do) because of the associated development suite.

So I think this shows that clearly the philosophy of DifferentialEquations.jl is not like differential equations suites in other languages/libraries. In other suites, like MATLAB or SciPy, most of the time the goal is to give you a few basic methods which are generally useful. Efficiency is not necessarily the goal (a good example of this is Shampine's conscious choice to not have high-order RK methods in MATLAB), and wrapping "standard" implementation is usually good enough (SciPy is only wrappers). Simplicity and standards make sense for these software suites.

But DifferentialEquations.jl's focus is on the development of new methods for handling modern computationally difficult equations, and using these new methods to solve problems which were previously infeasible. Because of this different focus, there are some choices that are made different. For example, having large pools of methods/implementations is good here because it allows users and methods researchers to compare between algorithms and find out how to continually improve them and the choices (I invite the users to test out tableau methods on their problems and see if some obscure RK method does well). I care about making sure the most efficient research tools are available, and handle the complexity through giving defaults and recommendations to users. If there's any part of the recommendations in the documentation which are not clear, please open an issue so we can discuss how to improve the documentation. But I see leading users to "correct" method choices as a documentation issue, not an issue which should constrain the design, capabilities, or efficiency.

I'm going to cover the philosophy of DifferentialEquations.jl in some near future blog posts, but that's the gist of it. I don't think this is much of a SO question because it's more of a personal question to me for why I keep known inefficient methods around (as not recommended though!), and I hope this gives an informative answer.