Can anyone explain why:
$$\int_{0}^{2\pi}\exp(\beta(2r_{1}r_{2}\cos(\theta)))\,\mathrm d\theta=2\left(\int_{0}^{\pi}\exp(\beta(2r_{1}r_{2}\cos(\theta)))\,\mathrm d\theta\right)$$
where the right hand side is the modified Bessel function of the first kind evaluated at $2r_{1}r_{2}\beta$ - here's the function for convenience:
$$I_{0}(x)=\frac{1}{\pi}\int_{0}^{\pi}\exp(x\cos\theta)\,\mathrm d\theta$$
Look at a plot of $\cos{\theta}$: $\cos{\theta}$ takes on exactly the same set of values in $[0,\pi]$ as it does in $[\pi,2 \pi]$. Therefore, the integrals over these regions are equal and the equality is true.