From patchwork Thu Nov 15 17:25:41 2012 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: Tobias Burnus X-Patchwork-Id: 199367 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 B11302C048A for ; Fri, 16 Nov 2012 04:26:02 +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=1353605163; 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=cFSX+OY E9Wori0Ump55wivc/rAE=; b=x0kmbXe7+DjIYs2U0X5j7961ikyGYFjmLBi24CZ LvZsh/r77a8qIa07RFaO61vI8Ki37LNWBphv4L0rPWTWQ6QzZav1yTxGuj7wWglk k930Bn699GsayDGewlg5VCmLnzR7Jo3g9kHiVWyuLEvVWH5TSbtl0TNY9dgCIrdY vtGI= 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=hgn9QN+5XPyklFDKAnmEDRrQ1Guh9jC+9eHI+BQE4ewiyXKZ4V2RlkSAWR8AoF oxKRHdiHMaxluOXeQJVR78tdv/N2R/rAozk22TmlGrxmc6XTDpm7iEEgFIu3VvH9 8ypb3FQ4sEAs4lQa6tOKUrhfuHuDyK55wPXp7EOGj0hZk=; Received: (qmail 28720 invoked by alias); 15 Nov 2012 17:25:53 -0000 Received: (qmail 28708 invoked by uid 22791); 15 Nov 2012 17:25:51 -0000 X-SWARE-Spam-Status: No, hits=-2.2 required=5.0 tests=AWL, BAYES_00, RCVD_IN_DNSWL_NONE X-Spam-Check-By: sourceware.org Received: from mx01.qsc.de (HELO mx01.qsc.de) (213.148.129.14) by sourceware.org (qpsmtpd/0.43rc1) with ESMTP; Thu, 15 Nov 2012 17:25:45 +0000 Received: from [192.168.188.20] (port-92-195-100-57.dynamic.qsc.de [92.195.100.57]) by mx01.qsc.de (Postfix) with ESMTP id EFDCE3CDEF; Thu, 15 Nov 2012 18:25:41 +0100 (CET) Message-ID: <50A52595.8020600@net-b.de> Date: Thu, 15 Nov 2012 18:25:41 +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: Update math/fmaq.c 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 Dear all, I have committed (Rev. 193538) attached patch, which does an other update from GLIBC. Tobias PS: I still want to update libquadmath's strtod and printf. Index: libquadmath/ChangeLog =================================================================== --- libquadmath/ChangeLog (Revision 193537) +++ libquadmath/ChangeLog (Arbeitskopie) @@ -1,3 +1,11 @@ +2012-11-15 Tobias Burnus + Joseph Myers + + * math/fmaq.c (fmaq): Merge from GLIBC. Fix fma + underflows with small x * y; Fix overflow results + outside round-to-nearest mode; make use of Dekker + and Knuth algorithms use round-to-nearest. + 2012-11-01 Tobias Burnus * math/fmaq.c (fmaq): Fix build. Index: libquadmath/math/fmaq.c =================================================================== --- libquadmath/math/fmaq.c (Revision 193537) +++ libquadmath/math/fmaq.c (Arbeitskopie) @@ -14,9 +14,8 @@ Lesser General Public License for more details. You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, write to the Free - Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA - 02111-1307 USA. */ + License along with the GNU C Library; if not, see + . */ #include "quadmath-imp.h" #include @@ -62,17 +61,18 @@ fmaq (__float128 x, __float128 y, __float128 z) underflows to 0. */ if (z == 0 && x != 0 && y != 0) return x * y; - /* If x or y or z is Inf/NaN, or if fma will certainly overflow, - or if x * y is less than half of FLT128_DENORM_MIN, - compute as x * y + z. */ + /* If x or y or z is Inf/NaN, or if x * y is zero, compute as + x * y + z. */ if (u.ieee.exponent == 0x7fff || v.ieee.exponent == 0x7fff || w.ieee.exponent == 0x7fff - || u.ieee.exponent + v.ieee.exponent - > 0x7fff + IEEE854_FLOAT128_BIAS - || u.ieee.exponent + v.ieee.exponent - < IEEE854_FLOAT128_BIAS - FLT128_MANT_DIG - 2) + || x == 0 + || y == 0) return x * y + z; + /* If fma will certainly overflow, compute as x * y. */ + if (u.ieee.exponent + v.ieee.exponent + > 0x7fff + IEEE854_FLOAT128_BIAS) + return x * y; /* 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. */ @@ -121,9 +121,18 @@ fmaq (__float128 x, __float128 y, __float128 z) { /* Similarly. If z exponent is very large and x and y exponents are - very small, it doesn't matter if we don't adjust it. */ - if (u.ieee.exponent > v.ieee.exponent) + very small, adjust them up to avoid spurious underflows, + rather than down. */ + if (u.ieee.exponent + v.ieee.exponent + <= IEEE854_FLOAT128_BIAS + FLT128_MANT_DIG) { + if (u.ieee.exponent > v.ieee.exponent) + u.ieee.exponent += 2 * FLT128_MANT_DIG + 2; + else + v.ieee.exponent += 2 * FLT128_MANT_DIG + 2; + } + else if (u.ieee.exponent > v.ieee.exponent) + { if (u.ieee.exponent > FLT128_MANT_DIG) u.ieee.exponent -= FLT128_MANT_DIG; } @@ -175,6 +184,12 @@ fmaq (__float128 x, __float128 y, __float128 z) if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0)) return x * y + z; +#ifdef USE_FENV_H + fenv_t env; + feholdexcept (&env); + fesetround (FE_TONEAREST); +#endif + /* Multiplication m1 + m2 = x * y using Dekker's algorithm. */ #define C ((1LL << (FLT128_MANT_DIG + 1) / 2) + 1) __float128 x1 = x * C; @@ -193,10 +208,25 @@ fmaq (__float128 x, __float128 y, __float128 z) t1 = m1 - t1; t2 = z - t2; __float128 a2 = t1 + t2; +#ifdef USE_FENV_H + feclearexcept (FE_INEXACT); +#endif + /* If the result is an exact zero, ensure it has the correct + sign. */ + if (a1 == 0 && m2 == 0) + { #ifdef USE_FENV_H - fenv_t env; - feholdexcept (&env); + feupdateenv (&env); +#endif + /* Ensure that round-to-nearest value of z + m1 is not + reused. */ + asm volatile ("" : "=m" (z) : "m" (z)); + return z + m1; + } + + +#ifdef USE_FENV_H fesetround (FE_TOWARDZERO); #endif /* Perform m2 + a2 addition with round to odd. */