openlibm icon indicating copy to clipboard operation
openlibm copied to clipboard

Test failures when using icc

Open certik opened this issue 6 years ago • 7 comments

I first tested with gcc, and all tests pass. Then I applied the following patch:

--- a/Make.inc
+++ b/Make.inc
@@ -43,7 +43,7 @@ CFLAGS_add += -fno-builtin -fno-strict-aliasing
 endif
 
 ifeq ($(USEGCC),1)
-CC = gcc
+CC = icc
 CFLAGS_add += -fno-gnu89-inline -fno-builtin
 endif

And recompiled with icc 16.0.2. Everything compiled and here are the test errors that I got:

test/test-double
testing double (without inline functions)
Failure: acos (NaN) == NaN: Exception "Invalid operation" set
Failure: asin (NaN) == NaN: Exception "Invalid operation" set
Failure: cos (inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: cos (inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: cos (-inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: cos (-inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: cos (NaN) == NaN
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: sin (inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: sin (inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: sin (-inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: sin (-inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sin (NaN) == NaN
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: sincos (inf, &sin_res, &cos_res) puts NaN in sin_res plus invalid exception: Exception "Invalid operation" not set
Failure: Test: sincos (inf, &sin_res, &cos_res) puts NaN in sin_res plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sincos (inf, &sin_res, &cos_res) puts NaN in cos_res
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: sincos (-inf, &sin_res, &cos_res) puts NaN in sin_res plus invalid exception: Exception "Invalid operation" not set
Failure: Test: sincos (-inf, &sin_res, &cos_res) puts NaN in sin_res plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sincos (-inf, &sin_res, &cos_res) puts NaN in cos_res
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sincos (NaN, &sin_res, &cos_res) puts NaN in sin_res
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sincos (NaN, &sin_res, &cos_res) puts NaN in cos_res
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: sincos (pi/6, &sin_res, &cos_res) puts 0.86602540378443864676372317075293616 in cos_res
Result:
 is:          8.66025403784438485566e-01   0x1.bb67ae8584ca90000000p-1
 should be:   8.66025403784438596588e-01   0x1.bb67ae8584caa0000000p-1
 difference:  1.11022302462515654042e-16   0x1.00000000000000000000p-53
 ulp       :  1.0000
 max.ulp   :  0.0000
Failure: tan (inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: tan (inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: tan (-inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: tan (-inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: tan (NaN) == NaN
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: acosh (7) == 2.633915793849633417250092694615937
Result:
 is:          2.63391579384963314681e+00   0x1.51242719804340000000p+1
 should be:   2.63391579384963359090e+00   0x1.51242719804350000000p+1
 difference:  4.44089209850062616169e-16   0x1.00000000000000000000p-51
 ulp       :  1.0000
 max.ulp   :  0.0000
Maximal error of `acosh'
 is      :  1.0000 ulp
 accepted:  0.0000 ulp
Failure: atanh (NaN) == NaN: Exception "Invalid operation" set
Failure: Test: cosh (0.7) == 1.255169005630943018
Result:
 is:          1.25516900563094280230e+00   0x1.4152c1862342e0000000p+0
 should be:   1.25516900563094302434e+00   0x1.4152c1862342f0000000p+0
 difference:  2.22044604925031308085e-16   0x1.00000000000000000000p-52
 ulp       :  1.0000
 max.ulp   :  0.0000
Maximal error of `cosh'
 is      :  1.0000 ulp
 accepted:  0.0000 ulp
Failure: Test: exp (0.7) == 2.0137527074704765216
Result:
 is:          2.01375270747047618869e+00   0x1.01c2a612689860000000p+1
 should be:   2.01375270747047663278e+00   0x1.01c2a612689870000000p+1
 difference:  4.44089209850062616169e-16   0x1.00000000000000000000p-51
 ulp       :  1.0000
 max.ulp   :  0.0000
Maximal error of `exp'
 is      :  1.0000 ulp
 accepted:  0.0000 ulp
Failure: Test: exp2 (0.7) == 1.6245047927124710452
Result:
 is:          1.62362532701732886764e+00   0x1.9fa5e8d07f29e0000000p+0
 should be:   1.62450479271247094637e+00   0x1.9fdf8bcce533d0000000p+0
 difference:  8.79465695142078729418e-04   0x1.cd17e3304f8000000000p-11
 ulp       :  3960761376927.0000
 max.ulp   :  0.0000
Maximal error of `exp2'
 is      :  3960761376927.0000 ulp
 accepted:  0.0000 ulp
Failure: Test: expm1 (-0) == -0
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -0.00000000000000000000e+00  -0x0.00000000000000000000p+0
 difference:  0.00000000000000000000e+00   0x0.00000000000000000000p+0
 ulp       :  0.0000
 max.ulp   :  0.0000
Failure: Test: log (10) == M_LN10l
Result:
 is:          2.30258509299404545700e+00   0x1.26bb1bbb555150000000p+1
 should be:   2.30258509299404590109e+00   0x1.26bb1bbb555160000000p+1
 difference:  4.44089209850062616169e-16   0x1.00000000000000000000p-51
 ulp       :  1.0000
 max.ulp   :  0.0000
Failure: log10 (NaN) == NaN: Exception "Invalid operation" set
Failure: Test: log10 (e) == log10(e)
Result:
 is:          4.34294481903251761157e-01   0x1.bcb7b1526e50d0000000p-2
 should be:   4.34294481903251816668e-01   0x1.bcb7b1526e50e0000000p-2
 difference:  5.55111512312578270212e-17   0x1.00000000000000000000p-54
 ulp       :  1.0000
 max.ulp   :  0.0000
Failure: log1p (NaN) == NaN: Exception "Invalid operation" set
Failure: log2 (NaN) == NaN: Exception "Invalid operation" set
Failure: Test: log2 (e) == M_LOG2El
Result:
 is:          1.44269504088896316496e+00   0x1.71547652b82fd0000000p+0
 should be:   1.44269504088896338700e+00   0x1.71547652b82fe0000000p+0
 difference:  2.22044604925031308085e-16   0x1.00000000000000000000p-52
 ulp       :  1.0000
 max.ulp   :  0.0000
Failure: tgamma (-inf) == NaN plus invalid exception: Exception "Invalid operation" not set
Failure: Test: tgamma (-inf) == NaN plus invalid exception
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: tgamma (NaN) == NaN
Result:
 is:          0.00000000000000000000e+00   0x0.00000000000000000000p+0
 should be:  -nan  -nan
Failure: Test: nearbyint (0.5) == 0.0
Result:
 is:          5.00000000000000000000e-01   0x1.00000000000000000000p-1
 should be:   0.00000000000000000000e+00   0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (1.5) == 2.0
Result:
 is:          1.50000000000000000000e+00   0x1.80000000000000000000p+0
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (-0.5) == -0
Result:
 is:         -5.00000000000000000000e-01  -0x1.00000000000000000000p-1
 should be:  -0.00000000000000000000e+00  -0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (-1.5) == -2.0
Result:
 is:         -1.50000000000000000000e+00  -0x1.80000000000000000000p+0
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Maximal error of `nearbyint'
 is      :  2251799813685248.0000 ulp
 accepted:  0.0000 ulp
Failure: Test: rint (0.5) == 0.0
Result:
 is:          5.00000000000000000000e-01   0x1.00000000000000000000p-1
 should be:   0.00000000000000000000e+00   0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: rint (1.5) == 2.0
Result:
 is:          1.50000000000000000000e+00   0x1.80000000000000000000p+0
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (2.5) == 2.0
Result:
 is:          2.50000000000000000000e+00   0x1.40000000000000000000p+1
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (3.5) == 4.0
Result:
 is:          3.50000000000000000000e+00   0x1.c0000000000000000000p+1
 should be:   4.00000000000000000000e+00   0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (4.5) == 4.0
Result:
 is:          4.50000000000000000000e+00   0x1.20000000000000000000p+2
 should be:   4.00000000000000000000e+00   0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (-0.5) == -0.0
Result:
 is:         -5.00000000000000000000e-01  -0x1.00000000000000000000p-1
 should be:  -0.00000000000000000000e+00  -0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: rint (-1.5) == -2.0
Result:
 is:         -1.50000000000000000000e+00  -0x1.80000000000000000000p+0
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (-2.5) == -2.0
Result:
 is:         -2.50000000000000000000e+00  -0x1.40000000000000000000p+1
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (-3.5) == -4.0
Result:
 is:         -3.50000000000000000000e+00  -0x1.c0000000000000000000p+1
 should be:  -4.00000000000000000000e+00  -0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (-4.5) == -4.0
Result:
 is:         -4.50000000000000000000e+00  -0x1.20000000000000000000p+2
 should be:  -4.00000000000000000000e+00  -0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Maximal error of `rint'
 is      :  2251799813685248.0000 ulp
 accepted:  0.0000 ulp
Failure: Test: j0 (8.0) == 0.17165080713755390609
Result:
 is:          2.02346525843941249523e-01   0x1.9e67daf7abec20000000p-3
 should be:   1.71650807137553901294e-01   0x1.5f8a7557e70250000000p-3
 difference:  3.06957187063873482291e-02   0x1.f6eb2cfe274e80000000p-6
 ulp       :  1105929818623645.0000
 max.ulp   :  0.0000
Failure: Test: j0 (10.0) == -0.24593576445134833520
Result:
 is:         -1.76071407662647311243e-01  -0x1.6898204d5fc090000000p-3
 should be:  -2.45935764451348348736e-01  -0x1.f7ad2b89e1e540000000p-3
 difference:  6.98643567887010374928e-02   0x1.1e2a1679044960000000p-4
 ulp       :  2517128729600587.0000
 max.ulp   :  2.0000
Maximal error of `j0'
 is      :  2517128729600587.0000 ulp
 accepted:  2.0000 ulp
Failure: Test: j1 (8.0) == 0.23463634685391462438
Result:
 is:          2.09141602191802888910e-01   0x1.ac526ead2cb800000000p-3
 should be:   2.34636346853914629085e-01   0x1.e0890561861830000000p-3
 difference:  2.54947446621117401744e-02   0x1.a1b4b5a2cb0180000000p-6
 ulp       :  918544980481539.0000
 max.ulp   :  1.0000
Failure: Test: j1 (10.0) == 0.043472746168861436670
Result:
 is:          1.71946391306713852742e-01   0x1.60256dfa9ed500000000p-3
 should be:   4.34727461688614383317e-02   0x1.6420f4e2009110000000p-5
 difference:  1.28473645137852421350e-01   0x1.071d30c21eb0c0000000p-3
 ulp       :  18515003531832368.0000
 max.ulp   :  2.0000
Maximal error of `j1'
 is      :  18515003531832368.0000 ulp
 accepted:  2.0000 ulp
Failure: Test: jn (0, 8.0) == 0.17165080713755390609
Result:
 is:          2.02346525843941249523e-01   0x1.9e67daf7abec20000000p-3
 should be:   1.71650807137553901294e-01   0x1.5f8a7557e70250000000p-3
 difference:  3.06957187063873482291e-02   0x1.f6eb2cfe274e80000000p-6
 ulp       :  1105929818623645.0000
 max.ulp   :  0.0000
Failure: Test: jn (0, 10.0) == -0.24593576445134833520
Result:
 is:         -1.76071407662647311243e-01  -0x1.6898204d5fc090000000p-3
 should be:  -2.45935764451348348736e-01  -0x1.f7ad2b89e1e540000000p-3
 difference:  6.98643567887010374928e-02   0x1.1e2a1679044960000000p-4
 ulp       :  2517128729600587.0000
 max.ulp   :  2.0000
Failure: Test: jn (1, 8.0) == 0.23463634685391462438
Result:
 is:          2.09141602191802888910e-01   0x1.ac526ead2cb800000000p-3
 should be:   2.34636346853914629085e-01   0x1.e0890561861830000000p-3
 difference:  2.54947446621117401744e-02   0x1.a1b4b5a2cb0180000000p-6
 ulp       :  918544980481539.0000
 max.ulp   :  1.0000
Failure: Test: jn (1, 10.0) == 0.043472746168861436670
Result:
 is:          1.71946391306713852742e-01   0x1.60256dfa9ed500000000p-3
 should be:   4.34727461688614383317e-02   0x1.6420f4e2009110000000p-5
 difference:  1.28473645137852421350e-01   0x1.071d30c21eb0c0000000p-3
 ulp       :  18515003531832368.0000
 max.ulp   :  2.0000
Failure: Test: jn (3, 10.0) == 0.058379379305186812343
Result:
 is:         -8.77621169371178150298e-02  -0x1.67793fe12991b0000000p-4
 should be:   5.83793793051868153965e-02   0x1.de3e6e9c0cc060000000p-5
 difference:  1.46141496242304630426e-01   0x1.2b4c3b9797f8f0000000p-3
 ulp       :  21061209216646716.0000
 max.ulp   :  3.0000
Failure: Test: jn (10, 10.0) == 0.20748610663335885770
Result:
 is:         -2.16508100295594030360e-02  -0x1.62ba1440914480000000p-6
 should be:   2.07486106633358868834e-01   0x1.a8ee79d2eacb00000000p-3
 difference:  2.29136916662918271870e-01   0x1.d545bc5afcf390000000p-3
 ulp       :  8255527459999545.0000
 max.ulp   :  4.0000
Maximal error of `jn'
 is      :  21061209216646716.0000 ulp
 accepted:  6.0000 ulp
Failure: Test: y0 (8.0) == 0.22352148938756622053
Result:
 is:          1.96167630581339291762e-01   0x1.91c055af0b9860000000p-3
 should be:   2.23521489387566218987e-01   0x1.c9c5a276fb82c0000000p-3
 difference:  2.73538588062269272250e-02   0x1.c02a663f7f5300000000p-6
 ulp       :  985526626614950.0000
 max.ulp   :  1.0000
Failure: Test: y0 (10.0) == 0.055671167283599391424
Result:
 is:          1.80506339242369695519e-01   0x1.71ad4ebe221e40000000p-3
 should be:   5.56711672835993945374e-02   0x1.c80ee65a05c590000000p-5
 difference:  1.24835171958770307921e-01   0x1.ff532a4f4159c0000000p-4
 ulp       :  17990644285319992.0000
 max.ulp   :  2.0000
Maximal error of `y0'
 is      :  17990644285319992.0000 ulp
 accepted:  3.0000 ulp
Failure: Test: y1 (8.0) == -0.15806046173124749426
Result:
 is:         -1.90518017678883433419e-01  -0x1.862e4f79d62080000000p-3
 should be:  -1.58060461731247492434e-01  -0x1.43b5340f693360000000p-3
 difference:  3.24575559476359409850e-02   0x1.09e46da9b3b480000000p-5
 ulp       :  1169406694969042.0000
 max.ulp   :  1.0000
Failure: Test: y1 (10.0) == 0.24901542420695388392
Result:
 is:          1.85291661088754555653e-01   0x1.7b7a31c4c81fa0000000p-3
 should be:   2.49015424206953883690e-01   0x1.fdfbcc7958f3c0000000p-3
 difference:  6.37237631181993280372e-02   0x1.0503356921a840000000p-4
 ulp       :  2295890526670146.0000
 max.ulp   :  3.0000
Maximal error of `y1'
 is      :  2295890526670146.0000 ulp
 accepted:  3.0000 ulp
Failure: Test: yn (0, 8.0) == 0.22352148938756622053
Result:
 is:          1.96167630581339291762e-01   0x1.91c055af0b9860000000p-3
 should be:   2.23521489387566218987e-01   0x1.c9c5a276fb82c0000000p-3
 difference:  2.73538588062269272250e-02   0x1.c02a663f7f5300000000p-6
 ulp       :  985526626614950.0000
 max.ulp   :  1.0000
Failure: Test: yn (0, 10.0) == 0.055671167283599391424
Result:
 is:          1.80506339242369695519e-01   0x1.71ad4ebe221e40000000p-3
 should be:   5.56711672835993945374e-02   0x1.c80ee65a05c590000000p-5
 difference:  1.24835171958770307921e-01   0x1.ff532a4f4159c0000000p-4
 ulp       :  17990644285319992.0000
 max.ulp   :  2.0000
Failure: Test: yn (1, 8.0) == -0.15806046173124749426
Result:
 is:         -1.90518017678883433419e-01  -0x1.862e4f79d62080000000p-3
 should be:  -1.58060461731247492434e-01  -0x1.43b5340f693360000000p-3
 difference:  3.24575559476359409850e-02   0x1.09e46da9b3b480000000p-5
 ulp       :  1169406694969042.0000
 max.ulp   :  1.0000
Failure: Test: yn (1, 10.0) == 0.24901542420695388392
Result:
 is:          1.85291661088754555653e-01   0x1.7b7a31c4c81fa0000000p-3
 should be:   2.49015424206953883690e-01   0x1.fdfbcc7958f3c0000000p-3
 difference:  6.37237631181993280372e-02   0x1.0503356921a840000000p-4
 ulp       :  2295890526670146.0000
 max.ulp   :  3.0000
Failure: Test: yn (3, 10.0) == -0.25136265718383732978
Result:
 is:         -2.42670863898602068298e-01  -0x1.f0fd6c0117ae50000000p-3
 should be:  -2.51362657183837323593e-01  -0x1.016536602934e0000000p-2
 difference:  8.69179328523525529526e-03   0x1.1cd00bf3abb700000000p-7
 ulp       :  156577428002267.5000
 max.ulp   :  1.0000
Failure: Test: yn (10, 10.0) == -0.35981415218340272205
Result:
 is:         -4.14786633080665156914e-01  -0x1.a8bdd3bf98fe90000000p-2
 should be:  -3.59814152183402735652e-01  -0x1.70731f01100ee0000000p-2
 difference:  5.49724808972624212622e-02   0x1.c255a5f4477d80000000p-5
 ulp       :  990296177938171.0000
 max.ulp   :  2.0000
Maximal error of `yn'
 is      :  17990644285319992.0000 ulp
 accepted:  6.0000 ulp

Test suite completed:
  1118 test cases plus 932 tests for exception flags executed.
  87 errors occurred.
Makefile:55: recipe for target 'test' failed
make: *** [test] Error 1

certik avatar Nov 14 '17 21:11 certik

Since icc by default does not follow IEEE behavior, I think we can ignore the NaN errors, as those are to be expected. Then I am going to ignore the 1 ulp errors, are those are comparable to Clang (#10). This leaves the following errors that point to either bugs in the Intel compiler or in the code (or in tests):

The exp2 seems to give incorrect results:

Failure: Test: exp2 (0.7) == 1.6245047927124710452
Result:
 is:          1.62362532701732886764e+00   0x1.9fa5e8d07f29e0000000p+0
 should be:   1.62450479271247094637e+00   0x1.9fdf8bcce533d0000000p+0
 difference:  8.79465695142078729418e-04   0x1.cd17e3304f8000000000p-11
 ulp       :  3960761376927.0000
 max.ulp   :  0.0000
Maximal error of `exp2'
 is      :  3960761376927.0000 ulp
 accepted:  0.0000 ulp

The nearbyint seems to give wrong answers:

Failure: Test: nearbyint (0.5) == 0.0
Result:
 is:          5.00000000000000000000e-01   0x1.00000000000000000000p-1
 should be:   0.00000000000000000000e+00   0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (1.5) == 2.0
Result:
 is:          1.50000000000000000000e+00   0x1.80000000000000000000p+0
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (-0.5) == -0
Result:
 is:         -5.00000000000000000000e-01  -0x1.00000000000000000000p-1
 should be:  -0.00000000000000000000e+00  -0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: nearbyint (-1.5) == -2.0
Result:
 is:         -1.50000000000000000000e+00  -0x1.80000000000000000000p+0
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Maximal error of `nearbyint'
 is      :  2251799813685248.0000 ulp
 accepted:  0.0000 ulp

So does rint:

Failure: Test: rint (0.5) == 0.0
Result:
 is:          5.00000000000000000000e-01   0x1.00000000000000000000p-1
 should be:   0.00000000000000000000e+00   0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: rint (1.5) == 2.0
Result:
 is:          1.50000000000000000000e+00   0x1.80000000000000000000p+0
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (2.5) == 2.0
Result:
 is:          2.50000000000000000000e+00   0x1.40000000000000000000p+1
 should be:   2.00000000000000000000e+00   0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (3.5) == 4.0
Result:
 is:          3.50000000000000000000e+00   0x1.c0000000000000000000p+1
 should be:   4.00000000000000000000e+00   0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (4.5) == 4.0
Result:
 is:          4.50000000000000000000e+00   0x1.20000000000000000000p+2
 should be:   4.00000000000000000000e+00   0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (-0.5) == -0.0
Result:
 is:         -5.00000000000000000000e-01  -0x1.00000000000000000000p-1
 should be:  -0.00000000000000000000e+00  -0x0.00000000000000000000p+0
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  2251799813685248.0000
 max.ulp   :  0.0000
Failure: Test: rint (-1.5) == -2.0
Result:
 is:         -1.50000000000000000000e+00  -0x1.80000000000000000000p+0
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (-2.5) == -2.0
Result:
 is:         -2.50000000000000000000e+00  -0x1.40000000000000000000p+1
 should be:  -2.00000000000000000000e+00  -0x1.00000000000000000000p+1
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  1125899906842624.0000
 max.ulp   :  0.0000
Failure: Test: rint (-3.5) == -4.0
Result:
 is:         -3.50000000000000000000e+00  -0x1.c0000000000000000000p+1
 should be:  -4.00000000000000000000e+00  -0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Failure: Test: rint (-4.5) == -4.0
Result:
 is:         -4.50000000000000000000e+00  -0x1.20000000000000000000p+2
 should be:  -4.00000000000000000000e+00  -0x1.00000000000000000000p+2
 difference:  5.00000000000000000000e-01   0x1.00000000000000000000p-1
 ulp       :  562949953421312.0000
 max.ulp   :  0.0000
Maximal error of `rint'
 is      :  2251799813685248.0000 ulp
 accepted:  0.0000 ulp

And finally the Bessel functions seem to be completely wrong:

Failure: Test: j0 (8.0) == 0.17165080713755390609
Result:
 is:          2.02346525843941249523e-01   0x1.9e67daf7abec20000000p-3
 should be:   1.71650807137553901294e-01   0x1.5f8a7557e70250000000p-3
 difference:  3.06957187063873482291e-02   0x1.f6eb2cfe274e80000000p-6
 ulp       :  1105929818623645.0000
 max.ulp   :  0.0000
Failure: Test: j0 (10.0) == -0.24593576445134833520
Result:
 is:         -1.76071407662647311243e-01  -0x1.6898204d5fc090000000p-3
 should be:  -2.45935764451348348736e-01  -0x1.f7ad2b89e1e540000000p-3
 difference:  6.98643567887010374928e-02   0x1.1e2a1679044960000000p-4
 ulp       :  2517128729600587.0000
 max.ulp   :  2.0000
Maximal error of `j0'
 is      :  2517128729600587.0000 ulp
 accepted:  2.0000 ulp
Failure: Test: j1 (8.0) == 0.23463634685391462438
Result:
 is:          2.09141602191802888910e-01   0x1.ac526ead2cb800000000p-3
 should be:   2.34636346853914629085e-01   0x1.e0890561861830000000p-3
 difference:  2.54947446621117401744e-02   0x1.a1b4b5a2cb0180000000p-6
 ulp       :  918544980481539.0000
 max.ulp   :  1.0000
Failure: Test: j1 (10.0) == 0.043472746168861436670
Result:
 is:          1.71946391306713852742e-01   0x1.60256dfa9ed500000000p-3
 should be:   4.34727461688614383317e-02   0x1.6420f4e2009110000000p-5
 difference:  1.28473645137852421350e-01   0x1.071d30c21eb0c0000000p-3
 ulp       :  18515003531832368.0000
 max.ulp   :  2.0000
Maximal error of `j1'
 is      :  18515003531832368.0000 ulp
 accepted:  2.0000 ulp
Failure: Test: jn (0, 8.0) == 0.17165080713755390609
Result:
 is:          2.02346525843941249523e-01   0x1.9e67daf7abec20000000p-3
 should be:   1.71650807137553901294e-01   0x1.5f8a7557e70250000000p-3
 difference:  3.06957187063873482291e-02   0x1.f6eb2cfe274e80000000p-6
 ulp       :  1105929818623645.0000
 max.ulp   :  0.0000
Failure: Test: jn (0, 10.0) == -0.24593576445134833520
Result:
 is:         -1.76071407662647311243e-01  -0x1.6898204d5fc090000000p-3
 should be:  -2.45935764451348348736e-01  -0x1.f7ad2b89e1e540000000p-3
 difference:  6.98643567887010374928e-02   0x1.1e2a1679044960000000p-4
 ulp       :  2517128729600587.0000
 max.ulp   :  2.0000
Failure: Test: jn (1, 8.0) == 0.23463634685391462438
Result:
 is:          2.09141602191802888910e-01   0x1.ac526ead2cb800000000p-3
 should be:   2.34636346853914629085e-01   0x1.e0890561861830000000p-3
 difference:  2.54947446621117401744e-02   0x1.a1b4b5a2cb0180000000p-6
 ulp       :  918544980481539.0000
 max.ulp   :  1.0000
Failure: Test: jn (1, 10.0) == 0.043472746168861436670
Result:
 is:          1.71946391306713852742e-01   0x1.60256dfa9ed500000000p-3
 should be:   4.34727461688614383317e-02   0x1.6420f4e2009110000000p-5
 difference:  1.28473645137852421350e-01   0x1.071d30c21eb0c0000000p-3
 ulp       :  18515003531832368.0000
 max.ulp   :  2.0000
Failure: Test: jn (3, 10.0) == 0.058379379305186812343
Result:
 is:         -8.77621169371178150298e-02  -0x1.67793fe12991b0000000p-4
 should be:   5.83793793051868153965e-02   0x1.de3e6e9c0cc060000000p-5
 difference:  1.46141496242304630426e-01   0x1.2b4c3b9797f8f0000000p-3
 ulp       :  21061209216646716.0000
 max.ulp   :  3.0000
Failure: Test: jn (10, 10.0) == 0.20748610663335885770
Result:
 is:         -2.16508100295594030360e-02  -0x1.62ba1440914480000000p-6
 should be:   2.07486106633358868834e-01   0x1.a8ee79d2eacb00000000p-3
 difference:  2.29136916662918271870e-01   0x1.d545bc5afcf390000000p-3
 ulp       :  8255527459999545.0000
 max.ulp   :  4.0000
Maximal error of `jn'
 is      :  21061209216646716.0000 ulp
 accepted:  6.0000 ulp
Failure: Test: y0 (8.0) == 0.22352148938756622053
Result:
 is:          1.96167630581339291762e-01   0x1.91c055af0b9860000000p-3
 should be:   2.23521489387566218987e-01   0x1.c9c5a276fb82c0000000p-3
 difference:  2.73538588062269272250e-02   0x1.c02a663f7f5300000000p-6
 ulp       :  985526626614950.0000
 max.ulp   :  1.0000
Failure: Test: y0 (10.0) == 0.055671167283599391424
Result:
 is:          1.80506339242369695519e-01   0x1.71ad4ebe221e40000000p-3
 should be:   5.56711672835993945374e-02   0x1.c80ee65a05c590000000p-5
 difference:  1.24835171958770307921e-01   0x1.ff532a4f4159c0000000p-4
 ulp       :  17990644285319992.0000
 max.ulp   :  2.0000
Maximal error of `y0'
 is      :  17990644285319992.0000 ulp
 accepted:  3.0000 ulp
Failure: Test: y1 (8.0) == -0.15806046173124749426
Result:
 is:         -1.90518017678883433419e-01  -0x1.862e4f79d62080000000p-3
 should be:  -1.58060461731247492434e-01  -0x1.43b5340f693360000000p-3
 difference:  3.24575559476359409850e-02   0x1.09e46da9b3b480000000p-5
 ulp       :  1169406694969042.0000
 max.ulp   :  1.0000
Failure: Test: y1 (10.0) == 0.24901542420695388392
Result:
 is:          1.85291661088754555653e-01   0x1.7b7a31c4c81fa0000000p-3
 should be:   2.49015424206953883690e-01   0x1.fdfbcc7958f3c0000000p-3
 difference:  6.37237631181993280372e-02   0x1.0503356921a840000000p-4
 ulp       :  2295890526670146.0000
 max.ulp   :  3.0000
Maximal error of `y1'
 is      :  2295890526670146.0000 ulp
 accepted:  3.0000 ulp
Failure: Test: yn (0, 8.0) == 0.22352148938756622053
Result:
 is:          1.96167630581339291762e-01   0x1.91c055af0b9860000000p-3
 should be:   2.23521489387566218987e-01   0x1.c9c5a276fb82c0000000p-3
 difference:  2.73538588062269272250e-02   0x1.c02a663f7f5300000000p-6
 ulp       :  985526626614950.0000
 max.ulp   :  1.0000
Failure: Test: yn (0, 10.0) == 0.055671167283599391424
Result:
 is:          1.80506339242369695519e-01   0x1.71ad4ebe221e40000000p-3
 should be:   5.56711672835993945374e-02   0x1.c80ee65a05c590000000p-5
 difference:  1.24835171958770307921e-01   0x1.ff532a4f4159c0000000p-4
 ulp       :  17990644285319992.0000
 max.ulp   :  2.0000
Failure: Test: yn (1, 8.0) == -0.15806046173124749426
Result:
 is:         -1.90518017678883433419e-01  -0x1.862e4f79d62080000000p-3
 should be:  -1.58060461731247492434e-01  -0x1.43b5340f693360000000p-3
 difference:  3.24575559476359409850e-02   0x1.09e46da9b3b480000000p-5
 ulp       :  1169406694969042.0000
 max.ulp   :  1.0000
Failure: Test: yn (1, 10.0) == 0.24901542420695388392
Result:
 is:          1.85291661088754555653e-01   0x1.7b7a31c4c81fa0000000p-3
 should be:   2.49015424206953883690e-01   0x1.fdfbcc7958f3c0000000p-3
 difference:  6.37237631181993280372e-02   0x1.0503356921a840000000p-4
 ulp       :  2295890526670146.0000
 max.ulp   :  3.0000
Failure: Test: yn (3, 10.0) == -0.25136265718383732978
Result:
 is:         -2.42670863898602068298e-01  -0x1.f0fd6c0117ae50000000p-3
 should be:  -2.51362657183837323593e-01  -0x1.016536602934e0000000p-2
 difference:  8.69179328523525529526e-03   0x1.1cd00bf3abb700000000p-7
 ulp       :  156577428002267.5000
 max.ulp   :  1.0000
Failure: Test: yn (10, 10.0) == -0.35981415218340272205
Result:
 is:         -4.14786633080665156914e-01  -0x1.a8bdd3bf98fe90000000p-2
 should be:  -3.59814152183402735652e-01  -0x1.70731f01100ee0000000p-2
 difference:  5.49724808972624212622e-02   0x1.c255a5f4477d80000000p-5
 ulp       :  990296177938171.0000
 max.ulp   :  2.0000
Maximal error of `yn'
 is      :  17990644285319992.0000 ulp
 accepted:  6.0000 ulp

certik avatar Nov 14 '17 21:11 certik

To make icc work, once can use the -fp-model=precise option:

--- a/Make.inc
+++ b/Make.inc
@@ -43,8 +43,8 @@ CFLAGS_add += -fno-builtin -fno-strict-aliasing
 endif
 
 ifeq ($(USEGCC),1)
-CC = gcc
-CFLAGS_add += -fno-gnu89-inline -fno-builtin
+CC = icc
+CFLAGS_add += -fno-gnu89-inline -fno-builtin -fp-model=precise
 endif
 
 ARCH ?= $(shell $(CC) -dumpmachine | sed "s/\([^-]*\).*$$/\1/")

Then all tests pass:

test/test-double
testing double (without inline functions)

Test suite completed:
  1118 test cases plus 932 tests for exception flags executed.
  All tests passed successfully.
test/test-float
testing float (without inline functions)

Test suite completed:
  1101 test cases plus 923 tests for exception flags executed.
  All tests passed successfully.

certik avatar Nov 14 '17 21:11 certik

The exp2 failure is caused by icc rewriting z = x - ((x+redux)-redux) into 0. One can prevent that with the following patch:

--- a/src/s_exp2.c
+++ b/src/s_exp2.c
@@ -342,6 +342,7 @@ exp2(double x)
 {
        double r, t, twopk, twopkp1000, z;
        u_int32_t hx, ix, lx, i0;
+       volatile double t_;
        int k;

        /* Filter out exceptional cases. */
@@ -370,7 +371,8 @@ exp2(double x)
        k = (i0 >> TBLBITS) << 20;
        i0 = (i0 & (TBLSIZE - 1)) << 1;
        t -= redux;
-       z = x - t;
+       t_ = t;
+       z = x - t_;

        /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */
        t = tbl[i0];            /* exp2t[i0] */

Then the exp2 test passes, with the default icc options.

It would be nice to figure out a cleaner option for the above.

certik avatar Nov 15 '17 04:11 certik

Ok, I finally figured out what it is doing. Here is a much better patch:

--- a/src/s_exp2.c
+++ b/src/s_exp2.c
@@ -369,7 +369,7 @@ exp2(double x)
        i0 += TBLSIZE / 2;
        k = (i0 >> TBLBITS) << 20;
        i0 = (i0 & (TBLSIZE - 1)) << 1;
-       t -= redux;
+       t = floor(x*TBLSIZE) / (double)TBLSIZE;
        z = x - t;

        /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */

Now at least the z computation is mathematically correct, without using IEEE conventions. The i0 computation still relies on IEEE, but the icc compiler does the right thing there, so I didn't work on it.

The new code is obviously slower, but perhaps there is a way to speed it up.

certik avatar Nov 15 '17 04:11 certik

Got rid of the multiplication and division:

--- a/src/s_exp2.c
+++ b/src/s_exp2.c
@@ -369,7 +369,7 @@ exp2(double x)
        i0 += TBLSIZE / 2;
        k = (i0 >> TBLBITS) << 20;
        i0 = (i0 & (TBLSIZE - 1)) << 1;
-       t -= redux;
+       t = ldexp(floor(ldexp(x, TBLBITS)), -TBLBITS);
        z = x - t;

        /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */

I think this is probably as simple as it gets, while still being mathematically correct. A good compiler like icc should be able to optimize this very well, since ldexp is basically just an exponent shift in the binary floating point, and floor can also be implemented efficiently. At the end of the day, the previous code was doing exactly that, directly. I think it might be possible for the compiler to figure this out, and essentially generate an equivalent code. I haven't checked.

I did check however that icc now compiles the exp2 function properly and this particular test fully passes.

certik avatar Nov 15 '17 05:11 certik

Discovered a slight mistake:

--- a/src/s_exp2.c
+++ b/src/s_exp2.c
@@ -369,7 +369,7 @@ exp2(double x)
        i0 += TBLSIZE / 2;
        k = (i0 >> TBLBITS) << 20;
        i0 = (i0 & (TBLSIZE - 1)) << 1;
-       t -= redux;
+       t = ldexp(floor(ldexp(x, TBLBITS)+0.5), -TBLBITS);
        z = x - t;

        /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */

Both versions work for the test suite, but only this latter actually works for all values of x.

certik avatar Nov 15 '17 05:11 certik

I submitted a fix in #170.

certik avatar Nov 15 '17 08:11 certik