Creating B-series

We have already seen some ways of creating B-series in the basic tutorial. In this tutorial, we look at ways to obtain the B-series of different time integration methods.

B-series for Runge-Kutta methods

BSeries.jl and RootedTrees.jl provide the type RungeKuttaMethod as wrapper of Butcher coefficients A, b, c of Runge-Kutta methods. For example, you can create the classical explicit, fourth-order Runge-Kutta method as follows.

using BSeries

# classical RK4
A = [0 0 0 0
      1//2 0 0 0
      0 1//2 0 0
      0 0 1 0]
b = [1 // 6, 1 // 3, 1 // 3, 1 // 6]
rk = RungeKuttaMethod(A, b)
RungeKuttaMethod{Rational{Int64}} with
A: 4×4 Matrix{Rational{Int64}}:
  0     0    0  0
 1//2   0    0  0
  0    1//2  0  0
  0     0    1  0
b: 4-element Vector{Rational{Int64}}:
 1//6
 1//3
 1//3
 1//6
c: 4-element Vector{Rational{Int64}}:
  0
 1//2
 1//2
  1

Instead of passing the Butcher coefficients explicitly to bseries, you can also pass the wrapper struct rk. In fact, this is the preferred method.

series = bseries(rk, 5)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 18 entries:
  RootedTree{Int64}: Int64[]         => 1
  RootedTree{Int64}: [1]             => 1
  RootedTree{Int64}: [1, 2]          => 1//2
  RootedTree{Int64}: [1, 2, 3]       => 1//6
  RootedTree{Int64}: [1, 2, 2]       => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]    => 1//24
  RootedTree{Int64}: [1, 2, 3, 3]    => 1//12
  RootedTree{Int64}: [1, 2, 3, 2]    => 1//8
  RootedTree{Int64}: [1, 2, 2, 2]    => 1//4
  RootedTree{Int64}: [1, 2, 3, 4, 5] => 0
  RootedTree{Int64}: [1, 2, 3, 4, 4] => 1//48
  RootedTree{Int64}: [1, 2, 3, 4, 3] => 1//48
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 1//24
  RootedTree{Int64}: [1, 2, 3, 3, 3] => 1//24
  RootedTree{Int64}: [1, 2, 3, 3, 2] => 1//16
  RootedTree{Int64}: [1, 2, 3, 2, 3] => 1//16
  RootedTree{Int64}: [1, 2, 3, 2, 2] => 5//48
  RootedTree{Int64}: [1, 2, 2, 2, 2] => 5//24

We can check that the classical Runge-Kutta method is indeed fourth-order accurate.

series - ExactSolution(series)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 18 entries:
  RootedTree{Int64}: Int64[]         => 0
  RootedTree{Int64}: [1]             => 0
  RootedTree{Int64}: [1, 2]          => 0
  RootedTree{Int64}: [1, 2, 3]       => 0
  RootedTree{Int64}: [1, 2, 2]       => 0
  RootedTree{Int64}: [1, 2, 3, 4]    => 0
  RootedTree{Int64}: [1, 2, 3, 3]    => 0
  RootedTree{Int64}: [1, 2, 3, 2]    => 0
  RootedTree{Int64}: [1, 2, 2, 2]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 5] => -1//120
  RootedTree{Int64}: [1, 2, 3, 4, 4] => 1//240
  RootedTree{Int64}: [1, 2, 3, 4, 3] => -1//240
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 1//120
  RootedTree{Int64}: [1, 2, 3, 3, 3] => -1//120
  RootedTree{Int64}: [1, 2, 3, 3, 2] => -1//240
  RootedTree{Int64}: [1, 2, 3, 2, 3] => 1//80
  RootedTree{Int64}: [1, 2, 3, 2, 2] => 1//240
  RootedTree{Int64}: [1, 2, 2, 2, 2] => 1//120
order_of_accuracy(series)
4

B-series for additive Runge-Kutta methods

BSeries.jl and RootedTrees.jl also support additive Runge-Kutta methods via the wrapper AdditiveRungeKuttaMethod. For example, we can write the Störmer-Verlet method as additive Runge-Kutta method following Table II.2.1 of Hairer, Lubich, and Wanner (2002).

