Simpson's rule (source code)

= Simpson's rule
{c}
{title2=$\int_{a}^{a+2h}f\approx\frac h3[f(a)+4f(a+h)+f(a+2h)]$}
{wiki}

<Simpson's rule> approximates a <definite integral> by integrating the <quadratic polynomial> obtained by <polynomial interpolation> at two endpoints and their <midpoint>. For an <interval> $[a,b]$ with $a<b$, write $m=(a+b)/2$ and $h=(b-a)/2$. The three <quadrature nodes> are $a,m,b$, and the <quadrature rule> is
$$
S(f)=\frac h3\bigl[f(a)+4f(m)+f(b)\bigr].
$$
Thus its <quadrature weights> are $h/3,4h/3,h/3$. Their sum is $2h=b-a$, as required for exactness on <constant functions>. This is the three-node <Newton-Cotes closed quadrature>; its <degree of precision of a quadrature rule> is three, even though the interpolating <polynomial> has <polynomial degree> at most two.<NIST Digital Library of Mathematical Functions/3.5 Quadrature, subsection (ii) Simpson’s Rule>[NIST Digital Library of Mathematical Functions; §3.5(ii), equation 3.5.6: elementary Simpson rule and its fourth-derivative remainder]{ref}

For an original worked example, take the <cubic polynomial> $f(x)=x^3-2x+4$ on the <interval> $[-1,2]$. Here $m=1/2$ and $h=3/2$, so the <function evaluations> give
$$
f(-1)=5,\qquad f(1/2)=\frac{25}{8}=3.125,\qquad f(2)=8.
$$
<Simpson's rule> gives
$$
S(f)=\frac12\left(5+4\cdot\frac{25}{8}+8\right)=\frac{51}{4}=12.75.
$$
An <antiderivative> supplies an independent exact <definite integral>:
$$
\int_{-1}^{2}(x^3-2x+4)\,dx
=\left[\frac{x^4}{4}-x^2+4x\right]_{-1}^{2}
=\frac{51}{4}.
$$
The <quadratic polynomial> through the same <quadrature nodes> is $p_2(x)=\tfrac32x^2-\tfrac12x+3$. It differs from $f$ between the <quadrature nodes>, but its <definite integral> agrees exactly. The picture shows why: the signed <polynomial interpolation error> has two lobes with opposite <integrals>.

\Image[/simpson-cubic-error-cancellation.png]
{title=An original cubic and its quadratic interpolant have the same integral because the two interpolation-error lobes have signed integrals plus and minus 81 over 64}

Three <function evaluations> therefore capture more information about an <integral> than about the whole <function>. When an <antiderivative> is unavailable, <composite Simpson's rule> repeats this inexpensive local construction over a <partition of an interval>. For a sufficiently smooth <function>, halving the node spacing reduces its <truncation error> by approximately a factor of sixteen.