This is the mail archive of the
libc-alpha@sourceware.org
mailing list for the glibc project.
Fix acos (-1) in round-downwards mode on x86 (bug 14034)
- From: "Joseph S. Myers" <joseph at codesourcery dot com>
- To: libc-alpha at sourceware dot org
- Date: Sun, 29 Apr 2012 15:01:00 +0000 (UTC)
- Subject: Fix acos (-1) in round-downwards mode on x86 (bug 14034)
Bug 14034 is the problem that acos (-1) returns -pi instead of pi in
round-downwards mode, for x86 for all floating-point formats and for
x86_64 for long double (x86_64 #includes the x86 implementation of
acosl). The ISO C requirement is that the result is always in the
range [0, pi].
The problem is that the fsqrt instruction follows the IEEE 754
requirement that sqrt(-0) is -0, and a previous subtraction resulted
in -0 (x-x is -0 in round-downwards mode), so the fpatan instruction
gets a -0 value where it needs a value without the sign bit set to
produce results in the correct quadrant. I propose this patch fixing
this in the obvious way by taking the absolute values of the square
roots. I've added tests for acos in all rounding modes - and for asin
although the similar square roots there do not cause any problems.
Tested x86 and x86_64 and ulps updated accordingly.
2012-04-29 Joseph Myers <joseph@codesourcery.com>
[BZ #14034]
* sysdeps/i386/fpu/e_acos.S (__ieee754_acos): Take absolute value
of square root.
* sysdeps/i386/fpu/e_acosf.S (__ieee754_acosf): Likewise.
* sysdeps/i386/fpu/e_acosl.c (__ieee754_acosl): Likewise.
* math/libm-test.inc (acos_test_tonearest): New function.
(acos_test_towardzero): Likewise.
(acos_test_downward): Likewise.
(acos_test_upward): Likewise.
(asin_test_tonearest): Likewise.
(asin_test_towardzero): Likewise.
(asin_test_downward): Likewise.
(asin_test_upward): Likewise.
(main): Call the new functions.
* sysdeps/i386/fpu/libm-test-ulps: Update.
* sysdeps/x86_64/fpu/libm-test-ulps: Update.
diff --git a/math/libm-test.inc b/math/libm-test.inc
index 1b5d1c7..0ee1d55 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -761,6 +761,126 @@ acos_test (void)
END (acos);
}
+
+static void
+acos_test_tonearest (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(acos) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (acos_tonearest);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_TONEAREST))
+ {
+ TEST_f_f (acos, 0, M_PI_2l);
+ TEST_f_f (acos, minus_zero, M_PI_2l);
+ TEST_f_f (acos, 1, 0);
+ TEST_f_f (acos, -1, M_PIl);
+ TEST_f_f (acos, 0.5, M_PI_6l*2.0);
+ TEST_f_f (acos, -0.5, M_PI_6l*4.0);
+ }
+
+ fesetround (save_round_mode);
+
+ END (acos_tonearest);
+}
+
+
+static void
+acos_test_towardzero (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(acos) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (acos_towardzero);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_TOWARDZERO))
+ {
+ TEST_f_f (acos, 0, M_PI_2l);
+ TEST_f_f (acos, minus_zero, M_PI_2l);
+ TEST_f_f (acos, 1, 0);
+ TEST_f_f (acos, -1, M_PIl);
+ TEST_f_f (acos, 0.5, M_PI_6l*2.0);
+ TEST_f_f (acos, -0.5, M_PI_6l*4.0);
+ }
+
+ fesetround (save_round_mode);
+
+ END (acos_towardzero);
+}
+
+
+static void
+acos_test_downward (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(acos) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (acos_downward);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_DOWNWARD))
+ {
+ TEST_f_f (acos, 0, M_PI_2l);
+ TEST_f_f (acos, minus_zero, M_PI_2l);
+ TEST_f_f (acos, 1, 0);
+ TEST_f_f (acos, -1, M_PIl);
+ TEST_f_f (acos, 0.5, M_PI_6l*2.0);
+ TEST_f_f (acos, -0.5, M_PI_6l*4.0);
+ }
+
+ fesetround (save_round_mode);
+
+ END (acos_downward);
+}
+
+
+static void
+acos_test_upward (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(acos) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (acos_upward);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_UPWARD))
+ {
+ TEST_f_f (acos, 0, M_PI_2l);
+ TEST_f_f (acos, minus_zero, M_PI_2l);
+ TEST_f_f (acos, 1, 0);
+ TEST_f_f (acos, -1, M_PIl);
+ TEST_f_f (acos, 0.5, M_PI_6l*2.0);
+ TEST_f_f (acos, -0.5, M_PI_6l*4.0);
+ }
+
+ fesetround (save_round_mode);
+
+ END (acos_upward);
+}
+
static void
acosh_test (void)
{
@@ -817,6 +937,126 @@ asin_test (void)
END (asin);
}
+
+static void
+asin_test_tonearest (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(asin) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (asin_tonearest);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_TONEAREST))
+ {
+ TEST_f_f (asin, 0, 0);
+ TEST_f_f (asin, minus_zero, minus_zero);
+ TEST_f_f (asin, 0.5, M_PI_6l);
+ TEST_f_f (asin, -0.5, -M_PI_6l);
+ TEST_f_f (asin, 1.0, M_PI_2l);
+ TEST_f_f (asin, -1.0, -M_PI_2l);
+ }
+
+ fesetround (save_round_mode);
+
+ END (asin_tonearest);
+}
+
+
+static void
+asin_test_towardzero (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(asin) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (asin_towardzero);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_TOWARDZERO))
+ {
+ TEST_f_f (asin, 0, 0);
+ TEST_f_f (asin, minus_zero, minus_zero);
+ TEST_f_f (asin, 0.5, M_PI_6l);
+ TEST_f_f (asin, -0.5, -M_PI_6l);
+ TEST_f_f (asin, 1.0, M_PI_2l);
+ TEST_f_f (asin, -1.0, -M_PI_2l);
+ }
+
+ fesetround (save_round_mode);
+
+ END (asin_towardzero);
+}
+
+
+static void
+asin_test_downward (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(asin) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (asin_downward);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_DOWNWARD))
+ {
+ TEST_f_f (asin, 0, 0);
+ TEST_f_f (asin, minus_zero, minus_zero);
+ TEST_f_f (asin, 0.5, M_PI_6l);
+ TEST_f_f (asin, -0.5, -M_PI_6l);
+ TEST_f_f (asin, 1.0, M_PI_2l);
+ TEST_f_f (asin, -1.0, -M_PI_2l);
+ }
+
+ fesetround (save_round_mode);
+
+ END (asin_downward);
+}
+
+
+static void
+asin_test_upward (void)
+{
+ int save_round_mode;
+ errno = 0;
+ FUNC(asin) (0);
+ if (errno == ENOSYS)
+ /* Function not implemented. */
+ return;
+
+ START (asin_upward);
+
+ save_round_mode = fegetround ();
+
+ if (!fesetround (FE_UPWARD))
+ {
+ TEST_f_f (asin, 0, 0);
+ TEST_f_f (asin, minus_zero, minus_zero);
+ TEST_f_f (asin, 0.5, M_PI_6l);
+ TEST_f_f (asin, -0.5, -M_PI_6l);
+ TEST_f_f (asin, 1.0, M_PI_2l);
+ TEST_f_f (asin, -1.0, -M_PI_2l);
+ }
+
+ fesetround (save_round_mode);
+
+ END (asin_upward);
+}
+
static void
asinh_test (void)
{
@@ -8192,7 +8432,15 @@ main (int argc, char **argv)
/* Trigonometric functions: */
acos_test ();
+ acos_test_tonearest ();
+ acos_test_towardzero ();
+ acos_test_downward ();
+ acos_test_upward ();
asin_test ();
+ asin_test_tonearest ();
+ asin_test_towardzero ();
+ asin_test_downward ();
+ asin_test_upward ();
atan_test ();
atan2_test ();
cos_test ();
diff --git a/sysdeps/i386/fpu/e_acos.S b/sysdeps/i386/fpu/e_acos.S
index d3505ba..d10a054 100644
--- a/sysdeps/i386/fpu/e_acos.S
+++ b/sysdeps/i386/fpu/e_acos.S
@@ -15,6 +15,7 @@ ENTRY(__ieee754_acos)
fld1 /* 1 : x^2 : x */
fsubp /* 1 - x^2 : x */
fsqrt /* sqrt (1 - x^2) : x */
+ fabs
fxch %st(1) /* x : sqrt (1 - x^2) */
fpatan /* atan (sqrt(1 - x^2) / x) */
ret
diff --git a/sysdeps/i386/fpu/e_acosf.S b/sysdeps/i386/fpu/e_acosf.S
index 6a843a5..54930af 100644
--- a/sysdeps/i386/fpu/e_acosf.S
+++ b/sysdeps/i386/fpu/e_acosf.S
@@ -16,6 +16,7 @@ ENTRY(__ieee754_acosf)
fld1
fsubp /* 1 - x^2 */
fsqrt /* sqrt (1 - x^2) */
+ fabs
fxch %st(1)
fpatan
ret
diff --git a/sysdeps/i386/fpu/e_acosl.c b/sysdeps/i386/fpu/e_acosl.c
index ec516ff..d249d5a 100644
--- a/sysdeps/i386/fpu/e_acosl.c
+++ b/sysdeps/i386/fpu/e_acosl.c
@@ -18,6 +18,7 @@ __ieee754_acosl (long double x)
"fld1\n"
"fsubp\n" /* 1 - x^2 */
"fsqrt\n" /* sqrtl (1 - x^2) */
+ "fabs\n"
"fxch %%st(1)\n"
"fpatan"
: "=t" (res) : "0" (x) : "st(1)");
diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
index c3a3ce0..fdaff35 100644
--- a/sysdeps/i386/fpu/libm-test-ulps
+++ b/sysdeps/i386/fpu/libm-test-ulps
@@ -5,6 +5,77 @@ Test "acos (0.75) == 0.722734247813415611178377352641333362":
ildouble: 1
ldouble: 1
+# acos_downward
+Test "acos_downward (-0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_downward (-0.5) == M_PI_6l*4.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+Test "acos_downward (-1) == pi":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_downward (0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_downward (0.5) == M_PI_6l*2.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+
+# acos_towardzero
+Test "acos_towardzero (-0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (-0.5) == M_PI_6l*4.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+Test "acos_towardzero (-1) == pi":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (0.5) == M_PI_6l*2.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+
+# acos_upward
+Test "acos_upward (-0) == pi/2":
+double: 1
+idouble: 1
+Test "acos_upward (-0.5) == M_PI_6l*4.0":
+ildouble: 1
+ldouble: 1
+Test "acos_upward (-1) == pi":
+double: 1
+idouble: 1
+Test "acos_upward (0) == pi/2":
+double: 1
+idouble: 1
+Test "acos_upward (0.5) == M_PI_6l*2.0":
+ildouble: 1
+ldouble: 1
+
# asin
Test "asin (-0.5) == -pi/6":
ildouble: 1
@@ -22,6 +93,76 @@ Test "asin (1.0) == pi/2":
ildouble: 1
ldouble: 1
+# asin_downward
+Test "asin_downward (-0.5) == -pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_downward (-1.0) == -pi/2":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "asin_downward (0.5) == pi/6":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_downward (1.0) == pi/2":
+float: 1
+ifloat: 1
+
+# asin_tonearest
+Test "asin_tonearest (-0.5) == -pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_tonearest (0.5) == pi/6":
+ildouble: 1
+ldouble: 1
+
+# asin_towardzero
+Test "asin_towardzero (-0.5) == -pi/6":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_towardzero (-1.0) == -pi/2":
+float: 1
+ifloat: 1
+Test "asin_towardzero (0.5) == pi/6":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_towardzero (1.0) == pi/2":
+float: 1
+ifloat: 1
+
+# asin_upward
+Test "asin_upward (-0.5) == -pi/6":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_upward (-1.0) == -pi/2":
+float: 1
+ifloat: 1
+Test "asin_upward (0.5) == pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_upward (1.0) == pi/2":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+
# atanh
Test "atanh (0.75) == 0.972955074527656652552676371721589865":
ildouble: 2
@@ -2132,10 +2273,60 @@ Function: "acos":
ildouble: 1
ldouble: 1
+Function: "acos_downward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "acos_towardzero":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "acos_upward":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+
Function: "asin":
ildouble: 1
ldouble: 1
+Function: "asin_downward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "asin_tonearest":
+ildouble: 1
+ldouble: 1
+
+Function: "asin_towardzero":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "asin_upward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
Function: "atanh":
ildouble: 2
ldouble: 1
diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps
index c700bdf..5068fe6 100644
--- a/sysdeps/x86_64/fpu/libm-test-ulps
+++ b/sysdeps/x86_64/fpu/libm-test-ulps
@@ -5,6 +5,52 @@ Test "acos (0.75) == 0.722734247813415611178377352641333362":
ildouble: 1
ldouble: 1
+# acos_downward
+Test "acos_downward (-0) == pi/2":
+ildouble: 1
+ldouble: 1
+Test "acos_downward (-0.5) == M_PI_6l*4.0":
+double: 1
+idouble: 1
+Test "acos_downward (-1) == pi":
+ildouble: 1
+ldouble: 1
+Test "acos_downward (0) == pi/2":
+ildouble: 1
+ldouble: 1
+Test "acos_downward (0.5) == M_PI_6l*2.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+
+# acos_towardzero
+Test "acos_towardzero (-0) == pi/2":
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (-0.5) == M_PI_6l*4.0":
+double: 1
+idouble: 1
+Test "acos_towardzero (-1) == pi":
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (0) == pi/2":
+ildouble: 1
+ldouble: 1
+Test "acos_towardzero (0.5) == M_PI_6l*2.0":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+
+# acos_upward
+Test "acos_upward (-0.5) == M_PI_6l*4.0":
+ildouble: 1
+ldouble: 1
+Test "acos_upward (0.5) == M_PI_6l*2.0":
+ildouble: 1
+ldouble: 1
+
# asin
Test "asin (-0.5) == -pi/6":
ildouble: 1
@@ -22,6 +68,72 @@ Test "asin (1.0) == pi/2":
ildouble: 1
ldouble: 1
+# asin_downward
+Test "asin_downward (-0.5) == -pi/6":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "asin_downward (0.5) == pi/6":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "asin_downward (1.0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+# asin_tonearest
+Test "asin_tonearest (-0.5) == -pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_tonearest (-1.0) == -pi/2":
+ildouble: 1
+ldouble: 1
+Test "asin_tonearest (0.5) == pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_tonearest (1.0) == pi/2":
+ildouble: 1
+ldouble: 1
+
+# asin_towardzero
+Test "asin_towardzero (-0.5) == -pi/6":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "asin_towardzero (-1.0) == -pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_towardzero (0.5) == pi/6":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "asin_towardzero (1.0) == pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+# asin_upward
+Test "asin_upward (-0.5) == -pi/6":
+ildouble: 1
+ldouble: 1
+Test "asin_upward (-1.0) == -pi/2":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "asin_upward (0.5) == pi/6":
+ildouble: 1
+ldouble: 1
+
# atan2
Test "atan2 (-0.75, -1.0) == -2.49809154479650885165983415456218025":
float: 1
@@ -2078,10 +2190,56 @@ Function: "acos":
ildouble: 1
ldouble: 1
+Function: "acos_downward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "acos_towardzero":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "acos_upward":
+ildouble: 1
+ldouble: 1
+
Function: "asin":
ildouble: 1
ldouble: 1
+Function: "asin_downward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "asin_tonearest":
+ildouble: 1
+ldouble: 1
+
+Function: "asin_towardzero":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "asin_upward":
+float: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
Function: "atan2":
float: 1
ifloat: 1
--
Joseph S. Myers
joseph@codesourcery.com