using BSeries

As = [
    [0 0; 1//2 1//2],
    [1//2 0; 1//2 0],
]
bs = [
    [1 // 2, 1 // 2],
    [1 // 2, 1 // 2],
]
ark = AdditiveRungeKuttaMethod(As, bs)
AdditiveRungeKuttaMethod{Rational{Int64}} with methods
1. RungeKuttaMethod{Rational{Int64}} with
A: 2×2 Matrix{Rational{Int64}}:
  0     0
 1//2  1//2
b: 2-element Vector{Rational{Int64}}:
 1//2
 1//2
c: 2-element Vector{Rational{Int64}}:
 0
 1
2. RungeKuttaMethod{Rational{Int64}} with
A: 2×2 Matrix{Rational{Int64}}:
 1//2  0
 1//2  0
b: 2-element Vector{Rational{Int64}}:
 1//2
 1//2
c: 2-element Vector{Rational{Int64}}:
 1//2
 1//2

We can create the B-series as usual, truncated to order 3.

series = bseries(ark, 3)
TruncatedBSeries{BicoloredRootedTree{Int64, Vector{Int64}, Vector{Bool}}, Rational{Int64}} with 21 entries:
  ColoredRootedTree{Int64}: (Int64[], Bool[])          => 1
  ColoredRootedTree{Int64}: ([1], Bool[0])             => 1
  ColoredRootedTree{Int64}: ([1], Bool[1])             => 1
  ColoredRootedTree{Int64}: ([1, 2], Bool[0, 0])       => 1//2
  ColoredRootedTree{Int64}: ([1, 2], Bool[1, 0])       => 1//2
  ColoredRootedTree{Int64}: ([1, 2], Bool[0, 1])       => 1//2
  ColoredRootedTree{Int64}: ([1, 2], Bool[1, 1])       => 1//2
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 0, 0]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 0, 0]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 1, 0]) => 0
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 1, 0]) => 0
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 0, 1]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 0, 1]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 1, 1]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 1, 1]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[0, 0, 0]) => 1//2
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[1, 0, 0]) => 1//2
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[0, 1, 0]) => 1//4
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[1, 1, 0]) => 1//4
  ⋮                                                    => ⋮

For additive Runge-Kutta methods like this, we use colored rooted trees. Again, we can check the order of accuracy by comparing the coefficients to the exact solution.

series - ExactSolution(series)
TruncatedBSeries{BicoloredRootedTree{Int64, Vector{Int64}, Vector{Bool}}, Rational{Int64}} with 21 entries:
  ColoredRootedTree{Int64}: (Int64[], Bool[])          => 0
  ColoredRootedTree{Int64}: ([1], Bool[0])             => 0
  ColoredRootedTree{Int64}: ([1], Bool[1])             => 0
  ColoredRootedTree{Int64}: ([1, 2], Bool[0, 0])       => 0
  ColoredRootedTree{Int64}: ([1, 2], Bool[1, 0])       => 0
  ColoredRootedTree{Int64}: ([1, 2], Bool[0, 1])       => 0
  ColoredRootedTree{Int64}: ([1, 2], Bool[1, 1])       => 0
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 0, 0]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 0, 0]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 1, 0]) => -1//6
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 1, 0]) => -1//6
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 0, 1]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 0, 1]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[0, 1, 1]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 3], Bool[1, 1, 1]) => 1//12
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[0, 0, 0]) => 1//6
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[1, 0, 0]) => 1//6
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[0, 1, 0]) => -1//12
  ColoredRootedTree{Int64}: ([1, 2, 2], Bool[1, 1, 0]) => -1//12
  ⋮                                                    => ⋮
order_of_accuracy(series)
2

We can also create LaTeX code for this B-series as follows.

using Latexify
latexify(series)

