mjo/polynomial.py

   1 from sage.all import *
   2
   3 def multidiv(f, gs):
   4     r"""
   5     Divide the multivariate polynomial ``f`` by the ordered list of
   6     polynomials ``gs`` from the same ring.
   7
   8     INPUT:
   9
  10     - ``f`` -- A multivariate polynomial, the "numerator," that we try
  11       to express as a combination of the ``gs``.
  12
  13     - ``gs`` -- An ordered list of multipolynomials, the "denominators,"
  14       in the same ring as ``f``.
  15
  16     OUTPUT:
  17
  18     A pair, the first component of which is the list of "quotients," and
  19     the second of which is the remainder. All quotients and the
  20     remainder will live in the same ring as ``f`` and the ``gs``. If the
  21     ordered list of "quotients" is ``qs``, then ``f`` is the remainder
  22     plus the sum of all ``qs[i]*gs[i]``.
  23
  24     The remainder is either zero, or a linear combination of monomials
  25     none of which are divisible by any of the leading terms of the ``gs``.
  26     Moreover if ``qs[i]*gs[i]`` is nonzero, then the multidegree of ``f``
  27     is greater than or equal to the multidegree of ``qs[i]*gs[i]``.
  28
  29     SETUP::
  30
  31         sage: from mjo.polynomial import multidiv
  32
  33     ALGORITHM:
  34
  35     The algorithm from Theorem 3 in  Section 2.3 of Cox, Little, and O'Shea
  36     is used almost verbatim.
  37
  38     REFERENCES:
  39
  40     - David A. Cox, John Little, Donal O'Shea. Ideals, Varieties, and
  41       Algorithms: An Introduction to Computational Algebraic Geometry
  42       and Commutative Algebra (Fourth Edition). Springer Undergraduate
  43       Texts in Mathematics. Springer International Publishing, Switzerland,
  44       2015. :doi:`10.1007/978-3-319-16721-3`.
  45
  46     EXAMPLES:
  47
  48     Example 1 in Section 2.3 of Cox, Little, and O'Shea::
  49
  50         sage: R = PolynomialRing(QQ, 'x,y', order='lex')
  51         sage: x,y = R.gens()
  52         sage: f = x*y^2 + 1
  53         sage: gs = [ x*y + 1, y + 1 ]
  54         sage: (qs, r) = multidiv(f, gs)
  55         sage: (qs, r)
  56         ([y, -1], 2)
  57         sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f
  58         True
  59         sage: not any( g.lt().divides(m) for m in r.monomials()
  60         ....:          for g in gs )
  61         True
  62
  63     Example 2 in Section 2.3 of Cox, Little, and O'Shea::
  64
  65         sage: R = PolynomialRing(QQ, 'x,y', order='lex')
  66         sage: x,y = R.gens()
  67         sage: f = x^2*y + x*y^2 + y^2
  68         sage: gs = [ x*y - 1, y^2 - 1 ]
  69         sage: (qs, r) = multidiv(f, gs)
  70         sage: (qs, r)
  71         ([x + y, 1], x + y + 1)
  72         sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f
  73         True
  74         sage: not any( g.lt().divides(m) for m in r.monomials()
  75         ....:          for g in gs )
  76         True
  77
  78     TESTS:
  79
  80     If we get a zero remainder, then the numerator should belong to
  81     the ideal generated by the denominators::
  82
  83         sage: set_random_seed()
  84         sage: R = PolynomialRing(QQ, 'x,y,z')
  85         sage: x,y,z = R.gens()
  86         sage: s = ZZ.random_element(1,5).abs()
  87         sage: gs = [ R.random_element() for idx in range(s) ]
  88         sage: f = R.random_element(ZZ.random_element(10).abs())
  89         sage: (qs, r) = multidiv(f,gs)
  90         sage: r != 0 or f in R.ideal(gs)
  91         True
  92
  93     The numerator is always the sum of the remainder and the quotients
  94     times the denominators, and the remainder's monomials aren't divisible
  95     by the leading term of any denominator::
  96
  97         sage: set_random_seed()
  98         sage: R = PolynomialRing(QQ, 'x,y,z')
  99         sage: x,y,z = R.gens()
 100         sage: s = ZZ.random_element(1,5).abs()
 101         sage: gs = [ R.random_element() for idx in range(s) ]
 102         sage: f = R.random_element(ZZ.random_element(10).abs())
 103         sage: (qs, r) = multidiv(f,gs)
 104         sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f
 105         True
 106         sage: r == 0 or (not any( g.lt().divides(m) for m in r.monomials()
 107         ....:                     for g in gs ))
 108         True
 109
 110     """
 111     R = f.parent()
 112     s = len(gs)
 113
 114     p = f
 115     r = R.zero()
 116     qs = [R.zero()]*s
 117
 118     while p != R.zero():
 119         for i in range(0,s):
 120             division_occurred = false
 121             # If gs[i].lt() divides p.lt(), then this remainder will
 122             # be zero and the quotient will be in R (and not the
 123             # fraction ring, which is important).
 124             (factor, lt_r) = p.lt().quo_rem(gs[i].lt())
 125             if lt_r.is_zero():
 126                 qs[i] += factor
 127                 p -= factor*gs[i]
 128                 division_occurred = true
 129                 break
 130
 131         if not division_occurred:
 132             r += p.lt()
 133             p -= p.lt()
 134
 135     return (qs,r)