数值积分中的辛普森方法及其误差估计

$f$在$(x_0,x_2)$上四阶可导，且在$[x_0,x_2]$上三阶导函数连续.则
\begin{equation}
\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{equation}其中$h=x_1-x_0=x_2-x_1$.

证明:我们仍然进行牛顿插值.设牛顿插值的插值点为
\begin{equation}
x_0,x_1,x_2
\end{equation}
则经过这三个插值点的牛顿插值多项式为
\begin{equation}
f(x_0)+(x-x_0)f[x_0,x_1]+(x-x_0)(x-x_1)f[x_0,x_1,x_2]
\end{equation}
插值余项为
\begin{equation}
(x-x_0)(x-x_1)(x-x_2)f[x,x_0,x_1,x_2]
\end{equation}
我们先来看
\begin{equation}
\int_{x_0}^{x_2}\{f(x_0)+(x-x_0)f[x_0,x_1]+(x-x_0)(x-x_1)f[x_0,x_1,x_2]\}dx
\end{equation}
我们先来算
\begin{equation}
\int_{x_0}^{x_2}f(x_0)dx=2hf(x_0)
\end{equation}
\begin{align*}
\int_{x_0}^{x_2}(x-x_0)f[x_0,x_1]dx=\frac{1}{2}(x_2-x_0)^2f[x_0,x_1]=\frac{1}{2}(2h)^2 \frac{f(x_1)-f(x_0)}{h}=2h(f(x_1)-f(x_0))
\end{align*}
\begin{align*}
\int_{x_0}^{x_2}(x-x_0)(x-x_1)f[x_0,x_1,x_2]dx&=f[x_0,x_1,x_2]\int_{x_0}^{x_2}(x-x_0)(x-x_1)dx\\&=f[x_0,x_1,x_2]\int_{x_0}^{x_2}(x^2-(x_0+x_1)x+x_0x_1)dx\\&=f[x_0,x_1,x_2][(\frac{1}{3}x_2^3-\frac{1}{2}(x_0+x_1)x_2^2+x_0x_1x_2)-(\frac{1}{3}x_0^3-\frac{1}{2}(x_0+x_1)x_0^2+x_0^2x_1)]\\&=\frac{1}{6}f[x_0,x_1,x_2](x_2-x_0)^2(2x_2+x_0-3x_1)
\end{align*}
而
\begin{align*}
f[x_0,x_1,x_2]=\frac{f(x_0)}{2h^2}+\frac{f(x_1)}{-h^2}+\frac{f(x_2)}{2h^2}
\end{align*}
可见，
\begin{align*}
\int_{x_0}^{x_2}(x-x_0)(x-x_1)f[x_0,x_1,x_2]dx=\frac{[f(x_0)-2f(x_1)+f(x_2)]h}{3}
\end{align*}
可见，
\begin{align*}
\int_{x_0}^{x_2}\{f(x_0)+(x-x_0)f[x_0,x_1]+(x-x_0)(x-x_1)f[x_0,x_1,x_2]\}dx=\frac{(f(x_0)+4f(x_1)+f(x_2))h}{3}
\end{align*}
下面我来证明
\begin{equation}\label{eq:12345}
\int_{x_0}^{x_2} (x-x_0)(x-x_1)(x-x_2)f[x,x_0,x_1,x_2]dx=\frac{-h^5}{90}f^{(4)}(\xi)
\end{equation}
令
\begin{align*}
g(t)=\int_{x_0}^t(x-x_0)(x-x_1)(x-x_2)f[x,x_0,x_1,x_2]dx
\end{align*}

则
\begin{align*}
g'(t)=(t-x_0)(t-x_1)(t-x_2)f[t,x_0,x_1,x_2]
\end{align*}
则
\begin{align*}
g'(x_0)=g'(x_1)=g'(x_2)=0,g(x_0)=0
\end{align*}
我们进行Hermite插值.由于Hermite插值是牛顿插值的极限情形，为此我们先进行牛顿插值.我们设立插值点
\begin{align*}
x_0,x_0',x_1,x_1',x_2,x_2'
\end{align*}
经过这几个插值点的牛顿插值多项式为
\begin{align*}
p(x)= g(x_0)&+(x-x_0)g[x_0,x_0']\\&+(x-x_0)(x-x_0')g[x_0,x_0',x_1]\\&+(x-x_0)(x-x_0')(x-x_1)g[x_0,x_0',x_1,x_1']\\&+(x-x_0)(x-x_0')(x-x_1)(x-x_1')g[x_0,x_0',x_1,x_1',x_2]\\&+(x-x_0)(x-x_0')(x-x_1)(x-x_1')(x-x_2)g[x_0,x_0',x_1,x_1',x_2,x_2']
\end{align