c5-2 (779486)
Текст из файла
5.2 Evaluation of Continued Fractions169CITED REFERENCES AND FURTHER READING:Goodwin, E.T. (ed.) 1961, Modern Computing Methods, 2nd ed. (New York: Philosophical Library), Chapter 13 [van Wijngaarden’s transformations]. [1]Dahlquist, G., and Bjorck, A. 1974, Numerical Methods (Englewood Cliffs, NJ: Prentice-Hall),Chapter 3.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), §3.6.Mathews, J., and Walker, R.L. 1970, Mathematical Methods of Physics, 2nd ed. (Reading, MA:W.A. Benjamin/Addison-Wesley), §2.3.
[2]5.2 Evaluation of Continued FractionsContinued fractions are often powerful ways of evaluating functions that occurin scientific applications. A continued fraction looks like this:f(x) = b0 +a1b1 +a2b2 +b3 +(5.2.1)a3a4a5b4 +b5 +···Printers prefer to write this asf(x) = b0 +a1a2a3a4a5···b1 + b2 + b3 + b4 + b5 +(5.2.2)In either (5.2.1) or (5.2.2), the a’s and b’s can themselves be functions of x, usuallylinear or quadratic monomials at worst (i.e., constants times x or times x2 ). Forexample, the continued fraction representation of the tangent function istan x =x x2 x2 x2···1− 3− 5− 7−(5.2.3)Continued fractions frequently converge much more rapidly than power seriesexpansions, and in a much larger domain in the complex plane (not necessarilyincluding the domain of convergence of the series, however).
Sometimes thecontinued fraction converges best where the series does worst, although this is notSample 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).into equation (5.1.11), and then setting z = 1.Sometimes you will want to compute a function from a series representationeven when the computation is not efficient.
For example, you may be using the valuesobtained to fit the function to an approximating form that you will use subsequently(cf. §5.8). If you are summing very large numbers of slowly convergent terms, payattention to roundoff errors! In floating-point representation it is more accurate tosum a list of numbers in the order starting with the smallest one, rather than startingwith the largest one. It is even better to group terms pairwise, then in pairs of pairs,etc., so that all additions involve operands of comparable magnitude.170Chapter 5.Evaluation of Functionsfn =AnBn(5.2.4)where An and Bn are given by the following recurrence:A−1 ≡ 1A0 ≡ b0Aj = bj Aj−1 + aj Aj−2B−1 ≡ 0B0 ≡ 1Bj = bj Bj−1 + aj Bj−2j = 1, 2, .
. . , n(5.2.5)This method was invented by J. Wallis in 1655 (!), and is discussed in his ArithmeticaInfinitorum [4]. You can easily prove it by induction.In practice, this algorithm has some unattractive features: The recurrence (5.2.5)frequently generates very large or very small values for the partial numerators anddenominators Aj and Bj . There is thus the danger of overflow or underflow of thefloating-point representation. However, the recurrence (5.2.5) is linear in the A’s andB’s.
At any point you can rescale the currently saved two levels of the recurrence,e.g., divide Aj , Bj , Aj−1, and Bj−1 all by Bj . This incidentally makes Aj = fjand is convenient for testing whether you have gone far enough: See if fj and fj−1from the last iteration are as close as you would like them to be. (If Bj happens tobe zero, which can happen, just skip the renormalization for this cycle. A fancierlevel of optimization is to renormalize only when an overflow is imminent, savingthe unnecessary divides.
All this complicates the program logic.)Two newer algorithms have been proposed for evaluating continued fractions.Steed’s method does not use Aj and Bj explicitly, but only the ratio Dj = Bj−1 /Bj .One calculates Dj and ∆fj = fj − fj−1 recursively usingDj = 1/(bj + aj Dj−1 )(5.2.6)∆fj = (bj Dj − 1)∆fj−1(5.2.7)Steed’s method (see, e.g., [5]) avoids the need for rescaling of intermediate results.However, for certain continued fractions you can occasionally run into a situationSample 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).a general rule. Blanch [1] gives a good review of the most useful convergence testsfor continued fractions.There are standard techniques, including the important quotient-difference algorithm, for going back and forth between continued fraction approximations, powerseries approximations, and rational function approximations.
Consult Acton [2] foran introduction to this subject, and Fike [3] for further details and references.How do you tell how far to go when evaluating a continued fraction? Unlikea series, you can’t just evaluate equation (5.2.1) from left to right, stopping whenthe change is small. Written in the form of (5.2.1), the only way to evaluate thecontinued fraction is from right to left, first (blindly!) guessing how far out tostart. This is not the right way.The right way is to use a result that relates continued fractions to rationalapproximations, and that gives a means of evaluating (5.2.1) or (5.2.2) from leftto right. Let fn denote the result of evaluating (5.2.2) with coefficients throughan and bn .
Then5.2 Evaluation of Continued Fractions171Cj = Aj /Aj−1 ,Dj = Bj−1 /Bj(5.2.8)and calculating fj byfj = fj−1 Cj Dj(5.2.9)From equation (5.2.5), one easily shows that the ratios satisfy the recurrence relationsDj = 1/(bj + aj Dj−1 ),Cj = bj + aj /Cj−1(5.2.10)In this algorithm there is the danger that the denominator in the expression for Dj ,or the quantity Cj itself, might approach zero.
Either of these conditions invalidates(5.2.10). However, Thompson and Barnett [5] show how to modify Lentz’s algorithmto fix this: Just shift the offending term by a small amount, e.g., 10−30 . If youwork through a cycle of the algorithm with this prescription, you will see that fj+1is accurately calculated.In detail, the modified Lentz’s algorithm is this:• Set f0 = b0 ; if b0 = 0 set f0 = tiny.• Set C0 = f0 .• Set D0 = 0.• For j = 1, 2, . . .Set Dj = bj + aj Dj−1 .If Dj = 0, set Dj = tiny.Set Cj = bj + aj /Cj−1 .If Cj = 0 set Cj = tiny.Set Dj = 1/Dj .Set ∆j = Cj Dj .Set fj = fj−1 ∆j .If |∆j − 1| < eps then exit.Here eps is your floating-point precision, say 10−7 or 10−15 .
The parameter tinyshould be less than typical values of eps|bj |, say 10−30 .The above algorithm assumes that you can terminate the evaluation of thecontinued fraction when |fj − fj−1 | is sufficiently small. This is usually the case,but by no means guaranteed.
Jones [7] gives a list of theorems that can be used tojustify this termination criterion for various kinds of continued fractions.There is at present no rigorous analysis of error propagation in Lentz’s algorithm.However, empirical tests suggest that it is at least as good as other methods.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).where the denominator in (5.2.6) approaches zero, so that Dj and ∆fj are verylarge.
The next ∆fj+1 will typically cancel this large change, but with loss ofaccuracy in the numerical running sum of the fj ’s. It is awkward to program aroundthis, so Steed’s method can be recommended only for cases where you know inadvance that no denominator can vanish. We will use it for a special purpose inthe routine bessik (§6.7).The best general method for evaluating continued fractions seems to be themodified Lentz’s method [6]. The need for rescaling intermediate results is avoidedby using both the ratios172Chapter 5.Evaluation of FunctionsManipulating Continued FractionsSeveral important properties of continued fractions can be used to rewrite themin forms that can speed up numerical computation.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.
















