about summary refs log tree commit diff
path: root/math/auto-libm-test-in
diff options
context:
space:
mode:
authorPaul Zimmermann <Paul.Zimmermann@inria.fr>2021-04-01 08:14:10 +0200
committerPaul Zimmermann <Paul.Zimmermann@inria.fr>2021-04-02 06:15:48 +0200
commit9acda61d94acc5348c2330f2519a14d1a4a37e73 (patch)
treedcad90e95870279c37b5be7c646b3a3f6edc15cb /math/auto-libm-test-in
parent595c22ecd8e87a27fd19270ed30fdbae9ad25426 (diff)
downloadglibc-9acda61d94acc5348c2330f2519a14d1a4a37e73.tar.gz
glibc-9acda61d94acc5348c2330f2519a14d1a4a37e73.tar.xz
glibc-9acda61d94acc5348c2330f2519a14d1a4a37e73.zip
Fix the inaccuracy of j0f/j1f/y0f/y1f [BZ #14469, #14470, #14471, #14472]
For j0f/j1f/y0f/y1f, the largest error for all binary32
inputs is reduced to at most 9 ulps for all rounding modes.

The new code is enabled only when there is a cancellation at the very end of
the j0f/j1f/y0f/y1f computation, or for very large inputs, thus should not
give any visible slowdown on average.  Two different algorithms are used:

* around the first 64 zeros of j0/j1/y0/y1, approximation polynomials of
  degree 3 are used, computed using the Sollya tool (https://www.sollya.org/)

* for large inputs, an asymptotic formula from [1] is used

[1] Fast and Accurate Bessel Function Computation,
    John Harrison, Proceedings of Arith 19, 2009.

Inputs yielding the new largest errors are added to auto-libm-test-in,
and ulps are regenerated for various targets (thanks Adhemerval Zanella).

Tested on x86_64 with --disable-multi-arch and on powerpc64le-linux-gnu.
Reviewed-by: Adhemerval Zanella  <adhemerval.zanella@linaro.org>
Diffstat (limited to 'math/auto-libm-test-in')
-rw-r--r--math/auto-libm-test-in20
1 files changed, 17 insertions, 3 deletions
diff --git a/math/auto-libm-test-in b/math/auto-libm-test-in
index 4edaaa8ee1..9fbd0c626b 100644
--- a/math/auto-libm-test-in
+++ b/math/auto-libm-test-in
@@ -5790,8 +5790,11 @@ j0 -0x1.001000001p+593
 j0 0x1p1023
 j0 0x1p16382
 j0 0x1p16383
-# the next value generates larger error bounds on x86_64 (binary32)
+# the next values yield large errors for binary32
+# (cf BZ #27670 for the xfail entry)
 j0 0x2.602774p+0 xfail-rounding:ibm128-libgcc
+j0 0x1.04c39cp+6
+j0 0x1.4b7066p+7
 # the next value exercises the flt-32 code path for x >= 2^127
 j0 0x8.2f4ecp+124
 
@@ -5825,8 +5828,11 @@ j1 0x1p-60
 j1 0x1p-100
 j1 0x1p-600
 j1 0x1p-10000
-# the next value generates larger error bounds on x86_64 (binary32)
+# the next values yield large errors in the binary32 format
+# (cf BZ #27670 for the xfail entries)
 j1 0x3.ae4b2p+0 xfail-rounding:ibm128-libgcc
+j1 0x1.2f28eap+7 xfail-rounding:binary64 xfail-rounding:binary128 xfail-rounding:intel96 xfail-rounding:ibm128-libgcc
+j1 0x1.a1d20ap+6 xfail-rounding:binary128 xfail-rounding:intel96 xfail-rounding:ibm128-libgcc
 j1 min
 j1 -min
 j1 min_subnorm
@@ -8273,8 +8279,11 @@ y0 0x1p-100
 y0 0x1p-110
 y0 0x1p-600
 y0 0x1p-10000
-# the next value generates larger error bounds on x86_64 (binary32)
+# the next values yield large errors for binary32
+# (cf BZ #16492 for the xfail entries)
 y0 0xd.3432bp-4
+y0 0x1.33eaacp+5 xfail:binary64 xfail:intel96 xfail-rounding:ibm128-libgcc
+y0 0x1.a681cep-1 xfail-rounding:ibm128-libgcc
 y0 min
 y0 min_subnorm
 
@@ -8303,6 +8312,11 @@ y1 0x1p-100
 y1 0x1p-110
 y1 0x1p-600
 y1 0x1p-10000
+# the next three values yield the largest error in the binary32 format
+# (cf BZ #27670 for the xfail entries)
+y1 0x1.065194p+7 xfail-rounding:binary64 xfail-rounding:intel96 xfail-rounding:ibm128-libgcc
+y1 0x1.c1badep+0 xfail-rounding:ibm128-libgcc
+y1 0x1.c1bc2ep+0
 y1 min
 y1 min_subnorm