On Wed, 17 Jun 2009, Klein, Christoph wrote:
I'm working on the bessel code right now. The attached patch
implements a bessel function which can be called like this:
bessel(<type>, <order>, <value>) example bessel("K",0.5,1).
Wow! How great to get such a well worked-out patch! I have
added to CVS most of what you submitted, in a slightly
trimmed-down version that (for now) only accepts a scalar value
for the order of the Bessel function.
the function works for J,I,K,Y bessel functions. For K bessel
there was no support for fractional orders in cephes, thats why
there is additonal code for these cases. the fractional order
for k and negative orders for I work reasonably well in the
range lambda<10 and value<20. The attached test script checks
against the implementation of R.
Nice! I'll check against libgsl too. On an initial scan it seems
that most things are fine but the Bessel 'K' function goes a bit
strange for large values of X (as you say, values < 20 are OK).
Allin Cottrell