Is there a closed form for the sum of the reciprocal cubes of all coordinate distances? I have a problem where I need to sum reciprocal cubed distances between the origin, $(0, 0)$, and all points on a coordinate grid. For a distance $d$, by reciprocal cube I mean $\frac{1}{d^3}$. All points are specified by $(m, n)$ for $m,n \in \mathbb{Z}$ (excluding $(0, 0)$), thus all the distances are $\sqrt{m^2 + n^2}$ (i.e. all hypotenuses of integer-legged right triangles).
I desire a closed form of the following:
$\sum_{m,n \in \mathbb{Z}}^{'} \frac{1}{\sqrt{m^2 + n^2}^3}$
This should be equal to $4\cdot\sum_{m,n \in \mathbb{N}}^{'} \frac{1}{\sqrt{m^2 + n^2}^3}$, as the distances in all quadrants are considered, where the sum over $\mathbb{N}$ ($0\in\mathbb{N}$ here) gives just the first quadrant. I welcome a closed form of either sum.
Supplemental:
If no closed form is possible, I would also be interested in other sums that would converge faster and would be easier to compute, in order to estimate a decimal value for the sum.
A less useful, but equivalent sum may be $\zeta(3)\cdot\sum_{coprime(m,n)}^{'} \frac{1}{\sqrt{m^2 + n^2}^3}$, which I have partially computed to be $\approx 1.6$ thus far, but which may change as the sum increases and converges. In english, this is the sum of all coprime distances (a line from the $(0, 0)$ to $(m, n)$ intersecting no other points), times Apérty's constant, which extends each coprime distance line to all other lines (e.g. $(3, 4) \rightarrow (3, 4), (6, 8), (9, 12), \dots$).
Disclaimer: I have little experience with sums and the $\zeta$ function. If I am blatantly wrong in an assumption (e.g. the series converges), please illustrate where.
 A: Such lattice sums were considered from an analytic point of view in this answer.
Other (possibly simpler) methods were used here and here.
The more general case $\,\Re(s)>1$ was considered :
$$\tag{1}S(s):=\sum_{m,n \in \mathbb{Z}}^{'}\frac{1}{\left ( m^2+n^2 \right )^s}=4\,\sum_{m=0}^{\infty}\;\sum_{n=1}^{\infty}\frac{1}{\left ( m^2+n^2 \right )^s}$$
(notice the initial values of $\,m,n$)
The different derivations proposed will all give you the wished :
$$\tag{2}\sum_{m,n \in \mathbb{Z}}^{'}\frac{1}{\left ( m^2+n^2 \right )^s}=4\,\zeta(s)\,\beta(s),\quad\Re(s)>1$$
with $\beta$ the Dirichlet beta function : $\;\displaystyle\beta(s) = \sum_{n=0}^\infty \frac{(-1)^n} {(2n+1)^s}$.
That is for $s=\dfrac 32$ :
$$\tag{3}S\left(\frac 32\right)=4\;\zeta\left(\frac 32\right)\beta\left(\frac 32\right)\approx 9.0336216831009503$$

Let's consider the modified sum 
$$\tag{4}S_C(s):=\sum_{m,n \in \mathbb{Z}}^{'}\frac{1}{\left ( m^2+n^2+C \right )^s}=4\,\sum_{m=0}^{\infty}\;\sum_{n=1}^{\infty}\frac{1}{\left ( m^2+n^2+C \right )^s}$$
then $(4)$ in the analytical method from the first link should become 
$$\tag{5}\Gamma(s)\,S_C(s)=\int_0^\infty t^{s-1}e^{-Ct}\left(\theta_3(0,e^{-t})^2-1\right)\;dt$$
while the inverse transformation between $(6)$ and $(7)$ should be :
$$\tag{6}S_C(s)=4\,\sum_{n=1}^\infty \sum_{m=0}^\infty (-1)^m \frac{1}{(C+n(2m+1))^s}$$
The convergence of this sum is better than using $(4)$ (as you may experiment $(*)$) but I don't see how to obtain a nicer factorization...
To conclude in a more positive way let's use $(6)$ to obtain at least an expansion in the case $|C|<1$ :
\begin{align}
S_C(s)&=4\,\sum_{k=0}^\infty\frac{(-C)^k}{k!}\sum_{n=1}^\infty \sum_{m=0}^\infty (-1)^m \frac{s(s+1)\cdots(s+k-1)}{(n(2m+1))^{s+k}}\\
&=4\,\sum_{k=0}^\infty\frac{(-C)^k}{k!}s(s+1)\cdots(s+k-1)\,\zeta(s+k)\,\beta(s+k)\\
\tag{7}S_C(s)&=4\,\sum_{k=0}^\infty(-C)^k\binom{s+k-1}{k}\,\zeta(s+k)\,\beta(s+k),\quad |C|<1\\
\end{align}
where I expressed the rising factorial $s^{(k)}=s(s+1)\cdots (s+k-1)\,$ divided by $\,k!\,$ as a binomial polynomial.
$(*)$ pari/gp scripts used :
S(C,s)=4*sumpos(k=0,sumpos(j=1,1/(k^2+j^2+C)^s))
L(C,s)=4*sumpos(n=1,sumalt(m=0,(-1)^m/(C+n*(2*m+1))^s))
beta(s)=sumalt(n=0,(-1)^n/(2*n+1)^s)
A(C,s,m)=4*sum(k=0,m,(-C)^k*binomial(s+k-1,k)*zeta(s+k)*beta(s+k))

