Conte, de Boor - Elementary Numerical Analysis. An Algorithmic Approach (523140), страница 60
Текст из файла (страница 60)
A particular solution of (a), obtained by tryingynPin (a), is found to be= -12. The characteristic equation of the homogeneous equation of (a) isβ 2 - (2 + h2 )β + 1 = 0By the quadratic formula the roots areOn expanding (1 + t)1/2 around t = 0 into a Taylor series and substituting h 2/4 fort, we obtainHence the general solution of the homogeneous system is3. The solution of (a) is thereforey n = yn P + y n Gwhich establishes the solution in the form (b).EXERCISES8.2-1 Find the general solution of the difference equations(a) y n + l - 3yn = 5(b) y n+2 - 4yn+1 + 4yn = n(Hint: To find a particular solution, try ynP = an + b.
)(c) y n+2 + 2y n + l + 2yn = 0(d) y n + 3 - y n + 2 + 2yn+1 - 2yn = 0(e) y n + 2 - y n + 1 - y n = 08.2-2 Find the solution of the initial-value difference equationsy0 = 0(a) y n+2 - 4yn+1 + 3yn = 2 ny1 = 1y0 = 0(b) y n + 2 - y n + 1 - y n = 0y1 = 1[Hint: To find a particular solution of (a), try ynP = A2 n .]354THE SOLUTION OF DIFFERENTIAL EQUATIONS8.2-3 Show that the general solution of the difference equationy n + 2 + 4hyn+1 - yn = 2hwhere h is a positive constant, can be expressed in the form8.2-4 Show that if y0 = 1, y1 = X, then the nth term, yn = yn(x), of the solution ofy n + 2 - 2xyn+1 + yn = 0is a polynomial of degree n in x with leading coefficient 2 n-1 .
[ Note: The yn(x) are theChebyshev polynomials considered in Sec. 6.1.]8.3 NUMERICAL INTEGRATION BY TAYLOR SERIESWe are now prepared to consider numerical methods for integratingdifferential equations. We shall first consider a first-order initial-valuedifferential equation of the form(8.19)y´ = f(x,y)y(x 0 ) = y 0The function f may be linear or nonlinear, but we assume that f issufficiently differentiable with respect to both x and y. It is known that(8.19) possesses a unique solution ifis continuous on the domain ofinterest.
If y(x) is the exact solution of (8.19), we can expand y(x) into aTaylor series about the point x = x0:(8.20)The derivatives in this expansion are not known explicitly since thesolution is not known. However, if f is sufficiently differentiable, they canbe obtained by taking the total derivative of (8.19) with respect to x,keeping in mind that y is itself a function of x (see Sec.
1.7). Thus weobtain for the first few derivatives:y´ = f(x,y)y´´ = f´ = fx + fyy´ = fx + fyfy´´´ = f´´ = fxx + fxyf + fyxf + fyyf2 + fyfx + fy2f= fxx + 2fxyf + fyyf2 + fxfy + fy2f(8.21)Continuing in this manner, we can express any derivative of y in terms off(x,y) and its partial derivatives. It is already clear, however, that unlessf(x,y) is a very simple function, the higher total derivatives becomeincreasingly complex. For practical reasons then, one must limit thenumber of terms in the expansion (8.20) to a reasonable number, and thisrestriction leads to a restriction on the value of x for which (8.20) is areasonable approximation. If we assume that the truncated series (8.20)8.3NUMERICAL INTEGRATION BY TAYLOR SERIES355yields a good approximation for a step of length h, that is, for x - x0 = h,we can then evaluate y at x0 + h; reevaluate the derivatives y´, y´´, etc., atx = x0 + h; and then use (8.20) to proceed to the next step.
If we continuein this manner, we will obtain a discrete set of values y n which areapproximations to the true solution at the points xn = x0 + nh ( n =0, 1, 2, . . . ). In this chapter we shall always denote the value of the exactsolution at a point xn by y(xn) and of an approximate solution by yn.In order to formalize this procedure, we first introduce the operatork = 1, 2, . . .(8.22)where we assume that a fixed step size h is being used, and where f (j)denotes the jth total derivative of the function f( x,y(x)) with respect to x.We can then state Algorithm 8.1.Algorithm 8.1: Taylor’s algorithm of order k To find an approximatesolution of the differential equationy´ = f(x,y)y(a) = y0over an interval [a, b]:1.
Choose a step h = (b - a)/ N. Setn = 0, 1, . . . , Nxn = a + nh2. Generate approximations yn to y(xn) from the recursiony n+1 = y n + hT k (x n ,y n )where Tk(x, y) is defined by (8.22).n = 0, l, . . . , N - 1Taylor’s algorithm, and other methods based on this algorithm, whichcalculate y at x = xn+1 by using only information about y and y´ at asingle point x = xn, are frequently called one-step methods.Taylor’s theorem with remainder shows that the local error of Taylor’salgorithm of order k isThe Taylor algorithm is said to be of order k if the local error E as definedabove is356THE SOLUTION OF DIFFERENTIAL EQUATIONSOn setting k = 1 in Algorithm 8.1 we obtain Euler’s method and itslocal error,(8.23)To illustrate Euler’s method, consider the initial-value problemy´ = yy(0) = 1On applying (8.23) with h = 0.01 and retaining six decimal places, weobtainy(0.01)y1 = 1 + 0.01 = 1.01y(0.02)y2 = 1.01 + 0.01(1.01) = 1.0201y(0.03)y3 = 1.0201 + 0.01(1.0201) = 1.030301y(0.04)y4 = 1.030301 + 0.0l( 1.030301) = 1.040606Since the exact solution of this equation is y = ex, the correct value atx = 0.04 is 1.0408.
It is clear that, to obtain more accuracy with Euler’smethod, we must take a considerably smaller value for h.If we take h = 0.005, we obtain the valuesy(0.005)y1 = 1.0050y(0.010)y2 = 1.0100y(0.015)y3 = 1.0151y(0.020)y4 = 1.0202y(0.025)y5 = 1.0253y(0.030)y6 = 1.0304y(0.035)y7 = 1.0356y(0.040)y8 = 1.0408These results are correct to four decimal places after the decimal point.Because of the relatively small step size required, Euler’s method is notcommonly used for integrating differential equations.We could, of course, apply Taylor’s algorithm of higher order toobtain better accuracy, and in general, we would expect that the higher theorder of the algorithm, the greater the accuracy for a given step size.
Iff(x,y) is a relatively simple function of x and y, then it is often possible togenerate the required derivatives relatively cheaply on a computer byemploying symbolic differentiation, or else by taking advantage of anyparticular properties the function f(x,y) may have (see Exercise 8.3-4).However, the necessity of calculating the higher derivatives makes Taylor’salgorithm completely unsuitable on high-speed computers for general8.3NUMERICAL INTEGRATION BY TAYLOR SERIES357integration purposes. Nevertheless, it is of great theoretical interest becausemost of the practical methods attempt to achieve the same accuracy as aTaylor algorithm of a given order without the disadvantage of having tocalculate the higher derivatives.
Although the general Taylor algorithm ishardly ever used for practical purposes, the special case of Euler’s methodwill be considered in more detail for its theoretical implications.Example 8.1 Using Taylor’s series, find the solution of the differential equationxy´ = x - yy (2) = 2at x = 2.1 correct to five decimal places.The first few derivatives and their values at x = 2, y = 2 areThe Taylor series expansion about x0 = 2 isy(x) - y0 + (x - 2)y´0 + 1/2 (x - 2)2 y´´0 + 1/6(x - 2)3 y´´´0 + 1/24 (x - 2)4 yiv0 + · · ·= 2 + (x - 2)0 + 1/4(x - 2)2 - 1/8(x - 2)3 + 1/16(x - 2)4 + · · ·At x = 2.1 we obtainy (2.1) = 2 + 0.0025 - 0.000125 + 0.0000062 - · · ·2.00238Since the terms in this Taylor series decrease in magnitude and alternate (see Exercise8.34) in sign, this result is correct to five decimal places. If we now wished to find y (2.2)to the same accuracy, we would have to carry the series through two additional terms.Alternatively, we could now make a new expansion about x = 2.1, reevaluate the firstfour derivatives at x = 2.1, and then compute y(2.2).Example 8.2 Solve the equationfrom x = 1 to x = 2.
Use Taylor’s algorithm of order 2. Solve the problem with h = 1/16,and estimate the accuracy of the results.S OLUTION Since358THE SOLUTION OF DIFFERENTIAL EQUATIONSthenandThe results as computed on the IBM 7094 are given below. The step size h is given in thefirst column, and the values of y(1.5), y´(1.5), y(2.0) y´(2.0) respectively, are given in thenext four columns. The exact solution of this equation is y = -1/x, so that the exactvalue of y(1.5) is -2/3, and the exact value of y (2.0) is -1/2.
We may estimate the totaldiscretization error as follows: The local error of Taylor’s algorithm of order 2 is(h 3 /6)y´´´. Since y´´´ = 6/x4, its maximum value on the interval [1, 2] is 6, and hence thelocal error is for each step, at most, h3. With h = 1/128, we will take 128 integration steps sothat the accumulated error will be, at most, 128h3 = (1/128)2 0.0006. The actual errorat x = 2.0 appears to be 0.00003, in close agreement with this estimate. In general, wewill not know the solution to check against.
Even without knowing the solution,however, we can estimate from the number of places of agreement as h0, theaccuracy of the solution. Since each halving of h appears to produce almost oneadditional digit of accuracy, it appears that in the absence of round-off error, a step of1/1,024 should produce at least seven places of accuracy. This same problem will besolved later by two other methods. For comparison purposes, the results for all threemethods are included here.COMPUTER RESULTS FOR EXAMPLE 8.2Method 1—Taylor expansion method of order 2HY(l.5)0.62500000E-010.31250000E-010.15625000E-010.78125000E-02-066787238E0.66696430E0.66674034E0.66668454E00000000YPRM(l.5)Y(2.)0.44363917E-000.44424593E-000.44439532E-000.44443253E-00-0.50187737E0.50046334E0.50011456E0.50002744EYPRM(2.)000000000.24905779E-000.24976812E-000.24994271E-000.24998628E-00Method 2—Simplified Runge-Kutta order 2HY(l.5)0.62500000E-010.31250000E-010.l5625000E-010.78125000E-02-0.66552725E066637699E066659356E0.66664808E00000000YPRM(l.5)Y(2.)YPRM(2.)0.44520275E-000.44463748E-000.44449317E-000.44445683E-00- 0.49822412E-00- 0.49954852E-00- 0.49988601E-0- 0.49997083E-000.25088478E-000.25022554E-000.25005698E-000.25001458E-00Y(2.)YPRM(2.)Method 3-Classical Runge-Kutta order 4HY(l.5)0.62500000E-010.31250000E-010.15625000E-010.78125000E-02-0.66666625E-066666664E-0.66666666E-066666667EYPRM(l.5)000000000.44444472E-00 -0.49999941E-00 0.250000129E-000.44444446E-00 -0.49999997E-00 0.25000001E-00-0.44444444E-00 -0.50000000E 00 0.25000000E-000.44444444E-00 -0.50000001E 00 0.24999999E-008.4ERROR ESTIMATES AND CONVERGENCE OF EULER’S METHOD359EXERCISES8.3-l For the equationy(1) = 1derive the difference equation corresponding to Taylor’s algorithm of order 3.
Carry out byhand one step of the integration with h = 0.0 1. Write a program for solving this problem, andcarry out the integration from x = 1 to x = 2, using h = 1/64 and h = 1/128.8.3-2 For the equationy´ = 2yy(0) = 1obtain the exact solution of the difference equation obtained from Euler’s method. Estimate avalue of h small enough to guarantee four-place accuracy in the solution over the interval[0, 1]. Carry out the solution with an appropriate value of h for 10 steps.8.3-3 From the Taylor series for y(x), find y (0.1) correct to six decimal places if y(x) satisfiesy´ = xy + ly (0) = 18.3-4 Prove that, for the function f(x,y) = 1 - y/x of Example 8.1, y´´ = (1 - 2 y ´ )/x,y(k) = -ky(k-1)/x, k = 3, 4 .
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
