Press, Teukolsly, Vetterling, Flannery - Numerical Recipes in C (523184), страница 13

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 13)

(The square root comes from arandom-walk.) However, this estimate can be very badly off the mark for two reasons:(i) It very frequently happens that the regularities of your calculation, or thepeculiarities of your computer, cause the roundoff errors to accumulate preferentiallyin one direction. In this case the total will be of order N m .(ii) Some especially unfavorable occurrences can vastly increase the roundofferror of single operations. Generally these can be traced to the subtraction of twovery nearly equal numbers, giving a result whose only significant bits are those(few) low-order ones in which the operands differed.

You might think that such a“coincidental” subtraction is unlikely to occur. Not always so. Some mathematicalexpressions magnify its probability of occurrence tremendously. For example, in thefamiliar formula for the solution of a quadratic equation,x=−b +√b2 − 4ac2a(1.3.2)the addition becomes delicate and roundoff-prone whenever ac b2 . (In §5.6 wewill learn how to avoid the problem in this particular case.)30Chapter 1.PreliminariesRoundoff error is a characteristic of computer hardware. There is another,different, kind of error that is a characteristic of the program or algorithm used,independent of the hardware on which the program is executed.

Many numericalalgorithms compute “discrete” approximations to some desired “continuous” quantity. For example, an integral is evaluated numerically by computing a functionat a discrete set of points, rather than at “every” point. Or, a function may beevaluated by summing a finite number of leading terms in its infinite series, ratherthan all infinity terms. In cases like this, there is an adjustable parameter, e.g., thenumber of points or of terms, such that the “true” answer is obtained only whenthat parameter goes to infinity. Any practical calculation is done with a finite, butsufficiently large, choice of that parameter.The discrepancy between the true answer and the answer obtained in a practicalcalculation is called the truncation error.

Truncation error would persist even on ahypothetical, “perfect” computer that had an infinitely accurate representation and noroundoff error. As a general rule there is not much that a programmer can do aboutroundoff error, other than to choose algorithms that do not magnify it unnecessarily(see discussion of “stability” below). Truncation error, on the other hand, is entirelyunder the programmer’s control. In fact, it is only a slight exaggeration to saythat clever minimization of truncation error is practically the entire content of thefield of numerical analysis!Most of the time, truncation error and roundoff error do not strongly interactwith one another.

A calculation can be imagined as having, first, the truncation errorthat it would have if run on an infinite-precision computer, “plus” the roundoff errorassociated with the number of operations performed.Sometimes, however, an otherwise attractive method can be unstable. Thismeans that any roundoff error that becomes “mixed into” the calculation at an earlystage is successively magnified until it comes to swamp the true answer. An unstablemethod would be useful on a hypothetical, perfect computer; but in this imperfectworld it is necessary for us to require that algorithms be stable — or if unstablethat we use them with great caution.Here is a simple, if somewhat artificial, example of an unstable algorithm:Suppose that it is desired to calculate all integer powers of the so-called “GoldenMean,” the number given by√5−1≈ 0.61803398(1.3.3)φ≡2It turns out (you can easily verify) that the powers φn satisfy a simple recursionrelation,φn+1 = φn−1 − φn(1.3.4)Thus, knowing the first two values φ0 = 1 and φ1 = 0.61803398, we cansuccessively apply (1.3.4) performing only a single subtraction, rather than a slowermultiplication by φ, at each stage.Unfortunately, the recurrence (1.3.4) also has another solution, namely the value√− 12 ( 5 + 1).

Since the recurrence is linear, and since this undesired solution hasmagnitude greater than unity, any small admixture of it introduced by roundoff errorswill grow exponentially. On a typical machine with 32-bit wordlength, (1.3.4) starts1.3 Error, Accuracy, and Stability31to give completely wrong answers by about n = 16, at which point φn is down to only10−4 . The recurrence (1.3.4) is unstable, and cannot be used for the purpose stated.We will encounter the question of stability in many more sophisticated guises,later in this book.CITED REFERENCES AND FURTHER READING:Stoer, J., and Bulirsch, R. 1980, Introduction to Numerical Analysis (New York: Springer-Verlag),Chapter 1.Kahaner, D., Moler, C., and Nash, S.

