Multiderivative multistep method

A multiderivative multistep method uses higher time derivatives as well as values and first derivatives from several time levels. For an autonomous ordinary differential equation $y^{\prime }=f(y)$, the second time derivative is $y^{\prime \prime }=f^{\prime }(y)f(y)$ by the chain rule. Its local truncation error must include the higher-derivative terms. Consequently the usual second-order barrier for A-stable linear multistep methods, which only use first derivatives, does not apply.

For $y^{\prime }=f(y)$ and $g=f^{\prime }f$, inserting a smooth exact solution into the displayed recurrence gives defect $h^7{y}^{(7)}(t_n)/4725+O(h^9)$. This follows by expanding $2\sinh z-(14/15)z\cosh z-(16/15)z+(2/15)z^2\sinh z$. The terms through $z^6$ vanish and the $z^7$ coefficient is nonzero. Its zero-step polynomial is ${\zeta }^2-1$, whose roots are simple and of unit modulus. Thus it is zero-stable and has order six, provided starting values and derivative evaluations have corresponding accuracy.

For the multiderivative multistep method $y_{n+1}-\alpha h{f}_{n+1}+\beta h^2{g}_{n+1}=\gamma h{f}_n+y_{n-1}+\alpha h{f}_{n-1}+\beta h^2{g}_{n-1}$, with $g=f^{\prime }f$, centered expansion gives defect coefficients $2-2\alpha -\gamma$, $1/3-\alpha +2\beta$, and $1/60-\alpha /12+\beta /3$ at powers $h,h^3,h^5$. Their simultaneous vanishing gives the displayed parameters. The next coefficient is $1/4725$ at $h^7$, so the formal order is six. Reflection symmetry removes every even defect power. The method satisfies fifth-order conditions but cannot have exact order five.

For the symmetric two-step two-derivative formula, the closed-disk root condition for a multistep method fails off the imaginary axis. Positive real part makes the product of root moduli exceed one; negative real part fails the complex quadratic Schur criterion, as in the empty strict decay domain of the symmetric two-derivative formula proof. On $z=iy$, write $D=15-y^2-7iy=|D|e^{i\delta }$ and set $\xi =i{e}^{-i\delta }\eta$. The amplification polynomial becomes ${\eta }^2-16y\eta /|D|+1=0$. Its two roots are distinct and on the unit circle exactly when $|16y/|D||<2$, or $y^4-45y^2+225>0$. Thus the bounded-root set consists of the imaginary intervals $|y|<\sqrt{(45-15\sqrt{5})/2}$ and $|y|>\sqrt{(45+15\sqrt{5})/2}$, including $y=0$ and excluding the double-root endpoints. In particular, the set is unbounded. Unit-modulus oscillations are bounded but do not tend to zero.

On the Dahlquist test equation, write the amplification polynomial as $D{\xi }^2-16z\xi -A$, where $D=z^2-7z+15$ and $A=z^2+7z+15$. The identity $|A{|}^2-|D{|}^2=28\mathrm{Re}z(|z{|}^2+15)$ shows that both roots of a polynomial can lie inside the open unit disk only if $\mathrm{Re}z<0$. But the complex quadratic Schur criterion would then require $32|\mathrm{Re}z|(|z{|}^2+15)<28|\mathrm{Re}z|(|z{|}^2+15)$, which is impossible. Thus the strict linear stability domain is empty. This concerns robust damping of all modes, rather than convergence on fixed time intervals or the bounded-root stability set of the symmetric two-derivative formula.

On the Dahlquist test equation, the symmetric two-step two-derivative formula has amplification polynomial $F_z(w)=(z^2-7z+15)w^2-16zw-(z^2+7z+15)$. For every real $z<0$, its leading coefficient is positive and $F_z(-1)=2z<0$, whereas $F_z(w)\to +\infty$ as $w\to -\infty$. The intermediate value theorem therefore gives a real root below minus one. Thus the method is not A-stable, without any appeal to a multistep order barrier.

The endpoint expansion of this multiderivative multistep method has unscaled defect $h^4{y}^{(4)}/3+O(h^5)$ and zero-stability roots $1,1/7$, giving third-order convergence with suitably accurate starting values. On the Dahlquist test equation its amplification roots obey $(7-6z+2z^2){\zeta }^2-8\zeta +1=0$. The complex quadratic Schur criterion puts both roots in the open unit disk for $\mathrm{Re}z<0$, and only the simple root one reaches the boundary at $z=0$. Thus it is A-stable; its higher derivative places it outside the hypotheses of the Second Dahlquist barrier.

For a fixed-step multiderivative multistep method, the ordinary zero-stability root condition still controls propagation of the starting errors when all derivative evaluation maps are uniformly Lipschitz continuous on the relevant bounded region. A local defect $O(h^{p+1})$ then yields global error $O(h^p)$ over a fixed time interval, provided the starting errors are $O(h^p)$ and the implicit updates use the nearby solution branch.

For example, writing the numerical error equation as $\rho (E)e_n=h{F}_h(e_{n+s})+d_n$, with $F_h$ uniformly Lipschitz continuous, a bounded impulse response for the root condition for a multistep method givesAbsorbing the last current-step term for small $h$ and applying the discrete Gronwall inequality gives the stated order. The higher derivatives enter through $F_h$, whose bound stays uniform as $h\to 0$.