Rounding errors using doubles?

[Sam Lantinga]

> >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.

> Take -ffloat-store off that list: it doesn't do the above.  It affects
> only explicit *stores* to *variables*, not intermediate coputations.  egcs
> provides no way to force 64-bit FP arithmetic.  In fact, I'd be surprised
> if -ffloat-store "improved" the behavior of the (buggy) program originally
> submitted in this thread.

Interestingly enough, the original project file (MS VC++ 5.0) did have the
/Op option, and -ffloat-store didn't force 64-bit FP arithmetic.  The
following code does force 64-bit arithmetic on Pentium style CPUs:

// Initialize the FPU at 53-bit precision to mimic Windows /Op behavior
// This code doesn't parse with the G++ parser because of the "::", but
// it works fine as a C file.
void Setup_FPU(void)
{
#ifdef i386
    unsigned short flags;

    asm volatile ("fnstcw %0":"=m" (flags));
    flags &= ~0x0300;       /* Clear precision flags */
    flags |=  0x0200;       /* Set 53-bit precision */
    asm volatile ("fldcw %0"::"m" (flags));
#endif
}

If this doesn't work on other x86 CPUs, I'd like to know.
Linux saves and restores the FPU flags on context switch, so this is
safe to call in an individual program.

BTW, I am on the list. :)

	-Sam Lantinga				(slouken@devolution.com)

Lead Programmer, Loki Entertainment Software
--
Author of Simple DirectMedia Layer -
	http://www.devolution.com/~slouken/SDL/
--