Asymptotic expansion of double integral Define
$$G(\theta) = \int\limits_0^\infty \int\limits_0^{2\pi} \frac{r\,e^{-r^2/2t}}{\sqrt{1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2}} \mathrm{d} r \,\mathrm{d} \varphi$$
Clearly, for $\theta=0$, this does not converge. However, I would like to obtain an asymptotic expansion for $\theta\searrow 0$.
How could I approach this problem? The function
$$F(r, \theta) = \int_0^{2\pi} \frac{\mathrm{d}\varphi}{\sqrt{1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2}}$$
looks somewhat like an elliptic integral, but I have no experience in dealing with these and I have no idea how this could help.
Any thoughts on how to approach this problem?
 A: Part 1: Expression of inner integral as an elliptic integral
You are right about the inner integral which you name $F{(r,\theta)}$ resembling an elliptic integral, and indeed it can expressed (with some difficulty) in terms of an elliptic integral of the first kind. For this reason I would avoid naming your function with a capital $F$, since this is usually reserved for a specific elliptic integral. In my work below, I will use a lowercase $f$ instead (and I've also elected to change $r$ to $\rho$ out of stylistic preferences). Argument conventions for defining elliptic integrals vary, but I will here be adopting the one used by Gradshteyn & Ryzhik since I will be invoking identities from their tome as part of my answer; they define:
$$F{\left(\varphi,k\right)}:=\int_{0}^{\varphi}\frac{\mathrm{d}\alpha}{\sqrt{1-k^2\sin^2{\alpha}}}=\int_{0}^{\sin{\varphi}}\frac{\mathrm{d}x}{\sqrt{(1-x^2)(1-k^2x^2)}},~~~\text{where }k^2<1;\\
K{\left(k\right)}:=F{\left(\frac{\pi}{2},k\right)}.$$
With that long-winded preamble out of the way, let's get to work solving the problem at hand. For real parameters $\theta,\rho\in\mathbb{R}\setminus\{n\pi|n\in\mathbb{Z}\}$ such that $\theta\neq\rho$,
$$\begin{align}
f{(\theta,\rho)}
&:=\int_{0}^{2\pi}\frac{\mathrm{d}\varphi}{\sqrt{1-\left(\cos{\theta}\cos{\rho}+\sin{\theta}\sin{\rho}\cos{\varphi}\right)^2}}\\
&=\int_{0}^{\pi}\frac{2\,\mathrm{d}\varphi}{\sqrt{1-\left(\cos{\theta}\cos{\rho}+\sin{\theta}\sin{\rho}\cos{\varphi}\right)^2}}\\
&=\int_{1}^{-1}\frac{-2\,\mathrm{d}x}{\sqrt{1-x^2}\sqrt{1-\left(\cos{\theta}\cos{\rho}+x\sin{\theta}\sin{\rho}\right)^2}};~~~\cos{\varphi}\rightarrow x\\
&=\int_{-1}^{1}\frac{2\,\mathrm{d}x}{\sqrt{1-x^2}\sqrt{1-\left(\cos{\theta}\cos{\rho}+x\sin{\theta}\sin{\rho}\right)^2}}\\
&=\int_{-1}^{1}\frac{2\left|\csc{\theta}\right|\left|\csc{\rho}\right|\,\mathrm{d}x}{\sqrt{1-x^2}\sqrt{\csc{\theta}\csc{\rho}-\left(\cot{\theta}\cot{\rho}+x\right)^2}}\\
&=\int_{-1}^{1}\frac{2\left|\csc{\theta}\right|\left|\csc{\rho}\right|\,\mathrm{d}x}{\sqrt{1-x^2}\sqrt{\left(\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}-x\right)\left(\csc{\theta}\csc{\rho}+\cot{\theta}\cot{\rho}+x\right)}}\\
&=\small{2\left|\csc{\theta}\right|\left|\csc{\rho}\right|\int_{-1}^{1}\frac{\mathrm{d}x}{\sqrt{1-x^2}\sqrt{\left(\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}-x\right)\left(\csc{\theta}\csc{\rho}+\cot{\theta}\cot{\rho}+x\right)}}}.\\
\end{align}$$
Since the plan is to ultimately consider the limit as $\theta\searrow0$, we may go ahead and restrict $\theta$ to the interval $0<\theta<\pi$. Because we shall need to integrate over $\rho$ from $\rho=0$ to $\rho\to+\infty$, it's not necessarily obvious at first that we can likewise restrict $\rho$ to the same interval, but since it is true that $f{(\theta,\rho+\pi)}=f{(\theta,\rho)}$ it follows that we can suppose $0<\rho<\pi$ without loss of generality.
Note that with the further restrictions $0<\theta,\rho<\pi$, we have the inequalities
$$\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}>1\\
\iff 1-\cos{\theta}\cos{\rho}>\sin{\theta}\sin{\rho}\\
\iff 1>\cos{\theta}\cos{\rho}+\sin{\theta}\sin{\rho}\\
\iff 1>\cos{\left(\theta-\rho\right)},$$
and
$$-1>-\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}\\
\iff -\sin{\theta}\sin{\rho}>-1-\cos{\theta}\cos{\rho}\\
\iff \cos{\theta}\cos{\rho}-\sin{\theta}\sin{\rho}>-1\\
\iff \cos{\left(\theta+\rho\right)}>-1.$$
In order to continue, I shall merely quote without proof a definite integral from Gradshteyn's Table of Integrals, Series, and Products, specifically proposition 3.147(4) on p.324 which states:

