$$p := x^4 + x + 1 = \begin{bmatrix} x^2\\ x\\ 1\end{bmatrix}^\top \begin{bmatrix} 1 & 0 & -t\\ 0 & 2 t & 0.5\\ -t & 0.5 & 1\end{bmatrix} \begin{bmatrix} x^2\\ x\\ 1\end{bmatrix}$$
where $t \in \Bbb R$. Using Sylvester's criterion, we learn that the (symmetric) matrix above is positive semidefinite for $t = 0.5$. Using the Cholesky decomposition,
$$\begin{bmatrix} 1 & 0 & -0.5\\ 0 & 1 & 0.5\\ -0.5 & 0.5 & 1\end{bmatrix} = {\rm L} {\rm L}^\top$$
where
$${\rm L} = \begin{bmatrix} \color{blue}{1} & 0 & 0\\ 0 & \color{magenta}{1} & 0\\ \color{blue}{-\frac{1}{2}} & \color{magenta}{\frac{1}{2}} & \color{red}{\frac{\sqrt{2}}{2}}\end{bmatrix}$$
and, thus,
$$p = \left( \color{blue}{x^2 - \frac12} \right)^2 + \left( \color{magenta}{x + \frac12} \right)^2 + \left( \color{red}{\frac{\sqrt{2}}{2}} \right)^2 > 0$$
which is the exact same sum of squares (SOS) decomposition in this answer.
>>> from sympy import *
>>> t = symbols('t', real=True)
>>> M = Matrix([[ 1, 0, -t],
[ 0, 2*t, 1/2],
[-t, 1/2, 1]])
>>> L = M.subs(t,1/2).cholesky()
>>> L
Matrix([
[ 1, 0, 0],
[ 0, 1, 0],
[-1/2, 1/2, sqrt(2)/2]])
Related
polynomials sum-of-squares-method matrices matrix-decomposition cholesky-decomposition