from sage.all import * def multidiv(f, gs): r""" Divide the multivariate polynomial ``f`` by the ordered list of polynomials ``gs`` from the same ring. INPUT: - ``f`` -- A multivariate polynomial, the "numerator," that we try to express as a combination of the ``gs``. - ``gs`` -- An ordered list of multipolynomials, the "denominators," in the same ring as ``f``. OUTPUT: A pair, the first component of which is the list of "quotients," and the second of which is the remainder. All quotients and the remainder will live in the same ring as ``f`` and the ``gs``. If the ordered list of "quotients" is ``qs``, then ``f`` is the remainder plus the sum of all ``qs[i]*gs[i]``. The remainder is either zero, or a linear combination of monomials none of which are divisible by any of the leading terms of the ``gs``. Moreover if ``qs[i]*gs[i]`` is nonzero, then the multidegree of ``f`` is greater than or equal to the multidegree of ``qs[i]*gs[i]``. SETUP:: sage: from mjo.polynomial import multidiv ALGORITHM: The algorithm from Theorem 3 in Section 2.3 of Cox, Little, and O'Shea is used almost verbatim. REFERENCES: - David A. Cox, John Little, Donal O'Shea. Ideals, Varieties, and Algorithms: An Introduction to Computational Algebraic Geometry and Commutative Algebra (Fourth Edition). Springer Undergraduate Texts in Mathematics. Springer International Publishing, Switzerland, 2015. :doi:`10.1007/978-3-319-16721-3`. EXAMPLES: Example 1 in Section 2.3 of Cox, Little, and O'Shea:: sage: R = PolynomialRing(QQ, 'x,y', order='lex') sage: x,y = R.gens() sage: f = x*y^2 + 1 sage: gs = [ x*y + 1, y + 1 ] sage: (qs, r) = multidiv(f, gs) sage: (qs, r) ([y, -1], 2) sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f True sage: not any( g.lt().divides(m) for m in r.monomials() ....: for g in gs ) True Example 2 in Section 2.3 of Cox, Little, and O'Shea:: sage: R = PolynomialRing(QQ, 'x,y', order='lex') sage: x,y = R.gens() sage: f = x^2*y + x*y^2 + y^2 sage: gs = [ x*y - 1, y^2 - 1 ] sage: (qs, r) = multidiv(f, gs) sage: (qs, r) ([x + y, 1], x + y + 1) sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f True sage: not any( g.lt().divides(m) for m in r.monomials() ....: for g in gs ) True A solution ``g`` to Exercise 6 in Section 2.3 of Cox, Little, and O'Shea that lives in the ideal generated by ``f1`` and ``f2`` but which has nonzero remainder after division:: sage: R = PolynomialRing(QQ, 'x,y', order='deglex') sage: x,y = R.gens() sage: f1 = 2*x*y^2 - x sage: f2 = 3*x^2*y - y - 1 sage: I = R.ideal(f1,f2) sage: g = 2*y*f2 sage: g in I True sage: (qs,r) = multidiv(g,[f1,f2]) sage: r.is_zero() False Two solutions ``g`` to Exercise 7 in Section 2.3 of Cox, Little, and O'Shea that live in the ideal generated by ``f1``, ``f2``, and ``f3`` but which have nonzero remainders after division:: sage: R = PolynomialRing(QQ, 'x,y,z', order='deglex') sage: x,y,z = R.gens() sage: f1 = x^4*y^2 - z sage: f2 = x^3*y^3 - 1 sage: f3 = x^2*y^4 - 2*z sage: I = R.ideal(f1,f2,f3) sage: g = x^2*f3 sage: g in I True sage: (qs, r) = multidiv(g, [f1,f2,f3]) sage: r.is_zero() False sage: g = x*f2 sage: g in I True sage: (qs, r) = multidiv(g, [f1,f2,f3]) sage: r.is_zero() False Example 4 in Section 2.3 of Cox, Little, and O'Shea. This is the same as Example 2, except with the order of ``gs`` reversed:: sage: R = PolynomialRing(QQ, 'x,y', order='lex') sage: x,y = R.gens() sage: f = x^2*y + x*y^2 + y^2 sage: gs = [ y^2 - 1, x*y - 1 ] sage: (qs, r) = multidiv(f, gs) sage: (qs, r) ([x + 1, x], 2*x + 1) sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f True sage: not any( g.lt().divides(m) for m in r.monomials() ....: for g in gs ) True TESTS: If we get a zero remainder, then the numerator should belong to the ideal generated by the denominators:: sage: R = PolynomialRing(QQ, 'x,y,z') sage: x,y,z = R.gens() sage: s = ZZ.random_element(1,5).abs() sage: gs = [ R.random_element() for idx in range(s) ] sage: # hack for SageMath Trac #28855 sage: f = R(R.random_element(ZZ.random_element(10).abs())) sage: (qs, r) = multidiv(f,gs) sage: r != 0 or f in R.ideal(gs) True The numerator is always the sum of the remainder and the quotients times the denominators, and the remainder's monomials aren't divisible by the leading term of any denominator:: sage: R = PolynomialRing(QQ, 'x,y,z') sage: s = ZZ.random_element(1,5).abs() sage: gs = [ R.random_element() for idx in range(s) ] sage: # hack for SageMath Trac #28855 sage: f = R(R.random_element(ZZ.random_element(10).abs())) sage: (qs, r) = multidiv(f,gs) sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f True sage: r == 0 or (not any( g.lt().divides(m) for m in r.monomials() ....: for g in gs )) True Exercise 8 in Section 2.4 of Cox, Little, and O'Shea says that we should always get a zero remainder if we divide an element of a monomial ideal by its generators:: sage: R = PolynomialRing(QQ,'x,y,z') sage: gs = R.random_element().monomials() sage: I = R.ideal(gs) sage: # hack for SageMath Trac #28855 sage: f = R(I.random_element(ZZ.random_element(5).abs())) sage: (qs, r) = multidiv(f, gs) sage: r.is_zero() True """ R = f.parent() s = len(gs) p = f r = R.zero() qs = [R.zero()]*s while p != R.zero(): for i in range(0,s): division_occurred = false # If gs[i].lt() divides p.lt(), then this remainder will # be zero and the quotient will be in R (and not the # fraction ring, which is important). (factor, lt_r) = p.lt().quo_rem(gs[i].lt()) if lt_r.is_zero(): qs[i] += factor p -= factor*gs[i] division_occurred = true break if not division_occurred: r += p.lt() p -= p.lt() return (qs,r)