Is the $d$-dimensional harmonic series is proportional to the surface-area of a $d$-sphere?

Question

TLDR:

How to Prove:

$$ \sum_{\sqrt{a_0^2 + ...+ a_k^2} \le n, (a_0 ,...,a_k) \ne (0)^{k+1} } \frac{1}{(a_0^2 + ...+ a_k^2)^{\frac{k}{2}}} = \frac{k\pi^{\frac{k}{2}}}{\Gamma(\frac{k}{2} + 1)} \ln(n) + O(1) $$

Context:

This originates as an offshoot of this question.

I have been able to experimentally verify (up to $n=30$) that the following seem true:

$$ \begin{matrix} \sum_{\sqrt{a_0^2} \le n, a_0 \ne 0} \frac{1}{(a_0^2)^{\frac{1}{2}}} = 2\ln(n)+O(1)\\ \sum_{\sqrt{a_0^2+a_1^2} \le n, (a_0,a_1) \ne (0,0)} \frac{1}{(a_0^2 + a_1^2)^{\frac{2}{2}}} = 2\pi \ln(n) +O(1) \\ \sum_{\sqrt{a_0^2+a_1^2+a_2^2} \le n, (a_0,a_1,a_2) \ne (0,0,0)} \frac{1}{(a_0^2 + a_1^2+a_2^2)^{\frac{3}{2}}} = 4\pi \ln(n) + O(1) \\ \sum_{\sqrt{a_0^2+a_1^2+a_2^2+a_3^2} \le n, (a_0,a_1,a_2,a_3) \ne (0,0,0,0)} \frac{1}{(a_0^2 + a_1^2+a_2^2+a_3^2)^{\frac{4}{2}}} = 2\pi^2 \ln(n) +O(1) \\ \end{matrix} $$

The first sum is just a verbose way of stating the harmonic series identity and the later sums are the $d$-dimensional generalizations of this.

In all the cases we see some mysterious coefficients in front of the $\ln(n)$ and those coefficients are precisely the $(n-1)$-content of an $n$-sphere.

In particular there are

1. $2$ endpoints of a unit line segment.
2. The perimeters of a unit circle is $2\pi$
3. The area of a unit spheres is $4\pi$
4. The boundary volume of a $3$-sphere is $2\pi^2$

So this leads to a natural conjecture which doesn't feel too surprising.

$$ \sum_{\sqrt{a_0^2 + ...+ a_k^2} \le n, (a_0 ,...,a_k) \ne (0)^{k+1} } \frac{1}{(a_0^2 + ...+ a_k^2)^{\frac{k}{2}}} = \frac{k\pi^{\frac{k}{2}}}{\Gamma(\frac{k}{2} + 1)} \ln(n) + O(1) $$

From using the closed form formula of an $n$-ball here and differentiating.

What would it take to prove this statement? Is there ideally a nice way to "see" why this ought to be true.

Some Observations:

The case of $n=1$ is easy to see in any Calc 1 textbook. And @Gerry Myerson found some very helpful links that all discuss $n=2$ (in the linked question)

1. Computing $\sum{\frac{1}{m^2+n^2}}$
2. Convergence of multiple series $\sum\frac{1}{n^2+m^2}$
3. How to prove $\sum_{m=1}^\infty \sum_{n=1}^\infty \frac{1}{m^2+n^2}=+\infty$

But the general case seems like it might be a bit harder to prove. It's clear these sums take the form of evaluating a harmonic like thing over the integer points of the $d$-dimensional sphere, so there should be some weird way to look at that sum as an integral over the sphere and then try some stokes-theorem like idea of transforming that integral over the interior to something involving just the boundary. Invoking techniques similar to proving the general formula for the volume of an $n$-ball might be relevant here.