Analytical evaluation of double integral I'm looking to evaluate the following integral:
$$\int^a_{-a}\int^a_{-a}\frac{\sqrt{x^2+y^2}}{4a^2} \,\mathrm{d}x \mathrm{d}y,$$
where $a> 0$.
Mathematica gives the following answer:
$$ \frac{a}{3} (\sqrt{2} + \sinh^{-1}(1)).$$
This goes well beyond my basic calculus 101 training.
Is anybody able to give me a step-by-step analytical solution?
Regards
 A: HINTS:
Exploiting the even symmetry and making use of 
$$\int \sqrt{x^2+y^2}\,dx=\frac x2\sqrt{x^2+y^2}+\frac{x^2}{2} \log(x+\sqrt{x^2+y^2})+C$$
evaluation of the inner integral yields
$$\begin{align}
\int_{-a}^a\int_{-a}^a \frac{\sqrt{x^2+y^2}}{4a^2}\,dx\,dy &= \frac1{a^2} \int_0^a \int_0^a \sqrt{x^2+y^2}\,dx\,dy \\\\
&=\frac1{a^2}\int_0^a\left(\frac a2\sqrt{a^2+y^2}+\frac{a^2}{2} \log(a+\sqrt{a^2+y^2})-\frac{a^2}{2}\log(y)\right) \,dy
\tag 1
\end{align}$$
Use integration by parts for the middle term in the integral on the right-hand side of $(1)$


EDIT:
The OP has requested to see a presentation on evaluation of $\displaystyle \int_0^a \log(a+\sqrt{y^2+a^2})\,dy$ by integration of parts. To that end, we proceed.

