Fix ldbl-128ibm ilogbl near powers of 2 (bug 18029).

The ldbl-128ibm implementation of ilogbl produces incorrect results
when the high part of the argument is a power of 2 and the low part a
nonzero number with the opposite sign (and so the returned exponent
should be 1 less than that of the high part).  For example, ilogbl
(0x1.ffffffffffffffp1L) returns 2 but should return 1.  (This is
similar to (fixed) bug 16740 for frexpl, and bug 18030 for logbl.)
This patch adds checks for that case.

Tested for powerpc.

	[BZ #18029]
	* sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c (__ieee754_ilogbl):
	Adjust exponent of power of 2 down when low part has opposite
	sign.
	* math/libm-test.inc (ilogb_test_data): Add more tests.

4 files changed, 33 insertions, 5 deletions

diff --git a/ChangeLog b/ChangeLog
index 4fbd340764..432c35d5aa 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,11 @@
+2015-02-26  Joseph Myers  <joseph@codesourcery.com>
+
+	[BZ #18029]
+	* sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c (__ieee754_ilogbl):
+	Adjust exponent of power of 2 down when low part has opposite
+	sign.
+	* math/libm-test.inc (ilogb_test_data): Add more tests.
+
 2015-02-26  Alexandre Oliva  <aoliva@redhat.com>
 
 	[BZ #15969]
diff --git a/NEWS b/NEWS
index db6c5f122c..75e83e0a17 100644
--- a/NEWS
+++ b/NEWS
@@ -12,7 +12,7 @@ Version 2.22
   4719, 14841, 13064, 14094, 15319, 15467, 15790, 15969, 16560, 16783,
   17269, 17523, 17569, 17588, 17792, 17836, 17912, 17916, 17932, 17944,
   17949, 17964, 17965, 17967, 17969, 17978, 17987, 17991, 17996, 17998,
-  17999, 18019, 18020.
+  17999, 18019, 18020, 18029.
 
 * Character encoding and ctype tables were updated to Unicode 7.0.0, using
   new generator scripts contributed by Pravin Satpute and Mike FABIAN (Red
diff --git a/math/libm-test.inc b/math/libm-test.inc
index b3396b6a79..6045764fd3 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -7219,6 +7219,11 @@ static const struct test_f_i_data ilogb_test_data[] =
     TEST_f_i (ilogb, 1024, 10, NO_INEXACT_EXCEPTION),
     TEST_f_i (ilogb, -2000, 10, NO_INEXACT_EXCEPTION),
 
+#if defined TEST_LDOUBLE && LDBL_MANT_DIG >= 57
+    TEST_f_i (ilogb, 0x1.ffffffffffffffp1L, 1, NO_INEXACT_EXCEPTION),
+    TEST_f_i (ilogb, -0x1.ffffffffffffffp1L, 1, NO_INEXACT_EXCEPTION),
+#endif
+
     /* ilogb (0.0) == FP_ILOGB0 plus invalid exception  */
     TEST_f_i (ilogb, 0.0, FP_ILOGB0, NO_INEXACT_EXCEPTION|INVALID_EXCEPTION|ERRNO_EDOM),
     /* ilogb (qNaN) == FP_ILOGBNAN plus invalid exception  */
diff --git a/sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c b/sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c
index aeace7c977..4088238f30 100644
--- a/sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c
+++ b/sysdeps/ieee754/ldbl-128ibm/e_ilogbl.c
@@ -31,12 +31,13 @@ static char rcsid[] = "$NetBSD: $";
 
 int __ieee754_ilogbl(long double x)
 {
-	int64_t hx;
+	int64_t hx, hxs;
 	int ix;
-	double xhi;
+	double xhi, xlo;
 
-	xhi = ldbl_high (x);
+	ldbl_unpack (x, &xhi, &xlo);
 	EXTRACT_WORDS64 (hx, xhi);
+	hxs = hx;
 	hx &= 0x7fffffffffffffffLL;
 	if(hx <= 0x0010000000000000LL) {
 	    if(hx==0)
@@ -45,7 +46,21 @@ int __ieee754_ilogbl(long double x)
 		for (ix = -1022, hx<<=11; hx>0; hx<<=1) ix -=1;
 	    return ix;
 	}
-	else if (hx<0x7ff0000000000000LL) return (hx>>52)-0x3ff;
+	else if (hx < 0x7ff0000000000000LL)
+	  {
+	    int hexp = (hx >> 52) - 0x3ff;
+	    /* If the high part is a power of 2, and the low part is
+	       nonzero with the opposite sign, the low part affects
+	       the exponent.  */
+	    if ((hx & 0x000fffffffffffffLL) == 0)
+	      {
+		int64_t lx;
+		EXTRACT_WORDS64 (lx, xlo);
+		if ((hxs ^ lx) < 0 && (lx & 0x7fffffffffffffffLL) != 0)
+		  hexp--;
+	      }
+	    return hexp;
+	  }
 	else if (FP_ILOGBNAN != INT_MAX) {
 	    /* ISO C99 requires ilogbl(+-Inf) == INT_MAX.  */
 	    if (hx==0x7ff0000000000000LL)