# Derive an $O(h^4)$ five point formula to approximate $f'(x_0)$ that uses $f(x_0 - h), f(x_0), f(x_0 + h), f(x_0 + 2h), f(x_0 + 3h)$

This is a Numerical Analysis textbook question on Numerical Differentiation (like https://en.wikipedia.org/wiki/Numerical_differentiation):

## Official Textbook Question

Q: Derive an $O(h^4)$ five point formula to approximate $f'(x_0)$ that uses $f(x_0 - h), f(x_0), f(x_0 + h), f(x_0 + 2h), f(x_0 + 3h)$

Hint: Consider the expression $A f(x_0 - h) + B f(x_0 + h) + C f(x_0 + 2h) + D f(x_0 + 3h)$. Expand in fourth Taylor polynomials, and choose $A, B, C, D$ appropriately.

## My Work

I can follow the full textbook derivation of the "Three Point" Formulas, the Three Point Endpoint Formula:

\begin{align*} f'(x_0) &= \frac{1}{2h} \left[ -3 f(x_0) + 4 f(x_0 + h) - f(x_0 + 2h) \right] + \frac{h^2}{3} f^{(3)}(\xi_0) \\ \end{align*}

where $\xi_0$ is between $x_0$ and $x_0 + 2h$

and the Three-Point Midpoint Formula:

\begin{align*} f'(x_0) &= \frac{1}{2h} \left[ f(x_0 + h) - f(x_0 - h) \right] - \frac{h^2}{6} f^{(3)}(\xi_1) \\ \end{align*}

where $\xi_1$ is between $x_0 - h$ and $x_0 + h$

These derivations involve calculating the Lagrange interpolating polynomial of three evenly spaced points, calculating the derivative, and solving for the derivative at each point. Using either the first or last point gives you the three point endpoint formula and using the midpoint gives you the three point midpoint formula.

Obviously, the same technique would work for five points, but the equations get very tedious and hard to work with, and the problem hint obviously doesn't want us to solve the problem that way.

I don't quite understand how the hint wants me to proceed. Is a fourth Taylor polynomial something like this:

\begin{align*} f(x) &= f(x_0) + \frac{f'(x_0)}{1!} (x - x_0) + \frac{f^{(2)}(x_0)}{2!} (x - x_0)^2 + \frac{f^{(3)}(x_0)}{3!} (x - x_0)^3 \end{align*}

I don't see how that would help solve this problem. Thanks!

• Use $$f(x+h) = f(x) + f'(x)h + \frac{f''(x)h^{2}}{2!} + \dots$$ and you want to take it to fourth order i.e the highest degree of the polynomial for $f(x+h)$ is $4$. Nov 25 '17 at 18:19

Consider that your formula must give the correct values for $f (x) = 1$, $f (x) = x$, $f (x) = x^2$, $f (x) = x^3$ and $f (x) = x^4$. The hint is weird, since it assumes that you are not using f (x).

So you get the equations:

$A + B + C + D + E = 0$, $A(x-h) + B(x) + C(x+h) + D(x + 2h) + E(x + 3h) = 1$, $A(x-h)^2 + B(x)^2 + C(x+h)^2 + D(x + 2h)^2 + E(x + 3h)^2 = 2x$, $A(x-h)^3 + B(x)^3 + C(x+h)^3 + D(x + 2h)^3 + E(x + 3h)^3 = 3x^2$, $A(x-h)^4 + B(x)^4 + C(x+h)^4 + D(x + 2h)^4 + E(x + 3h)^4 = 4x^3$.

Consider that your formula must give the correct values for $f (x) = 1$, $f (x) = x$, $f (x) = x^2$, $f (x) = x^3$ and $f (x) = x^4$. The hint is weird, since it assumes that you are not using f (x).

These equations look a bit difficult. So we use different polynomials:

$f (x) = 1$, $f (x) = (x-x_0)/h$, $f (x) = ((x-x_0)/h)^2$, $f (x) = ((x-x_0)/h)^3$ and $f (x) = ((x-x_0)/h)^4$.

(1) $A + B + C + D + E = 0$

(2) $-A + C + 2D + 3E = 1/h$

(3) $A + C + 4D + 9E = 2/h^2$

(4) $-A + C + 8D + 27E = 3/h^3$

(5) $A + C + 16D + 81E = 4/h^4$.

This is easy enough to solve for A, B, C, D and E by subtracting (4) - (2) and (5) - (3), calculating E and D, subtracting (3) - (2) and substituting D, E to get C, adding (2) + (3) and substituting D, E to get A, and substituting all in (1) to get B. To follow the original "hint" set B = 0 instead.

What the hint is hiding is that the final formula is $$f'(x)=A\frac{f(x-h)-f(x)}h+B\frac{f(x+h)-f(x)}h+C\frac{f(x+2h)-f(x)}h+D\frac{f(x+3h)-f(x)}h.$$ Replacing $f$ by its Taylor expansion $$\frac{f(x+kh)-f(x)}h=f'(x)k+\frac12f''(x)k^2h+\frac16f'''(x)k^3h^2+\frac1{24}f^{(4)}(x)k^4h^3+O(h^4)$$ leads to the system $$-A+B+\ 2C+\ 3D=1\\ \ \ A+B+\ 4C+9D=0\\ -A+B+\ 8C+27D=0\\ \ \ A+B+16C+81D=0$$ to get the first derivative with factor $1$ and cancel the coefficients of the 2nd, 3th and 4th derivatives. It does not seem that complicated to solve, starting with the combination of I and III as well as II and IV to $$6C+24D=-1\\ 12C+72D=0\\$$ leading to $C=-6D$, $C=-\frac12$, $D=\frac1{12}$, $-A+B=2-\frac14$, $A+B=2-\frac34$, $B=\frac32$, $A=-\frac14$.

• This isn't close to what I got. How did you get the 4x4 linear system? FYI, I used [1 1 1 1; -1 1 2 3; 1 1 4 9; -1 1 8 27] \ [1; 1; 0; 0] to calculate $A=-5/24, B=7/4, C=-2/3, D=1/8$ in Julia.
– clay
Nov 26 '17 at 4:46
• You were right, it was the wrong system. However, the correction is in the other direction, it needs the equation to cancel the 4th derivatives. Nov 26 '17 at 8:03

Here's my answer to my own question, that I came to after help tips from here and retrying. This doesn't match the two other answers given, so if it's wrong please let me know why.

First basic four term Taylor expansion and remainder term:

\begin{align*} f(x) &= f(x_0) + \frac{f'(x_0)}{1!} (x - x_0) + \frac{f''(x_0)}{2!} (x - x_0)^2 + \frac{f'''(x_0)}{3!} (x - x_0)^3 + \frac{f''''(\xi(x))}{4!} (x - x_0)^4 \\ f(x) &= f(x_0) + f'(x_0) (x - x_0) + \frac{1}{2} f''(x_0) (x - x_0)^2 + \frac{1}{6} f'''(x_0) (x - x_0)^3 + \frac{1}{24} f''''(\xi(x)) (x - x_0)^4 \\ \end{align*}

Apply for $x_0 \pm h, x_0 + 2h, x_0 + 3h$:

\begin{align*} f(x_0 - h) &= f(x_0) - f'(x_0) h + \frac{1}{2} f''(x_0) h^2 - \frac{1}{6} f'''(x_0) h^3 + \frac{1}{24} f''''(\xi(x_0 - h)) h^4 \\ f(x_0 + h) &= f(x_0) + f'(x_0) h + \frac{1}{2} f''(x_0) h^2 + \frac{1}{6} f'''(x_0) h^3 + \frac{1}{24} f''''(\xi(x_0 + h)) h^4 \\ f(x_0 + 2h) &= f(x_0) + 2 f'(x_0) h + 2 f''(x_0) h^2 + \frac{4}{3} f'''(x_0) h^3 + \frac{2}{3} f''''(\xi(x_0 + 2h)) h^4 \\ f(x_0 + 3h) &= f(x_0) + 3 f'(x_0) h + \frac{9}{2} f''(x_0) h^2 + \frac{9}{2} f'''(x_0) h^3 + \frac{27}{8} f''''(\xi(x_0 + 3h)) h^4 \\ \end{align*}

Converting to matrix form and dropping remainder terms:

\begin{align*} \begin{bmatrix} 1 & -1 & \frac{1}{2} & -\frac{1}{6} \\ 1 & 1 & \frac{1}{2} & \frac{1}{6} \\ 1 & 2 & 2 & \frac{4}{3} \\ 1 & 3 & \frac{9}{2} & \frac{9}{2} \\ \end{bmatrix} \begin{bmatrix} f(x_0) \\ h \cdot f'(x_0) \\ h^2 \cdot f''(x_0) \\ h^3 \cdot f'''(x_0) \\ \end{bmatrix} &\approx \begin{bmatrix} f(x_0 - h) \\ f(x_0 + h) \\ f(x_0 + 2h) \\ f(x_0 + 3h) \\ \end{bmatrix} \\ \end{align*}

Inverting the matrix:

Wolfram Alpha:

inverse {{1, -1, 1/2, -1/6}, {1, 1, 1/2, 1/6}, {1, 2, 2, 4/3}, {1, 3, 9/2, 9/2}}


\begin{align*} \begin{bmatrix} \frac{1}{4} & \frac{3}{2} & -1 & \frac{1}{4} \\ -\frac{11}{24} & \frac{1}{4} & \frac{1}{3} & -\frac{1}{8} \\ \frac{1}{2} & -2 & 2 & -\frac{1}{2} \\ -\frac{1}{4} & \frac{3}{2} & -2 & \frac{3}{4} \\ \end{bmatrix} \begin{bmatrix} f(x_0 - h) \\ f(x_0 + h) \\ f(x_0 + 2h) \\ f(x_0 + 3h) \\ \end{bmatrix} &\approx \begin{bmatrix} f(x_0) \\ h \cdot f'(x_0) \\ h^2 \cdot f''(x_0) \\ h^3 \cdot f'''(x_0) \\ \end{bmatrix} \\ \end{align*}

This gives us:

\begin{align*} f'(x_0) &\approx \frac{1}{h} \left[ -\frac{11}{24} f(x_0 - h) + \frac{1}{4} f(x_0 + h) + \frac{1}{3} f(x_0 + 2h) - \frac{1}{8} f(x_0 + 3h) \right] \\ \end{align*}

• You are missing that $f(x_0)$ is also an input. Essentially, what you should be computing is $$f'(x_0)=A\frac{f(x_0-h)-f(x_0)}h+B\frac{f(x_0+h)-f(x_0)}h+C\frac{f(x_0+2h)-f(x_0)}h+D\frac{f(x_0+3h)-f(x_0)}h.$$ Nov 26 '17 at 7:45