\[F_{f}\mathopen{}\left( \varnothing \right)\mathclose{} + h \cdot F_{f}\mathopen{}\left( \rootedtree[.] \right)\mathclose{} + h \cdot F_{f}\mathopen{}\left( \rootedtree[o] \right)\mathclose{} + \frac{1}{2} \cdot h^{2} \cdot F_{f}\mathopen{}\left( \rootedtree[.[.]] \right)\mathclose{} + \frac{1}{2} \cdot h^{2} \cdot F_{f}\mathopen{}\left( \rootedtree[o[.]] \right)\mathclose{} + \frac{1}{2} \cdot h^{2} \cdot F_{f}\mathopen{}\left( \rootedtree[.[o]] \right)\mathclose{} + \frac{1}{2} \cdot h^{2} \cdot F_{f}\mathopen{}\left( \rootedtree[o[o]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[.[.]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[.[.]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[.[o]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[.[o]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[o[o]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[o[o]]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[.][.]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[.][.]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[o][.]] \right)\mathclose{} + \frac{1}{4} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[o][.]] \right)\mathclose{} + \frac{1}{8} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[.[o][o]] \right)\mathclose{} + \frac{1}{8} \cdot h^{3} \cdot F_{f}\mathopen{}\left( \rootedtree[o[o][o]] \right)\mathclose{}\]

When compiled using the preamble code shown in the docstring of RootedTrees.latexify, the output looks as follows.

bseries-SV

References

B-series for Rosenbrock methods

BSeries.jl and RootedTrees.jl also support Rosenbrock (Rosenbrock-Wanner, ROW) methods via the wrapper RosenbrockMethod. For example, a classical ROW method of Kaps and Rentrop (1979) can be parameterized as follows.

using BSeries

γ = [0.395 0 0 0;
     -0.767672395484 0.395 0 0;
     -0.851675323742 0.522967289188 0.395 0;
     0.288463109545 0.880214273381e-1 -0.337389840627 0.395]
A = [0 0 0 0;
     0.438 0 0 0;
     0.796920457938 0.730795420615e-1 0 0;
     0.796920457938 0.730795420615e-1 0 0]
b = [0.199293275701, 0.482645235674, 0.680614886256e-1, 0.25]
ros = RosenbrockMethod(γ, A, b)
RosenbrockMethod{Float64} with
γ: 4×4 Matrix{Float64}:
  0.395     0.0         0.0      0.0
 -0.767672  0.395       0.0      0.0
 -0.851675  0.522967    0.395    0.0
  0.288463  0.0880214  -0.33739  0.395
A: 4×4 Matrix{Float64}:
 0.0      0.0        0.0  0.0
 0.438    0.0        0.0  0.0
 0.79692  0.0730795  0.0  0.0
 0.79692  0.0730795  0.0  0.0
b: 4-element Vector{Float64}:
 0.199293275701
 0.482645235674
 0.0680614886256
 0.25
c: 4-element Vector{Float64}:
 0.0
 0.438
 0.8699999999995001
 0.8699999999995001

We can create the B-series as usual, truncated to order 5.

series = bseries(ros, 5)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Float64} with 18 entries:
  RootedTree{Int64}: Int64[]         => 1.0
  RootedTree{Int64}: [1]             => 1.0
  RootedTree{Int64}: [1, 2]          => 0.5
  RootedTree{Int64}: [1, 2, 3]       => 0.166667
  RootedTree{Int64}: [1, 2, 2]       => 0.333333
  RootedTree{Int64}: [1, 2, 3, 4]    => 0.0416667
  RootedTree{Int64}: [1, 2, 3, 3]    => 0.0833333
  RootedTree{Int64}: [1, 2, 3, 2]    => 0.125
  RootedTree{Int64}: [1, 2, 2, 2]    => 0.25
  RootedTree{Int64}: [1, 2, 3, 4, 5] => 0.00872428
  RootedTree{Int64}: [1, 2, 3, 4, 4] => 0.00418004
  RootedTree{Int64}: [1, 2, 3, 4, 3] => 0.0320513
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 0.0467417
  RootedTree{Int64}: [1, 2, 3, 3, 3] => 0.05
  RootedTree{Int64}: [1, 2, 3, 3, 2] => 0.00387949
  RootedTree{Int64}: [1, 2, 3, 2, 3] => 0.0469263
  RootedTree{Int64}: [1, 2, 3, 2, 2] => 0.09295
  RootedTree{Int64}: [1, 2, 2, 2, 2] => 0.19998

Again, we can check the order of accuracy by comparing the coefficients to the exact solution.

series - ExactSolution(series)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Float64} with 18 entries:
  RootedTree{Int64}: Int64[]         => 0.0
  RootedTree{Int64}: [1]             => 5.99965e-13
  RootedTree{Int64}: [1, 2]          => -1.52767e-13
  RootedTree{Int64}: [1, 2, 3]       => -2.12941e-13
  RootedTree{Int64}: [1, 2, 2]       => -2.50522e-13
  RootedTree{Int64}: [1, 2, 3, 4]    => -1.37459e-13
  RootedTree{Int64}: [1, 2, 3, 3]    => -8.44186e-14
  RootedTree{Int64}: [1, 2, 3, 2]    => -8.9137e-14
  RootedTree{Int64}: [1, 2, 2, 2]    => -3.5999e-13
  RootedTree{Int64}: [1, 2, 3, 4, 5] => 0.000390949
  RootedTree{Int64}: [1, 2, 3, 4, 4] => -0.0124866
  RootedTree{Int64}: [1, 2, 3, 4, 3] => 0.00705128
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 0.0134083
  RootedTree{Int64}: [1, 2, 3, 3, 3] => -9.4473e-14
  RootedTree{Int64}: [1, 2, 3, 3, 2] => -0.0627872
  RootedTree{Int64}: [1, 2, 3, 2, 3] => -0.00307372
  RootedTree{Int64}: [1, 2, 3, 2, 2] => -0.00705
  RootedTree{Int64}: [1, 2, 2, 2, 2] => -2.0e-5
