Английская Википедия:Five-point stencil
In numerical analysis, given a square grid in one or two dimensions, the five-point stencil of a point in the grid is a stencil made up of the point itself together with its four "neighbors". It is used to write finite difference approximations to derivatives at grid points. It is an example for numerical differentiation.
In one dimension
In one dimension, if the spacing between points in the grid is h, then the five-point stencil of a point x in the grid is
- <math display="block"> \{x-2h, x-h, x, x+h, x+2h\}.</math>
1D first derivative
The first derivative of a function f of a real variable at a point x can be approximated using a five-point stencil as:[1]
- <math display="block">f'(x) \approx \frac{-f(x+2 h)+8 f(x+h)-8 f(x-h)+f(x-2h)}{12 h}</math>
The center point f(x) itself is not involved, only the four neighboring points.
Derivation
This formula can be obtained by writing out the four Taylor series of f(x ± h) and f(x ± 2h) up to terms of h3 (or up to terms of h5 to get an error estimation as well) and solving this system of four equations to get f ′(x). Actually, we have at points x + h and x − h:
- <math display="block">f(x \pm h) = f(x) \pm h f'(x) + \frac{h^2}{2}f(x) \pm \frac{h^3}{6} f^{(3)}(x) + O_{1\pm}(h^4). \qquad (E_{1\pm}).</math>
Evaluating <math>(E_{1+})-(E_{1-})</math> gives us
- <math display="block">f(x+h) - f(x-h) = 2hf'(x) + \frac{h^3}{3}f^{(3)}(x) + O_1(h^4). \qquad (E_1).</math>
The residual term O1(h4) should be of the order of h5 instead of h4 because if the terms of h4 had been written out in (E1+) and (E1−), it can be seen that they would have canceled each other out by Шаблон:Math. But for this calculation, it is left like that since the order of error estimation is not treated here (cf below).
Similarly, we have
- <math display="block">f(x \pm 2h) = f(x) \pm 2h f'(x) + \frac{4h^2}{2!} f(x) \pm \frac{8h^3}{3!} f^{(3)}(x) + O_{2\pm}(h^4). \qquad (E_{2\pm})</math>
and <math>(E_{2+}) - (E_{2-})</math> gives us
- <math display="block">f(x+2h) - f(x-2h) = 4hf'(x) + \frac{8h^3}{3}f^{(3)}(x) + O_2(h^4). \qquad (E_2).</math>
In order to eliminate the terms of ƒ(3)(x), calculate 8 × (E1) − (E2)
- <math display="block">8f(x+h) - 8f(x-h) - f(x+2h) + f(x-2h) = 12h f'(x) + O(h^4)</math>
thus giving the formula as above. Note: the coefficients of f in this formula, (8, -8,-1,1), represent a specific example of the more general Savitzky–Golay filter.
Error estimate
The error in this approximation is of order h 4. That can be seen from the expansion[2]
- <math display="block"> \frac{-f(x+2 h)+8 f(x+h)-8 f(x-h)+f(x-2h)}{12 h}=f'(x)-\frac{1}{30} f^{(5)}(x) h^4+O(h^5)</math>
which can be obtained by expanding the left-hand side in a Taylor series. Alternatively, apply Richardson extrapolation to the central difference approximation to <math>f'(x)</math> on grids with spacing 2h and h.
1D higher-order derivatives
The centered difference formulas for five-point stencils approximating second, third, and fourth derivatives are
- <math display="block"> \begin{align}
f(x) &\approx \frac{-f(x+2 h)+16 f(x+h)-30 f(x) + 16 f(x-h) - f(x-2h)}{12 h^2} \\[1ex] f^{(3)}(x) &\approx \frac{f(x+2 h)-2 f(x+h) + 2 f(x-h) - f(x-2h)}{2 h^3} \\[1ex] f^{(4)}(x) &\approx \frac{f(x+2 h)-4 f(x+h)+6 f(x) - 4 f(x-h) + f(x-2h)}{h^4}
\end{align} </math>
The errors in these approximations are O(hШаблон:I sup), O(hШаблон:I sup) and O(hШаблон:I sup) respectively.[2]
Relationship to Lagrange interpolating polynomials
As an alternative to deriving the finite difference weights from the Taylor series, they may be obtained by differentiating the Lagrange polynomials
- <math display="block">\ell_j(\xi) = \prod_{i=0,\, i\neq j}^{k} \frac{\xi-x_i}{x_j-x_i}, </math>
where the interpolation points are
- <math display="block">
x_0 = x-2h,\quad x_1 = x-h,\quad x_2 = x,\quad x_3 = x+h,\quad x_4 = x+2h. </math>
Then, the quartic polynomial <math>p_4(x)</math> interpolating Шаблон:Math at these five points is
- <math display="block"> p_4(x) = \sum_{j=0}^4 f(x_j) \ell_j(x) </math>
and its derivative is
- <math display="block"> p_4'(x) = \sum_{j=0}^4 f(x_j) \ell'_j(x). </math>
So, the finite difference approximation of Шаблон:Math at the middle point Шаблон:Math is
- <math display="block">
f'(x_2) = \ell_0'(x_2) f(x_0) + \ell_1'(x_2) f(x_1) + \ell_2'(x_2) f(x_2) + \ell_3'(x_2) f(x_3) + \ell_4'(x_2) f(x_4) + O(h^4) </math>
Evaluating the derivatives of the five Lagrange polynomials at Шаблон:Math gives the same weights as above. This method can be more flexible as the extension to a non-uniform grid is quite straightforward.
In two dimensions
In two dimensions, if for example the size of the squares in the grid is h by h, the five point stencil of a point (x, y) in the grid is
- <math display="block">\{(x-h, y), (x, y), (x+h, y), (x, y-h), (x, y+h)\}, </math>
forming a pattern that is also called a quincunx. This stencil is often used to approximate the Laplacian of a function of two variables:
- <math display="block"> \nabla^2 f(x,y) \approx \frac{f(x-h,y) + f(x+h,y) + f(x,y-h) + f(x,y+h) - 4f(x,y)}{h^2}. </math>
The error in this approximation is O(h 2),[3] which may be explained as follows:
From the 3 point stencils for the second derivative of a function with respect to x and y:
- <math display="block">\begin{align}
\frac{\partial ^2 f}{\partial x^2} &= \frac{f\left(x + \Delta x,y\right) + f\left(x - \Delta x,y\right) - 2f(x,y)}{\Delta x^2} - 2\frac{f^{(4)}(x,y)}{4!} \Delta x^2 + \cdots \\[1ex] \frac{\partial ^2 f}{\partial y^2} &= \frac{f\left(x, y + \Delta y\right) + f\left(x, y - \Delta y\right) - 2f(x,y)}{\Delta y^2} - 2\frac{f^{(4)}(x,y)}{4!} \Delta y^2 + \cdots \end{align}</math>
If we assume <math>\Delta x = \Delta y = h</math>:
- <math display="block">\begin{align}
\nabla^2 f &= \frac{\partial ^2 f}{\partial x^2} + \frac{\partial ^2 f}{\partial y^2} \\[1ex]
&= \frac{f\left(x + h,y\right) + f\left(x - h,y\right) + f\left(x, y + h\right) + f\left(x, y - h\right) - 4f(x,y)}{h^2} - 4 \frac{f^{(4)}(x,y)}{4!}h^2 + \cdots \\[1ex] &= \frac{f\left(x + h,y\right) + f\left(x - h,y\right) + f\left(x, y + h\right) + f\left(x, y - h\right) - 4f(x,y)}{h^2} + O\left(h^2\right)\\
\end{align}</math>
See also
- Шаблон:Annotated link
- Шаблон:Annotated link
- Nine-point stencil
- Stencil (numerical analysis)
- Шаблон:Annotated link
References
- Шаблон:Citation. Ninth printing. Table 25.2.
- ↑ Шаблон:Cite book
- ↑ 2,0 2,1 Abramowitz & Stegun, Table 25.2
- ↑ Abramowitz & Stegun, 25.3.30