Conte, de Boor - Elementary Numerical Analysis. An Algorithmic Approach (523140), страница 64
Текст из файла (страница 64)
. . , 4) in the Adams-Bashforth formula (8.45).8.7-2 Set m = 4 in (8.45) and derive the corresponding Adams-Bashforth formula in terms ofordinates as in formula (8.47). Also derive the error term for this formula.8.7-3 Derive Milne’s formula (8.51) and its corresponding error term.8.7-4 Write a program using Milne’s formula for integrating a differential equation withequally spaced points. Assume that the first three starting values are known.8.7-5 Solve the equation of Example 8.5 using the Milne program with h = 1/32 and compareyour results with those given in Example 8.5.8.7-6 Solve the equation xy´ = x - y, y (2) = 2 from x = 2 to x = 3 with h = 0.05 using theAdams-Bashforth method (8.47).
Obtain the starting values from the exact solution8.7-7 Using the Adams-Bashforth method (8.47) solve the equation y´ + y = e-x from x = 0to x = 1 using h = 1/64 and h = 1/128. Estimate the accuracy of your results, Starting values canbe obtained from the exact solution y = xe -x.8.7-8 Derive the formulas in (8.51) using (8.49) with m = 2 (not 3) and the error (2.18) inpolynomial interpolation. In this the discussion at the beginning of Sec. 7.2 will be helpful.8.7-9 Verify (8.50) by expanding yn+1 , and yn-1 , about x = xn through third-order terms,assuming that the starting values are exact.8.8PREDICTOR CORRECTOR METHODS3798.8 PREDICTOR-CORRECTOR METHODSThe multistep methods of Sec.
8.7 were derived using polynomials whichinterpolated at the point xn and at points backward from xn . These aresometimes known as formulas of open type. Formulas of closed type arederived by basing the interpolating polynomial on the point xn+1, as wellas on xn and points backward from xn. The simplest formula of this type isobtained if we approximate the integral in (8.43) by the trapezoidalformula (7.26). This leads to the formulan = 0, l, . . .(8.52)3The error of this formula is -(h /12)y´´´ and thus represents an improvement over Euler’s method.
However, (8.52) is an implicit equation for yn+1since yn+1 , appears as an argument on the right-hand side.If f(x,y) is a nonlinear function, we will, in general, not be able tosolve (8.52) for yn+1 exactly. We can, however, attempt to obtain yn+1 bymeans of iteration. Thus, keeping xn fixed, we obtain a first approximationto yn+1 by means of Euler’s formula(8.53)and substitute in the right-hand side ofWe then evaluate(8.52) to obtain the approximationNext we evaluateand again use (8.52) to obtain a nextapproximation.
In general, the iteration is defined byk = 1, 2, . . .(8.54)The iteration is terminated when two successive iterates agree to thedesired accuracy. This iteration for obtaining improved values of yn+1 at afixed point xn+1 is sometimes called an inner iteration to distinguish it from(8.52), which is used to generate values of yn at n = 0, 1, .
. . . We shallsummarize this procedure in Algorithm 8.4.Algorithm 8.4: A second-order predictor-corrector method For the differential equation y´ = f(x,y), y(x0) = y0 with h given and xn = x0 +nh, for each fixed n = 0, 1, . . . :1. Computeusing (8.53).380THE SOLUTION OF DIFFERENTIAL EQUATIONS2.
Compute(k = 1, 2, . . . ), using (8.54), iterating on k untilfor a prescribed εIn specifying e in Algorithm 8.4, we must keep in mind that theaccuracy that can be expected on each step is limited by the error of thebasic formula (8.52) and by the step size h.To adapt this algorithm to the solution of a specific problem, wewould have to specify (a) the number N of steps desired; (b) a maximumnumber K of inner iterations; (c) what to do in case k exceeds K.It is customary to call an explicit formula such as Euler’s formula anopen-type formula, while an implicit formula such as (8.52) is said to be ofclosed type. When they are used as a pair of formulas, the open-typeformula is also called a predictor, while the closed-type formula is called acorrector.
A corrector formula is generally more accurate than a predictorformula, even when both have a discretization error of the same order,primarily because the coefficient in the error term is smaller. Two questions arise naturally in connection with corrector formulas. The first is,“Under what conditions will the inner iteration on k converge?,” and thesecond, “How many iterations will be needed to produce the requiredaccuracy?” The answer to the latter question will depend on many factors.However, if the predictor and corrector formulas are of the same order,experience has shown that only one or two applications of the correctorare sufficient, provided that the step size h has been properly selected.
Ifwe find that one or two corrections are not sufficient, it is better to reducethe step size h than to continue to iterate. The answer to the first questionis contained in Theorem 8.2.Theorem 8.2 If f(x,y) andare continuous in x and y on theclosed interval [a,b] the inner iteration defined by (8.54) will converge,provided h is chosen small enough so that, for x = xn, and all y with(8.55)To prove this, we first observe that in the iteration (8.54) xn is fixed.Hence, if we setwe can write (8.54) in the formY(k) = F(Y(k-1))whereand where C depends on n but not on Y.
This can be viewed as an8.8PREDICTOR CORRECTOR METHODS381instance of fixed-point iteration considered in Sec. 3.3. In a corollary toTheorem 3.1 we proved that such an iteration will converge provided thatF´(Y) is continuous and satisfies|F´(Y)| < 1for all Y with | Y - yn+1 | < | Y(0) - yn+1|, where yn+1 is the fixed point ofF(Y). Since F´(Y) =and sinceis bounded and nonvanishing by assumption, the iteration (8.54) will converge ifi.e., ifSince F´(Y) = (h/ 2 )this proves the theorem.Example 8.6 Solve the equationy´ = x - 1/yy(0) = 1from x = 0 to x = 0.2, using Algorithm 8.4 with h = 0.1.Since the error of (8.54) is -(h 3 /12) y´´´, and since by differentiating above we findthat y´´´(0)-2, the error will be approximately 0.0002. We cannot therefore expectmuch more than three decimal places of accuracy in the results.Step 1By Euler’s method:By (8.54):y 1 (0) = 0.9y 1 (1) = 0.8994y 1 (2) = 0.8994Since y 1 (1) and y1 (2) agree to four places, we accept this answer, and we computey´ 1 = f(x 1 ,y 1 ) = -1.0118.Step 2 By Euler’s method,y 2 (0) = 0.8994 + 0.1(-1.0118) = 0.7982By (8.54),y2(3) = 0.7960We accept y2 = 0.7960, compute y´2 and proceed to the next step.As the computation proceeds, we can expect a gradual loss of accuracy.
It appearshere that for h = 0.1 we need two or three applications of the corrector. This is primarilydue to the fact that we are using a predictor which is of lower order than the corrector.To verify that the inner iterations for this example will converge for h = 0.1, weand hence, from Theorem 8.2, we want h to be less that 2 y2. Wecomputedo not know the solution y, but it is clear from the above steps that y > 0.7 on theinterval [0, 0.2]. Hence the inner iterations will converge if h < 2(0.7)2 = 0.98.382THE SOLUTION OF DIFFERENTIAL EQUATIONSEXERCISES8.8-1 For the special equation y´ = Ay, y(0) = 1, show that the trapezoidal corrector formula(8.52) leads to a difference equation whose solution isprovided that |Ah/2| < 1.8.8-2 For the solution obtained in Exercise 8.8-l show thatfor a fixed value of x = xn = nh.8.8-3 Solve the equation y´ = x2 + y, y(0) = 1, from x = 0 to x = 0.5, using Euler’s methodas a predictor and (8.54) as a corrector.
Determine the step h so that four decimal places ofaccuracy are obtained at x = 0.5. Start with h = 0.05.8.9 THE ADAMS-MOULTON METHODCorrector formulas of higher order can be obtained by using a polynomialwhich interpolates at xn+1, xn . . . , xn-m for an integer m > 0. The Newton backward formula which interpolates at these m + 2 points in terms ofs = (x - xn ) /h is(8.56)These differences are based on the values f n+1 , fn , .
. . , f n-m . If we integrate (8.56) from xn to xn+1 and use (8.43), we obtain(8.57)whereThe first few values of γ´k areThe error of (8.57), based on the error of the interpolating polynomial, is(8.58)The case m = 2 is frequently used. If the differences in (8.57) are expressedin terms of ordinates for m = 2, we obtain(8.59)with the error(8.60)8.9THE ADAMS-MOULTON METHOD383The formula (8.57) is known as the Adams-Moulton formula. The fourthorder Adams-Moulton formula (8.59) is clearly a corrector formula ofclosed type since fn+1 = f(xn+1, yn+1 ) involves the unknown quantity yn+1.It must therefore be solved by iteration. It can be shown that the iterationbased on (8.59) will converge, provided that h is small enough so that theconditionis satisfied.
A convenient predictor to use withthis corrector is the Adams-Bashforth fourth-order formula (8.47). In thiscase the predictor is of the same order as the corrector. If h is properlychosen, then one application of the corrector will yield a significantimprovement in accuracy.Specifications for a fourth-order predictor-corrector method are givenin Algorithm 8.5.Algorithm 8.5: The Adams-Moulton predictor-corrector method Forthe differential equation y´ = f(x,y) with h fixed and xn = x0 + nhand with (y0,f0), (y1,f1), (y2,f2), (y3,f3) given, for each fixed n =3, 4, . . .
:using the formula1. Compute2. Compute3. Computek = 1, 2, . . .4. Iterate on k untilfor ε prescribedAgain this algorithm is not complete unless we specify what to do incase of nonconvergence in step 4. A subroutine like DVERK contains amore complete specification for a general-purpose subroutine to solvedifferential equations.Besides yielding improved accuracy, the corrector formula servesanother useful function. It provides an estimate of the local discretizationerror, which can then be used to decide whether the step h is adequate forthe required accuracy. To examine this error estimation procedure for thepredictor-corrector pair consisting of the Adams-Bashforth and AdamsMoulton fourth-order formulas, we write the local-error estimate for each:(8.61)384THE SOLUTION OF DIFFERENTIAL EQUATIONSLetrepresent the value of yn+1 obtained from (8.47), andtheresult obtained with one application of Algorithm 8.5. If the values off areassumed to be exact at all points up to and including xn , and if y(xn+1 )represents the exact value of y at xn+1, then from (8.61) we obtain the errorestimates(8.62a)(8.62b)In general,However, if we assume that over the interval ofinterest yv(x) is approximately constant, then on subtracting (8.62b ) from(8.62a), we obtain the following estimate for yv:Substituting this into (8.62b ), we find that(8.63)Thus the error of the corrected value is approximately - 1/14 of thedifference between the corrected and predicted values.As mentioned before, it is advisable to use the corrector only once.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
