19.4 Numerical Integration in Detail

The really wonderful thing about this is how easy it is to implement it.

Suppose we want to integrate sin2(x) from 0 to 3. This is a random thing to want to do, and you can perform any similar task in the same way; and once you have done it you can change limits of functions almost immediately.

Here is what I would do, you of course may choose differently, and better.

First I would leave some space at the top for input data and use labels so that I can recognize what the spreadsheet is about at some future date.

Thus in A1 I would write Integration, and in A2 write a, A3 b, A4 n, A5 d, where a, b, n, and d refer respectively to the lower and upper limits of integration, n to the number of evaluation points -1, and d to their separation, which is .

Then in B2 through B5, I would insert the values of these that we want to choose, namely, 0, 3, 1000 (why not), and =(b3-b2)/b4.

So far we have only labeled our input data. Now we want to construct columns for, the values of the evaluation points, the function at them and the sum of same.

If we only want to do the trapezoid rule, that is enough. If we want Simpson's rule as well, we can add a column that gives the value of the function at every second point only. This will allow us to do the trapezoid rule with 2d separation in parallel with that for d separation, and combine the two to get Simpson's.

So in column a we label, say in a9, x; b9, f(x); c9, left rule, d9 right rule, e9 left by 2, f9 right by 2.

And now let us construct these columns. In a10 we put =b2. In a11 we put =a10+b$5 and copy a11 this down to a1010.
Then in b10 we put =(sin(a10))^2 (which you can later change to whatever you want to integrate) and copy this down to b1010.
Now in c11 we put =c10+b10*$b$5 and copy down to c1010. This produces the left hand rule computation.
In d11 put =d10+b11*$b$5 and copy down this will give the right hand rule computation.
In e12 put =e10+b10 *$b$5*2 and in f12 put = f10 +b12*$b$5*2, and copy these down to row 1010.

What you have done is to produce two left hand and two right hand rule computations.

And what is the answer? If you want the integral to the argument in a1010 (or mutatis mutandis, anywhere else), by the trapezoid rule for separation d it is a=(1010b+1010c)/2, and for the trapezoid rule for separation 2d it is =(1010d+1010e)/2.

So let us set b6=(1010b+1010c)/2, and b7 =(1010d+1010e)/2
Then we can put the Simpson rule answer in b8, and it will be =(4*b6-b7)/3, and we are done.

Notice that what we did involved three kinds of acts: first labeling. Second entering the successive arguments at which we evaluate f and its values there. Third, setting up a bunch of sums, to produce left and right rules for d and 2d separations of evaluation arguments.

The nice thing to notice is that the first and last steps need be done only once ever. Changing limits only involves changing b2 and b3. Changing functions only involves changing the entry in b10 and copying it down to b1010.

What do you do if you want to use infinite limits?

You can't get there from here. You can however, change variables so that the range of integration in the new variables becomes finite. (Or you can break your interval into finite pieces, and write your integral as a sum of integrals over these. For most integrals you will encounter, the contributions from all but the first few pieces will be negligible and you can get by summing over what is essentially a finite domain.)


19.1 Set up a spreadsheet as we have described it. Check your answer for accuracy.

19.2 Add two more columns (computing results for separation 4d) and use these to compute the super Simpson answer.

We computed the left and right hand rules for d and 2d in separate columns and then calculated the trapezoid and Simpson answers in B6, B7 and B8.

Instead we can devote a column to each of these, setting f10=(b10+c10)/2, g10=(d10+e10)/2 and h10=(4*f10-g10)/3, and copying these down and then these will represent the trapezoid and Simpson's rule integrals to all endpoints that are of the form a + 4kd. (You can ignore the other entries in the Simpson's rule column, and ignore every second entry in the trapezoid column)

If you do this you can use a chart to plot all columns and see what the various rules lead to in the integral as you go along. You can get an idea from this of the relative accuracy of solutions, and see visually how that changes as you change d.


19.3 Do this for integrands (sin(x)^2, and sin(x), and try a few others, for example try exp(-x^2) from 0 to say 10.

19.4 Do the same for the integrands (exp(-x)-1))/x, from 0 on.to say 10, and also sin(x)/x. For each of these you have to put in values at 0 (-1 and 1 respectively)separately by hand or the machine will accuse you of dividing by 0.) Extimate the accuracy of your answer.

By the way, if you want to check on an integral you have evaluated algebraically, you can add still one more column, putting in i10=your answer(a10) and copy it down, and have your chart show this column as well. If it is wrong, that will stand out like a sore thumb.

What happens when the integrand becomes infinite in the region of integration?

This, and using an infinite limit point creates what is called an "improper integral". Such a thing runs the risk of being infinite, or of being not well defined.

However, it is often possible to make perfect sense out of such "improper integrals".

For example, suppose we want to integrate x-1 from -1 to 1. This integrand behaves badly near x = 0, and the area between 0 and 1 can be shown to be unbounded.
However, if you think about it, there is a negative contribution between -1 and 0 that exactly matches the positive contributions, and it is reasonable to claim that this integral should be 0.

If an integrand is singular at x = a, we can define the principal part of an integral whose lower limit is below a and upper limit above a, by omitting a tiny interval (from a - d to a + d) and letting the size of the interval converge to 0. If the integral converges to something, that is its principal part. Thus the principal part of the integral of x - 1 from -1 to 1 is indeed 0.

To see that you have to be careful, consider the integral of integrand between -1 and 1. This function and the area under it are both unbounded, and both positive for all x in this interval. It is tempting to notice that is the derivative of , so that by the fundamental theorem of calculus, the integral should be or -2.

That an infinite positive area is -2 is a bit hard to swallow. One way to see what went wrong is to notice that the integral from 0 to 1 is exactly the same as the integral from -1 to 1 and each is . The answer -2 comes from ignoring the contributions to this integral from the lower limit at x = 0.

Strangely enough there is a sense in which this integral actually is -2, which you can discover by learning about integration in the complex plane. If you are allowed to deform your interval of integration to go into the complex plane and avoid 0, the answer -2 is always correct if you do so!