fix for scaled besseli at negative non-integer orders

John W. Eaton jwe at bevo.che.wisc.edu
Tue Oct 28 10:46:38 CDT 2008


On 28-Oct-2008, Brian Gough wrote:

| In the current tip there is a problem with the scaling factor for
| besseli with negative non-integer order in liboctave/lo-specfun.cc.  
| 
| Computing the ratio of scaled and unscaled besseli(alpha,2) should
| give the same value for all orders but gives a different result for
| orders like -3.5:
| 
|    $ ./run-octave
|    GNU Octave, version 3.1.51+
|    ....
|    octave:1> for alpha = [1,2,3,3.5,-1,-2,-3,-3.5] ; r=besseli(alpha,2,1)/besseli(alpha,2); disp(r); end
|     0.13534  # exp(-2) correct ratio
|     0.13534
|     0.13534
|     0.13534
|     0.13534
|     0.13534
|     0.13534
|     8.6239   # should also be 0.13534 = exp(-2)
| 
| The current code in zbesi() computes besseli(-alpha,z) as 
| 
|   zbesi(alpha,z) + (2/pi)*sin(pi*alpha)*zbesk(alpha,z) 
| 
| This is correct in the unscaled case, but needs an extra term to
| compensate for the different scaling factors of zbesi and zbesk in the
| scaled case.
| 
| The scaling factor of zbesi() is si=exp(-abs(real(z))) but for zbesk it
| is sk=exp(z).  Hence the formula in the scaled case needs to be
| 
|   zbesi(alpha,z) + (si/sk)*(2/pi)*sin(pi*alpha)*zbesk(alpha,z) 
| 
| i.e. an extra factor of exp(-z-abs(real(z))).  I've attached three
| patches which
| 
| 1) add all the different scaling factors for besselj,y,i,k,h to the
| documentation
| 
| 2) add tests for the bessel functions for combinations of
| odd,even,fractional order and real and complex argument
| 
| 3) correct the calculation of zbesi for negative order with scaling.
| I've put the same correction into cbesi as well,

I applied the changesets.

| although the single
| precision calls to besseli don't currently seem to work (they always
| return NaN for me).

I think we need the following changes because the argument lists for
the CBES{H,I,J,K,Y} subroutines from the Amos library are different
from the corresponding ZBES{H,I,J,K,Y} subroutines.

Thanks,

jwe

-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: diffs
Url: https://www-old.cae.wisc.edu/pipermail/bug-octave/attachments/20081028/3a20d32d/attachment-0001.ksh 


More information about the Bug-octave mailing list