Combining error terms in Simpson's rule

My numerical analysis textbook (Burden and Faires) derives Simpson's rule as

$$\begin{align} \int_{x_0}^{x_2}f(x)\,dx&=2hf(x_1)+\frac{h^3}{3}f''(x_1)+\frac{h^5}{60}f^{(4)}(\xi_1) \\&=\frac{h}{3}[f(x_0)+4f(x_1)+f(x_2)]-\frac{h^5}{12}\left[\frac{1}{3}f^{(4)}(\xi_2)-\frac{1}{5}f^{(4)}(\xi_1)\right] \end{align}$$ for some $\xi_1,\xi_2\in(x_0,x_2)$. Here $x_0=x_1-h$ and $x_2=x_1+h$.

Note that $\xi_1$ comes from the error term when integrating the Taylor series of $f(x)$ while $\xi_2$ comes from the error term in rewriting $f''(x_1)$.

How can we show that we can replace $\xi_1$ and $\xi_2$ in the above formula with some $\xi\in(x_0,x_2)$?

My textbook leaves this part as an exercise, but I am unconvinced the exercise demonstrates what we want. The following is the exercise verbatim:

"Derive Simpson's rule with error term by using $$\int_{x_0}^{x_2}f(x)\,dx=a_0f(x_0)+a_1f(x_1)+a_2f(x_2)+kf^{(4)}(\xi)$$ Find $a_0$, $a_1$, and $a_2$ from the fact that Simpson's rule is exact for $x^n$ when $n=1$, $2$, and $3$. Then find $k$ by applying the integration formula with $f(x)=x^4$."

Performing this exercise, we certainly get the correct coefficients, but I don't see how this tells us we can combine $\xi_1$ and $\xi_2$ - it seems that we just assumed it.

Solution 1:

