Re: Rounding errors using doubles?

["Ross Smith" <ross.s@ihug.co.nz>]

From: Sam Lantinga <slouken@devolution.com>
>
>I have to thank everyone here on this list for your responsiveness.
>The problem was caused by unexpected values in code using doubles
>that was ported from Windows.  Apparently, by default, Windows uses 
>53 bits of precision for it's double operations.  Linux uses 64 bits.
>
>Since the code relies on the exact behavior of doubles quite extensively,
>we are setting the precision of the FPU to 53 bits using the fldcw
>instruction. :)
>
>BTW, does anyone know what precision is used on the PPC and Alpha
>architectures?

It's not really a difference between operating systems or (to a first
approximation) CPUs.

There's a standard called IEC 559 (formerly IEEE 754), which specifies
three standard floating-point arithmetic formats. Two of them are
32-bit and 64-bit formats (with 24-bit and 53-bit precision,
respectively); modern C/C++ compilers almost universally equate these
to float and double. The third is an 80-bit (64-bit precision) format
used for intermediate results in multi-step calculations.

The Intel X86 processors have can perform arithmetic in all three
modes, although it always works internally in 80-bit mode and
inserts automatic conversions for the other two. I don't know
anything about PPCs or Alphas, but the IEC standard is pretty much
universal now, and I'd be mildly amazed if either of them differed
in any important ways.

Exactly how floating-point arithmetic is done is a function of the
compiler, not the operating system. You didn't say which compiler you
were using on Windows, but Microsoft Visual C++ is the most likely.
Both MSVC and EGCS perform intermediate calculations in 80-bit internal
registers wherever possible, but allow values to spill into 64-bit
memory when the compiler can't manage to fit the entire calculation
into registers (a common problem on the register-poor X86
architecture). (The details of exactly how EGCS should handle this
were the subject of heated debate on this list not so long ago.)

It looks like what happened was that, at some critical point in your
program, MSVC allowed a value to spill to 64 bits (thus truncating it
to 53-bit precision) while EGCS was able to keep it in an 80-bit
register (retaining 64-bit precision). (Incidentally, because people
tend to attach far too much importance to this sort of observation, I
hasten to add that this implies nothing about the relative performance
of the two compilers in general; perhaps they were using different
optimisation settings, or perhaps they simple made different decisions
about which intermediate results should be kept in registers and which
sacrificed.)

Both compilers have options to force pure 64-bit arithmetic (53-bit
precision) throughout (at some cost in speed): -ffloat-store on EGCS,
/Op on MSVC.

--
Ross Smith ................................... mailto:ross.s@ihug.co.nz
.............. The Internet Group, Auckland, New Zealand ..............
         "The award for the Most Effective Promotion of Linux
         goes to Microsoft."             -- Nicholas Petreley