Rename some of the time steppers

[wthrowe]

Proposed changes

This renames

• Cerk2 to Heun as the standard name for the method,
• RungeKutta4 to ClassicalRungeKutta4 to distinguish it from Cerk4, and
• AdamsBashforthN to AdamsBashforth because the N always seemed unnecessary and this seemed like a good opportunity.

Of the remaining classes, DormandPrince5 seems fine, and the others I don't know better names for:

• RungeKutta3 is based on a strong-stability-preserving RK3 method in Hesthaven, but I think I broke the SSP property in Stop storing old values in time-stepper history #3253 and even before that the requirements to obtain SSP were undocumented (step size must be below about 0.8 of the stable step).
• Cerk3, Cerk4, Cerk5 are from https://epubs.siam.org/doi/10.1137/0913084 derived by minimizing some sort of truncation error, but my knowledge of RK theory isn't yet good enough to understand the details. That paper and https://www.sciencedirect.com/science/article/pii/S0021999110005802 (both cited by the classes) disagree on what CERK stands for, but as far as I can tell it just means an explicit RK method with dense output, which describes all our time steppers except AB.

Note to reviewers: The changes that are not straight find-replace/formatting are the removal of two unnecessary "\file" comments and rewording the help string for Heun.

Upgrade instructions

Some of the time steppers have been renamed. Update your input files accordingly:

• AdamsBashforthN -> AdamsBashforth
• Cerk2 -> Heun
• RungeKutta4 -> ClassicalRungeKutta4

Code review checklist

• The code is documented and the documentation renders correctly. Run
  make doc to generate the documentation locally into BUILD_DIR/docs/html.
  Then open index.html.
• The code follows the stylistic and code quality guidelines listed in the
  code review guide.
• The PR lists upgrade instructions and is labeled bugfix or
  new feature if appropriate.

Further comments

[wthrowe]

How do I fix Unit.Options.Python.ExtractInputSourceYAMLFromH5? I assume I have to regenerate the h5 file it's reading, but I don't see instructions for how to do that.

[wthrowe]

@geoffrey4444 ^

[knelli2]

@wthrowe ?

import h5py
h5file = h5py.File("SurfaceTestData0.h5", "a")
input_file_str = h5file.attrs['InputSource.yaml'][0]
h5file.attrs['InputSource.yaml'] = [input_file_str.replace("AdamsBashforthN", "AdamsBashforth")]
h5file.close()

[wthrowe]

Reading a bit more on the CERK terminology:

• Gassner et al. (2011) define CERK as "Continuous Extension Runge-Kutta" and refer to two Owren and Zennaro papers (1991 and 1992) for details. This is an application of CERK methods, and it is not immediately clear to me what order dense output they require.
• Owren and Zennaro (1991 and 1992) define CERK as "continuous, explicit Runge-Kutta" and define the order of the method as the local truncation order of the interpolant (u(t0+x dt) - y0(t0+x dt) = O(dt^{p+1}) with u the interpolant and y0 a solution of the ODE with y0(t0) = u(t0)), but comment that "An important issue in many of these [referenced] papers has been whether or not the continuous approximation should yield the same order of consistency in the interior of the step as at the end-points."
• Hairer et al. (Solving Ordinary Differential Equations I) never use anything that acronyms to CERK, but define a continuous Runge-Kutta method as a Runge-Kutta method along with a dense output polynomial. (They do refer to a "continuous extension" of an RK method in a few places, but to my reading this is just standard English use rather than a technical term.) Such a method is defined to have order p^* as the local truncation order of the dense output (u(t0+x dt) - y0(t0+x dt) = O(dt^{p^*+1})), which is considered distinct from the order p of the discrete method (u1 - y0(t0 + dt) = O(dt^{p+1}), where u1 is the result of the RK step, equal to u(t0 + dt) if the dense output is continuous). They then consider p^* = p, p^* = p - 1, and in a few cases smaller values.
• Skimming a few random papers referenced by Owren as examples of CERK methods, the only one that was straightforward enough that I could quickly understand the definitions was Dormand and Prince (1986). They take the desired accuracy of dense output as the global truncation error of the discrete method (u(t) - y(t) = O(dt^p) with y(0) = u(0) with the full integration starting at 0). This is equivalent to p^* = p - 1 in the notation above. (It is not clear what they consider the order of the dense output itself to be.) They never use the term "continuous", but just talk about dense output.

