Numerical Integration (Quadrature)

NRiC Chapter 4.

Already seen Monte Carlo integration.

Can cast problem as a differential equation (DE):

is equivalent to solving for I ≡ y(b) the DE dy/dx = f(x) with the boundary condition (BC) y(a) = 0.

Will learn about ODE solution methods next class.

Trapezoidal & Simpson's Rule

Have abscissas xi = x0 + ih, i = 0, 1, ..., N + 1.

A function f(x) has known values f(xi) = fi.

Want to integrate f(x) between endpoints a & b.

Trapezoidal rule (2-point closed formula):

Simpson's rule (3-point closed formula):

Extended Trapezoidal Rule

If we apply the Trapezoidal rule N -1 times and add the results, we get:

Big advantage is it builds on previous work:

Coarsest step: average f at endpoints a and b.

Next refinement: add value at midpoint to average.

Next: add values at ¼ and ¾ points.

And so on. This is implemented as trapzd() in NRiC.

More Sophistication

Usually don't know N in advance, so iterate to a desired accuracy: qtrap().

Higher-order method by cleverly adding refinements to cancel error terms: qsimp().

Generalization to order 2k (Richardson's deferred approach to the limit): qromb().

For improper integrals, generally need open formulae (not evaluated at endpoints).

For multi-D, use nested 1-D techniques.

2001-02-15T09:13:37

Derek Richardson