[SciPy-user] special.jn(3,x) became inexact in recent versions

Robert Kern robert.kern at gmail.com
Sat Jan 20 19:26:52 CST 2007


Amit Aronovitch wrote:

> Possible Explanation:
>   Older releases of cephes used the original names (j0,jn,...), while
> the newer release contains macros converting these names to "cephes_jn" etc.
>   It is possible that the older cephes module did not really use the
> cephes functions, but the ones from glibc instead (glibc contains many
> mathematical functions, including a "jn", which does work for large x).

Excellent detective work! This seems quite likely to me.

>  I wonder how do cephes functions compare to glibc ones
> performance-wise. Might be a good idea to generate ufunc wrappers to
> those (on platforms where glibc is available).

And several other C runtimes as well. For some reason I was under the impression
that at least j0 was part of the C99 standard, but I do not see any of the
Bessel functions in the standard itself (section 7.12):

  http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1124.pdf

It would be possible to check for each of these and construct a configuration
header with HAVE_LIBM_JN #defines. The #define hackery that renames jn to
cephes_jn would be conditional on the HAVE_LIBM settings.

I hesitate to add more to the configuration steps, though. Better all around
would be for us to find (or write) a suitably licensed implementation that's better.

-- 
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless enigma
 that is made terrible by our own mad attempt to interpret it as though it had
 an underlying truth."
  -- Umberto Eco


More information about the SciPy-user mailing list