Let $u=\log(a+\sqrt{y^2+a^2})$ and $v=y$.  Then, we can write
$$\begin{align}
\int_0^a \log(a+\sqrt{y^2+a^2})\,dy&=a\log(a(1+\sqrt{2}))-\int_0^a \frac{y^2}{a+\sqrt{y^2+a^2}}\frac{1}{\sqrt{y^2+a^2}}\,dy\\\\
&=a\log(a(1+\sqrt{2}))-\int_0^a \frac{\sqrt{y^2+a^2}-a}{\sqrt{y^2+a^2}}\,dy\\\\
&=a\log(a(1+\sqrt{2}))-a+a\int_0^a \frac{1}{\sqrt{y^2+a^2}}\,dy\\\\
&=a\log(a(1+\sqrt{2}))-a+a\log(a(1+\sqrt{2}))-a\log(a)\\\\
&=a(\log(a)-1+2\log(1+\sqrt{2}))
\end{align}$$
A: Use polar coordinates (I obtain it by dividing each quarter of the square along the diagonal starting at the origin and noticing that the integral of the considered function is the same when calculated on each of these two domains) :
$$I=\int^a_0\int^a_0 \frac{\sqrt{x^2+y^2}}{4a^2} \,dx \, dy=\frac{2}{4a^2} \int^{\frac\pi4}_0 \int^{a\sec \theta}_0 r^2 \, dr \, d\theta$$
So by symmetry :
$$\int^a_{-a} \int^a_{-a}\frac{\sqrt{x^2+y^2}}{4a^2} \, dx \, dy = 4I = \frac2{a^2} \int^{\frac\pi4}_0\int^{a\sec \theta}_0 r^2 \, dr \, d\theta \\
$$$$= \frac2{a^2}\int^{\frac\pi4}_0 \frac{a^3}3\sec^3\theta \, d\theta = \frac {2a} 3 \int^{\frac\pi4}_0 \sec^3\theta \, d\theta\\
$$$$=\frac {2a}3\left[\frac{\sec(\theta)\tan(\theta)+\ln(\sec(\theta)+\tan(\theta)}{2}\right]_0^{\frac\pi4}\\=\frac{a}{3} (\sqrt{2} + \sinh^{-1}(1)).$$
A: $\newcommand{\bbx}[1]{\,\bbox[8px,border:1px groove navy]{\displaystyle{#1}}\,}
 \newcommand{\braces}[1]{\left\lbrace\,{#1}\,\right\rbrace}
 \newcommand{\bracks}[1]{\left\lbrack\,{#1}\,\right\rbrack}
 \newcommand{\dd}{\mathrm{d}}
 \newcommand{\ds}[1]{\displaystyle{#1}}
 \newcommand{\expo}[1]{\,\mathrm{e}^{#1}\,}
 \newcommand{\ic}{\mathrm{i}}
 \newcommand{\mc}[1]{\mathcal{#1}}
 \newcommand{\mrm}[1]{\mathrm{#1}}
 \newcommand{\pars}[1]{\left(\,{#1}\,\right)}
 \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}}
 \newcommand{\root}[2][]{\,\sqrt[#1]{\,{#2}\,}\,}
 \newcommand{\totald}[3][]{\frac{\mathrm{d}^{#1} #2}{\mathrm{d} #3^{#1}}}
 \newcommand{\verts}[1]{\left\vert\,{#1}\,\right\vert}$

Note that
  $\ds{\left.\int_{-a}^{a}\int_{-a}^{a}{\root{x^{2} + y^{2}} \over 4a^{2}}
\,\dd x\,\dd y\,\right\vert_{\ a\ >\ 0} =
a\int_{0}^{1}\int_{0}^{1}\root{x^{2} + y^{2}}\,\dd x\,\dd y}$

With $\ds{\vec{r} = x\,\hat{x} + y\,\hat{y}}$, note that in $\ds{2}$D:
\begin{align}
r & = {1 \over 3}\nabla\cdot\pars{r\vec{r}\,} =
{1 \over 3}
\bracks{\partiald{\pars{rx}}{x} - \partiald{\pars{-ry}}{y}} =
{1 \over 3}\,\bracks{\nabla \times \pars{-ry\,\hat{x} + rx\,\hat{y}\,}}_{\ z}
\end{align}

With the Stokes Theorem:

\begin{align}
&\int_{0}^{1}\int_{0}^{1}\root{x^{2} + y^{2}}\,\dd x\,\dd y =
{1 \over 3}\oint_{\mrm{square}}\pars{-ry\,\hat{x} + rx\,\hat{y}\,}\cdot\dd\vec{r}
\\[5mm] = &\
{1 \over 3}\bracks{
\int_{0}^{1}\root{1 + y^{2}} \times 1\,\dd y +
\int_{1}^{0}\pars{-\root{x^{2} + 1} \times 1}\,\dd x} =
{2 \over 3}\int_{0}^{1}\root{1 + x^{2}}\,\dd x
\\[5mm] = &\
\bbx{\ds{\root{2} + \,\mrm{arcsinh}\pars{1} \over 3}}
\end{align}

The integral is evaluated with the change of variables
  $\ds{x = \sinh\pars{\theta}}$:

\begin{align}
\int_{0}^{1}\root{1 + x^{2}}\,\dd x & =
\int_{0}^{\mrm{arcsinh\pars{1}}}\cosh^{2}\pars{\theta}\,\dd\theta =
\int_{0}^{\mrm{arcsinh\pars{1}}}{1 + \cosh\pars{2\theta} \over 2}\,\dd\theta
\\[5mm] & =
{1 \over 2}\,\mrm{arcsinh}\pars{1} +
\bracks{{1 \over 2}\,\sinh\pars{\theta}
\root{1 + \sinh^{2}\pars{\theta}}}_{\ 0}^{\ \,\mrm{arcsin}\pars{1}}
\\[5mm] & =
{\,\mrm{arcsinh}\pars{1} + \root{2} \over 2}
\end{align}
A: $$
\frac 1 {4a^2} \int^a_{-a}\int^a_{-a} \sqrt{x^2+y^2} \, dx \, dy = \text{what?}
$$
$$
\int_{-a}^a (\text{an even function}) = 2\times\int_0^a (\text{the same function}), \text{ so}
$$
\begin{align}
& \int_{-a}^a \sqrt{x^2+y^2} \  dx = 2\int_0^a \sqrt{x^2+y^2}\  dx = 2\int_0^{\arctan a} \sqrt{y^2\tan^2\theta + y^2} (\sec^2\theta\,d\theta) \\[10pt]
= {} & 2|y| \int_0^{\arctan a} \Big(\sec\theta\Big) (\sec^2\theta\,d\theta) = 2|y| \int_0^{\arctan a} \sec^3\theta\,d\theta.
\end{align}
This is one of the more challenging among freshman-level integrals. This article tells us how to find that
$$
\int \sec^3\theta\,d\theta = \frac 1 2 \sec\theta\tan\theta + \frac 1 2 \log|\sec\theta+\tan\theta| + C.
$$
So we need
$$
\left[ \frac 1 2 \sec\theta\tan\theta + \frac 1 2 \log|\sec\theta+\tan\theta| \right]_0^{\arctan\theta}.
$$
From trigonometry we have
\begin{align}
\tan\arctan a & = a, \\[10pt]
\sec\arctan a & = \sqrt{1+a^2}.
\end{align}
Maybe you can do the rest. (?)