order_of_accuracy(series)
4

References

  • Peter Kaps and Peter Rentrop. "Generalized Runge-Kutta methods of order four with stepsize control for stiff ordinary differential equations." Numerische Mathematik 33, no. 1 (1979): 55-68. DOI: 10.1007/BF01396495

B-series for the average vector field method

Consider the autonomous ODE

\[u'(t) = f\bigl( u(t) \bigr).\]

The average vector field method

\[u^{n+1} = u^{n} + \Delta t \int_0^1 f\bigl(\xi u^{n+1} + (1 - \xi) u^{n}\bigr) \mathrm{d} \xi\]

introduced by McLachlan, Quispel, and Robidoux (1999) is a second-order accurate method. Quispel and McLaren (2008) discovered that it is indeed a B-series method. Its coefficients are given explicitly as

\[\begin{align*} b(.) &= 1, \\ b([t_1, ..., t_n]) &= b(t_1)...b(t_n) / (n + 1). \end{align*}\]

by Celledoni, McLachlan, McLaren, Owren, Quispel, and Wright (2009). We can implement this up to order 5 in BSeries.jl as follows.

using BSeries

series = bseries(5) do t, series
    if order(t) in (0, 1)
        return 1 // 1
    else
        v = 1 // 1
        n = 0
        for subtree in SubtreeIterator(t)
            v *= series[subtree]
            n += 1
        end
        return v / (n + 1)
    end
end
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 18 entries:
  RootedTree{Int64}: Int64[]         => 1
  RootedTree{Int64}: [1]             => 1
  RootedTree{Int64}: [1, 2]          => 1//2
  RootedTree{Int64}: [1, 2, 3]       => 1//4
  RootedTree{Int64}: [1, 2, 2]       => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]    => 1//8
  RootedTree{Int64}: [1, 2, 3, 3]    => 1//6
  RootedTree{Int64}: [1, 2, 3, 2]    => 1//6
  RootedTree{Int64}: [1, 2, 2, 2]    => 1//4
  RootedTree{Int64}: [1, 2, 3, 4, 5] => 1//16
  RootedTree{Int64}: [1, 2, 3, 4, 4] => 1//12
  RootedTree{Int64}: [1, 2, 3, 4, 3] => 1//12
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 1//12
  RootedTree{Int64}: [1, 2, 3, 3, 3] => 1//8
  RootedTree{Int64}: [1, 2, 3, 3, 2] => 1//9
  RootedTree{Int64}: [1, 2, 3, 2, 3] => 1//12
  RootedTree{Int64}: [1, 2, 3, 2, 2] => 1//8
  RootedTree{Int64}: [1, 2, 2, 2, 2] => 1//5