Given $a>b\ge u>c>d$,
  $$\int_{c}^{u}\frac{\mathrm{d}x}{\sqrt{(a-x)(b-x)(x-c)(x-d)}}=\frac{2}{\sqrt{(a-c)(b-d)}}F{(\delta,q)},$$
  where $\delta=\arcsin{\sqrt{\frac{(b-d)(u-c)}{(b-c)(u-d)}}},~~q=\sqrt{\frac{(b-c)(a-d)}{(a-c)(b-d)}}$.

Letting $b=u=1$, $c=-1$, $a=\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}$, and $d=-\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}$, we have $\delta=\frac{\pi}{2}$ and 
$$\begin{align}
q
&=\sqrt{\frac{2(a-d)}{(1+a)(1-d)}}\\
&=\sqrt{\frac{2((\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho})-(-\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}))}{(1+(\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}))(1-(-\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}))}}\\
&=\sqrt{\frac{4\csc{\theta}\csc{\rho}}{(1+\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho})(1+\csc{\theta}\csc{\rho}+\cot{\theta}\cot{\rho})}}\\
&=\sqrt{\frac{4\sin{\theta}\sin{\rho}}{\left(\sin{\theta}+\sin{\rho}\right)^2}}\\
&=\frac{2\sqrt{\sin{\theta}\sin{\rho}}}{\left|\sin{\theta}+\sin{\rho}\right|}.\\
\end{align}$$
Thus, for $0<\theta\neq\rho<\pi$, we have:
$$\begin{align}
f{(\theta,\rho)}
&=\small{2\left|\csc{\theta}\right|\left|\csc{\rho}\right|\int_{-1}^{1}\frac{\mathrm{d}x}{\sqrt{1-x^2}\sqrt{\left(\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho}-x\right)\left(\csc{\theta}\csc{\rho}+\cot{\theta}\cot{\rho}+x\right)}}}\\
&=2\left|\csc{\theta}\right|\left|\csc{\rho}\right|\cdot\frac{2\,F{\left(\frac{\pi}{2},\frac{2\sqrt{\sin{\theta}\sin{\rho}}}{\left|\sin{\theta}+\sin{\rho}\right|}\right)}}{\sqrt{(1+\csc{\theta}\csc{\rho}-\cot{\theta}\cot{\rho})(1+\csc{\theta}\csc{\rho}+\cot{\theta}\cot{\rho})}}\\
&=\frac{4\left|\csc{\theta}\right|\left|\csc{\rho}\right|}{\sqrt{(\csc{\theta}+\csc{\rho})^2}}\,K{\left(\frac{2\sqrt{\sin{\theta}\sin{\rho}}}{\left|\sin{\theta}+\sin{\rho}\right|}\right)}\\
&=\frac{4}{\sin{\theta}+\sin{\rho}}\,K{\left(\frac{2\sqrt{\sin{\theta}\sin{\rho}}}{\sin{\theta}+\sin{\rho}}\right)}\\
&=4\csc{\theta}\,\frac{1}{1+\frac{\sin{\rho}}{\sin{\theta}}}\,K{\left(\frac{2\sqrt{\frac{\sin{\rho}}{\sin{\theta}}}}{1+\frac{\sin{\rho}}{\sin{\theta}}}\right)}\\
&=\begin{cases}
4\csc{\theta}\,K{\left(\frac{\sin{\rho}}{\sin{\theta}}\right)};~~~\rho<\theta\\
4\csc{\rho}\,K{\left(\frac{\sin{\theta}}{\sin{\rho}}\right)};~~~\theta<\rho.\\
\end{cases}\\
\end{align}$$

