Thoughts on Numerical Integration (Part 21): Simpson’s rule and global rate of convergence

Numerical integration is a standard topic in first-semester calculus. From time to time, I have received questions from students on various aspects of this topic, including:
  • Why is numerical integration necessary in the first place?
  • Where do these formulas come from (especially Simpson’s Rule)?
  • How can I do all of these formulas quickly?
  • Is there a reason why the Midpoint Rule is better than the Trapezoid Rule?
  • Is there a reason why both the Midpoint Rule and the Trapezoid Rule converge quadratically?
  • Is there a reason why Simpson’s Rule converges like the fourth power of the number of subintervals?
In this series, I hope to answer these questions. While these are standard questions in a introductory college course in numerical analysis, and full and rigorous proofs can be found on Wikipedia and Mathworld, I will approach these questions from the point of view of a bright student who is currently enrolled in calculus and hasn’t yet taken real analysis or numerical analysis.
In this previous post in this series, we showed that the Simpson’s Rule approximation of \displaystyle \int_{x_i}^{x_i+2h} x^k \, dx has an error of 

-\displaystyle \frac{k(k-1)(k-2)(k-3)}{90} x_i^{k-4} h^5 + O(h^6).

In this post, we consider the global error when integrating on the interval [a,b] instead of a subinterval [x_i,x_i+2h]. The total error when approximating \displaystyle \int_a^b x^k \, dx = \int_{x_0}^{x_n} x^k \, dx will be the sum of the errors for the integrals over [x_0,x_2], [x_2,x_4], through [x_{n-2},x_n]. Therefore, the total error will be

$latex E \approx \displaystyle \frac{k(k-1)(k-2)(k-3)}{90} \left(x_0^{k-4} + x_2^{k-4} + \dots + x_{n-2}^{k-4} \right) h^5.

So that this formula doesn’t appear completely mystical, this actually matches the numerical observations that we made earlier. The figure below shows the left-endpoint approximations to \displaystyle \int_1^2 x^9 \, dx for different numbers of subintervals. If we take n = 100 and h = 0.01, then the error should be approximately equal to

\displaystyle \frac{9 \times 8 \times 7 \times 6}{90} \left(1^5 + 1.02^5 + 1.04^5 + \dots + 1.98^5 \right) (0.01)^5 \approx 0.0000017,

which, as expected, is close to the observed error of 102.3000018 - 102.3 \approx 0.0000018.
Let y_i = x_i^{k-4}, so that the error becomes

E \approx \displaystyle \frac{k(k-1)(k-2)(k-3)}{90} \left(y_0 + y_2 + \dots + y_{n-2} \right) h^5 = \displaystyle \frac{k(k-1)(k-2)(k-3)}{90} \overline{y} \frac{n}{2} h^5,

where \overline{y} = (y_0 + y_2 + \dots + y_{n-2})/(n/2) is the average of the y_i. (We notice that there are only n/2 terms in this sum since we’re adding only the even terms.) Clearly, this average is somewhere between the smallest and the largest of the y_i. Since y = x^{k-4} is a continuous function, that means that there must be some value of x_* between x_0 and x_{k-2} — and therefore between a and b — so that x_*^{k-4} = \overline{y} by the Intermediate Value Theorem. We conclude that the error can be written as

E \approx \displaystyle \frac{k(k-1)(k-2)(k-3)}{180} x_*^{k-4} nh^5,

Finally, since h is the length of one subinterval, we see that nh = b-a is the total length of the interval [a,b]. Therefore,

E \approx \displaystyle \frac{k(k-1)(k-2)(k-3)}{180} x_*^{k-4} (b-a)h^4 \equiv ch^4,

where the constant c is determined by a, b, and k. In other words, for the special case f(x) = x^k, we have established that the error from Simpson’s Rule is approximately quartic in h — without resorting to the generalized mean-value theorem and confirming the numerical observations we made earlier.

One thought on “Thoughts on Numerical Integration (Part 21): Simpson’s rule and global rate of convergence

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

This site uses Akismet to reduce spam. Learn how your comment data is processed.