Esta respuesta es sólo una sugerencia, la fórmula puede simplificarse.
Es dado $ \enspace\displaystyle S{k,n}|{gcd(k,n)=1}=\sum\limits{{\substack{j=1\n~\nmid~j}}}^\infty\frac{1}{j}\cot\left(\frac{k\pi}{n}j\right) \enspace $ $\enspace\displaystyle S{k,n}=\sum\limits_{m=1}^\infty\frac{am}{m}\enspace $ $\enspace\displaystyle a{m+ln}=am:=\cot\left(\frac{k\pi}{n}m\right)\enspace$ y $\enspace a{ln}=a_n:=0\enspace $, $\enspace l\in\mathbb{N}_0$.
Uso de $ \enspace\displaystyle E1(x):=\sum\limits{k=1}^\infty \frac{1}{k}e^{i2\pi kx}=i\pi\left(\frac{1}{2}-x\right)-\ln(2\sin(\pi x)) \enspace $ $\enspace 0<x href="https://www.fernuni-hagen.de/analysis/download/bachelorarbeit_aschauer.pdf" rel="nofollow noreferrer" y="">https://www.fernuni-hagen.de/analysis/download/bachelorarbeit_aschauer.pdf, página 4 (2.1),</x>
uno consigue $\enspace\displaystyle bm |{n\nmid m}=\frac{1}{n}\sum\limits_{v=1}^{n-1}a_v e^{-i2\pi \frac{m}{n}v}\enspace$ y % adicional $\enspace\displaystyle bn=\frac{1}{n}\sum\limits{v=1}^{n-1}av =0\enspace$y $\enspace\displaystyle \sum\limits{m=1}^\infty\frac{am}{m}=\sum\limits{m=1}^{n-1}b_m E_1\left(\frac{m}{n}\right)\enspace $.
Resultado:
$$S{k,n}=\frac{1}{n}\sum\limits{m=1}^{n-1} \left(i\pi\left(\frac{1}{2}-\frac{m}{n}\right)-\ln\left(2\sin\frac{m\pi}{n}\right)\right)\sum\limits_{v=1}^{n-1}\cot\left(\frac{k\pi}{n}v\right) e^{-i2\pi \frac{m}{n}v}$$