Double Gaussian definite integral with one variable limit

I am interested in solving a definite double integral of the following form:

\begin{align} f(a,b) &= \int_0^\infty \exp\Big(\frac{-x^2}{2a}\Big)\int_{x}^{\infty} \exp\Big(\frac{-y^2}{2b}\Big) dydx\\ &= \int_0^\infty \exp\Big(\frac{-x^2}{2a}\Big)\text{erfc}\Big(\frac{x}{\sqrt{2b}}\Big)dx, \end{align}

for $a,b>0$, where erfc is the complementary error function. One potential way to go would be to use a power-series expansion (e.g., see answer by robjohn to this question or this paper), but I'm finding that a bit difficult to follow. I'm wondering if anybody has any ideas about ways to get an approximate answer. For now, I'm just trying to see if I can fit the numerical solution with a function of $a$ and $b$.

Assuming $a,b$ to be real and positive $$I=\int_0^\infty e^{-\frac{x^2}{2 a}}\, \text{erfc}\left(\frac{x}{\sqrt{2b} }\right)\,dx=\sqrt{2a} \int_0^\infty e^{-y^2} \text{erfc}\left({\sqrt{\frac a b}\, y}\right)\,dy$$ Now, using series $$\text{erfc}(z)=1-\frac{2}{\sqrt{\pi}}\sum_{n=0}^\infty\frac{(-1)^n }{n!\, (2n+1)}z^{2n+1}$$ integrate termwise using the fact that $$\int_0^\infty e^{-z^2}z^{2n+1}\,dz=\frac{n!}2$$ to get $$I=\sqrt{2a}\left(\sqrt{\frac \pi 2}-\frac 1{\sqrt \pi}\sum_{n=0}^\infty\frac{(-1)^n }{ (2n+1)}\left(\sqrt{\frac a b}\right)^{2n+1}\right)$$ $$I=\sqrt{2a}\left(\sqrt{\frac \pi 2}-\frac 1{\sqrt \pi}\tan ^{-1}\left(\sqrt{\frac{a}{b}}\right)\right)=\sqrt{\frac{2a}{\pi }} \tan ^{-1}\left(\sqrt{\frac{b}{a}}\right)$$
My answer has been disclaimed by @Did; however, I still think it is correct for the last integral written in the post. Were there is a mistake is at the very start of the post since $$\int_{x}^{\infty} \exp\Big(\frac{-y^2}{2b}\Big) dy=\color{blue}{\sqrt{\frac{b \pi} 2 }}\text{erfc}\left(\frac{x}{\sqrt{2b} }\right)$$
So, $$\color{blue}{\sqrt{\frac{b \pi} 2 }} \sqrt{\frac{2a}{\pi }} \tan ^{-1}\left(\sqrt{\frac{b}{a}}\right)=\sqrt{a b}\tan ^{-1}\left(\sqrt{\frac{b}{a}}\right)$$ which is Did's answer.
Consider $(X,Y)$ i.i.d. standard normal, then $$f(a,b)=2\pi\sqrt{ab}P(\sqrt bY>\sqrt aX>0)$$ The distribution of $(X,Y)$ is rotationally invariant hence, for every $c>0$, $$P(Y>cX>0)=\frac{\pi/2-\vartheta}{2\pi}$$ where $\vartheta$ in $(0,\pi/2)$ solves $$\tan\vartheta=c$$ Thus, $$f(a,b)=\sqrt{ab}(\pi/2-\arctan\sqrt{a/b})=\sqrt{ab}\arctan\sqrt{b/a}$$
• @funtoast No.  – Did Sep 11 '18 at 10:19