, 19 min read

Comparing BDF vs. Tendler vs. Tischer formulas

We tested the BDF, Tendler's formulas, new Tendler-like formulas, and Tischer's formulas on the classical test equation

$$ y´(t) = \lambda y \in \mathbb{C}, \quad y(t_0) = e^{\lambda t_0}, \quad t\in[t_0, t_1]. $$

We used

$$ \lambda = r e^{i \varphi}, \quad t_0=0, \quad t1=-40, \quad h\lt 0. $$

The radius is fixed to $r=100$, and $\varphi$ varies from 5° to 90° in steps of 5°. $\varphi$ directly tests the Widlund wedge angle of the formula.

The stepsize $h$ of each formula is varied repeatedly from -0.1, -0.01, and -0.001.

Above differential equation is now solved with the following formulas:

  1. BDF order 1–6
  2. Tendler order 3–7
  3. new Tendler-like formulas order 3–9, see Tendler-like formulas
  4. Tischer's formulas of order 2–8 using $S=0$

This creates 1350 data points.

We report the summed global error computed as

$$ g_\hbox{err} = \sum_{i=0}^n |y(t_i) - y_i|, \quad n={t_1-t_0\over h}. $$

Computations were done in double precision (double complex). Results are here. As the global error $|y(t_i) - y_i|$ per step varies wildly between $10^{-324}$ to $10^{+304}$, we computed

$$ \hat g_\hbox{err} = \begin{cases} 5, & \hbox{if } g_\hbox{err} \gt 30 \hbox{ or NaN}\\ \log_{10} g_\hbox{err}, & \hbox{else} \end{cases} $$

That's what is shown in below graphic.

Below 3D chart can be rotated and zoomed in or out. Clicking on any entry shows the log global error $\hat g_\hbox{err}$.

Once can clearly see that the higher order methods quickly lose precision when a higher Widlund wedge angle is required. This is the reason why GEAR, EPISODE, LSODE and CVODE all do not use BDF6. However, in reality, this is mainly a problem for the stepsize and order control segment to properly switch order.

For the new Tendler-like formulas we intend to use even the higher order methods of order 7, 8, and 9 as we aim to provide a type-insensitive code, which just switches between fixed point iteration and modified Newton method. It is well known, see Montenbruck/Gill (2000) §4.1.6, that higher order methods are indeed required for certain precision.

As Gaffney (1984) noted a sensitivity regarding the used precision, we repeated above test in single precision (float complex). The qualitative results didn't differ in any way. Results are here.