Evaluation of a triple sum does not finish in reasonable time

Question

I'm trying to compute the following triple sum, but no result is produced within a reasonable amount of time. What to do?

Sum[1/( i j k (i + j + k + 1)), {i, 1, Infinity}, {j, 1, Infinity}, {k, 1, Infinity}]

I use Mathematica 8.0.

Answer 1

@whuber gave me another idea : first change the summand to an exponential and then do the sums using the symmetry.

Integrate[Exp[-t i j k (i + j + k + 1)], {t, 0, Infinity}, Assumptions -> {i > 0, j > 0, k > 0}]
(* 1/(i j k (1 + i + j + k)) *)

so we need to do the sums of the exponential and then integrate at the end. Now we can make a change of variables t i j k = y which will simplify the exponent and bring a factor 1/(i j k) from the Jacobian. This makes it easier to see that

$$\sum_{i,j,k=1}^{\infty} \frac{1}{i j k} \exp(-y\ (i+j+k+1)) = \exp(-y) \left( \sum_{i=1}^{\infty} \frac{1}{i}\exp(-y\ i)\right)^3$$

Simplify[Sum[Exp[-y i]/i, {i, 1, Infinity}, Assumptions -> {y > 0}], 
  Assumptions -> {y > 0}]
(* -Log[1 - E^-y] *)

and we are left with the integration :

Integrate[Exp[-y] (-Log[1 - E^-y])^3, {y, 0, Infinity}]
(* 6 *)

It's also straightforward to generalize to the sum over n integers and the result is :

Integrate[Exp[-y] (-Log[1 - E^-y])^n, {y, 0, Infinity}, 
  Assumptions -> {n \[Element] Integers, n >= 1}] 
(* Gamma[1 + n] *)

Answer 2

Here is another variant along the same lines as @whuber proposed. Consider function $$ f(x,y,z,t)=\sum _{i=1}^{\infty } \sum _{j=1}^{\infty } \sum _{k=1}^{\infty } x^{i-1} y^{j-1} z^{k-1} t^{i+j+k}. $$ Mma gives

f[x_, y_, z_, t_] = 
 Sum[x^(i - 1) y^(j - 1) z^(k - 1) t^(i + j + k), {i, 
   1, \[Infinity]}, {j, 1, \[Infinity]}, {k, 1, \[Infinity]}]

$$ -\frac{t^3}{(t x-1) (t y-1) (t z-1)}. $$ Then the required sum is equal to $$ \int _0^1\int _0^1\int _0^1\int _0^1f(x,y,z,t)\,dzdydxdt, $$

Integrate[f[x, y, z, t], {t, 0, 1}, {x, 0, 1}, {y, 0, 1}, {z, 0, 1}]

which gives 6 as in the answer of @b.gatessucks.

Answer 3

If I take the challenge to be to get Mathematica to do the sum in the form given (as distinct from applying some mathematics insightfully), then here is a way to get the answer. The only trick I use is to break the triple sum into a double and single sum.

s2 = Sum[1/(i j (i + j + k + 1)), {i, 1, Infinity}, {j, 1, Infinity}]; // Timing
(s2 = FullSimplify[s2]) // Timing
Sum[s2/k, {k, 1, Infinity}] // Timing

{93.696694, Null}
{0.234598, (\[Pi]^2 + 6 HarmonicNumber[1 + k]^2 - 6 PolyGamma[1, 2 + k])/(6 + 6 k)}
{12.027743, 6}

The intermediate step of simplifying s2 seems to be necessary.

Interestingly Mathematica caches auxiliary expressions in computing the sum. The second time around is faster. (Almost all the time is spent in FullSimplify)

(s2 = Sum[1/(i j (i + j + k + 1)), {i, 1, Infinity}, {j, 1, Infinity}];
  s2 = FullSimplify[s2];
  Sum[s2/k, {k, 1, Infinity}]) // Timing

{0.228539, 6}

Answer 4

Not an answer, but it might help. You can get one summation in closed form :

sum1 = Sum[1/(i j k (i + j + k + 1)), {i, 1, Infinity}]
(* (EulerGamma + PolyGamma[0, 2 + j + k])/(j k (1 + j + k)) *)

You can also do the easy bit :

sum2 = Sum[Sum[EulerGamma/(j k (1 + j + k)), {j, 1, Infinity}], {k, 1, Infinity}]
(* 2 EulerGamma *)

so you're left with a double sum

sum = 2 EulerGamma + 
      Sum[PolyGamma[0, 2 + j + k]/(j k (1 + j + k)), {j, 1, Infinity}, {k, 1, Infinity}]

You can make some further simplification by using the symmetry of the summand or proceed numerically.