From patchwork Wed Jan 11 17:22:32 2017 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: Joseph Myers X-Patchwork-Id: 18857 Received: (qmail 125287 invoked by alias); 11 Jan 2017 17:22:56 -0000 Mailing-List: contact libc-alpha-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Unsubscribe: List-Subscribe: List-Archive: List-Post: List-Help: , Sender: libc-alpha-owner@sourceware.org Delivered-To: mailing list libc-alpha@sourceware.org Received: (qmail 125157 invoked by uid 89); 11 Jan 2017 17:22:55 -0000 Authentication-Results: sourceware.org; auth=none X-Virus-Found: No X-Spam-SWARE-Status: No, score=-1.9 required=5.0 tests=AWL, BAYES_00, RCVD_IN_DNSWL_NONE, SPF_PASS, URIBL_RED autolearn=ham version=3.3.2 spammy=sigma, Tulio, libm-test.inc, libmtestinc X-HELO: relay1.mentorg.com Date: Wed, 11 Jan 2017 17:22:32 +0000 From: Joseph Myers To: Tulio Magno Quites Machado Filho CC: Subject: Re: XFAIL libm-test.inc tests as needed for ibm128 [committed] In-Reply-To: <87eg09k2ou.fsf@linux.vnet.ibm.com> Message-ID: References: <87eg09k2ou.fsf@linux.vnet.ibm.com> User-Agent: Alpine 2.20 (DEB 67 2015-01-07) MIME-Version: 1.0 X-ClientProxiedBy: svr-ies-mbx-01.mgc.mentorg.com (139.181.222.1) To svr-ies-mbx-01.mgc.mentorg.com (139.181.222.1) On Wed, 11 Jan 2017, Tulio Magno Quites Machado Filho wrote: > > (The case of patched libgcc was already clean up to ulps > > Which patch are you referring to? This one (only suitable for building GCC specifically for testing glibc, not for more general use, since the changes are unconditional and make performance significantly worse). This patch doesn't fully support soft float (that would require an additional interface exported from glibc) and as a minimal patch doesn't include testcases for the various bugs fixed (although I have such testcases in an earlier patch version). With this patch it should be possible to define TEST_COND_ibm128_libgcc to 0 in libm-test.inc, and also possible to do random test generation to find cases of large errors in the ldbl-128ibm functions without running into problems with libgcc issues (being able to use random test generation to find bugs in glibc is probably the clearest use of this patch). Index: libgcc/config/rs6000/ibm-ldouble.c =================================================================== --- libgcc/config/rs6000/ibm-ldouble.c (revision 228280) +++ libgcc/config/rs6000/ibm-ldouble.c (working copy) @@ -47,6 +47,10 @@ see the files COPYING3 and COPYING.RUNTIME respect #if defined (__MACH__) || defined (__powerpc__) || defined (_AIX) +typedef _Bool bool; +#define false 0 +#define true 1 + #define fabs(x) __builtin_fabs(x) #define isless(x, y) __builtin_isless (x, y) #define inf() __builtin_inf() @@ -87,6 +91,58 @@ __asm__ (".symver __gcc_qadd,_xlqadd@GCC_3.4\n\t" ".symver .__gcc_qdiv,._xlqdiv@GCC_3.4"); #endif +#define FE_TONEAREST 0 + +/* Set the rounding mode to MODE, returning the previous mode. */ + +static int +swapround (int mode) +{ + int old_mode; +#ifdef __NO_FPRS__ +# ifdef _SOFT_FLOAT + /* Not implemented. */ + old_mode = FE_TONEAREST; +# else + /* SPE floating point. */ + int spefscr; + asm volatile ("mfspefscr %0" : "=r" (spefscr)); + old_mode = spefscr & 3; + spefscr = (spefscr & ~3) | mode; + asm volatile ("mtspefscr %0" : : "r" (spefscr)); +# endif +#else + /* Logic as in glibc. */ + asm volatile ("mcrfs 7,7\n\t" + "mfcr %0" : "=r" (old_mode) : : "cr7"); + old_mode &= 3; + if ((mode ^ old_mode) & 2) + { + if (mode & 2) + asm volatile ("mtfsb1 30"); + else + asm volatile ("mtfsb0 30"); + } + if ((mode ^ old_mode) & 1) + { + if (mode & 1) + asm volatile ("mtfsb1 31"); + else + asm volatile ("mtfsb0 31"); + } +#endif + return old_mode; +} + +/* Overflowing values in each rounding mode (FE_TONEAREST, + FE_TOWARDZERO, FE_UPWARD, FE_DOWNWARD). */ + +static const long double overflow_values[2][4] = + { + { __builtin_infl (), __LDBL_MAX__, __builtin_infl (), __LDBL_MAX__ }, + { -__builtin_infl (), -__LDBL_MAX__, -__LDBL_MAX__, -__builtin_infl () } + }; + /* Combine two 'double' values into one 'long double' and return the result. */ static inline long double pack_ldouble (double dh, double dl) @@ -110,40 +166,74 @@ pack_ldouble (double dh, double dl) long double __gcc_qadd (double a, double aa, double c, double cc) { - double xh, xl, z, q, zz; + double xh, xhs, xl, z, q, zz; + bool scale = false; + int old_mode; + if (nonfinite (a) || nonfinite (c)) + return a + c; + + old_mode = swapround (FE_TONEAREST); + + if (fabs (a) >= 0x1p1022 || fabs (c) >= 0x1p1022) + { + scale = true; + if (fabs (a) >= 0x1p914) + { + a *= 0.5; + aa *= 0.5; + } + if (fabs (c) >= 0x1p914) + { + c *= 0.5; + cc *= 0.5; + } + } + z = a + c; + q = a - z; + zz = q + c + (a - (q + z)) + aa + cc; - if (nonfinite (z)) + /* Keep -0 result. */ + if (zz == 0.0) { - if (fabs (z) != inf()) - return z; - z = cc + aa + c + a; - if (nonfinite (z)) - return z; - xh = z; /* Will always be DBL_MAX. */ - zz = aa + cc; - if (fabs(a) > fabs(c)) - xl = a - z + c + zz; - else - xl = c - z + a + zz; + double ret = scale ? 2.0 * z : z; + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (ret) == inf ()) + return overflow_values[ret < 0][old_mode]; + if (ret == 0) + { + asm volatile ("" : "+m" (a)); + ret = a + c; + } + } + return ret; } - else + xh = z + zz; + xhs = xh; + if (scale) + xhs *= 2.0; + if (nonfinite (xhs)) { - q = a - z; - zz = q + c + (a - (q + z)) + aa + cc; + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (xhs) == inf ()) + return overflow_values[xhs < 0][old_mode]; + } + return xhs; + } - /* Keep -0 result. */ - if (zz == 0.0) - return z; + xl = z - xh + zz; + if (scale) + xl *= 2.0; - xh = z + zz; - if (nonfinite (xh)) - return xh; + if (unlikely (old_mode != FE_TONEAREST)) + swapround (old_mode); - xl = z - xh + zz; - } - return pack_ldouble (xh, xl); + return pack_ldouble (xhs, xl); } long double @@ -159,13 +249,41 @@ static double fmsub (double, double, double); long double __gcc_qmul (double a, double b, double c, double d) { - double xh, xl, t, tau, u, v, w; + double xh, xl, t, tau, u, us, v, w; + bool scale = false; + int old_mode; + + if (nonfinite (a) || nonfinite (c)) + return a * c; + + old_mode = swapround (FE_TONEAREST); + + if (fabs (a) >= 0x1p512) + { + scale = true; + a *= 0.5; + b *= 0.5; + } + else if (fabs (c) >= 0x1p512) + { + scale = true; + c *= 0.5; + d *= 0.5; + } t = a * c; /* Highest order double term. */ if (unlikely (t == 0) /* Preserve -0. */ || nonfinite (t)) - return t; + { + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (t) == inf ()) + return overflow_values[t < 0][old_mode]; + } + return t; + } /* Sum terms of two highest orders. */ @@ -179,12 +297,29 @@ __gcc_qmul (double a, double b, double c, double d w = b*c; tau += v + w; /* Add in other second-order terms. */ u = t + tau; + us = u; + if (scale) + us *= 2.0; /* Construct long double result. */ - if (nonfinite (u)) - return u; - xh = u; + if (nonfinite (us)) + { + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (us) == inf ()) + return overflow_values[us < 0][old_mode]; + } + return us; + } + xh = us; xl = (t - u) + tau; + if (scale) + xl *= 2.0; + + if (unlikely (old_mode != FE_TONEAREST)) + swapround (old_mode); + return pack_ldouble (xh, xl); } @@ -191,13 +326,41 @@ __gcc_qmul (double a, double b, double c, double d long double __gcc_qdiv (double a, double b, double c, double d) { - double xh, xl, s, sigma, t, tau, u, v, w; + double xh, xl, s, sigma, t, tau, u, us, v, w; + bool scale = false; + int old_mode; + + if (nonfinite (a) || nonfinite (c) || c == 0) + return a / c; + + old_mode = swapround (FE_TONEAREST); + + if (fabs (a) >= 0x1p512) + { + scale = true; + a *= 0.5; + b *= 0.5; + } + else if (fabs (c) <= 0x1p-511) + { + scale = true; + c *= 2.0; + d *= 2.0; + } t = a / c; /* highest order double term */ if (unlikely (t == 0) /* Preserve -0. */ || nonfinite (t)) - return t; + { + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (t) == inf ()) + return overflow_values[t < 0][old_mode]; + } + return t; + } /* Finite nonzero result requires corrections to the highest order term. These corrections require the low part of c * t to be @@ -224,12 +387,29 @@ __gcc_qdiv (double a, double b, double c, double d tau = ((v-sigma)+w)/c; /* Correction to t. */ u = t + tau; + us = u; + if (scale) + us *= 2.0; /* Construct long double result. */ - if (nonfinite (u)) - return u; - xh = u; + if (nonfinite (us)) + { + if (unlikely (old_mode != FE_TONEAREST)) + { + swapround (old_mode); + if (fabs (us) == inf ()) + return overflow_values[us < 0][old_mode]; + } + return us; + } + xh = us; xl = (t - u) + tau; + if (scale) + xl *= 2.0; + + if (unlikely (old_mode != FE_TONEAREST)) + swapround (old_mode); + return pack_ldouble (xh, xl); }