From 8a9fa0086dfe1b2b9d828712ef41fc7dd4b8f94c Mon Sep 17 00:00:00 2001
From: Joseph Myers <joseph@codesourcery.com>
Date: Thu, 18 Feb 2016 21:40:39 +0000
Subject: Fix ldbl-128ibm ceill for non-default rounding modes (bug 19592).

The ldbl-128ibm implementation of ceill is only correct in
round-to-nearest mode (in other modes, there are incorrect results and
overflow exceptions in some cases).  It is also unnecessarily
complicated, rounding both high and low parts to the nearest integer
and then adjusting for the semantics of ceil, when it seems more
natural to take the ceiling of the high part (__ceil optimized inline
versions can be used), and that of the low part if the high part is an
integer, as was done for floorl.  This patch makes it use that simpler
approach.

Tested for powerpc.

	[BZ #19592]
	* sysdeps/ieee754/ldbl-128ibm/s_ceill.c (__ceill): Use __ceil on
	high and low parts then use ldbl_canonicalize_int if needed.
---
 sysdeps/ieee754/ldbl-128ibm/s_ceill.c | 52 +++++++++++------------------------
 1 file changed, 16 insertions(+), 36 deletions(-)

(limited to 'sysdeps/ieee754/ldbl-128ibm/s_ceill.c')

diff --git a/sysdeps/ieee754/ldbl-128ibm/s_ceill.c b/sysdeps/ieee754/ldbl-128ibm/s_ceill.c
index ac649b7..635fddc 100644
--- a/sysdeps/ieee754/ldbl-128ibm/s_ceill.c
+++ b/sysdeps/ieee754/ldbl-128ibm/s_ceill.c
@@ -35,42 +35,22 @@ __ceill (long double x)
 			&& __builtin_isless (__builtin_fabs (xh),
 					     __builtin_inf ()), 1))
     {
-      double orig_xh;
-
-      /* Long double arithmetic, including the canonicalisation below,
-	 only works in round-to-nearest mode.  */
-
-      /* Convert the high double to integer.  */
-      orig_xh = xh;
-      hi = ldbl_nearbyint (xh);
-
-      /* Subtract integral high part from the value.  */
-      xh -= hi;
-      ldbl_canonicalize (&xh, &xl);
-
-      /* Now convert the low double, adjusted for any remainder from the
-         high double.  */
-      lo = ldbl_nearbyint (xh);
-
-      /* Adjust the result when the remainder is non-zero.  nearbyint
-         rounds values to the nearest integer, and values halfway
-         between integers to the nearest even integer.  ceill must
-         round towards +Inf.  */
-      xh -= lo;
-      ldbl_canonicalize (&xh, &xl);
-
-      if (xh > 0.0 || (xh == 0.0 && xl > 0.0))
-	lo += 1.0;
-
-      /* Ensure the final value is canonical.  In certain cases,
-         rounding causes hi,lo calculated so far to be non-canonical.  */
-      xh = hi;
-      xl = lo;
-      ldbl_canonicalize (&xh, &xl);
-
-      /* Ensure we return -0 rather than +0 when appropriate.  */
-      if (orig_xh < 0.0)
-	xh = -__builtin_fabs (xh);
+      hi = __ceil (xh);
+      if (hi != xh)
+	{
+	  /* The high part is not an integer; the low part does not
+	     affect the result.  */
+	  xh = hi;
+	  xl = 0;
+	}
+      else
+	{
+	  /* The high part is a nonzero integer.  */
+	  lo = __ceil (xl);
+	  xh = hi;
+	  xl = lo;
+	  ldbl_canonicalize_int (&xh, &xl);
+	}
     }
 
   return ldbl_pack (xh, xl);
-- 
cgit v1.1