So I see a definition of CERK requiring the truncation error to be at least the local truncation error (Orwen), an example of a CERK method with it the global truncation error (DP), and a definition of "continuous" placing no specific requirements on accuracy (Hairer). In light of this, I think it would be good to somehow explicitly include the higher-order property in the names of the classes, but I'm not coming up with anything reasonable.

[nilsdeppe]

Ooof, that's confusing... What if we do something like RkXDoY where X gives the truncation error at the end of the step and Y of the dense output. We could reasonably add Continuous to the front and document the TimeSteppers namespace that by this we mean it does not require additional RHS evaluations to get to the specified dense output. This probably still leaves a large possible overload set because of the many different ways of choosing the coefficients but I'm also not coming up with anything elegant

[wthrowe]

Here's what I'm going with:

• Methods with well established names get those names: AdamsBashforth, ClassicalRungeKutta4, DormandPrince5, Heun. I may change Heun to Heun2. Haven't decided.
• Other Runge-Kutta methods are named Rk${N}${Ref}${Prop}, where ${N} is the order of the method, ${Ref} is a reference for the method, and ${Prop} is an optional important property or disambiguation. This gives us Rk3Hesthaven (potentially Rk3HesthavenSsp if that were still true), and Rk${N}Owren for N in {3, 4, 5}. In the case that a method takes options that determine its order, N can be omitted.
• In the theoretical case where there is no clear reference to cite, either because of a well-known method without a standard name or because we derived the method ourselves to have some particular property, the Ref above can be left out, but Prop is then mandatory.
• Non-Runge-Kutta methods can be handled later if necessary, but could presumably follow similar rules with "Rk" replaced with a different family. The only potential case I have thought about would be reasonably named AdamsPredictorCorrector or AdamsMoultonPredictorCorrector. (Aside: My understanding is that, unlike in the Adams-Bashforth case, which is entirely due to Adams, this method actually is due to Moulton, but resulted in his name being attached to a different "Adams-Moulton" method, which is solely due to Adams. Naming is complicated.)
• IMEX methods are more complicated and I haven't thought them out yet.

It is still unclear to me how important dense-output order is. (I think I'd need to go through the Gassner reference in detail.) The methods that have dense output of lower order than the method order are ClassicalRungeKutta4, DormandPrince5, and Rk3Hesthaven, and the last of those has a third-order formula we just haven't coded.

As for using Continuous as a disambiguator for not needing additional evaluations for dense output:

• Hesthaven gives it a different meaning,
• I think all our current methods have that property,
• We don't support optimizing based on that anyway so it's not really important to us.

(I'm also a bit confused by interest in this property, since it is not obviously desirable. As far as I can tell, those methods avoid requiring extra evaluations by doing more evaluations for the non-dense output, rather than fewer for the dense output.)

I will also add something to the docs for Rk3Hesthaven explaining the whole SSP situation.

[nilsdeppe]

That all sounds good to me. A couple thoughts:

• We should have some way of doing SSP since this is actually sometimes important for hydro. SSP allows high-order while being able to do things like guarantee positivity with flux limiters. While I doubt we will need that in practice (we could drop time step size or something to achieve a similar effect), it would be useful to be able to show by testing that we don't need SSP in a situation where the code is crashing. It would also be good for code comparison since SSP is very common, especially in the disk community.
• The CERK methods are used for providing a high-order initial guess for the ADER-DG (or CG) "implicit" solve in time. Since you need to interpolate to Gauss (-Lobatto) nodes in time, you always need dense output. Using a CERK method of one order lower than the corrector allows the solver to converge in one iteratio. I believe this is provable in the linear equation case because ADER-DG is a fixed point iteration scheme that is also guaranteed to converge in N iterations for an N-th order method.

[wthrowe]

Updated with all the new names.

For RungeKutta3, I went with Rk3HesthavenSsp because I have a local fix for the problem that I will submit a few PRs down the line and renaming twice is just extra churn. It's not any worse than the current situation.