Example of a Package in a Notebook: Trapezoid Integration of Data Points

Integrate by the Trapezoid Rule

This example of writing your own numerical integration function is from Bahder's excellent book, Mathematica for Scientists and Engineers. I use it to show how to write a Package in a Notebook, as opposed to in a .m file.

While it always pays to see if Mathematica has a built-in function to do what you want, I don't see that NIntegrate can integrate over data points, although it has this nice introduction for the TrapezoidRule: "The trapezoidal rule for integral estimation is one of the simplest and oldest rules (possibly used by the Babylonians and certainly by the ancient Greek mathematicians)."

I needed a function to integrate over data points from an irregular, piecewise, sawtooth-type function. Bahder notes that integrating via the trapezoid rule approximates the integral over a list of (x,y) data points that may be irregularly spaced. Here is the formula. Keyboard shortcuts are faster than a palette; here they are Ctrl+- for subscript and Ctrl+spacebar to exit a subscript. Remember to use a double equal sign for 'equals'.

Here is the numerical integral in a Package. The Package is in one Cell and I set the Cell to Initialization (right-click on the Cell and select Initialization Cell) to automatically load the Package when I evaluate any function in the Notebook.

Note Bahder's use of a named Pattern and Repeated to type-check the data argument. It compares a List of repeated two-element Lists to the pairs that we give it. He then makes sure that there is more than one data pair using Condition. Some say the Wolfram Language is weakly-typed. That is incorrect; it has powerful type-creation and type-checking capabilities.

BeginPackage@"BahderMSE`";

IntegrateTrapezoid::usage="IntegrateTrapezoid[data:{{_,_}..}/;Length@data>1] takes a list of (x,y) pairs, x1 < x2 <x3 ..., that are not necessarily evenly spaced and returns a numerical approximation to their integral using the trapezoid rule.";

Begin@"`Private`";

IntegrateTrapezoid[data:{{_,_}..}/;Length@data>1]:=1/2 Sum[(data[[n+1,1]]-data[[n,1]])(data[[n,2]]+data[[n+1,2]]),{n,Length[data]-1}]

End[];
EndPackage[]

Let's check the usage statement.

In[7]:= ?IntegrateTrapezoid

IntegrateTrapezoid[data:{{_,_}..}/;Length@data>1] takes a list of (x,y) pairs, x1 < x2 <x3 ..., that are not necessarily evenly spaced and returns a numerical approximation to their integral using the trapezoid rule.

Now test it on the same data Bahder used...

In[8]:= data1={{1.2,51.3},{2.3,72.8},{3.0,87.1},{3.7,53.2}};

And we do get the same result.

In[9]:= IntegrateTrapezoid@data1

Out[9]= 173.325

Let's try a different example using Mathematica's SawtoothWave function. This is a typical use of Thread to combine the x and y points in their own Lists.

In[10]:= Thread@{Range[0,3,0.1],SawtoothWave@#&/@Range[0,3,0.1]}

Out[10]= {{0.,0.},{0.1,0.1},{0.2,0.2},{0.3,0.3},{0.4,0.4},{0.5,0.5},{0.6,0.6},{0.7,0.7},{0.8,0.8},{0.9,0.9},{1.,0.},{1.1,0.1},{1.2,0.2},{1.3,0.3},{1.4,0.4},{1.5,0.5},{1.6,0.6},{1.7,0.7},{1.8,0.8},{1.9,0.9},{2.,0.},{2.1,0.1},{2.2,0.2},{2.3,0.3},{2.4,0.4},{2.5,0.5},{2.6,0.6},{2.7,0.7},{2.8,0.8},{2.9,0.9},{3.,0.}}

We see that it is a numerical approximation of the full SawtoothWave where the peaks fall short of 1.

There are three triangles of unit length and height = 0.9, so the total area by a = 1/2 b * h is 1.35. The function works for this example but would need to be validated and its accuracy evaluated for others.

IntegrateTrapezoid@%%

Out[12]= 1.35

Thomas Bahder, Mathematica for Scientists and Engineers