Fix ldbl-128 expm1l (-min_subnorm) result sign (bug 18619) [committed]

[Joseph Myers]

In the ldbl-128 implementation of expm1l, when expm1l's result should
underflow to 0 (argument minus the least subnormal, in some rounding
modes), it can be a zero of the wrong sign.  This patch fixes this in
the same way previously used for the x86 / x86_64 versions.

Tested for mips64.  Committed.

2015-07-01  Joseph Myers  <joseph@codesourcery.com>

	[BZ #18619]
	* sysdeps/ieee754/ldbl-128/s_expm1l.c (__expm1l): Force underflow
	and return argument in case of subnormal argument.

diff --git a/sysdeps/ieee754/ldbl-128/s_expm1l.c b/sysdeps/ieee754/ldbl-128/s_expm1l.c
index f708af5..573d00b 100644
--- a/sysdeps/ieee754/ldbl-128/s_expm1l.c
+++ b/sysdeps/ieee754/ldbl-128/s_expm1l.c
@@ -137,9 +137,18 @@ __expm1l (long double x)
   if (x < minarg)
     return (4.0/big - 1.0L);
 
-  /* Avoid internal underflow when result does not underflow.  */
-  if (fabsl (x) < 0x1p-113L && fabsl (x) >= LDBL_MIN)
-    return x;
+  /* Avoid internal underflow when result does not underflow, while
+     ensuring underflow (without returning a zero of the wrong sign)
+     when the result does underflow.  */
+  if (fabsl (x) < 0x1p-113L)
+    {
+      if (fabsl (x) < LDBL_MIN)
+	{
+	  long double force_underflow = x * x;
+	  math_force_eval (force_underflow);
+	}
+      return x;
+    }
 
   /* Express x = ln 2 (k + remainder), remainder not exceeding 1/2. */
   xx = C1 + C2;			/* ln 2. */

-- 
Joseph S. Myers
joseph@codesourcery.com