Xjyuk

I am seeking a closed-form solution for this double sum:

begineqnarray*
sum_n=1^infty sum_m=1^infty frac1nm(colorblue3n+m)= ?.
endeqnarray*

I will turn it into $3$ tough integrals in a moment. But first I will state some similar results:

begineqnarray*
sum_n=1^infty sum_m=1^infty frac1nm(n+m) &=& 2 zeta(3) \
sum_n=1^infty sum_m=1^infty frac1nm(colorblue2n+m) &=& frac118 zeta(3) \
sum_n=1^infty sum_m=1^infty frac1nm(colorblue4n+m) &=& frac6732 zeta(3) -fracG pi2. \
endeqnarray*

where $G$ is the Catalan constant. The last result took some effort ...

Now I know most of you folks prefer integrals to sums, so lets turn this into an integral. Using

begineqnarray*
frac1n &=& int_0^1 x^n-1 dx\
frac1m &=& int_0^1 y^m-1 dy\
frac13n+m &=& int_0^1 z^3n+m-1 dz \
endeqnarray*
and summing the geometric series, we have the following triple integral
begineqnarray*
int_0^1 int_0^1 int_0^1 fracz^3 dx dy dz(1-xz^3)(1-yz).
endeqnarray*

Now doing the $x$ and $y$ integrations we have
begineqnarray*
I=int_0^1 fracln(1-z) ln(1-z^3)z dz.
endeqnarray*

Factorize the argument of the second logarithm ...

begineqnarray*
I= underbraceint_0^1 fracln(1-z) ln(1-z)z dz_=2zeta(2) + int_0^1 fracln(1-z) ln(1+z+z^2)z dz.
endeqnarray*

So if you prefer my question is ... find a closed form for:

begineqnarray*
I_1 = - int_0^1 fracln(1-z) ln(1+z+z^2)z dz.
endeqnarray*

Integrating by parts gives:

begineqnarray*
I_1 = - int_0^1 fracln(z) ln(1+z+z^2)1-z dz + int_0^1 frac(1+2z)ln(z) ln(1-z)1+z+z^2 dz.
endeqnarray*

and let us call these integrals $I_2$ and $I_3$ respectively.

All $3$ of these integrals are not easy for me to evaluate and any help with their resolution will be gratefully received.

$$boxedI=int_0^1 fracln(1-x) ln(1-x^3)xdx=frac53zeta(3) +frac2pi^327sqrt 3 -fracpi9sqrt 3psi_1left(frac13right)$$
Where the value of the integral can be extracted from here, since using $2ab=a^2+b^2-(a-b)^2$ we have:
$$2I=underbraceint_0^1 fracln^2(1-x)+ln^2(1-x^3)xdx_=frac83zeta(3)-int_0^1 fracln^2(1+x+x^2)xdx$$

Alternatively we can make use of the following series:
$$ln(1+x+x^2)=-2sum_n=1^infty fraccosleft(frac2n pi3right)n x^n $$
$$mathcal J =int_0^1 fracln(1-x)ln(1+x+x^2)xdx=-2sum_n=1^infty fraccosleft(frac2n pi3right)nint_0^1 ln(1-x) x^n-1dx $$
$$=-2sum_n=1^infty fraccosleft(frac2n pi3right)n^2H_n=-2Re left(sum_n=1^infty fracz^nn^2H_nright),quad z=e^frac2pi i3$$
Now using the following generating function:
$$sum_n=1^infty fracx^nn^2H_n=operatornameLi_3(x)-operatornameLi_3(1-x)+operatornameLi_2(1-x)ln(1-x)+frac12ln x ln^2(1-x)+zeta(3)$$
And simplifying should give the value for the integral, but I don't really know how to obtain the values for the real part of the polylogarithms.
Some can be found here.

Yet another idea: $$mathcal J=int_0^1 fracln (1-x)ln(1+x+x^2)xdxoversetxto 1-x=int_0^1 fracln xln(3-3x+x^2)1-xdx$$
Now consider the following integral (quite a weird one) and derivate under the integral sign:
$$I(t)=int_0^1 fracln xln[(1-x)colorblue(2(1+cos t))+x^2]1-xdxRightarrow I'(t)=-2sin tint_0^1 fracln x2(1-x)(1+cos t)+x^2dx$$
$$overset1-xto x=-2int_0^1 fracsin tln(1-x)x^2+(2cos t )x+1dx=2sum_n=1^infty (-1)^n sin (nt)int_0^1 x^n-1ln(1-x)dx$$
$$=2sum_n=1^infty (-1)^n sin (nt)fracH_nn=2Imsum_n=1^infty frac(-1)^n e^in tnH_n=2Im sum_n=1^infty fracz^nnH_n,quad z=-e^it $$
Also, using the generating function:
$$sum_n =1^infty fracH_nnz^n=operatornameLi_2(z)+frac12ln^2(1-z)Rightarrow I'(t)=Imleft(2 operatornameLi_2(z)+ln^2(1-z)right)$$
We are looking to find $I(0)$ since $cos 0=1$, but we have:
$$Ileft(piright)=int_0^1 fracln xln(x^2)1-xdx=2 sum_n=0^infty int_0^1 x^n ln ^2 xdx=8sum_n=1^infty frac1n^3=8zeta(3)$$
$$mathcal J=-(I(0)-I(pi)+8zeta(3))=-int_0^pi Im(2operatornameLi_2(z)+ln^2(1-z))dt+8zeta(3)$$
Maybe noticing that: $Im operatornameLi(e^ix)=operatornameCl_2(x)$ can help here, although I have a strong feeling it's going to the approach from above.

