What is the probability distribution for the squared distance between random points in an $n$-dimensional hypercube? I choose random points $X,\,Y$ in $[0,\,1]^n$ (I assume all $2n$ Cartesian coordinates are $U(0,\,1)$ iids). What is the probability distribution of $\left\Vert X-Y\right\Vert _{2}^{2}$? Even the $n=1$ case requires some care, since if I first fix $X$ then $Y$ has a uniform distribution with extrema of opposite sign.
 A: For an approximate cdf, one possibility is the saddlepoint approximation. For that the mgf (moment generating function), or its logarithm, the cgf (cumulant generating function) is needed. For background see How does saddlepoint approximation work?.   So let $X_1, \dotsc, X_n, Y_1, \dotsc, Y_n$ be iid $\mathcal{U}(0,1)$. Then $D^2 = \sum_i (X_i-Y_i)^2$ have mgf
$$ \DeclareMathOperator{\E}{\mathbb{E}} \DeclareMathOperator{\P}{\mathbb{P}}  
M_n(t)=M_{\sum_i (X_i-Y_i)^2}(t)=M_{(X_1-Y_1)^2}(t)^n=\left(\E e^{(X_1-Y_1)^2 t}\right)^n
$$ so the cgf is
$$
   K_n(t)=n\cdot \log\left\{ \E e^{(X_1-Y_1)^2 t} \right\}
$$ which in this case is defined for all real $t$.  
Without going into details, and with some help from maple, we find that
$$
 \E e^{(X_1-Y_1)^2 t} = {\frac {\sqrt {\pi}{\rm erf} \left(\sqrt {-t}\right){t}^{2}+{{\rm e}^{
t}} \left( -t \right) ^{3/2}- \left( -t \right) ^{3/2}}{ \left( -t
 \right) ^{5/2}}}
$$  or alternatively
$$
 \E e^{(X_1-Y_1)^2 t} = {\frac {-2\,i{{\rm e}^{t}}{\it dawson} \left( i\sqrt {-t} \right) {t}^
{2}+{{\rm e}^{t}} \left( -t \right) ^{3/2}- \left( -t \right) ^{3/2}}{
 \left( -t \right) ^{5/2}}}
$$ but before using this numerically, more work will be needed to find a useful form. 
I will come back later to add the resulting approximation. 
