c5-5 (779489), страница 2
Текст из файла (страница 2)
For example,the relation (5.5.2) becomes5.5 Recurrence Relations and Clenshaw’s Recurrence Formula1811 = J0 (x) + 2J2 (x) + 2J4 (x) + 2J6 (x) + · · ·(5.5.16)Incidentally, there is an important relation between three-term recurrencerelations and continued fractions. Rewrite the recurrence relation (5.5.9) asynyn−1=−bnan + yn+1 /yn(5.5.17)Iterating this equation, starting with n, givesbnbn+1yn···=−yn−1an − an+1 −(5.5.18)Pincherle’s Theorem [2] tells us that (5.5.18) converges if and only if (5.5.9) has aminimal solution fn , in which case it converges to fn /fn−1 .
This result, usually forthe case n = 1 and combined with some way to determine f0 , underlies many of thepractical methods for computing special functions that we give in the next chapter.Clenshaw’s Recurrence FormulaClenshaw’s recurrence formula [5] is an elegant and efficient way to evaluate asum of coefficients times functions that obey a recurrence formula, e.g.,f(θ) =NXck cos kθorf(x) =k=0NXck Pk (x)k=0Here is how it works: Suppose that the desired sum isf(x) =NXck Fk (x)(5.5.19)k=0and that Fk obeys the recurrence relationFn+1 (x) = α(n, x)Fn (x) + β(n, x)Fn−1 (x)(5.5.20)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).and fn is again the minimal solution to (5.5.9).
Cases other than those in thesetwo theorems are inconclusive for the existence of minimal solutions. (For moreon the stability of recurrences, see [3].)How do you proceed if the solution that you desire is the minimal solution?The answer lies in that old aphorism, that every cloud has a silver lining: If arecurrence relation is catastrophically unstable in one direction, then that (undesired)solution will decrease very rapidly in the reverse direction.
This means that youcan start with any seed values for the consecutive fj and fj+1 and (when you havegone enough steps in the stable direction) you will converge to the sequence offunctions that you want, times an unknown normalization factor. If there is someother way to normalize the sequence (e.g., by a formula for the sum of the fn ’s),then this can be a practical means of function evaluation. The method is calledMiller’s algorithm. An example often given [1,4] uses equation (5.5.2) in just thisway, along with the normalization formula182Chapter 5.Evaluation of Functionsfor some functions α(n, x) and β(n, x).
Now define the quantities yk (k =N, N − 1, . . . , 1) by the following recurrence:yN+2 = yN+1 = 0yk = α(k, x)yk+1 + β(k + 1, x)yk+2 + ck(k = N, N − 1, . . . , 1)(5.5.21)f(x) = · · ·+ [y8 − α(8, x)y9 − β(9, x)y10 ]F8 (x)+ [y7 − α(7, x)y8 − β(8, x)y9 ]F7 (x)+ [y6 − α(6, x)y7 − β(7, x)y8 ]F6 (x)+ [y5 − α(5, x)y6 − β(6, x)y7 ]F5 (x)(5.5.22)+···+ [y2 − α(2, x)y3 − β(3, x)y4 ]F2 (x)+ [y1 − α(1, x)y2 − β(2, x)y3 ]F1 (x)+ [c0 + β(1, x)y2 − β(1, x)y2 ]F0 (x)Notice that we have added and subtracted β(1, x)y2 in the last line.
If you examinethe terms containing a factor of y8 in (5.5.22), you will find that they sum to zero asa consequence of the recurrence relation (5.5.20); similarly all the other yk ’s downthrough y2 . The only surviving terms in (5.5.22) aref(x) = β(1, x)F0 (x)y2 + F1 (x)y1 + F0 (x)c0(5.5.23)Equations (5.5.21) and (5.5.23) are Clenshaw’s recurrence formula for doing the sum(5.5.19): You make one pass down through the yk ’s using (5.5.21); when you havereached y2 and y1 you apply (5.5.23) to get the desired answer.Clenshaw’s recurrence as written above incorporates the coefficients ck in adownward order, with k decreasing.
At each stage, the effect of all previous ck ’sis “remembered” as two coefficients which multiply the functions Fk+1 and Fk(ultimately F0 and F1 ). If the functions Fk are small when k is large, and if thecoefficients ck are small when k is small, then the sum can be dominated by smallFk ’s. In this case the remembered coefficients will involve a delicate cancellationand there can be a catastrophic loss of significance.
An example would be to sumthe trivial seriesJ15 (1) = 0 × J0 (1) + 0 × J1 (1) + . . . + 0 × J14 (1) + 1 × J15 (1)(5.5.24)Here J15 , which is tiny, ends up represented as a canceling linear combination ofJ0 and J1 , which are of order unity.Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).If you solve equation (5.5.21) for ck on the left, and then write out explicitly thesum (5.5.19), it will look (in part) like this:5.6 Quadratic and Cubic Equations183The solution in such cases is to use an alternative Clenshaw recurrence thatincorporates ck ’s in an upward direction.
The relevant equations aref(x) = cN FN (x) − β(N, x)FN−1 (x)yN−1 − FN (x)yN−2(5.5.25)(5.5.26)(5.5.27)The rare case where equations (5.5.25)–(5.5.27) should be used instead ofequations (5.5.21) and (5.5.23) can be detected automatically by testing whetherthe operands in the first sum in (5.5.23) are opposite in sign and nearly equal inmagnitude. Other than in this special case, Clenshaw’s recurrence is always stable,independent of whether the recurrence for the functions Fk is stable in the upwardor downward direction.CITED REFERENCES AND FURTHER READING:Abramowitz, M., and Stegun, I.A.
1964, Handbook of Mathematical Functions, Applied Mathematics Series, Volume 55 (Washington: National Bureau of Standards; reprinted 1968 byDover Publications, New York), pp. xiii, 697. [1]Gautschi, W. 1967, SIAM Review, vol. 9, pp. 24–82. [2]Lakshmikantham, V., and Trigiante, D. 1988, Theory of Difference Equations: Numerical Methodsand Applications (San Diego: Academic Press). [3]Acton, F.S.
1970, Numerical Methods That Work; 1990, corrected edition (Washington: Mathematical Association of America), pp. 20ff. [4]Clenshaw, C.W. 1962, Mathematical Tables, vol. 5, National Physical Laboratory (London: H.M.Stationery Office). [5]Dahlquist, G., and Bjorck, A. 1974, Numerical Methods (Englewood Cliffs, NJ: Prentice-Hall),§4.4.3, p.
111.Goodwin, E.T. (ed.) 1961, Modern Computing Methods, 2nd ed. (New York: Philosophical Library), p. 76.5.6 Quadratic and Cubic EquationsThe roots of simple algebraic equations can be viewed as being functions of theequations’ coefficients. We are taught these functions in elementary algebra. Yet,surprisingly many people don’t know the right way to solve a quadratic equationwith two real roots, or to obtain the roots of a cubic equation.There are two ways to write the solution of the quadratic equationax2 + bx + c = 0(5.6.1)with real coefficients a, b, c, namelyx=−b ±√b2 − 4ac2a(5.6.2)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).y−2 = y−1 = 01[yk−2 − α(k, x)yk−1 − ck ],yk =β(k + 1, x)(k = 0, 1, . .
. , N − 1).
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
