Problem definition
I'm finding some difficulties in the computation of the following integral \begin{equation*} I_C(y)\triangleq \int_{C(a)} \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right] \exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\text{ d}\nu_1\text{d}\nu_2 \end{equation*} where $y=[y_1\,\,y_2]'\in\mathbb{R}^2$ and the domain $C$ is a disc with radius $a$ centered in the origin, i.e. \begin{equation*}C(a)\triangleq \{\nu\in\mathbb{R}^2: \lVert \nu\rVert\leq a\}\end{equation*} I have in hand the solution of a simplified version of this problem and I'm wondering if I can exploit so solve the current problem. I'm also want to point out that a reasonable approximation (better than the one given by the simplified problem) would be acceptable - the important thing is to find an analytic expression that does not involve numerical procedures such as Monte Carlo.
Simplified problem
The simplified problem consists in the same Gaussian integral but defined over a simplified domain, where $C(a)$ is replaced by a square with edge $2a$ centered in the origin, i.e. \begin{equation*} I_S(y)\triangleq \int_{S(a)} \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right] \exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\text{ d}\nu_1\text{d}\nu_2\\ \end{equation*} and \begin{equation*}S(a)\triangleq[-a,a]^2 \end{equation*} in this case the solution is pretty straightforward and is given as follows. First, the integral can be splitted into two independent (with the same form!) integrals, \begin{equation*} I_S(y)\triangleq \left(\int_{-a}^a \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right]\text{ d}\nu_1\right)\left( \int_{-a}^a\exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\text{d}\nu_2\right)\\ \end{equation*} Now, focus on the first factor, say $I_1(y)$. This can be written as \begin{equation*}\begin{aligned} I_1(y)&\triangleq \int_{-a}^a \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right]\text{ d}\nu_1\\ &=\sqrt{\frac{\pi}{2}}\sigma_v \left[\frac{2}{\sqrt{\pi}} \int_{t(-a,y_1)}^{t(a,y_1)}\exp\left(-t^2\right)\text{ d}t\right] \end{aligned}\end{equation*} where I've used the change of integration variable \begin{equation*} \nu_1\mapsto t(\nu_1,y_1)\triangleq \frac{1}{\sqrt{2}\sigma_v}(\nu_1-y_1) \end{equation*} By introducing the error function \begin{equation*} \textrm{erf}(\alpha)\triangleq \frac{2}{\sqrt{\pi}} \int_0^\alpha \exp\left(-t^2\right)\text{ d}t \end{equation*} it turns out that \begin{equation*} I_1(y)= \sqrt{\frac{\pi}{2}} \sigma_v \left[\textrm{erf}(t(a,y_1))-\textrm{erf}(t(-a,y_1))\right] \end{equation*} and so, in conclusion, \begin{equation*} \boxed{ I_S(y)=\frac{\pi}{2} \sigma_v^2 \left[\textrm{erf}(t(a,y_1))-\textrm{erf}(t(-a,y_1))\right] \left[\textrm{erf}(t(a,y_2))-\textrm{erf}(t(-a,y_2))\right]} \end{equation*}
Current problem
Despite the integrand is still the same, the change of the integration domain causes a substancial difference in the difficulty of the problem. Indeed, it seems not possible to split the double integral into the product of two single integrals due to the fact that the integration extrema are not indepedent \begin{equation*}\begin{aligned} I_C(y)&= \int_{-a}^a \int_0^{\pm\sqrt{a^2-\nu_2^2}} \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right] \exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\text{ d}\nu_1\text{d}\nu_2\\ &=\int_{-a}^a \exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\left(\int_0^{\pm\sqrt{a^2-\nu_2^2}} \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right] \text{ d}\nu_1\right)\text{d}\nu_2\\ \end{aligned} \end{equation*} where I'm using the symbol $\pm$ on one extrema to denote the fact that there are two different integral to be evaluated, more precisely I mean \begin{equation*}\int_0^{\pm x}\triangleq\int_0^x+\int_0^{-x}\end{equation*}
The inner integral is pretty easy because we can use the previous solution found in the simplified problem, obtaining \begin{equation*} \begin{aligned} \int_0^{\pm\sqrt{a^2-\nu_2^2}} \exp\left[-\frac{1}{2}\left(\frac{\nu_1-y_1}{\sigma_v}\right)^2\right] \text{ d}\nu_1 &= \sqrt{\frac{\pi}{2}} \sigma_v \left[\textrm{erf}\left(t_+(\nu_2,y_1)\right)-\textrm{erf}\left({t}(0,y_1)\right)\right]+\\ &\qquad\qquad \sqrt{\frac{\pi}{2}} \sigma_v \left[\textrm{erf}\left(t_-(\nu_2,y_1)\right)-\textrm{erf}\left({t}(0,y_1)\right)\right] \end{aligned} \end{equation*} where \begin{equation*} t_\star(\nu_2,y_1)\triangleq t\left(\star\sqrt{a^2-\nu_2^2},y_1\right) \qquad \star=+,- \end{equation*}
Now let's define the ugliest function in the world, i.e. \begin{equation*}\begin{aligned} u(\nu_2,y_1)&\triangleq \left[\textrm{erf}\left(t_+(\nu_2,y_1)\right)-\textrm{erf}\left({t}(0,y_1)\right)\right]+ \left[\textrm{erf}\left(t_-(\nu_2,y_1)\right)-\textrm{erf}\left({t}(0,y_1)\right)\right]\\ &=\textrm{erf}\left(t_+(\nu_2,y_1)\right)+\textrm{erf}\left(t_-(\nu_2,y_1)\right)-2\,\textrm{erf}\left({t}(0,y_1)\right) \end{aligned} \end{equation*}
so that the searched integral gets the form \begin{equation*} I_C(y)= \sqrt{\frac{\pi}{2}} \sigma_v\int_{-a}^a u(\nu_2,y_1)\,\exp\left[-\frac{1}{2}\left(\frac{\nu_2-y_2}{\sigma_v}\right)^2\right]\text{ d}\nu_2 \end{equation*}
At this point I get stuck because I have no idea how to treat this bad boy here.
Observation
Let's reduce the expectations and let's try to find an approximation for the integral above. A trivial observation is that, written in that form, the integral can be seen as a weighted sum of $u$ where the weight function is a Gaussian kernel centered in $y_2$ with spread $\sigma_v$. Hence, in such sum the more important values $\nu_2$ are the one near $y_2$. Now, if we consider the limit case where $\sigma_v=0$, the Gaussian kernel collapses to a delta concentred in $\nu_2=y_2$, and so only the point $\nu_2=y_2$ matters in the total sum \begin{equation*} \lim_{\sigma_v \to 0}I_C(y)= \sqrt{\frac{\pi}{2}} \sigma_v\int_{-a}^a u(\nu_2,y_1)\,\delta(\nu_2-y_2)\text{ d}\nu_2= \begin{cases} \sqrt{\frac{\pi}{2}}\,\sigma_\, u(y_2,y_1) & \text{if } -a\leq y_2 \leq a\\ 0 & \text{otherwise} \end{cases} \end{equation*} I think that this result can be considered as a "zero-order" approximation of what I'm trying to compute. The problem is that in my applications $\sigma_v$ can be large (with respect to the domain of integration given here by $a$), so that my "zero-order" approximation is pretty rough.
Questions
- Is it possible to solve my problem exactly? Maybe by employing some exotic function?
- Is it possible to improve my approximation by taking into account a better approximation of the Gaussian kernel? Maybe it can be useful to consider a Talyor approximation of such term.
- Any further ideas or suggestions?
Thank you so much in advance.