Bessel functions of the first kind, \(J_l\) where \(l\) is an integer, are found in scipy.special.jn. This is a function that takes two arguments: the first one is the value of \(l\) and the second one is the value of the point at which to evaluate the functions. If one is interested in \(J_0\) or \(J_1\) only, then there are faster implementations of these functions as j0 and j1.
If the goal is to find the zeros of the Bessel function, then rather than using a root finding routine (such as scipy.optimize.newton) it is more efficient to call scipy.special.jn_zeros. This function also takes two arguments, the first one is the value of \(l\) and the second one is the number of zeros to display (it starts at the first positive zero).
For Bessel functions of half-integer order, so called spherical Bessel function, the situation is a little bit more complicated. The appropriate function is scipy.special.sph_jn. To compute \(J_{l + 1/2}\), the first argument to sph_jn is the value of \(l\). The second argument is the value of the point at which to evaluate the function. Where it gets complicated is in the return value of the function. sph_jn returns a pair of arrays. The first array contains the values of the function, and the second array contains the values of the derivative. Each array has length \(l + 1\) and contains a value for all spherical Bessel functions with \(l\) smaller than or equal to the supplied value.
It is useful to wrap this function into something friendlier:
def j3o2(x): """Compute the spherical Bessel function J_{3/2}.""" return sph_jn(1, x)[0][1] def j5o2(x): """Compute the spherical Bessel function J_{5/2}.""" return sph_jn(2, x)[0][2]
To obtain the zeros of the spherical Bessel function, one has to call a root finding procedure. For instance, scipy.optimize.newton does the job quite well.
John D. Cook has a great page about Bessel functions and another one about how to work with the implementation in Scipy.