Section 1.11 Numerical Integration
By now the reader will have come to appreciate that integration is generally quite a bit more difficult than differentiation. There are a great many simple-looking integrals, such as that are either very difficult or even impossible to express in terms of standard functions . Such integrals are not merely mathematical curiosities, but arise very naturally in many contexts. For example, the error function
β1β
We apologise for being a little sloppy here β but we just want to say that it can be very hard or even impossible to write some integrals as some finite sized expression involving polynomials, exponentials, logarithms and trigonometric functions. We donβt want to get into a discussion of computability, though that is a very interesting topic.
is extremely important in many areas of mathematics, and also in many practical applications of statistics.
In such applications we need to be able to evaluate this integral (and many others) at a given numerical value of In this section we turn to the problem of how to find (approximate) numerical values for integrals, without having to evaluate them algebraically. To develop these methods we return to Riemann sums and our geometric interpretation of the definite integral as the signed area.
We start by describing (and applying) three simple algorithms for generating, numerically, approximate values for the definite integral In each algorithm, we begin in much the same way as we approached Riemann sums.
- We first select an integer
called the βnumber of stepsβ. - We then divide the interval of integration,
into equal subintervals, each of length The first subinterval runs from to The second runs from to and so on. The last runs from to
This splits the original integral into pieces:
Each subintegral is approximated by the area of a simple geometric figure. The three algorithms we consider approximate the area by rectangles, trapezoids and parabolas (respectively).
We will explain these rules in detail below, but we give a brief overview here:
- The midpoint rule approximates each subintegral by the area of a rectangle of height given by the value of the function at the midpoint of the subinterval
- The trapezoidal rule approximates each subintegral by the area of a trapezoid with vertices at
- Simpsonβs rule approximates two adjacent subintegrals by the area under a parabola that passes through the points
and
Subsection 1.11.1 The midpoint rule
The integral represents the area between the curve and the -axis with running from to The width of this region is The height varies over the different values that takes as runs from to
The midpoint rule approximates this area by the area of a rectangle of width and height which is the exact height at the midpoint of the range covered by
The area of the approximating rectangle is and the midpoint rule approximates each subintegral by
Applying this approximation to each subinterval and summing gives us the following approximation of the full integral:
So notice that the approximation is the sum of the function evaluated at the midpoint of each interval and then multiplied by Our other approximations will have similar forms.
In summary:
Equation 1.11.2. The midpoint rule.
Example 1.11.3. .
We approximate the above integral using the midpoint rule with step.
Solution:
- First we set up all the
-values that we will need. Note that and - In this case we can compute the integral exactly (which is one of the reasons it was chosen as a first example):
- So the error in the approximation generated by eight steps of the midpoint rule is
- The relative error is then
times the actual value of the integral. - We can write this as a percentage error by multiplying it by 100
of the exact value.
The midpoint rule gives us quite good estimates of the integral without too much work β though it is perhaps a little tedious to do by hand . Of course, it would be very helpful to quantify what we mean by βgoodβ in this context and that requires us to discuss errors.
β2β
Thankfully it is very easy to write a program to apply the midpoint rule.
Definition 1.11.4.
Example 1.11.5. .
As a second example, we apply the midpoint rule with steps to the above integral.
- We again start by setting up all the
-values that we will need. So and - Again, we have chosen this example so that we can compare it against the exact value:
- So with eight steps of the midpoint rule we achieved
of its true value.
Subsection 1.11.2 The trapezoidal rule
Consider again the area represented by the integral The trapezoidal rule (unsurprisingly) approximates this area by a trapezoid whose vertices lie at
β3β
This method is also called the βtrapezoid ruleβ and βtrapezium ruleβ.
β4β
A trapezoid is a four sided polygon, like a rectangle. But, unlike a rectangle, the top and bottom of a trapezoid need not be parallel.
The trapezoidal approximation of the integral is the shaded region in the figure on the right above. It has width Its left hand side has height and its right hand side has height
As the figure below shows, the area of a trapezoid is its width times its average height.
So the trapezoidal rule approximates each subintegral by
Applying this approximation to each subinterval and then summing the result gives us the following approximation of the full integral
So notice that the approximation has a very similar form to the midpoint rule, excepting that
- we evaluate the function at the
βs rather than at the midpoints, and - we multiply the value of the function at the endpoints
by
In summary:
Equation 1.11.6. The trapezoidal rule.
To compare and contrast we apply the trapezoidal rule to the examples we did above with the midpoint rule.
Example 1.11.7. β using the trapezoidal rule.
Solution: We proceed very similarly to Example 1.11.3 and again use steps.
- We again have
and - Applying the trapezoidal rule, Equation 1.11.6, gives
- The exact value of the integral is still
So the error in the approximation generated by eight steps of the trapezoidal rule is which is of the exact answer. Notice that this is roughly twice the error that we achieved using the midpoint rule in Example 1.11.3.
Example 1.11.8. β using the trapezoidal rule.
Solution: We proceed very similarly to Example 1.11.5 and again use steps.
Subsection 1.11.3 Simpsonβs Rule
When we use the trapezoidal rule we approximate the area by the area between the -axis and a straight line that runs from to β that is, we approximate the function on this interval by a linear function that agrees with the function at each endpoint. An obvious way to extend this β just as we did when extending linear approximations to quadratic approximations in our differential calculus course β is to approximate the function with a quadratic. This is precisely what Simpsonβs rule does.
β5β
Simpsonβs rule is named after the 18th century English mathematician Thomas Simpson, despite its use a century earlier by the German mathematician and astronomer Johannes Kepler. In many German texts the rule is often called Keplerβs rule.
Simpsonβs rule approximates the integral over two neighbouring subintervals by the area between a parabola and the -axis. In order to describe this parabola we need 3 distinct points (which is why we approximate two subintegrals at a time). That is, we approximate
by the area bounded by the parabola that passes through the three points and the -axis and the vertical lines and
We repeat this on the next pair of subintervals and approximate by the area between the -axis and the part of a parabola with This parabola passes through the three points and And so on. Because Simpsonβs rule does the approximation two slices at a time, must be even.
To derive Simpsonβs rule formula, we first find the equation of the parabola that passes through the three points and Then we find the area between the -axis and the part of that parabola with To simplify this computation consider a parabola passing through the points and
Write the equation of the parabola as
Adding the first and third equations together gives us
To this we add four times the middle equation
This means that
Note that here
is one half of the length of the -interval under consideration is the height of the parabola at the left hand end of the interval under consideration is the height of the parabola at the middle point of the interval under consideration is the height of the parabola at the right hand end of the interval under consideration
In summary
Equation 1.11.9. Simpsonβs rule.
Notice that Simpsonβs rule requires essentially no more work than the trapezoidal rule. In both rules we must evaluate at but we add those terms multiplied by different constants .
β6β
There is an easy generalisation of Simpsonβs rule that uses cubics instead of parabolas. It is known as Simpsonβs second rule and Simpsonβs rule. While one can push this approach further (using quartics, quintics etc), it can sometimes lead to larger errors β the interested reader should look up Rungeβs phenomenon.
Letβs put it to work on our two running examples.
Example 1.11.10. β using Simpsonβs rule.
Solution: We proceed almost identically to Example 1.11.7 and again use steps.
It is striking that the absolute error approximating with Simpsonβs rule is so much smaller than the error from the midpoint and trapezoidal rules.
Example 1.11.11. β Simpsonβs rule.
This completes our derivation of the midpoint, trapezoidal and Simpsonβs rules for approximating the values of definite integrals. So far we have not attempted to see how efficient and how accurate the algorithms are in general. Thatβs our next task.
Subsection 1.11.4 Three Simple Numerical Integrators β Error Behaviour
Now we are armed with our three (relatively simple) methods for numerical integration we should give thought to how practical they might be in the real world . Two obvious considerations when deciding whether or not a given algorithm is of any practical value are
β7β
Indeed, even beyond the βreal worldβ of many applications in first year calculus texts, some of the methods we have described are used by actual people (such as ship builders, engineers and surveyors) to estimate areas and volumes of actual objects!
- the amount of computational effort required to execute the algorithm and
- the accuracy that this computational effort yields.
For algorithms like our simple integrators, the bulk of the computational effort usually goes into evaluating the function The number of evaluations of required for steps of the midpoint rule is while the number required for steps of the trapezoidal and Simpsonβs rules is So all three of our rules require essentially the same amount of effort β one evaluation of per step.
To get a first impression of the error behaviour of these methods, we apply them to a problem whose answer we know exactly:
To be a little more precise, we would like to understand how the errors of the three methods change as we increase the effort we put in (as measured by the number of steps ). The following table lists the error in the approximate value for this number generated by our three rules applied with three different choices of It also lists the number of evaluations of required to compute the approximation.
Midpoint | Trapezoidal | Simpsonβs | ||||
n | error | # evals | error | # evals | error | # evals |
10 | 10 | 11 | 11 | |||
100 | 100 | 101 | 101 | |||
1000 | 1000 | 1001 | 1001 |
Observe that
- Using 101 evaluations of
worth of Simpsonβs rule gives an error 75 times smaller than 1000 evaluations of worth of the midpoint rule. - The trapezoidal rule error with
steps is about twice the midpoint rule error with steps. - With the midpoint rule, increasing the number of steps by a factor of 10 appears to reduce the error by about a factor of
- With the trapezoidal rule, increasing the number of steps by a factor of 10 appears to reduce the error by about a factor of
- With Simpsonβs rule, increasing the number of steps by a factor of 10 appears to reduce the error by about a factor of
To test these conjectures for the behaviour of the errors we apply our three rules with about ten different choices of of the form with integer. Figure 1.11.12 contains two graphs of the results. The left-hand plot shows the results for the midpoint and trapezoidal rules and the right-hand plot shows the results for Simpsonβs rule.
For each rule we are expecting (based on our conjectures above) that the error
with steps is (roughly) of the form
for some constants and We would like to test if this is really the case, by graphing against and seeing if the graph βlooks rightβ. But it is not easy to tell whether or not a given curve really is for some specific by just looking at it. However, your eye is pretty good at determining whether or not a graph is a straight line. Fortunately, there is a little trick that turns the curve into a straight line β no matter what is.
Instead of plotting against we plot against This transformation works because when
β8β
There is a variant of this trick that works even when you donβt know the answer to the integral ahead of time. Suppose that you suspect that the approximation satisfies where is the exact value of the integral and suppose that you donβt know the values of and Then so plotting against gives the straight line
The three graphs in Figure 1.11.12 plot against for our three rules. Note that we have chosen to use logarithms with this βunusual baseβ because it makes it very clear how much the error is improved if we double the number of steps used. To be more precise β one unit step along the -axis represents changing For example, applying Simpsonβs rule with steps results in an error of so the point has been included on the graph. Doubling the effort used β that is, doubling the number of steps to β results in an error of So, the data point lies on the graph. Note that the -coordinates of these points differ by 1 unit.
β9β
Now is a good time for a quick revision of logarithms β see βWhirlwind review of logarithmsβ in Section 2.7 of the CLP-1 text.
For each of the three sets of data points, a straight line has also been plotted βthroughβ the data points. A procedure called linear regression has been used to decide precisely which straight line to plot. It provides a formula for the slope and -intercept of the straight line which βbest fitsβ any given set of data points. From the three lines, it sure looks like for the midpoint and trapezoidal rules and for Simpsonβs rule. It also looks like the ratio between the value of for the trapezoidal rule, namely and the value of for the midpoint rule, namely is pretty close to 2:
β10β
Linear regression is not part of this course as its derivation requires some multivariable calculus. It is a very standard technique in statistics.
The intuition, about the error behaviour, that we have just developed is in fact correct β provided the integrand is reasonably smooth. To be more precise
Theorem 1.11.13. Numerical integration errors.
The first of these error bounds in proven in the following (optional) section. Here are some examples which illustrate how they are used. First let us check that the above result is consistent with our data in Figure 1.11.12
Example 1.11.14. Midpoint rule error approximating .
- The integral
has - The second derivative of the integrand satisfies
- So the error,
introduced when steps are used is bounded by
In a typical application we would be asked to evaluate a given integral to some specified accuracy. For example, if you are manufacturer and your machinery can only cut materials to an accuracy of of a millimeter, there is no point in making design specifications more accurate than of a millimeter.
Example 1.11.15. How many steps for a given accuracy?
Suppose, for example, that we wish to use the midpoint rule to evaluate
β11β
This is our favourite running example of an integral that cannot be evaluated algebraically β we need to use numerical methods.
to within an accuracy of
Solution:
- The integral has
and - The first two derivatives of the integrand are
- As
runs from 0 to 1, increases from to so that - The error introduced by the
step midpoint rule is at most - We need this error to be smaller than
so steps of the midpoint rule will do the job. - In fact
results in an error of about
That seems like far too much work, and the trapezoidal rule will have twice the error. So we should look at Simpsonβs rule.
Example 1.11.16. How many steps using Simpsonβs rule?
Suppose now that we wish evaluate to within an accuracy of β but now using Simpsonβs rule. How many steps should we use?
Solution:
- Again we have
- We then need to bound
on the domain of integration, - Now, for any
Also, for is bounded by and so - The error introduced by the
step Simpsonβs rule is at most - In order for this error to be no more than
we require to satisfy steps of Simpsonβs rule will do the job. steps actually results in an error of The reason that we get an error so much smaller than we need is that we have overestimated the number of steps required. This, in turn, occurred because we made quite a rough bound of If we are more careful then we will get a slightly smaller It actually turns outthat you only needβ12β
The authors tested this empirically. to approximate within
Subsection 1.11.5 Optional β An error bound for the midpoint rule
We now try develop some understanding as to why we got the above experimental results. We start with the error generated by a single step of the midpoint rule. That is, the error introduced by the approximation
To do this we are going to need to apply integration by parts in a sneaky way. Let us start by considering a subinterval and letβs call the width of the subinterval so that If we were to now apply the midpoint rule to this subinterval, then we would write
β13β
We chose this interval so that we didnβt have lots of subscripts floating around in the algebra.
The sneaky trick we will employ is to write
and then examine each of the integrals on the right-hand side (using integration by parts) and show that they are each of the form
Let us apply integration by parts to β with so and we will make the slightly non-standard choice of
Notice that the first term on the right-hand side is the term we need, and that our non-standard choice of allowed us to avoid introducing an term.
To obtain a similar expression for the other integral, we repeat the above steps and obtain:
Now add together these two expressions
Then since
Rearrange this expression a little and take absolute values
where we have also made use of the triangle inequality . By assumption on the interval so
β14β
The triangle inequality says that for any real numbers
where we have used in the last step.
Putting everything together we see that the error using the midpoint rule is bounded by
as required.
Exercises 1.11.6 Exercises
Recall that we are using to denote the logarithm of with base In other courses it is often denoted
Exercises β Stage 1 .
1.
Suppose we approximate an object to have volume when its exact volume is Give the relative error, absolute error, and percent error of our approximation.
2.
3.
4.
5.
6.
7.
Suppose my mother is under 100 years old, and I am under 200 years old. Who is older?
β15β
Weβre going somewhere with this.
8.
- True or False: for fixed positive constants
and with - True or False: for a function
and fixed constants and with the -interval midpoint approximation of is more accurate than the -interval trapezoidal approximation.
9. (β³).
Decide whether the following statement is true or false. If false, provide a counterexample. If true, provide a brief justification.
Whenis positive and concave up, any trapezoidal rule approximation for will be an upper estimate for
10.
Give a polynomial with the property that the Simpsonβs rule approximation of is exact for all and
Exercises β Stage 2 .
11.
Write out all three approximations of with (That is: midpoint, trapezoidal, and Simpsonβs.) You do not need to simplify your answers.
12. (β³).
13. (β³).
The solid is 40 cm high and the horizontal cross sections are circular disks. The table below gives the diameters of the cross sections in centimeters at 10 cm intervals. Use the trapezoidal rule to estimate the volume of
height | 0 | 10 | 20 | 30 | 40 |
diameter | 24 | 16 | 10 | 6 | 4 |
14. (β³).
A metre long cedar log has cross sections that are approximately circular. The diameters of the log, measured at one metre intervals, are given below:
metres from left end of log | 0 | 1 | 2 | 3 | 4 | 5 | 6 |
diameter in metres | 1.2 | 1 | 0.8 | 0.8 | 1 | 1 | 1.2 |
Use Simpsonβs Rule to estimate the volume of the log.
15. (β³).
The circumference of an 8 metre high tree at different heights above the ground is given in the table below. Assume that all horizontal cross-sections of the tree are circular disks.
height (metres) | 0 | 2 | 4 | 6 | 8 |
circumference (metres) | 1.2 | 1.1 | 1.3 | 0.9 | 0.2 |
Use Simpsonβs rule to approximate the volume of the tree.
16. (β³).
By measuring the areas enclosed by contours on a topographic map, a geologist determines the cross sectional areas in of a m high hill. The table below gives the cross sectional area at various heights The volume of the hill is
0 | 10 | 20 | 30 | 40 | 50 | 60 | |
10,200 | 9,200 | 8,000 | 7,100 | 4,500 | 2,400 | 100 |
- If the geologist uses the Trapezoidal Rule to estimate the volume of the hill, what will be their estimate, to the nearest 1,000
- What will be the geologistβs estimate of the volume of the hill if they use Simpsonβs Rule instead of the Trapezoidal Rule?
17. (β³).
The graph below applies to both parts (a) and (b).
- Use the Trapezoidal Rule, with
to estimate the area under the graph between and Simplify your answer completely. - Use Simpsonβs Rule, with
to estimate the area under the graph between and
18. (β³).
The integral is estimated using the Midpoint Rule with intervals. Show that the absolute error in this approximation is at most
You may use the fact that when approximating with the Midpoint Rule using points, the absolute value of the error is at most when for all
19. (β³).
The total error using the midpoint rule with subintervals to approximate the integral of over is bounded by if for all
Using this bound, if the integral is approximated using the midpoint rule with subintervals, what is the largest possible error between the approximation and the true value of the integral?
20. (β³).
Both parts of this question concern the integral
- Write down the Simpsonβs Rule approximation to
with Leave your answer in calculator-ready form. - Which method of approximating
results in a smaller error bound: the Midpoint Rule with intervals, or Simpsonβs Rule with intervals? You may use the formulas is an upper bound for and is an upper bound for and and are the absolute errors arising from the midpoint rule and Simpsonβs rule, respectively.
21. (β³).
Find a bound for the error in approximating using Simpsonβs rule with Do not write down the Simpsonβs rule approximation
In general the error in approximating using Simpsonβs rule with steps is bounded by where and for all
22. (β³).
Find a bound for the error in approximating
In general, the error in approximating using Simpsonβs rule with steps is bounded by where and for all
23. (β³).
Let
- Write down the trapezoidal approximation
for You do not need to simplify your answer. - Write down the Simpsonβs approximation
for You do not need to simplify your answer. - Without computing
find an upper bound for You may use the fact that if on the interval then the error in using to approximate has absolute value less than or equal to
24. (β³).
- Find the best Trapezoidal Rule and Simpsonβs Rule approximations that you can for
- Determine the maximum possible sizes of errors in the approximations you gave in part (a). Recall that if a function
satisfies on then
25. (β³).
Consider the trapezoidal rule for making numerical approximations to The error for the trapezoidal rule satisfies , where for If for find a value of to guarantee the trapezoidal rule will give an approximation for with absolute error, less than
Exercises β Stage 3 .
26. (β³).
A swimming pool has the shape shown in the figure below. The vertical cross-sections of the pool are semi-circular disks. The distances in feet across the pool are given in the figure at 2--foot intervals along the sixteen--foot length of the pool. Use Simpsonβs Rule to estimate the volume of the pool.
27. (β³).
A piece of wire 1m long with radius 1mm is made in such a way that the density varies in its cross-section, but is radially symmetric (that is, the local density in depends only on the distance in mm from the centre of the wire). Take as given that the total mass of the wire in kg is given by
Data from the manufacturer is given below:
0 | 1/4 | 1/2 | 3/4 | 1 | |
8051 | 8100 | 8144 | 8170 | 8190 |
- Find the best Trapezoidal Rule approximation that you can for
based on the data in the table. - Suppose that it is known that
and for all values of Determine the maximum possible size of the error in the approximation you gave in part (a). Recall that if a function satisfies on then and is the Trapezoidal Rule approximation to using subintervals.
28. (β³).
- Use Simpsonβs rule with 6 subintervals to approximate
-
How many subintervals are required in order to guarantee that the absolute error is less thanNote that if
is the error using subintervals, then where is the maximum absolute value of the fourth derivative of the function being integrated and and are the end points of the interval.
29. (β³).
- Estimate the maximum absolute error in using
to approximate - How large should
be in order to ensure that
Note: The graph of where is shown below. The absolute error in the Simpsonβs rule approximation is bounded by when on the interval
30. (β³).
Estimate how many subdivisions are needed to calculate to five decimal places of accuracy using the trapezoidal rule.
Note that if is the error using subintervals, then where is the maximum absolute value of the second derivative of the function being integrated and and are the limits of integration.
31.
Let be a function with
β18β
For example, will do, but you donβt need to know what is for this problem.
- Show that
whenever is in the interval - Find the maximum value of
over the interval - Assuming
how many intervals should you use to approximate to within - Using the value of
you found in (b), how many intervals should you use to approximate to within
32.
Approximate the function with a rational function by approximating the integral using Simpsonβs rule. Your rational function should approximate with an error of not more than 0.1 for any in the interval
33.
You may assume use without proof that You may use a calculator, but only to add, subtract, multiply, and divide.