Define $I(s)$ by
$$ I(s) = \int_{0}^{\infty} \frac{e^{-su^2}}{1+u^2} \, du. $$
Then $I(s)$ solves the following equation:
$$ I(s) - I'(s)
= \int_{0}^{\infty} e^{-su^2} \, du
= \frac{1}{2}\sqrt{\frac{\pi}{s}}. $$
This is a 1st-order differential equation, which can be solved systematically by means of integration factor. The result is that
$$ I(s) = e^s \left( \mathsf{C}-\int \frac{1}{2}\sqrt{\frac{\pi}{s}}e^{-s} \, ds \right) $$
for some appropriate choice of constant $\mathsf{C}$. Together with the boundary condition $I(\infty) = 0$, it turns out that
$$ I(s)
= e^s \int_{s}^{\infty} \frac{1}{2}\sqrt{\frac{\pi}{s'}}e^{-s'} \, ds'
= \frac{\pi}{2}e^s \operatorname{erfc}(\sqrt{s}). $$
Plugging $s = 1$ gives the value
$$ I(1) = \frac{e\pi}{2}\operatorname{erfc}(1) \approx 0.67164671082336758522\cdots. $$