Shampine, Allen, Pruess - Fundamentals of Numerical Computing (523185), страница 35
Текст из файла (страница 35)
Generally we have inmind {xi } that lie in [a,b], but in some important applications this is not the case. Forexample, the very important Adams formulas for the solution of ordinary differentialequations are based on quadrature rules that use the end points a and b as nodes, butall other xi lie outside the interval.
The same is true of a method for integrating tabulardata that is considered later in this chapter.The following theorem develops some bounds on the error of a formula with degree of precision d. It is stated using the notation ||f|| for the maximum over [a, b] of|f(x)|.
Also, as in Chapter 3, Mq is used for ||f(q)||. Finally, the absolute error of thequadrature formula is denoted by E(f), that is,Theorem 5.1.If the quadrature formula (5.2) has degree of precision d, thenfor any polynomial p(x) of degree q < d,(5.4)If each Ai > 0, then(5.5)ProofFor p(x) any polynomial of degree q < d,5.1 BASIC QUADRATURE RULES173where we have used the fact that E(p) = 0. (Why?) This is (5.4). For (5.5), when eachAi > 0 the absolute values in (5.4) can be dropped. Because the quadrature formulaintegrates constants exactly, applying it to f(x) = 1 shows thatwhich results in (5.5).Corollary 5.1.If f(x) has d + 1 continuous derivatives on [a, b], then(5.6)If each Ai > 0, then(5.7)Proof Since the bounds of Theorem 5.1 hold for any p(x), we can use the p(x)coming from Taylor’s theorem (see the appendix) with x0 = (a + b)/2 and n = q:whereandfor some z between x0 and x.
This implies that(5.8)Substituting this with q = d into (5.4) or (5.5) yields (5.6) or (5.7).nWhen we studied polynomial interpolation, we learned that interpolants of highdegree are likely to oscillate and provide unsatisfactory fits. The situation now is different because it is the area under the curve that is being approximated, and it seems at174CHAPTER 5NUMERICAL INTEGRATIONleast possible that the oscillations will average out. This is the importance of the special case of formulas with all Ai > 0.
At least as far as the bound of the Theorem 5.1goes, increasing the degree of precision with such formulas can only help. Unfortunately, the interpolatory quadrature formulas forbased on {xi } equallyspaced in [a, b], which are called Newton-Cotes formulas, have some Ai that are negative for even modest degrees. The results of these formulas may not converge to thevalue of the integral as the degree is increased.
However, we shall take up another family of formulas of arbitrarily high degree of precision for which all the Ai are positive.In the bounds (5.4), (5.5) we can choose any polynomial p(x) of any degree q < d.For finite a, b there is a polynomial p*(x) of degree at most d that is as close as possibleto f in the sense thatThe code that accompanies this chapter, called Adapt, is based on two formulas withAi > 0 for all i. One has d = 5 and the other, d = 11. In the bound (5.5), in the one casewe have ||f - p*|| for the best possible approximation by a polynomial p*5 of degree 5and in the other, p*11 of degree 11.
According to the bound, the formula of degree 11cannot be worse because the polynomial of degree 5 can be considered a polynomialof degree 11 with some zero coefficients. It would be remarkable if it were not the casethat ||f - p*11|| is quite a bit smaller than ||f - p*5||, hence that the formula of degree 11is quite a bit more accurate than the formula of degree 5.A more detailed analysis of the error shows that for many formulas with w(x) = 1,including all those considered in this book, the error E(f) can be expressed as(5.9)for some ξ in (a,b). Note that this is an equality rather than a bound.The result (5.9) is a traditional one, but when it involves a derivative of high order,it causes people to doubt whether the formula is practical.
For instance, the formulaof degree 11 mentioned earlier satisfies (5.9) with f(12). It is hard to come to anyunderstanding of such a high order derivative, and a natural question is, What happensif you use the formula and the derivative does not exist? We appreciate now thatthe form (5.9) is just a consequence of the method of analysis. The inequality (5.8)provides bounds when f has only q + 1 continuous derivatives, and the bound basedon best approximation does not directly assume any continuous derivatives. There isno reason to fear a formula of high degree of precision because of an expression like(5.9); other expressions for its error are applicable when the function is less smooth.If a quadrature formula has the degree of precision d, then(5.10)and(5.11)5.1 BASIC QUADRATURE RULES175If we assume that the error has the form (5.9), it is easy to find c fromThe equations (5.10), (5.11) furnish another way to generate quadrature rules.
Theapproach is known as the method of undetermined coeficients. In this approach thecoefficients {Ai } are regarded as unknowns that are found by satisfying the system oflinear equations (5.10) for d as large as possible. Before giving examples, we note thatIt is often convenient to apply the method of undetermined coefficients to a standardinterval [-1,l] and then transform to a general interval [a,b] by a simple change ofvariable.
If we haveletThen dt = (b - a)dx/2 andSinceit follows thatso that the change of variable yieldsExample 5.1.We seek a quadrature formula of the formIn the method of undetermined coefficientsf(x) = 1 implies 2 = A1 + A2f(x) = x implies 0 = -A1 + A2.176CHAPTER 5NUMERICAL INTEGRATIONFigure 5.1 trapezoid rule.Hence, A1 = A2 = 1. We also observe that, by construction, d > 1. Then f(x) = x2yieldsSince E(x2)0 this tells us that d = 1 and c = E(x2)/2 ! = -2/3, which is to say thatfor some ξ in (-1,l).For a general interval [a, b] we apply the result of the change of variable formula(5.12).
We have (in braces) thetrapezoid rule:for some ξt in (a,b). The name of the quadrature rule comes from the fact that itamounts to approximating the area under the curve f(x) by the area of a trapezoid.See Figure 5.1 for an illustration.nExample 5.2.Let us find the most accurate formula of the form5.1 BASIC QUADRATURE RULES177In the method of undetermined coefficients, we try to integrate powers of x exactly toas high a degree as possible.This set of equations determines the coefficients to bethe degree of precision d, we check the next higher power,To findf(x) = x3 implies 0 = -Al + A3 + E(x3) = E(x3)and find that d is at least 3, higher than we might have expected.
For the next power,so d = 3 and E(x4 ) = (2/5) - (2/3) = -(4/15). Then c = -(4/15)/4! = -l/90 andfor some ξ. As in (5.12), for a general [a,b] this gives usSimpson’s rule:See Figure 5.2 for an illustration.nBoth these formulas are Newton-Cotes formulas because the nodes are equallyspaced in [a,b]. The procedure was to select the nodes {xi } and then solve a systemof linear equations for the {Ai }. This is typical when the nodes are specified in advance.
But what if the {xi } are allowed to be unknown as well? With twice as manyunknowns, 2N, at our disposal, we might hope to find formulas with a much higherdegree of precision, perhaps even 2N - 1; that is, we might hope to get a (2N - 1)stdegree formula that uses only N evaluations of f. Unfortunately, the system of equations for {Ai } and {xi } is then not linear. It is not obvious that there are real solutions,and if there are, how to obtain them. Gauss elegantly solved the problem for generalN, even with rather general weight functions and infinite intervals.
The resulting formulas are known as Gaussian quadrature formulas. Special cases can be worked outin an elementary way.Example 5.3.For N = 1 the Gaussian formula has the form178CHAPTER 5NUMERICAL INTEGRATIONFigure 5.2 Simpson’s rule.In the method of undetermined coefficients,f(x)f(x)==1ximpliesimplies20==A1A1x1,hence, A1 = 2 and x1 = 0. To determine the error, we tryf(x) = x2 impliesand find d = 1, c == 2 × 0 + E(x2),andOn [a,b] this becomes(5.14)This formula is known as the midpoint rule.Example 5.4.nFor N = 3, the form isBecause of the symmetry of the interval [-1,1], it is plausible that A1 = A3, x2 = 0,and x1 = -x3, so let us try5.1 BASIC QUADRATURE RULES179Now,f(x)f(x)f(x)f(x)f(x)==1ximpliesimplies20==2A1A1x1+-A2A1x1=0==2xx3impliesimplies0=A1x31-A1x31=0=4impliesx(automatic)(automatic)At this point we have three equations in the three unknowns A1, A2, and x1. The lasttwo equations require that x21 = 3/5, A1 = 5/9 and the first that A2 = 8/9.
To find theerror, we tryf(x) = x5 implies 0 = A1x51 - A1x51 + E(x5) = E(x5).This implies that d > 5. Finallyf(x) = x6 impliesThis says that d = 5 and c =Collecting the results,On [a, b] the resulting quadrature rule is called thethree-point Gaussian quadrature formula,(5.15)and its error is(5.16)See Figure 5.3 for an illustration.nFor larger N the method of undetermined coefficients is impractical for derivingGaussian quadrature rules. Besides, the questions of the existence of formulas and thebest possible degree of precision are left open in this approach. Gauss used the theoryof orthogonal polynomials to answer these questions. We cannot develop the theoryhere (see [8, pp. 327-331]), but it is possible to see how high degrees of precision can180CHAPTER 5NUMERICAL, INTEGRATIONFigure 5.3 Three-point Gaussian quadrature.be achieved. With reasonable conditions on w(x) and [a,b], it is known that there is asequence of polynomials θN+1 (x), N = 0, 1, .
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
