From owner-freebsd-standards@FreeBSD.ORG Thu Jan 14 01:50:02 2010 Return-Path: Delivered-To: freebsd-standards@hub.freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2001:4f8:fff6::34]) by hub.freebsd.org (Postfix) with ESMTP id 0C7F6106566B for ; Thu, 14 Jan 2010 01:50:02 +0000 (UTC) (envelope-from gnats@FreeBSD.org) Received: from freefall.freebsd.org (freefall.freebsd.org [IPv6:2001:4f8:fff6::28]) by mx1.freebsd.org (Postfix) with ESMTP id BF8168FC16 for ; Thu, 14 Jan 2010 01:50:01 +0000 (UTC) Received: from freefall.freebsd.org (localhost [127.0.0.1]) by freefall.freebsd.org (8.14.3/8.14.3) with ESMTP id o0E1o1qR080601 for ; Thu, 14 Jan 2010 01:50:01 GMT (envelope-from gnats@freefall.freebsd.org) Received: (from gnats@localhost) by freefall.freebsd.org (8.14.3/8.14.3/Submit) id o0E1o1Z6080600; Thu, 14 Jan 2010 01:50:01 GMT (envelope-from gnats) Resent-Date: Thu, 14 Jan 2010 01:50:01 GMT Resent-Message-Id: <201001140150.o0E1o1Z6080600@freefall.freebsd.org> Resent-From: FreeBSD-gnats-submit@FreeBSD.org (GNATS Filer) Resent-To: freebsd-standards@FreeBSD.org Resent-Reply-To: FreeBSD-gnats-submit@FreeBSD.org, "Steven G. Kargl" Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2001:4f8:fff6::34]) by hub.freebsd.org (Postfix) with ESMTP id 1D1F3106568B for ; Thu, 14 Jan 2010 01:43:18 +0000 (UTC) (envelope-from kargl@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (troutmask.apl.washington.edu [128.208.78.105]) by mx1.freebsd.org (Postfix) with ESMTP id D507B8FC08 for ; Thu, 14 Jan 2010 01:43:17 +0000 (UTC) Received: from troutmask.apl.washington.edu (localhost.apl.washington.edu [127.0.0.1]) by troutmask.apl.washington.edu (8.14.3/8.14.3) with ESMTP id o0E1hOlt040084 for ; Wed, 13 Jan 2010 17:43:24 -0800 (PST) (envelope-from kargl@troutmask.apl.washington.edu) Received: (from kargl@localhost) by troutmask.apl.washington.edu (8.14.3/8.14.3/Submit) id o0E1hOm8040083; Wed, 13 Jan 2010 17:43:24 -0800 (PST) (envelope-from kargl) Message-Id: <201001140143.o0E1hOm8040083@troutmask.apl.washington.edu> Date: Wed, 13 Jan 2010 17:43:24 -0800 (PST) From: "Steven G. Kargl" To: FreeBSD-gnats-submit@FreeBSD.org X-Send-Pr-Version: 3.113 Cc: Subject: standards/142803: j0 Bessel function inaccurate near zeros of the function X-BeenThere: freebsd-standards@freebsd.org X-Mailman-Version: 2.1.5 Precedence: list Reply-To: "Steven G. Kargl" List-Id: Standards compliance List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Thu, 14 Jan 2010 01:50:02 -0000 >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: