Fix powl inaccuracy for ldbl-128ibm (bug 14914)

["Joseph S. Myers" <joseph at codesourcery dot com>]

Bug 14914 is inaccuracy of powl for ldbl-128ibm, which showed up in
validation for 2.17.

In various places the code, based on that for ldbl-128, tries to split
numbers into high and low parts in order to carry out operations on
those parts separately, with at least some of the high-part operations
needing to be exact.  However, the way the code extracts high parts is
to clear the low 47 bits of the mantissa.  I've no idea how that
figure was chosen - the ldbl-128 code clears 59 bits, i.e. more than
half the mantissa - but as shown by the large errors for powl in the
testsuite, it doesn't work well.  This patch changes the code to clear
the whole low double of the long double value.  This eliminates all
the ulps seen for the tests recently added or enabled for long double;
no ulps at all are needed for pow for powerpc test-ldouble after this
patch.  (Tested powerpc32, hard float; if this goes in then I'll
retest powerpc-nofpu to see if any pow ulps updates are needed there.)

2012-12-03  Joseph Myers  <joseph@codesourcery.com>

	[BZ #14914]
	* sysdeps/ieee754/ldbl-128ibm/e_powl.c (__ieee754_powl): Clear
	whole low double instead of just low 47 bits when splitting values
	into two parts.

diff --git a/sysdeps/ieee754/ldbl-128ibm/e_powl.c b/sysdeps/ieee754/ldbl-128ibm/e_powl.c
index 8216c49..8bd35d0 100644
--- a/sysdeps/ieee754/ldbl-128ibm/e_powl.c
+++ b/sysdeps/ieee754/ldbl-128ibm/e_powl.c
@@ -324,13 +324,13 @@ __ieee754_powl (long double x, long double y)
 
   o.value = s_h;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   s_h = o.value;
   /* t_h=ax+bp[k] High */
   t_h = ax + bp[k];
   o.value = t_h;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   t_h = o.value;
   t_l = ax - (t_h - bp[k]);
   s_l = v * ((u - s_h * t_h) - s_h * t_l);
@@ -344,7 +344,7 @@ __ieee754_powl (long double x, long double y)
   t_h = 3.0 + s2 + r;
   o.value = t_h;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   t_h = o.value;
   t_l = r - ((t_h - 3.0) - s2);
   /* u+v = s*(1+...) */
@@ -354,7 +354,7 @@ __ieee754_powl (long double x, long double y)
   p_h = u + v;
   o.value = p_h;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   p_h = o.value;
   p_l = v - (p_h - u);
   z_h = cp_h * p_h;		/* cp_h+cp_l = 2/(3*log2) */
@@ -364,7 +364,7 @@ __ieee754_powl (long double x, long double y)
   t1 = (((z_h + z_l) + dp_h[k]) + t);
   o.value = t1;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   t1 = o.value;
   t2 = z_l - (((t1 - t) - dp_h[k]) - z_h);
 
@@ -377,7 +377,7 @@ __ieee754_powl (long double x, long double y)
   y1 = y;
   o.value = y1;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   y1 = o.value;
   p_l = (y - y1) * t1 + y * t2;
   p_h = y1 * t1;
@@ -421,7 +421,7 @@ __ieee754_powl (long double x, long double y)
   t = p_l + p_h;
   o.value = t;
   o.parts32.w3 = 0;
-  o.parts32.w2 &= 0xffff8000;
+  o.parts32.w2 = 0;
   t = o.value;
   u = t * lg2_h;
   v = (p_l - (t - p_h)) * lg2 + t * lg2_l;

-- 
Joseph S. Myers
joseph@codesourcery.com