# Numerov's method

> Mediated Wiki article. Canonical URL: https://mediated.wiki/source/Numerov's_method
> Markdown URL: https://mediated.wiki/source/Numerov's_method.md
> Source: https://en.wikipedia.org/wiki/Numerov's_method
> Source revision: 1320532658
> License: Creative Commons Attribution-ShareAlike 4.0 International (https://creativecommons.org/licenses/by-sa/4.0/)

**Numerov's method** (also called **Cowell's method**) is a [numerical method](/source/Numerical_method) to solve [ordinary differential equations](/source/Ordinary_differential_equation) of second order in which the first-order term does not appear. It is a fourth-order [linear multistep method](/source/Linear_multistep_method). The method is implicit, but can be made explicit if the differential equation is linear.

Numerov's method was developed by the Russian astronomer [Boris Vasil'evich Numerov](/source/Boris_Vasil'evich_Numerov).

## The method

The Numerov method can be used to solve differential equations of the form

- \frac{d^2 y}{dx^2} = - g(x) y(x) + s(x).

In it, three values of y_{n-1}, y_n, y_{n+1} taken at three equidistant points x_{n-1}, x_n, x_{n+1} are related as follows:

- y_{n+1} \left(1 + \frac{h^2}{12} g_{n+1}\right) = 2 y_n \left(1 - \frac{5 h^2}{12} g_n\right) - y_{n-1} \left(1 + \frac{h^2}{12} g_{n-1}\right) + \frac{h^2}{12} (s_{n+1} + 10 s_n + s_{n-1}) + \mathcal{O}(h^6),

where y_n = y(x_n), g_n = g(x_n), s_n = s(x_n), and h = x_{n+1} - x_n.

### Nonlinear equations

For nonlinear equations of the form

- \frac{d^2 y}{dx^2} = f(x,y),

the method gives

- y_{n+1} - 2 y_n + y_{n-1} = \frac{h^2}{12} (f_{n+1} + 10 f_n + f_{n-1}) + \mathcal{O} (h^6).

This is an implicit [linear multistep method](/source/Linear_multistep_method), which reduces to the explicit method given above if f is linear in y by setting f(x,y) = - g(x) y(x) + s(x). It achieves order-4 accuracy (Hairer, Nørsett & Wanner 1993, §III.10).

## Application

In numerical physics the method is used to find solutions of the unidimensional [Schrödinger equation](/source/Schr%C3%B6dinger_equation) for arbitrary potentials. An example of which is solving the radial equation for a spherically symmetric potential. In this example, after separating the variables and analytically solving the angular equation, we are left with the following equation of the radial function R(r):

- \frac{d}{dr} \left(r^2 \frac{dR}{dr} \right) - \frac{2 m r^2}{\hbar^2} (V(r) - E) R(r) = l (l + 1) R(r).

This equation can be reduced to the form necessary for the application of Numerov's method with the following substitution:

- u(r) = r R(r) \Rightarrow R(r) = \frac{u(r)}{r},
- \frac{dR}{dr} = \frac{1}{r} \frac{du}{dr} - \frac{u(r)}{r^2} = \frac{1}{r^2} \left (r \frac{du}{dr} - u(r) \right) \Rightarrow \frac{d}{dr} \left(r^2 \frac{dR}{dr} \right) = \frac{du}{dr} + r \frac{d^2 u}{dr^2} - \frac{du}{dr} = r \frac{d^2 u}{dr^2}.

And when we make the substitution, the radial equation becomes

- r \frac{d^2 u}{dr^2} - \frac{2 m r}{\hbar^2} (V(r) - E) u(r) = \frac{l (l + 1)}{r} u(r),

or

- -\frac{\hbar^2}{2 m} \frac{d^2 u}{dr^2} + \left(V(r) + \frac{\hbar^2}{2 m} \frac{l (l + 1)}{r^2} \right) u(r) = E u(r),

which is equivalent to the one-dimensional Schrödinger equation, but with the modified effective potential

- V_\text{eff}(r) = V(r) + \frac{\hbar^2}{2 m} \frac{l (l + 1)}{r^2} = V(r) + \frac{L^2}{2 m r^2}, \quad L^2 = l (l + 1) \hbar^2.

This equation we can proceed to solve the same way we would have solved the one-dimensional Schrödinger equation. We can rewrite the equation a little bit differently and thus see the possible application of Numerov's method more clearly:

- \frac{d^2 u}{dr^2} = - \frac{2 m}{\hbar^2} (E - V_\text{eff}(r)) u(r),

- g(r) = \frac{2 m}{\hbar^2} (E - V_\text{eff}(r)),

- s(r) = 0.

## Derivation

We are given the differential equation

- y''(x) = - g(x) y(x) + s(x).

To derive the Numerov's method for solving this equation, we begin with the [Taylor expansion](/source/Taylor_expansion) of the function we want to solve, y(x), around the point x_0:

- y(x) = y(x_0) + (x-x_0)y'(x_0) + \frac{(x-x_0)^2}{2!}y''(x_0) + \frac{(x-x_0)^3}{3!}y'''(x_0) + \frac{(x-x_0)^4}{4!}y''''(x_0) + \frac{(x-x_0)^5}{5!}y'''''(x_0) + \mathcal{O}(h^6).

Denoting the distance from x to x_0 by h = x - x_0, we can write the above equation as

- y(x_0 + h) = y(x_0) + hy'(x_0) + \frac{h^2}{2!}y''(x_0) + \frac{h^3}{3!}y'''(x_0) + \frac{h^4}{4!}y''''(x_0) + \frac{h^5}{5!}y'''''(x_0) + \mathcal{O}(h^6).

If we evenly discretize the space, we get a grid of x points, where h = x_{n+1} - x_n. By applying the above equations to this discrete space, we get a relation between the y_n and y_{n+1}:

- y_{n+1} = y_n + hy'(x_n) + \frac{h^2}{2!}y''(x_n) + \frac{h^3}{3!}y'''(x_n) + \frac{h^4}{4!}y''''(x_n) + \frac{h^5}{5!}y'''''(x_n) + \mathcal{O}(h^6).

Computationally, this amounts to taking a step *forward* by an amount h. If we want to take a step *backwards*, we replace every h with - h and get the expression for y_{n-1}:

- y_{n-1} = y_n - hy'(x_n) + \frac{h^2}{2!}y''(x_n) - \frac{h^3}{3!}y'''(x_n) + \frac{h^4}{4!}y''''(x_n) - \frac{h^5}{5!}y'''''(x_n) + \mathcal{O}(h^6).

Note that only the odd powers of h experienced a sign change. By summing the two equations, we derive that

- y_{n+1} - 2 y_n + y_{n-1} = h^2 y''_n + \frac{h^4}{12}y''''_n + \mathcal{O}(h^6).

We can solve this equation for y_{n+1} by substituting the expression given at the beginning, that is y''_n = - g_n y_n + s_n. To get an expression for the y''''_n factor, we simply have to differentiate y''_n = - g_n y_n + s_n twice and approximate it again in the same way we did this above:

- y''''_n = \frac{d^2}{d x^2} (-g_n y_n + s_n),
- h^2 y''''_n = -g_{n+1} y_{n+1} + s_{n+1} + 2 g_n y_n - 2 s_n - g_{n-1} y_{n-1} + s_{n-1} + \mathcal{O}(h^4).

If we now substitute this to the preceding equation, we get

- y_{n+1} - 2 y_n + y_{n-1} = {h^2} (- g_n y_n + s_n) + \frac{h^2}{12} (- g_{n+1} y_{n+1} + s_{n+1} + 2 g_n y_n - 2 s_n - g_{n-1} y_{n-1} + s_{n-1}) + \mathcal{O}(h^6),

or

- y_{n+1} \left(1 + \frac{h^2}{12} g_{n+1} \right) - 2 y_n \left(1 - \frac{5 h^2}{12} g_n \right) + y_{n-1} \left(1 + \frac{h^2}{12} g_{n-1} \right) = \frac{h^2}{12} (s_{n+1} + 10 s_n + s_{n-1}) + \mathcal{O}(h^6).

This yields the Numerov's method if we ignore the term of order h^6. It follows that the order of convergence (assuming stability) is 4.

## References

- Hairer, Ernst; Nørsett, Syvert Paul; Wanner, Gerhard (1993), *Solving ordinary differential equations I: Nonstiff problems*, Berlin, New York: [Springer-Verlag](/source/Springer-Verlag), ISBN 978-3-540-56670-0. This book includes the following references:
- Numerov, Boris Vasil'evich (1924), "A method of extrapolation of perturbations", *[Monthly Notices of the Royal Astronomical Society](/source/Monthly_Notices_of_the_Royal_Astronomical_Society)*. **84** (8): 592–601, [Bibcode:1924MNRAS..84..592N](https://ui.adsabs.harvard.edu/abs/1924MNRAS..84..592N). [doi:10.1093/mnras/84.8.592](https://doi.org/10.1093/mnras/84.8.592).
- Numerov, Boris Vasil'evich (1927), "Note on the numerical integration of d2*x*/d*t*2 = *f*(*x*,*t*)", *[Astronomische Nachrichten](/source/Astronomische_Nachrichten)*. **230** (19): 359–364, [Bibcode:1927AN....230..359N](https://ui.adsabs.harvard.edu/abs/1927AN....230..359N). [doi:10.1002/asna.19272301903](https://doi.org/10.1002/asna.19272301903).

---
Adapted from the Wikipedia article [Numerov's method](https://en.wikipedia.org/wiki/Numerov's_method) by Wikipedia contributors ([contributor history](https://en.wikipedia.org/wiki/Numerov's_method?action=history)). Available under [Creative Commons Attribution-ShareAlike 4.0 International](https://creativecommons.org/licenses/by-sa/4.0/). Changes may have been made.