diff options
| author | Paul Zimmermann <Paul.Zimmermann@inria.fr> | 2021-10-05 10:32:36 +0200 |
|---|---|---|
| committer | Paul Zimmermann <Paul.Zimmermann@inria.fr> | 2021-10-05 13:45:37 +0200 |
| commit | 6bbf7298323bf31bc43494b2201465a449778e10 (patch) | |
| tree | 4dabd416c1e18864b18b09549d494f101df36420 /sysdeps | |
| parent | a312e8fe6d89f5eae6a4583d5db577121e61c0b5 (diff) | |
| download | glibc-6bbf7298323bf31bc43494b2201465a449778e10.tar.gz glibc-6bbf7298323bf31bc43494b2201465a449778e10.tar.xz glibc-6bbf7298323bf31bc43494b2201465a449778e10.zip | |
Fixed inaccuracy of j0f (BZ #28185)
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 Inputs that were yielding huge errors have been added to "make check". Reviewed-by: Adhemeral Zanella <adhemerval.zanella@linaro.org>
Diffstat (limited to 'sysdeps')
| -rw-r--r-- | sysdeps/ieee754/flt-32/e_j0f.c | 6 |
1 files changed, 3 insertions, 3 deletions
diff --git a/sysdeps/ieee754/flt-32/e_j0f.c b/sysdeps/ieee754/flt-32/e_j0f.c index 2c7f062c4a..0453a30109 100644 --- a/sysdeps/ieee754/flt-32/e_j0f.c +++ b/sysdeps/ieee754/flt-32/e_j0f.c @@ -39,7 +39,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 @@ -211,7 +211,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 */ + 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); @@ -285,7 +285,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; + float threshold = 0x1.579b26p-4f; if (fabsf (cc) > threshold) return z; else |