I want to know if we can simplify the following double integral $$ I = \int_0^{\pi}dx\int_0^{2\pi-x}dy\sqrt{1+4\cos{x}\cos{\frac{y}{3}} + 4\cos^2{\frac{y}{3}}}. $$
I have tried to simplify it, but failed. I do not expect to find a closed-form of it (its numerical value is $\approx$23.3110) because it seems to be impossible for me (how about you?). However, it would nice if someone could point out whether $I$ is a transcendental number or not. I also wonder that if it is possible to transfer it to a single integral, like $I = \int_a^b dx\rho(x)$ where $a$ and $b$ are real numbers, and $\rho(x)$ is the kernel function.
Let $\mathcal{I}$ denote the following double integral:
$$\mathcal{I}:=\int_{0}^{\pi}\mathrm{d}x\int_{0}^{2\pi-x}\mathrm{d}y\,\sqrt{1+4\cos{\left(x\right)}\cos{\left(\frac{y}{3}\right)}+4\cos^{2}{\left(\frac{y}{3}\right)}}.$$
The trigonometric expression under the radical is easily seen to be nonnegative for all real arguments by rewriting it in the form
$$\begin{align} 1+4\cos{\left(x\right)}\cos{\left(\frac{y}{3}\right)}+4\cos^{2}{\left(\frac{y}{3}\right)} &=\sin^{2}{\left(x\right)}+\left[\cos{\left(x\right)}+2\cos{\left(\frac{y}{3}\right)}\right]^{2}\ge0,\\ \end{align}$$
and so the value of $\mathcal{I}$ is a well-defined positive real number, and according to WolframAlpha it has an approximate numerical value of
$$\mathcal{I}\approx23.311.$$
We can rewrite $\mathcal{I}$ as a double integral of an algebraic function over a rectangular region after the following manipulations:
$$\begin{align} \mathcal{I} &=\int_{0}^{\pi}\mathrm{d}x\int_{0}^{2\pi-x}\mathrm{d}y\,\sqrt{1+4\cos{\left(x\right)}\cos{\left(\frac{y}{3}\right)}+4\cos^{2}{\left(\frac{y}{3}\right)}}\\ &=9\int_{0}^{\frac{\pi}{3}}\mathrm{d}x\int_{0}^{\frac{2\pi}{3}-x}\mathrm{d}y\,\sqrt{1+4\cos{\left(3x\right)}\cos{\left(y\right)}+4\cos^{2}{\left(y\right)}};~~~\small{\left[\left(x,y\right)\mapsto\left(3x,3y\right)\right]}\\ &=9\int_{0}^{\frac{\pi}{3}}\mathrm{d}x\int_{x-\frac{\pi}{6}}^{\frac{\pi}{2}}\mathrm{d}\varphi\,\sqrt{1+4\cos{\left(3x\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}};~~~\small{\left[y=\frac{\pi}{2}-\varphi\right]}\\ &=9\int_{-\frac{\pi}{6}}^{\frac{\pi}{6}}\mathrm{d}\theta\int_{\theta}^{\frac{\pi}{2}}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(3\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}};~~~\small{\left[x=\theta+\frac{\pi}{6}\right]}\\ &=9\int_{-\frac{\pi}{6}}^{\frac{\pi}{6}}\mathrm{d}\theta\,\int_{-\frac{\pi}{2}}^{\theta}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(3\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}};~~~\small{\left[\left(\theta,\varphi\right)\mapsto\left(-\theta,-\varphi\right)\right]}\\ &=\frac{9}{2}\int_{-\frac{\pi}{6}}^{\frac{\pi}{6}}\mathrm{d}\theta\,\int_{\theta}^{\frac{\pi}{2}}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(3\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}}\\ &~~~~~+\frac{9}{2}\int_{-\frac{\pi}{6}}^{\frac{\pi}{6}}\mathrm{d}\theta\,\int_{-\frac{\pi}{2}}^{\theta}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(3\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}}\\ &=\frac{9}{2}\int_{-\frac{\pi}{6}}^{\frac{\pi}{6}}\mathrm{d}\theta\,\int_{-\frac{\pi}{2}}^{\frac{\pi}{2}}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(3\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}}\\ &=\frac32\int_{-\frac{\pi}{2}}^{\frac{\pi}{2}}\mathrm{d}\theta\int_{-\frac{\pi}{2}}^{\frac{\pi}{2}}\mathrm{d}\varphi\,\sqrt{1-4\sin{\left(\theta\right)}\sin{\left(\varphi\right)}+4\sin^{2}{\left(\varphi\right)}};~~~\small{\left[\theta\mapsto\frac13\theta\right]}\\ &=\frac32\int_{-1}^{1}\mathrm{d}t\int_{-1}^{1}\mathrm{d}u\,\frac{\sqrt{1-4tu+4u^{2}}}{\sqrt{1-t^{2}}\sqrt{1-u^{2}}};~~~\small{\left[\left(\theta,\varphi\right)=\left(\arcsin{\left(t\right)},\arcsin{\left(u\right)}\right)\right]}.\\ \end{align}$$
Finally, the double integral $\mathcal{I}$ can be reduced to the one single-variable integral,
$$\begin{align} \mathcal{I} &=\frac32\int_{-1}^{1}\mathrm{d}t\int_{-1}^{1}\mathrm{d}u\,\frac{\sqrt{1-4tu+4u^{2}}}{\sqrt{1-t^{2}}\sqrt{1-u^{2}}}\\ &=\frac32\int_{-1}^{1}\mathrm{d}u\int_{-1}^{1}\mathrm{d}t\,\frac{\sqrt{1-4tu+4u^{2}}}{\sqrt{1-t^{2}}\sqrt{1-u^{2}}}\\ &=3\int_{0}^{1}\mathrm{d}u\int_{-1}^{1}\mathrm{d}t\,\frac{\sqrt{1-4tu+4u^{2}}}{\sqrt{1-t^{2}}\sqrt{1-u^{2}}};~~~\small{\left[\left(t,u\right)\mapsto\left(-t,-u\right)\,\text{symmetry}\right]}\\ &=\int_{0}^{1}\mathrm{d}u\,\frac{3}{\sqrt{1-u^{2}}}\int_{-1}^{1}\mathrm{d}t\,\frac{\sqrt{4u^{2}+1-4ut}}{\sqrt{1-t^{2}}}\\ &=\int_{0}^{1}\mathrm{d}u\,\frac{3}{\sqrt{1-u^{2}}}\int_{0}^{1}\mathrm{d}w\,\frac{\sqrt{\left(2u+1\right)^{2}-8uw^{2}}}{\sqrt{1-w^{2}}};~~~\small{\left[t=2w^{2}-1\right]}\\ &=\int_{0}^{2}\mathrm{d}x\,\frac{6}{\sqrt{4-x^{2}}}\int_{0}^{1}\mathrm{d}w\,\frac{\sqrt{\left(x+1\right)^{2}-4xw^{2}}}{\sqrt{1-w^{2}}};~~~\small{\left[u=\frac12x\right]}\\ &=\int_{0}^{2}\mathrm{d}x\,\frac{6\left(1+x\right)}{\sqrt{4-x^{2}}}\int_{0}^{1}\mathrm{d}w\,\frac{\sqrt{1-\left[\frac{4x}{\left(1+x\right)^{2}}\right]w^{2}}}{\sqrt{1-w^{2}}}\\ &=\int_{0}^{2}\mathrm{d}x\,\frac{6\left(1+x\right)}{\sqrt{4-x^{2}}}\,E{\left(\frac{2\sqrt{x}}{1+x}\right)},\\ \end{align}$$
where the complete elliptic integral of the second kind $E$ is defined here by the integral representation,
$$E{\left(\kappa\right)}:=\int_{0}^{1}\mathrm{d}t\,\frac{\sqrt{1-\kappa^{2}t^{2}}}{\sqrt{1-t^{2}}};~~~\small{-1\le\kappa\le1}.$$
I suspect there's more progress to be made here if we keep at it, but that's all I got for now at least...
Cheers!