BSeries.jl also offers a convenience constructor using the type AverageVectorFieldMethod as follows.

series == bseries(AverageVectorFieldMethod(), 5)
true

We can check that this method is second-order accurate by comparing it to the B-series of the exact solution, truncated at the same order.

series - ExactSolution(series)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 18 entries:
  RootedTree{Int64}: Int64[]         => 0
  RootedTree{Int64}: [1]             => 0
  RootedTree{Int64}: [1, 2]          => 0
  RootedTree{Int64}: [1, 2, 3]       => 1//12
  RootedTree{Int64}: [1, 2, 2]       => 0
  RootedTree{Int64}: [1, 2, 3, 4]    => 1//12
  RootedTree{Int64}: [1, 2, 3, 3]    => 1//12
  RootedTree{Int64}: [1, 2, 3, 2]    => 1//24
  RootedTree{Int64}: [1, 2, 2, 2]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 5] => 13//240
  RootedTree{Int64}: [1, 2, 3, 4, 4] => 1//15
  RootedTree{Int64}: [1, 2, 3, 4, 3] => 7//120
  RootedTree{Int64}: [1, 2, 3, 4, 2] => 1//20
  RootedTree{Int64}: [1, 2, 3, 3, 3] => 3//40
  RootedTree{Int64}: [1, 2, 3, 3, 2] => 2//45
  RootedTree{Int64}: [1, 2, 3, 2, 3] => 1//30
  RootedTree{Int64}: [1, 2, 3, 2, 2] => 1//40
  RootedTree{Int64}: [1, 2, 2, 2, 2] => 0
order_of_accuracy(series)
2

References

  • Robert I. McLachlan, G. Reinout W. Quispel, and Nicolas Robidoux. "Geometric integration using discrete gradients." Philosophical Transactions of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences 357, no. 1754 (1999): 1021-1045. DOI: 10.1098/rsta.1999.0363
  • G. Reinout W. Quispel, and David Ian McLaren. "A new class of energy-preserving numerical integration methods." Journal of Physics A: Mathematical and Theoretical 41, no. 4 (2008): 045206. DOI: 10.1088/1751-8113/41/4/045206
  • Elena Celledoni, Robert I. McLachlan, David I. McLaren, Brynjulf Owren, G. Reinout W. Quispel, and William M. Wright. "Energy-preserving Runge-Kutta methods." ESAIM: Mathematical Modelling and Numerical Analysis 43, no. 4 (2009): 645-649. DOI: 10.1051/m2an/2009020

Composing B-series

B-series representing a mapping such as a time integration method can be composed. This operation is equivalent to performing a step with one method after another. Here, we present the example of Butcher's method of effective order 5. This is a fourth-order Runge-Kutta method resulting in a fifth-order method when composed with a special starting and finishing procedure.

