Compute integral $\int_0^1\int_0^1\int_0^1 \sqrt{x^2+y^2+z^2} \,\mathrm{d}x\mathrm{d}y\mathrm{d}z$

Consider the following figure:

Due to symmetry one sixth of the rays emanating from $(0,0,0)$ leave the cube through the red triangle $T$. We therefore use the parametrization $$\psi:\quad \left.\eqalign{x&=t\> u\cr y&=t\>v \cr z&=t\cr}\right\}\qquad(0\leq v\leq u, \ 0\leq u\leq 1, \ 0\leq t\leq 1)$$ with Jacobian $J_\psi=t^2$. It follows that the integral $(=:Q)$ in question is given by $$Q=6\int_0^1\int_0^1\int_0^u t\sqrt{1+u^2+v^2}\>t^2 dv\ du\ dt={3\over2}\int_0^1\int_0^u \sqrt{1+u^2+v^2}\ dv\ du\ .$$ Herewith the geometrical difficulties are out of the way and it remains to compute the nested integral using standard calculus techniques. A first step could be introducing polar coordinates $$u=r\cos\phi, \quad v=r\sin\phi\qquad\left(0\leq r\leq{1\over\cos\phi},\quad 0\leq\phi\leq{\pi\over4}\right)\ ,$$ which leads to $$Q={3\over2}\int_0^{\pi/4}\int_0^{1/\cos\phi}r\>\sqrt{1+r^2}\ dr\ d\phi\ .$$ The inner integral is easy, and we are left with $$Q={1\over2}\int_0^{\pi/4}\left((1+\cos^{-2}\phi)^{3/2}-1\right)\ d\phi\ ,$$ which has the quoted value $\ldots$

You are computing the average distance from a vertex in a unit hypercube. Assuming that $X,Y,Z$ are three independent random variables, uniformly distributed over $(0,1)$, and computing the probability density function of $W=X^2+Y^2+Z^2$ (that is a piecewise-simple function), it follows that

$$\mathbb{E}[\sqrt{W}] = \color{red}{\frac{\sqrt{3}}{4}+\frac{\log(2+\sqrt{3})}{2}-\frac{\pi}{24}}=0.960591956455\ldots $$ As an alternative, you may just integrate $\sqrt{x^2+y^2+z^2}$ with respect to $x$ first, getting
$\frac{1}{2} x \sqrt{x^2+y^2+z^2}+\frac{1}{2} \left(y^2+z^2\right) \log\left(x+\sqrt{x^2+y^2+z^2}\right)$, then integrate the resulting expression with respect to $y$ and $z$. Lengthy but perfectly doable.