I have started with the following relation:
$ 2\sum_{j=1}^{N}\frac{1}{\tan\left(\frac{\pi j}{2N+1}\right)}=\frac{4N+2}{\pi}H_N+2\sum_{j=1}^{N}\left[\frac{1}{\tan\left(\frac{\pi j}{2N+1}\right)}-\frac{1}{(\frac{\pi j}{2N+1})}\right]\tag1$
According to WolframAlpha, $\frac{1}{\tan(x)}-\frac{1}{x}$ can be integrated over interval $\left(0,\frac{\pi}{2}\right)$ to equal $\log \left(\frac{2}{\pi }\right)$.
So using $\log \left(\frac{2}{\pi }\right)$ in the right hand side of (1) gives:
$(\frac{4 \text{N}+2}{\pi}) (H_{\text{N}}+\log (\frac{2}{\pi }))\tag2$ Empirically (1) seems to be related to (2).
My questions are
How do I show (1) is related to (2)?
How can I prove this?
How do I find bounds for any error? or does the harmonic number and nearer the $\tan$ part approximates simply reduce this error the larger N gets? If so then how do I prove this?
The question is equivalent to an asymptotic expansion for the following, as $n \to \infty,$ $$ S_n:=\sum_{k=1}^n \cot{(\pi\ k/(2n+1))} $$ The proof will include the terms proposed, plus an extra one. The key formula to use can be derived from Gradshteyn 4.131.4 $$ \frac{1}{2a} - \frac{\pi}{2b} \cot{(\pi \ a/b)} = \int_0^\infty \frac{e^{-bx}}{\sinh{bx}} \ \sinh{(2ax)} \ dx $$ The integral converges for 0<a<b. Let $a=k$, $b=2n+1,$ then sum over $k.$ We get $$ \frac{\pi}{2(2n+1)}S_n=\frac{1}{2}\sum_{k=1}^n \frac{1}{n} -\int_0^\infty dx \ \frac{e^{-(2n+1)x}}{\sinh{((2n+1)x)}} \ \sum_{k=0}^n \sinh{(2kx)} $$ $$ = \frac{1}{2}H_n -\int_0^\infty dx \ \frac{e^{-(2n+1)x}}{\sinh{((2n+1)x)}} \frac{\sinh{(nx)}\sinh{((n+1)x)}}{\sinh{x}}$$ The Harmonic sum notation has been used, and the finite sum withing the integral is well-known (really, it's just a variation of the geometric series). We now want to find how the integral behaves as $n \to \infty$. Scale $x \to x/(2n+1)$.
$$ I_n:=\frac{1}{2n+1}\int_0^\infty dx \ \frac{e^{-x}}{\sinh{(x)}} \frac{\sinh{(n/(2n+1)x)}\sinh{((n+1)/(2n+1)x)}}{\sinh{(x/(2n+1))}} $$ The pair of $\sinh$'s to a first approx. will behave like $\sinh{x/2}.$ Thus we'll expand to second order in $x^2$ as
$$ \sinh{(n/(2n+1)x)}\sinh{((n+1)/(2n+1)x)}= $$ $$= (\sinh{(x/2)})^2 \Big(\frac{4n(n+1)}{(2n+1)^2}+x^2\frac{n(n+1)}{3(2n+1)^4}+O((x/n)^4)\Big)$$ Likewise, $$ \frac{1}{\sinh{(x/(2n+1))} }=\frac{2n+1}{x}\Big(1- \frac{x^2}{6(2n+1)^2}+O((x/n)^4) \Big)$$
Therefore, to order $O(n^{-2}),$
$$I_n \sim \frac{4n(n+1)}{(2n+1)^2}\int_0^\infty \frac{e^{-x}(\sinh{(x/2)})^2} {x\ \sinh{x}} \Big(1-\frac{x^2}{12(2n+1)^2} \Big) \ dx$$ The integrals can be performed in Mathematica $$ \int_0^\infty \frac{e^{-x}(\sinh{(x/2)})^2}{x\ \sinh{x}} dx = \frac{1}{2} \log{(\pi/2)} $$
$$ \int_0^\infty \frac{x \ e^{-x}(\sinh{(x/2)})^2}{ \sinh{x}} dx = \frac{1}{2}(\frac{\pi^2}{6} - 1) $$
Collecting results,
$$ \boxed{ \ \ S_n \sim \frac{2n+1}{\pi}H_n + \frac{4n(n+1)}{\pi(2n+1)}\Big(\log{(2/\pi)}+ \frac{\pi^2/6-1}{12(2n+1)^2} \Big) \ \ } $$
The factor involving the ratios of polys in $n$ reduces to $2n,$ which is the same reduction as what the problem states.
For a numerical check with $n=10,$ the absolute error when using only the first term in the expansion (in big parentheses) is about 0.006%; for both terms, about 0.001%.