higher-order taylor methods for odes
Higher-order Taylor methods for ODEs are advanced numerical techniques for solving initial value problems. These methods use Taylor series expansions to approximate solutions, offering improved accuracy and efficiency compared to simpler methods like Euler's method. The order of the method determines the degree of the Taylor polynomial used, with higher orders providing better accuracy but requiring more complex computations. Taylor methods excel in problems with smooth solutions, balancing accuracy and computational cost for optimal results in various scientific and engineering applications.
Consider an initial value problem $y'(t) = f(t, y(t))$ with $y(t_0) = y_0$
The Taylor series expansion of the solution $y(t)$ around $t_0$ is given by:
$y(t) = y(t_0) + y'(t_0)(t - t_0) + \frac{y''(t_0)}{2!}(t - t_0)^2 + \frac{y'''(t_0)}{3!}(t - t_0)^3 + \cdots$
The derivatives of $y(t)$ can be expressed in terms of $f(t, y(t))$ using the chain rule:
$y'(t) = f(t, y(t))$ $y''(t) = \frac{\partial f}{\partial t} + \frac{\partial f}{\partial y}f(t, y(t))$ $y'''(t) = \frac{\partial^2 f}{\partial t^2} + 2\frac{\partial^2 f}{\partial t \partial y}f(t, y(t)) + \frac{\partial^2 f}{\partial y^2}(f(t, y(t)))^2 + \frac{\partial f}{\partial y}\left(\frac{\partial f}{\partial t} + \frac{\partial f}{\partial y}f(t, y(t))\right)$
The $n$-th order Taylor method approximates the solution by truncating the Taylor series at the $n$-th term:
$y(t_0 + h) \approx y(t_0) + hf(t_0, y(t_0)) + \frac{h^2}{2!}y''(t_0) + \cdots + \frac{h^n}{n!}y^{(n)}(t_0)$
The local truncation error is of order $O(h^{n+1})$, while the global error is of order $O(h^n)$
The coefficients of the Taylor series can be computed efficiently using automatic differentiation or symbolic computation
Consider the initial value problem $y'(t) = t^2 + y^2$ with $y(0) = 1$. Implement a 4th-order Taylor method to approximate the solution at $t = 0.5$ with a step size of $h = 0.1$. Compare the result with the exact solution $y(t) = \tan(t^3/3 + \pi/4)$.
Solution:
Solve the Van der Pol oscillator equation $y''(t) = \mu(1 - y^2)y'(t) - y(t)$ with $\mu = 1$ and initial conditions $y(0) = 2$, $y'(0) = 0$ using a 6th-order Taylor method. Compute the solution for $t \in [0, 10]$ with a step size of $h = 0.1$ and plot the result.
Solution:
Implement an adaptive step size control scheme for a 5th-order Taylor method and apply it to the Lorenz system: $x'(t) = \sigma(y - x)$, $y'(t) = x(\rho - z) - y$, $z'(t) = xy - \beta z$ with parameters $\sigma = 10$, $\rho = 28$, $\beta = 8/3$ and initial conditions $x(0) = 1$, $y(0) = 1$, $z(0) = 1$. Compute the solution for $t \in [0, 50]$ and plot the result in 3D.
Solution: