Integrals are Easy: Visualized Riemann Integration in Python
The integral is not so complicated as it seems to be. It is one of the fundamental and universal tools in mathematics allowing us to calculate the area or the volume of any arbitrary body. It is one of the cornerstones of mathematics having a multitude of applications in many disciplines.
The Riemann Integral is the simplest form of integration, yet it lays down the foundation of all other types of integrals. It offers a rigorous method for approximating the area under the curve of some function \( f \) over some interval \( [a, b] \). This fact assigns to it an intuitive geometrical interpretation.
In this blog post we will introduce and elaborate more on the Riemann Integration. We will start with intuitive reasoning on the process of integration in order to have a smooth transition towards the mathematical foundations. Then, we will see how to transform the theory into an easy Python implementation. Finally, we visually enhance and complement our understanding with an animated visualization of the Riemann Sums using the Matplotlib's Animation API.
Intuition
We can easily calculate the area of any regular-shaped bodies like the rectangle because it consists of only straight lines. Thus, for a rectangle with width \( m \) and height \( n \), the area is calculated as simple as \(m \times n \). It is straightforward to deduct this because for every small unit on the side \( m \) there is still a regular rectangle with height \( n \).
However, if we only change the top side from a straight line to some arbitrary line that can be described with some non-linear function, the circumstances get complicated. We can still divide the surface in small rectangles, but now, they have varying height and on top of this, they do not entirely fit inside the body. This is illustrated in the image below:
Now, to reduce the calculation error we would need to fit infinitely many rectangles inside the irregular body. This leads us to the definition of the Riemann Integral which has exactly the same geometrical motivation and interpretation. We only need to give a more rigorous definition of this procedure, which we do in the next section.
Riemann Integration Definition
To formally define the Riemann Integral, we start with some real function \( f: [a, b] \rightarrow \mathbb{R} \) which is non-negative (it includes zero values) and continuous over the interval \( [a, b] \). This is our arbitrary top line in the example above depicted in Figure 1. To complete the missing parts, we need to define the widths and heights of the mini rectangles.
To define the rectangles' widths, we make a partition of the interval \( [a, b] \). That means we divide the interval \( [a, b] \) into \( N \) sub-intervals for some \( N \in \mathbb{N} \), i.e.
For simplicity reasons we can consider all sub-intervals equidistant, although in the general case they can take any length. Thus, the width of any rectangle would be \( w = (b - a) \div N \), such that for any \( i \in [1, N] \) it holds that \( x_{i} - x_{i-1} = w = (b - a) \div N \).
To define the rectangles' heights, we simply chose a random point from each sub-interval, i.e. \( x_{i}^{*} \in [x_{i - 1}, x_{i}] \) for any \( i \in [1, N] \). To make things easier, we select this point to be the left-most one, i.e. \( x_{i}^{*} = x_{i-1} \). Thus, the height of each mini rectangle is \( f(x_{i - 1}) \) for \( i \in [1, N] \).
Having said all of this, the formal and simplified definition of the Riemann Integral is as follows:
However, this simple definition is too lossy and it would need a large \( N \) to converge properly. To escape this limitation, we make a simple trick: transforming the mini rectangles to mini right trapezoids. The right trapezoids fit better under the curve, accounting for the loss, as depicted in the figure below:
With this in mind, we only need to slightly modify our formal Riemann Integral definition. We switch from an area of a rectangle to an area of a right trapezoid. In literature, this is referred to as Trapezoidal Rule. In this case, the sum would be:
Numerical Integration in Python
We only need to translate the last equation into a Python set of instructions. Thus, the Python implementation is a piece of cake as given below:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 | def calculate_integral(f, a, b, n): '''Calculates the integral based on the composite trapezoidal rule relying on the Riemann Sums. :param function f: the integrand function :param int a: lower bound of the integral :param int b: upper bound of theintergal :param int n: number of trapezoids of equal width :return float: the integral of the function f between a and b ''' w = (b - a)/n result = 0.5*f(a) + sum([f(a + i*w) for i in range(1, n)]) + 0.5*f(b) result *= w return result |
Once we have the implementation, it is necessary to test it against some universal mathematical truth. For instance, it is well known and we can mathematically calculate that:
To test the convergence of our numerical integration implementation, we calculate the absolute difference between the exact and approximated value of \( \pi \). Therefore, at the same time we try to approximate \( \pi \) and test our implementation.
Moreover, to enhance the perception of this approximation it is necessary to show a geometrical and visual interpretation of the process. For this reason, we make an animated visualization using Matplotlib's Animation API. We make the following observation: as the number of trapezoids \( N \) increase, the approximation error decrease. On top of that, we see how the number of trapezoids geometrically reflects in the calculation of the integral. The animation is shown below:
We can notice that for a fairly small number of trapezoids, i.e. 200 in total, the approximation error is already in an order of magnitude of \( 10^{-5} \). That means our implementation is correct, although we can apply additional error bound analysis.
The full source code related to all we have discussed during this blog can be found on GitHub. For more information, please follow me on Twitter.
If you liked what you just saw, it would be really helpful to subscribe to the mailing list below. You will not get spammed that's a promise! You will get updates for the newest blog posts and visualizations from time to time.
Conclusion
In this short blog post, we extended our general geometry knowledge to calculate the area of some irregular-shaped bodies. We achieved this with the simplest form of integration, the Riemann Sums, for which we gave a formal definition. Later on, we provided a straightforward Python implementation and an animated visualization of the integration process using Matplotlib's Animation API.
The Riemann Integral is one simple but yet powerful tool to calculate the area under the curve. However, the fact that we fit mini rectangles or trapezoids inside the area is quite limiting. More generally, the body can have any irregular shape for which we need other methods like the Stieltjes or Lebesgue integrals.
References
[1] Svein Linge, Hans Petter Langtangen,
"Programming for
Computations - Python" (2016), Springer Open
Leave a comment