From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: (qmail 3796 invoked by alias); 2 Nov 2012 02:15:06 -0000 Received: (qmail 3730 invoked by uid 22791); 2 Nov 2012 02:15:03 -0000 X-SWARE-Spam-Status: No, hits=-1.2 required=5.0 tests=AWL,BAYES_00,KAM_STOCKTIP,KHOP_RCVD_UNTRUST,RCVD_IN_HOSTKARMA_W,RCVD_IN_HOSTKARMA_WL,TW_MX,TW_XC X-Spam-Check-By: sourceware.org Received: from relay1.mentorg.com (HELO relay1.mentorg.com) (192.94.38.131) by sourceware.org (qpsmtpd/0.43rc1) with ESMTP; Fri, 02 Nov 2012 02:14:54 +0000 Received: from svr-orw-exc-10.mgc.mentorg.com ([147.34.98.58]) by relay1.mentorg.com with esmtp id 1TU6ma-0006cA-KO from joseph_myers@mentor.com ; Thu, 01 Nov 2012 19:14:52 -0700 Received: from SVR-IES-FEM-01.mgc.mentorg.com ([137.202.0.104]) by SVR-ORW-EXC-10.mgc.mentorg.com with Microsoft SMTPSVC(6.0.3790.4675); Thu, 1 Nov 2012 19:14:51 -0700 Received: from digraph.polyomino.org.uk (137.202.0.76) by SVR-IES-FEM-01.mgc.mentorg.com (137.202.0.104) with Microsoft SMTP Server id 14.1.289.1; Fri, 2 Nov 2012 02:14:49 +0000 Received: from jsm28 (helo=localhost) by digraph.polyomino.org.uk with local-esmtp (Exim 4.76) (envelope-from ) id 1TU6mX-0005qq-0Y; Fri, 02 Nov 2012 02:14:49 +0000 Date: Fri, 02 Nov 2012 02:15:00 -0000 From: "Joseph S. Myers" To: CC: David Miller , Subject: Make fma use of Dekker and Knuth algorithms use round-to-nearest (bug 14796) Message-ID: MIME-Version: 1.0 Content-Type: text/plain; charset="US-ASCII" Mailing-List: contact libc-ports-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Subscribe: List-Post: List-Help: , Sender: libc-ports-owner@sourceware.org X-SW-Source: 2012-11/txt/msg00002.txt.bz2 [This patch includes architecture changes, see below; most are pretty mechanical, but the SPARC changes could perhaps do with review / testing.] Another instance of fma code assuming round-to-nearest is bug 14796, the use of Dekker multiplication and Knuth addition. Dekker multiplication definitely only works in round-to-nearest in the case of an odd number of mantissa bits (in other rounding modes, some of the intermediate multplications can end up multiplying two numbers with (M+1)/2 bits, and so being inexact). I don't have testcases for it not working in the case of an even number of mantissa bits for ldbl-96 - but in any case, Knuth addition cannot work reliably in modes other than round-to-nearest (because if x and y are sufficiently different in size, the difference between the exact value of x+y and a not-to-nearest rounded value may not itself be exactly representable). And while the proofs given by Boldo and Melquiond apply more or less unchanged to the case of the desired rounding for the final result not being round-to-nearest, they certainly don't allow for the addition step not forming an exact low part. So this patch makes all relevant fma implementations set round-to-nearest for the Dekker and Knuth algorithms (even in the absence of testcases showing wrong results from the present ldbl-96 code, it doesn't seem wise to go outside the conditions of the proofs). In turn this requires other changes: * Exact zero returns from nonzero arguments previously managed to get the correct sign of zero through the Knuth addition being carried out in the rounding mode of the final computation. Now that addition is carried out in round-to-nearest, getting -0 from x - x in round-downward mode requires new special cases. Such exact zero results can only occur when the multiplication is exact (so m2 = 0) and m1 = -z (so a1 = 0, which implies a2 = 0). Thus appropriate checks are introduced. * Inexact exceptions from the Dekker/Knuth code must be cleared before the round-to-odd code (previously the feholdexcept call was located in a position to do that); clearing them is a good idea anyway, since they are irrelevant to whether the final fma is exact, and as a fully specified IEEE arithmetic operation, the fma function should be raising inexact exceptins if and only if the final fma is inexact. To avoid check-localplt failures, this requires a hidden alias for feclearexcept, so requiring changes to each file implementing that function. * As the dbl-64 code will on x86_64 only be setting the SSE rounding mode with libc_feholdexcept_setround, and only restoring that mode with libc_feupdateenv, a corresponding libc_fesetround is needed, now we have a function that sets the rounding mode twice to two different modes before restoring it for the final calculation and return. This results in the fenv_private / math_private.h changes. Tested x86_64 (multi-arch both enabled and disabled) and x86. 2012-11-02 Joseph Myers [BZ #14796] * sysdeps/ieee754/dbl-64/s_fma.c (__fma): Set rounding mode to FE_TONEAREST before applying Dekker multiplication and Knuth addition. Clear inexact exceptions and check for exact zero results afterwards. * sysdeps/ieee754/ldbl-128/s_fmal.c (__fmal): Likewise. * sysdeps/ieee754/ldbl-96/s_fma.c (__fma): Likewise. * sysdeps/ieee754/ldbl-96/s_fmal.c (__fmal): Likewise. * math/libm-test.inc (fma_test): Add more tests. (fma_test_towardzero): Likewise. (fma_test_downward): Likewise. (fma_test_upward): Likewise. * sysdeps/generic/math_private.h (default_libc_fesetround): New function. (libc_fesetround): New macro. (libc_fesetroundf): Likewise. (libc_fesetroundl): Likewise. * sysdeps/i386/fpu/fenv_private.h (libc_fesetround_sse): New function. (libc_fesetround_387): Likewise. (libc_fesetroundf): New macro. (libc_fesetround): Likewise. (libc_fesetroundl): Likewise. * sysdeps/sparc/fpu/fenv_private.h (libc_fesetround): New function. (libc_fesetroundf): New macro. (libc_fesetround): Likewise. (libc_fesetroundl): Likewise. * include/fenv.h (feclearexcept): Add libm_hidden_proto. * math/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. * sysdeps/i386/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. * sysdeps/powerpc/fpu/fclrexcpt.c (feclearexcept): Likewise. * sysdeps/s390/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_def. * sysdeps/sh/sh4/fpu/fclrexcpt.c (feclearexcept): Likewise. * sysdeps/sparc/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. * sysdeps/x86_64/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_def. ports/ChangeLog.alpha 2012-11-02 Joseph Myers * sysdeps/alpha/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. ports/ChangeLog.am33 2012-11-02 Joseph Myers * sysdeps/am33/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. ports/ChangeLog.arm 2012-11-02 Joseph Myers * sysdeps/arm/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. ports/ChangeLog.hppa 2012-11-02 Joseph Myers * sysdeps/hppa/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_def. ports/ChangeLog.ia64 2012-11-02 Joseph Myers * sysdeps/ia64/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_def. ports/ChangeLog.m68k 2012-11-02 Joseph Myers * sysdeps/m68k/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. ports/ChangeLog.mips 2012-11-02 Joseph Myers * sysdeps/mips/fpu/fclrexcpt.c (feclearexcept): Add libm_hidden_def. ports/ChangeLog.powerpc 2012-11-02 Joseph Myers * sysdeps/powerpc/nofpu/fclrexcpt.c (feclearexcept): Add libm_hidden_ver. diff --git a/include/fenv.h b/include/fenv.h index 59d4c3f..8e36fdb 100644 --- a/include/fenv.h +++ b/include/fenv.h @@ -19,5 +19,6 @@ libm_hidden_proto (fesetround) libm_hidden_proto (feholdexcept) libm_hidden_proto (feupdateenv) libm_hidden_proto (fetestexcept) +libm_hidden_proto (feclearexcept) #endif diff --git a/math/fclrexcpt.c b/math/fclrexcpt.c index dcdcfbb..bbf1011 100644 --- a/math/fclrexcpt.c +++ b/math/fclrexcpt.c @@ -30,6 +30,7 @@ __feclearexcept (int excepts) strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); stub_warning (feclearexcept) diff --git a/math/libm-test.inc b/math/libm-test.inc index 258c960..48241e0 100644 --- a/math/libm-test.inc +++ b/math/libm-test.inc @@ -4649,6 +4649,10 @@ fma_test (void) TEST_fff_f (fma, 0x1p-149, -0x1p-149, 0x1p-149, 0x1p-149, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, 0x1p-149, -0x1p-149, -0x1p-149, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, -0x1p-149, -0x1p-149, -0x1p-149, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffp0, 0x0.fffp0, -0x0.ffep0, 0x1p-24); + TEST_fff_f (fma, 0x0.fffp0, -0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, 0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, -0x0.fffp0, -0x0.ffep0, 0x1p-24); #endif #if defined (TEST_DOUBLE) && DBL_MANT_DIG == 53 TEST_fff_f (fma, 0x1.7fp+13, 0x1.0000000000001p+0, 0x1.ffep-48, 0x1.7f00000000001p+13); @@ -4695,6 +4699,10 @@ fma_test (void) TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, 0x1p-1074, 0x1p-1074, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, 0x1p-1074, -0x1p-1074, -0x1p-1074, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, -0x1p-1074, -0x1p-1074, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 64 TEST_fff_f (fma, -0x8.03fcp+3696L, 0xf.fffffffffffffffp-6140L, 0x8.3ffffffffffffffp-2450L, -0x8.01ecp-2440L); @@ -4727,6 +4735,10 @@ fma_test (void) TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, 0x1p-16445L, 0x1p-16445L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, -0x1p-16445L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 113 TEST_fff_f (fma, 0x1.bb2de33e02ccbbfa6e245a7c1f71p-2584L, -0x1.6b500daf0580d987f1bc0cadfcddp-13777L, 0x1.613cd91d9fed34b33820e5ab9d8dp-16378L, -0x1.3a79fb50eb9ce887cffa0f09bd9fp-16360L); @@ -4766,6 +4778,10 @@ fma_test (void) TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, 0x1p-16494L, 0x1p-16494L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, -0x1p-16494L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); #endif END (fma); @@ -4846,6 +4862,10 @@ fma_test_towardzero (void) TEST_fff_f (fma, 0x1p-149, -0x1p-149, 0x1p-149, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, 0x1p-149, -0x1p-149, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, -0x1p-149, -0x1p-149, -0x1p-149, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffp0, 0x0.fffp0, -0x0.ffep0, 0x1p-24); + TEST_fff_f (fma, 0x0.fffp0, -0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, 0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, -0x0.fffp0, -0x0.ffep0, 0x1p-24); #endif #if defined (TEST_DOUBLE) && DBL_MANT_DIG == 53 TEST_fff_f (fma, 0x1.4p-1022, 0x1.0000000000002p-1, 0x1p-1024, 0x1.c000000000002p-1023, UNDERFLOW_EXCEPTION); @@ -4872,6 +4892,10 @@ fma_test_towardzero (void) TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, 0x1p-1074, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, 0x1p-1074, -0x1p-1074, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, -0x1p-1074, -0x1p-1074, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 64 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000004p-1L, 0x1p-16384L, 0x1.c000000000000004p-16383L, UNDERFLOW_EXCEPTION); @@ -4898,6 +4922,10 @@ fma_test_towardzero (void) TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, 0x1p-16445L, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, 0x1p-16445L, -0x1p-16445L, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, -0x1p-16445L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 113 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000000000000000002p-1L, 0x1p-16384L, 0x1.c000000000000000000000000002p-16383L, UNDERFLOW_EXCEPTION); @@ -4924,6 +4952,10 @@ fma_test_towardzero (void) TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, 0x1p-16494L, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, 0x1p-16494L, -0x1p-16494L, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, -0x1p-16494L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); #endif } @@ -5007,6 +5039,10 @@ fma_test_downward (void) TEST_fff_f (fma, 0x1p-149, -0x1p-149, 0x1p-149, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, 0x1p-149, -0x1p-149, -0x1p-149, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, -0x1p-149, -0x1p-149, -0x1p-148, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffp0, 0x0.fffp0, -0x0.ffep0, 0x1p-24); + TEST_fff_f (fma, 0x0.fffp0, -0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, 0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, -0x0.fffp0, -0x0.ffep0, 0x1p-24); #endif #if defined (TEST_DOUBLE) && DBL_MANT_DIG == 53 TEST_fff_f (fma, 0x1.4p-1022, 0x1.0000000000002p-1, 0x1p-1024, 0x1.c000000000002p-1023, UNDERFLOW_EXCEPTION); @@ -5033,6 +5069,10 @@ fma_test_downward (void) TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, 0x1p-1074, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, 0x1p-1074, -0x1p-1074, -0x1p-1074, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, -0x1p-1074, -0x1p-1073, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 64 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000004p-1L, 0x1p-16384L, 0x1.c000000000000004p-16383L, UNDERFLOW_EXCEPTION); @@ -5059,6 +5099,10 @@ fma_test_downward (void) TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, 0x1p-16445L, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, -0x1p-16444L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 113 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000000000000000002p-1L, 0x1p-16384L, 0x1.c000000000000000000000000002p-16383L, UNDERFLOW_EXCEPTION); @@ -5085,6 +5129,10 @@ fma_test_downward (void) TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, 0x1p-16494L, plus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, -0x1p-16493L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); #endif } @@ -5168,6 +5216,10 @@ fma_test_upward (void) TEST_fff_f (fma, 0x1p-149, -0x1p-149, 0x1p-149, 0x1p-149, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, 0x1p-149, -0x1p-149, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-149, -0x1p-149, -0x1p-149, -0x1p-149, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffp0, 0x0.fffp0, -0x0.ffep0, 0x1p-24); + TEST_fff_f (fma, 0x0.fffp0, -0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, 0x0.fffp0, 0x0.ffep0, -0x1p-24); + TEST_fff_f (fma, -0x0.fffp0, -0x0.fffp0, -0x0.ffep0, 0x1p-24); #endif #if defined (TEST_DOUBLE) && DBL_MANT_DIG == 53 TEST_fff_f (fma, 0x1.4p-1022, 0x1.0000000000002p-1, 0x1p-1024, 0x1.c000000000004p-1023, UNDERFLOW_EXCEPTION); @@ -5194,6 +5246,10 @@ fma_test_upward (void) TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, 0x1p-1074, 0x1p-1074, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, 0x1p-1074, -0x1p-1074, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-1074, -0x1p-1074, -0x1p-1074, -0x1p-1074, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); + TEST_fff_f (fma, 0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0, 0x0.fffffffffffffp0, -0x1p-106); + TEST_fff_f (fma, -0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0, -0x0.fffffffffffffp0, 0x1p-106); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 64 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000004p-1L, 0x1p-16384L, 0x1.c000000000000008p-16383L, UNDERFLOW_EXCEPTION); @@ -5220,6 +5276,10 @@ fma_test_upward (void) TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, 0x1p-16445L, 0x1p-16445L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, 0x1p-16445L, -0x1p-16445L, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16445L, -0x1p-16445L, -0x1p-16445L, -0x1p-16445L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); + TEST_fff_f (fma, 0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, 0x0.ffffffffffffffffp0L, 0x0.fffffffffffffffep0L, -0x1p-128L); + TEST_fff_f (fma, -0x0.ffffffffffffffffp0L, -0x0.ffffffffffffffffp0L, -0x0.fffffffffffffffep0L, 0x1p-128L); #endif #if defined (TEST_LDOUBLE) && LDBL_MANT_DIG == 113 TEST_fff_f (fma, 0x1.4p-16382L, 0x1.0000000000000000000000000002p-1L, 0x1p-16384L, 0x1.c000000000000000000000000004p-16383L, UNDERFLOW_EXCEPTION); @@ -5246,6 +5306,10 @@ fma_test_upward (void) TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, 0x1p-16494L, 0x1p-16494L, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, 0x1p-16494L, -0x1p-16494L, minus_zero, UNDERFLOW_EXCEPTION); TEST_fff_f (fma, 0x1p-16494L, -0x1p-16494L, -0x1p-16494L, -0x1p-16494L, UNDERFLOW_EXCEPTION); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); + TEST_fff_f (fma, 0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffffffffffffffffffffffffp0L, -0x1p-226L); + TEST_fff_f (fma, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffffffffffffffffffffffffp0L, 0x1p-226L); #endif } diff --git a/ports/sysdeps/alpha/fpu/fclrexcpt.c b/ports/sysdeps/alpha/fpu/fclrexcpt.c index 2d2bd94..8459652 100644 --- a/ports/sysdeps/alpha/fpu/fclrexcpt.c +++ b/ports/sysdeps/alpha/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997,99,2000,01 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Richard Henderson , 1997. @@ -43,4 +43,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/ports/sysdeps/am33/fpu/fclrexcpt.c b/ports/sysdeps/am33/fpu/fclrexcpt.c index 2b15f45..492ea38 100644 --- a/ports/sysdeps/am33/fpu/fclrexcpt.c +++ b/ports/sysdeps/am33/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1998, 1999, 2000, 2002, 2004 Free Software Foundation, Inc. + Copyright (C) 1998-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Alexandre Oliva based on corresponding file in the MIPS port. @@ -48,4 +48,5 @@ __feclearexcept (int excepts) return 0; } +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/ports/sysdeps/arm/fclrexcpt.c b/ports/sysdeps/arm/fclrexcpt.c index 23435fb..fb106a5 100644 --- a/ports/sysdeps/arm/fclrexcpt.c +++ b/ports/sysdeps/arm/fclrexcpt.c @@ -54,4 +54,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/ports/sysdeps/hppa/fpu/fclrexcpt.c b/ports/sysdeps/hppa/fpu/fclrexcpt.c index 5d1e593..bbc2c70 100644 --- a/ports/sysdeps/hppa/fpu/fclrexcpt.c +++ b/ports/sysdeps/hppa/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 2000 Free Software Foundation, Inc. + Copyright (C) 2000-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by David Huggins-Daines , 2000 @@ -33,3 +33,4 @@ feclearexcept (int excepts) /* Success. */ return 0; } +libm_hidden_def (feclearexcept) diff --git a/ports/sysdeps/ia64/fpu/fclrexcpt.c b/ports/sysdeps/ia64/fpu/fclrexcpt.c index 84f8327..c099fd3 100644 --- a/ports/sysdeps/ia64/fpu/fclrexcpt.c +++ b/ports/sysdeps/ia64/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997, 1999, 2000 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Christian Boissat , 1999 and Jes Sorensen , 2000 @@ -36,3 +36,4 @@ feclearexcept (int excepts) /* success */ return 0; } +libm_hidden_def (feclearexcept) diff --git a/ports/sysdeps/m68k/fpu/fclrexcpt.c b/ports/sysdeps/m68k/fpu/fclrexcpt.c index ceda99c..cfc8d9e 100644 --- a/ports/sysdeps/m68k/fpu/fclrexcpt.c +++ b/ports/sysdeps/m68k/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997,99,2000,01 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Andreas Schwab @@ -46,4 +46,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/ports/sysdeps/mips/fpu/fclrexcpt.c b/ports/sysdeps/mips/fpu/fclrexcpt.c index f97d892..f4709b4 100644 --- a/ports/sysdeps/mips/fpu/fclrexcpt.c +++ b/ports/sysdeps/mips/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1998, 1999, 2000, 2002 Free Software Foundation, Inc. + Copyright (C) 1998-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Andreas Jaeger , 1998. @@ -44,3 +44,4 @@ feclearexcept (int excepts) /* Success. */ return 0; } +libm_hidden_def (feclearexcept) diff --git a/ports/sysdeps/powerpc/nofpu/fclrexcpt.c b/ports/sysdeps/powerpc/nofpu/fclrexcpt.c index 768fd8f..f4b9016 100644 --- a/ports/sysdeps/powerpc/nofpu/fclrexcpt.c +++ b/ports/sysdeps/powerpc/nofpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear floating-point exceptions (soft-float edition). - Copyright (C) 2002 Free Software Foundation, Inc. + Copyright (C) 2002-2012 Free Software Foundation, Inc. Contributed by Aldy Hernandez , 2002. This file is part of the GNU C Library. @@ -33,4 +33,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/sysdeps/generic/math_private.h b/sysdeps/generic/math_private.h index b375bc0..7661788 100644 --- a/sysdeps/generic/math_private.h +++ b/sysdeps/generic/math_private.h @@ -402,6 +402,22 @@ default_libc_feholdexcept (fenv_t *e) #endif static __always_inline void +default_libc_fesetround (int r) +{ + (void) fesetround (r); +} + +#ifndef libc_fesetround +# define libc_fesetround default_libc_fesetround +#endif +#ifndef libc_fesetroundf +# define libc_fesetroundf default_libc_fesetround +#endif +#ifndef libc_fesetroundl +# define libc_fesetroundl default_libc_fesetround +#endif + +static __always_inline void default_libc_feholdexcept_setround (fenv_t *e, int r) { feholdexcept (e); diff --git a/sysdeps/i386/fpu/fclrexcpt.c b/sysdeps/i386/fpu/fclrexcpt.c index f24d07f..f28ef6b 100644 --- a/sysdeps/i386/fpu/fclrexcpt.c +++ b/sysdeps/i386/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997,99,2000, 2001, 2003, 2004 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper , 1997. @@ -65,4 +65,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/sysdeps/i386/fpu/fenv_private.h b/sysdeps/i386/fpu/fenv_private.h index f33f57c..03f4c97 100644 --- a/sysdeps/i386/fpu/fenv_private.h +++ b/sysdeps/i386/fpu/fenv_private.h @@ -77,6 +77,24 @@ libc_feholdexcept_387 (fenv_t *e) } static __always_inline void +libc_fesetround_sse (int r) +{ + unsigned int mxcsr; + asm (STMXCSR " %0" : "=m" (*&mxcsr)); + mxcsr = (mxcsr & ~0x6000) | (r << 3); + asm volatile (LDMXCSR " %0" : : "m" (*&mxcsr)); +} + +static __always_inline void +libc_fesetround_387 (int r) +{ + fpu_control_t cw; + _FPU_GETCW (cw); + cw = (cw & ~0xc00) | r; + _FPU_SETCW (cw); +} + +static __always_inline void libc_feholdexcept_setround_sse (fenv_t *e, int r) { unsigned int mxcsr; @@ -247,6 +265,7 @@ libc_feresetround_387 (fenv_t *e) #ifdef __SSE_MATH__ # define libc_feholdexceptf libc_feholdexcept_sse +# define libc_fesetroundf libc_fesetround_sse # define libc_feholdexcept_setroundf libc_feholdexcept_setround_sse # define libc_fetestexceptf libc_fetestexcept_sse # define libc_fesetenvf libc_fesetenv_sse @@ -256,6 +275,7 @@ libc_feresetround_387 (fenv_t *e) # define libc_feresetroundf libc_feresetround_sse #else # define libc_feholdexceptf libc_feholdexcept_387 +# define libc_fesetroundf libc_fesetround_387 # define libc_feholdexcept_setroundf libc_feholdexcept_setround_387 # define libc_fetestexceptf libc_fetestexcept_387 # define libc_fesetenvf libc_fesetenv_387 @@ -267,6 +287,7 @@ libc_feresetround_387 (fenv_t *e) #ifdef __SSE2_MATH__ # define libc_feholdexcept libc_feholdexcept_sse +# define libc_fesetround libc_fesetround_sse # define libc_feholdexcept_setround libc_feholdexcept_setround_sse # define libc_fetestexcept libc_fetestexcept_sse # define libc_fesetenv libc_fesetenv_sse @@ -276,6 +297,7 @@ libc_feresetround_387 (fenv_t *e) # define libc_feresetround libc_feresetround_sse #else # define libc_feholdexcept libc_feholdexcept_387 +# define libc_fesetround libc_fesetround_387 # define libc_feholdexcept_setround libc_feholdexcept_setround_387 # define libc_fetestexcept libc_fetestexcept_387 # define libc_fesetenv libc_fesetenv_387 @@ -286,6 +308,7 @@ libc_feresetround_387 (fenv_t *e) #endif /* __SSE2_MATH__ */ #define libc_feholdexceptl libc_feholdexcept_387 +#define libc_fesetroundl libc_fesetround_387 #define libc_feholdexcept_setroundl libc_feholdexcept_setround_387 #define libc_fetestexceptl libc_fetestexcept_387 #define libc_fesetenvl libc_fesetenv_387 diff --git a/sysdeps/ieee754/dbl-64/s_fma.c b/sysdeps/ieee754/dbl-64/s_fma.c index 68c63a1..07fe715 100644 --- a/sysdeps/ieee754/dbl-64/s_fma.c +++ b/sysdeps/ieee754/dbl-64/s_fma.c @@ -167,6 +167,9 @@ __fma (double x, double y, double z) if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0)) return x * y + z; + fenv_t env; + libc_feholdexcept_setround (&env, FE_TONEAREST); + /* Multiplication m1 + m2 = x * y using Dekker's algorithm. */ #define C ((1 << (DBL_MANT_DIG + 1) / 2) + 1) double x1 = x * C; @@ -185,9 +188,20 @@ __fma (double x, double y, double z) t1 = m1 - t1; t2 = z - t2; double a2 = t1 + t2; + feclearexcept (FE_INEXACT); - fenv_t env; - libc_feholdexcept_setround (&env, FE_TOWARDZERO); + /* If the result is an exact zero, ensure it has the correct + sign. */ + if (a1 == 0 && m2 == 0) + { + libc_feupdateenv (&env); + /* Ensure that round-to-nearest value of z + m1 is not + reused. */ + asm volatile ("" : "=m" (z) : "m" (z)); + return z + m1; + } + + libc_fesetround (FE_TOWARDZERO); /* Perform m2 + a2 addition with round to odd. */ u.d = a2 + m2; diff --git a/sysdeps/ieee754/ldbl-128/s_fmal.c b/sysdeps/ieee754/ldbl-128/s_fmal.c index c6a3d71..d576403 100644 --- a/sysdeps/ieee754/ldbl-128/s_fmal.c +++ b/sysdeps/ieee754/ldbl-128/s_fmal.c @@ -170,6 +170,10 @@ __fmal (long double x, long double y, long double z) if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0)) return x * y + z; + fenv_t env; + feholdexcept (&env); + fesetround (FE_TONEAREST); + /* Multiplication m1 + m2 = x * y using Dekker's algorithm. */ #define C ((1LL << (LDBL_MANT_DIG + 1) / 2) + 1) long double x1 = x * C; @@ -188,9 +192,19 @@ __fmal (long double x, long double y, long double z) t1 = m1 - t1; t2 = z - t2; long double a2 = t1 + t2; + feclearexcept (FE_INEXACT); + + /* If the result is an exact zero, ensure it has the correct + sign. */ + if (a1 == 0 && m2 == 0) + { + feupdateenv (&env); + /* Ensure that round-to-nearest value of z + m1 is not + reused. */ + asm volatile ("" : "=m" (z) : "m" (z)); + return z + m1; + } - fenv_t env; - feholdexcept (&env); fesetround (FE_TOWARDZERO); /* Perform m2 + a2 addition with round to odd. */ u.d = a2 + m2; diff --git a/sysdeps/ieee754/ldbl-96/s_fma.c b/sysdeps/ieee754/ldbl-96/s_fma.c index 001d806..bf2d794 100644 --- a/sysdeps/ieee754/ldbl-96/s_fma.c +++ b/sysdeps/ieee754/ldbl-96/s_fma.c @@ -42,6 +42,10 @@ __fma (double x, double y, double z) if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0)) return x * y + z; + fenv_t env; + feholdexcept (&env); + fesetround (FE_TONEAREST); + /* Multiplication m1 + m2 = x * y using Dekker's algorithm. */ #define C ((1ULL << (LDBL_MANT_DIG + 1) / 2) + 1) long double x1 = (long double) x * C; @@ -60,9 +64,19 @@ __fma (double x, double y, double z) t1 = m1 - t1; t2 = z - t2; long double a2 = t1 + t2; + feclearexcept (FE_INEXACT); + + /* If the result is an exact zero, ensure it has the correct + sign. */ + if (a1 == 0 && m2 == 0) + { + feupdateenv (&env); + /* Ensure that round-to-nearest value of z + m1 is not + reused. */ + asm volatile ("" : "=m" (z) : "m" (z)); + return z + m1; + } - fenv_t env; - feholdexcept (&env); fesetround (FE_TOWARDZERO); /* Perform m2 + a2 addition with round to odd. */ a2 = a2 + m2; diff --git a/sysdeps/ieee754/ldbl-96/s_fmal.c b/sysdeps/ieee754/ldbl-96/s_fmal.c index b5fdcba..32e71a1 100644 --- a/sysdeps/ieee754/ldbl-96/s_fmal.c +++ b/sysdeps/ieee754/ldbl-96/s_fmal.c @@ -168,6 +168,10 @@ __fmal (long double x, long double y, long double z) if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0)) return x * y + z; + fenv_t env; + feholdexcept (&env); + fesetround (FE_TONEAREST); + /* Multiplication m1 + m2 = x * y using Dekker's algorithm. */ #define C ((1LL << (LDBL_MANT_DIG + 1) / 2) + 1) long double x1 = x * C; @@ -186,9 +190,19 @@ __fmal (long double x, long double y, long double z) t1 = m1 - t1; t2 = z - t2; long double a2 = t1 + t2; + feclearexcept (FE_INEXACT); + + /* If the result is an exact zero, ensure it has the correct + sign. */ + if (a1 == 0 && m2 == 0) + { + feupdateenv (&env); + /* Ensure that round-to-nearest value of z + m1 is not + reused. */ + asm volatile ("" : "=m" (z) : "m" (z)); + return z + m1; + } - fenv_t env; - feholdexcept (&env); fesetround (FE_TOWARDZERO); /* Perform m2 + a2 addition with round to odd. */ u.d = a2 + m2; diff --git a/sysdeps/powerpc/fpu/fclrexcpt.c b/sysdeps/powerpc/fpu/fclrexcpt.c index 87376bf..e99cc58 100644 --- a/sysdeps/powerpc/fpu/fclrexcpt.c +++ b/sysdeps/powerpc/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997,99,2000,01 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -44,4 +44,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/sysdeps/s390/fpu/fclrexcpt.c b/sysdeps/s390/fpu/fclrexcpt.c index 2352d74..3e8d9bb 100644 --- a/sysdeps/s390/fpu/fclrexcpt.c +++ b/sysdeps/s390/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 2000 Free Software Foundation, Inc. + Copyright (C) 2000-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -37,3 +37,4 @@ feclearexcept (int excepts) /* Success. */ return 0; } +libm_hidden_def (feclearexcept) diff --git a/sysdeps/sh/sh4/fpu/fclrexcpt.c b/sysdeps/sh/sh4/fpu/fclrexcpt.c index b4b2ead..2d31fa6 100644 --- a/sysdeps/sh/sh4/fpu/fclrexcpt.c +++ b/sysdeps/sh/sh4/fpu/fclrexcpt.c @@ -39,3 +39,4 @@ feclearexcept (int excepts) return 0; } +libm_hidden_def (feclearexcept) diff --git a/sysdeps/sparc/fpu/fclrexcpt.c b/sysdeps/sparc/fpu/fclrexcpt.c index 82b6ac4..50ec4f4 100644 --- a/sysdeps/sparc/fpu/fclrexcpt.c +++ b/sysdeps/sparc/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 1997, 1999, 2000 Free Software Foundation, Inc. + Copyright (C) 1997-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -39,4 +39,5 @@ strong_alias (__feclearexcept, __old_feclearexcept) compat_symbol (libm, __old_feclearexcept, feclearexcept, GLIBC_2_1); #endif +libm_hidden_ver (__feclearexcept, feclearexcept) versioned_symbol (libm, __feclearexcept, feclearexcept, GLIBC_2_2); diff --git a/sysdeps/sparc/fpu/fenv_private.h b/sysdeps/sparc/fpu/fenv_private.h index a6e8e95..8690879 100644 --- a/sysdeps/sparc/fpu/fenv_private.h +++ b/sysdeps/sparc/fpu/fenv_private.h @@ -14,6 +14,15 @@ libc_feholdexcept (fenv_t *e) } static __always_inline void +libc_fesetround (int r) +{ + fenv_t etmp; + __fenv_stfsr(etmp); + etmp = (etmp & ~__FE_ROUND_MASK) | (r); + __fenv_ldfsr(etmp); +} + +static __always_inline void libc_feholdexcept_setround (fenv_t *e, int r) { fenv_t etmp; @@ -79,6 +88,7 @@ libc_feresetround (fenv_t *e) } #define libc_feholdexceptf libc_feholdexcept +#define libc_fesetroundf libc_fesetround #define libc_feholdexcept_setroundf libc_feholdexcept_setround #define libc_fetestexceptf libc_fetestexcept #define libc_fesetenvf libc_fesetenv @@ -87,6 +97,7 @@ libc_feresetround (fenv_t *e) #define libc_feholdsetroundf libc_feholdsetround #define libc_feresetroundf libc_feresetround #define libc_feholdexcept libc_feholdexcept +#define libc_fesetround libc_fesetround #define libc_feholdexcept_setround libc_feholdexcept_setround #define libc_fetestexcept libc_fetestexcept #define libc_fesetenv libc_fesetenv @@ -95,6 +106,7 @@ libc_feresetround (fenv_t *e) #define libc_feholdsetround libc_feholdsetround #define libc_feresetround libc_feresetround #define libc_feholdexceptl libc_feholdexcept +#define libc_fesetroundl libc_fesetround #define libc_feholdexcept_setroundl libc_feholdexcept_setround #define libc_fetestexceptl libc_fetestexcept #define libc_fesetenvl libc_fesetenv diff --git a/sysdeps/x86_64/fpu/fclrexcpt.c b/sysdeps/x86_64/fpu/fclrexcpt.c index 5ef3162..2d40313 100644 --- a/sysdeps/x86_64/fpu/fclrexcpt.c +++ b/sysdeps/x86_64/fpu/fclrexcpt.c @@ -1,5 +1,5 @@ /* Clear given exceptions in current floating-point environment. - Copyright (C) 2001 Free Software Foundation, Inc. + Copyright (C) 2001-2012 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -49,3 +49,4 @@ feclearexcept (int excepts) /* Success. */ return 0; } +libm_hidden_def (feclearexcept) -- Joseph S. Myers joseph@codesourcery.com