Heath - Scientific Computing (523150), страница 72
Текст из файла (страница 72)
8.5, because the lowest-orderterm in h is linear.F...F...1.0.......................................... .......... ............................................ ............................................................................................ ....................................................................................................................................................................0.5 •01.0 •extrapolated valuecomputed values••0.250.5.....................
.......................................................................................................................................................................................................................................................................................................................................h0.5extrapolated value••computed values0π/4π/2hFigure 8.5: Richardson extrapolation in Examples 8.7 (left) and 8.8 (right).Example 8.8 Romberg Integration. As another example of Richardson extrapolation,we evaluate the integralZ π/2sin(x) dx.0If we use the composite trapezoid rule, we recall from Section 8.4.1 thatF (h) = a0 + a1 h2 + O(h4 ),which means that p = 2 and r = 4 in this case. With h = π/4, we obtain the valueF (h) = 0.948059.
With h = π/2 (i.e., q = 2), we obtain the value F (h) = 0.785398. Theextrapolated value is then given byF (π/4) − F (π/2)22 − 10.948059 − 0.785398= 0.948059 += 1.00228,4−1which is substantially more accurate than either value previously computed (the exactanswer is 1). In this example the extrapolation is quadratic, as can be seen on the right inFig. 8.5, because the lowest-order term in h is quadratic.Evaluation of the trapezoid rule for additional values of h would permit further extrapolations to attain even higher accuracy, up to the limit imposed by the arithmetic precision.Continued use of Richardson extrapolation in this manner, using the trapezoid quadraturerule with various stepsizes, is called Romberg integration.
It is capable of producing veryhigh accuracy for well-behaved problems.F (0) = a0 = F (π/4) +2668.9CHAPTER 8. NUMERICAL INTEGRATION AND DIFFERENTIATIONSoftware for Numerical Integration and DifferentiationTable 8.1 is a list of some of the software available for numerical quadrature. Most of theone-dimensional quadrature routines listed are adaptive routines based on Gauss-Kronrodquadrature rules. We note that software for solving initial value problems for ordinarydifferential equations, which will be covered in Chapter 9, can also be used for computingdefinite integrals (see Computer Problem 9.5). Several routines are available for generating the nodes and weights for various Gaussian and other quadrature rules, includinggaussq from netlib; and iqpack(#655), extend(#672), and gauss (the latter is part ofthe orthpol(#726) package), all from TOMS.SourceFMMHSLIMSLMATLABKMNNAGNUMALNRQUADPACKSLATECTOMSTable 8.1: Software for numerical integration and differentiationOne dimensionTwo dimensions n dimensionsquanc8qa02/qa04/qa05qb01/qm01qdag/qdagstwodqqandquad/quad8q1dad01ajfd01dafd01fcfquadrattricubvegas/miserdfridrqag/qagsqag/qnc/qng/gaus8squank(#379)/quad(#468) dcutri(#706)dcuhre(#698)Differentiationtd01derivdiffd04aafSoftware for numerical integration typically requires the user to supply the name ofa routine that computes the value of the integrand function for any argument.
The usermust also supply the endpoints of the interval of integration, as well as absolute or relativeerror tolerances. In addition to the approximate value of the integral, the output usuallyincludes an estimate of the error, a status flag indicating any warnings or error conditions,and possibly a count of the number of function evaluations that were required.Although adaptive quadrature routines can often be used as black boxes, they can beineffective for integrals having discontinuities, singularities, or other such difficulties.
Insuch cases, it may be advantageous to transform the problem to enable the automaticroutine to arrive at an accurate result more efficiently. For practical advice on handlingsuch problematic integrals, see [2, 3].In the last column of Table 8.1 are listed some routines for numerical differentiation. Inaddition, a number of packages are available that implement automatic differentiation (seeSection 8.7.2), including ADIC, ADIFOR, ADOL-C, ADOL-F, AMC, GRESS, Odyssée, and PADRE2.See URL http://www.mcs.anl.gov/adifor/ for further information.8.10. HISTORICAL NOTES AND FURTHER READING8.10267Historical Notes and Further ReadingAs mentioned earlier, quadrature is an ancient technique.
Most of the methods we discusseddate from the nineteenth century or earlier, as the names associated with them suggest—Simpson, Newton, Cotes, Gauss, and others. Kronrod published the quadrature rules thatbear his name in 1964. One of the earliest adaptive quadrature routines was published byMcKeeman in 1962. Many others have followed, most notably squank, cadre, qnc7, andquanc8, culminating with the quadpack package, which represents the current state of theart (see also TOMS #691).Comprehensive general references on numerical integration are [52, 70, 74, 153]. Thecomputation of multiple integrals is discussed in [113, 169, 232, 250].
The quadpack packageis documented in [203]. Cautionary advice on using automatic quadrature routines can befound in [168, 170]. For a comprehensive survey of extrapolation techniques, see [141]. Formore details on the numerical solution of integral equations, see [54, 277].Review Questions8.1 True or false: Since the midpointquadrature rule is based on interpolation by aconstant, whereas the trapezoid rule is basedon linear interpolation, the trapezoid rule isgenerally more accurate than the midpointrule.8.2 True or false: The polynomial degree ofa quadrature rule is the degree of the interpolating polynomial on which the rule is based.8.3 True or false: An n-point Newton-Cotesquadrature rule is always of polynomial degreen − 1.8.4 True or false: Gaussian quadrature rulesof different orders never have any points incommon.8.5 How can you estimate the error in aquadrature formula without computing thederivatives of the integrand function thatwould be required by a Taylor series approximation?8.6 (a) If a quadrature rule for an interval[a, b] is based on polynomial interpolation at nequally spaced points in the interval, what isthe highest degree such that the rule integratesall polynomials of that degree exactly?(b) How would your answer change if thepoints were optimally placed to integrate thehighest possible degree polynomials exactly?8.7 Would you expect an n-point NewtonCotes quadrature rule to Rwork well for inte1grating Runge’s function, −1 (1 + 25x2 )−1 dx,if n is very large? Why?8.8 (a) What is the polynomial degree ofSimpson’s rule for numerical quadrature?(b) What is the polynomial degree of an npoint Gaussian quadrature rule?8.9 Newton-Cotes and Gaussian quadraturerules are both based on polynomial interpolation.(a) What specific property characterizes aNewton-Cotes quadrature rule for a givennumber of nodes?(b) What specific property characterizes aGaussian quadrature rule for a given numberof nodes?8.10 (a) Explain how the midpoint rule,which is based on interpolation by a polynomial of degree zero, can nevertheless integratepolynomials of degree one exactly.(b) Is the midpoint rule a Gaussian quadraturerule? Explain your answer.8.11 Suppose that the quadrature ruleZ bnXf (x) dx ≈wi f (xi )ai=1is exact for all constant functions.
What doesthis imply about the weights wi or the nodesxi ?268CHAPTER 8. NUMERICAL INTEGRATION AND DIFFERENTIATION8.12 If the integrand has an integrable singularity at one endpoint of the interval of integration, which type of quadrature rule wouldbe better to use, a closed Newton-Cotes ruleor a Gaussian rule? Why?8.13 What is the polynomial degree of eachof the following types of numerical quadraturerules?(a) An n-point Newton-Cotes rule, where n isodd(b) An n-point Newton-Cotes rule, where n iseven8.18 (a) What is a composite quadraturerule?(b) Why is a composite quadrature rule preferable to an ordinary quadrature rule for achieving high accuracy in numerically computing adefinite integral on a given interval?(c) In using the composite trapezoid quadrature rule to approximate a definite integral onan interval [a, b], by what factor is the overall error reduced if the mesh size (i.e. panelwidth) h is halved?(c) An n-point Gaussian rule8.19 (a) Describe in general terms how adaptive quadrature works.(d ) What accounts for the difference betweenthe answers to parts a and b?(b) How can the necessary error estimate beobtained?(e) What accounts for the difference betweenthe answers to parts b and c?(c) Under what circumstances might such aprocedure produce a result that is seriously inerror?8.14 For each of the following properties,state which type of quadrature, Newton-Cotesor Gaussian, more accurately fits the description:(a) Easier to compute nodes and weights(b) Easier to apply for a general interval [a, b](c) More accurate for the same number ofnodes(d ) Has maximal polynomial degree for thenumber of nodes(e) Nodes easy to reuse as order of rule changes8.15 What is the relationship between Gaussian quadrature and orthogonal polynomials?8.16 (a) What is the advantage of using aGauss-Kronrod pair of quadrature rules, suchas G7 and K15 , compared with using twoGaussian rules, such as G7 and G15 , to obtainan approximate integral with error estimate?(b) How many evaluations of the integrandfunction are required to evaluate both of therules G7 and K15 in a given panel?8.17 Rank the following types of quadraturerules in order of their polynomial degree forthe same number of nodes (1 for highest polynomial degree, etc.):(a) Newton-Cotes(b) Gaussian(c) Kronrod(d ) Under what circumstances might such aprocedure be very inefficient?8.20 What is the most efficient way to usean adaptive quadrature routine for computing a definite integral whose integrand has aknown discontinuity within the interval of integration?8.21 What is a good way to integrate tabulardata (i.e., an integrand whose value is knownonly at a discrete set of points)?8.22 (a) How might one use a standardquadrature routine, designed for integratingover a finite interval, to integrate a functionover an infinite interval?(b) What precautions would need to be takento ensure a good result?8.23 How might one use a standard onedimensional quadrature routine to computethe value of a double integral over a rectangular region?8.24 Why is Monte Carlo not a practicalmethod for computing one-dimensional integrals?8.25 Relative to other methods for numericalquadrature, why is the Monte Carlo methodmore effective in higher dimensions than in lowdimensions?EXERCISES8.26 Explain why integral equations of thefirst kind with smooth kernels are always illconditioned.8.27 Explain how a quadrature rule can beused to solve an integral equation numerically.What type of computational problem results?8.28 In solving an integral equation of thefirst kind by numerical quadrature, does thesolution always improve if the order of thequadrature rule is increased or the mesh sizeis decreased? Why?8.29 List three approaches for obtaining ameaningful solution to an ill-conditioned linearsystem approximating an integral equation ofthe first kind.8.30 Consider the problem of approximatingthe derivative of a function that is measuredor sampled at only a finite number of points.(a) One way to obtain an approximate derivative is to interpolate the discrete data pointsand then differentiate the interpolant.
Is this269a good method for approximating the derivative? Why?(b) Similarly, one can approximate the integral of a function given by such discrete databy integrating the interpolant. Is this a goodmethod for computing the integral? Why?8.31 Comparing integration and differentiation, which problem is inherently better conditioned? Why?8.32 (a) Suggest a good method for numerically approximating the derivative of a function whose value is given only at a discrete setof data points.(b) For this problem, what would be the effectof noisy data, and how would you cope with itin your numerical method?8.33 (a) Explain the basic idea of Richardsonextrapolation.(b) Does it give a more accurate answer thanthe values on which it is based?8.34 What is meant by Romberg integration?Exercises8.1 (a) Compute the approximate value ofR1the integral 0 x3 dx, first by the midpoint ruleand then by the trapezoid rule.(b) Use the difference between these two results to estimate the error in each of them.(c) Combine the two results to obtain theSimpson’s rule approximation to the integral.(d ) Would you expect the latter to be exactfor this problem? Why?Pn8.3 If Q(f ) = i=1 wi f (xi ) is an interpolatory quadrature rule (i.e., based on polynomialinterpolation)Pn on the interval [0, 1], then is ittrue that i=1 wi = 1? Prove your answer.8.4 Fill in the details of the derivation of theerror estimates for the midpoint and trapezoidquadrature rules given in Section 8.2.3.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
