§3.7 Ordinary Differential Equations

Referenced by:: §11.13(iv), §14.32, §28.34(iv), §29.20(i), Ch.3, §33.23(iii), §4.45(iii)
Permalink:: http://dlmf.nist.gov/3.7

§3.7(i) Introduction
§3.7(ii) Taylor-Series Method: Initial-Value Problems
§3.7(iii) Taylor-Series Method: Boundary-Value Problems
§3.7(iv) Sturm–Liouville Eigenvalue Problems
§3.7(v) Runge–Kutta Method

§3.7(i) Introduction

Keywords:: differential equations
Permalink:: http://dlmf.nist.gov/3.7.i

Consideration will be limited to ordinary linear second-order differential equations

3.7.1 $\frac{{d}^{2}w}{{dz}^{2}}+f(z)\frac{dw}{dz}+g(z)w=h(z),$

Symbols:: $\frac{df}{dx}$ : derivative of $f$ with respect to $x$ , $w(z)$ : function, $f(z)$ : function, $g(z)$ : function and $h(z)$ : function
Referenced by:: §3.7(i), §3.7(ii), §3.7(ii), §3.7(ii)
Permalink:: http://dlmf.nist.gov/3.7.E1
Encodings:: TeX, pMML, png

where $f$ , $g$ , and $h$ are analytic functions in a domain $D\subset\Complex$ . If $h=0$ the differential equation is homogeneous, otherwise it is inhomogeneous. For applications to special functions $f$ , $g$ , and $h$ are often simple rational functions.

For general information on solutions of equation (3.7.1) see §1.13. For classification of singularities of (3.7.1) and expansions of solutions in the neighborhoods of singularities, see §2.7. For an introduction to numerical methods for ordinary differential equations, see Ascher and Petzold (1998), Hairer et al. (1993), and Iserles (1996).

§3.7(ii) Taylor-Series Method: Initial-Value Problems

Keywords:: differential equations, initial-value problems
Referenced by:: §10.74(ii), §13.29(ii), §15.19(ii), §28.34(i), §31.18, §9.17(ii)
Permalink:: http://dlmf.nist.gov/3.7.ii

Assume that we wish to integrate (3.7.1) along a finite path $\mathscr{P}$ from $z=a$ to $z=b$ in a domain $D$ . The path is partitioned at $P+1$ points labeled successively $z_{0},z_{1},\dots,z_{P}$ , with $z_{0}=a$ , $z_{P}=b$ .

By repeated differentiation of (3.7.1) all derivatives of $w(z)$ can be expressed in terms of $w(z)$ and $w^{{\prime}}(z)$ as follows. Write

3.7.2 $w^{{(s)}}(z)=f_{s}(z)w(z)+g_{s}(z)w^{{\prime}}(z)+h_{s}(z),$ $s=0,1,2,\dots$ ,

Symbols:: $w(z)$ : function, $g_{s}(z)$ : coefficient, $h(z)$ : function and $f_{s}(z)$ : coefficient
Referenced by:: §3.7(ii)
Permalink:: http://dlmf.nist.gov/3.7.E2
Encodings:: TeX, pMML, png

with

3.7.3

$f_{0}(z)=1,$

$g_{0}(z)=0,$

$h_{0}(z)=0,$

$f_{1}(z)=0,$

$g_{1}(z)=1,$

$h_{1}(z)=0.$

Symbols:: $g_{s}(z)$ : coefficient, $h(z)$ : function and $f_{s}(z)$ : coefficient
Permalink:: http://dlmf.nist.gov/3.7.E3
Encodings:: TeX, TeX, TeX, TeX, TeX, TeX, pMML, pMML, pMML, pMML, pMML, pMML, png, png, png, png, png, png

Then for $s=2,3,\dots$ ,

3.7.4

$f_{s}(z)=f^{{\prime}}_{{s-1}}(z)-g(z)g_{{s-1}}(z),$

$g_{s}(z)=f_{{s-1}}(z)-f(z)g_{{s-1}}(z)+g^{{\prime}}_{{s-1}}(z),$

