I'm looking for an open-source (to use and learn from) software which computes Bessel functions of integer order of real argument to double precision the fastest among all such implementations. Currently I've tried Boost.Math and GSL. From these GSL appeared to be faster for smaller arguments and much slower for larger ones.

Are there any implementations specially designed to be very fast for the whole working ranges of arguments and orders, but still not neglecting precision?

  • $\begingroup$ The NIST Digital Library of Mathematical Functions has a section on computation of Bessel Functions with references to various methods that can be used. It doesn't have benchmarks of performance and accuracy though. See dlmf.nist.gov $\endgroup$ Commented Jun 10, 2015 at 21:07
  • $\begingroup$ Have you profiled to determine that calculating Bessel functions is your performance bottleneck? $\endgroup$ Commented Jun 10, 2015 at 21:16
  • $\begingroup$ @GeoffOxberry I think switching from one implementation to another and seeing overall performance increase is a good enough indicator of this?.. $\endgroup$
    – Ruslan
    Commented Jun 11, 2015 at 3:03

3 Answers 3


The appropriate and fastest library depends on several things. Which Bessel functions (only J, Y & Hankel or modified Bessel functions I & K too), for which types of arguments (real or complex, integer, fractional or general order)?

Amos's libraries are written in Fortran-77 (there are Fortran-90 coverted versions of TOMS 644 on a mirror of Alan Miller's fortran page: http://jblevins.org/mirror/amiller/) and are numerically accurate and fast for complex argument of all Bessel functions, but they are difficult to understand the source code. He produced several versions of his library, but only the one included in SLATEC (http://www.netlib.org/slatec/) is actually open source, the other versions have a copyright associated with the transactions on mathematical software journal.


mpmath (http://mpmath.org/) uses a totally different and more general approach, by expressing all Bessel functions as special cases of hypergeometric functions. Mpmath is also arbitrary precision, so this is probably overkill if you just need double precision. This is more general, but slower. The author of mpmath has a C-based library called arb (http://fredrikj.net/arb/) that implements some of the functions in mpmath but much more quickly (since it is C-based rather than Python-based). Since arb and mpmath are arbitrary precision, they might be good for benchmarking, but probably not for speed testing.

I haven't personally worked much with the Gnu Scientific Library, but I believe it is fast and accurate for what functions it covers. It may not have all the functions you need. It would probably be a better example of well-written modern code, compared to the Amos libraries.


This is really going to require some testing with regards to speed, but here are some examples I remember coming across that aren't Boost or GSL:

Also, a related SO question here with links to other C implementations

  • 1
    $\begingroup$ I would add mpmath, that is used in Sympy $\endgroup$
    – nicoguaro
    Commented Jun 10, 2015 at 21:02

You can try the Colt library for Java:


In spite of having had trouble setting them up, I have been able to run the GNU GMP and MPFR libraries. MPFR has some Bessel functions as well. See the documentation below in the link:


One nice thing about GSL is that you can use functions to get the error of your calculation. If you're calculating functions in Fortran, you might want to look at FGSL, the Fortran version of GSL.

  • $\begingroup$ You haven't said anything about performance, and this was the main subject of the question. $\endgroup$
    – Ruslan
    Commented Jul 26, 2015 at 20:00
  • $\begingroup$ I have done some calculations concerning the Riemann Zeta function using the power-series expansion for positive arguments. I performed the calculations using some Java code I developed with BigDecimal and using MPFR. MPFR instantaneously gave me the answer for the Riemann Zeta function to an arbitrary precision. The Java code I developed didn't work too well. Maybe an unfair comparison, as I wasn't using optimized routines for calculating pi. $\endgroup$ Commented Jul 26, 2015 at 20:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.