Accurate recursive generation of spherical Bessel and Neumann functions for a large range of indices

Partial wave expansion is encountered in a wide range of problems in theoretical physics. Typically, spherical Bessel functions are encountered in scattering of nuclei, atoms, and molecules. A routine is described for computing spherical Bessel and Neumann functions j_l( x) and n_l( x), respectively, by backward and forward recurrence relations. The method is fast, very accurate, not code intensive, and easily adopted for a personal computer. These features make the code well suited for the physics student research environment.