Part 2: Representation of double integral as a single integral over bounded domain
For $t,\theta\in\mathbb{R}^{+}$ with $0<\theta<\pi$,
$$\begin{align}
G{\left(t,\theta\right)}
&:=\int_{0}^{\infty}\mathrm{d}r\int_{0}^{2\pi}\mathrm{d}\varphi\,\frac{r\,e^{-\frac{r^2}{2t}}}{\sqrt{1-\left(\cos{\theta}\cos{r}+\sin{\theta}\sin{r}\cos{\varphi}\right)^2}}\\
&=\int_{0}^{\infty}\mathrm{d}r\,r\,e^{-\frac{r^2}{2t}}\int_{0}^{2\pi}\frac{\mathrm{d}\varphi}{\sqrt{1-\left(\cos{\theta}\cos{r}+\sin{\theta}\sin{r}\cos{\varphi}\right)^2}}\\
&=\int_{0}^{\infty}\rho\,e^{-\frac{\rho^2}{2t}}\,f{\left(\theta,\rho\right)}\,\mathrm{d}\rho\\
&=\sum_{n=0}^{\infty}\int_{n\pi}^{(n+1)\pi}r\,e^{-\frac{r^2}{2t}}\,f{\left(\theta,r\right)}\,\mathrm{d}r\\
&=\sum_{n=0}^{\infty}\int_{0}^{\pi}\left(\rho+n\pi\right)\,e^{-\frac{\left(\rho+n\pi\right)^2}{2t}}\,f{\left(\theta,\rho+n\pi\right)}\,\mathrm{d}\rho;~~~r\rightarrow\rho+n\pi\\
&=\sum_{n=0}^{\infty}\int_{0}^{\pi}\left(\rho+n\pi\right)\,e^{-\frac{\left(\rho+n\pi\right)^2}{2t}}\,f{\left(\theta,\rho\right)}\,\mathrm{d}\rho\\
&=\int_{0}^{\pi}\sum_{n=0}^{\infty}\left[\left(\rho+n\pi\right)\,e^{-\frac{\left(\rho+n\pi\right)^2}{2t}}\right]\,f{\left(\theta,\rho\right)}\,\mathrm{d}\rho\\
&=:\int_{0}^{\pi}h{\left(t,\rho\right)}\,f{\left(\theta,\rho\right)}\,\mathrm{d}\rho,\\
\end{align}$$
where in the last line we've defined the new auxiliary function $h{\left(t,\rho\right)}$ by the infinite series
$$h{\left(t,\rho\right)}:=\sum_{n=0}^{\infty}\left(\rho+n\pi\right)\,e^{-\frac{\left(\rho+n\pi\right)^2}{2t}}.$$
I'm not completely sure, but I believe the function $h{\left(t,\rho\right)}$ can be written in terms of elliptic theta functions somehow.
A: The following is not a proof, however numerics shows that the expressions seem to be correct:
In your integral for $G(\theta)$, the denominator $\sqrt{1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2}$ is responsible for the divergence. 
In fact it becomes minimal for $r \to r_n = n\pi + \theta \cos(\varphi)$ (to first order in $\theta$), where it approaches the value $\theta |\sin(\varphi)|$. We expand the denominator close to this minimal value and have the approximation
$$ 1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2 \approx \theta^2 \sin^2(\varphi) + (r- r_n)^2.$$ As the denominator is only large close to the points $r=r_n$, we can put for a first approximation the corresponding value in the numerator.
It is thus useful to write the integral $G(\theta)$ as $G(\theta) =\int\limits_0^{2\pi}\! d\varphi \sum_n g_n$ with 
$$ \begin{align}
 g_0 &= \int_0^{\pi/2}\!dr  \frac{r\,e^{-r^2/2t}}{\sqrt{1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2}}\\
 g_{n\geq1} &= \int_{(n-1/2) \pi}^{(n+1/2)\pi}\!dr  \frac{r\,e^{-r^2/2t}}{\sqrt{1-(\sin\theta\sin r \cos\varphi + \cos \theta \cos r)^2}}.
\end{align} $$
The individual integrals, we evaluate with logarithmic accuracy
$$
\begin{align}
 g_{n\geq 1} &\approx \int_{-\pi/2}^{\pi/2}\!d\delta \frac{r_n\,e^{-r_n^2/2t}}{ \sqrt{ \theta^2 \sin^2(\varphi) + \delta^2} } 
 \approx 2
 r_n\,e^{-r_n^2/2t} \int_{ \theta \sin(\varphi)}^1 \frac{1}{\delta}=
2
r_n\,e^{-r_n^2/2t} \left|\ln\theta |\sin(\varphi)| \right|\\
&\approx 2
r_n\,e^{-r_n^2/2t} (|\ln\theta| + |\ln \sin(\varphi)|)
\end{align}
$$
For $n=0$, the value is only half of the one above as the range of integration is halved.
In total, we obtain the estimate
$$ G(\theta) \sim 2\pi \sum_{n\geq1} n e^{-\pi^2 n^2/2 t} \int_0^{2\pi}\!d\varphi
\, (|\ln\theta| + |\ln \sin(\varphi)|)
\sim 4\pi^2 |\ln\theta| \sum_{n\geq 1} n e^{-\pi^2 n^2/2 t}. $$
The sum can only be evaluated in terms of elementary functions in the limit $t\ll1$ and $t\gg1$. In the former case, we have (taking only the term with $n=1$)
$$ G(\theta) \sim 4 \pi^2 e^{-\pi^2/2 t} |\ln \theta|, \qquad t \ll 1.$$
In the later case, we have
$$ G(\theta) \sim  4 \pi^2 |\ln\theta| \int_0^\infty\!dn\,ne^{-\pi^2 n^2/2 t}
=   4 t |\ln\theta|, \qquad t \gg 1.$$
