niansa-misc/glibc
1
0
Fork 0
mirror of git://sourceware.org/git/glibc.git synced 2025-03-06 20:58:33 +01:00
Code Issues Projects Releases Packages Wiki Activity

Fix acosh (1) in round-downward mode (bug 16927).

Browse source 
According to C99 and C11 Annex F, acosh (1) should be +0 in all
rounding modes.  However, some implementations in glibc wrongly return
-0 in round-downward mode (which is what you get if you end up
computing log1p (-0), via 1 - 1 being -0 in round-downward mode).
This patch fixes the problem implementations, by correcting the test
for an exact 1 value in the ldbl-96 implementation to allow for the
explicit high bit of the mantissa, and by inserting fabs instructions
in the i386 implementations; tests of acosh are duly converted to
ALL_RM_TEST.  I believe all the other sysdeps/ieee754 implementations
are already OK (I haven't checked the ia64 versions, but if buggy then
that will be obvious from the results of test runs after this patch is
in).

Tested x86_64 and x86 and ulps updated accordingly.

	[BZ #16927]
	* sysdeps/i386/fpu/e_acosh.S (__ieee754_acosh): Use fabs on x-1
	value.
	* sysdeps/i386/fpu/e_acoshf.S (__ieee754_acoshf): Likewise.
	* sysdeps/i386/fpu/e_acoshl.S (__ieee754_acoshl): Likewise.
	* sysdeps/ieee754/ldbl-96/e_acoshl.c (__ieee754_acoshl): Correct
	for explicit high bit of mantissa when testing for argument equal
	to 1.
	* math/libm-test.inc (acosh_test): Use ALL_RM_TEST.
	* sysdeps/i386/fpu/libm-test-ulps: Update.
	* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
This commit is contained in:
Joseph Myers 2014-05-14 12:35:40 +00:00
parent 0bf061d3e3
commit 913d03c864
9 changed files with 45 additions and 5 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

12
ChangeLog
View file

@ -1,5 +1,17 @@
2014-05-14 Joseph Myers <joseph@codesourcery.com>
[BZ #16927]
* sysdeps/i386/fpu/e_acosh.S (__ieee754_acosh): Use fabs on x-1
value.
* sysdeps/i386/fpu/e_acoshf.S (__ieee754_acoshf): Likewise.
* sysdeps/i386/fpu/e_acoshl.S (__ieee754_acoshl): Likewise.
* sysdeps/ieee754/ldbl-96/e_acoshl.c (__ieee754_acoshl): Correct
for explicit high bit of mantissa when testing for argument equal
to 1.
* math/libm-test.inc (acosh_test): Use ALL_RM_TEST.
* sysdeps/i386/fpu/libm-test-ulps: Update.
* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
[BZ #16516]
* sysdeps/ieee754/dbl-64/s_erf.c (efx8): Remove variable.
(__erf): Scale by 16 instead of 8 in potentially underflowing

2
NEWS
View file

@ -17,7 +17,7 @@ Version 2.20
16713, 16714, 16731, 16739, 16740, 16743, 16754, 16758, 16759, 16760,
16770, 16786, 16789, 16791, 16799, 16800, 16815, 16823, 16824, 16831,
16838, 16854, 16876, 16877, 16885, 16888, 16890, 16912, 16916, 16922,
16932.
16927, 16932.
* The minimum Linux kernel version that this version of the GNU C Library
can be used with is 2.6.32.

4
math/libm-test.inc
View file

@ -1792,9 +1792,7 @@ static const struct test_f_f_data acosh_test_data[] =
static void
acosh_test (void)
{
START (acosh, 0);
RUN_TEST_LOOP_f_f (acosh, acosh_test_data, );
END;
ALL_RM_TEST (acosh, 0, acosh_test_data, RUN_TEST_LOOP_f_f, END);
}
static const struct test_f_f_data asin_test_data[] =

1
sysdeps/i386/fpu/e_acosh.S
View file

@ -52,6 +52,7 @@ ENTRY(__ieee754_acosh)
// 1 <= x <= 2 => y = log1p(x-1+sqrt(2*(x-1)+(x-1)^2))
fsubl MO(one) // x-1 : log(2)
fabs // acosh(1) is +0 in all rounding modes
fld %st // x-1 : x-1 : log(2)
fmul %st(1) // (x-1)^2 : x-1 : log(2)
fadd %st(1) // x-1+(x-1)^2 : x-1 : log(2)

1
sysdeps/i386/fpu/e_acoshf.S
View file

@ -52,6 +52,7 @@ ENTRY(__ieee754_acoshf)
// 1 <= x <= 2 => y = log1p(x-1+sqrt(2*(x-1)+(x-1)^2))
fsubl MO(one) // x-1 : log(2)
fabs // acosh(1) is +0 in all rounding modes
fld %st // x-1 : x-1 : log(2)
fmul %st(1) // (x-1)^2 : x-1 : log(2)
fadd %st(1) // x-1+(x-1)^2 : x-1 : log(2)

1
sysdeps/i386/fpu/e_acoshl.S
View file

@ -59,6 +59,7 @@ ENTRY(__ieee754_acoshl)
// 1 <= x <= 2 => y = log1p(x-1+sqrt(2*(x-1)+(x-1)^2))
fsubl MO(one) // x-1 : log(2)
fabs // acosh(1) is +0 in all rounding modes
fld %st // x-1 : x-1 : log(2)
fmul %st(1) // (x-1)^2 : x-1 : log(2)
fadd %st(1) // x-1+(x-1)^2 : x-1 : log(2)

12
sysdeps/i386/fpu/libm-test-ulps
View file

@ -21,6 +21,18 @@ Function: "acos_upward":
ildouble: 1
ldouble: 1
Function: "acosh_downward":
ildouble: 2
ldouble: 2
Function: "acosh_towardzero":
ildouble: 2
ldouble: 2
Function: "acosh_upward":
ildouble: 1
ldouble: 1
Function: "asin_downward":
double: 1
float: 1

2
sysdeps/ieee754/ldbl-96/e_acoshl.c
View file

@ -48,7 +48,7 @@ __ieee754_acoshl(long double x)
return x+x;
} else
return __ieee754_logl(x)+ln2; /* acoshl(huge)=logl(2x) */
} else if(((se-0x3fff)|i0|i1)==0) {
} else if(((se-0x3fff)|(i0^0x80000000)|i1)==0) {
return 0.0; /* acosh(1) = 0 */
} else if (se > 0x4000) { /* 2**28 > x > 2 */
t=x*x;

15
sysdeps/x86_64/fpu/libm-test-ulps
View file

@ -31,6 +31,21 @@ Function: "acosh":
double: 1
idouble: 1
Function: "acosh_downward":
float: 1
ildouble: 1
ldouble: 2
Function: "acosh_towardzero":
float: 1
ildouble: 1
ldouble: 2
Function: "acosh_upward":
double: 1
ildouble: 1
ldouble: 1
Function: "asin_downward":
double: 1
float: 1