First, we set up the B-series of the main method (method "a" in Butcher's paper):

using BSeries

A = [0 0 0 0 0;
     1//5 0 0 0 0;
     0 2//5 0 0 0;
     3//16 0 5//16 0 0;
     1//4 0 -5//4 2 0]
b = [1 // 6, 0, 0, 2 // 3, 1 // 6]
rk_a = RungeKuttaMethod(A, b)
series_a = bseries(rk_a, 6)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 38 entries:
  RootedTree{Int64}: Int64[]            => 1
  RootedTree{Int64}: [1]                => 1
  RootedTree{Int64}: [1, 2]             => 1//2
  RootedTree{Int64}: [1, 2, 3]          => 1//6
  RootedTree{Int64}: [1, 2, 2]          => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]       => 1//24
  RootedTree{Int64}: [1, 2, 3, 3]       => 1//12
  RootedTree{Int64}: [1, 2, 3, 2]       => 1//8
  RootedTree{Int64}: [1, 2, 2, 2]       => 1//4
  RootedTree{Int64}: [1, 2, 3, 4, 5]    => 1//120
  RootedTree{Int64}: [1, 2, 3, 4, 4]    => 1//60
  RootedTree{Int64}: [1, 2, 3, 4, 3]    => 1//48
  RootedTree{Int64}: [1, 2, 3, 4, 2]    => 1//30
  RootedTree{Int64}: [1, 2, 3, 3, 3]    => 1//24
  RootedTree{Int64}: [1, 2, 3, 3, 2]    => 1//15
  RootedTree{Int64}: [1, 2, 3, 2, 3]    => 5//96
  RootedTree{Int64}: [1, 2, 3, 2, 2]    => 5//48
  RootedTree{Int64}: [1, 2, 2, 2, 2]    => 5//24
  RootedTree{Int64}: [1, 2, 3, 4, 5, 6] => 0
  ⋮                                     => ⋮

Note that this method is fourth-order accurate:

order_of_accuracy(series_a)
4

Next, we set up the starting procedure (method "b" in Butcher's paper):

A = [0 0 0 0 0;
     1//5 0 0 0 0;
     0 2//5 0 0 0;
     75//64 -9//4 117//64 0 0;
     -37//36 7//3 -3//4 4//9 0]
b = [19 // 144, 0, 25 // 48, 2 // 9, 1 // 8]
rk_b = RungeKuttaMethod(A, b)
series_b = bseries(rk_b, 6)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 38 entries:
  RootedTree{Int64}: Int64[]            => 1
  RootedTree{Int64}: [1]                => 1
  RootedTree{Int64}: [1, 2]             => 1//2
  RootedTree{Int64}: [1, 2, 3]          => 1//6
  RootedTree{Int64}: [1, 2, 2]          => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]       => 13//320
  RootedTree{Int64}: [1, 2, 3, 3]       => 13//160
  RootedTree{Int64}: [1, 2, 3, 2]       => 121//960
  RootedTree{Int64}: [1, 2, 2, 2]       => 121//480
  RootedTree{Int64}: [1, 2, 3, 4, 5]    => 13//1600
  RootedTree{Int64}: [1, 2, 3, 4, 4]    => 13//800
  RootedTree{Int64}: [1, 2, 3, 4, 3]    => 139//6400
  RootedTree{Int64}: [1, 2, 3, 4, 2]    => 13//400
  RootedTree{Int64}: [1, 2, 3, 3, 3]    => 139//3200
  RootedTree{Int64}: [1, 2, 3, 3, 2]    => 13//200
  RootedTree{Int64}: [1, 2, 3, 2, 3]    => 2003//38400
  RootedTree{Int64}: [1, 2, 3, 2, 2]    => 2003//19200
  RootedTree{Int64}: [1, 2, 2, 2, 2]    => 2003//9600
  RootedTree{Int64}: [1, 2, 3, 4, 5, 6] => 0
  ⋮                                     => ⋮
order_of_accuracy(series_b)
3

Note that this method is only third-order accurate - as is the finishing procedure given by

A = [0 0 0 0 0;
     1//5 0 0 0 0;
     0 2//5 0 0 0;
     161//192 -19//12 287//192 0 0;
     -27//28 19//7 -291//196 36//49 0]
b = [7 // 48, 0, 475 // 1008, 2 // 7, 7 // 72]
rk_c = RungeKuttaMethod(A, b)
series_c = bseries(rk_c, 6)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 38 entries:
  RootedTree{Int64}: Int64[]            => 1
  RootedTree{Int64}: [1]                => 1
  RootedTree{Int64}: [1, 2]             => 1//2
  RootedTree{Int64}: [1, 2, 3]          => 1//6
  RootedTree{Int64}: [1, 2, 2]          => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]       => 41//960
  RootedTree{Int64}: [1, 2, 3, 3]       => 41//480
  RootedTree{Int64}: [1, 2, 3, 2]       => 119//960
  RootedTree{Int64}: [1, 2, 2, 2]       => 119//480
  RootedTree{Int64}: [1, 2, 3, 4, 5]    => 41//4800
  RootedTree{Int64}: [1, 2, 3, 4, 4]    => 41//2400
  RootedTree{Int64}: [1, 2, 3, 4, 3]    => 463//19200
  RootedTree{Int64}: [1, 2, 3, 4, 2]    => 41//1200
  RootedTree{Int64}: [1, 2, 3, 3, 3]    => 463//9600
  RootedTree{Int64}: [1, 2, 3, 3, 2]    => 41//600
  RootedTree{Int64}: [1, 2, 3, 2, 3]    => 639//12800
  RootedTree{Int64}: [1, 2, 3, 2, 2]    => 639//6400
  RootedTree{Int64}: [1, 2, 2, 2, 2]    => 639//3200
  RootedTree{Int64}: [1, 2, 3, 4, 5, 6] => 0
  ⋮                                     => ⋮
order_of_accuracy(series_c)
3

Finally, we can compose the three methods to obtain

series = compose(series_b, series_a, series_c, normalize_stepsize = true)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 38 entries:
  RootedTree{Int64}: Int64[]            => 1
  RootedTree{Int64}: [1]                => 1
  RootedTree{Int64}: [1, 2]             => 1//2
  RootedTree{Int64}: [1, 2, 3]          => 1//6
  RootedTree{Int64}: [1, 2, 2]          => 1//3
  RootedTree{Int64}: [1, 2, 3, 4]       => 1//24
  RootedTree{Int64}: [1, 2, 3, 3]       => 1//12
  RootedTree{Int64}: [1, 2, 3, 2]       => 1//8
  RootedTree{Int64}: [1, 2, 2, 2]       => 1//4
  RootedTree{Int64}: [1, 2, 3, 4, 5]    => 1//120
  RootedTree{Int64}: [1, 2, 3, 4, 4]    => 1//60
  RootedTree{Int64}: [1, 2, 3, 4, 3]    => 1//40
  RootedTree{Int64}: [1, 2, 3, 4, 2]    => 1//30
  RootedTree{Int64}: [1, 2, 3, 3, 3]    => 1//20
  RootedTree{Int64}: [1, 2, 3, 3, 2]    => 1//15
  RootedTree{Int64}: [1, 2, 3, 2, 3]    => 1//20
  RootedTree{Int64}: [1, 2, 3, 2, 2]    => 1//10
  RootedTree{Int64}: [1, 2, 2, 2, 2]    => 1//5
  RootedTree{Int64}: [1, 2, 3, 4, 5, 6] => 121//87480
  ⋮                                     => ⋮

Note that this composition has to be read from left to right. Finally, we check that the resulting series is indeed fifth-order accurate:

order_of_accuracy(series)
5
series - ExactSolution(series)
TruncatedBSeries{RootedTree{Int64, Vector{Int64}}, Rational{Int64}} with 38 entries:
  RootedTree{Int64}: Int64[]            => 0
  RootedTree{Int64}: [1]                => 0
  RootedTree{Int64}: [1, 2]             => 0
  RootedTree{Int64}: [1, 2, 3]          => 0
  RootedTree{Int64}: [1, 2, 2]          => 0
  RootedTree{Int64}: [1, 2, 3, 4]       => 0
  RootedTree{Int64}: [1, 2, 3, 3]       => 0
  RootedTree{Int64}: [1, 2, 3, 2]       => 0
  RootedTree{Int64}: [1, 2, 2, 2]       => 0
  RootedTree{Int64}: [1, 2, 3, 4, 5]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 4]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 3]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 2]    => 0
  RootedTree{Int64}: [1, 2, 3, 3, 3]    => 0
  RootedTree{Int64}: [1, 2, 3, 3, 2]    => 0
  RootedTree{Int64}: [1, 2, 3, 2, 3]    => 0
  RootedTree{Int64}: [1, 2, 3, 2, 2]    => 0
  RootedTree{Int64}: [1, 2, 2, 2, 2]    => 0
  RootedTree{Int64}: [1, 2, 3, 4, 5, 6] => -1//174960
  ⋮                                     => ⋮

References

  • Butcher, J. C. "The effective order of Runge-Kutta methods." In Conference on the numerical solution of differential equations, pp. 133-139. Springer, Berlin, Heidelberg, 1969. DOI: 10.1007/BFb0060019