How would I proceed in computing:
$$\int_D e^{-(x_1^2 + ... + x_n^2)/2} dx \qquad \qquad \text{ with } \quad D=\left\{x \in \mathbb R^n: \left|\sum_{i=1}^n x_i\right| \le \sqrt{n}\right\}.$$
I need to somehow get rid of the implicitly given domain. Perhaps, by using properties of the Gaussian?
This is best understood in probabilistic terms.
Since the pdf of the standard normal distribution is given by $$ f(x)=\frac{1}{\sqrt{2\pi}}e^{-x^2/2} \tag{1}$$ we have that $$ P = \frac{1}{(2\pi)^{n/2}}\int_{D}e^{-(x_1^2+\ldots+x_n^2)/2}\,d\mu \tag{2}$$ is the probability that the sum of $n$ independent standard normal variables lies in $[-\sqrt{n},\sqrt{n}]$.
This sum is just a normal variable with mean zero and $\sigma^2=n$, hence:
$$ P = \frac{1}{\sqrt{2n\pi}}\int_{-\sqrt{n}}^{\sqrt{n}}\exp\left(-\frac{x^2}{2n}\right)\,dx =\frac{1}{\sqrt{2\pi}}\int_{-1}^{1}e^{-y^2/2}\,dy=\operatorname{Erf}\left(\frac{1}{\sqrt{2}}\right)\tag{3}$$ and: