Average distance of the mean of $n$ random complex numbers in a unit disc Let $z_1,z_2,\dots,z_n$ be $n$ complex numbers distributed uniformly and randomly over the unit disc $x^2+y^2 \leq 1$. Let $z$ be the complex number  defined by the mean of the of these numbers,that is,
$$ z=\frac1n \sum_{i=1}^{n} z_i.$$
Numerical simulations reveal that $E(|z|)$ exists. I want to know if there is some analytical way of finding a closed form expression for the expected radial distance of $z$ in terms of $n$? Any hints/responses will be greatly appreciated . The graph shows the variation of $E(|z|)$ with $n$ 
 A: Here is another method.
Since the uniform distribution on the unit disk is rotationally symmetric (invariant under orthogonal transformations), the problem can be reduced to a random walk problem in $\mathbb{R}$.
Let $Z_1,Z_2,\ldots$ be independent, and identically uniformly distributed on the unit disk, with $Z_i=(X_i,Y_i)$.
Clearly $X_1,X_2\ldots$ are i.i..d with density $p(x)=\frac{2}{\pi}\sqrt{1-x^2}1_{[-1,1]}(x)$.
For any rotationally symmetric
random vector $V=(V_1,V_2)$ with finite expectation of $|V|:=\sqrt{V_1^2+V_2^2}$  it holds that
$$\mathbb{E}(|V|)=\frac{\pi}{2}\mathbb{E}(|V_1|)$$
Applying this to $Z^{(n)}:=                                 |Z_1+\ldots+Z_n|$ we get that
\begin{align*}
\mathbb{E}\big(\frac{Z^{(n)}}{n}\big)&=\frac{1}{n}\frac{\pi}{2}\mathbb{E}|X_1+\ldots + X_n|\\
                                 &=\frac{1}{n}\big(\frac{2}{\pi}\big)^{n-1}\,\int_{-1}^1\ldots\int_{-1}^1|x_1+\ldots+x_n|\prod_{i=1}^n\sqrt{1-x_i^2}\,dx_1\,\ldots dx_n\\
                                 \end{align*}
The integrals can be solved explicitly for $n=1$ and $n=2$, but that seems to be as far as it goes. But (using the central limit theorem and uniform integrability) it is easy to see that $$\mathbb{E}\big(\frac{Z^{(n)}}{\sqrt{n}}\big)\longrightarrow \sqrt{\frac{\pi}{8}}$$
Remarks:
(1) a related problem appeared here An interesting triple integral
(2) for similar considerations see Feller II (1971), p. 30 ff
A: This does not answer the question, but reformulates it as an integration problem and addresses a question asked in the comments.
Conditionally on $|z_1| = a_1$, ..., $|z_n| = a_n$, the probability that $x<c$ was computed in Kluyver, "A local probability problem" (1906) to be
$$
\mathbb{P}\Big( \big|\sum_{i=1}^n z_i\big| < c \,\Big|\, |z_i|=a_i \Big) = c \int_0^\infty J_1(uc)\prod_{i=1}^n J_0(u a_i)\,\mathrm{d}u.
$$
Using that $\int_0^1 2a\,J_0(u a) \mathrm{d}a = 2 J_1(u) /u$ and $\frac{\partial}{\partial c}( c J_1(u c)) = uc J_0(u c)$, this gives the density as
$$
\rho(x) = \frac{\partial}{\partial x}\mathbb{P}\Big( \big|\sum_{i=1}^n z_i\big| < x\Big) = \int_0^\infty xu J_0(ux) \,\left(\tfrac{2}{u}J_1(u)\right)^n\mathrm{d}u.
$$
(Note: there was a small mistake in the formula given in the comments.)
By construction $\int_0^n\rho(x)\,\mathrm{d}x = 1$ as $\rho(x)$ vanishes for $x > n$.
Hence
\begin{align}
\mathbb{E}\Big( \frac{1}{n}\big|\sum_{i=1}^n z_i\big|\Big) &= \frac{1}{n} \int_0^n \mathrm{d}x \int_0^\infty \mathrm{d}u\, x^2u\, J_0(ux) \left(\tfrac{2}{u}J_1(u)\right)^n \\
&= \frac{n^2}{3} \int_0^\infty u\,{_1F_2}(\tfrac{3}{2};1,\tfrac{5}{2}; -\tfrac14 n^2u^2) \left(\tfrac{2}{u}J_1(u)\right)^n,
\end{align}
where the second integral was kindly provided by mathematica.
One could hope to extract some exact values for $n\geq 3$ from Borwein, Jonathan M. "A short walk can be beautiful." Journal of Humanistic Mathematics 6, no. 1 (2016): 86-109.
