7. Numerical Differentiation and Integration
7.3. Newton–Cotes Formulas
is called Richardson’s extrapolation. A similar procedure can be applied also in the case when the Taylor expansion of the truncation error contains all powers ofh(see Exercises 2 and 3), but later we will use the case presented in this section.
Example 7.6. In the previous section we saw that the central difference formula (7.9) is second-order in h. Using Taylor’s method we get a more precise form of the truncation error.
Suppose that f ∈C2m+3, and consider the following Taylor’s expansion:
f(x0+h) =f(x0) +f′(x0)h+· · ·+f(2m+2)(x0)
(2m+ 2)! h2m+2+f(2m+3)(ξ1)
(2m+ 3)! h2m+3.
We apply the previous relation with−hinstead ofh, subtracting the two equations, and solving it forf′(x0) we get:
f′(x0) = f(x0+h)−f(x0−h)
2h −f′′′(x0)
3! h2−f(5)(x0) 5! h4
− · · · − f(2m+1)(x0)
(2m+ 1)! h2m−f(2m+3)(ξ1) +f(2m+3)(ξ2)
(2m+ 3)! h2m+2.
Hence we have that the central difference satisfies relation (7.21). Therefore, we get a higher order formula using Richardson’s extrapolation. We have that formula
K(1)(h) =
4f(x0+h/2)−f(x0+h/2)
h −f(x0+h)−f(x0−h) 2h
3
= f(x0−h)−8f(x0−h/2) + 8f(x0+h/2)−f(x0+h) 6h
has fourth-order error in h. We note that this formula is equivalent to (7.11).
□
Exercises
1. Derive a sixth-order approximation formula for the first derivative of a function starting from the central difference formula (7.9) using the Richardson’s extrapolation. Apply the formula for approximating the first derivative of f(x) = exsinx atx = 0 using step size h= 0.25.
2. Reformulate the Richardson’s extrapolation for the case when the Taylor expansion of the truncation error contains all powers ofh, i.e.,
M =K(h) +a1h+a2h2+· · ·+amhm+b(x), where|b(h)| ≤Bhm+1 with someB >0.
3. Reformulate the Richardson’s extrapolation for the general case when M =K(h) +a1hα1 +a2hα2+· · ·+amhαm+b(x),
where 1≤α1< α2 <· · ·< αm are integers, and|b(h)| ≤Bhαm+1 with someB >0.
4. Derive a third-order approximation of the first derivative using Richardson’s extrapolation starting from the first-order difference formula.
7.3. Newton–Cotes Formulas
Let f ∈ C[a, b]. The definite integral, similarly to the derivative, is defined by a limit.
The definition using Riemann’s sum is the following: consider a finite partition of the interval [a, b] using the mesh pointsa =x0 < x1 <· · ·< xn =b, and in each subinterval [xi−1, xi] select a point ξi. Then the integral ∫︁b
a f(x)dx is a limit of the Riemann’s sum
∑︁n
i=1f(ξi)(xi −xi−1) as the norm of the partition, max{xi−xi−1: i = 1, . . . , n} goes to zero. Such a Riemann’s sum is for example
∫︂ b rule. (See Exercises 5 and 6.)
Similarly to the numerical differentiation, we can use the Lagrange’s method to derive approximation formulas for definite integrals. Consider a partition of the interval [a, b]
(typically with equidistant mesh points), and let Ln be the Lagrange interpolating poly-nomial of the function f corresponding to the given mesh. Consider ∫︁b
aLn(x)dx as an approximation of ∫︁b
a f(x)dx. We suppose that f ∈Cn+1[a, b]. Then Theorem 6.5 yields the error of the approximation:
∫︂ b
wherelk(x) (corresponding to the mesh points) is the Lagrange basis polynomial of degree n defined by (6.2). Here we get an approximation formula of the form
∫︂ b
where the weights ck are defined by ck =
∫︂ b a
lk(x)dx. (7.30)
Approximation formulas of the form (7.29) are called quadrature formulas. Those quadra-ture formulas when the weights ck are defined by the integrals (7.30) are calledNewton–
Cotes formulas. If the end points of the interval a and b belong to the mesh points, then formulas (7.29)–(7.30) are called closed Newton–Cotes formulas, and if all mesh points belong to the open interval (a, b), then they are called open Newton–Cotes formulas.
We say that the degree of precision of a quadrature formula is n if the formula gives back the exact value of the definite integral for all polynomials with degree at mostn, and there exists a polynomial of degree n+ 1 for which the quadrature formula is not exact.
Therefore, the degree of precision of the (n+1)-point Newton–Cotes formula (7.29)–(7.30)
7.3. Newton–Cotes Formulas 147 is at least n, since in this case the Lagrange polynomialLn is identical to the function f.
It is possible to show that for evenn the (n+ 1)-point Newton–Cotes formulas are exact for polynomials with degree n+ 1 too.
Next we consider the closed Newton–Cotes formula for n = 1. Letx0 =a,x1 =b and
The error of this formula, according to (7.28), is
∫︂ x1
We obtained the so-called trapezoidal rule:
∫︂ b a
f(x)dx= h
2(f(a) +f(b))− h3
12f′′(ξ), ξ∈(a, b). (7.31) The name of the formula comes from the fact that h2(f(a) +f(b)) gives back the area of the region bounded by the secant line of the function corresponding to the points a and b, the x-axis, and the vertical lines x=a and x=b.
The trapezoidal rule gives a good approximation of the integral if the length of the interval is small. If we have a large interval, then we divide it inton subintervals of equal
length by the mesh pointsxi =a+ih(i= 0,1, . . . , n), whereh= (b−a)/n, and we apply the trapezoidal rule for each subintervals:
∫︂ b
We suppose that f ∈ C2[a, b]. Then it follows from Theorem 2.2 that the average value
1 This formula is called composite trapezoidal rule.
Example 7.7. We compute approximate values of the integral ∫︁1
0 x2exdxusing the basic or composite trapezoidal rule with h = 1, h = 0.5 and h = 0.25, respectively. It can be checked that the exact value of the integral is ∫︁1
0 x2exdx=e−2 = 0.7182818 (with 7 digits precision).
For the first case we have
∫︂ 1
0
x2exdx≈ 1
2(0 +e) = 1.3591409,
where we computed the numerical values with 7 digits precision. The error in this case is 0.6408591. With h= 0.5 the composite trapezoidal rule gives
∫︂ 1 0
x2exdx≈ 0.5
2 (0 + 0.52e0.5+e) = 0.8856606.
Hence its error is 0.1673788. Finally, for h= 0.25 we get
∫︂ 1 0
x2exdx≈ 0.25
2 (0 + 0.252e0.25+ 0.52e0.5+ 0.752e0.75+e) = 0.7605963,
so its error is 0.0423145. We can observe that if the step size reduces to its half, then the corresponding error in the approximation reduces to its quarter, which indicates that the error
in his quadratic.
□
Consider formula (7.28) for n = 2 and using equidistant mesh points, i.e., x0 = a,
7.3. Newton–Cotes Formulas 149
Now there is a difference compared to the previous case: the function (x − x0)(x − x1)(x−x2) has opposite signs on the intervals (x0, x1) and (x1, x2), so Theorem 2.6 is not applicable on (x0, x2). We have a different method to simplify the formula for the error term. Let
p is a nonnegative function, hence applying Theorems 2.6 and 6.8, we get
∫︂ x2
We have proved the relation which is called Simpson’s rule.
This error formula yields that the Simpson’s rule is precise for third-order polynomials, since then f(4) is identically equal to 0. On the other hand, the order of approximation inh is five. Similar higher order of precision can be shown for all Newton–Cotes formulas with even n.
Similarly to the composite trapezoidal rule, we can derive the composite Simpson’s rule: We divide the interval [a, b] into 2n equal parts, so let h= (b−a)/2n. Then
Example 7.8. Compute the approximate values of ∫︁1
0 x2exdx using (composite) Simpson’s formula with h= 0.5,h= 0.25 and h= 0.125. First we get
∫︂ 1 0
x2exdx≈ 0.5
3 (0 + 4·0.52e0.5+e) = 0.7278339.
The error is 0.0095520. For h= 0.25 we apply the composite Simpson’s formula:
∫︂ 1 0
x2exdx≈ 0.25
3 (0 + 4·0.252e0.25+ 2·0.52e0.5+ 4·0.752e0.75+e) = 0.7189082.
Its error is 0.0006264. Finally, for h= 0.125 we get
∫︂ 1
which has the error 0.0000396.
□
Next we present some other closed Newton–Cotes formulas.
Simpson’s 38 formula:
7.3. Newton–Cotes Formulas 151 Finally, we present some open Newton–Cotes formulas:
∫︂ x1 We close this section with the investigation of the numerical stability of the integration.
Theorem 7.9. Let ∑︁n
i=1cif(xi) be a quadrature formula which is exact for constant functions and each coefficient ci is positive. Let yi be an approximate value of the exact function value f(xi), and suppose|yi−f(xi)| ≤ε. Then
Proof. According to the assumptions, (b−a) = ∫︁b
a 1dx=∑︁n
We note that all quadrature formulas we presented in this section were exact for constant functions, and most of them had positive weights. Therefore, all such formulas are stable for the rounding error.
Exercises
1. Compute approximate values of the integrals using the trapezoidal rule with step sizes h= 0.5,0.25,0.125, respectively:
(a) ∫︁1
2. Repeat Exercise 1 using the Simpson’s rule.
3. Repeat Exercise 1 using formulas (7.35)-(7.36).
4. Repeat Exercise 1 using formulas Newton–Cotes Formulas (7.37)-(7.40).
5. Prove that the midpoint formula (7.27) gives back the sum of the areas under tangent lines at the midpoints of the intervals [xi, xi+1].
6. Show that the midpoint formula is a Newton–Cotes formula, and derive its error term.
7. Derive formulas (7.35)-(7.36) (without computing the error terms).
8. Derive formulas (7.37)-(7.40) (without computing the error terms).