Calculating in closed form $\int _0^1\int _0^1\frac{1}{1+x y (x+y)} \ dx \ dy$

Solution 1:

By symmetry (i.e. by exploiting the fact that our integral is twice the integral over the sub-region $0\leq y\leq x\leq 1$) we just have to compute:

$$ I=2 \int_{0}^{1}\int_{0}^{1}\frac{x}{1+x^3 y(1+y)}\,dy\,dx = 2 \int_{0}^{1}\int_{0}^{2}\frac{x}{(1+x^3y)\sqrt{1+4y}}\,dy\,dx$$ Integrating with respect to $x$ first,

$$\begin{eqnarray*} I &=& 2\int_{0}^{2}\frac{2 \sqrt{3}\arctan\left(\frac{\sqrt{3}}{1-2y^{1/3}}\right)-2\log\left(1+y^{1/3}\right)+\log\left(1-y^{1/3}+y^{2/3}\right)}{6 y^{2/3}\sqrt{1+4y}}\,dy\\&=&\int_{0}^{2^{1/3}}\frac{2\sqrt{3}\arctan\left(\frac{\sqrt{3}}{1-2y}\right)-3\log(1+y)+\log(1+y^3)}{\sqrt{1+4y^3}}\,dy\end{eqnarray*}$$ but the resulting integrals in just one variable do not look so appealing.

Am I missing some crucial simplification that follows from replacing $y$ with a Jacobi elliptic function (maybe $\text{dn}$) or with the Weierstrass elliptic function $\wp(z)$ (corresponding to $g_2=0,g_3=-1$) then exploiting some weird/mystical product formulas?

Solution 2:

Another single integral form, $$ 2 \int_0^1 \frac{\arctan{\big(\sqrt{y(4-y^3)}/(y^2+2)\big)}}{\sqrt{y(4-y^3)}}\,dy. $$ It doesn't look simpler than the first 1D integral given, however, so I'll skip the proof.

Solution 3:

Starting in the polar coordinates, \begin{align} &J=\int\limits_0^1\int\limits_0^1 \dfrac{dx\,dy}{1+xy(x+y)} = 2\int\limits_0^{\dfrac\pi4}\int\limits_0^\dfrac1{\cos\varphi}\dfrac{\rho\,d\rho\,d\varphi}{1+\rho^3\cos\varphi\sin\varphi(\cos\varphi+\sin\varphi)} = 2\int\limits_0^{\pi/4}F(\rho,\varphi)\,d\varphi,\\ \end{align} and substitution \begin{align} &\rho=\dfrac1{v \cos\varphi}\tag1 \end{align} gives \begin{align} &F(\rho,\varphi) = \int\limits_1^\infty\dfrac1 {1+\dfrac1{v^3\cos^3\varphi}\cos\varphi\sin\varphi(\cos\varphi+\sin\varphi)}\,\dfrac{dv}{v^3\cos^2\varphi} =\\ &\dfrac1{\cos^2\varphi}\int\limits_1^\infty\dfrac{dv}{\tan\varphi(\tan\varphi+1)+ v^3},\quad\text{or}\\ &J=2\int\limits_1^\infty G(v,\varphi)\, dv,\quad\text{where}\\ &G(v,\varphi) = \int\limits_0^{\pi/4}\dfrac{dv}{\tan\varphi(\tan\varphi+1)+v^3}\cdot\dfrac{d\varphi}{\cos^2\varphi}=\int\limits_0^1\dfrac{dt}{t^2+t+v^3} = \\ &\dfrac1{\sqrt{v^3-\dfrac14}}\left.\arctan\dfrac{t+\dfrac12}{\sqrt{v^3-\dfrac14}}\right|_{t=0}^1 = \dfrac1{\sqrt{v^3-\dfrac14}}\left(\arctan\dfrac{\dfrac32}{\sqrt{v^3-\dfrac14}}-\arctan\dfrac{\dfrac12}{\sqrt{v^3-\dfrac14}}\right),\tag2\\[4pt] &G(v,\varphi) = \dfrac1{\sqrt{v^3-\dfrac14}}\arctan\dfrac{\dfrac1{\sqrt{v^3-\dfrac14}}}{1+\dfrac3{4\left(v^3-\dfrac14\right)}} = \dfrac2{\sqrt{4v^3-1}}\arctan\dfrac{2\sqrt{4v^3-1}}{4v^3+2}.\tag3 \end{align} Therefore, \begin{align} &J=\int\limits_1^\infty \dfrac4{\sqrt{4v^3-1}}\arctan\dfrac{2\sqrt{4v^3-1}}{4v^3+2}\, dv\approx 0.798965\tag4 \end{align} (see also Wolfram Alpha).

