Page 113 -

Simpson’s Algorithm: We shall assume that the interval of integration is sam-
pled at an odd number of points (2N + 1), so that we have an even number of
intervals. The algorithm groups the intervals in pairs.
Referring to Figure 4.3, the points A, H, and G are the first three points in
the sampled x-interval. The assumption underlying Simpson’s rule is that the
curve passing through the points B, D, and F, on the curve of the integrand,
can have their locations approximated by a parabola. The line CDE is tangent
to this parabola at the point D.
Under the above approximation, the value of the integral of the y-function
between the points A and G is then simply the sum of the area of the trape-
zoid ABFG plus 2/3 the area of the parallelogram BCEF, namely:

4 ∆x   y()1 + y( )3 
Area of the first two slices = ∆xy(( )1 + y( ))3 +  y( )2 − 
3   2 
(4.21)
= ∆x y (( )1 + y( )4 2 + y( ))3
3

In a similar fashion, we can find the area of the third and fourth slices,

Area of the third and fourth slices = ∆x (( + y4 4( ) + y 5( )) (4.22)
y 3)
3
Continuing for each successive pair of slices, we obtain for the total integral,
or total area of all slices, the expression:

Total area of all slices = ∆x  y 1() + y4 2( ) + y2 3( ) + y4 4( ) + y2 5( ) +…  (4.23)

1)
3  … +…+ yN42( ) + yN2( + 
that is, the weights are equal to 1 for the first and last elements, equal to 4 for
even elements, and equal to 2 for odd elements.

Example 4.7
Using Simpson’s rule, compute the integral of sin(x) over the interval 0 ≤ x ≤π.

Solution: Edit and execute the following script M-file:

a=0;b=pi;N=4;
x=linspace(a,b,2*N+1);
y=sin(x);
for k=1:2*N+1
if k==1 | k==2*N+1
w(k)=1;

© 2001 by CRC Press LLC