$h_{s}(z)=h(z)g_{{s-1}}(z)+h^{{\prime}}_{{s-1}}(z).$

Symbols:: $g_{s}(z)$ : coefficient, $f(z)$ : function, $g(z)$ : function, $h(z)$ : function and $f_{s}(z)$ : coefficient
Permalink:: http://dlmf.nist.gov/3.7.E4
Encodings:: TeX, TeX, TeX, pMML, pMML, pMML, png, png, png

Write $\tau _{j}=z_{{j+1}}-z_{j}$ , $j=0,1,\dots,P$ , expand $w(z)$ and $w^{{\prime}}(z)$ in Taylor series (§1.10(i)) centered at $z=z_{j}$ , and apply (3.7.2). Then

3.7.5 $\begin{bmatrix}w(z_{{j+1}})\\ w^{{\prime}}(z_{{j+1}})\end{bmatrix}=\mathbf{A}(\tau _{j},z_{j})\begin{bmatrix}w(z_{{j}})\\ w^{{\prime}}(z_{{j}})\end{bmatrix}+\mathbf{b}(\tau _{j},z_{j}),$

Symbols:: $(a,b)$ : open interval, $w(z)$ : function and $\mathbf{A}$ : matrix
Referenced by:: §3.7(ii), §3.7(iii)
Permalink:: http://dlmf.nist.gov/3.7.E5
Encodings:: TeX, pMML, png

where $\mathbf{A}(\tau,z)$ is the matrix

3.7.6 $\mathbf{A}(\tau,z)=\begin{bmatrix}A_{{11}}(\tau,z)&A_{{12}}(\tau,z)\\ A_{{21}}(\tau,z)&A_{{22}}(\tau,z)\end{bmatrix},$

Symbols:: $(a,b)$ : open interval, $\mathbf{A}$ : matrix and $\mathbf{A}$ : matrix
Permalink:: http://dlmf.nist.gov/3.7.E6
Encodings:: TeX, pMML, png

and $\mathbf{b}(\tau,z)$ is the vector

3.7.7 $\mathbf{b}(\tau,z)=\begin{bmatrix}b_{1}(\tau,z)\\ b_{2}(\tau,z)\end{bmatrix},$

Symbols:: $(a,b)$ : open interval and $\mathbf{b}(\tau,z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E7
Encodings:: TeX, pMML, png

with

3.7.8

$A_{{11}}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}f_{s}(z),$

$A_{{12}}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}g_{s}(z),$

$A_{{21}}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}f_{{s+1}}(z),$

$A_{{22}}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}g_{{s+1}}(z),$

Symbols:: $!$ : $n!$ : factorial, $(a,b)$ : open interval, $g_{s}(z)$ : coefficient, $\mathbf{A}$ : matrix and $f_{s}(z)$ : coefficient
Permalink:: http://dlmf.nist.gov/3.7.E8
Encodings:: TeX, TeX, TeX, TeX, pMML, pMML, pMML, pMML, png, png, png, png

3.7.9

$b_{1}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}h_{s}(z),$

$b_{2}(\tau,z)=\sum _{{s=0}}^{\infty}\frac{\tau^{s}}{s!}h_{{s+1}}(z).$

