I need to compute (or best approximate?) the following integral
$$\int_0^\infty \int_0^\infty (1 + \alpha u)^{-1}(1 + v)^{-1} \Phi\left(\frac{\beta}{\sqrt{\gamma + uv}}\right) \text{d}u \text{d}v,\qquad \gamma > 0, \quad \beta \in \mathbb{R}, \quad \alpha \in \mathbb{N}^*,$$ where $\Phi(\cdot)$ is the standard normal cumulative distribution function and $\alpha$ is typically a large natural number. What strategy would you advise for this?
I don't know if this is helpful but I can also reformulate my problem (up to multiplicative constants) as computing the following expectation:
$$\text{E}\left\{ \Phi\left(\frac{\beta}{\sqrt{\gamma + X^2Y^2}}\right)\right\},\qquad \gamma > 0, \quad \beta \in \mathbb{R},$$ where $$X \sim \text{Half-Cauchy}(0,1), \qquad Y \sim \text{Half-Cauchy}(0,\alpha^{-1/2}), \quad \alpha \in \mathbb{N}^* \, (\text{large}).$$
EDIT: I made a mistake in the first integral, it should be:
$$\int_0^\infty \int_0^\infty (1 + \alpha x^2)^{-1}(1 + y^2)^{-1} \Phi\left(\frac{\beta}{\sqrt{\gamma + x^2y^2}}\right) \text{d}x \text{d}y,\quad \gamma > 0, \beta \in \mathbb{R}, \alpha \in \mathbb{N}^*\, \text{(large)},$$ which now corresponds to the above expectation.
EDIT 2: Another way to solve my problem would be to characterise the difference (hopefully small?) $$\text{E}\left\{ \Phi\left(\frac{\beta}{\sqrt{\gamma + X^2Y^2}}\right)\right\} - \Phi\left(\frac{\beta}{\sqrt{\gamma}}\right),\qquad \gamma > 0, \quad \beta \in \mathbb{R},$$ would it be sensible to approximate the difference $$\Phi\left(\frac{\beta}{\sqrt{\gamma + X^2Y^2}}\right) - \Phi\left(\frac{\beta}{\sqrt{\gamma}}\right),$$ with a Taylor series in zero or is this unreasonable given that the half-Cauchy distributions are heavy-tailed?