On Mon, Feb 6, 2012 at 1:29 PM, Vincent Lefevre <vincent+...@vinc17.org> wrote:
> On 2012-02-06 12:54:09 +0100, Richard Guenther wrote:
>> Note that you are comparing a constant folded sin() result against
>> sincos() (or sin() and cos()). Use
>>
>> #include <stdio.h>
>> #include <math.h>
>>
>> int main (void)
>> {
>> double x, c, s;
>> volatile double v;
>>
>> x = 1.0e22;
>> v = x;
>> x = v;
>> s = sin (x);
>> printf ("sin(%.17g) = %.17g\n", x, s);
>>
>> v = x;
>> x = v;
>> c = cos (x);
>> s = sin (x);
>> printf ("sin(%.17g) = %.17g\n", x, s);
>>
>> return c == 0;
>> }
>>
>> to compare libm sin against sincos. Works for me, btw. Note that I remember
>> that at some point sincos() was just fsincos even on x86_64 (ugh).
>
> This is unfortunately still true under Debian/unstable, which has
> glibc 2.13, and that's why I get (with -O, so that glibc's sincos
> is used):
>
> sin(1e+22) = -0.85220084976718879
> sin(1e+22) = 0.46261304076460175
>
> What do you mean by "Works for me"? i.e. you get the same value or
> you can reproduce the bug? If the former, has this changed in later
> glibc versions or has glibc been built with some specific option on
> your machine?
glibc for openSUSE or SLE carries patches to replace all x86_64 libm
routines and does not suffer from the above problem. Note that the
issue in glibc is that sin/cos use a software implementation while sincos
uses the hardware fsincos. sinl/cosl use hardware fsin/fcos, similar to
sincosl so those should be not prone to the error (and sinf/cosf/sincosf
consistently use a software implementation).
Richard.
