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.
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 :