Higham - Accuracy and Stability of Numerical Algorithms (523152), страница 16
Текст из файла (страница 16)
This approach is superficially similar to the earlierCESTAC (permutation-perturbation) method of La Porte and Vignes [153,1986], [682, 1974], [1054, 1986], but differs from it in several respects. CESTAC deals with the arithmetic reliability of algorithms, whereas PRECISEis designed as a tool to explore the robustness of numerical algorithms as afunction of parameters such as mesh size, time step, and nonnormality.Several authors have investigated the distribution of rounding errors underthe assumption that the mantissas of the operands are from a logarithmicdistribution, and for different modes of rounding; see Barlow and Bareiss [62,1985] and the references therein.Other work concerned with statistical modelling of rounding errors includes that of Tienari [1002, 1970] and Linnainmaa [704, 1975].2.7.
Alternative Number SystemsThe floating point format is not the only means for representing numbers infinite precision arithmetic. Various alternatives have been proposed, thoughnone has achieved widespread use.A particularly elegant system is the “level index arithmetic” of Clenshaw,Olver, and Turner, in which a number x > 1 is represented by= l + f,where f[0,1] andorf = ln(ln(. . .
(lnx) . . .)),where the exponentiation or logarithm is performed l times (l is the “level”).If 0 < x < 1, then x is represented by the reciprocal of the representationfor l/x. An obvious feature of the level index system is that it can representmuch larger and smaller numbers than the floating point system, for similarword lengths. A price to be paid is that addition and subtraction are morecomplicated (and more costly) than in floating point arithmetic.
For veryreadable introductions to level index arithmetic see Clenshaw and Olver [213,1984] and Turner [1033, 1991], and for more details see Clenshaw, Olver, andTurner [214, 1989]. Level index arithmetic is somewhat controversial in thatthere is disagreement about its advantages and disadvantages with respect tofloating point arithmetic; see Demmel [282, 1987]. A number system involvinglevels has also been proposed by Matsui and Iri [736, 1981]; in their system,54F LOATING P OINT A RITHMETICTable 2.3.
Test arithmetics.HardwareCasio fx-140 ( 1979)Casio fx-992VB ( 1990)HP 48G (1993)Sharp EL-5020 (1994)486DX486DX486DX486DXSoftwareMATLAB 4.2 (1994)WATFOR-77c V3.0 (1988)FTN 90d (1993)MS-DOS QBasic 1.1|3 × (4/3 - 1) -1|a1 × 10-91 × 10-131 × 10-110.0 b2.2. . . × 10-162.2.
. . × 10-162.2. . . × 10-161.1. . . × 10-19 eaIntegers in the test expression are typed as real constants 3.0, etc., for the Fortrantests.b1 × 10-9 if 4/3 is stored and recalled from memory.cWATCOM Systems Inc.dSalford Software/Numerical Algorithms Group, Version 1.2.e2.2. . . . × 10-16 if 4/3 is stored and recalled from a variable.the number of bits allocated to the mantissa and exponent is allowed to vary(within a fixed word size).Other number systems include those of Swartzlander and Alexopolous [981,1975], Matula and Kornerup [741, 1985], and Hamada [495, 1987].
For summaries of alternatives to floating point arithmetic see the section “Alternativesto Floating-Point Some Candidates” in [214, 198 9], and Knuth [668, 198 1,Chap. 4].2.8. Accuracy TestsHow can you test the accuracy of the floating point arithmetic on a computeror pocket calculator? There is no easy way, though a few software packages areavailable to help with the tasks in specific programming languages (see §25.6).There are, however, a few quick and easy tests that may reveal weaknesses.The following list is far from comprehensive and good performance on thesetests does not imply that an arithmetic is correct.
Results from the tests aregiven in Tables 2.4-2.5 for the selected floating point arithmetics described inTable 2.3. Double precision was used for the compiled languages. The lastcolumn of Table 2.3 gives an estimate of the unit roundoff (see Problem 2.14).The estimate produced by QBasic indicates that the compiler used extendedprecision in evaluating the estimate.1. (Cody [221, 1982]) Evaluate sin(22) = -8.8513 0929 0403 8759 2169 ×10-3 (shown correct to 21 digits).
This is a difficult test for the range2.8 A CCURACY T ESTS55Table 2.4. Sine test.Machine |sin(22)Exact -8.8513 0929 0403 8759 × 10-3Casio fx-140 -8.8513 62 × 10 -3Casio fx-992VB -8.8513 0929 096 × 10-3HP 48G -8.8513 0929 040 × 10-3Sharp EL-5020 -8.8513 0915 4 × 10-3MATLAB 4.2 -8.8513 0929 0403 876 × 10-3WATFOR-77 -8.8513 0929 0403 880 × 10-3FTN 90 -8.8513 0929 0403 876 × 10-3QBasic -8.8513 0929 0403 876 × 10-3Table 2.5. Exponentation test. No entry for last column means same value asprevious column.MachineExactCasio fx-140Casio fx-992VBHP 48GSharp EL-5020MATLAB 4.2WATFOR-77FTN 90QBasic5.52715.52715.52715.52715.52715.52715.52715.52715.52712.5 1254787 5260 4446 × 1049477 × 10494787 526 × 10494787 526 × 10494787 3 × 10494787 5260 445 × 10494787 5260 450 × 10494787 5260 445 × 10494787 5260 444 × 1049exp(125log(2.5))5.5271 4787 5260 4446 × 10495.5271 463 × 10495.52715.52715.52715.52715.527147874796478747874787377 × 10492 × 10495260 459 × 10495260 460 × 10495260 459 × 1049reduction used in the sine evaluation (which brings the argument withint h e r a n g e [ -π/2,π/2], and which necessarily uses an approximate valueof π), since 22 is close to an integer multiple of π.2.
(Cody [221, 1982]) Evaluate 2.5125 = 5.5271 4787 5260 4445 6025 × 1049(shown correct to 21 digits). One way to evaluate z = xy is as z =exp(ylogx). But to obtain z correct to within a few ulps it is not sufficient to compute exp and log correct to within a few ulps; in otherwords, the composition of two functions evaluated to high accuracy isnot necessarily obtained to the same accuracy.
To examine this particular case, writew : = ylogx,z = exp(w).F LOATING P OINT A RITHMETIC56If w w + ∆w then zz + ∆z, where z + ∆z = exp(w+∆w) =exp(w) exp(∆w)exp(w)(l+∆w), so ∆ z / z∆ w. In other words,the relative error of z depends on the absolute error of w and hence onthe size of w. To obtain z correct to within a few ulps it is necessaryto use extra precision in calculating the logarithm and exponential [228,1980, Chap. 7].3.
(Karpinski [645, 1985]) A simple test for the presence of a guard digiton a pocket calculator is to evaluate the expressions9/27 * 3 - 1,9/27 * 3 - 0.5 - 0.5,which are given in a form that can be typed directly into most fourfunction calculators. If the results are equal then a guard digit is present.Otherwise there is probably no guard digit (we cannot be completelysure from this simple test).
To test for a guard digit on a computer itis best to run one of the diagnostic codes described in §25.5.2.9. Notes and ReferencesThe classic reference on floating point arithmetic, and on all aspects of rounding error analysis, is Wilkinson’s Rounding Errors in Algebraic Processes(REAP) [1088, 1963]. Wilkinson was uniquely qualified to write such a book,for not only was he the leading expert in rounding error analysis, but he wasone of the architects and builders of the Automatic Computing Engine (ACE)at the National Physical Laboratory [1082, 1954].
The Pilot (prototype) ACEfirst operated in May 1950, and an engineered version was later sold commercially as the DEUCE Computer by the English Electric Company. Wilkinsonand his colleagues were probably the first to write subroutines for floatingpoint arithmetic. and this enabled them to accumulate practical experienceof floating point arithmetic much earlier than anyone else [357, 1976], [1099,1980].In REAP, Wilkinson gives equal consideration to fixed point and floatingpoint arithmetic. In fixed point arithmetic, all numbers are constrained to liein a range such as [-1,1], as if the exponent were frozen in the floating pointrepresentation (2.1).
Preliminary analysis and the introduction of scale factorsduring the computation is needed to keep numbers in the permitted range.We consider only floating point arithmetic in this book. REAP, together withWilkinson’s second book, The Algebraic Eigenvalue Problem (AEP) [1089,1965], has been immensely influential in the areas of floating point arithmeticand rounding error analysis.Wilkinson’s books were preceded by the paper Error Analysis of FloatingPoint Computation [1084, 196 0], in which he presents the model (2.4) forfloating point arithmetic and applies the model to several algorithms for the2.9 NOTESANDREFERENCES57eigenvalue problem. This paper has hardly dated and is still well worth reading.Another classic book devoted entirely to floating point arithmetic is Sterbenz’s Floating-Point Computation [938, 1974]. It contains a thorough treatment of low-level details of floating point arithmetic, with particular referenceto IBM 360 and IBM 7090 machines.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
