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

>> 
Advertisements