Compute $\int\limits_{x_1}^{x_2} \int\limits_{y_1}^{y_2} \int\limits_{z_1}^{z_2} {1 \over \sqrt{x^2+y^2+z^2}} dx\ dy\ dz$

Here is what I get

Used Rubi for the indefinite integration part, since it gives simpler anti-derivatives, (no complex numbers show up). Then took the limits (assumed proper integrals)

 solR = Int[Int[Int[1/Sqrt[x^2 + y^2 + z^2], x], y], z]

Gives

(-(1/2))*z^2*ArcTan[(x*y)/(z*Sqrt[x^2 + y^2 + z^2])] - 
   (1/2)*y^2*ArcTan[(x*z)/(y*Sqrt[x^2 + y^2 + z^2])] - 
  (1/2)*x^2*ArcTan[(y*z)/(x*Sqrt[x^2 + y^2 + z^2])] + 
  y*z*ArcTanh[x/Sqrt[x^2 + y^2 + z^2]] + 
  x*z*ArcTanh[y/Sqrt[x^2 + y^2 + z^2]] + 
  x*y*ArcTanh[z/Sqrt[x^2 + y^2 + z^2]]

Now just took the limits, one at a time

solR2 = Assuming[Element[{z,y},Reals]&&x2>x1,
           Simplify[Limit[solR,x->x2]-Limit[solR,x->x1]]]
solR3 = Assuming[Element[{z,x1,x2},Reals]&&y2>y1&&x2>x1,
           Simplify[Limit[solR2,y->y2]-Limit[solR2,y->y1]]]
solR4 = Assuming[Element[{x1,x2,y1,y2},Reals]&&y2>y1&&x2>x1&&z2>z1,
           Simplify[Limit[solR3,z->z2]-Limit[solR3,z->z1]]]

And the result is 

(1/2)*(z1^2*ArcTan[(x1*y1)/(z1*Sqrt[x1^2 + y1^2 + z1^2])] + y1^2*ArcTan[(x1*z1)/(y1*Sqrt[x1^2 + y1^2 + z1^2])] + x1^2*ArcTan[(y1*z1)/(x1*Sqrt[x1^2 + y1^2 + z1^2])] - 
   z1^2*ArcTan[(x2*y1)/(z1*Sqrt[x2^2 + y1^2 + z1^2])] - y1^2*ArcTan[(x2*z1)/(y1*Sqrt[x2^2 + y1^2 + z1^2])] - x2^2*ArcTan[(y1*z1)/(x2*Sqrt[x2^2 + y1^2 + z1^2])] - 
   z1^2*ArcTan[(x1*y2)/(z1*Sqrt[x1^2 + y2^2 + z1^2])] - y2^2*ArcTan[(x1*z1)/(y2*Sqrt[x1^2 + y2^2 + z1^2])] - x1^2*ArcTan[(y2*z1)/(x1*Sqrt[x1^2 + y2^2 + z1^2])] + 
   z1^2*ArcTan[(x2*y2)/(z1*Sqrt[x2^2 + y2^2 + z1^2])] + y2^2*ArcTan[(x2*z1)/(y2*Sqrt[x2^2 + y2^2 + z1^2])] + x2^2*ArcTan[(y2*z1)/(x2*Sqrt[x2^2 + y2^2 + z1^2])] - 
   z2^2*ArcTan[(x1*y1)/(z2*Sqrt[x1^2 + y1^2 + z2^2])] - y1^2*ArcTan[(x1*z2)/(y1*Sqrt[x1^2 + y1^2 + z2^2])] - x1^2*ArcTan[(y1*z2)/(x1*Sqrt[x1^2 + y1^2 + z2^2])] + 
   z2^2*ArcTan[(x2*y1)/(z2*Sqrt[x2^2 + y1^2 + z2^2])] + y1^2*ArcTan[(x2*z2)/(y1*Sqrt[x2^2 + y1^2 + z2^2])] + x2^2*ArcTan[(y1*z2)/(x2*Sqrt[x2^2 + y1^2 + z2^2])] + 
   z2^2*ArcTan[(x1*y2)/(z2*Sqrt[x1^2 + y2^2 + z2^2])] + y2^2*ArcTan[(x1*z2)/(y2*Sqrt[x1^2 + y2^2 + z2^2])] + x1^2*ArcTan[(y2*z2)/(x1*Sqrt[x1^2 + y2^2 + z2^2])] - 
   z2^2*ArcTan[(x2*y2)/(z2*Sqrt[x2^2 + y2^2 + z2^2])] - y2^2*ArcTan[(x2*z2)/(y2*Sqrt[x2^2 + y2^2 + z2^2])] - x2^2*ArcTan[(y2*z2)/(x2*Sqrt[x2^2 + y2^2 + z2^2])] - 
   2*y1*z1*ArcTanh[x1/Sqrt[x1^2 + y1^2 + z1^2]] - 2*x1*z1*ArcTanh[y1/Sqrt[x1^2 + y1^2 + z1^2]] - 2*x1*y1*ArcTanh[z1/Sqrt[x1^2 + y1^2 + z1^2]] + 2*y1*z1*ArcTanh[x2/Sqrt[x2^2 + y1^2 + z1^2]] + 
   2*x2*z1*ArcTanh[y1/Sqrt[x2^2 + y1^2 + z1^2]] + 2*x2*y1*ArcTanh[z1/Sqrt[x2^2 + y1^2 + z1^2]] + 2*y2*z1*ArcTanh[x1/Sqrt[x1^2 + y2^2 + z1^2]] + 2*x1*z1*ArcTanh[y2/Sqrt[x1^2 + y2^2 + z1^2]] + 
   2*x1*y2*ArcTanh[z1/Sqrt[x1^2 + y2^2 + z1^2]] - 2*y2*z1*ArcTanh[x2/Sqrt[x2^2 + y2^2 + z1^2]] - 2*x2*z1*ArcTanh[y2/Sqrt[x2^2 + y2^2 + z1^2]] - 2*x2*y2*ArcTanh[z1/Sqrt[x2^2 + y2^2 + z1^2]] + 
   2*y1*z2*ArcTanh[x1/Sqrt[x1^2 + y1^2 + z2^2]] + 2*x1*z2*ArcTanh[y1/Sqrt[x1^2 + y1^2 + z2^2]] + 2*x1*y1*ArcTanh[z2/Sqrt[x1^2 + y1^2 + z2^2]] - 2*y1*z2*ArcTanh[x2/Sqrt[x2^2 + y1^2 + z2^2]] - 
   2*x2*z2*ArcTanh[y1/Sqrt[x2^2 + y1^2 + z2^2]] - 2*x2*y1*ArcTanh[z2/Sqrt[x2^2 + y1^2 + z2^2]] - 2*y2*z2*ArcTanh[x1/Sqrt[x1^2 + y2^2 + z2^2]] - 2*x1*z2*ArcTanh[y2/Sqrt[x1^2 + y2^2 + z2^2]] - 
   2*x1*y2*ArcTanh[z2/Sqrt[x1^2 + y2^2 + z2^2]] + 2*y2*z2*ArcTanh[x2/Sqrt[x2^2 + y2^2 + z2^2]] + 2*x2*z2*ArcTanh[y2/Sqrt[x2^2 + y2^2 + z2^2]] + 2*x2*y2*ArcTanh[z2/Sqrt[x2^2 + y2^2 + z2^2]])

