floating point question
On Fri, 31 Jan 2003, Thomas Lumley wrote:
On 31 Jan 2003, Peter Dalgaard BSA wrote:
Deepayan Sarkar <deepayan at stat.wisc.edu> writes:
Might have something to do with .Machine$double.eps on the respective
machines.
From help(.Machine),
double.eps: the smallest positive floating-point number `x' such that
`1 + x != 1'. It equals `base^ulp.digits' if either `base'
is 2 or `rounding' is 0; otherwise, it is `(base^ulp.digits)
/ 2'.
No, it's trickier than that. I think it's due to the guard digits that intel FPUs use. Both intel and Sun have 53bit mantissas in double pr. but intel stores intermediate results to 64 bit precision, so you get some double rounding effects. Essentially, during alignment, before adding to 1, 2^-53(1+2^-52) is getting rounded up to binary
(Sparc FPUs have some guard digits, as far as I recall, and the fine details do vary by hardware, or at least they used to and we still have some of the machines that varied running.)
There's some discussion of double rounding in a Appendix to David Goldberg's "What every computer scientist should know about floating point arithmetic", widely available on the Web. However, as Brian Ripley points out, we have no control and there isn't any point worrying about it. Incidentally, this may well be a reason for the irritating habit of DLLs in setting the floating point precision to 53 bits (well, the irritating bit is not changing it back). This makes the results of computations more nearly independent of where the numbers end up being stored.
I used to think that! But when it became possible to build R for Windows not chopping in-FPU calculations to 53 bits, the results were more consistent, and we changed over. At that point most of the inconsistencies were in printing, hence my remark earlier. The print routines in Windows `libc' (really msvcrt.dll) does lose a digit or two quite often, quite possibly to avoid giving meaningless ones (which glibc is happy to if asked).
Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595