I was using `sympy.sturm()`

to calculate the Sturm functions of a polynomial and came across an error for a 4th degree polynomial

```
using SymPy
begin
@vars x
# Polynomial coefficients
c4 = Sym[1; -4; 8; -8; 4] # Problem with Sympy sturm() ?
c3 = Sym[1; 2; 1; -2]
p4 = sympy.Poly(c4, x, domain="QQ")
p3 = sympy.Poly(c3, x, domain="QQ")
print("p3.sturm() == my_sturm(p3): ", p3.sturm() == my_sturm(p3),"\n")
print("p4.sturm() == my_sturm(p4): ", p4.sturm() == my_sturm(p4),"\n")
print("\n")
print("p4.sturm() :", p4.sturm(),"\n\n")
print("my_sturm(p4):", my_sturm(p4),"\n")
end
```

Output:

p3.sturm() == my_sturm(p3): true

p4.sturm() == my_sturm(p4): false

p4.sturm() :SymPy.Sym[Poly(x^2 - 2x + 2, x, domain=‘QQ’), Poly(2x - 2, x, domain=‘QQ’), Poly(-1, x, domain=‘QQ’)]

my_sturm(p4):SymPy.Sym[Poly(x^4 - 4x^3 + 8x^2 - 8x + 4, x, domain=‘QQ’), Poly(4x^3 - 12x^2 + 16x - 8, x, domain=‘QQ’), Poly(-x^2 + 2*x - 2, x, domain=‘QQ’)]

I believe my implementation

```
function my_sturm(p)
@vars x
n = p.degree()
pₛ = zeros(Sym,n+1)
pₛ[1] = p
pₛ[2] = p.diff()
i = 3
imax = 0
while sympy.degree(pₛ[i-1]) != 0 && i <= n+1
if pₛ[i-1].args[1] != 0
# quotient and remainder
q, r = sympy.div(pₛ[i-2], pₛ[i-1]) #, domain="QQ")
if r.args[1] != 0
pₛ[i] = -(pₛ[i-2] - pₛ[i-1]*q) # -r
imax = i
else
pₛ[i] = sympy.Poly([0],x)
end
end
i += 1
end
return pₛ[1:imax]
end
```

is giving the correct answer.

Does SymPy use a different definition of Sturm functions or am I missing something here?