Aloha :)
Sorry für die späte Antwort, ich hatte gestern keine Zeit für Mathe ;)
Die erste Idee ist natürlich die Verwendung von Kugelkoordinaten. Damit würde der Integrand zu \(r^2e^{-r^2}\sin\vartheta\) und wir müssten das Integral von \(r^2e^{-r^2}\) über \(dr\) bestimmen. Das ist nicht einfach zu lösen, weil die Ableitung des Exponeten (\(-2r\)) nicht als Faktor im Integranden auftaucht. In 3 Dimensionen kriegen wir das also nicht gelöst.
Das Integral zerfällt jedoch in 3 unabhängig voneinander berechnbare Integrale:$$I=\int\limits_{\mathbb R^3}e^{-(x^2+y^2+z^2)}\,dx\,dy\,dz=\int\limits_{\mathbb R}e^{-x^2}dx\int\limits_{\mathbb R}e^{-y^2}dy\int\limits_{\mathbb R}e^{-z^2}dz=\left[\int\limits_{\mathbb R}e^{-x^2}dx\right]^3$$Bei dem verbliebenen 1-dimensionalen Integral fehlt uns wieder die Ableitung des Exponenten (\(-2x\)) als Faktor im Integranden, damit die Integration einfach wird.
Bleibt noch der Versuch mit 2 Dimensionen:$$I^{2/3}=\left[\int\limits_{\mathbb R}e^{-x^2}dx\right]\cdot\left[\int\limits_{\mathbb R}e^{-x^2}dx\right]=\left[\int\limits_{\mathbb R}e^{-x^2}dx\right]\cdot\left[\int\limits_{\mathbb R}e^{-y^2}dy\right]$$$$\phantom{I^{2/3}}=\int\limits_{\mathbb R}\int\limits_{\mathbb R}e^{-(x^2+y^2)}dx\,dy$$Gehen wir nun zu Polarkoordinaten über$$\binom{x}{y}=\binom{r\cos\varphi}{r\sin\varphi}\quad;\quad r=[0;\infty[\;\;;\;\;\phi\in[0;2\pi]\;\;;\;\;dx\,dy=r\,dr\,d\varphi$$erhalten wir die Ableitung des Exponenten als Faktor im Integranden und die Integration wird sehr einfach:$$I^{2/3}=\int\limits_0^\infty\int\limits_0^{2\pi}e^{-(r^2\cos^2+r^2\sin^2\varphi)}r\,dr\,d\varphi=\int\limits_0^\infty\int\limits_0^{2\pi}e^{-r^2}r\,dr\,d\varphi=\int\limits_0^\infty e^{-r^2}r\,dr\int\limits_0^{2\pi}d\varphi$$$$\phantom{I^{2/3}}=\left[-\frac{1}{2}e^{-r^2}\right]_0^\infty\,\left[\varphi\right]_0^{2\pi}=\left(0-\left(-\frac{1}{2}\right)\right)\cdot\left(2\pi-0\right)=\pi$$Damit ist also \(I^{1/3}=\sqrt{\pi}\) und das ursrünglich zu bestimmende Integral ist:$$\boxed{I=(\sqrt{\pi})^3}$$