Substitution $$u=\sqrt{4v^3-1}$$ leads to \begin{align} &v=\dfrac{\sqrt[3]2}2\sqrt[3]{u^2+1},\quad dv=\dfrac{2\sqrt[3]2}3\dfrac1{\sqrt[3]{(u^2+1)^2}}u\,du,\\ &J=\dfrac{4\sqrt[3]2}3\int\limits_\sqrt3^\infty\arctan\dfrac {2u}{3+u^2}\dfrac{du}{\sqrt[3]{(u^2+1)^2}}\tag5\\ &J=\dfrac{4\sqrt[3]2}3\int\limits_\sqrt3^\infty\arctan\dfrac {u-\dfrac13u}{1+u\cdot\dfrac13u}\dfrac{du}{\sqrt[3]{(u^2+1)^2}},\\ &J=\dfrac{4\sqrt[3]2}3\int\limits_\sqrt3^\infty\left(\arctan u - \arctan\dfrac u3\right)\dfrac{du}{\sqrt[3]{(u^2+1)^2}} = J_1-J_2, \tag6\\ \end{align} (see also Wolfram Alpha), wherein $$J_1=\dfrac{4\sqrt[3]2}3\int\limits_{\pi/3}^{\pi/2}\dfrac{z\,dz}{\sqrt[3]{\cos^2z}}\\$$ has the closed form.

But the further attempts to get a closed form weren't successful.

Solution 4:

$\def\l{\ell}$I have found (many) integral forms no nicer than those in the other answers. Here I give a series expansion.

We have \begin{align*} \frac{1}{1+xy(x+y)} &= \frac{1}{2\left(1+\frac{xy(x+y)-1}{2}\right)} \\ &= \sum_{k=0}^\infty \frac{(-1)^k}{2^{k+1}} [xy(x+y)-1)]^k \\ &= \sum_{k=0}^\infty \frac{(-1)^k}{2^{k+1}} \sum_{\l=0}^k {k\choose \l} [xy(x+y)]^\l (-1)^{k-\l} \\ &= \sum_{k=0}^\infty \frac{(-1)^k}{2^{k+1}} \sum_{\l=0}^k (-1)^{k-\l}{k\choose \l} x^\l y^\l \sum_{m=0}^\l {\l\choose m} x^{\l-m}y^m \\ &= \sum_{k=0}^\infty \sum_{\l=0}^k \sum_{m=0}^\l \frac{(-1)^\l}{2^{k+1}}{k\choose \l}{\l\choose m} x^{2\l-m}y^{\l+m}. \end{align*} Thus, \begin{align}\tag{1} \int_0^1\int_0^1 \frac{1}{1+xy(x+y)} dx\, dy &= \sum_{k=0}^\infty \sum_{\l=0}^k \sum_{m=0}^\l \frac{(-1)^\l}{2^{k+1}}{k\choose \l}{\l\choose m} \frac{1}{(2\l-m+1)(\l+m+1)}. \end{align} Similarly, \begin{align*} \int_0^1 \cdots \int_0^1 & \frac{1}{1+x_1\cdots x_n(x_1+\ldots +x_n)} dx_1\cdots dx_n \\ &= \sum_{k=0}^\infty \sum_{\l=0}^k \sum_{{m_1,\ldots,m_n}\atop{m_1+\ldots +m_n=\l}} (-1)^\l\frac{(n-1)^{k-\l}}{n^{k+1}} {k\choose \l} \frac{\l!}{m_1!\cdots m_n!} \frac{1}{(\l+m_1+1)\cdots (\l+m_n+1)}. \end{align*} It is possible to write the sum in equation (1) in terms of a single (infinite) sum over special functions, but it is not particularly illuminating. I am curious to see if there is some nice closed form at all!

Solution 5:

Let $a$ be real and $|a|\leq1$. Then we set $$ I_n(a):=\int^{a}_{0}\int^{a}_{0}x^{n}y^n(x-y)^ndxdy\textrm{, }n=0,1,2,\ldots. $$ and $$ f(a,t)=\int^{a}_{0}\int^{a}_{0}\frac{1}{1-txy(x-y)}dxdy $$ Hence $$ f(a,t)=\sum^{\infty}_{l=0}t^l\int^{a}_{0}\int^{a}_{0}x^{n}y^n(x-y)^ndxdy=\sum^{\infty}_{n=0}t^nI_n(a). $$ Under the change of coordinates $x\rightarrow r\sin(\theta)$, $y\rightarrow r\cos(\theta)$, we get $$ I_n(a)=2\int^{\pi/4}_{0}\int^{a/\cos(\theta)}_{0}(r\sin(\theta))^{n}(r\cos(\theta))^{n}r^n(\sin(\theta)-\cos(\theta))^nrdrd\theta= $$ $$ 2\int^{\pi/4}_{0}\sin^n(\theta)\cos^n(\theta)(\sin(\theta)-\cos(\theta))^n\left(\int^{a/\cos(\theta)}_{0}r^{3n+1}dr\right)d\theta= $$ $$ 2\frac{a^{3n+2}}{3n+2}\int^{\pi/4}_{0}\sin^n(\theta)\cos^{-2n-2}(\theta)(\sin(\theta)-\cos(\theta))^nd\theta $$ But (see Wolfram alpha) it holds $$ \int^{\pi/4}_{0}\sin^n(\theta)\cos^{-2n-2}(\theta)(\sin(\theta)-\cos(\theta))^nd\theta=\frac{2(-1)^n\Gamma(n+1)\Gamma(n+2)}{\Gamma(2n+3)}. $$ Hence when $n=0,2,4,\ldots$ $$ I_n(a)=\frac{4(-1)^na^{3n+2}}{3n+2}\frac{\Gamma(n+1)\Gamma(n+2)}{\Gamma(2n+3)} $$ and $$ I_n(a)=0 $$ when $n=1,3,5,\ldots$.

