floating point inconsistency

Andrew Haley aph@redhat.com
Tue Feb 16 17:20:00 GMT 2010


On 02/16/2010 02:39 PM, Christoph Groth wrote:
> Andrew Haley <aph@redhat.com> writes:
> 
>> I've just checked, and the 64-bit glibc version of sincos does indeed
>> do
>>
>> ENTRY (BP_SYM (__sincos))
>> 	ENTER
>>
>> 	movsd	%xmm0, -8(%rsp)
>> 	fldl	-8(%rsp)
>> 	fsincos
>>         ...
> 
> Thank you all very much for clarifying this issue.  I am very much
> surprised that x87 instructions behave differently on AMD and Intel
> processors but this seems indeed to be the case.

I don't know why you're surprised; I'm not.  ;-)

AMD and Intel processors use different algorithms for transcendental
functions, and different Intel processors may not return the same
results.  The usual guarantee is not correctly rounded results, but to
specify some accuracy, often 1ulp (Unit in the Last Place).

Correctly rounded sines and cosines require extended intermediate
precision, and in a few cases enormously extended intermediate
precision.  So, accurate routines can be slow.

There is active research in this area, in particular
Jean-Michel Muller's group at Lyon
(http://perso.ens-lyon.fr/jean-michel.muller/publications.html)
They have been working on making correctly-rounded algorithms that are
fast enough to be used everywhere, so in the future correct rounding
may be the norm.

> I found the following article (especially section 4.1) helpful:
> 
> D. Monniaux, "The Pitfalls of Verifying Floating-Point Computations",
> http://dx.doi.org/10.1145/1353445%2E1353446

cf. http://docs.sun.com/source/806-3568/ncg_goldberg.html

BTW, if your Monte Carlo simulation is affected by such tiny
differences in transcendentals, then your simulation really is broken!

Andrew.



More information about the Gcc-help mailing list