Expected distance between two points inside unit circle 
Question: What is the expected distance between two points inside a unit circle? 

I know the expected distance on the unit circle is $\dfrac 4\pi$.
The expected distance $E$ must equal the integral
$$E=\frac1{\pi^2} I=\frac1{\pi^2}\iint_{x^2+y^2\le 1}\iint_{z^2+w^2\le 1} \sqrt{(x-w)^2+(y-z)^2}\,\mathrm d x\,\mathrm d y\,\mathrm d z\,\mathrm d w.$$
Someone gives the beautiful answer $E=\dfrac1{\pi^2}I=\dfrac{128}{45\pi}$.
Polar coordinates seems complicated here. Any ideas to simplify the integral?
 A: Polar coordinates are the best way because of the circular geometry, you just have to survive the integration. Start with the integral over position of both points (as suggested in the comments by @Zang), and you can already integrate over one of the angles due to symmetry:
$$2\pi \int_0^1\int_0^1\int_0^{2\pi}\sqrt{R^2+r^2-2rR\cos\theta} r \,\mathrm d r\, R \mathrm d R\, \mathrm d \theta$$
This is tricky, and can be integrated in many ways, but I prefer this one:
Introduce new variables
$$u=rR, \quad v=\frac{r}{R}$$
The Jacobian of the transformation is $|J|=\frac{1}{2v}$. We have yet to determine the limits:
$$2\pi \iiint_{?}\sqrt{uv+u/v+2u \cos\theta}\frac{u}{2v} \, {\rm d}u\,{\rm d}v\, {\rm d}\theta$$
The integrand simplifies a lot for integration over $u$, so we want to do that first (innermost integral should be over $u$). The previous domain in $R,r$ coordinates was a unit square. Due to symmetry, we can fold it in half, and integrate over the lower triangle. $v=r/R$ will run from $0$ to $1$, and $u$ will run from $0$ to the maximum value of $rR$, which happens when $R=1$ and $r=r/1=r/R=v$. We have:
$$4\pi \int_0^{2\pi}\int_0^1\int_0^v\sqrt{v+1/v+2 \cos\theta}\frac{1}{2v} u^{3/2}\,{\rm d}u\,{\rm d}v\, {\rm d}\theta=$$
$$=2\pi \int_0^{2\pi}\int_0^1\sqrt{v+1/v+2 \cos\theta}\frac{1}{v}\frac25 v^{5/2} \, {\rm d}v\, {\rm d}\theta=$$
$$=\frac{4}{5}\pi \int_0^{2\pi}\int_0^1\sqrt{v^2-2 v\cos\theta+1}v{\rm d}v\, {\rm d}\theta=\frac{4}{5}\pi \int_0^{2\pi}\int_0^1\underbrace{\sqrt{(v-\cos\theta)^2+\sin^2\theta}}_{\text{distance}}\underbrace{v \, {\rm d}v\, {\rm d} \theta}_{dS}$$
The remaining integral is the expected value of distance from an interior point to the point on the boundary (my transformation is not the only one that does this, you could simply take out the $R$ and integrate over $r$, or something like that - I just wanted a symmetric substitution).
You can solve this one by moving the center to the single point at the edge (e.g. (1,0)) so the distance to the "current" point becomes the new radius. This can be done "graphically", so instead of doing a complicated change of variables by hand, move the origin, which transfers all the complicated geometry into the limits of the integral and you get:
$$=\frac{4}{5}\pi \int_{-\pi/2}^{\pi/2}\int_0^{\sqrt{2(1+\cos2\theta')}} v'^2\,{\rm d}v'\,{\rm d}\theta'$$
The limits are determined just by the length of the chord at angle $\theta'$ with respect to the diameter.
$$=\frac{4}{5}\pi \int_{-\pi/2}^{\pi/2}\frac13(2(1+\cos2\theta'))^{3/2} \, {\rm d} \theta'=\frac{4}{5}\pi \int_{-\pi/2}^{\pi/2}\frac13(4\cos^2\theta')^{3/2} \, {\rm d}\theta'=$$
$$=\frac{4}{15}\pi \int_{-\pi/2}^{\pi/2}8|\cos\theta'|^3 \, {\rm d}\theta' = \frac{64}{15} \pi\int_0^{\pi/2}\cos\theta'^3 \, {\rm d}\theta' = \frac{64\pi}{15}\frac{2}{3}=\frac{128}{45}\pi$$

TLDR: Calculation can be done in two steps:

*

*Solve interior-point-to-edge-point problem by centering the polar coordinates at the edge.

*Average this solution over possible radii.

I did it in reverse and moved the coordinate origin by hand mid-calculation, but now you know the trick, you can start with that.
