How to compute this difficult looking integral

$x^3+x+1$ has a single real root $\alpha$ (close to $-0.7$) and two complex conjugate roots $\beta,\overline{\beta}$.
By Vieta's formulas, $\alpha+\beta+\overline{\beta}=0$, $\alpha\beta\overline{\beta}=-1$ and $\alpha(\beta+\overline{\beta})+\beta\overline{\beta}=1$.
By partial fraction decomposition, $$ \frac{1}{x^3+x+1}=\frac{1}{(x-\alpha)(x-\beta)(x-\overline{\beta})} = \sum_{\xi\in\{\alpha,\beta,\overline{\beta}\}}\frac{\text{Res}\left(\frac{1}{x^3+x+1},x=\xi\right)}{x-\xi} \tag{1}$$ or: $$ \frac{1}{x^3+x+1} = \sum_{\xi\in\{\alpha,\beta,\overline{\beta}\}}\frac{1}{(3\xi^2+1)(x-\xi)} \tag{2} $$ so that: $$ \int_{0}^{M}\frac{dx}{x^3+x+1}=\sum_{\xi\in\{\alpha,\beta,\overline{\beta}\}}\frac{\log(M-\xi)-\log(-\xi)}{(3\xi^2+1)} \tag{3}$$ and since the sum of the residues at $\alpha,\beta,\overline{\beta}$ is zero, $$ \int_{0}^{+\infty}\frac{dx}{x^3+x+1} = \color{red}{\frac{-\log(-\alpha)}{3\alpha^2+1}-2\log\|\beta\|\cdot\text{Re}\left(\frac{1}{3\beta^2+1}\right)+2\text{Arg}(-\beta)\cdot \text{Im}\left(\frac{1}{3\beta^2+1}\right)}. \tag{4} $$ Numerically, such integral is $0.921763372185057543420173329\ldots$.
It has to be real, obviously.

Using the methods of residues, you can prove that:

$$\int_0^{\infty}\frac{x^{-p}dx}{1+x} = \frac{\pi}{\sin(\pi p)}\tag{1}$$

This then allows you to evaluate the integral of any rational function from zero to infinity, you start with writing down the partial fraction decomposition of the integrand. The idea is then that after multiplying the expression by $x^{-p}$ we obtain an analytic function of $p$, the integrals of the terms of the partial fraction expansion of the form $\frac{a}{b + x} x^{-p}$ which for $b$ not on the negative real axis can be easily evaluated using (1), possibly using an analytic continuation to allow for complex $b$. If $b$ is on the negative real axis, then this can also be handled via analytic continuation; as long as the original integral is convergent, this method will always work.

Adding up the integrals of all the terms as a function of $p$ and taking the limit of $p$ to zero will yield the desired integral.