IEEE floating point support for guile?

[Chris Hanson]

   Date: Fri, 3 Nov 2000 16:17:11 -0600
   From: "John W. Eaton" <address@hidden>

   On  3-Nov-2000, Jim Blandy <address@hidden> wrote:

   Things like (/ 1 0) or (* 1e200 1e200) should produce Inf, (/ 0 0)
   should produce NaNs, these sorts of operations should be (optionally)
   raise exceptions, we should be able to test for Inf (isinf x) and NaN
   (isnan x), etc. (though maybe the spelling of isinf and isnan should
   be inf? and nan?).  Control over what exceptions are raised is system
   dependent, but that's what configure is for.  Octave has some tests
   for simple things like isinf and isnan, and the GSL provides some code
   for controlling exceptions.

Actually, I think you have it backwards.  These operations should
produce exceptions by default, and optionally not do so.

The default should be safe in that either the computation proceeds
correctly, or it signals an exception.  This is the right behavior for
naive users, because it guarantees that any problems such as overflow
or underflow are noticed and not silently ignored.

Sophisticated users can turn this off to take advantage of the
extended representation.  Such users presumably understand how to deal
with this added complexity.

I have never understood why hardware manufacturers insist on disabling
traps by default.  MIT Scheme enables traps by default so that you
will get an error if something like this happens.

A related problem on the x86 architecture, and perhaps on others, is
that the 80-bit extended representation is enabled by default.  This
is incorrect, because it means that the results of a computation
depend on details of the compiled code.  For example, a computation in
which a floating-point register is spilled to memory can get a
different answer than an identical computation in which the register
is not spilled.  It's far better to switch the hardware into the mode
where all results are stored internally in the standard 64-bit format.
That makes the computations completely repeatable.  It also makes them
portable to other architectures.

FYI, here is some gcc code that initializes the x86 architecture in
what I consider a sane manner:

/*

Code to initialize the ix87 FP coprocessor control word.
This code must be run once before starting a computation.

Bit(s)  Description
------  -----------
0       invalid operation (FP stack overflow or IEEE invalid arithmetic op)
1       denormalized operand
2       zero divide
3       overflow
4       underflow
5       precision (indicates that precision was lost; happens frequently)

        The first 6 bits control how the chip responds to various
        exceptions.  If a given mask bit is 0, then that exception
        will generate a processor trap.  If the mask bit is 1, then
        that exception will not trap but is handled by a "default
        action", usually substituting an infinity or NaN.

        Default is all masks set to 1.

8/9     precision control

        00 IEEE single precision
        01 (reserved)
        10 IEEE double precision
        11 non-IEEE extended precision

        Default is non-IEEE extended precision.

10/11   rounding control

        00 round to nearest or even
        01 round toward negative infinity
        10 round toward positive infinity
        11 truncate toward zero

        Default is round to nearest or even.

This code (0x0220) sets these bits as follows:

1. Precision mask 1, all others 0.
2. Precision control: IEEE double precision.
3. Rounding control: round to nearest or even.

*/

#ifdef __GNUC__
#if #cpu (i386)

void
initialize_387_to_ieee (void)
{
  unsigned short control_word;
  asm ("fclex" : : );
  asm ("fnstcw %0" : "=m" (control_word) : );
  asm ("andw %2,%0" : "=m" (control_word) : "0" (control_word), "n" (0xf0e0));
  asm ("orw %2,%0" : "=m" (control_word) : "0" (control_word), "n" (0x0220));
  asm ("fldcw %0" : : "m" (control_word));
}

#endif /* #cpu (i386) */
#endif /* __GNUC__ */