float compares on x86 (i.e. 32 bit) bug?

[Andrew Haley]

On 04/13/2010 01:22 PM, Matthias Kretz wrote:
> Hi,
> 
> I was wondering why the following code failed:
> 
> int main() {
>   const float_v x = ...;
>   const float_v y = ...;
>   float_v r = Vc::sqrt(x * x + y * y);
> 
>   for (int i = 0; i < Count; ++i) {
>     const float a = r[0];
>     const float b = std::sqrt(x[0] * x[0] + y[0] * y[0]);
>     if (a != b) {
>       std::cout << a << " != " << b << ": " << a - b << std::endl;
>       return 1;
>     }
>   }
>   return 0;
> }
> 
> Output is:
> 1.41421 != 1.41421: 0
> 
> float_v is a class from my Vc library (http://gitorious.org/vc/) that wraps 
> __m128. Vc::sqrt thus uses the SSE sqrt with single precision. std::sqrt OTOH 
> is implemented as fsqrt and the comparison "a != b" is implemented with 
> fucomip using the result from fsqrt directly and the result from Vc::sqrt 
> loaded from memory. After fucomip the following flags are set:
> (gdb) p $eflags
> $6 = [ CF IF ]
> 
> Obviously I expected ZF to be set. Here's the problem:
> st6 1.41421353816986083984375            (raw 0x3fffb504f30000000000)
> st7 1.4142135623730950487637880730318329 (raw 0x3fffb504f333f9de6484)
> 
> Looking at the binary representation the single precision values are equal. 
> Just st7 (the result from fsqrt) is not stored as single.
> 
> So now I understand why it failed. Am I right in my understanding that this is 
> a bug? According to the C++ code the comparison should not have returned false 
> because std::sqrt returns a single precision value.

Isn't this http://gcc.gnu.org/bugzilla/show_bug.cgi?id=323 ?

We won't know until you give us a real test case, one that we can compile
and run.

Andrew.