1989, Numerical Methods and Software (Englewood Cliffs,NJ: Prentice Hall), Chapter 2.Johnson, L.W., and Riess, R.D. 1982, Numerical Analysis, 2nd ed. (Reading, MA: AddisonWesley), §1.3.Wilkinson, J.H. 1964, Rounding Errors in Algebraic Processes (Englewood Cliffs, NJ: PrenticeHall).Chapter 2.

Solution of LinearAlgebraic Equations2.0 IntroductionA set of linear algebraic equations looks like this:a11 x1 + a12 x2 + a13 x3 + · · · + a1N xN = b1a21 x1 + a22 x2 + a23 x3 + · · · + a2N xN = b2a31 x1 + a32 x2 + a33 x3 + · · · + a3N xN = b3···(2.0.1)···a M 1 x 1 + a M 2 x 2 + aM 3 x 3 + · · · + aM N x N = b MHere the N unknowns xj , j = 1, 2, . . . , N are related by M equations. Thecoefficients aij with i = 1, 2, .

. . , M and j = 1, 2, . . ., N are known numbers, asare the right-hand side quantities bi , i = 1, 2, . . . , M .Nonsingular versus Singular Sets of EquationsIf N = M then there are as many equations as unknowns, and there is a goodchance of solving for a unique solution set of xj ’s. Analytically, there can fail tobe a unique solution if one or more of the M equations is a linear combination ofthe others, a condition called row degeneracy, or if all equations contain certainvariables only in exactly the same linear combination, called column degeneracy.(For square matrices, a row degeneracy implies a column degeneracy, and viceversa.) A set of equations that is degenerate is called singular.

We will considersingular matrices in some detail in §2.6.Numerically, at least two additional things can go wrong:• While not exact linear combinations of each other, some of the equationsmay be so close to linearly dependent that roundoff errors in the machinerender them linearly dependent at some stage in the solution process. Inthis case your numerical procedure will fail, and it can tell you that ithas failed.32332.0 Introduction• Accumulated roundoff errors in the solution process can swamp the truesolution. This problem particularly emerges if N is too large. Thenumerical procedure does not fail algorithmically.

However, it returns aset of x’s that are wrong, as can be discovered by direct substitution backinto the original equations. The closer a set of equations is to being singular,the more likely this is to happen, since increasingly close cancellationswill occur during the solution. In fact, the preceding item can be viewedas the special case where the loss of significance is unfortunately total.Much of the sophistication of complicated “linear equation-solving packages”is devoted to the detection and/or correction of these two pathologies. As youwork with large linear sets of equations, you will develop a feeling for when suchsophistication is needed.

It is difficult to give any firm guidelines, since there is nosuch thing as a “typical” linear problem. But here is a rough idea: Linear sets withN as large as 20 or 50 can be routinely solved in single precision (32 bit floatingrepresentations) without resorting to sophisticated methods, if the equations are notclose to singular. With double precision (60 or 64 bits), this number can readilybe extended to N as large as several hundred, after which point the limiting factoris generally machine time, not accuracy.Even larger linear sets, N in the thousands or greater, can be solved when thecoefficients are sparse (that is, mostly zero), by methods that take advantage of thesparseness.

We discuss this further in §2.7.At the other end of the spectrum, one seems just as often to encounter linearproblems which, by their underlying nature, are close to singular. In this case, youmight need to resort to sophisticated methods even for the case of N = 10 (thoughrarely for N = 5). Singular value decomposition (§2.6) is a technique that cansometimes turn singular problems into nonsingular ones, in which case additionalsophistication becomes unnecessary.MatricesEquation (2.0.1) can be written in matrix form asA·x=b(2.0.2)Here the raised dot denotes matrix multiplication, A is the matrix of coefficients, andb is the right-hand side written as a column vector,a11 a21A=aM 1a12a22···aM 2......a1Na2N . .

Характеристики

Список файлов книги