Integration

LA home
Computing
 Algorithms
 Bioinformatics
 FP,  λ
 Logic,  π
 MML
 Prog.Langs

Algorithms
 glossary
 Numerical
  Num'Errors
  Polynomials
  Stirling
  Mean&S.D.
  Integration

Integration is the operation of finding the area under a curve, f(x), between given limits, from x=lo to x=hi. If the form of f(x) is suitable, and if we are smart enough, it may be possible to solve the integration analytically to get an exact answer. If this is not possible, f(x) can be integrated numerically to give an approximate answer. Such a computation is known as quadrature. The errors introduced in the process depend upon the properties of f(x) and on the method of integration. Common methods require that f(x) and its derivatives are continuous and that its high-order derivatives are "small".

Some simple methods of numerical integration are illustrated below.

Rectangle Rule

The range [lo,hi] is divided up into N equally sized intervals of width (hi-lo)/N. In a given interval, f(x) is approximated by its value at the centre of an interval, hence rectangles:

function rectangle(f, lo, hi, N)
 { var width = (hi-lo)/N;
   var sum = 0;

   var i;
   for(i=0; i < N; i++)
      sum += f( lo+(i+0.5)*width );
      // f() at centre of i-th interval

   return sum*width;
 }

Trapezoidal Rule

The range [lo,hi] is divided up into N intervals, as before. A straight-line approximation for f(x) is used in each interval, i.e. the area under f(x) is approximated by a series of trapeziums. Note that the area of a trapezium is its width multiplied by the average of the two parallel sides.

C
o
m
p
.
S
c
i
.
function trapezoidal(f, lo, hi, N)
 { var width = (hi-lo)/N;
   var sum = (f(lo)+f(hi))/2;

   var i;
   for(i=1; i < N; i++)
      sum += f( (lo*(N-i) + hi*i)/N );

   return sum*width;
 }

Simpson's Rule

A straight line can be fitted through any two points and a quadratic can be fitted through three points.

  |         |         |
  |      **q***       |
  |  ****   |  ****   |      * ax2+bx+c
  |*        |      ****r
 p*         |         |
  |         |         |
  |         |         |
  |         |         |
  |         |         |
  |         |         |
 -.---------------------------> x
  |         |         |
  0         w        2w

Fitting ax2+bx+c through (0,p), (w,q) & (2w,r) gives us some constraints:

         c = p             -- (0,p)
 aw2+ bw+c = q             -- (w,q)
4aw2+2bw+c = r             -- (2w,r)
We could solve the constraints to find a, b & c, but there is no need to do so directly...

The area under a quadratic is easy to work out:

x=0..2w {ax2+bx+c}   --i.e. area under quadratic
= [x(ax2/3 + bx/2 +c)]0..2w
= 2w(4aw2/3 + bw +c)
= 2w(8aw2 + 6bw + 6c)/6   --Note r+4q = 8aw+6bw+5c
= w*(r + 4q + c)/3
= w(r + 4q + p)/3
So, if the range [lo,hi] is divided into an even number of intervals, giving an odd number of boundaries, the area under f(x) can approximated by a sequence of quadratics. Code for the method can be optimised by realising that one needs a sum of f(x) evaluated at the boundaries of the intervals, weighted by four at odd-numbered internal boundaries and by two at even-numbered internal boundaries.

function Simpson(f, lo, hi, N)
// PRE: N is even
 { var width = (hi-lo)/N;
   var sum = f(lo)+f(hi);

   var i, odd=true;
   for(i=1; i < N; i++)
    { sum += f(lo+i*width) * (odd ? 4 : 2);
      odd = !odd;
    }

   return sum*width/3;
 }

Simpson's method is much more accurate than the rectangle and trapezoidal rules, for a given number of intervals.

Complexity

The rectangle rule evaluates f(x) N times. The trapezoidal rule and Simpson's method evaluate f(x) N+1 times.

Demonstration

The HTML FORM demonstrates the integration of a function, f(x), between x=lo and x=hi in N intervals. In this particular case f(x) is a polynomial so it is easy to do the integration analytically to find the correct answer for comparison with the numerical methods; this is not possible in general. Change f(x), lo, hi and the number of intervals, and experiment:

f(x)=x4+ x3+ x2+ x1+
lo=[  ], hi=[  ], N=[  ]
op=
-- 1999 L.A.
window on the wide world:

Computer Science Education Week

Linux
 Ubuntu
free op. sys.
OpenOffice
free office suite,
ver 3.4+

The GIMP
~ free photoshop
Firefox
web browser
FlashBlock
like it says!

© L. Allison   http://www.allisons.org/ll/   (or as otherwise indicated),
Faculty of Information Technology (Clayton), Monash University, Australia 3800 (6/'05 was School of Computer Science and Software Engineering, Fac. Info. Tech., Monash University,
was Department of Computer Science, Fac. Comp. & Info. Tech., '89 was Department of Computer Science, Fac. Sci., '68-'71 was Department of Information Science, Fac. Sci.)
Created with "vi (Linux + Solaris)",  charset=iso-8859-1,  fetched Saturday, 30-Aug-2014 08:09:25 EST.