6.4. Higher-Order Taylor Methods#

Remark

(Forward) Euler’s method is Taylor’s method of order one.

Assume that \(y(t)\) is the solution of the IVP (6.11) such that it has \((n+1)\) continuous derivatives. If expand \(y(t)\) in terms of its \(n\)th Taylor polynomial about \(t_{i}\) and then evaluate \(t_{i+1}\), we have [Burden and Faires, 2005],

(6.37)#\[y(t_{i+1}) = y(t_{i}) + \frac{h}{1!}y'(t_{i}) + \frac{h^2}{2!}y''(t_{i})+ \ldots + \frac{h^n}{n!}y^{(n)}(t_{i}) + \frac{h^{n+1}}{(n+1)!}y^{(n+1)}(\xi_{i}),\]

where \(\xi_{i}\) is a number between \(t_{i}\) and \(t_{i+1}\). Now since,

(6.38)#\[\begin{split}\begin{cases} y'(t) = f(t,~y(t)),\\ y''(t) = f'(t,~y(t)),\\ ~~ \vdots\\ y^{(k)}(t) = f^{(k-1)}(t,~y(t)). \end{cases}\end{split}\]

It follows from (6.37) and (6.38) that,

(6.39)#\[\begin{align}\label{eq.5.17} y(t_{i+1}) &= y(t_{i}) + h\,f(t_{i},~y(t_{i})) + \frac{h^2}{2}f'(t_{i},~y(t_{i})) + \ldots \notag \\ & + \frac{h^n}{n!}f^{(n-1)}(t_{i},~y(t_{i})) + \frac{h^{n+1}}{(n+1)!}f^{(n)}(\xi_{i},~y(\xi_{i})). \end{align}\]

The difference-equation method can be achieved by deleting the remainder term involving \(\xi_{i}\) in the above equation.

The following definition is Definition 5.11 from [Burden and Faires, 2005].

The following theorem isTheorem 5.12 from [Burden and Faires, 2005].

Proof: Observe that the equation \eqref{eq.5.17} can be expressed as,

(6.40)#\[\begin{equation} y(t_{i+1}) - y(t_{i}) - h\,f(t_{i},~y_{i}) - \frac{h^2}{2}f'(t_{i},~y_{i})- \ldots - \frac{h^n}{n!}f^{(n-1)}(t_{i},~y_{i}) = \frac{h^{n+1}}{(n+1)!}f^{(n)}(\xi_{i},~y(\xi_{i})), \end{equation}\]

where \(\xi_{i}\) is a number between \(t_{i}\) and \(t_{i+1}\). Thus, the local truncation error is

(6.41)#\[\begin{equation} \tau_{i+1}(h) = \frac{y_{i+1} - y_{i}}{h} - T^{(n)}(t_{i},~y_{i}) = \frac{h^{n}}{(n+1)!}f^{(n)}(\xi_{i},~y(\xi_{i}))\qquad i = 0, 1,\ldots, n - 1. \end{equation}\]

As \(y\in C^{n+1}[a,~b]\), we have \(y^{(n+1)}(t) = f^{(n)} (t,~y(t))\) bounded on \([a,b]\) and \(\tau_{i}(h) = \mathcal{O}(h^{n})\) for each \(i = 0, 1,\ldots, n\).

Example: Consider the following IVP, Apply Taylor’s method of orders two with \(n = 10\) to the initial-value problem problem

\[\begin{align*} \begin{cases} y' = t^3+3\,t^2-y+1,& 0 \leq t \leq 1,\\ y\left(0\right)=1, \end{cases} \end{align*}\]

with exact solution

\[\begin{align*} y\left(t \right) = t^3+1. \end{align*}\]

Solution:

For the method of order two we need the 1st, the 2nd and 3rd derivatives of

\[\begin{align*} f(t, y) = t^3+3\,t^2-y+1. \end{align*}\]

with respect to the variable \(t\). Thus,

\[\begin{align*} f'(t, y(t))= -y'+3\,t^2+6\,t. \end{align*}\]

Since \(y' = t^3+3\,t^2-y+1\),

\[\begin{align*} % 6*t + y(t) - t^3 - 1 f'(t, y(t)) = -(t^3+3\,t^2-y+1)+3\,t^2+6\,t = -t^3+6\,t+y-1. \end{align*}\]

Therefore,

(6.42)#\[\begin{align} T^{(2)}(t_{i}, w_{i}) = f(t_{i}, w_{i}) + \frac{h}{2}\,f'(t_{i}, w_{i}) = t_{i}^3+3\,t_{i}^2-w_{i}+1 + \frac{h}{2}\,\left(- t_{i}^3+6\,t_{i}+w_{i}-1 \right) \end{align}\]

and

(6.43)#\[\begin{align} \begin{cases} w_{0} = \alpha,&\\ w_{j+1} = w_{j} + h\,T^{(2)}(t_{i}, w_{i}),&i = 0, 1,\ldots, 9. \end{cases} \end{align}\]
import numpy as np
import pandas as pd
a = 0; b = 1
def Example1(N, a = a, b = b):
    y_exact = lambda t:  t**3 + 1;
    alpha = 1
    t = np.linspace(a, b, N+1)
    h = (b-a)/(N)
    f = lambda t, y: t**3 + 3*(t**2) - y + 1;
    f1 = lambda t, y: -t**3 + 6*t + y -1;
    fun = lambda t, y, h: y + h*f(t, y) + (h**2)*f1(t, y)/2
    w = np.zeros(t.shape, dtype=float)
    w[0] = alpha
    for j in range(N):
        w[j+1] = fun(t[j], w[j], h)
    Table = pd.DataFrame({'t':t, 'y': y_exact(t),
                          'w':w, 'Absolute Error':  abs(y_exact(t) - w)})
    return Table

N = [2**(i) for i in range(3, 11)]
Eh_list = []
for n in N:
    TB = Example1(N = n)
    Eh_list.append(np.max(np.abs(TB['Absolute Error'])[1:]))
Table = pd.DataFrame(data = {'h':[(b-a)/n for n in N], 'n':N, 'Eh': Eh_list})
display(Table.style.set_properties(subset=['h'], **{'background-color': 'Lavender', 'color': 'Navy',
                                                'border-color': 'DarkGreen'}).format({'Eh': "{:.6e}"}))
  h n Eh
0 0.125000 8 1.051778e-02
1 0.062500 16 2.547861e-03
2 0.031250 32 6.270429e-04
3 0.015625 64 1.555377e-04
4 0.007812 128 3.873265e-05
5 0.003906 256 9.664255e-06
6 0.001953 512 2.413705e-06
7 0.000977 1024 6.031316e-07