Decomposing a quartic polynomial into a difference of squares Given quartic polynomial $1 + x + x^2 +x^3 + x^4$, I want to find polynomials $p, q \in \mathbb R[x]$ such that 
$$1 + x + x^2 +x^3 + x^4 = p^2 (x) -  q^2 (x)$$
 A: With $$f(x)=x^4+x^3+x^2+x+1$$ choose $$g(x):=\frac{f(x)+1}{2}$$ $$h(x):=\frac{f(x)-1}{2}$$ then you have $$g^2-h^2=f$$
A: $$1 + x + x^2 +x^3 + x^4 = \frac 12 \begin{bmatrix} 1\\ x\\ x^2\end{bmatrix}^\top \begin{bmatrix} 2 & 1 & t\\ 1 & 2-2t & 1\\ t & 1 & 2\end{bmatrix} \begin{bmatrix} 1\\ x\\ x^2\end{bmatrix}$$
Since we want the matrix to be rank-$2$,
$$\det \begin{bmatrix} 2 & 1 & t\\ 1 & 2-2t & 1\\ t & 1 & 2\end{bmatrix} = 2 (t - 2) (t^{2} + t - 1) = 0$$
For $\color{blue}{t = 2}$, the matrix has one positive eigenvalue, one negative eigenvalue, and one zero eigenvalue. Using the spectral decomposition, eventually we obtain
$$\begin{aligned} p (x) &:= \left(\frac{\sqrt{\sqrt{11} + 1}}{4 \sqrt{11 + 3 \sqrt{11}}}\right) \left( \left(\sqrt{22} + 3 \sqrt{2}\right) x^{2} + 2 \sqrt{2} x + \left( \sqrt{22} + 3 \sqrt{2} \right)\right)\\ q (x) &:= \left(\frac{\sqrt{\sqrt{11} - 1}}{4 \sqrt{11 - 3 \sqrt{11}}} \right) \left( \left(\sqrt{22} - 3 \sqrt{2}\right) x^{2} - 2 \sqrt{2} x + \left( \sqrt{22} - 3 \sqrt{2} \right)\right)\end{aligned}$$
such that
$$1 + x + x^2 +x^3 + x^4 = p^2 (x) -  q^2 (x)$$

SymPy code
from sympy import *

t = Symbol('t', real=True)
x = Symbol('x', real=True)

M0 = Matrix([[ 2, 1, 0],
             [ 1, 2, 1],
             [ 0, 1, 2]])

M1 = Matrix([[ 0, 0, 1],
             [ 0,-2, 0],
             [ 1, 0, 0]])

M = M0 + t * M1

# find roots of the determinant
roots = solve(M.det(), t)

for r in roots:
    print "Root:       ", r.evalf()
    print "Determinant:", M.subs(t,r).det()
    print "Eigenvalues:", M.subs(t,r).eigenvals()
    print "\n"

# compute eigendecomposition    
V, D = M.subs(t,roots[0]).diagonalize()

# normalize eigenvectors
Q = V * ( diag(V[:,0].norm(), V[:,1].norm(), V[:,2].norm()) )**-1

p = collect(simplify(sqrt( D[1,1] / 2) * (Q[:,1].T * Matrix([1,x,x**2]))[0,0]), x)
q = collect(simplify(sqrt(-D[2,2] / 2) * (Q[:,2].T * Matrix([1,x,x**2]))[0,0]), x)

print "p: ", p
print "q: ", q

print "Difference of squares:", simplify(p**2 - q**2)

