Note that we can try if we have a purely imaginary root, and indeed we have this particular factorisation: $$(x^2+a)(x^4-ax^2+a^2+1)=x^6+x^2+a^3+a$$
Which fits the orignal polynomial for $\ a^3+a=1$, and this cubic has at least one real root that we can calculate by Cardan's method (I skip to the result though...)
Le set $\alpha=\sqrt[3]{108+12\sqrt{93}}\ $ then $\ a=\dfrac{\alpha}6-\dfrac 2{\alpha}$
Now we get to factorize the remaining quartic, there are formulas for this that should lead you to a factorization in $\mathbb R[x]$ in product of quadratics, but considering there are no odd powers we can try to factorize under the form:
$$(x^2+bx+c)(x^2-bx+c)=x^4+(2c-b^2)x^2+c^2$$
This gives us by identification of the coefficients (note that we are interested in one possible factorization, so I just take the positive square root each time) $\begin{cases}c=\sqrt{a^2+1}\\b=\sqrt{2c+a}\end{cases}$
And there you have it:
$$x^6+x^2+1=(x^2+a)(x^2+bx+c)(x^2-bx+c)$$