Numerical integration

When to Use Numerical Integration?

Numerical integration becomes essential in various scenarios, predominantly when analytical methods are impractical, insufficient, or entirely non-existent. Below are some key circumstances where numerical integration is highly applicable:

  • Non-integrable Functions: In many real-world applications, you'll encounter functions that are not easily integrable in closed-form, making numerical methods the go-to strategy.
  • Computational Efficiency: Even when an analytical solution exists, its complexity might render it computationally inefficient. In such cases, a numerical approximation can often provide sufficient accuracy more quickly.
  • Data-Driven Models: When dealing with empirical data that doesn't conform to an analytical model, numerical integration is an effective way to derive value from such data.
    • Electrical Systems: For instance, you may need to integrate the power consumption curve of a machine over time to calculate total energy usage.
    • Computational Geometry: When you have a complex object represented as a 3D mesh, numerical integration can help in calculating its volume, surface area, or other geometric properties.

It's worth noting that the choice of the numerical integration method may depend on factors such as the required accuracy, the smoothness of the function, and computational resources available.

Method

Numerical integration is a technique used to approximate definite integrals of functions for which analytical solutions are either impractical or impossible to obtain. We start by segmenting the interval [a, b] into $N$ subintervals or 'panels'.

The mathematical framework for this is:

$$ \int_a^b f(x) dx = \sum_{i=0}^{N-1} \int\limits_{x_i}^{x_{i+1}} f(x) dx $$

Here, $x_0 = a$ and $x_N = b$. The width of each panel is given by:

$$ \Delta x = \frac{b-a}{N} $$

Visual representation of panels

If the panels are sufficiently narrow, we can approximate the function $f(x)$ within each panel $i$ using a Taylor series expansion about the point $x_i$:

$$ f(x) = f(x_i) + f'(x_i) \cdot h + \frac{1}{2} f''(x_i) \cdot h^2 + \ldots $$

Here, $h = x - x_i$.

Different numerical integration methods can be derived based on how many terms we include from this Taylor series expansion, thus affecting the degree of precision.

Rectangle Method

The Rectangle Method, also known as the Riemann sum method, is one of the simplest forms of numerical integration. In this approach, we approximate the integral by taking only the first term of the Taylor series expansion:

$$ \int_a^b f(x) dx \approx \Delta x \sum_{i=0}^{N-1} f(x_i) $$

Visual representation of the Rectangle Method

This formula uses the left endpoint of each panel to approximate the integral. We can also use the right endpoint, which would yield:

$$ \int_a^b f(x) dx \approx \Delta x \sum_{i=1}^{N} f(x_i) $$

The accuracy of this method increases with the number of panels $N$, but it is generally less accurate than methods that use higher-order Taylor series approximations.

right rectangles left rectangles

Or we can expand around the middle of the interval:

$$ \int_a^b f(x) dx = \sum_{i=0}^{N-1} \int\limits_{x_i}^{x_{i+1}} f(m_i) dx \approx \frac{b-a}{N}\sum_{i=1}^{N}f(m_i) = \Delta x \sum_{i=1}^{N}f(m_i) $$

where

$$m_i = \frac{x_i+x_{i-1}}{2}$$

Error estimate

In “normal” conditions the error scales with increasing number of panels $N$ as

  • $N^{-1}$ for the left- or right-rectangle method
  • $N^{-2}$ for the midpoint method

Trapezium rule

We can improve the situation by increasing the order of the approximation. If we include the first derivative information we obtain the trapezium rule, which approximates the curve with linear segments:

Derivation

$$ \int_a^b f(x) dx \approx \sum_{i=0}^{N-1} \int\limits_{x_i}^{x_{i+1}} f(x_i)+ (x-x_i)\frac{f(x_{i+1})-f(x_i)}{x_{i+1}-x_1} $$

$$ = \sum_{i=0}^{N-1} \left( \Delta x f(x_i) + \frac{1}{2}\Delta x^2\frac{f(x_{i+1})-f(x_i)}{\Delta x}\right) $$

$$ = \sum_{i=0}^{N-1} \Delta x \frac12 \left( f(x_i)+f(x_{i+1}) \right) =\Delta x\left[\sum_{i=1}^{N-1}f(x_i) +\frac12\left(f(a)+f(b)\right)\right]$$

Cost

  • requires $N+1$ function evaluation
  • one more than the rectangle method
  • error scales like $N^{-2}$

Error estimate

Simpson Rule

  • For each panel we can approximate the function by a quadratic function.
  • we need three points for each panel

Derivation

Ansatz:

$$g(x)=a x^2 + bx +c$$

Solve equations

$ g(0) = f_1 $,

$ g(h) = f_2 $,

$ g(2h) = f_3 $ for $a$, $b$, $c$

Integrate $g(x)$ over the interval $[0,2h]$ gives:

$$ \frac{h}{3}\left(f_1+4f_2+f_3\right) = \Delta x\frac{1}{6}\left(f_1+4f_2+f_3\right). $$

import sympy
x, h, f1, f2, f3 = sympy.symbols("x h f1 f2 f3")
a, b, c = sympy.symbols("a b c")

def g(x):
    return a*x**2 + b*x + c

abcRule = sympy.solve([
    g(0) - f1,
    g(h) - f2,
    g(2 * h) - f3,
],a,b,c)

integ = sympy.integrate(g(x) , (x , 0, 2*h ) )

integ.subs(abcRule).simplify()

$\displaystyle \frac{h \left(f_{1} + 4 f_{2} + f_{3}\right)}{3}$

Now combine all panels:

$$ \int_a^b f(x) dx \approx \frac{ \Delta x}{6} \bigg( f(x_0) + 4 f(m_1) + 2 f(x_1) + 4 f(m_2) + 2 f(x_2) + … + 4 f(m_{N-1})+ 2 f(x_{N-1})+ 4 f(m_{N})+ f(x_N) \bigg)$$

Error estimate

Higher order approximation

  • we can use higher and higher order approximations
  • faster convergence for smooth functions
  • numerical round-off issues

Validity

  • Higher order only helps if the function and its derivatives are smooth

$$f(x) = |x - \pi |\qquad I=\int\limits_0^6 f(x),dx$$

  • Only help if there are no divergences

$$ f(x) = \sqrt{1-x^2},\qquad I=\int\limits_0^1 f(x),dx$$

  • a too high order approximation can be less optimal if the function being integrated is not of high order

$$ f(x)=x^2,,\qquad I=\int\limits_{-1}^2 f(x),dx$$

Comparison table

This table shows the number of function evaluations and the asymptotic behaviour of the integration error (under appropriate smoothness condition) for the methods we have seen in this lecture.

method # function evaluation error
Rectangle $N$ $N^{-1}$
Mid-point $N$ $N^{-2}$
Trapezium $N+1$ $N^{-2}$
Simpson $2N+1$ $N^{-4}$
Quartic $4N+1$ $N^{-6}$
6th order $6N+1$ $N^{-8}$

N.B. This is for single variable integrals, it gets worse as the number of variables increases.

Assignment 2

You should now be able to start the second assignment which investigates numerical integration using Simpson’s rule.