Hi, Your code tries to to evaluate z = 1263309.3633394379 + 101064.74910119522j jv(536, z) # -> (inf+inf*j) In reality, this number is not infinite, but jv(536, z) == -2.3955170861527422e+43888 + 9.6910119847300024e+43887 These numbers (~ 10^43888) are too large for the floating point numbers that computers use (maximum ~ 10^308). This is why you get infinities and NaNs in the result. The same is true for the spherical Bessel functions. You will not be able to do this calculation using any software that uses only floating point numbers (Scipy, Matlab, ...). You need to use analytical properties of your problem to get rid of such large numbers. Alternatively, you can use arbitrary precision numbers. Python has libraries for that: http://code.google.com/p/mpmath/ By the way, the proper place for this discussion is the following mailing list: http://mail.scipy.org/mailman/listinfo/scipy-user -- Pauli Virtanen