pylib.numerical.integration module

Numerical integration, numerical quadrature.

de: numerische Integration, numerische Quadratur.

Date

2015-10-15

trapez(f, a=0, b=1, N=10, x=None, verbose=False, save_values=False)

Integration of \(f(x)\) using the trapezoidal rule (Simpson’s rule, Kepler’s rule).

de: Trapezregel, Simpsonregel (Thomas Simpson), Keplersche Fassregel (Johannes Kepler)

Parameters

• f (function or list) – function to integrate.

• a (float) – lower limit of integration (default = 0).

• b (float) – upper limit of integration (default = 1).

• N (int) – specify the number of subintervals.

• x (list) – variable of integration, necessary if f is a list (default = None).

• verbose (bool) – print information (default = False)

Returns

the definite integral as approximated by trapezoidal rule.

Return type

float

The trapezoidal rule approximates the integral by the area of a trapezoid with base h=b-a and sides equal to the values of the integrand at the two end points.

\[f_n(x) = f(a)+\frac{f(b)-f(a)}{b-a}(x-a)\]

\[\begin{split}I &= \int\limits_a^b f(x) \,\mathrm{d}x \\ I &\approx \int\limits_a^b f_n(x) \,\mathrm{d}x \\ &= \int\limits_a^b \left( f(a)+\frac{f(b)-f(a)}{b-a}(x-a) \right) \mathrm{d}x \\ &= \left.\left( f(a)-a\frac{f(b)-f(a)}{b-a} \right) x \right\vert_a^b + \left. \frac{f(b)-f(a)}{b-a} \frac{x^2}{2} \right\vert_a^b \\ &= \frac{b-a}{2}\left[f(a)+f(b)\right]\end{split}\]

The composite trapezium rule. If the interval is divided into n segments (not necessarily equal)

\[a = x_0 \leq x_1 \leq x_2 \leq \ldots \leq x_n = b\]

\[\begin{split}I &\approx \sum\limits_{i=0}^{n-1} \frac{1}{2} (x_{i+1}-x_i) \left[f(x_{i+1})+f(x_i)\right] \\\end{split}\]

Special Case (Equaliy spaced base points)

\[x_{i+1}-x_i = h \quad \forall i\]

\[I \approx h \left\{ \frac{1}{2} \left[f(x_0)+f(x_n)\right] + \sum\limits_{i=1}^{n-1} f(x_i) \right\}\]

Example

\[\begin{split}I &= \int\limits_a^b f(x) \,\mathrm{d}x \\ f(x) &= x^2 \\ a &= 0 \\ b &= 1\end{split}\]

analytical solution

\[I = \int\limits_{0}^{1} x^2 \,\mathrm{d}x = \left. \frac{1}{3} x^3 \right\vert_0^1 = \frac{1}{3}\]

numerical solution

>>> f = lambda(x): x**2
>>> trapez(f, 0, 1, 1)
0.5
>>> trapez(f, 0, 1, 10)
0.3350000000000001
>>> trapez(f, 0, 1, 100)
0.33335000000000004