Numerical Methods: The Trapezium Rule and Simpson's Rule

Integrals don't have to get very complicated before symbolic methods fail to work. Consider, for example, the integral

01cos(x3+x)dx: 

there are no know symbolic methods, based on indefinite integration, that can be brought to bear on this problem.

Instead, we can sample the integrand at regular integrals and carry out an estimate based on this. One way of doing that is to approximate the function by a sequence of straight line segments, as in the figure. The area between each segment and the x-axis is a trapezium, meaning that if the width of the interval is h, and the y-values at each end of the interval are yi and yi+1, then the area of the trapezium is

2h(yi+yi+1)

The entire area between the curve and the x-axis, which is to say the integral, can be approximated by adding together several such trapezia. If there are n trapezia, and n+1 y-values (ordinates) running from y0 to yn, then the integral is approximately

Tn=2h(y0+2y2+2y2++2yn−2+2yn−1+yn)

This is called the trapezium rule.

The trapezium rule with four ordinates

For example, let us estimate

01cos(x3+x)dx 

using five ordinates (four intervals). The interval width, h, is 025. From the table,

T4=2025(058385+2216432)=06141

An alternative, which is often (though not always) more accurate, is to approximate the curve with a sequence of quadratic parabolic segments instead of straight lines. Each parabola requires three points to specify it, so each parabola spans two intervals. This method can only be used, therefore, if the number of intervals is even (and the number of ordinates, therefore, odd).

Figure 2: Approximating a curve (red) with a quadratic (blue): the basis of Simpson's rule

Now, consider the quadratic passing through the points (−hy0), (0y1) and (hy2). It's not hard to show that this has the equation

y=2h2y0−2y1+y2x2+2hy2−y0x+y1 

If we integrate this between −h and h we get

3h(y0+4y1+y2)

This would also have been the area if we'd situated the middle of the interval somewhere other than zero (though the calculations would have been harder). So from our quadratic approximation we get the estimate

Sn=3h(y0+4y2+2y2+4y3++4yn−3+2yn−2+4yn−1+yn)

where n must be even. This is called Simpson's rule.

From the table,

S4=3025(058395+4135335+2081096)=06349