The poles of the meromorphic function $f(z)=z^2/(e^z-1)$ are all on the imaginary axis $z_k = 2\pi k i$. If we try to integrate along the perimeter of a circular sector, formed two straight lines ($L_1$ and $L_2$) and a circle arc $C_R$, given by:
$$L_1 = \{x\in\mathbb{R}| 0 \le x \le R\}\\ L_2 = \{z\in\mathbb{C}| z=\rho e^{(\frac{\pi}{2}-\alpha)i}, 0 \le \rho \le R\}\\ C_R = \{z\in\mathbb{C}| z=R e^{\theta i}, 0 \le \theta \le \frac{\pi}{2}-\alpha, 0< \alpha < \frac{\pi}{4} \}$$
Then, the Cauchy–Goursat theorem or the Residue theorem give:
$$\int_{L_1} \frac{x^2}{e^x-1}\text{d}x + \int_{C_R} \frac{z^2}{e^z-1}\text{d}z - e^{2\alpha i}\int_{L_2} \frac{\rho^2}{e^{\rho(\sin\alpha +i\cos \alpha)}-1}\text{d}z = 0$$
It can be proven that the first and second integrals are not zero, but the second one tends to zero for $R\to \infty$. So, the first and third integrals cancel each other out. In addition, I know the final result:
$$\int_0^\infty \frac{x^2}{e^x-1}\text{d}x = 2\zeta(3)$$
What is the "missing ingredient" here (assuming it is possible to compute the integral using complex integration)?