Excess precision problem on IA-64

[Eric Botcazou <ebotcazou at adacore dot com>]

Hi,

We have run into an excess precision problem on IA-64, similar to the infamous 
problem on IA-32.  A C testcase is attached.

botcazou@taff:~/EA26-004> gcc -v
Using built-in specs.
Target: ia64-sgi-linux-gnu
Configured with: /home01/botcazou/cvs/gcc/configure ia64-sgi-linux-gnu 
--prefix=/nile.build/botcazou/fsf/install_ia64 --disable-nls 
--disable-libmudflap --disable-checking --enable-languages=c,c++
Thread model: posix
gcc version 4.1.0 20051024 (experimental)
botcazou@taff:~/EA26-004> gcc -o t t.c
botcazou@taff:~/EA26-004> ./t
botcazou@taff:~/EA26-004> gcc -o t t.c -O
botcazou@taff:~/EA26-004> ./t
Aborted

The bottom line is that:

  fmpy.s f8 = f8, f15
  fma.s f8 = f11, f12, f8

doesn't always give 0.0 in f8 when f8 = f15 = f11 = -f12, because the
computations are done in "infinite precision" and only rounded at the end.


What should we do about that?  Conditionalize the combined FP operations on 
-ffast-math or something along these lines?

Thanks in advance.

-- 
Eric Botcazou

extern void abort (void);

typedef struct
{
  float X;
  float Y;
  float Z;
  float S;
} Unit_Quaternion_Type;

Unit_Quaternion_Type Mult (Unit_Quaternion_Type Left, Unit_Quaternion_Type Right)
{
  Unit_Quaternion_Type ret
    = {Left.S * Right.X + Left.X * Right.S + Left.Y * Right.Z - Left.Z * Right.Y,
       Left.S * Right.Y - Left.X * Right.Z + Left.Y * Right.S + Left.Z * Right.X,
       Left.S * Right.Z + Left.X * Right.Y - Left.Y * Right.X + Left.Z * Right.S,
       Left.S * Right.S - Left.X * Right.X - Left.Y * Right.Y - Left.Z * Right.Z};
  return ret;
}

int main(void)
{
  Unit_Quaternion_Type A = { 0.707106769, 0.0, 0.0, 0.707106769};
  Unit_Quaternion_Type B = {-0.707106769, 0.0, 0.0, 0.707106769};
  Unit_Quaternion_Type C;

  C = Mult (A, B);
  if (C.X != 0.0)
    abort ();

  return 0;
}