From 1248c1c41508387ff282b208636737e8cdc9b5b0 Mon Sep 17 00:00:00 2001 From: Petr Baudis Date: Fri, 9 Sep 2011 22:16:10 -0400 Subject: [PATCH] Fix jn precision --- ChangeLog | 20 +++++ NEWS | 5 +- math/libm-test.inc | 12 ++- sysdeps/i386/fpu/libm-test-ulps | 54 +++++++++++++- sysdeps/ieee754/dbl-64/e_jn.c | 11 ++- sysdeps/ieee754/flt-32/e_jnf.c | 11 ++- sysdeps/ieee754/ldbl-128/e_jnl.c | 11 ++- sysdeps/ieee754/ldbl-128ibm/e_jnl.c | 11 ++- sysdeps/ieee754/ldbl-96/e_jnl.c | 11 ++- sysdeps/x86_64/fpu/libm-test-ulps | 142 +++++++++++++++++++++++++----------- 10 files changed, 232 insertions(+), 56 deletions(-) diff --git a/ChangeLog b/ChangeLog index 6ea3e11..663765f 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,5 +1,25 @@ 2011-09-09 Ulrich Drepper + * sysdeps/i386/fpu/libm-test-ulps: Adjust ULPs for jn tests. + +2011-07-03 Andreas Jaeger + + * math/libm-test.inc (jn_test): Add tests for BZ#11589. + * sysdeps/x86_64/fpu/libm-test-ulps: Add new ULPs for jn_test, + regenerate with gen-libm-tests.pl. + +2010-05-12 Petr Baudis + + [BZ #11589] + * sysdeps/ieee754/dbl-64/e_jn.c: Compensate major precision loss + around j0() zero points by switching to j1(). + * sysdeps/ieee754/flt-32/e_jnf.c: Likewise. + * sysdeps/ieee754/ldbl-128/e_jnl.c: Likewise. + * sysdeps/ieee754/ldbl-128ibm/e_jnl.c: Likewise. + * sysdeps/ieee754/ldbl-96/e_jnl.c: Likewise. + +2011-09-09 Ulrich Drepper + * sysdeps/unix/bsd/bsd4.4/bits/socket.h (__cmsg_nxthdr): Use NULL instead of 0. * sysdeps/unix/sysv/linux/bits/socket.h (__cmsg_nxthdr): Use (void*)0 diff --git a/NEWS b/NEWS index f431c9d..d122a73 100644 --- a/NEWS +++ b/NEWS @@ -9,8 +9,9 @@ Version 2.15 * The following bugs are resolved with this release: - 9696, 12403, 12847, 12868, 12852, 12874, 12885, 12907, 12922, 12935, - 13007, 13021, 13068, 13092, 13114, 13118, 13123, 13134, 13138, 13150 + 9696, 11589, 12403, 12847, 12868, 12852, 12874, 12885, 12907, 12922, + 12935, 13007, 13021, 13068, 13092, 13114, 13118, 13123, 13134, 13138, + 13150 * New program pldd to list loaded object of a process Implemented by Ulrich Drepper. diff --git a/math/libm-test.inc b/math/libm-test.inc index 301d4a8..70d936c 100644 --- a/math/libm-test.inc +++ b/math/libm-test.inc @@ -1,4 +1,4 @@ -/* Copyright (C) 1997-2006, 2007, 2009, 2010 Free Software Foundation, Inc. +/* Copyright (C) 1997-2006, 2007, 2009, 2010, 2011 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Andreas Jaeger , 1997. @@ -3285,6 +3285,16 @@ jn_test (void) TEST_ff_f (jn, 10, 2.0, 0.251538628271673670963516093751820639e-6L); TEST_ff_f (jn, 10, 10.0, 0.207486106633358857697278723518753428L); + /* BZ #11589 .*/ + TEST_ff_f (jn, 2, 2.4048255576957729, 0.43175480701968038399746111312430703L); + TEST_ff_f (jn, 3, 2.4048255576957729, 0.19899990535769083404042146764530813L); + TEST_ff_f (jn, 4, 2.4048255576957729, 0.647466661641779720084932282551219891E-1L); + TEST_ff_f (jn, 5, 2.4048255576957729, 0.163892432048058525099230549946147698E-1L); + TEST_ff_f (jn, 6, 2.4048255576957729, 0.34048184720278336646673682895929161E-2L); + TEST_ff_f (jn, 7, 2.4048255576957729, 0.60068836573295394221291569249883076E-3L); + TEST_ff_f (jn, 8, 2.4048255576957729, 0.92165786705344923232879022467054148E-4L); + TEST_ff_f (jn, 9, 2.4048255576957729, 0.12517270977961513005428966643852564E-4L) + END (jn); } diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps index 4b1a9e7..ebd46b0 100644 --- a/sysdeps/i386/fpu/libm-test-ulps +++ b/sysdeps/i386/fpu/libm-test-ulps @@ -640,6 +640,52 @@ double: 1 idouble: 1 ildouble: 1 ldouble: 1 +Test "jn (2, 2.4048255576957729) == 0.43175480701968038399746111312430703": +float: 1 +ifloat: 1 +double: 1 +idouble: 1 +ldouble: 82 +ildouble: 82 +Test "jn (3, 2.4048255576957729) == 0.19899990535769083404042146764530813": +ldouble: 186 +ildouble: 186 +Test "jn (4, 2.4048255576957729) == 0.647466661641779720084932282551219891E-1": +ldouble: 185 +ildouble: 185 +Test: "jn (5, 2.4048255576957729) == 0.163892432048058525099230549946147698E-1": +float: 1 +ifloat: 1 +double: 1 +idouble: 1 +ldouble: 249 +ildouble: 249 +Test "jn (6, 2.4048255576957729) == 0.34048184720278336646673682895929161E-2": +float: 2 +ifloat: 2 +double: 1 +idouble: 1 +ldouble: 511 +ildouble: 511 +Test "jn (7, 2.4048255576957729) == 0.60068836573295394221291569249883076E-3": +float: 2 +ifloat: 2 +double: 1 +idouble: 1 +ldouble: 428 +ildouble: 428 +Test "jn (8, 2.4048255576957729) == 0.92165786705344923232879022467054148E-4": +float: 3 +ifloat: 3 +double: 1 +idouble: 1 +ldouble: 609 +ildouble: 609 +Test "jn (9, 2.4048255576957729) == 0.12517270977961513005428966643852564E-4": +float: 4 +ifloat: 4 +ldouble: 750 +ildouble: 750 # lgamma Test "lgamma (-0.5) == log(2*sqrt(pi))": @@ -1168,11 +1214,11 @@ ldouble: 1 Function: "jn": double: 5 -float: 2 +float: 4 idouble: 5 -ifloat: 2 -ildouble: 2 -ldouble: 2 +ifloat: 4 +ildouble: 750 +ldouble: 750 Function: "lgamma": double: 1 diff --git a/sysdeps/ieee754/dbl-64/e_jn.c b/sysdeps/ieee754/dbl-64/e_jn.c index bf4a13d..d9d6f91 100644 --- a/sysdeps/ieee754/dbl-64/e_jn.c +++ b/sysdeps/ieee754/dbl-64/e_jn.c @@ -215,7 +215,16 @@ static double zero = 0.00000000000000000000e+00; } } } - b = (t*__ieee754_j0(x)/b); + /* j0() and j1() suffer enormous loss of precision at and + * near zero; however, we know that their zero points never + * coincide, so just choose the one further away from zero. + */ + z = __ieee754_j0 (x); + w = __ieee754_j1 (x); + if (fabs (z) >= fabs (w)) + b = (t * z / b); + else + b = (t * w / a); } } if(sgn==1) return -b; else return b; diff --git a/sysdeps/ieee754/flt-32/e_jnf.c b/sysdeps/ieee754/flt-32/e_jnf.c index de2e53d..dd3d551 100644 --- a/sysdeps/ieee754/flt-32/e_jnf.c +++ b/sysdeps/ieee754/flt-32/e_jnf.c @@ -165,7 +165,16 @@ static float zero = 0.0000000000e+00; } } } - b = (t*__ieee754_j0f(x)/b); + /* j0() and j1() suffer enormous loss of precision at and + * near zero; however, we know that their zero points never + * coincide, so just choose the one further away from zero. + */ + z = __ieee754_j0f (x); + w = __ieee754_j1f (x); + if (fabsf (z) >= fabsf (w)) + b = (t * z / b); + else + b = (t * w / a); } } if(sgn==1) return -b; else return b; diff --git a/sysdeps/ieee754/ldbl-128/e_jnl.c b/sysdeps/ieee754/ldbl-128/e_jnl.c index a4a4e24..a7f6772 100644 --- a/sysdeps/ieee754/ldbl-128/e_jnl.c +++ b/sysdeps/ieee754/ldbl-128/e_jnl.c @@ -285,7 +285,16 @@ __ieee754_jnl (n, x) } } } - b = (t * __ieee754_j0l (x) / b); + /* j0() and j1() suffer enormous loss of precision at and + * near zero; however, we know that their zero points never + * coincide, so just choose the one further away from zero. + */ + z = __ieee754_j0l (x); + w = __ieee754_j1l (x); + if (fabsl (z) >= fabsl (w)) + b = (t * z / b); + else + b = (t * w / a); } } if (sgn == 1) diff --git a/sysdeps/ieee754/ldbl-128ibm/e_jnl.c b/sysdeps/ieee754/ldbl-128ibm/e_jnl.c index 0eea745..372f942 100644 --- a/sysdeps/ieee754/ldbl-128ibm/e_jnl.c +++ b/sysdeps/ieee754/ldbl-128ibm/e_jnl.c @@ -286,7 +286,16 @@ __ieee754_jnl (n, x) } } } - b = (t * __ieee754_j0l (x) / b); + /* j0() and j1() suffer enormous loss of precision at and + * near zero; however, we know that their zero points never + * coincide, so just choose the one further away from zero. + */ + z = __ieee754_j0l (x); + w = __ieee754_j1l (x); + if (fabsl (z) >= fabsl (w)) + b = (t * z / b); + else + b = (t * w / a); } } if (sgn == 1) diff --git a/sysdeps/ieee754/ldbl-96/e_jnl.c b/sysdeps/ieee754/ldbl-96/e_jnl.c index 3d715d3..bedff7d 100644 --- a/sysdeps/ieee754/ldbl-96/e_jnl.c +++ b/sysdeps/ieee754/ldbl-96/e_jnl.c @@ -281,7 +281,16 @@ __ieee754_jnl (n, x) } } } - b = (t * __ieee754_j0l (x) / b); + /* j0() and j1() suffer enormous loss of precision at and + * near zero; however, we know that their zero points never + * coincide, so just choose the one further away from zero. + */ + z = __ieee754_j0l (x); + w = __ieee754_j1l (x); + if (fabsl (z) >= fabsl (w)) + b = (t * z / b); + else + b = (t * w / a); } } if (sgn == 1) diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps index 0ced4be..25135e3 100644 --- a/sysdeps/x86_64/fpu/libm-test-ulps +++ b/sysdeps/x86_64/fpu/libm-test-ulps @@ -103,8 +103,8 @@ double: 1 float: 1 idouble: 1 ifloat: 1 -ldouble: 1 ildouble: 1 +ldouble: 1 # catan Test "Real part of: catan (-2 - 3 i) == -1.4099210495965755225306193844604208 - 0.22907268296853876629588180294200276 i": @@ -168,10 +168,10 @@ ifloat: 1 ildouble: 1 ldouble: 1 Test "Imaginary part of: ccos (0.75 + 1.25 i) == 1.38173873063425888530729933139078645 - 1.09193013555397466170919531722024128 i": -ildouble: 1 -ldouble: 1 float: 1 ifloat: 1 +ildouble: 1 +ldouble: 1 # ccosh Test "Real part of: ccosh (-2 - 3 i) == -3.72454550491532256547397070325597253 + 0.511822569987384608834463849801875634 i": @@ -304,6 +304,9 @@ idouble: 1 ifloat: 1 # cos +Test "cos (0.80190127184058835) == 0.69534156199418473": +double: 1 +idouble: 1 Test "cos (M_PI_6l * 2.0) == 0.5": double: 1 float: 1 @@ -323,16 +326,13 @@ idouble: 1 ifloat: 1 ildouble: 1 ldouble: 1 -Test "cos (0.80190127184058835) == 0.69534156199418473": -double: 1 -idouble: 1 # cpow Test "Real part of: cpow (0.75 + 1.25 i, 0.0 + 1.0 i) == 0.331825439177608832276067945276730566 + 0.131338600281188544930936345230903032 i": float: 1 ifloat: 1 -ldouble: 1 ildouble: 1 +ldouble: 1 Test "Imaginary part of: cpow (0.75 + 1.25 i, 0.0 + 1.0 i) == 0.331825439177608832276067945276730566 + 0.131338600281188544930936345230903032 i": float: 1 ifloat: 1 @@ -380,15 +380,15 @@ ildouble: 1 ldouble: 1 # csin +Test "Imaginary part of: csin (-2 - 3 i) == -9.15449914691142957346729954460983256 + 4.16890695996656435075481305885375484 i": +double: 1 +idouble: 1 Test "Real part of: csin (0.75 + 1.25 i) == 1.28722291002649188575873510790565441 + 1.17210635989270256101081285116138863 i": ildouble: 1 ldouble: 1 Test "Imaginary part of: csin (0.75 + 1.25 i) == 1.28722291002649188575873510790565441 + 1.17210635989270256101081285116138863 i": float: 1 ifloat: 1 -Test "Imaginary part of: csin (-2 - 3 i) == -9.15449914691142957346729954460983256 + 4.16890695996656435075481305885375484 i": -double: 1 -idouble: 1 # csinh Test "Real part of: csinh (-2 - 3 i) == 3.59056458998577995201256544779481679 - 0.530921086248519805267040090660676560 i": @@ -440,12 +440,12 @@ ldouble: 3 # ctanh Test "Real part of: ctanh (-2 - 3 i) == -0.965385879022133124278480269394560686 + 0.988437503832249372031403430350121098e-2 i": +double: 1 float: 2 +idouble: 1 ifloat: 2 ildouble: 5 ldouble: 5 -double: 1 -idouble: 1 Test "Imaginary part of: ctanh (-2 - 3 i) == -0.965385879022133124278480269394560686 + 0.988437503832249372031403430350121098e-2 i": ildouble: 25 ldouble: 25 @@ -456,10 +456,10 @@ Test "Real part of: ctanh (0.75 + 1.25 i) == 1.372607570533783202580486065712268 double: 1 idouble: 1 Test "Imaginary part of: ctanh (0.75 + 1.25 i) == 1.37260757053378320258048606571226857 + 0.385795952609750664177596760720790220 i": -ildouble: 1 -ldouble: 1 double: 1 idouble: 1 +ildouble: 1 +ldouble: 1 # erf Test "erf (1.25) == 0.922900128256458230136523481197281140": @@ -481,26 +481,26 @@ ldouble: 1 # exp10 Test "exp10 (-1) == 0.1": +double: 2 +float: 1 +idouble: 2 +ifloat: 1 ildouble: 1 ldouble: 1 +Test "exp10 (0.75) == 5.62341325190349080394951039776481231": +double: 1 float: 1 +idouble: 1 ifloat: 1 -double: 2 -idouble: 2 -Test "exp10 (0.75) == 5.62341325190349080394951039776481231": ildouble: 2 ldouble: 2 -float: 1 -ifloat: 1 -double: 1 -idouble: 1 Test "exp10 (3) == 1000": -ildouble: 8 -ldouble: 8 -float: 2 -ifloat: 2 double: 6 +float: 2 idouble: 6 +ifloat: 2 +ildouble: 8 +ldouble: 8 # expm1 Test "expm1 (0.75) == 1.11700001661267466854536981983709561": @@ -660,10 +660,19 @@ ifloat: 3 ildouble: 2 ldouble: 2 Test "jn (10, 2.0) == 0.251538628271673670963516093751820639e-6": +double: 1 float: 4 +idouble: 1 ifloat: 4 ildouble: 1 ldouble: 1 +Test "jn (2, 2.4048255576957729) == 0.43175480701968038399746111312430703": +double: 2 +float: 1 +idouble: 2 +ifloat: 1 +ildouble: 82 +ldouble: 82 Test "jn (3, -1.0) == -0.0195633539826684059189053216217515083": ildouble: 1 ldouble: 1 @@ -694,6 +703,51 @@ idouble: 1 ifloat: 2 ildouble: 1 ldouble: 1 +Test "jn (3, 2.4048255576957729) == 0.19899990535769083404042146764530813": +double: 3 +idouble: 3 +ildouble: 186 +ldouble: 186 +Test "jn (4, 2.4048255576957729) == 0.647466661641779720084932282551219891E-1": +double: 1 +idouble: 1 +ildouble: 185 +ldouble: 185 +Test "jn (5, 2.4048255576957729) == 0.163892432048058525099230549946147698E-1": +double: 3 +float: 1 +idouble: 3 +ifloat: 1 +ildouble: 249 +ldouble: 249 +Test "jn (6, 2.4048255576957729) == 0.34048184720278336646673682895929161E-2": +double: 4 +float: 3 +idouble: 4 +ifloat: 3 +ildouble: 511 +ldouble: 511 +Test "jn (7, 2.4048255576957729) == 0.60068836573295394221291569249883076E-3": +double: 3 +float: 5 +idouble: 3 +ifloat: 5 +ildouble: 428 +ldouble: 428 +Test "jn (8, 2.4048255576957729) == 0.92165786705344923232879022467054148E-4": +double: 3 +float: 2 +idouble: 3 +ifloat: 2 +ildouble: 609 +ldouble: 609 +Test "jn (9, 2.4048255576957729) == 0.12517270977961513005428966643852564E-4": +double: 1 +float: 2 +idouble: 1 +ifloat: 2 +ildouble: 750 +ldouble: 750 # lgamma Test "lgamma (-0.5) == log(2*sqrt(pi))": @@ -714,12 +768,12 @@ ldouble: 1 # log10 Test "log10 (0.75) == -0.124938736608299953132449886193870744": -ildouble: 1 -ldouble: 1 -float: 2 -ifloat: 2 double: 1 +float: 2 idouble: 1 +ifloat: 2 +ildouble: 1 +ldouble: 1 Test "log10 (e) == log10(e)": float: 1 ifloat: 1 @@ -732,6 +786,9 @@ float: 1 ifloat: 1 # sincos +Test "sincos (0.80190127184058835, &sin_res, &cos_res) puts 0.69534156199418473 in cos_res": +double: 1 +idouble: 1 Test "sincos (M_PI_6l*2.0, &sin_res, &cos_res) puts 0.5 in cos_res": double: 1 float: 1 @@ -754,9 +811,6 @@ ldouble: 1 Test "sincos (pi/6, &sin_res, &cos_res) puts 0.86602540378443864676372317075293616 in cos_res": float: 1 ifloat: 1 -Test "sincos (0.80190127184058835, &sin_res, &cos_res) puts 0.69534156199418473 in cos_res": -double: 1 -idouble: 1 # tan Test "tan (pi/4) == 1": @@ -1178,12 +1232,12 @@ ildouble: 5 ldouble: 5 Function: Imaginary part of "ctanh": +double: 1 float: 1 +idouble: 1 ifloat: 1 ildouble: 25 ldouble: 25 -double: 1 -idouble: 1 Function: "erf": double: 1 @@ -1196,12 +1250,12 @@ ildouble: 1 ldouble: 1 Function: "exp10": -ildouble: 8 -ldouble: 8 -float: 2 -ifloat: 2 double: 6 +float: 2 idouble: 6 +ifloat: 2 +ildouble: 8 +ldouble: 8 Function: "expm1": double: 1 @@ -1235,11 +1289,11 @@ ldouble: 1 Function: "jn": double: 4 -float: 4 +float: 5 idouble: 4 -ifloat: 4 -ildouble: 2 -ldouble: 2 +ifloat: 5 +ildouble: 750 +ldouble: 750 Function: "lgamma": double: 1 @@ -1250,12 +1304,12 @@ ildouble: 1 ldouble: 1 Function: "log10": +double: 1 float: 2 +idouble: 1 ifloat: 2 ildouble: 1 ldouble: 1 -double: 1 -idouble: 1 Function: "log1p": float: 1 -- 2.7.4