diff options
author | Joseph Myers <joseph@codesourcery.com> | 2014-06-29 11:49:08 +0000 |
---|---|---|
committer | Joseph Myers <joseph@codesourcery.com> | 2014-06-29 11:49:08 +0000 |
commit | edea402804bce917cfd7cd1af76212e6364c23db (patch) | |
tree | d7688292b5b41de2f5707affccea4b992db30df2 | |
parent | dd0ba018122e88937a5f14b6594b9a40693b2e58 (diff) | |
download | glibc-edea402804bce917cfd7cd1af76212e6364c23db.tar glibc-edea402804bce917cfd7cd1af76212e6364c23db.tar.gz glibc-edea402804bce917cfd7cd1af76212e6364c23db.tar.bz2 glibc-edea402804bce917cfd7cd1af76212e6364c23db.zip |
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.
-rw-r--r-- | ChangeLog | 7 | ||||
-rw-r--r-- | NEWS | 2 | ||||
-rw-r--r-- | sysdeps/ieee754/ldbl-128/e_powl.c | 26 |
3 files changed, 21 insertions, 14 deletions
@@ -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 @@ -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 d131750718..f531385232 100644 --- 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) |