Calculating squared reciprocals of arithmetic series

What follows is not a closed form but a different representation. The series part of the sum is $$S(x) = \sum_{k\ge 1} \frac{1}{(1+kx)^2}$$ evaluated at $x=n.$

The sum term is harmonic and may be evaluated by inverting its Mellin transform.

Recall the harmonic sum identity $$\mathfrak{M}\left(\sum_{k\ge 1} \lambda_k g(\mu_k x);s\right) = \left(\sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} \right) g^*(s)$$ where $g^*(s)$ is the Mellin transform of $g(x).$

In the present case we have $$\lambda_k = 1, \quad \mu_k = k \quad \text{and} \quad g(x) = \frac{1}{(1+x)^2}.$$ We need the Mellin transform $g^*(s)$ of $g(x)$ which is $$\int_0^\infty \frac{1}{(1+x)^2} x^{s-1} dx = \left[- \frac{1}{1+x} x^{s-1} \right]_0^\infty + (s-1) \int_0^\infty \frac{1}{1+x} x^{s-2} dx.$$ The bracketed term vanishes for $1<\Re(s)<2.$ Now the Mellin transform $h^*(s)$ of $h(x) = \frac{1}{1+x}$ is easily seen to obey (contour integration with a keyhole contour, slot on the positive real axis) $$h^*(s) (1-e^{2\pi i (s-1)}) = 2\pi i\times \mathrm{Res}\left(\frac{1}{1+x} x^{s-1} ; x=-1 \right)$$ which gives $$h^*(s) = \frac{2\pi i}{1-e^{2\pi i s}} e^{i\pi (s-1)} = - \frac{2\pi i}{1-e^{2\pi i s}} e^{i\pi s} = - \pi \frac{2i}{e^{-i\pi s}-e^{i\pi s}} = \frac{\pi}{\sin(\pi s)}.$$ Hence the Mellin transform $g^*(s)$ of $g(x)$ has the form $$g^*(s) = (s-1) \frac{\pi}{\sin(\pi (s-1))} = (1-s) \frac{\pi}{\sin(\pi s)}.$$ The harmonic sum identity then implies that the Mellin transform $Q(s)$ of $S(x)$ is given by $$(1-s) \frac{\pi}{\sin(\pi s)} \zeta(s) \quad\text{because}\quad \sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} = \zeta(s).$$ The Mellin inversion integral for this case is $$\frac{1}{2\pi i}\int_{3/2-i\infty}^{3/2+i\infty} Q(s)/x^s ds.$$ Now inverting the transform in the right half plane to obtain an expansion about infinity, we have that from the poles at $q\ge 2$ the sum of their residues is $$S(x) = -\sum_{q\ge 2} (1-q) (-1)^q \zeta(q) / x^q.$$ This produces for the original sum the representation $$1+ \sum_{q\ge 2} (q-1) (-1)^q\zeta(q) \frac{1}{n^q},$$ convergent for $n>1.$

For $n > 0$ (whether integer or not is irrelevant) $$\sum_{j=0}^\infty \dfrac{1}{(1+jn)^2} = \dfrac{1}{n^2} \Psi'(1/n)$$ where $\Psi$ is the Digamma function. Known values include $$ \eqalign{n = 1, & sum = \dfrac{\pi^2}{6}\cr n = 2, & sum = \dfrac{\pi^2}{8}\cr n = 4, & sum = \dfrac{\pi^2}{16} + \dfrac{{\rm Catalan}}{2}\cr}$$