Polynomial identity

edit

Here is some code which "linearizes" your problem:

def polynomial_identity(polys):
    monomials = set()
    for p in polys:
        monomials.update(p.monomials())
        base_ring = p.base_ring()
    monomial_order = list(monomials)
    m = []
    for p in polys:
        row = []
        for monomial in monomial_order:
            try:
                row.append(base_ring(p.coefficient(monomial)))
            except TypeError:
                row.append(base_ring(0))
        m.append(row)
    m = matrix(m)
    return m.kernel()

The result is a vector space which gives you all of your solutions. Your example

sage: polynomial_identity([x^2 + y^2, x*y + x*z, x*y*z])
Vector space of degree 3 and dimension 0 over Rational Field
Basis matrix:
[]

shows that (0,0,0) is the only answer. Another example:

sage: polynomial_identity([x+y, 2*x+2*y, 3*(x+y)])        
Vector space of degree 3 and dimension 2 over Rational Field
Basis matrix:
[   1    0 -1/3]
[   0    1 -2/3]