I'm trying to carry out the following integral numerically

$$\int_{r_\mathrm{in}}^{r_\mathrm{out}} \Sigma\left(r'\right) \frac{r'}{r} \left( \frac{1}{r-r'}\, E(L) + \frac{1}{r+r'}\, K(L) \right) \mathrm{d}r',$$


$$L = \frac{4 r r'}{(r+r')^2},$$

and $K$ and $E$ are the complete elliptical integrals of the first and second kind. $r$, $r_\mathrm{in}$, and $r_\mathrm{out}$ are all $>0$, $r_\mathrm{in}\leq r \leq r_\mathrm{out}$ and $\Sigma(r)$ is some positive, well behaved function (no singularities, or discontinuities, ...).

In the end, I need to evaluate the result numerically for a given $\Sigma(r')$, but my usual python integrators scipy.integrate.quad or scipy.integrate.romberg fail to converge or give unreliable results due to the two singularities:

  • $1/(r-r')$ in the first summand
  • $K(L)\rightarrow \infty$ as $r'\rightarrow r$ in the second summand ($L \rightarrow 1$)

Can anyone help me to evaluate this numerically or to rewrite it such that the usual methods don't fail?

  • $\begingroup$ It looks like even in the benign case of $\Sigma(r')=1$, $K(L)=0$, $E(L)=1$, the value of the integral over $(r-\epsilon, r+\epsilon)$ isn't well-defined (compare with $\int_{-\epsilon}^{\epsilon} 1/x\,\mathrm{d}x$). $\endgroup$ – Nico Schlömer Oct 6 '13 at 17:09

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy