From patchwork Fri Nov 2 16:59:47 2012 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: Tobias Burnus X-Patchwork-Id: 196604 Return-Path: X-Original-To: incoming@patchwork.ozlabs.org Delivered-To: patchwork-incoming@bilbo.ozlabs.org Received: from sourceware.org (server1.sourceware.org [209.132.180.131]) by ozlabs.org (Postfix) with SMTP id 735812C00D0 for ; Sat, 3 Nov 2012 04:00:18 +1100 (EST) Comment: DKIM? See http://www.dkim.org DKIM-Signature: v=1; a=rsa-sha1; c=relaxed/relaxed; d=gcc.gnu.org; s=default; x=1352480419; h=Comment: DomainKey-Signature:Received:Received:Received:Received: Message-ID:Date:From:User-Agent:MIME-Version:To:Subject: Content-Type:Mailing-List:Precedence:List-Id:List-Unsubscribe: List-Archive:List-Post:List-Help:Sender:Delivered-To; bh=XRXgtf3 zApAvETq9mp4iSYZg2mM=; b=KdQ+qRUMc+DKWxztFCBrMAfCykCtWqLlI89FS+R PeOV4pMoPPwORGQQ0qU8umDdYFoiXxAu+Zx1alUUUjtEc1SHPbXbvtTROiU95lMP Kmu12QAwcO+Lc/yjyU4PS4y7EDRf1VGy3MV44NFnanvNWRSdbxDNICHu1NZqWQGF d4ew= Comment: DomainKeys? See http://antispam.yahoo.com/domainkeys DomainKey-Signature: a=rsa-sha1; q=dns; c=nofws; s=default; d=gcc.gnu.org; h=Received:Received:X-SWARE-Spam-Status:X-Spam-Check-By:Received:Received:Message-ID:Date:From:User-Agent:MIME-Version:To:Subject:Content-Type:Mailing-List:Precedence:List-Id:List-Unsubscribe:List-Archive:List-Post:List-Help:Sender:Delivered-To; b=KiJtBhm5v8WoyAiW0llFZaao6kvioC8Q/oCaCfKlfGp3kPFNfMkm2q1j+0rTHJ Mofm+bEyeWK/w4sO01gFKJPE+RJ8EmF4oGV3on+Nvffrstc27td8Ou+ZCuJ/taqp TWqvceRpaMdzePTrjCU4gwlWqd9/ocm0O4ucjKsGSee1c=; Received: (qmail 31013 invoked by alias); 2 Nov 2012 16:59:58 -0000 Received: (qmail 31005 invoked by uid 22791); 2 Nov 2012 16:59:56 -0000 X-SWARE-Spam-Status: No, hits=-2.0 required=5.0 tests=AWL, BAYES_00, KHOP_SPAMHAUS_DROP, RCVD_IN_DNSWL_NONE, TW_IB, TW_OG X-Spam-Check-By: sourceware.org Received: from mx02.qsc.de (HELO mx02.qsc.de) (213.148.130.14) by sourceware.org (qpsmtpd/0.43rc1) with ESMTP; Fri, 02 Nov 2012 16:59:49 +0000 Received: from [192.168.178.22] (port-92-195-110-241.dynamic.qsc.de [92.195.110.241]) by mx02.qsc.de (Postfix) with ESMTP id F16F825439; Fri, 2 Nov 2012 17:59:47 +0100 (CET) Message-ID: <5093FC03.7000306@net-b.de> Date: Fri, 02 Nov 2012 17:59:47 +0100 From: Tobias Burnus User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:16.0) Gecko/20121025 Thunderbird/16.0.2 MIME-Version: 1.0 To: gcc patches Subject: [patch, committed] libquadmath: Small fmaq and lgamma update Mailing-List: contact gcc-patches-help@gcc.gnu.org; run by ezmlm Precedence: bulk List-Id: List-Unsubscribe: List-Archive: List-Post: List-Help: Sender: gcc-patches-owner@gcc.gnu.org Delivered-To: mailing list gcc-patches@gcc.gnu.org The attached patch ports Joseph's GLIBC fmal patch from today; additionally, I made a small adjustment to lgamma's signgam handling. Build and tested on x86-64-gnu-linux and applied as Rev. 193099 Otherwise, it still applies what I wrote at http://gcc.gnu.org/ml/gcc-patches/2012-11/msg00042.html: Besides those isses [i.e. the test failures], one should - as mentioned in the yesterday's thread - update strtod and printf; the GLIBC test should be modified in such a way that they work properly. I would also like to see the test automatically be run (with "make test") and possibly to have some additional DejaGNU support for libquadmath-specific tests. (One could also add support for __complex__ __float128 as Andreas suggested.) Tobias 2012-11-01 Tobias Burnus Joseph Myers * math/fmaq.c (fmaq): Merge from GLIBC. Handle cases with small x * y using scaling, not as x * y + z. * math/lgammaq.c (lgammaq): Fix to signgam handling. diff --git a/libquadmath/math/fmaq.c b/libquadmath/math/fmaq.c index 5616c1a..e5a9d37 100644 --- a/libquadmath/math/fmaq.c +++ b/libquadmath/math/fmaq.c @@ -73,6 +73,37 @@ fmaq (__float128 x, __float128 y, __float128 z) || u.ieee.exponent + v.ieee.exponent < IEEE854_FLOAT128_BIAS - FLT128_MANT_DIG - 2) return x * y + z; + /* If x * y is less than 1/4 of FLT128_DENORM_MIN, neither the + result nor whether there is underflow depends on its exact + value, only on its sign. */ + if (u.ieee.exponent + v.ieee.exponent + < IEEE854_FLT128_DOUBLE_BIAS - FLT128_MANT_DIG - 2) + { + int neg = u.ieee.negative ^ v.ieee.negative; + __float128 tiny = neg ? -0x1p-16494L : 0x1p-16494L; + if (w.ieee.exponent >= 3) + return tiny + z; + /* Scaling up, adding TINY and scaling down produces the + correct result, because in round-to-nearest mode adding + TINY has no effect and in other modes double rounding is + harmless. But it may not produce required underflow + exceptions. */ + v.value = z * 0x1p114L + tiny; + if (TININESS_AFTER_ROUNDING + ? v.ieee.exponent < 115 + : (w.ieee.exponent == 0 + || (w.ieee.exponent == 1 + && w.ieee.negative != neg + && w.ieee.mantissa3 == 0 + && w.ieee.mantissa2 == 0 + && w.ieee.mantissa1 == 0 + && w.ieee.mantissa0 == 0))) + { + volatile __float128 force_underflow = x * y; + (void) force_underflow; + } + return v.value * 0x1p-114L; + } if (u.ieee.exponent + v.ieee.exponent >= 0x7fff + IEEE854_FLOAT128_BIAS - FLT128_MANT_DIG) { @@ -228,17 +259,17 @@ fmaq (__float128 x, __float128 y, __float128 z) for proper rounding. */ if (v.ieee.exponent == 226) { - /* If the exponent would be in the normal range when - rounding to normal precision with unbounded exponent - range, the exact result is known and spurious underflows - must be avoided on systems detecting tininess after - rounding. */ - if (TININESS_AFTER_ROUNDING) - { - w.value = a1 + u.value; - if (w.ieee.exponent == 227) - return w.value * 0x1p-226L; - } + /* If the exponent would be in the normal range when + rounding to normal precision with unbounded exponent + range, the exact result is known and spurious underflows + must be avoided on systems detecting tininess after + rounding. */ + if (TININESS_AFTER_ROUNDING) + { + w.value = a1 + u.value; + if (w.ieee.exponent == 227) + return w.value * 0x1p-226L; + } /* v.ieee.mant_low & 2 is LSB bit of the result before rounding, v.ieee.mant_low & 1 is the round bit and j is our sticky bit. */ diff --git a/libquadmath/math/lgammaq.c b/libquadmath/math/lgammaq.c index 361f703..485939a 100644 --- a/libquadmath/math/lgammaq.c +++ b/libquadmath/math/lgammaq.c @@ -759,7 +759,8 @@ lgammaq (__float128 x) { __float128 p, q, w, z, nx; int i, nn; - int sign; + + signgam = 1; if (! finiteq (x)) return x * x; @@ -767,11 +768,9 @@ lgammaq (__float128 x) if (x == 0.0Q) { if (signbitq (x)) - sign = -1; + signgam = -1; } - signgam = sign; - if (x < 0.0Q) { q = -x; @@ -780,9 +779,9 @@ lgammaq (__float128 x) return (one / (p - p)); i = p; if ((i & 1) == 0) - sign = -1; + signgam = -1; else - sign = 1; + signgam = 1; z = q - p; if (z > 0.5Q) { @@ -790,10 +789,8 @@ lgammaq (__float128 x) z = p - q; } z = q * sinq (PIQ * z); - signgam = sign; - if (z == 0.0Q) - return (sign * huge * huge); + return (signgam * huge * huge); w = lgammaq (q); z = logq (PIQ / z) - w; return (z); @@ -1025,7 +1022,7 @@ lgammaq (__float128 x) } if (x > MAXLGM) - return (sign * huge * huge); + return (signgam * huge * huge); q = ls2pi - x; q = (x - 0.5Q) * logq (x) + q;