Date: Fri, 28 Sep 2007 08:22:27 -0700 From: Steve Kargl <sgk@troutmask.apl.washington.edu> To: freebsd-standards@freebsd.org Subject: long double broken on i386? Message-ID: <20070928152227.GA39233@troutmask.apl.washington.edu>
next in thread | raw e-mail | index | archive | help
So, in my never ending crusade to implement missing C99 math functions, I decided to tackle sinl, cosl, and tanl. I have implementations that use argument reduction from k_rem_pio2.c and Chebyshev interpolation for the trig functions over the range [0,pi/4]. On amd64, these routines give <= 1/2 ULP for the various values I've tested in the reduced range. Now, if I test these on i386, I see between 400 and 1200 ULP. I spent a week or so trying to understand the descrepancy and "fix" the problem using a double-double precision arithemetic. I finally found the problem! /usr/include/float.h on i386 is lying about its numerical representation of long double. In particular, at least these 3 values appear to be wrong #define LDBL_MANT_DIG 64 #define LDBL_EPSILON 1.0842021724855044340E-19L #define LDBL_DIG 18 #include <stdio.h> int main (void) { int i; long double a; a = 1.L; for (i = 1; i < 64; i++) { a /= 2; if (1.L - a == 1.L) break; } printf("%d %Le\n", i, a); return 0; } amd64:sgk[206] ./z 64 1.084202e-19 i386:kargl[205] ./z 54 5.551115e-17 The above results for i386 are consistent with DBL_MANT_DIG and DBL_EPSILON. Is this intentional, and should float.h be fixed? -- Steve
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20070928152227.GA39233>