Shampine, Allen, Pruess - Fundamentals of Numerical Computing (523185), страница 52
Текст из файла (страница 52)
The expressionis the geometric mean (average value) of the factors(σ - ri )/ri, and the bound merely states that theleast value of a factor is no more than the averagevalue. If some factor is much larger than thesmallest one, the average will be substantiallylarger than the smallest factor and the bound willbe rather poor. If σ is approximately a large root,say rj, then the factor (σ - rj)/ rj is small. But ifsome root, say rk, is much smaller than σ, that is,|rk| >> |σ|, then the factorExercise Set 4.4 (Page 160)4.23 m = 2 and α = 1/3, so the right side of (4.14) is10 - 610 -7.is quite large and the bound is poor.Exercise Set 5.1 (Page 183)From f(x) = x4Exercise Set 4.5 (Page 162)5.1 A1 = A4 =4.25we get d = 3 and c = -2/405.
(The name -rulecomes from the coefficient values for the interval[0,1] instead of [-1,1].)5.3 With h = π/N we haveN-panel trapezoidN0.811148922853600.8111557334221000.811155735194140n012345x0.50000.28080.21380.21130.21140.2114y1.OOOO0.62780.46350.43050.42930.4293z0.00000.03650.05650.06220.06230.0623and A2 = A3 =ANSWERS TO SELECTED EXERCISESNote: the routine Adapt from Section 5.3 requires441 function evaluations for comparable accuracy.Exercise Set 5.3 (Page 189)5.4 (a) This integrand is very smooth, so Adapt has notrouble getting an accurate answer in only 21function calls. The output is Answer= 3.14159122, Flag = 0, Errest = 1.431E-5. Sincethe exact answer is π, the actual error is 1.431E-5,which was what was estimated by Adapt.(d) This integrand is highly oscillatory, yet Adapthas little trouble getting an accurate answer in 49function calls.
The output is ANSWER= 1.5000000, FLAG = 0, ERREST =-1036E-16. Since the exact answer is 3/2,Adapt’s error estimate is quite accurate.(f) This integrand has a vertical asymptote at 0.25but it still is integrable with exact answer of1+Adapt has great difficulty with thisproblem, even when f is defined to be finite, say1.E+8, at 0.25 to avoid the division by zero. Theoutput is ANSWER = 2.7319539, FLAG = 2,ERREST = 9.264E-5.
The actual error is9.264E-5, which is close to the actual error of9.691E-5. The Flag of 2 indicates that Adapt wasunable to compute a sufficiently accurate answer inthe 3577 function calls allowed.5.5 ANSWER = 1.500000 and ERREST = -2.4E-18are comparable to the results from Exercise 5.4d.The number of integrand evaluations was 21,which is faster than for Exercise 5.4d.5.7xAnswerF callsFlag70.00.0000000070.10.0993360070.20.194751000.30.28263170770.40.359943500.4244364070.55.9 T (0.50) = 3.736767, T(0.75) = 3.904747,T (1.00) = 4.187184.2615.11 For a machine with u = 1.4 × 10-17, the first n = 4and the second is 12.
Exercise 1.3 suggests that forthe first n, 4 digits of accuracy (out of 17) might belost due to the noise. For the second n, 12 digitsmight be lost. This idea was used with integrand(b) of Exercise 5.4. The output ANSWER= 0.909090912, with ERREST = -7.82E-9 thatrequired 49 function calls. The noise level was notmagnified by Adapt, so the code appears to bestable. When n = 12 we got ANSWER= 0.909099105 in 2807 function calls withERREST = -7.82E-6. The accuracy wascontaminated by the noise after the expected 5digits; the number of function evaluations requiredalso increased due to the lack of smoothness.Exercise Set 5.4 (Page 198)5.13PartAnswer(a)(b)(c)0.62053700.62053620.6205366Flag000F calls1473549For the first method in (b) we usedwith ε = 0.23; the first integral was doneanalytically and the second by Adapt.
For (c), weused the change of variable t = x2.5.15 (a) For ABSERR = 10-9 Adapt returnedANSWER = 0.2865295 and an error estimate of-0.34E-9 with 1897 function calls.(c) Method (a) is faster at the larger tolerances andis more accurate at ABSERR = 10-9.5.17 The singularity at x = 0 causes no problems.t1.02.03.0E1(t) estimate0.21938380.04890040.0130475262ANSWERS TO SELECTED EXERCISESMiscellaneous Exercises for Chapter 5 (Page208)5.195.25 For convenience in notation, letSincewe have, after some algebraThis function has the formit is mostefficiently evaluated as a reverse iteration. Letfor i = 5, 4, .
. . , 1 with the integral for F(T6) takenover the interval [0, 1900/T 6] to get the iterationstarted. Then Cv(T) = g(T)F(T).Exercise Set 5.5 (Page 201)5.21MethodGlucoseInsulinGlucagonCr. Hormone(a)6930.91905.192319.4256.5Exercise Set 5.6 (Page 203)5.24 (a) Al = A2 = A3 = 1/6; hence,(b)6927.11904.692688.3256.7which should have better numerical properties for rboth thenear 1. This is because asnumerator f(θ´) - f(θ) and the denominator1 - 2rcos(θ − θ´) + r2 become small. Since theintegral acts as a small correction to the term f(θ)and since the integrand is of one sign, thenumerator must balance out the effect of the smalldenominator. Moreover, even if the integral isgrossly in error, it has little effect on the result ofthe sum asTo illustrate this, let θ = π/2with f(θ) = sin θ so that the solution becomesφ(r,θ) = r.
The function φ was evaluated usingboth forms, the original and the modified, andfunction counts were used as a measure of thecomputational effort. For example, atr = 0.9921875 we got for the original φ, a value of0.99218012 requiring 1057 function evaluations.The modified φ value was 0.99218750 requiring329 function evaluations. The modified version isclearly faster and a little more accurate.5.27 For n = 3 we have qn = 7.08638 so Adapt returned0.000049746 with Errest =-7.21E-9; hence, An = 0.011046.5.29 The main program is the usual driver to producethe answer T0.
The function f for Zero hasindependent variable T0 and outputwhere ANSWER is the output from Adapt that iscalled inside F since the upper limit of integrationANSWERS TO SELECTED EXERCISESis T0. The result is T0 = 246.73.Exercise Set 6.1 (Page 216)6.1 Clearly y(0) = 0 and y´(x) = 0 for 0 < x < c andy´ = (x - c)/2 for c < x < b. The limits asare the same in both cases, so y´ is continuous on[0, b]. Also, for 0 < x < c,and for c < x < b,so that y(x) satisfies the differential equation.6.3 Recall that Dawson’s integral is the expressionNow, y(0) = 0 and differentiation shows thatObviously, y´ is continuous on [0,b] for any finite b.6.5 (a) fy = 2y; the partial derivative is not bounded forall y and so does not satisfy a Lipschitz conditionfor 0 < x < π/2.(c) fy = 1/x; since |fy| < 1/1 = 1 for x in [1,2], thefunction f does satisfy a Lipschitz condition withconstant L = 1.(e) fy = cos x cos y; since |fy| < 1.1 for any x andy, the function f does satisfy a Lipschitz conditionwith constant L = 1.6.7 (a) One solution is Y1 = u, Y2 = u´, Y3 = u´´, andY4 = u´´´ ; then Y´1 = Y2, Y´2 = Y3, Y´3 = Y4, andY´4 = cos αt + tY1 - etY2.(c) One solution is Y1 = u, Y2 = u´, and Y3 = v; thensolving the system for u´´ and v´, we have263(e) One solution is Yl = u, Y2 = u´, Y3 = u´´,Y4 = u(3), Y5 = u(4), and Y6 = u(5) ; then Y´1 = Y2,Y´2 = Y3, Y´3 = Y4, Y´4 = Y5, Y´5 = Y6, andY´6 = et - Y1Y2.Exercise Set 6.2 (Page 221)6.9 For h = l/40 the error at 0.5 is 0.0021; at 1.0 theerror is 0.0010.
For h = l/80 the error at 0.5 is0.0023; at 1.0 the error is 0.00116. For both xvalues the error does drop by a factor of 2 as h ishalved. To achieve an absolute error of 10-6 inmagnitude would require h 0.000011.Exercise Set 6.3 (Page 227)6.11 For h = l/40 the error at 0.5 is -0.l434E-4; at1.0 the error is -0.l395E-4. For h = l/80 theerror at 0.5 is -0.3564E-5; at 1.0 the error is-0.3471. For both x values the error does drop bya factor of 4 as h is halved. To achieve an absoluteerror of 10-6 in magnitude would requireh0.0067.6.13 To show that y(r) can be computed via the formulay(r) = Pry + Qr we use induction. The result is truefor r = 2 sincey´´= P´1y + P1y´ + Q´1= P´1y + P1(P1y + Q1) + Q´1= (P´1 + P1 · P1)y + (Q´1 + Q1P1)= P2 y + Q 2 .Assume it to be true for r = k; theny(k+1) = P´ky + Pky´ + Q´k= P´ky + Pk(P1y + Ql) + Q´k= (P´k + P1Pk)y + (Q´k + Q1Pk)=Pk+ly + Qk+l.A fifth order formula is obtained by dropping theremainder term in the expressionHence Y´1 = Y2, Y´2 = 2t + (3 cos t)/4 - 7Y 1/4 - Y3,and Y´3 = 2t - (cos t)/4 - 3Y 1 /4.264ANSWERS TO SELECTED EXERCISESwhere hn = xn+l - xn,y´ = 1 - 2 xy, andY(r) = Pry + Qr.
The polynomials Pr and Qr aregiven recursively byP r = P´ r - 1 + P 1 P r - 1 ,Q r = Q´ r - 1 , + Q 1 P r - 1 ,with Pl = -2x and Ql = l.Exercise Set 6.6 (Page 238)6.16 (a) The solution is quite smooth on [0,3] and Rkecomputes an accurate solution very efficiently. Theoutput y = 0.49999991 with 8 calls to Rke hasabsolute error 0.00000009 well within therequested tolerance.(d) The large Lipschitz constant of 100 is a cluethat this problem may be difficult.
Rke producesy = 0.83598456 requiring 63 calls. Thge actualglobal error is -0.0000002, which is within therequested tolerance. This problem is “stiff” in thesense of Section 6.7.6.17 The equation x´´ + (x2 - 1)x´ + x = 0 is of the formtreated by Lienhard with f(x) = x2 - 1. Theindefinite integral G(x) = x3/3 - x, so the Liénardvariables areY1(t) = x(t)Y 2 (t) = x´(t) + G(x(t))6.19 The results are comparable to those ofExercise 5.15. For example, P(12) = 778.030,P(21) = 1112.06, and P(30) = 1246.86.6.21 Note that the equation can be written in the formIntegrate this to obtainy´´ - y´ sinx - y cos x = x ln x - x + c2,the first integral relation.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