Symbols:: $!$ : $n!$ : factorial, $\mathbf{b}(\tau,z)$ : function and $h(z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E9
Encodings:: TeX, TeX, pMML, pMML, png, png

If the solution $w(z)$ that we are seeking grows in magnitude at least as fast as all other solutions of (3.7.1) as we pass along $\mathscr{P}$ from $a$ to $b$ , then $w(z)$ and $w^{{\prime}}(z)$ may be computed in a stable manner for $z=z_{0},z_{1},\dots,z_{P}$ by successive application of (3.7.5) for $j=0,1,\dots,P-1$ , beginning with initial values $w(a)$ and $w^{{\prime}}(a)$ .

Similarly, if $w(z)$ is decaying at least as fast as all other solutions along $\mathscr{P}$ , then we may reverse the labeling of the $z_{j}$ along $\mathscr{P}$ and begin with initial values $w(b)$ and $w^{{\prime}}(b)$ .

§3.7(iii) Taylor-Series Method: Boundary-Value Problems

Keywords:: boundary-value methods or problems, differential equations
Referenced by:: §11.13(iv), §16.25, §28.34(iii), §9.17(ii)
Permalink:: http://dlmf.nist.gov/3.7.iii

Now suppose the path $\mathscr{P}$ is such that the rate of growth of $w(z)$ along $\mathscr{P}$ is intermediate to that of two other solutions. (This can happen only for inhomogeneous equations.) Then to compute $w(z)$ in a stable manner we solve the set of equations (3.7.5) simultaneously for $j=0,1,\dots,P$ , as follows. Let $\mathbf{A}$ be the $(2P)\times(2P+2)$ band matrix

3.7.10 $\mathbf{A}=\begin{bmatrix}-\mathbf{A}(\tau _{0},z_{0})&\mathbf{I}&\boldsymbol{{0}}&\cdots&\boldsymbol{{0}}&\boldsymbol{{0}}\\ \boldsymbol{{0}}&-\mathbf{A}(\tau _{1},z_{1})&\mathbf{I}&\cdots&\boldsymbol{{0}}&\boldsymbol{{0}}\\ \vdots&\vdots&\ddots&\ddots&\vdots&\vdots\\ \boldsymbol{{0}}&\boldsymbol{{0}}&\cdots&-\mathbf{A}(\tau _{{P-2}},z_{{P-2}})&\mathbf{I}&\boldsymbol{{0}}\\ \boldsymbol{{0}}&\boldsymbol{{0}}&\cdots&\boldsymbol{{0}}&-\mathbf{A}(\tau _{{P-1}},z_{{P-1}})&\mathbf{I}\end{bmatrix}$

Symbols:: $(a,b)$ : open interval, $\tau _{j}$ : change of variable and $P$ : partitioning point
Permalink:: http://dlmf.nist.gov/3.7.E10
Encodings:: TeX, pMML, png

( $\mathbf{I}$ and 0 being the identity and zero matrices of order $2\times 2$ .) Also let $\mathbf{w}$ denote the $(2P+2)\times 1$ vector

3.7.11 $\mathbf{w}=\left[w(z_{0}),w^{{\prime}}(z_{0}),w(z_{1}),w^{{\prime}}(z_{1}),\dots,w(z_{P}),w^{{\prime}}(z_{P})\right]^{{\rm T}},$

Symbols:: $w(z)$ : function and $P$ : partitioning point
Permalink:: http://dlmf.nist.gov/3.7.E11
Encodings:: TeX, pMML, png

and $\mathbf{b}$ the $(2P)\times 1$ vector

3.7.12 $\mathbf{b}=\left[b_{1}(\tau _{0},z_{0}),b_{2}(\tau _{0},z_{0}),b_{1}(\tau _{1},z_{1}),b_{2}(\tau _{1},z_{1}),\ldots,b_{1}(\tau _{{P-1}},z_{{P-1}}),b_{2}(\tau _{{P-1}},z_{{P-1}})\right]^{{\rm T}}.$

Symbols:: $(a,b)$ : open interval, $\tau _{j}$ : change of variable and $P$ : partitioning point
Permalink:: http://dlmf.nist.gov/3.7.E12
Encodings:: TeX, pMML, png

Then

3.7.13 $\mathbf{A}\mathbf{w}=\mathbf{b}.$

Referenced by:: §3.7(iii), §3.7(iv)
Permalink:: http://dlmf.nist.gov/3.7.E13
Encodings:: TeX, pMML, png

This is a set of $2P$ equations for the $2P+2$ unknowns, $w(z_{j})$ and $w^{{\prime}}(z_{j})$ , $j=0,1,\dots,P$ . The remaining two equations are supplied by boundary conditions of the form

3.7.14

$\alpha _{0}w(z_{0})+\beta _{0}w^{{\prime}}(z_{0})=\gamma _{0},$

$\alpha _{1}w(z_{P})+\beta _{1}w^{{\prime}}(z_{P})=\gamma _{1},$

Symbols:: $w(z)$ : function, $\alpha$ : constant, $\beta$ : constant, $\gamma$ : constant and $P$ : partitioning point
Permalink:: http://dlmf.nist.gov/3.7.E14
Encodings:: TeX, TeX, pMML, pMML, png, png

where the $\alpha$ ’s, $\beta$ ’s, and $\gamma$ ’s are constants.

If, for example, $\beta _{0}=\beta _{1}=0$ , then on moving the contributions of $w(z_{0})$ and $w(z_{P})$ to the right-hand side of (3.7.13) the resulting system of equations is not tridiagonal, but can readily be made tridiagonal by annihilating the elements of $\mathbf{A}$ that lie below the main diagonal and its two adjacent diagonals. The equations can then be solved by the method of §3.2(ii), if the differential equation is homogeneous, or by Olver’s algorithm (§3.6(v)). The latter is especially useful if the endpoint $b$ of $\mathscr{P}$ is at $\infty$ , or if the differential equation is inhomogeneous.

It will be observed that the present formulation of the Taylor-series method permits considerable parallelism in the computation, both for initial-value and boundary-value problems.

For further information and examples, see Olde Daalhuis and Olver (1998, §7) and Lozier and Olver (1993). General methods for boundary-value problems for ordinary differential equations are given in Ascher et al. (1995).

§3.7(iv) Sturm–Liouville Eigenvalue Problems

Keywords:: Sturm–Liouville eigenvalue problems, differential equations
Referenced by:: §28.34(ii)
Permalink:: http://dlmf.nist.gov/3.7.iv

Let $(a,b)$ be a finite or infinite interval and $q(x)$ be a real-valued continuous (or piecewise continuous) function on the closure of $(a,b)$ . The Sturm–Liouville eigenvalue problem is the construction of a nontrivial solution of the system

3.7.15 $\frac{{d}^{2}w_{k}}{{dx}^{2}}+(\lambda _{k}-q(x))w_{k}=0,$

Symbols:: $\frac{df}{dx}$ : derivative of $f$ with respect to $x$ , $w(z)$ : function, $q(x)$ : real-valued function and $\lambda _{k}$ : eigenvalues
Permalink:: http://dlmf.nist.gov/3.7.E15
Encodings:: TeX, pMML, png

3.7.16 $w_{k}(a)=w_{k}(b)=0,$

Symbols:: $w(z)$ : function
Referenced by:: §3.7(iv)
Permalink:: http://dlmf.nist.gov/3.7.E16
Encodings:: TeX, pMML, png

with limits taken in (3.7.16) when $a$ or $b$ , or both, are infinite. The values $\lambda _{k}$ are the eigenvalues and the corresponding solutions $w_{k}$ of the differential equation are the eigenfunctions. The eigenvalues $\lambda _{k}$ are simple, that is, there is only one corresponding eigenfunction (apart from a normalization factor), and when ordered increasingly the eigenvalues satisfy

3.7.17 $\lambda _{1}<\lambda _{2}<\lambda _{3}<\cdots,$ $\lim _{{k\to\infty}}\lambda _{k}=\infty$ .

Symbols:: $\lambda _{k}$ : eigenvalues
Permalink:: http://dlmf.nist.gov/3.7.E17
Encodings:: TeX, pMML, png

If $q(x)$ is $\mathop{C^{{\infty}}\/}\nolimits$ on the closure of $(a,b)$ , then the discretized form (3.7.13) of the differential equation can be used. This converts the problem into a tridiagonal matrix problem in which the elements of the matrix are polynomials in $\lambda$ ; compare §3.2(vi). The larger the absolute values of the eigenvalues $\lambda _{k}$ that are being sought, the smaller the integration steps $\left|\tau _{j}\right|$ need to be.

For further information, including other methods and examples, see Pryce (1993, §2.5.1).

§3.7(v) Runge–Kutta Method

Keywords:: Runge–Kutta methods, differential equations
Referenced by:: §28.34(i), §32.17
Permalink:: http://dlmf.nist.gov/3.7.v

The Runge–Kutta method applies to linear or nonlinear differential equations. The method consists of a set of rules each of which is equivalent to a truncated Taylor-series expansion, but the rules avoid the need for analytic differentiations of the differential equation.

¶ First-Order Equations

For $w^{{\prime}}=f(z,w)$ the standard fourth-order rule reads

3.7.18 $w_{{n+1}}=w_{n}+\tfrac{1}{6}(k_{1}+2k_{2}+2k_{3}+k_{4})+\mathop{O\/}\nolimits\!\left(h^{5}\right),$

Symbols:: $\mathop{O\/}\nolimits\!\left(x\right)$ : order not exceeding, $w(z)$ : function and $h(z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E18
Encodings:: TeX, pMML, png

where $h=z_{{n+1}}-z_{n}$ and

3.7.19

$k_{1}=hf(z_{n},w_{n}),$

$k_{2}=hf(z_{n}+\tfrac{1}{2}h,w_{n}+\tfrac{1}{2}k_{1}),$

$k_{3}=hf(z_{n}+\tfrac{1}{2}h,w_{n}+\tfrac{1}{2}k_{2}),$

$k_{4}=hf(z_{n}+h,w_{n}+k_{3}).$

Symbols:: $w(z)$ : function, $f(z)$ : function and $h(z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E19
Encodings:: TeX, TeX, TeX, TeX, pMML, pMML, pMML, pMML, png, png, png, png

The order estimate $\mathop{O\/}\nolimits\!\left(h^{5}\right)$ holds if the solution $w(z)$ has five continuous derivatives.

¶ Second-Order Equations

Keywords:: Runge–Kutta methods, differential equations

For $w^{{\prime\prime}}=f(z,w,w^{{\prime}})$ the standard fourth-order rule reads

3.7.20

$w_{{n+1}}=w_{n}+\tfrac{1}{6}h(6w^{{\prime}}_{n}+k_{1}+k_{2}+k_{3})+\mathop{O\/}\nolimits\!\left(h^{5}\right),$

$w^{{\prime}}_{{n+1}}=w^{{\prime}}_{n}+\tfrac{1}{6}(k_{1}+2k_{2}+2k_{3}+k_{4})+\mathop{O\/}\nolimits\!\left(h^{5}\right),$

Symbols:: $\mathop{O\/}\nolimits\!\left(x\right)$ : order not exceeding, $w(z)$ : function and $h(z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E20
Encodings:: TeX, TeX, pMML, pMML, png, png

where

3.7.21

$k_{1}=hf(z_{n},w_{n},w^{{\prime}}_{n}),$

$k_{2}=hf(z_{n}+\tfrac{1}{2}h,w_{n}+\tfrac{1}{2}hw^{{\prime}}_{n}+\tfrac{1}{8}hk_{1},w^{{\prime}}_{n}+\tfrac{1}{2}k_{1}),$

$k_{3}=hf(z_{n}+\tfrac{1}{2}h,w_{n}+\tfrac{1}{2}hw^{{\prime}}_{n}+\tfrac{1}{8}hk_{2},w^{{\prime}}_{n}+\tfrac{1}{2}k_{2}),$

$k_{4}=hf(z_{n}+h,w_{n}+hw^{{\prime}}_{n}+\tfrac{1}{2}hk_{3},w^{{\prime}}_{n}+k_{3}).$

Symbols:: $w(z)$ : function, $f(z)$ : function and $h(z)$ : function
Permalink:: http://dlmf.nist.gov/3.7.E21
Encodings:: TeX, TeX, TeX, TeX, pMML, pMML, pMML, pMML, png, png, png, png

The order estimates $\mathop{O\/}\nolimits\!\left(h^{5}\right)$ hold if the solution $w(z)$ has five continuous derivatives.

An extensive literature exists on the numerical solution of ordinary differential equations by Runge–Kutta, multistep, or other methods. See, for example, Butcher (1987), Dekker and Verwer (1984, Chapter 3), Hairer et al. (1993, Chapter 2), and Hairer and Wanner (1996, Chapter 4).