L_TRAPZ

=L_TRAPZ(y, [x], [dx], [cumulative])

Description

L_TRAPZ is a numerical integration function that uses the Trapezoidal Rule [1] to approximate a definite integral from a set of ordered (x,y) values. If the x values are not given, you can either specify a value for dx, or leave it blank to assume unit spacing, Δx=1.

For our examples we'll use the polynomial function \(p(x)=-0.5x^3+5x^2-10x+10\). The image below shows a graph of 4 points between x=0 and x=6, with uniform spacing Δx=2. We'll use the trapezoidal rule to approximate the definite integral \(\int_0^6p(x)dx\). The approximation of the integral is the area shown in purple.

The trapezoidal rule uses a linear approximation between each point which sometimes overestimates and sometimes underestimates the true area under the curve. The area (integral) of each separate interval \([x_i,x_{i+1}]\) is:

We can use L_POLYINT to calculate the exact integral for this polynomial example:

=LET(
    px, {-0.5, 5, -10, 10},
    pxint, L_POLYINT(px),
    L_DIFF(L_POLYVAL(pxint,{0;6}))
)
Result: 78

Now we'll use L_TRAPZ to approximate the integral using the 4 points shown in the image above.

=LET(
    x, {0;2;4;6},
    y, {10;6;18;22},
    L_TRAPZ(y,x)
)
Result: 80

L_TRAPZ(y,,2) would produce the same result because the x values have uniform spacing of Δx=2.

To get a more accurate approximation, we can use a smaller Δx interval. In the following example we're using 1000 intervals (1001 points) of size Δx=(6-0)/1000=0.006.

=LET(
    x, L_LINSPACE(0,6,1001),
    y, L_POLYVAL({-0.5,5,-10,10},x),
    L_TRAPZ(y,x)
)
Result: 78.000018

When Y is a Matrix

This function is similar to the Matlab and NumPy trapz functions but it behaves a bit differently. When Y is a matrix, each of the column vectors in Y is integrated separately and the function returns a row vector (not a column vector) where each element is the integration of a column of Y.

To integrate a matrix Y by row instead of by column, use TRANSPOSE(Y).

If x is a vector, each column of Y uses the same x vector for the Δx spacing. If X is a matrix, then each column of Y uses the corresponding column of X. Each column is an independent integration.

The example below uses L_MESHGRID to create a 3D surface and L_TRAPZ to estimate the volume under the surface.

Cumulative Trapezoidal Integration

The updated function includes a cumulative option which if TRUE will return the intermediate values of the integral (a cumulative sum of the area under the curve). Remember that Excel starts array indexing at 1, so we use \(F(x_k)\) to represent the cumulative integral, with \(F(x_1)=0\). For k = 2 to N (where N is the number of values in vector y), we have:

When the cumulative option is TRUE, the function returns the intermediate values as column vectors, so the size of the array that is returned is the same as the size of the input Y. This also works when Y and X are matrices - you just have to make sure to remember that L_TRAPZ integrates each column, so use TRANSPOSE as needed.

The following image is a screenshot from the Lambda Library template file showing the results of the second example below. Note that one of the benefits of trapezoidal integration is that the cumulative result is the same size as the x and y vectors, making it easy to set up for graphing.

Lambda Formula

This code for using L_TRAPZ in Excel is provided under the License as part of the LAMBDA Library, but to use just this function, you may copy the following code directly into your spreadsheet.

/**
* Numerical Integration using the Trapezoidal Rule
* Defaults: dx=1 if x and dx are blank, cumulative=FALSE
*/
L_TRAPZ = LAMBDA(y,[x],[dx],[cumulative],
LET(doc,"https://www.vertex42.com/lambda/trapz.html",
    dx,IF(ISBLANK(dx),1,dx),
    areas,IF(ISOMITTED(x),
        (DROP(y,-1)+DROP(y,1))/2*dx,
        ((DROP(y,-1)+DROP(y,1))/2)*(DROP(x,1)-DROP(x,-1))
    ),
    REDUCE("",SEQUENCE(COLUMNS(areas)),
        LAMBDA(acc,i,LET(
            col,CHOOSECOLS(areas,i),
            new,IF(cumulative=TRUE,
                SCAN(0,SEQUENCE(ROWS(col)+1),
                    LAMBDA(acc,i,IF(i=1,0,acc+INDEX(col,i-1,1)))
                ),
                SUM(col)
            ),
            IF(i=1,new,HSTACK(acc,new))
        ))
    )
));

Name: L_TRAPZ
Description: Numeric Integration using the Trapezoidal Rule
Arguments: y, x, dx 
Function:
=LET(doc,"https://www.vertex42.com/lambda/trapz.html",
    dx,IF(ISBLANK(dx),1,dx),
    areas,IF(ROWS(x)>1,
        ARRAYFORMULA(((L_DROP(y,-1,0)+L_DROP(y,1,0))/2)*(L_DROP(x,1,0)-L_DROP(x,-1,0))),
        ARRAYFORMULA((L_DROP(y,-1,0)+L_DROP(y,1,0))/2*dx)        
    ),
    REDUCE("",SEQUENCE(COLUMNS(areas)),
        LAMBDA(acc,i,LET(
            col,CHOOSECOLS(areas,i),
            new,IF(cumulative=TRUE,
                SCAN(0,SEQUENCE(ROWS(col)+1),
                    LAMBDA(acc,i,IF(i=1,0,acc+INDEX(col,i-1,1)))
                ),
                SUM(col)
            ),
            IF(i=1,new,HSTACK(acc,new))
        ))
    )
)

L_TRAPZ Examples

Test: Copy and Paste this LET function into a cell
=LET(
    xvec, L_LINSPACE(-2,2,101),
    yvec, L_LINSPACE(-3,3,201),
    X, L_MESHGRID(xvec),yvec,"X"),
    Y, L_MESHGRID(xvec),yvec,"Y"),
    Z, X^2-X*Y+Y^2,
    L_TRAPZ(TRANSPOSE(L_TRAPZ(TRANSPOSE(Z),xvec)),yvec)
)

Result: 104.01

=LET(
    ...
    L_TRAPZ(TRANSPOSE(L_TRAPZ(Z,yvec)),xvec)
)

Result: 104.01

=LET(
    x, L_LINSPACE(0,10,11),
    y, x^2,
    L_TRAPZ(y,x,,TRUE)
)

Result: {0; 0.5; 3; 9.5; 22; 42.5; 73; 115.5; 172; 244.5; 335}