Closed form for $\sum\limits_{n=2}^{\infty}\frac1{n^3-1}$ I am investigating (just for fun) the sum 
$$S=\sum_{n=2}^{\infty}\frac1{n^3-1}$$
Wolfram Alpha gives me the 'value'
$$S=-\frac13\sum_{\{\omega\,:\,\omega^3+6\omega^2+12\omega+7=0\}}\frac{\psi_{0}(-\omega)}{\omega^2+4\omega+4}$$
Where $\psi_{0}(s)$ is the di-gamma function. I would like to know how this is found.
My attempts: recall that $$x^3-1=\prod_{k=1}^{3}(x-\alpha_k)$$
Where $\alpha_k=\exp\frac{2i\pi(k-1)}{3}$. Hence 
$$S=\sum_{n\geq2}\prod_{k=1}^3\frac1{n-\alpha_k}$$
I have shown in other posts that given some non-zero sequence $\{a_k:k=1,2,..,m\}$ where $j\neq k\iff a_j\neq a_k$ then 
$$\prod_{k=1}^{m}\frac1{x-a_k}=\sum_{k=1}^{m}\frac1{x-a_k}\prod_{j=1\\j\neq k}^{m}\frac1{a_k-a_j}$$
So 
$$S=\sum_{n\geq2}\bigg[\sum_{k=1}^{3}\frac1{n-\alpha_k}\prod_{j=1\\j\neq k}^{3}\frac1{\alpha_k-\alpha_j}\bigg]$$
Setting $b(k)=\prod_{j=1\\j\neq k}^{3}\frac1{\alpha_k-\alpha_j}$,
$$S=\sum_{n\geq2}\sum_{k=1}^3\frac{b(k)}{n-\alpha_k}$$
$$S=b(1)\sum_{n\geq2}\frac1{n-\alpha_1}+b(2)\sum_{n\geq2}\frac1{n-\alpha_2}+b(3)\sum_{n\geq2}\frac1{n-\alpha_3}$$
Then focusing on 
$$\begin{align}
S_k=&\sum_{n\geq2}\frac1{n-\alpha_k}\\
=&\sum_{n\geq2}\int_0^1x^{n-1}\frac{\mathrm dx}{x^{\alpha_k}}\\
=&\int_0^1\frac{1}{x^{\alpha_k-1}}\sum_{n\geq0}x^{n}\mathrm dx\\
=&\int_0^1x^{1-\alpha_k}(1-x)^{-1}\mathrm dx\\
=&\text{???}
\end{align}$$
for this final integral I considered using the beta function but that wouldn't work because $\Gamma(0)$ is undefined. Plus I'm not sure if $S_k$ even converges. 
As you can see, I'm stuck. Could I have a little help?
Thanks
Edit:
It can be seen here that there seems to be a pattern to evaluations of the function $$S(k)=\sum_{n\geq2}\frac1{n^{2k+1}-1},\qquad k\in\Bbb N$$
Which seem to be in the form 
$$S(k)=-\frac1{1+2k}\sum_{\{\omega\,:\,R_k(\omega)=0\}}\frac{\psi_0(-\omega)}{P_k(\omega)}$$
Where $R_k(x)$ is a polynomial of degree $k$ and $P_k(x)$ is a polynomial of degree $k-1$. 
 A: $$\sum_{n\geq 0}\frac{1}{(n+a)(n+b)}=\frac{\psi(a)-\psi(b)}{a-b} $$
leads by partial fraction decomposition to
$$\sum_{n\geq 0}\frac{1}{(n+a)(n+b)(n+c)} = \frac{1}{(a-b)(b-c)(c-a)}\left[(b-c)\psi(a)+(c-a)\psi(b)+(a-b)\psi(c)\right] $$
In your case $a,b,c$ are the opposite of the roots of $(x+2)^3-1$, and I'll let you finish the computations: $\psi(1)=-\gamma$, but $\frac{1}{2}(5\pm i\sqrt{3})$ are not special arguments for $\psi(x)=\frac{d}{dx}\log\Gamma(x)$.
A: Partial fraction decomposition tells us that 
$$\frac{1}{n^k-1}=\sum_{\omega^k=1} \frac{1}{n-\omega}\lim_{z\to\omega}\frac{z-\omega}{z^k-1}\overset{\text{L'H}}{=}\frac{1}{k}\sum_{\omega^k=1} \frac{\omega}{n-\omega}$$
and by the property that 
$$\psi(z+1)=\psi(z)+\frac{1}{z}$$
we see that
$$\frac{1}{n^k-1}=\frac{1}{k}\sum_{\omega^k=1} \omega\Big(\psi(n-\omega+1)-\psi(n-\omega)\Big)$$
Then we may sum both sides from 2 to infinity: 
$$\sum_{n=2}^{\infty}\frac{1}{n^k-1}=\frac{1}{k} \lim_{h \to \infty}\sum_{\omega^k=1} \omega\Big(\psi(h-\omega)-\psi(2-\omega)\Big)=-\frac{1}{k}\sum_{\omega^k=1} \omega \,\psi(2-\omega)$$
where the second equality is because $\psi(h)=\log(h)+O(1/h)$.

Arriving at the WolframAlpha result: 
Plug in $k=3$, and thus 
$$\begin{align*}&\quad\sum_{n=2}^{\infty}\frac{1}{n^3-1}\\
&=-\frac{1}{3}\sum_{\omega^3=1}\omega\,\psi(2-\omega)\\
&=-\frac{1}{3}\sum_{\omega^3=1}\frac{\psi(2-\omega)}{\omega^2}\quad\big(\omega^3=1\big)\\
&=-\frac{1}{3}\sum_{(\omega+2)^3=1}\frac{\psi(-\omega)}{(\omega+2)^2}\quad\big(\omega\mapsto\omega+2\big)\\
&=-\frac{1}{3}\sum_{\omega^3+6\omega^2+12\omega+7=0}\frac{\psi(-\omega)}{\omega^2+4\omega+4}
\end{align*}$$
which is what we wanted. 
