Another approach. Since:

\[ \int_{0}^{+\infty}\sin(ax)e^{-bx}\,dx = \frac{a}{a^2+b^2} \tag{1}\] we have:

\[ \sum_{j=1}^{n^2}\frac{n}{n^2+j^2} = \int_{0}^{+\infty}\frac{1-e^{-n^2 x}}{e^x-1}\sin(nx)\,dx \tag{2}\] where: \[ \int_{0}^{+\infty}\frac{\sin(nx)}{e^x-1}\,dx = \text{Im}\int_{0}^{+\infty}\frac{e^{inx}}{e^x-1}\,dx=\sum_{k\geq 1}\frac{n}{n^2+k^2}=\frac{-1+\pi n \coth(\pi n)}{2n}\tag{3}\] by Frullani's theorem and/or the logarithmic derivative of the Weierstrass product for the \(\sinh\) function. It is quite trivial that the limit of the RHS of \((3)\) as \(n\to +\infty\) is \(\frac{\pi}{2}\), hence it is enough to prove that the contribute given by \[ \int_{0}^{+\infty}\frac{\sin(nx) e^{-n^2 x}}{e^x-1}\,dx\quad\text{or}\quad\sum_{k>n^2}\frac{n}{n^2+k^2}\tag{4}\] is negligible. But that is quite easy.

Here is a method in the same spirit as what you have done. The main difference is that we explicitly estimate the difference between your sum and \(\int_0^n\frac{{\rm d}x}{1+x^2}\).

Since \(f\) is monotonely decreasing we have 

\[\left|\sum_{j=1}^n \frac{n}{n^2+(j+nm)^2} - \int_{m}^{m+1}\frac{{\rm d}x}{1+x^2}\right| < \frac{1}{n(1 + (m+1)^2)}\]

Using this we can estimate the difference between your sum and the integral of \(\frac{1}{1+x^2}\) over \([0,n]\), which can be written \(\sum_{m=0}^{n-1}\int_0^1\frac{{\rm d}x}{1+(x+m)^2}\), as

\[\left|\sum_{j=1}^{n^2}\frac{n}{n^2 + j^2} - \int_0^n\frac{{\rm d}x}{1+x^2}\right| \leq \sum_{m=0}^{n-1} \left|\frac{1}{n}\sum_{j=1}^n \frac{1}{1+\left(\frac{j}{n}+m\right)^2} - \int_{0}^{1}\frac{{\rm d}x}{1+(x+m)^2}\right| \\\leq \frac{1}{n}\sum_{m=1}^{n} \frac{1}{1+m^2} \leq \frac{1}{n}\sum_{m=1}^\infty \frac{1}{1+m^2}\]

This last sum can be evaluated as \(\frac{\pi\coth(\pi)-1}{2}\), however all we need for this argument is that it's finite. This can be proved for example by comparison to \(\sum_{m=1}^\infty \frac{1}{m^2} = \frac{\pi^2}{6}\). Taking \(n\to\infty\) the result follows

\[\lim_{n\to\infty}\sum_{j=1}^{n^2}\frac{n}{n^2 + j^2} = \lim_{n\to\infty}\int_0^n\frac{{\rm d}x}{1+x^2} = \frac{\pi}{2}\]