Numerical integration

22. Numerical integration#

Just as we demonstrated numerical differentiation (with np.gradient()) a few lessons ago, now we will work with numerical integration, specifically the trapezoidal method for its efficiency.

Summary of commands#

In this exercise, we will demonstrate the following:

This function comes from the SciPy package and has the signature trapezoid(y, x), where:

  • y is an array of \(y\) values corresponding to \(f(x)\). It is not the analytical expression for \(f(x)\) the way it was for fmin().

  • x is an array of \(x\) values corresponding to each of the values in the y array. This argument is technically optional, but then the default assumes x = [0, 1, 2, ...] (or \(\Delta x = 1\)) which probably will not result in the correct answer!

Integrating functions#

Compute the integral of the function \(f(x) = x^2 + x + 1\) for \(0 \le x \le 1\) using the trapezoidal method. Compare the result to the exact analytical solution.

import numpy as np
from scipy.integrate import trapezoid

# define domain and function
x = np.linspace(0, 1, 1000)
fx = x ** 2 + x + 1

# integrate
res = trapezoid(fx, x)
print(res)
1.833333500333834

The analytical result is

\[ \int_0^1 x^2 + x + 1\ \mathrm{d}x = \dfrac{1}{3}x^3 + \dfrac{1}{2}x^2 + x\ \bigg|_0^1 = \dfrac{11}{6} \approx 1.833 \]