bugGNU Octave - Bugs: bug #48316, besselj gives NaN for moderately...

 
 

bug #48316: besselj gives NaN for moderately large inputs (1e10)

Submitted by:  Colin Macdonald <cbm>
Submitted on:  Sun 26 Jun 2016 05:55:18 AM UTC  
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  Open/Closed:  Closed
Release:  dev Operating System:  Any

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

Tue 27 Mar 2018 12:18:55 AM UTC, comment #3: 

Fixes from GSOC2018 special function project were incorporated during OctConf 2018.  Closing report.

Rik <rik5>
Project Administrator
Thu 30 Jun 2016 07:44:31 PM UTC, comment #2: 

Scipy appears to use a combination of the amos library and the cephes bessel functions (http://www.netlib.org/cephes/).

Mike Miller <mtmiller>
Project Administrator
Thu 30 Jun 2016 07:35:48 PM UTC, comment #1: 

Changing status to confirmed, and adding some details from the mailing list (http://thread.gmane.org/gmane.comp.gnu.octave.maintainers/46617).

The bessel functions come from the amos library (http://netlib.org/amos/) which is embedded in liboctave, so any numerical problems are inherited from that implementation.

Mike Miller <mtmiller>
Project Administrator
Sun 26 Jun 2016 05:55:18 AM UTC, original submission:  

Octave gives "besselj(1, 1e10)" as "NaN + NaNI".

Octave:

>> O = besselj(1, 10.^[5:12]')
O =
   0.00184675756288257 +                    0i
  -0.000725968356813763 +                    0i
  0.000236899240312058 +                    0i
  7.30639118155185e-05 +                    0i
  -5.21042264155388e-06 +                    0i
                 NaN +               NaNi
                 NaN +               NaNi
                 NaN +               NaNi

Symbolic pkg, converted to double:

>> S = double(besselj(1, sym(10).^[5:12]'))
S =
   0.00184675756288257
  -0.000725968356813763
  0.000236899240312058
  7.30639118155185e-05
  -5.21042264155388e-06
  -7.67650817568416e-06
  9.95266266795786e-07
  -7.91380268385044e-07

The relative accuracy is fine before blow-up:

>> abs(O-S)./abs(S)
ans =
  1.17416838493212e-16
  1.49345651543824e-16
  1.14415961302507e-16
                     0
  1.62565113336314e-16
                                NaN
                                NaN
                                NaN

(Ironically, I found this when I went looking for a function with good relative error behaviour even though it decays to zero!  viz., #48307)

Colin Macdonald <cbm>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by mtmiller (Posted a comment)
  • -email is unavailable- added by cbm (Submitted the item)
  • -email is unavailable- added by cbm
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

     

     

    Follow 5 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2018-03-27 rik5 StatusConfirmed => Fixed
        Open/ClosedOpen => Closed
    2016-11-16 mtmiller Operating SystemOther => Any
    2016-06-30 mtmiller StatusNone => Confirmed
    2016-06-26 cbm Carbon-Copy- => Added sebald

    Back to the top


    Powered by Savane 3.5