Not a full answer, but another interesting expression for the series.

Let's introduce a function:

$$S(x,y)=sum_n=1^infty sum_m=1^infty fracx^n y^m n m (3n +m)$$

Assuming $|x|<1$ and $|y|<1$ we avoid any convergence problems, and can use partial fractions:

$$S(x,y)=sum_n=1^infty sum_m=1^infty fracx^n y^m n m^2-sum_n=1^infty sum_m=1^infty fracx^n y^m m^2(n+ frac13 m)$$

$$S(x,y)=-log(1-x) textLi_2(y) -x sum_m=1^infty fracy^m m^2 Phi left(x,1,frac13 m+1 right)$$

Let's use the integral representation of the Lerch transcendent:

$$Phi left(x,1,frac13 m+1 right)= int_0^infty frace^-(1+frac13 m)t ~dt1-x e^-t$$

Summation under the integral gives us:

$$S(x,y)=-log(1-x) textLi_2(y) -x int_0^infty textLi_2 left(y e^-t/3 right) fracdte^t-x$$

So we can assume:

$$S=lim_x to 1 left[-log(1-x) textLi_2(x) -x int_0^infty textLi_2 left(x e^-t/3 right) fracdte^t-x right]$$

Which does seem to work numerically, though of course is quite hard to evaluate symbolically.

Numerical check:

In[22]:= x=9999999/10000000;
y=9999999/10000000;
N[-Log[1-x]PolyLog[2,y],10]-x NIntegrate[PolyLog[2,y Exp[-t/3]]/(Exp[t]-x),t,0,Infinity,WorkingPrecision->10]
Out[24]= 1.29484017

Compare to exact expression:

In[25]:= N[(1/216)*(-15*Pi^2*Log[3] + 9*Log[3]^3 + 4*Sqrt[3]*Pi*(-PolyGamma[1, 1/3] + 
 PolyGamma[1, 2/3]) - 216*(PolyLog[3, (-1)^(1/6)/Sqrt[3]] + 
 PolyLog[3, -((-1)^(5/6)/Sqrt[3])]) + 672*Zeta[3]), 10]
Out[25]= 1.2948652620+0.*10^-11 I

It may be that $x=y$ is not the best choice for the limit. For example, we can assume $x=y^a$ where $a$ is some real number. A good choice may lead to a better numerical convergence or even a closed form.

Using the dilogarithm properties, we have:

$$F(x,y)=-x int_0^infty textLi_2 left(y e^-t/3 right) fracdte^t-x=x int_0^infty int_0^1 fraclog(1-e^-t/3 y u) du dtu (e^t-x)$$

Let's change the variable:

$$e^-t=v \ t=- log v$$

$$F(x,y)=x int_0^1 int_0^1 fraclog(1- y u v^1/3) du dvu (1-x v)$$

Let's take:

$$y=x^1/3$$

We have:

$$F(x,y)=x int_0^1 int_0^1 fraclog(1- u (xv)^1/3) du dvu (1-x v)$$

$$v=w/x$$

$$F(x,y)=int_0^x int_0^1 fraclog(1- u w^1/3) du dwu (1-w)$$

$$F(x,y)=-int_0^x fractextLi_2 (w^1/3) dw1-w$$

So, there's a neater expression for the limit:

$$ colorblueS=lim_x to 1 left[-log(1-x) textLi_2(x^1/3) -int_0^x fractextLi_2 (w^1/3) dw1-w right]$$

This allows a simple generalization:

$$sum_n=1^infty sum_m=1^infty fracx^n y^m n m (an +m)=lim_x to 1 left[-log(1-x) textLi_2(x^1/a) -int_0^x fractextLi_2 (w^1/a) dw1-w right]$$

Which checks out numerically with the examples from the OP.

I wonder if we can somehow use L'Hospital here to deal with the integral and get a closed form for the limit.

Integration by parts could also work.