Multivariate polynomial division

[kuzand]

x, y = symbols('x y')
p = Poly(x**2+y, x, y, domain=QQ)
p0 = Poly(y, x, y, domain=QQ)
div(p, p0)  # or p.div(p0)

gives me 0 and x**2 + y as quotient and remainder respectively, whereas I was expecting 1 and x**2. I can get the correct result using e.g. reduced, but was wondering why it doesn't work using div on polynomials.

I have sympy 1.7.1.

[jksuom]

Division of multivariate polynomials depends on the ordering of the generators.

In [1]: p = Poly(x**2+y, y, x, domain=QQ)

In [2]: p0 = Poly(y, y, x, domain=QQ)

In [3]: div(p, p0)
Out[3]: (Poly(1, y, x, domain='QQ'), Poly(x**2, y, x, domain='QQ'))

[oscarbenjamin]

Should anything be changed here?

[kuzand]

Division of multivariate polynomials depends on the ordering of the generators.

In both cases the remainder should be x**2.

Have a look at the following:

_, x, y = ring('x y', domain=QQ, order='lex')  # or ilex
p =  x**2+y
p0 = y
p.div(p0)

(1, x**2)

or by using reduced we also get (1, x**2).

[kuzand]

_, x, y = ring('x y', domain=QQ, order='lex')
f =  x**2+y
f0 = y
f.div(f0)

(1, x**2)

_, x, y = ring('x y', domain=QQ, order='ilex')
f =  x**2+y
f0 = y
f.div(f0)

(1, x**2)

Both give the same result. Whereas

x, y = symbols('x y')
p = Poly(x**2+y, x, y)
p0 = Poly(y, x, y)
p.div(p0)

(0, x**2 + y)

x, y = symbols('x y')
p = Poly(x**2+y, y, x)
p0 = Poly(y, y, x)
p.div(p0)

(1, x**2)

However if we consider e.g. p = x**2+y and p0 = x - y, the results are consistent:

_, x, y = ring('x y', domain=QQ, order='lex')
p =  x**2+y
p0 = x - y
p.div(p0)

(x + y, y**2 + y)

x, y = symbols('x y')
p = Poly(x**2+y, x, y)
p0 = Poly(x - y, x, y)
p.div(p0)

(x + y, y**2 + y)

And for ilex:

_, x, y = ring('x y', domain=QQ, order='ilex')
p =  x**2+y
p0 = x - y
p.div(p0)

(-1, x + x**2)

x, y = symbols('x y')
p = Poly(x**2+y, y, x)
p0 = Poly(x - y, y, x)
p.div(p0)

(-1, x + x**2)