Answer
$$\iint \limits_{R} \frac{x y}{\sqrt{x^{2}+y^{2}+1}} =\sqrt{3}- \frac{4\sqrt{2}}{3}+\frac{1}{3}$$
Work Step by Step
Given $$\iint \limits_{R} \frac{x y}{\sqrt{x^{2}+y^{2}+1}} d A, \quad R=[0,1] \times[0,1]$$
So, we have
\begin{aligned}
I&=\iint \limits_{R} \frac{x y}{\sqrt{x^{2}+y^{2}+1}} d A\\
&=\int_{0}^{1} \int_{0}^{1} \frac{x y}{\sqrt{x^{2}+y^{2}+1}} d y d x\\
&=\frac {1}{2}\int_{0}^{1} \int_{0}^{1} 2x y (x^{2}+y^{2}+1)^{-\frac {1}{2}} d y d x\\
&=\int_{0}^{1} \left( x (x^{2}+y^{2}+1)^{\frac {1}{2}} \right) _{0}^{1} d x\\
&= \int_{0}^{1} \left( x (x^{2}+1+1)^{\frac {1}{2}} \right) d x-\ \int_{0}^{1} \left( x (x^{2}+0+1)^{\frac {1}{2}} \right) d x\\
&=\frac {1}{2}\int_{0}^{1} \left(2 x (x^{2}+2)^{\frac {1}{2}} \right) d x-\frac {1}{2}\int_{0}^{1} \left( 2 x (x^{2}+ 1)^{\frac {1}{2}} \right) d x\\
&=\frac {1}{3} \left( (x^{2}+2)^{\frac {3}{2}} \right)_{0}^{1}-\frac {1}{3} \left( (x^{2}+ 1)^{\frac {3}{2}} \right)_{0}^{1}\\
&=\frac {1}{3} \left( (1^{2}+2)^{\frac {3}{2}}- (0+2)^{\frac {3}{2}} \right) -\frac {1}{3} \left( (1+ 1)^{\frac {3}{2}}-(0+ 1)^{\frac {3}{2}} \right) \\
&=\frac {1}{3} \left( (3^{\frac {3}{2}}- 2^{\frac {3}{2}} \right) -\frac {1}{3} \left( 2^{\frac {3}{2}}-1^{\frac {3}{2}} \right)\\
&=\frac {1}{3} \left( 3^{\frac {3}{2}}- 2^{\frac {3}{2}}-2^{\frac {3}{2}}+1^{\frac {3}{2}} \right) \\
&=\sqrt{3}- \frac{4\sqrt{2}}{3}+\frac{1}{3}
\end{aligned}