The

Trapezoidal rule is a mathematical method for doing

numerical
integration ("calculating the

area under a curve"). This is often
useful when an

*exact* integral does not exist, can not easily
be obtained, or is mathematically too time consuming for repetitious
automated calculations.

One approach to obtain a numerical solution of an integral is to
approximate the function with an `n`^{th}
order polynomial, since these are relatively simple to integrate. The choice of the
order of the polynomial depends on the required accuracy (higher order
generally results in a higher precision), and the number data points over
the selected interval (higher order requires more data points).

Consider the function f(x). We want to calculate the area under the
curve over the interval `a` ≤ x ≤ b:

`b`
I = ∫ y(x)dx
`a`

The *simplest* polynomial approximation, a 0^{th} order polynomial
would draw a straight horizontal line, halfway between y(`a`) and y(`b`),
marking a rectangle. The area of this rectangle is simply the length
times the height. But that is not a very accurate method; the
approximation of the area yields a relatively large error.

A *better* approximation is to approximate the curve with a
1^{st} order polynomial; a line that goes trough the points
`a` and `b`. This can be seen in the following figure:

y(x)
y(b)_| ..
| ./|C
| ./ |
| ./ |
| / | ..... y(x)
| /. |
| / . |
| / . |
| . /. |
y(a)_| . D| |
| . | |
| .. | |
|. | |
| | |
|_______A|________|B___ x
`a` `b`

Over the range

`a` ≤ x ≤

`b`, the function y(x) is
now approximated with a

trapezoid marked by ABCD. The area of the
trapezoid is:

ABCD = 1/2 AB(AD + BC) = 1/2 (b-a)(y(`a`) + y(`b`))

We can gain a higher precision if we divide the interval `a`
≤ x ≤ `b` up into several equal parts, `a`
= `x`_{0}, `x`_{1},
`x`_{2}, ..., `x`_{n-1},
`x`_{n} = `b`. We find the corresponding y
values, `y`_{0} = y(`x`_{0}),
`y`_{1} = y(`x`_{1}),
`y`_{2} = y(`x`_{2}), ...,
`y`_{n} = y(`x`_{n}). Each part has a
width equal to h = (`b` - `a`)/n

The Trapezoidal Rule for approximating the integral now becomes:

`b`
∫ y(x)dx = h/2(`y`_{0} + 2`y`_{1} + 2`y`_{2} + ... + 2`y`_{n-1} + `y`_{n})
`a`

**Example:** Let's approximate the integral of the function:

y = x^{2} + 2x -3

over the range `3` ≤ x ≤ `4`. The *exact* solution is straightforward:

4 4
∫ y(x)dx = [ 1/3 x^{3} + x^{2} - 3x] = 76/3 - 9 = 16.3333...
3 3

Now let's try the numerical solution. First generate a table with x-, and y-values for the selected interval:

n x_{n} y
0 3.0 12.00
1 3.2 13.64
2 3.4 15.36
3 3.6 17.16
4 3.8 19.04
5 4.0 21.00

The numerical solution is given by (h=1/5):

4
∫ y(x)dx = (1/10)(12.00 + 2×13.64 + 2×15.36 + 2×17.16 + 2×19.04 + 21.00) = 16.34
3

The example shows that the *trapezoidal rule* can obtain
relatively accurate results (in this case an error of approximately 0.04%) with only a
small number of simple calculations. The error between the exact integral (I) and the numerical solution I_{n} is given by:

I - I_{n} = -1/12 h^{3}y"(ε)

with `a` ≤ ε ≤ `b`

Of course, the integral can also be calculated with higher order polynomial approximations. The 2^{nd} order polynomial approximation is called Simpson's Rule.