Date: Wed, 13 Jan 2010 17:43:24 -0800 (PST) From: "Steven G. Kargl" <kargl@troutmask.apl.washington.edu> To: FreeBSD-gnats-submit@FreeBSD.org Subject: standards/142803: j0 Bessel function inaccurate near zeros of the function Message-ID: <201001140143.o0E1hOm8040083@troutmask.apl.washington.edu> Resent-Message-ID: <201001140150.o0E1o1Z6080600@freefall.freebsd.org>
next in thread | raw e-mail | index | archive | help
>Number: 142803
>Category: standards
>Synopsis: j0 Bessel function inaccurate near zeros of the function
>Confidential: no
>Severity: serious
>Priority: low
>Responsible: freebsd-standards
>State: open
>Quarter:
>Keywords:
>Date-Required:
>Class: sw-bug
>Submitter-Id: current-users
>Arrival-Date: Thu Jan 14 01:50:01 UTC 2010
>Closed-Date:
>Last-Modified:
>Originator: Steven G. Kargl
>Release: FreeBSD 9.0-CURRENT amd64
>Organization:
APL-UW
>Environment:
System: FreeBSD troutmask.apl.washington.edu 9.0-CURRENT FreeBSD 9.0-CURRENT #3 r201269M: Mon Jan 4 13:52:03 PST 2010 kargl@troutmask.apl.washington.edu:/usr/obj/usr/src/sys/SPEW amd64
>Description:
The j0 bessel function supplied by libm is fairly inaccurate at
arguments at and near a zero of the function. Here's the first
30 zeros computed by j0f, my implementation of j0f, a 4000-bit
significand computed via MPFR (the assumed exact value), and the
relative absolute error.
x my j0f(x) libm j0f(x) MPFR j0 my err libm err
2.404825 5.6434434E-08 5.9634296E-08 5.6434400E-08 1.64 152824.59
5.520078 2.4476664E-08 2.4153294E-08 2.4476659E-08 0.31 18878.52
8.653728 1.0355323E-07 1.0359805E-07 1.0355306E-07 6.36 1694.47
11.791534 -2.4511966E-09 -3.5193941E-09 -3.5301714E-09 78243.14 781.53
14.930918 -6.6019510E-09 -6.3911618E-09 -6.4815052E-09 8962.94 6722.88
18.071064 5.1734359E-09 5.3149818E-09 5.1532318E-09 1362.83 10910.50
21.211637 -1.5023797E-07 -1.5002509E-07 -1.5023348E-07 1213.07 56347.01
24.352472 1.2524691E-07 1.2516310E-07 1.2524570E-07 41.65 2834.53
27.493479 5.4330716E-08 5.4263626E-08 5.4331104E-08 18.77 3261.75
30.634607 1.2205560E-07 1.2203689E-07 1.2205546E-07 4.85 645.39
33.775822 -2.0213101E-07 -2.0206903E-07 -2.0213095E-07 5.48 6263.95
36.917099 8.4751605E-08 8.4749573E-08 8.4751581E-08 0.99 82.59
40.058426 -1.7484851E-08 -1.7475532E-08 -1.7484840E-08 0.91 767.56
43.199791 -9.2091398E-08 -9.2135146E-08 -9.2091406E-08 2.47 13530.51
46.341187 2.1663259E-07 2.1664336E-07 2.1663259E-07 0.16 268.90
49.482609 -1.2502527E-07 -1.2504157E-07 -1.2502526E-07 2.69 23512.60
52.624050 1.8706569E-07 1.8707487E-07 1.8706569E-07 0.01 251.43
55.765511 -2.0935554E-08 -2.0932896E-08 -2.0935556E-08 0.19 227.03
58.906982 1.5637660E-07 1.5634730E-07 1.5637661E-07 0.26 892.22
62.048470 3.5779891E-08 3.5787338E-08 3.5779888E-08 0.14 403.17
Note, my j0f(x) currently uses double precision to accumulate intermediate
values. Below x = 10, I use the ascending series to compute the value.
Above x = 10, I'm using an asymptotic approximation. I haven't investigated
whether additional terms in an asymptotic approximation would pull 'my err'
for x = 11, 14, 18, and 21 closer to the exact value.
>How-To-Repeat:
>Fix:
I'll get to it eventually.
>Release-Note:
>Audit-Trail:
>Unformatted:
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?201001140143.o0E1hOm8040083>