[PATCH 09/11] math: Use erfcf from CORE-MATH
DJ Delorie
dj@redhat.com
Thu Nov 21 01:35:41 GMT 2024
A few minor nits about comments, two about table entries. Ok if those
are addressed, even if unchanged.
Reviewed-by: DJ Delorie <dj@redhat.com>
Adhemerval Zanella <adhemerval.zanella@linaro.org> writes:
> The CORE-MATH implementation is correctly rounded (for any rounding mode)
> and shows better performance to the generic erfcf.
> diff --git a/SHARED-FILES b/SHARED-FILES
> index ccc5179f80..e6e29bcadc 100644
> --- a/SHARED-FILES
> +++ b/SHARED-FILES
> @@ -276,3 +276,7 @@ sysdeps/ieee754/flt-32/s_erff.c
> (file src/binary32/erf/erff.c in CORE-MATH)
> - The code was adapted to use glibc code style and internal
> functions to handle errno, overflow, and underflow.
> +sysdeps/ieee754/flt-32/s_erfcf.c
> + (file src/binary32/erfc/erfcf.c in CORE-MATH)
> + - The code was adapted to use glibc code style and internal
> + functions to handle errno, overflow, and underflow.
Ok.
> diff --git a/sysdeps/alpha/fpu/libm-test-ulps b/sysdeps/alpha/fpu/libm-test-ulps
> diff --git a/sysdeps/arc/fpu/libm-test-ulps b/sysdeps/arc/fpu/libm-test-ulps
> diff --git a/sysdeps/arc/nofpu/libm-test-ulps b/sysdeps/arc/nofpu/libm-test-ulps
> diff --git a/sysdeps/arm/libm-test-ulps b/sysdeps/arm/libm-test-ulps
> diff --git a/sysdeps/csky/fpu/libm-test-ulps b/sysdeps/csky/fpu/libm-test-ulps
> diff --git a/sysdeps/csky/nofpu/libm-test-ulps b/sysdeps/csky/nofpu/libm-test-ulps
> diff --git a/sysdeps/hppa/fpu/libm-test-ulps b/sysdeps/hppa/fpu/libm-test-ulps
> diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
> diff --git a/sysdeps/i386/i686/fpu/multiarch/libm-test-ulps b/sysdeps/i386/i686/fpu/multiarch/libm-test-ulps
Ok.
> diff --git a/sysdeps/ieee754/flt-32/s_erfcf.c b/sysdeps/ieee754/flt-32/s_erfcf.c
> index fdba278eff..6f8c7c6a57 100644
> --- a/sysdeps/ieee754/flt-32/s_erfcf.c
> +++ b/sysdeps/ieee754/flt-32/s_erfcf.c
> @@ -1,173 +1,185 @@
> -/* s_erfcf.c -- float version of s_erfc.c.
> - */
> +/* Correctly-rounded complementary error function for the binary32 format
>
> -/*
> - * ====================================================
> - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
> - *
> - * Developed at SunPro, a Sun Microsystems, Inc. business.
> - * Permission to use, copy, modify, and distribute this
> - * software is freely granted, provided that this notice
> - * is preserved.
> - * ====================================================
> - */
> +Copyright (c) 2023, 2024 Alexei Sibidanov.
>
> -#if defined(LIBM_SCCS) && !defined(lint)
> -static char rcsid[] = "$NetBSD: s_erff.c,v 1.4 1995/05/10 20:47:07 jtc Exp $";
> -#endif
> +This file is part of the CORE-MATH project
> +project (file src/binary32/erfc/erfcf.c revision bc385c2).
> +
> +Permission is hereby granted, free of charge, to any person obtaining a copy
> +of this software and associated documentation files (the "Software"), to deal
> +in the Software without restriction, including without limitation the rights
> +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
> +copies of the Software, and to permit persons to whom the Software is
> +furnished to do so, subject to the following conditions:
> +
> +The above copyright notice and this permission notice shall be included in all
> +copies or substantial portions of the Software.
> +
> +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
> +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
> +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
> +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
> +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
> +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
> +SOFTWARE.
> +*/
Ok.
> #include <errno.h>
> -#include <fix-int-fp-convert-zero.h>
> -#include <libm-alias-float.h>
> -#include <math-narrow-eval.h>
> #include <math.h>
> -#include <math_private.h>
> -
> +#include <stdint.h>
> +#include <libm-alias-float.h>
> +#include "math_config.h"
Ok.
> -static const float
> -tiny = 1e-30,
> -half= 5.0000000000e-01, /* 0x3F000000 */
> -one = 1.0000000000e+00, /* 0x3F800000 */
> -two = 2.0000000000e+00, /* 0x40000000 */
> - /* c = (subfloat)0.84506291151 */
> -erx = 8.4506291151e-01, /* 0x3f58560b */
> -/*
> - * Coefficients for approximation to erf on [0,0.84375]
> - */
> -pp0 = 1.2837916613e-01, /* 0x3e0375d4 */
> -pp1 = -3.2504209876e-01, /* 0xbea66beb */
> -pp2 = -2.8481749818e-02, /* 0xbce9528f */
> -pp3 = -5.7702702470e-03, /* 0xbbbd1489 */
> -pp4 = -2.3763017452e-05, /* 0xb7c756b1 */
> -qq1 = 3.9791721106e-01, /* 0x3ecbbbce */
> -qq2 = 6.5022252500e-02, /* 0x3d852a63 */
> -qq3 = 5.0813062117e-03, /* 0x3ba68116 */
> -qq4 = 1.3249473704e-04, /* 0x390aee49 */
> -qq5 = -3.9602282413e-06, /* 0xb684e21a */
> -/*
> - * Coefficients for approximation to erf in [0.84375,1.25]
> - */
> -pa0 = -2.3621185683e-03, /* 0xbb1acdc6 */
> -pa1 = 4.1485610604e-01, /* 0x3ed46805 */
> -pa2 = -3.7220788002e-01, /* 0xbebe9208 */
> -pa3 = 3.1834661961e-01, /* 0x3ea2fe54 */
> -pa4 = -1.1089469492e-01, /* 0xbde31cc2 */
> -pa5 = 3.5478305072e-02, /* 0x3d1151b3 */
> -pa6 = -2.1663755178e-03, /* 0xbb0df9c0 */
> -qa1 = 1.0642088205e-01, /* 0x3dd9f331 */
> -qa2 = 5.4039794207e-01, /* 0x3f0a5785 */
> -qa3 = 7.1828655899e-02, /* 0x3d931ae7 */
> -qa4 = 1.2617121637e-01, /* 0x3e013307 */
> -qa5 = 1.3637083583e-02, /* 0x3c5f6e13 */
> -qa6 = 1.1984500103e-02, /* 0x3c445aa3 */
> -/*
> - * Coefficients for approximation to erfc in [1.25,1/0.35]
> - */
> -ra0 = -9.8649440333e-03, /* 0xbc21a093 */
> -ra1 = -6.9385856390e-01, /* 0xbf31a0b7 */
> -ra2 = -1.0558626175e+01, /* 0xc128f022 */
> -ra3 = -6.2375331879e+01, /* 0xc2798057 */
> -ra4 = -1.6239666748e+02, /* 0xc322658c */
> -ra5 = -1.8460508728e+02, /* 0xc3389ae7 */
> -ra6 = -8.1287437439e+01, /* 0xc2a2932b */
> -ra7 = -9.8143291473e+00, /* 0xc11d077e */
> -sa1 = 1.9651271820e+01, /* 0x419d35ce */
> -sa2 = 1.3765776062e+02, /* 0x4309a863 */
> -sa3 = 4.3456588745e+02, /* 0x43d9486f */
> -sa4 = 6.4538726807e+02, /* 0x442158c9 */
> -sa5 = 4.2900814819e+02, /* 0x43d6810b */
> -sa6 = 1.0863500214e+02, /* 0x42d9451f */
> -sa7 = 6.5702495575e+00, /* 0x40d23f7c */
> -sa8 = -6.0424413532e-02, /* 0xbd777f97 */
> -/*
> - * Coefficients for approximation to erfc in [1/.35,28]
> - */
> -rb0 = -9.8649431020e-03, /* 0xbc21a092 */
> -rb1 = -7.9928326607e-01, /* 0xbf4c9dd4 */
> -rb2 = -1.7757955551e+01, /* 0xc18e104b */
> -rb3 = -1.6063638306e+02, /* 0xc320a2ea */
> -rb4 = -6.3756646729e+02, /* 0xc41f6441 */
> -rb5 = -1.0250950928e+03, /* 0xc480230b */
> -rb6 = -4.8351919556e+02, /* 0xc3f1c275 */
> -sb1 = 3.0338060379e+01, /* 0x41f2b459 */
> -sb2 = 3.2579251099e+02, /* 0x43a2e571 */
> -sb3 = 1.5367296143e+03, /* 0x44c01759 */
> -sb4 = 3.1998581543e+03, /* 0x4547fdbb */
> -sb5 = 2.5530502930e+03, /* 0x451f90ce */
> -sb6 = 4.7452853394e+02, /* 0x43ed43a7 */
> -sb7 = -2.2440952301e+01; /* 0xc1b38712 */
Ok.
> +static const double E[] =
> + {
> + 0x1p+0, 0x1.0163da9fb3335p+0, 0x1.02c9a3e778061p+0,
> + 0x1.04315e86e7f85p+0, 0x1.059b0d3158574p+0, 0x1.0706b29ddf6dep+0,
> + 0x1.0874518759bc8p+0, 0x1.09e3ecac6f383p+0, 0x1.0b5586cf9890fp+0,
> + 0x1.0cc922b7247f7p+0, 0x1.0e3ec32d3d1a2p+0, 0x1.0fb66affed31bp+0,
> + 0x1.11301d0125b51p+0, 0x1.12abdc06c31ccp+0, 0x1.1429aaea92dep+0,
> + 0x1.15a98c8a58e51p+0, 0x1.172b83c7d517bp+0, 0x1.18af9388c8deap+0,
> + 0x1.1a35beb6fcb75p+0, 0x1.1bbe084045cd4p+0, 0x1.1d4873168b9aap+0,
> + 0x1.1ed5022fcd91dp+0, 0x1.2063b88628cd6p+0, 0x1.21f49917ddc96p+0,
> + 0x1.2387a6e756238p+0, 0x1.251ce4fb2a63fp+0, 0x1.26b4565e27cddp+0,
> + 0x1.284dfe1f56381p+0, 0x1.29e9df51fdee1p+0, 0x1.2b87fd0dad99p+0,
> + 0x1.2d285a6e4030bp+0, 0x1.2ecafa93e2f56p+0, 0x1.306fe0a31b715p+0,
> + 0x1.32170fc4cd831p+0, 0x1.33c08b26416ffp+0, 0x1.356c55f929ff1p+0,
> + 0x1.371a7373aa9cbp+0, 0x1.38cae6d05d866p+0, 0x1.3a7db34e59ff7p+0,
> + 0x1.3c32dc313a8e5p+0, 0x1.3dea64c123422p+0, 0x1.3fa4504ac801cp+0,
> + 0x1.4160a21f72e2ap+0, 0x1.431f5d950a897p+0, 0x1.44e086061892dp+0,
> + 0x1.46a41ed1d0057p+0, 0x1.486a2b5c13cdp+0, 0x1.4a32af0d7d3dep+0,
> + 0x1.4bfdad5362a27p+0, 0x1.4dcb299fddd0dp+0, 0x1.4f9b2769d2ca7p+0,
> + 0x1.516daa2cf6642p+0, 0x1.5342b569d4f82p+0, 0x1.551a4ca5d920fp+0,
> + 0x1.56f4736b527dap+0, 0x1.58d12d497c7fdp+0, 0x1.5ab07dd485429p+0,
> + 0x1.5c9268a5946b7p+0, 0x1.5e76f15ad2148p+0, 0x1.605e1b976dc09p+0,
> + 0x1.6247eb03a5585p+0, 0x1.6434634ccc32p+0, 0x1.6623882552225p+0,
> + 0x1.68155d44ca973p+0, 0x1.6a09e667f3bcdp+0, 0x1.6c012750bdabfp+0,
> + 0x1.6dfb23c651a2fp+0, 0x1.6ff7df9519484p+0, 0x1.71f75e8ec5f74p+0,
> + 0x1.73f9a48a58174p+0, 0x1.75feb564267c9p+0, 0x1.780694fde5d3fp+0,
> + 0x1.7a11473eb0187p+0, 0x1.7c1ed0130c132p+0, 0x1.7e2f336cf4e62p+0,
> + 0x1.80427543e1a12p+0, 0x1.82589994cce13p+0, 0x1.8471a4623c7adp+0,
> + 0x1.868d99b4492edp+0, 0x1.88ac7d98a6699p+0, 0x1.8ace5422aa0dbp+0,
> + 0x1.8cf3216b5448cp+0, 0x1.8f1ae99157736p+0, 0x1.9145b0b91ffc6p+0,
> + 0x1.93737b0cdc5e5p+0, 0x1.95a44cbc8520fp+0, 0x1.97d829fde4e5p+0,
> + 0x1.9a0f170ca07bap+0, 0x1.9c49182a3f09p+0, 0x1.9e86319e32323p+0,
> + 0x1.a0c667b5de565p+0, 0x1.a309bec4a2d33p+0, 0x1.a5503b23e255dp+0,
> + 0x1.a799e1330b358p+0, 0x1.a9e6b5579fdbfp+0, 0x1.ac36bbfd3f37ap+0,
> + 0x1.ae89f995ad3adp+0, 0x1.b0e07298db666p+0, 0x1.b33a2b84f15fbp+0,
> + 0x1.b59728de5593ap+0, 0x1.b7f76f2fb5e47p+0, 0x1.ba5b030a1064ap+0,
> + 0x1.bcc1e904bc1d2p+0, 0x1.bf2c25bd71e09p+0, 0x1.c199bdd85529cp+0,
> + 0x1.c40ab5fffd07ap+0, 0x1.c67f12e57d14bp+0, 0x1.c8f6d9406e7b5p+0,
> + 0x1.cb720dcef9069p+0, 0x1.cdf0b555dc3fap+0, 0x1.d072d4a07897cp+0,
> + 0x1.d2f87080d89f2p+0, 0x1.d5818dcfba487p+0, 0x1.d80e316c98398p+0,
> + 0x1.da9e603db3285p+0, 0x1.dd321f301b46p+0, 0x1.dfc97337b9b5fp+0,
> + 0x1.e264614f5a129p+0, 0x1.e502ee78b3ff6p+0, 0x1.e7a51fbc74c83p+0,
> + 0x1.ea4afa2a490dap+0, 0x1.ecf482d8e67f1p+0, 0x1.efa1bee615a27p+0,
> + 0x1.f252b376bba97p+0, 0x1.f50765b6e454p+0, 0x1.f7bfdad9cbe14p+0,
> + 0x1.fa7c1819e90d8p+0, 0x1.fd3c22b8f71f1p+0
> + };
Ok.
> -float __erfcf(float x)
> +float
> +__erfcf (float xf)
> {
Ok.
> - int32_t hx,ix;
> - float R,S,P,Q,s,y,z,r;
> - GET_FLOAT_WORD(hx,x);
> - ix = hx&0x7fffffff;
> - if(ix>=0x7f800000) { /* erfc(nan)=nan */
> - /* erfc(+-inf)=0,2 */
> - float ret = (float)(((uint32_t)hx>>31)<<1)+one/x;
> - if (FIX_INT_FP_CONVERT_ZERO && ret == 0.0f)
> - return 0.0f;
> - return ret;
> - }
> -
> - if(ix < 0x3f580000) { /* |x|<0.84375 */
> - if(ix < 0x32800000) /* |x|<2**-26 */
> - return one-x;
> - z = x*x;
> - r = pp0+z*(pp1+z*(pp2+z*(pp3+z*pp4)));
> - s = one+z*(qq1+z*(qq2+z*(qq3+z*(qq4+z*qq5))));
> - y = r/s;
> - if(hx < 0x3e800000) { /* x<1/4 */
> - return one-(x+x*y);
> - } else {
> - r = x*y;
> - r += (x-half);
> - return half - r ;
> - }
Ok.
> + float axf = fabsf (xf);
> + double axd = axf;
> + double x2 = axd * axd;
> + uint32_t t = asuint (xf);
> + unsigned int at = t & (~0u >> 1);
> + unsigned int sgn = t >> 31;
> + int64_t i = at > 0x40051000;
Ok.
> + /* for x < -0x1.ea8f94p+1, erfc(x) rounds to 2 (to nearest) */
> + if (__glibc_unlikely (t > 0xc07547ca))
> + { /* xf < -0x1.ea8f94p+1 */
> + if (__glibc_unlikely (t >= 0xff800000))
> + { /* -Inf or NaN */
> + if (t == 0xff800000)
> + return 2.0f; /* -Inf */
> + return xf + xf; /* NaN */
> }
> + return 2.0f - 0x1p-25f; /* rounds to 2 or nextbelow(2) */
> + }
Ok.
> - if(ix < 0x3fa00000) { /* 0.84375 <= |x| < 1.25 */
> - s = fabsf(x)-one;
> - P = pa0+s*(pa1+s*(pa2+s*(pa3+s*(pa4+s*(pa5+s*pa6)))));
> - Q = one+s*(qa1+s*(qa2+s*(qa3+s*(qa4+s*(qa5+s*qa6)))));
> - if(hx>=0) {
> - z = one-erx; return z - P/Q;
> - } else {
> - z = erx+P/Q; return one+z;
> - }
Ok.
> + /* at is the absolute value of xf
> + for x >= 0x1.41bbf8p+3, erfc(x) < 2^-150, thus rounds to 0 or 2^-149 */
"round to 0 or 2^-149" is confusing. Does it mean "is less than 2^-149,
which is effectively zero" ? Or should we add ", depending on rounding
mode" to clarify?
> + if (__glibc_unlikely (at >= 0x4120ddfc))
> + { /* |xf| >= 0x1.41bbf8p+3 */
> + if (__glibc_unlikely (at >= 0x7f800000))
> + { /* +Inf or NaN */
> + if (at == 0x7f800000)
> + return 0.0f; /* +Inf */
> + return xf + xf; /* NaN */
> }
> + __set_errno (ERANGE);
> + return 0x1p-149f * 0.25f; /* 0 or 2^-149 wrt rounding */
> + }
So ", depending on rounding" then ;-)
> - if (ix < 0x41e00000) { /* |x|<28 */
> - x = fabsf(x);
> - s = one/(x*x);
> - if(ix< 0x4036DB6D) { /* |x| < 1/.35 ~ 2.857143*/
> - R=ra0+s*(ra1+s*(ra2+s*(ra3+s*(ra4+s*(
> - ra5+s*(ra6+s*ra7))))));
> - S=one+s*(sa1+s*(sa2+s*(sa3+s*(sa4+s*(
> - sa5+s*(sa6+s*(sa7+s*sa8)))))));
> - } else { /* |x| >= 1/.35 ~ 2.857143 */
> - if(hx<0&&ix>=0x40c00000) return two-tiny;/* x < -6 */
> - R=rb0+s*(rb1+s*(rb2+s*(rb3+s*(rb4+s*(
> - rb5+s*rb6)))));
> - S=one+s*(sb1+s*(sb2+s*(sb3+s*(sb4+s*(
> - sb5+s*(sb6+s*sb7))))));
> - }
> - GET_FLOAT_WORD(ix,x);
> - SET_FLOAT_WORD(z,ix&0xffffe000);
> - r = __ieee754_expf(-z*z-(float)0.5625)*
> - __ieee754_expf((z-x)*(z+x)+R/S);
> - if(hx>0) {
> - float ret = math_narrow_eval (r/x);
> - if (ret == 0)
> - __set_errno (ERANGE);
> - return ret;
> - } else
> - return two-r/x;
> - } else {
> - if(hx>0) {
> - __set_errno (ERANGE);
> - return tiny*tiny;
> - } else
> - return two-tiny;
Ok.
> + if (__glibc_unlikely (at <= 0x3db80000))
> + { /* |x| <= 0x1.7p-4 */
> + if (__glibc_unlikely (t == 0xb76c9f62))
> + return 0x1.00010ap+0f + 0x1p-25f; /* exceptional case */
> + /* for |x| <= 0x1.c5bf88p-26. erfc(x) rounds to 1 (to nearest) */
> + if (__glibc_unlikely (at <= 0x32e2dfc4))
> + { /* |x| <= 0x1.c5bf88p-26 */
> + if (__glibc_unlikely (at == 0))
> + return 1.0f;
> + static const float d[] = { -0x1p-26, 0x1p-25 };
> + return 1.0f + d[sgn];
> }
Ok.
> + /* around 0, erfc(x) behaves as 1 + (odd polynomial) */
> + static const double c[] =
> + {
> + 0x1.20dd750429b6dp+0, -0x1.812746b03610bp-2, 0x1.ce2f218831d2fp-4,
> + -0x1.b82c609607dcbp-6, 0x1.553af09b8008ep-8
> + };
> + double f0 = xf
> + * (c[0] + x2 * (c[1] + x2 * (c[2] + x2 * (c[3] + x2 * (c[4])))));
> + return 1.0 - f0;
> + }
> +
The comment says "1 + poly" but the code says "1 - poly". Which is correct?
(technically, both, I know, but they should be consistent)
> + /* now -0x1.ea8f94p+1 <= x <= 0x1.41bbf8p+3, with |x| > 0x1.7p-4 */
> + const double iln2 = 0x1.71547652b82fep+0;
> + const double ln2h = 0x1.62e42fefap-8;
> + const double ln2l = 0x1.cf79abd6f5dc8p-47;
> + uint64_t jt = asuint64 (fma (x2, iln2, -(1024 + 0x1p-8)));
> + int64_t j = (int64_t) (jt << 12) >> 48;
> + double S = asdouble (((j >> 7) + (0x3ff | sgn << 11)) << 52);
> + static const double ch[] =
> + {
> + -0x1.ffffffffff333p-2, 0x1.5555555556a14p-3, -0x1.55556666659b4p-5,
> + 0x1.1111074cc7b22p-7
> + };
> + double d = (x2 + ln2h * j) + ln2l * j;
> + double d2 = d * d;
> + double e0 = E[j & 127];
> + double f = d + d2 * ((ch[0] + d * ch[1]) + d2 * (ch[2] + d * ch[3]));
> + static const double ct[][16] =
> + {
> + {
> + 0x1.c162355429b28p-1, 3.70, 0x1.da951cece2b85p-2,
3.70 is not exactly representable in float32; the hex constant should be
used here.
> + -0x1.70ef6cff4bcc4p+0, 0x1.3d7f7b3d617dep+1, -0x1.9d0aa47537c51p+1,
> + 0x1.9754ea9a3fcb1p+1, -0x1.27a5453fcc015p+1, 0x1.1ef2e0531aebap+0,
> + -0x1.eca090f5a1c06p-3, -0x1.7a3cd173a063cp-4, 0x1.30fa68a68fdddp-4,
> + 0x1.55ad9a326993ap-10, -0x1.07e7b0bb39fbfp-6, 0x1.2328706c0e95p-10,
> + 0x1.d6aa0b7b19cfep-9
> + },
> + {
> + 0x1.137c8983f8516p+2, 2.95, 0x1.05b53aa241333p-3,
Same here.
> + -0x1.a3f53872bf87p-3, 0x1.de4c30742c9d5p-4, -0x1.cb24bfa591986p-5,
> + 0x1.666aec059ca5fp-6, -0x1.a61250eb26b0bp-8, 0x1.2b28b7924b34dp-10,
> + 0x1.41b13a9d45013p-15, -0x1.6dd5e8a273613p-14, 0x1.09ce8ea5e8da5p-16,
> + 0x1.33923b4102981p-18, -0x1.1dfd161e3f984p-19, -0x1.c87618fcae3b3p-23,
> + 0x1.e8a6ffa0ba2c7p-23
> + }
> + };
> + double z = (axd - ct[i][0]) / (axd + ct[i][1]);
> + double z2 = z * z, z4 = z2 * z2;
> + double z8 = z4 * z4;
> + const double *c = ct[i] + 3;
> + double s = (((c[0] + z * c[1]) + z2 * (c[2] + z * c[3]))
> + + z4 * ((c[4] + z * c[5]) + z2 * (c[6] + z * c[7])))
> + + z8 * (((c[8] + z * c[9]) + z2 * (c[10] + z * c[11])) + z4 * (c[12]));
> + s = ct[i][2] + z * s;
> + static const double off[] = { 0, 2 };
> + double r = (S * (e0 - f * e0)) * s;
> + double y = off[sgn] + r;
> + return y;
> }
> libm_alias_float (__erfc, erfc)
Ok.
> diff --git a/sysdeps/loongarch/lp64/libm-test-ulps b/sysdeps/loongarch/lp64/libm-test-ulps
> diff --git a/sysdeps/m68k/m680x0/fpu/libm-test-ulps b/sysdeps/m68k/m680x0/fpu/libm-test-ulps
> diff --git a/sysdeps/microblaze/libm-test-ulps b/sysdeps/microblaze/libm-test-ulps
> diff --git a/sysdeps/mips/mips32/libm-test-ulps b/sysdeps/mips/mips32/libm-test-ulps
> diff --git a/sysdeps/mips/mips64/libm-test-ulps b/sysdeps/mips/mips64/libm-test-ulps
> diff --git a/sysdeps/nios2/libm-test-ulps b/sysdeps/nios2/libm-test-ulps
> diff --git a/sysdeps/or1k/fpu/libm-test-ulps b/sysdeps/or1k/fpu/libm-test-ulps
> diff --git a/sysdeps/or1k/nofpu/libm-test-ulps b/sysdeps/or1k/nofpu/libm-test-ulps
> diff --git a/sysdeps/powerpc/fpu/libm-test-ulps b/sysdeps/powerpc/fpu/libm-test-ulps
> diff --git a/sysdeps/powerpc/nofpu/libm-test-ulps b/sysdeps/powerpc/nofpu/libm-test-ulps
> diff --git a/sysdeps/riscv/nofpu/libm-test-ulps b/sysdeps/riscv/nofpu/libm-test-ulps
> diff --git a/sysdeps/riscv/rvd/libm-test-ulps b/sysdeps/riscv/rvd/libm-test-ulps
> diff --git a/sysdeps/s390/fpu/libm-test-ulps b/sysdeps/s390/fpu/libm-test-ulps
> diff --git a/sysdeps/sh/libm-test-ulps b/sysdeps/sh/libm-test-ulps
> diff --git a/sysdeps/sparc/fpu/libm-test-ulps b/sysdeps/sparc/fpu/libm-test-ulps
> diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps
Ok.
More information about the Libc-alpha
mailing list