Fix ldbl-128 powl sign of result in overflow / underflow cases (bug 17097).

This patch fixes bug 17097, ldbl-128 powl producing overflowing /
underflowing results with positive sign when the result should have
been negative.  This was shown up by the tests in non-default rounding
modes added by my patch for bug 16315, but isn't actually limited to
non-default rounding modes: rather, when rounding to nearest the
wrappers produced a result with the correct sign and so always hid the
bug unless -lieee was used to disable the wrappers.  The problem is
that in the cases where Y is large enough that the result overflows or
underflows for X not very close to 1, but not large enough to overflow
or underflow for all X != +/- 1 (in the latter case Y is always an
even integer), a positive overflowing / underflowing result is always
returned, rather than one with the correct sign.  This patch moves the
relevant part of computation of the sign earlier and returns a result
of the correct sign.

Tested for mips64.

	[BZ #17097]
	* sysdeps/ieee754/ldbl-128/e_powl.c (__ieee754_powl): Return
	result with correct sign in case of exponents that produce
	overflow except for X very close to 1.

diff --git a/ChangeLog b/ChangeLog
index 33ce8eded702c01f9139daa413f757a8cb6762f0..070be3c6a5e18238eba280ecf0604ecfeaa7efcb 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,10 @@
+2014-06-29  Joseph Myers  <joseph@codesourcery.com>
+
+       [BZ #17097]
+       * sysdeps/ieee754/ldbl-128/e_powl.c (__ieee754_powl): Return
+       result with correct sign in case of exponents that produce
+       overflow except for X very close to 1.
+
 2014-06-28  Paul Eggert  <eggert@cs.ucla.edu>
 
        mktime: merge #if/#ifdef usage from glibc

diff --git a/NEWS b/NEWS
index 02e3cd83ff8c4d5ca7594438a065f06aeaf81e09..a07ea6615c34b6358ee70c24aa2d2bb32cf09f16 100644 (file)
--- a/NEWS
+++ b/NEWS
@@ -21,7 +21,7 @@ Version 2.20
   16882, 16885, 16888, 16890, 16912, 16915, 16916, 16917, 16918, 16922,
   16927, 16928, 16932, 16943, 16958, 16965, 16966, 16967, 16977, 16978,
   16984, 16990, 16996, 17009, 17022, 17031, 17042, 17048, 17050, 17058,
-  17061, 17062, 17069, 17075, 17079, 17084, 17086, 17092.
+  17061, 17062, 17069, 17075, 17079, 17084, 17086, 17092, 17097.
 
 * Optimized strchr implementation for AArch64.  Contributed by ARM Ltd.
 

diff --git a/sysdeps/ieee754/ldbl-128/e_powl.c b/sysdeps/ieee754/ldbl-128/e_powl.c
index d131750718d287ef435da9be05d2262a32df44d7..f531385232625bd4002a0a84caaf74015b4377f8 100644 (file)
--- a/sysdeps/ieee754/ldbl-128/e_powl.c
+++ b/sysdeps/ieee754/ldbl-128/e_powl.c
@@ -148,7 +148,7 @@ long double
 __ieee754_powl (long double x, long double y)
 {
   long double z, ax, z_h, z_l, p_h, p_l;
-  long double y1, t1, t2, r, s, t, u, v, w;
+  long double y1, t1, t2, r, s, sgn, t, u, v, w;
   long double s2, s_h, s_l, t_h, t_l, ay;
   int32_t i, j, k, yisint, n;
   u_int32_t ix, iy;
@@ -262,6 +262,11 @@ __ieee754_powl (long double x, long double y)
   if (((((u_int32_t) hx >> 31) - 1) | yisint) == 0)
     return (x - x) / (x - x);
 
+  /* sgn (sign of result -ve**odd) = -1 else = 1 */
+  sgn = one;
+  if (((((u_int32_t) hx >> 31) - 1) | (yisint - 1)) == 0)
+    sgn = -one;                        /* (-ve)**(odd int) */
+
   /* |y| is huge.
      2^-16495 = 1/2 of smallest representable value.
      If (1 - 1/131072)^y underflows, y > 1.4986e9 */
@@ -277,9 +282,9 @@ __ieee754_powl (long double x, long double y)
        }
       /* over/underflow if x is not close to one */
       if (ix < 0x3ffeffff)
-       return (hy < 0) ? huge * huge : tiny * tiny;
+       return (hy < 0) ? sgn * huge * huge : sgn * tiny * tiny;
       if (ix > 0x3fff0000)
-       return (hy > 0) ? huge * huge : tiny * tiny;
+       return (hy > 0) ? sgn * huge * huge : sgn * tiny * tiny;
     }
 
   ay = y > 0 ? y : -y;
@@ -366,11 +371,6 @@ __ieee754_powl (long double x, long double y)
   t1 = o.value;
   t2 = z_l - (((t1 - t) - dp_h[k]) - z_h);
 
-  /* s (sign of result -ve**odd) = -1 else = 1 */
-  s = one;
-  if (((((u_int32_t) hx >> 31) - 1) | (yisint - 1)) == 0)
-    s = -one;                  /* (-ve)**(odd int) */
-
   /* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
   y1 = y;
   o.value = y1;
@@ -386,11 +386,11 @@ __ieee754_powl (long double x, long double y)
     {
       /* if z > 16384 */
       if (((j - 0x400d0000) | o.parts32.w1 | o.parts32.w2 | o.parts32.w3) != 0)
-       return s * huge * huge; /* overflow */
+       return sgn * huge * huge;       /* overflow */
       else
        {
          if (p_l + ovt > z - p_h)
-           return s * huge * huge;     /* overflow */
+           return sgn * huge * huge;   /* overflow */
        }
     }
   else if ((j & 0x7fffffff) >= 0x400d01b9)     /* z <= -16495 */
@@ -398,11 +398,11 @@ __ieee754_powl (long double x, long double y)
       /* z < -16495 */
       if (((j - 0xc00d01bc) | o.parts32.w1 | o.parts32.w2 | o.parts32.w3)
          != 0)
-       return s * tiny * tiny; /* underflow */
+       return sgn * tiny * tiny;       /* underflow */
       else
        {
          if (p_l <= z - p_h)
-           return s * tiny * tiny;     /* underflow */
+           return sgn * tiny * tiny;   /* underflow */
        }
     }
   /* compute 2**(p_h+p_l) */
@@ -441,6 +441,6 @@ __ieee754_powl (long double x, long double y)
       o.parts32.w0 = j;
       z = o.value;
     }
-  return s * z;
+  return sgn * z;
 }
 strong_alias (__ieee754_powl, __powl_finite)