Bessel function integration for zero argument

[jmsull]

Not sure this will ever actually be an issue in practice, but right now the bessel function interpolator is constructed using a point at x =0 (were we have j_\ell(x)) when \eta = \eta_final.

At least for the polarization source function, there is a problem because j_2(x)/x^2 |x=0 is not zero, but the way the source function integration is done is by calling the bessel interpolator at x=0 and then multiplying it into the source function.

I caught a NaN in the source function for x=0, but the bessel interpolator returning j_2(0)=0 will still be wrong.

This might be solved by the new Bessel integration?