Heath - Scientific Computing (523150), страница 73
Текст из файла (страница 73)
In particular, show that the odd-order terms dropout in both cases, as claimed.8.2 (a) Using the composite midpointquadrature rule, computeR 1 3 the approximatevalue for the integral 0 x dx, using a meshsize (panel width) of h = 0.5 and also using amesh size of h = 1.8.5 Suppose that Lagrange interpolation ata given set of nodes x1 , . . . , xn is used to derive a quadrature rule. Prove that the corresponding weights are given by the integrals ofRbthe Lagrange basis functions, wi = a li (x) dx,i = 1, .
. . , n.(b) Based on the two approximate values computed in part a, use Richardson extrapolationto compute a more accurate approximation tothe integral.8.6 Let p be a real polynomial of degree nsuch thatZ bp(x)xk dx = 0, k = 0, . . . , n − 1.(c) Would you expect the extrapolated resultcomputed in part b to be exact in this case?Why?a(a) Show that the n zeros of p are real, simple, and lie in the open interval (a, b).
(Hint:270CHAPTER 8. NUMERICAL INTEGRATION AND DIFFERENTIATIONConsider the polynomial qk (x) = (x − x1 )(x −x2 ) · · · (x − xk ), where xi , i = 1, . . . k, are theroots of p in [a, b].)(b) Show that the n-point interpolatoryquadrature rule on [a, b] whose nodes are thezeros of p has polynomial degree 2n−1. (Hint:Consider the quotient and remainder polynomials when a given polynomial is divided byp.)8.7 Newton-Cotes quadrature rules are derived by fixing the nodes and then determiningthe corresponding weights by the method ofundetermined coefficients so that the polynomial degree is maximized for the given nodes.The opposite approach could also be taken,with the weights fixed and the nodes to be determined.
In a Chebyshev quadrature rule, forexample, all of the weights are taken to havethe same constant value, w, thereby eliminating n multiplications in evaluating the resulting quadrature formula, since the single weightcan be factored out of the summation.(a) Use the method of undetermined coefficients to derive a three-point Chebyshevquadrature rule on the interval [−1, 1].(b) What is the polynomial degree of the resulting rule?8.8 In approximating the first derivative ofa function f : R → R, the forward differenceformulaf (x + h) − f (x)f 0 (x) ≈hand the backward difference formulaf (x) − f (x − h)f 0 (x) ≈hare both first-order accurate, meaning thattheir dominant error terms are O(h).
Showhow these two formulas can be combined toproduce a difference approximation for thefirst derivative of f that is second-order accurate, i.e., whose dominant error term is O(h2 ).8.9 Given a sufficiently smooth functionf : R → R, use Taylor series to derive a secondorder accurate, one-sided difference approximation to f 0 (x) in terms of the values of f (x),f (x + h), and f (x + 2h).8.10 Consider the following two methods forapproximating the second derivative of a function f at a point x:1. Evaluate the finite difference quotientf (x + h) − 2f (x) + f (x − h).h22. Interpolate f at the points x − h, x, andx + h by a quadratic polynomial p(x) andthen evaluate p00 (x).Do these two methods produce the same result? Why?8.11 Suppose that the first-order accurate,forward difference approximation to thederivative of a function at a given point produces the value −0.8333 for h = 0.2 and thevalue −0.9091 for h = 0.1.
Use Richardsonextrapolation to obtain a better approximatevalue for the derivative.8.12 Archimedes approximated the value of πby computing the perimeter of a regular polygon inscribing or circumscribing a circle of diameter 1. The perimeter of an inscribed polygon with n sides is given bypn = n sin(π/n),and that of a circumscribed polygon byqn = n tan(π/n),and these values provide lower and upperbounds, respectively, on the value of π.(a) Using the Taylor series expansions for thesine and tangent functions, show that pn andqn can be expressed in the formpn = a0 + a1 h2 + a2 h4 + · · ·andq n = b 0 + b1 h 2 + b2 h 4 + · · · ,where h = 1/n. What are the true values ofa0 and b0 ?(b) Given the values p6 = 3.0000 and p12 =3.1058, use Richardson extrapolation to produce a better estimate for π.
Similarly, giventhe values q6 = 3.4641 and q12 = 3.2154, useRichardson extrapolation to produce a betterestimate for π.COMPUTER PROBLEMS271Computer Problems(c)8.1 Since1Z0one can compute an approximate value for πusing numerical integration of the given function.(a) Use the midpoint, trapezoid, and Simpsoncomposite quadrature rules to compute the approximate value for π in this manner for various stepsizes h. Try to characterize the erroras a function of h for each rule, and also compare the accuracy of the rules with each other(based on the known value of π). Is there anypoint beyond which decreasing h yields no further improvement? Why?(b) Implement Romberg integration and repeat part a using it.(c) Compute π again by the same method,this time using a library routine for adaptive quadrature and various error tolerances.How reliable is the error estimate it produces?Compare the work required (integrand evaluations and elapsed time) with that for parts aand b.(d ) Compute π again by the same method,this time using Monte Carlo integration withvarious numbers n of sample points.
Try tocharacterize the error as a function of n, andalso compare the work required with that forthe previous methods. For a suitable randomnumber generator, see Section 13.5.8.2 The integral in the previous problem israther easy. Repeat the problem, this timecomputing the more difficult integralZ1√0Z4dx = π,1 + x24x log(x) dx = − .91p|x| dx−1Try several composite quadrature rules for various fixed mesh sizes and compare their efficiency and accuracy. Also, try one or moreautomatic adaptive quadrature routines usingvarious error tolerances, and again compare efficiency for a given accuracy.8.4 Use numerical integration to verify or refute each of the following conjectures.(a)1Z√x3 dx = 0.40(b)1Z1dx = 0.41 + 10x20(c)Z0122e−9x + e−1024(x−1/4)√dx = 0.2π(d )Z01050dx = 0.5π(2500x2 + 1)(e)Z100Z10−91pdx = 26|x|(f )25e−25x dx = 10(g)Z1log(x) dx = −108.3 Evaluate each of the following integrals.(a)Z111dx1 + 100x2cos(x) dx−1(b)Z−18.5 Each of the following integrands is defined piecewise over the indicated interval.
Usean adaptive quadrature routine to evaluateeach integral over the given interval. For thesame overall accuracy requirement, comparethe cost of evaluating the integral using a single subroutine call over the whole interval withthe cost when the routine is called separately272CHAPTER 8. NUMERICAL INTEGRATION AND DIFFERENTIATIONin each appropriate subinterval. Experimentwith both loose and strict error tolerances.(a)f (x) =0 ≤ x < 0.30.3 ≤ x ≤ 1018.7 The intensity of diffracted light near astraight edge is determined by the values ofthe Fresnel integrals 2Z xπtdtcosC(x) =20and(b)f (x) =0≤x<e−2e−2≤x≤11/(x + 2)0(c)f (x) =exe1−x−1≤x<00≤x≤2(d )f (x) =e10xe10(1−x)− 1 ≤ x < 0.50.5 ≤ x ≤ 1.58.6 Evaluate the following quantities usingeach of the given methods:(a) Use an adaptive quadrature routine toevaluate each of the integralsIk = e−1Z1xk ex dx0for k = 0, 1, .
. . , 20.(b) Verify that the integrals just defined satisfythe recurrenceIk = 1 − kIk−1 ,and use it to generate the same quantities,starting with I0 = 1 − e−1 .(c) Generate the same quantities using thebackward recurrenceIk−1 = (1 − Ik )/k,beginning with In = 0 for some chosen valuen > 20. Experiment with different values of nto see the effect on the accuracy of the valuesgenerated.(d ) Compare the three methods with respectto accuracy, stability, and execution time. Canyou explain these results?S(x) =Z0xsinπt22dt.Use an adaptive quadrature routine to evaluate these integrals for enough values of x todraw a smooth plot of C(x) and S(x) over therange 0 ≤ x ≤ 5.
You may wish to check yourresults by obtaining a routine for computingFresnel integrals from a special function library(see Section 7.4.1).8.8 The period of a simple pendulum is determined by the complete elliptic integral ofthe first kindZ π/2dθpK(x) =.01 − x2 sin2 θUse an adaptive quadrature routine to evaluate this integral for enough values of x todraw a smooth plot of K(x) over the range0 ≤ x ≤ 1. You may wish to check your results by obtaining a routine for computing elliptic integrals from a special function library(see Section 7.4.1).8.9 The gamma function is defined byZ ∞Γ(x) =tx−1 e−t dt, x > 0.0Write a program to compute the value of thisfunction from the definition using each of thefollowing approaches:(a) Truncate the infinite interval of integration and use a composite quadrature rule, suchas trapezoid or Simpson.
You will need todo some experimentation or analysis to determine where to truncate the interval, based onthe usual trade-off between efficiency and accuracy.(b) Truncate the interval and use a standardadaptive quadrature routine. Again, explorethe trade-off between accuracy and efficiency.(c) Gauss-Laguerre quadrature is designed forthe interval [0, ∞] and the weight function e−t ,COMPUTER PROBLEMS273so it is ideal for approximating this integral.Look up the nodes and weights for GaussLaguerre quadrature rules of various orders(see [1, 251, 282], for example) and computethe resulting estimates for the integral.(d ) If available, use an adaptive quadratureroutine designed for an infinite interval of integration.For each method, compute the approximatevalue of the integral for several values of x inthe range 1 to 10.
Compare your results withthe values given by the built-in gamma function or with the known values for integer arguments,Γ(n) = (n − 1)! .How do the various methods compare in efficiency for a given level of accuracy?8.10 Planck’s theory of blackbody radiationleads to the integralZ ∞x3dx.ex − 10Evaluate this integral using each of the methods in the previous exercise, and compare theirefficiency and accuracy.8.11 In two dimensions, suppose that thereis a uniform charge distribution in the region−1 ≤ x ≤ 1, −1 ≤ y ≤ 1. Then, with suitablychosen units, the electrostatic potential at apoint (x̂, ŷ) outside the region is given by thedouble integralΦ(x̂, ŷ) =Z1−1Z1−1dx dyp(x̂ − x)2 + (ŷ − y)2.Evaluate this integral for enough points (x̂, ŷ)to plot the Φ(x̂, ŷ) surface over the region2 ≤ x̂ ≤ 10, 2 ≤ ŷ ≤ 10.8.12 Using any method you choose, evaluatethe double integralZZe−xy dx dy8.13 (a) Write an automatic quadrature routine using the composite Simpson rule.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