Let $F$ be an anti-derivative of $f$, $F'=f$. W.l.o.g. $x_1=0$, set $x=h$ then we are interested in the error expression $$ g(x)=F(x)-F(-x)-\frac{x}3(f(x)+4f(0)+f(-x)). $$ This has derivatives \begin{alignat}{2} g'(x)&=\frac23(f(x)-2f(0)+f(-x))&&-\frac x3 (f'(x)-f'(-x)) \\ g''(x)&=\frac13(f'(x)-f'(-x))&&-\frac x3(f''(x)+f''(-x)) \\ g'''(x)&=&&-\frac x3(f'''(x)-f'''(-x)). \end{alignat} Consequently, by the extended mean value theorem $$ \frac{g(x)}{x^m} =\frac{g'(x_1)}{mx_1^{m-1}} =\frac{g''(x_2)}{m(m-1)x_2^{m-2}} =\frac {g'''(x_3)}{m(m-1)(m-2)x_3^{m-3}} $$ with $0<x_3<x_2<x_1<x$. Using $m=5$ this gives $$ \frac{g(x)}{x^5}=\frac{g'''(x_3)}{60x_3^2}=-\frac1{90}·\frac{f'''(x_3)-f'''(-x_3)}{2x_3}=-\frac1{90}·f^{(4)}(x_4). $$ with $|x_4|<x_3<x$. This results in the error formula $$ g(x)=-\frac1{90}·f^{(4)}(x_4)·x^5. $$ or after translating the initial simplifications back, \begin{align} \int_{x_0}^{x_2}f(x)\,dx &=\frac{h}{3}[f(x_0)+4f(x_1)+f(x_2)]-\frac{h^5}{90}f^{(4)}(\xi) \end{align} with $\xi\in(x_0,x_2)$

Solution 2:

Oh well, it looks like my brilliant counterexample is in the shop for repairs. But here's the idea of how we try to make this sort of error estimate, and it's the same idea as goes into proving the error estimate for the Taylor series. Simpson's rule is based on a polynomial $p(x)$ that matches the target function $f(x)$ at three points. Because we can always make the interval of integration $[-1,1]$ by a linear transformation, it is conventional to assume this has been done so we have $$y(x_0)=y(-1)=p(x_0)$$ $$y(x_1)=y(0)=p(x_1)$$ $$y(x_2)=y(1)=p(x_2)$$ We can construct a polynomial $q(x)$ of degree at most $3$ which also matches $y(x)$ at these three points and also at $x_3$ where $x_3\in(x_0,x_1)\cup(x_1,x_2)$ as follows: $$q(x)=p(x)+(y(x_3)-p(x_3))\prod_{k=0}^2\frac{(x-x_k)}{(x_3-x_k)}$$ Just evaluate $q(x)$ at each of the four points to assure yourself that this is so. Now consider the error in that approximation $$e(x)=y(x)-q(x)$$ We just now checked that $e(x)=0$ at $x_0$, $x_2$, and at two more points in $(x_0,x_2)$. By $3$ applications of Rolle's theorem we find that there is at least one point in $(x_0,x_2)$ such that $e^{\prime\prime\prime}(\xi)=y^{\prime\prime\prime}(\xi)-q^{\prime\prime\prime}(\xi)=0$. That is $$0=y^{\prime\prime\prime}(\xi)-3!\frac{y(x_3)-p(x_3)}{\prod_{k=0}^2(x_3-x_k)}$$ Since $x_3$ could have been any point in $[x_0,x_2]$ we have shown that $$y(x)=p(x)+\frac1{3!}y^{\prime\prime\prime}(\xi(x))\prod_{k=0}^2(x-x_k)$$ although we have to double check for $x\in\{x_0,x_1,x_2\}$. Notice that in general the point $\xi$ varies depending on $x$. But if we try to use this estimate to also estimate the error in Simpson's rule there is a complication. $$\int_{-1}^1y(x)dx=\int_{-1}^1p(x)dx+\frac1{3!}\int_{-1}^1y^{\prime\prime\prime}(\xi(x))\prod_{k=0}^2(x-x_k)dx$$ And at this point we would like to say something like 'the smallest the error term could be is' $$\frac1{3!}\max\left(y^{\prime\prime\prime}(\xi)\right)\int_{-1}^1\prod_{k=0}^2(x-x_k)dx$$ But we can't do that because the integrand changes sign in the interval so that if the wiggles in $y^{\prime\prime\prime}(\xi(x))$ line up with the wiggles in the integrand that could make the integral even smaller. The problem of the integrand in the error estimate changing sign in the interval of integration can be fixed by revamping our derivation of Simpson's rule from the start. We will attempt $4$-point osculation with a cubic polynomial $q(x)$ such that $$y(-1)=y(x_0)=q(-1)$$ $$y(0)=y(x_1)=q(0)$$ $$y^{\prime}(0)=y^{\prime}(x_1)=q^{\prime}(0)$$ $$y(1)=y(x_2)=q(1)$$ This time we construct a polynomial $r(x)$ of degree at most $4$ such that it satisfies the above conditions and also $r(x_3)=y(x_3)$ at a point $x_3\in(x_0,x_1)\cup(x_1,x_2)$. $$r(x)=q(x)+(y(x_3)-q(x_3))\frac{(x-x_0)(x-x_1)^2(x-x_2)}{(x_3-x_x)(x_3-x_1)^2(x_3-x_2)}$$ Again it can be checked to satisfy the $5$ conditions. Now the error in this approximation is $e(x)=y(x)-r(x)$ and has $4$ zeros in $[x_0,x_2]$. but this time $e^{\prime}(x)$ also has $4$ zeros in $(x_0,x_2)$: $3$ it gets from Rolle's theorem, and one more for the double root at $x=x_1$. So after $3$ more applications of Rolle's theorem we arrive at $$e^{(4)}(\xi)=0=y^{(4)}(\xi)-\left(y(x_3)-q(x_3)\right)\frac{4!}{(x_3-x_0)(x_3-x_1)^2(x_3-x_2)}$$ For some $\xi\in(x_0,x_2)$. So we have established $$y(x)=q(x)+\frac1{4!}(x-x_0)(x-x_1)^2(x-x_2)y^{(4)}(\xi(x))$$ And when we progress to $$\int_{-1}^1y(x)dx=\int_{-1}^1q(x)dx+\frac1{4!}\int_{-1}^1(x-x_0)(x-x_1)^2(x-x_2)y^{(4)}(\xi(x))dx$$ We really can say that $$\begin{align}-\frac4{15}\frac1{4!}\max\left(y^{(4)}(\xi)\right)&\le\frac1{4!}\int_{-1}^1(x-x_0)(x-x_1)^2(x-x_2)y^{(4)}(\xi(x))dx\\ &\le-\frac4{15}\frac1{4!}\min\left(y^{(4)}(\xi)\right)\end{align}$$ Because the integrand now has the same sign in $(x_0,x_2)$ and $$\int_{-1}^1(x-x_0)(x-x_1)^2(x-x_2)dx=\int_{-1}^1(x^4-x^2)dx=\frac25-\frac23=-\frac4{15}$$ By the intermediate value theorem we can conclude that $$\frac1{4!}\int_{-1}^1(x-x_0)(x-x_1)^2(x-x_2)y^{(4)}(\xi(x))dx= -\frac4{15}\frac1{4!}y^{(4)}(\xi)$$ for some $\xi\in(x_0,x_2)$ and then conclude that $$\int_{-1}^1y(x)dx=\int_{-1}^1q(x)dx-\frac1{90}y^{(4)}(\xi)$$ We just have to find $q(x)$ and integrate. Let $q(x)=a+bx+c(1-x^2)+d(x-x^3)$. Then $$q(0)=y(0)=a+c$$ $$q(-1)=y(-1)=a-b$$ $$q(1)=y(1)=a+b$$ $$q^{\prime}(0)=y^{\prime}(0)=b+d$$ It's easy to solve these for $$a=\frac{y(1)+y(-1)}2$$ $$b=\frac{y(1)-y(-1)}2$$ $$c=\frac{-y(1)+2y(0)-y(-1)}2$$ $$d=\frac{-y(1)+2y^{\prime}(0)+y(-1)}2$$ And $$\int_{-1}^1q(x)dx=2a+\frac43c=\frac13y(-1)+\frac43y(0)+\frac13y(1)$$ So even though we included $y^{\prime}(0)$ in our derivation its value was ignored in the eventual formula. This made it easier to get the error estimate for Simpson's rule and relied on the fact that the rule is really a Gauss-Lobatto formula that can actually ignore the values of the first derivatives of the function being approximated at all interior points. The fact that the error is not a bound, but actually attained for some $\xi\in(x_0,x_2)$ is useful for estimating the error in the composite Simpson's rule.

Hoping that this discussion at least gives you some insight on the leap of faith you are being asked to make as described in your question.