From patchwork Wed Aug 4 09:46:09 2021 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: Paul Zimmermann X-Patchwork-Id: 44569 Return-Path: X-Original-To: patchwork@sourceware.org Delivered-To: patchwork@sourceware.org Received: from server2.sourceware.org (localhost [IPv6:::1]) by sourceware.org (Postfix) with ESMTP id 8D7D93972037 for ; Wed, 4 Aug 2021 09:46:24 +0000 (GMT) X-Original-To: libc-alpha@sourceware.org Delivered-To: libc-alpha@sourceware.org Received: from mail3-relais-sop.national.inria.fr (mail3-relais-sop.national.inria.fr [192.134.164.104]) by sourceware.org (Postfix) with ESMTPS id 4A36838515F4 for ; Wed, 4 Aug 2021 09:46:12 +0000 (GMT) DMARC-Filter: OpenDMARC Filter v1.4.1 sourceware.org 4A36838515F4 Authentication-Results: sourceware.org; dmarc=none (p=none dis=none) header.from=inria.fr Authentication-Results: sourceware.org; spf=pass smtp.mailfrom=loria.fr IronPort-HdrOrdr: A9a23:vnez6qDRDoIhVyLlHem955DYdb4zR+YMi2TDGXocdfUzSL37qynAppomPHPP4gr5O0tQ+uxoRpPgfZq0z/ccirX5Eo3SOTUO01HGEGgN1+bfKkXbexHDyg== X-IronPort-AV: E=Sophos;i="5.84,293,1620684000"; d="scan'208";a="389661031" Received: from tomate.loria.fr ([152.81.10.51]) by mail3-relais-sop.national.inria.fr with ESMTP/TLS/DHE-RSA-AES256-GCM-SHA384; 04 Aug 2021 11:46:10 +0200 Received: from zimmerma by tomate.loria.fr with local (Exim 4.94.2) (envelope-from ) id 1mBDTV-0011wA-VK; Wed, 04 Aug 2021 11:46:10 +0200 From: Paul Zimmermann To: libc-alpha@sourceware.org Subject: [PATCH] Fixed inaccuracy of j0f (BZ #28185) Date: Wed, 4 Aug 2021 11:46:09 +0200 Message-Id: <20210804094609.245719-1-Paul.Zimmermann@inria.fr> X-Mailer: git-send-email 2.30.2 MIME-Version: 1.0 X-Spam-Status: No, score=-9.9 required=5.0 tests=BAYES_00, GIT_PATCH_0, HEADER_FROM_DIFFERENT_DOMAINS, KAM_DMARC_STATUS, RCVD_IN_MSPIKE_H3, RCVD_IN_MSPIKE_WL, SPF_HELO_NONE, SPF_PASS, TXREP autolearn=ham autolearn_force=no version=3.4.4 X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on server2.sourceware.org X-BeenThere: libc-alpha@sourceware.org X-Mailman-Version: 2.1.29 Precedence: list List-Id: Libc-alpha mailing list List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , Errors-To: libc-alpha-bounces+patchwork=sourceware.org@sourceware.org Sender: "Libc-alpha" The largest errors over the full binary32 range are after this patch (on x86_64): RNDN: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6 RNDZ: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6 RNDU: libm wrong by up to 9.00e+00 ulp(s) [9] for x=0x1.04c39cp+6 RNDD: libm wrong by up to 8.98e+00 ulp(s) [9] for x=0x1.4b7066p+7 Some constants used in the code are also changed to "static const". --- sysdeps/ieee754/flt-32/e_j0f.c | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/sysdeps/ieee754/flt-32/e_j0f.c b/sysdeps/ieee754/flt-32/e_j0f.c index 462518c365..080395eab5 100644 --- a/sysdeps/ieee754/flt-32/e_j0f.c +++ b/sysdeps/ieee754/flt-32/e_j0f.c @@ -40,7 +40,7 @@ S04 = 1.1661400734e-09; /* 0x30a045e8 */ static const float zero = 0.0; /* This is the nearest approximation of the first zero of j0. */ -#define FIRST_ZERO_J0 0xf.26247p-28f +#define FIRST_ZERO_J0 0x1.33d152p+1f #define SMALL_SIZE 64 @@ -212,7 +212,7 @@ j0f_asympt (float x) /* Now x - pi/4 - alpha0 = h + n*pi/2 mod (2*pi). */ float xr = (float) h; n = n & 3; - float cst = 0xc.c422ap-4; /* sqrt(2/pi) rounded to nearest */ + static const float cst = 0xc.c422ap-4f; /* sqrt(2/pi) rounded to nearest */ float t = cst / sqrtf (x) * (float) beta0; if (n == 0) return t * __cosf (xr); @@ -286,7 +286,7 @@ __ieee754_j0f(float x) /* The following threshold is optimal: for x=0x1.3b58dep+1 and rounding upwards, |cc|=0x1.579b26p-4 and z is 10 ulps far from the correctly rounded value. */ - float threshold = 0x1.579b26p-4; + static const float threshold = 0x1.579b26p-4f; if (fabsf (cc) > threshold) return z; else @@ -493,7 +493,7 @@ y0f_asympt (float x) /* Now x - pi/4 - alpha0 = h + n*pi/2 mod (2*pi). */ float xr = (float) h; n = n & 3; - float cst = 0xc.c422ap-4; /* sqrt(2/pi) rounded to nearest */ + static const float cst = 0xc.c422ap-4; /* sqrt(2/pi) rounded to nearest */ float t = cst / sqrtf (x) * (float) beta0; if (n == 0) return t * __sinf (xr); @@ -584,7 +584,7 @@ __ieee754_y0f(float x) z = (invsqrtpi*ss)/sqrtf(x); /* The following threshold is optimal (determined on aarch64-linux-gnu). */ - float threshold = 0x1.be585ap-4; + static const float threshold = 0x1.be585ap-4; if (fabsf (ss) > threshold) return z; else