screen shot

Mathematica graphics

Quick verification

   solR4/.{x1->2,x2->3,y1->1,y2->5,z1->3,z2->9}//N

Mathematica graphics

    NIntegrate[1/Sqrt[x^2+y^2+z^2],{x,2,3},{y,1,5},{z,3,9}]

Mathematica graphics


There is no easy solution. See arXiv:chem-ph/9508002.

However, in the simplest case of a unit cube at the origin MA is able to compute the integral

 Integrate[1/Sqrt[x^2 + y^2 + z^2], {x, -1/2, 1/2}, {y, -1/2, 1/2}, {z, -1/2, 1/2}]
  (*-(\[Pi]/2) - 2 ArcSinh[1/Sqrt[2]] + Log[97 + 56 Sqrt[3]]*)
  FullSimplify[%]
  (*-(\[Pi]/2) + Log[26 + 15 Sqrt[3]]*)

Which is the same result as in the paper.

I would like to add that for a unit cube MA can produce results in many partial cases, for instance

 Integrate[1/Sqrt[x^2 + y^2 + z^2], {x, 0, 1}, {y, 0, 1}, {z, 0, 1}]
 (*-(\[Pi]/4) - ArcSinh[1/Sqrt[2]] + Log[7 + 4 Sqrt[3]]*)

 Integrate[1/Sqrt[x^2 + y^2 + z^2], {x, 1, 2}, {y, 1, 2}, {z, 1, 2}]
 (*(5 \[Pi])/4 - 8 ArcCot[2] + 8 ArcCot[3] - 8 ArcCot[2 Sqrt[6]] + 
  4 ArcSinh[1/(2 Sqrt[2])] - 3 ArcSinh[1/Sqrt[2]] + ArcSinh[Sqrt[2]] - 
  4 ArcTan[1/3] + 3/2 ArcTan[4/3] - 3 ArcTan[Sqrt[2/3]] + 
  2 ArcTan[1/(2 Sqrt[6])] + Log[2] - 2 Log[625] - Log[40960000] + 
  12 Log[1 + Sqrt[3]] + 12 Log[1 + Sqrt[6]] + 2 Log[2 + Sqrt[6]]*)

 Integrate[1/Sqrt[x^2 + y^2 + z^2], {x, 1, 2}, {y, -1/2, 1/2}, {z, -1/2, 1/2}]
  (* -8 ArcCot[12 Sqrt[2]] + 2 ArcCot[2 Sqrt[6]] - ArcSinh[Sqrt[2]] + 
 ArcSinh[2 Sqrt[2]] + ArcTan[Sqrt[2/3]] - ArcTan[(2 Sqrt[2])/3] + 
 8 Log[1 + 3 Sqrt[2]] - 2 Log[289/5 (7 + 2 Sqrt[6])]*)

Often, Mathematica is quicker with indefinite integrals than with definite ones. It's annoying but a result is also producible as follows

int1 = Integrate[1/Sqrt[x^2 + y^2 + z^2], x];
int2 = Integrate[(int1 /. x -> x2) - (int1 /. x -> xx1), y];
int3 = Integrate[(int2 /. y -> y2) - (int2 /. y -> y1), z];
(int3 /. z -> z2) - (int3 /. z -> z1)

However, I would not rely on this result without further testing. (We might jump from one branch of, e.g. the complex logarithm to another.)

Tags:

Calculus And Analysis

Recent Posts