MATLAB - Trapezoidal Rule
The trapezoidal rule is a numerical method for approximating the definite integral of a function. It works by dividing the area under the curve into a series of trapezoids rather than rectangles (as in the rectangular rule), and then summing their areas to get an estimate of the total integral.
Step-by-Step Details of the Trapezoidal Rule
Step 1 − Divide the Interval [a, b]
The first step in applying the Trapezoidal Rule is to divide the interval [a,b] over which you want to integrate the function f(x) into n subintervals of equal width x.
The width x of each subinterval is calculated as,
$$\mathrm{\Delta x\: = \: \frac{b \: - \: a}{n}}$$
This division results in n+1 points x0,x1,x2,,xn where
$$\mathrm{x_{i} \: = \: a \: + \: i\cdot\Delta x \:\:for\: i \: = \: 0,1,2,\dots, n }$$
Here, x0 = a and xn = b, representing the endpoints of the interval.
Step 2 − Forming the Trapezoids
$$\mathrm{\begin{bmatrix} x_{i - 1}, & x_{i} \end{bmatrix}}$$
For each subinterval a trapezoid is formed using the function values at the endpoints of the subinterval,
$$\mathrm{f(x_{i-1}) \: and \: f(x_{i})}$$
$$\mathrm{\begin{bmatrix} x_{i - 1}, & x_{i} \end{bmatrix}}$$
The area of a trapezoid with bases f(xi-1) and f(xi) and height x is given by
$$\mathrm{Area \: of \: Trapezoid_{i} \: = \: \frac{1}{2} \: \times \: \Delta x \: \times \: [f(x_{i-1}) \: + \: f(x_{i})]}$$
The total approximate area under the curve, which represents the integral, is the sum of the areas of all trapezoids.
Step 3 − Summing the Areas of Trapezoids
The integral of the function f(x) over the interval [a,b] is approximated by summing the areas of all trapezoids
$$\mathrm{\displaystyle\int_{a}^{b} f(x)\: dx\: \approx \: \displaystyle\sum\limits_{i=1}^n \: Area \: of \: Trapezoid_{i}}$$
Substituting the area formula, this becomes −
$$\mathrm{\displaystyle\int_{a}^{b} f(x)\: dx\: \approx \: \frac{\Delta x}{2} \left[f(x_{0}) \: + \: 2 \displaystyle\sum\limits_{i=1}^{n-1} f(x_{i}) \: + \: f(x_{n})\right]}$$
This equation gives a weighted sum of the function values, where the endpoints f(x0) and f(xn) are each counted once , while all intermediate function values f(x1), f(x2), ..f(xn-1) are counted twice.
Now let us see how we can implement trapezoidal rule in matlab.
The Trapezoidal Rule can be implemented in MATLAB using either a manual approach with loops or by utilizing the built-in trapz function.
Syntax
Q = trapz(Y) Q = trapz(X,Y) Q = trapz(___,dim)
Syntax Explanation
trapz(Y) − For vector Y, it calculates the approximate integral of the values in Y using the trapezoidal rule, assuming the points are evenly spaced.
For a matrix Y , it calculates the integral for each column separately, giving you a row vector with the results. For multidimensional array Y, it calculates the integral along the first dimension that has more than one element. After the calculation, this dimension becomes a single value, while the other dimensions stay the same.
trapz(X,Y) − It calculates the integral of Y using the trapezoidal rule, with respect to the values or spacing given by X.
If X is a vector of coordinates: The number of elements in X must match the size of the first dimension of Y that has more than one element.
If X is a single number (scalar), trapz(X,Y) is the same as multiplying that number X by the result of trapz(Y), which assumes uniform spacing.
trapz(___,dim) − For Q = trapz(___,dim): It calculates the integral along the specific dimension dim in Y using the trapezoidal rule.You must specify Y, and you can also optionally specify X. If you provide X, it can either be a single number (scalar) or a vector whose length matches the size of Y along the specified dimension dim.For example: If Y is a matrix and you use trapz(X,Y,2), it will integrate across each row of Y.
Example 1: Using trapz(Y) where Y is a vector
Y = [1, 4, 9, 16, 25]; Q = trapz(Y); disp(Q);
The vector represents a series of values at evenly spaced points (assumed to be 1 unit apart). The trapz(Y) function calculates the approximate integral of these values.
The output Q is the approximate integral of the values in Y, giving you the total area under the curve represented by these points.
On execution the output is −
>> Y = [1, 4, 9, 16, 25];
Q = trapz(Y);
disp(Q);
42
>>
Example 2: Using trapz(Y) where Y is a matrix
Y = [1 2 3; 4 5 6; 7 8 9]; Q = trapz(Y); disp(Q);
The matrix represents multiple sets of data, one in each column. The trapz(Y) function calculates the integral for each column separately.The output Q is a row vector containing the integral values for each column of Y. For this matrix, trapz(Y) gives a result like [8 10 12].
On execution the output is −
>> Y = [1 2 3; 4 5 6; 7 8 9];
Q = trapz(Y);
disp(Q);
8 10 12
>>
Example 3: Using trapz(X,Y) where X and Y are vectors
X = [1, 2, 3, 4, 5]; Y = [1, 4, 9, 16, 25]; Q = trapz(X, Y); disp(Q);
Here, X is a vector, and its length must match the size of Y's first dimension (both have 5 elements).
trapz(X,Y) calculates the integral of Y with respect to X, considering the actual spacing between points in X.
On code execution the output we get is as follows −
>> X = [1, 2, 3, 4, 5];
Y = [1, 4, 9, 16, 25];
Q = trapz(X, Y);
disp(Q);
42
>>
Example 4: Using trapz(X,Y) where X is a scalar number
X = 1; Y = [1, 4, 9, 16, 25]; Q = trapz(X, Y); disp(Q);
Here, X is a single number representing uniform spacing between the points.
trapz(X,Y) is equivalent to multiplying X by the result of trapz(Y). Since X=1, it simply returns the result of trapz(Y) assuming evenly spaced points.
The output Q will give the approximate integral of Y assuming a uniform spacing of 1 unit between the points.
When the code is executed the output we get is as follows.
>> X = 1;
Y = [1, 4, 9, 16, 25];
Q = trapz(X, Y);
disp(Q);
42
>>
Example 5: Integrating along columns with dim=1
The code we have is −
Y = [1 2 3;
4 5 6;
7 8 9];
Q = trapz(Y, 1);
disp(Q);
Dimension: dim = 1 means integrating along each column of the matrix.The result will be a row vector where each element represents the integral of the corresponding column of Y.Here, each value in Q is the integral of the values in each column of Y.
On execution the output we get is.
>> Y = [1 2 3;
4 5 6;
7 8 9];
Q = trapz(Y, 1);
disp(Q);
8 10 12
>>
Example 6: Integrating Along Rows (dim = 2)
The code we have is.
Y = [1 2 3;
4 5 6;
7 8 9];
Q = trapz(Y, 2);
disp(Q);
Dimension − dim = 2 means integrating along each row of the matrix.
The result will be a column vector where each element represents the integral of the corresponding row of Y.
When the code is executed the output is :
>> Y = [1 2 3;
4 5 6;
7 8 9];
Q = trapz(Y, 2);
disp(Q);
4
10
16
>>