Consequently $$ f(a,t)=2\sum^{\infty}_{n=0}\frac{a^{6n+2}}{6n+2}B(2n+1,2n+1)t^{2n}=a^2\cdot { }_4F_{3}\left(\frac{1}{2},\frac{1}{3},1,1;\frac{3}{4},\frac{5}{4},\frac{4}{3};\frac{a^6t^2}{16}\right) $$ and clearly if $|a|\leq1$, $|t|\leq 1$, then $$ \int^{a}_{0}\int^{a}_{0}\frac{1}{1-txy(x-y)}dxdy=a^2\cdot { }_4F_{3}\left(\frac{1}{2},\frac{1}{3},1,1;\frac{3}{4},\frac{5}{4},\frac{4}{3};\frac{a^6t^2}{16}\right)\tag 1 $$

Note that, an interesting formula is $$ \int^{\arctan(a)}_{0}\int^{b/\cos(\theta)}_{0}(r\sin(\theta))^n(r\cos(\theta))^n r^n (\sin(\theta)-\cos(\theta))^n rdrd\theta= $$ $$ =\frac{(-1)^n b^{3n+2}}{3n+2}B(a,n+1,n+1)\textrm{, }Re(n)>-1 $$ Also for $|a|\leq 1$ and $n=0,1,2,\ldots$ it holds (see Wolfram alpha): $$ \int^{a}_{0}\int^{a}_{0}(xy(x+y))^ndxdy= $$ $$ =\int^{\pi/4}_{0}\int^{a/\cos(\theta)}_{0}(r\sin(\theta))^n(r\cos(\theta))^n r^n (\sin(\theta)+\cos(\theta))^n rdrd\theta= $$ $$ =2\frac{a^{3n+2}}{3n+2}\cdot { }_2F_1\left(-n,n+1;n+2;-1\right) $$ Hence $$ J=\int^{a}_{0}\int^{a}_{0}\frac{1}{1-txy(x+y)}dxdy=2\sum^{\infty}_{n=0}\frac{a^{3n+2}}{3n+2}\cdot \frac{{ }_2F_{1}\left(-n,n+1;n+2;-1\right)}{n+1}t^n= $$ $$ =2a^2\sum^{\infty}_{n=0}\frac{(-1)^{n+1}}{3n+2}\cdot B\left(-1,n+1,n+1\right)(a^3t)^n\tag 2 $$ where $|a|<1$ and $|t|\leq1$ or if $|a|\leq 1$ and $|t|<1$.

Also is $$ B(-1,n+1,n+1)=\int^{-1}_{0}t^{n}(1-t)^ndt=(-1)^{n+1}\int^{1}_{0}t^n(1+t)^ndt $$ Having in mind the above integral we can write (easy) $$ \int^{a}_{0}\int^{a}_{0}\frac{1}{1-txy(x+y)}dxdy=a^2\int^{1}_{0}{ }_2F_1\left(1,\frac{2}{3};\frac{5}{3};a^3tx(1+x)\right)dx. $$ Setting $x(1+x)\rightarrow y$ and then $y\rightarrow 2w$, we arrive to $$ J=2a^2\int^{1}_{0}\frac{{ }_2F_1\left(1,\frac{2}{3};\frac{5}{3};2a^3tw\right)}{\sqrt{1+8w}}dw. $$ But $$ { }_2F_1\left(1,\frac{2}{3};\frac{5}{3};x\right)=-\frac{2\log(1-x^{1/3})}{3x^{2/3}}+\frac{2(-1)^{1/3}\log\left(1-e^{-2\pi i/3}x^{1/3}\right)}{3x^{2/3}}-\frac{2(-1)^{2/3}\log\left(1-e^{2\pi i/3}x^{1/3}\right)}{3x^{2/3}} $$ Hence if we set $$ F(z)=\int^{1}_{0}\frac{\log(1-zx)}{\sqrt{1+8x^3}}dx,\tag 3 $$ we get $$ J=-2\cdot 2^{1/3}t^{-2/3}F\left(2^{1/3}at^{1/3}\right)+2\cdot (-2)^{1/3}t^{-2/3}F\left(2^{1/3}e^{-2\pi i/3}at^{1/3}\right)-2\cdot (-1)^{2/3}2^{1/3}t^{-2/3}F\left(2^{1/3}e^{2\pi i/3}at^{1/3}\right)\tag 4 $$