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     A solution ``g`` to Exercise 6 in Section 2.3 of Cox, Little, and
  79     O'Shea that lives in the ideal generated by ``f1`` and ``f2`` but
  80     which has nonzero remainder after division::
  81
  82         sage: R = PolynomialRing(QQ, 'x,y', order='deglex')
  83         sage: x,y = R.gens()
  84         sage: f1 = 2*x*y^2 - x
  85         sage: f2 = 3*x^2*y - y - 1
  86         sage: I = R.ideal(f1,f2)
  87         sage: g = 2*y*f2
  88         sage: g in I
  89         True
  90         sage: (qs,r) = multidiv(g,[f1,f2])
  91         sage: r.is_zero()
  92         False
  93
  94     Two solutions ``g`` to Exercise 7 in Section 2.3 of Cox, Little, and
  95     O'Shea that live in the ideal generated by ``f1``, ``f2``, and ``f3``
  96     but which have nonzero remainders after division::
  97
  98         sage: R = PolynomialRing(QQ, 'x,y,z', order='deglex')
  99         sage: x,y,z = R.gens()
 100         sage: f1 = x^4*y^2 - z
 101         sage: f2 = x^3*y^3 - 1
 102         sage: f3 = x^2*y^4 - 2*z
 103         sage: I = R.ideal(f1,f2,f3)
 104         sage: g = x^2*f3
 105         sage: g in I
 106         True
 107         sage: (qs, r) = multidiv(g, [f1,f2,f3])
 108         sage: r.is_zero()
 109         False
 110         sage: g = x*f2
 111         sage: g in I
 112         True
 113         sage: (qs, r) = multidiv(g, [f1,f2,f3])
 114         sage: r.is_zero()
 115         False
 116
 117     Example 4 in Section 2.3 of Cox, Little, and O'Shea. This is the
 118     same as Example 2, except with the order of ``gs`` reversed::
 119
 120         sage: R = PolynomialRing(QQ, 'x,y', order='lex')
 121         sage: x,y = R.gens()
 122         sage: f = x^2*y + x*y^2 + y^2
 123         sage: gs = [ y^2 - 1, x*y - 1 ]
 124         sage: (qs, r) = multidiv(f, gs)
 125         sage: (qs, r)
 126         ([x + 1, x], 2*x + 1)
 127         sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f
 128         True
 129         sage: not any( g.lt().divides(m) for m in r.monomials()
 130         ....:          for g in gs )
 131         True
 132
 133     TESTS:
 134
 135     If we get a zero remainder, then the numerator should belong to
 136     the ideal generated by the denominators::
 137
 138         sage: set_random_seed()
 139         sage: R = PolynomialRing(QQ, 'x,y,z')
 140         sage: x,y,z = R.gens()
 141         sage: s = ZZ.random_element(1,5).abs()
 142         sage: gs = [ R.random_element() for idx in range(s) ]
 143         sage: # hack for SageMath Trac #28855
 144         sage: f = R(R.random_element(ZZ.random_element(10).abs()))
 145         sage: (qs, r) = multidiv(f,gs)
 146         sage: r != 0 or f in R.ideal(gs)
 147         True
 148
 149     The numerator is always the sum of the remainder and the quotients
 150     times the denominators, and the remainder's monomials aren't divisible
 151     by the leading term of any denominator::
 152
 153         sage: set_random_seed()
 154         sage: R = PolynomialRing(QQ, 'x,y,z')
 155         sage: s = ZZ.random_element(1,5).abs()
 156         sage: gs = [ R.random_element() for idx in range(s) ]
 157         sage: # hack for SageMath Trac #28855
 158         sage: f = R(R.random_element(ZZ.random_element(10).abs()))
 159         sage: (qs, r) = multidiv(f,gs)
 160         sage: r + sum( qs[i]*gs[i] for i in range(len(gs)) ) == f
 161         True
 162         sage: r == 0 or (not any( g.lt().divides(m) for m in r.monomials()
 163         ....:                     for g in gs ))
 164         True
 165
 166     Exercise 8 in Section 2.4 of Cox, Little, and O'Shea says that we
 167     should always get a zero remainder if we divide an element of a
 168     monomial ideal by its generators::
 169
 170         sage: set_random_seed()
 171         sage: R = PolynomialRing(QQ,'x,y,z')
 172         sage: gs = R.random_element().monomials()
 173         sage: I = R.ideal(gs)
 174         sage: # hack for SageMath Trac #28855
 175         sage: f = R(I.random_element(ZZ.random_element(5).abs()))
 176         sage: (qs, r) = multidiv(f, gs)
 177         sage: r.is_zero()
 178         True
 179
 180     """
 181     R = f.parent()
 182     s = len(gs)
 183
 184     p = f
 185     r = R.zero()
 186     qs = [R.zero()]*s
 187
 188     while p != R.zero():
 189         for i in range(0,s):
 190             division_occurred = false
 191             # If gs[i].lt() divides p.lt(), then this remainder will
 192             # be zero and the quotient will be in R (and not the
 193             # fraction ring, which is important).
 194             (factor, lt_r) = p.lt().quo_rem(gs[i].lt())
 195             if lt_r.is_zero():
 196                 qs[i] += factor
 197                 p -= factor*gs[i]
 198                 division_occurred = true
 199                 break
 200
 201         if not division_occurred:
 202             r += p.lt()
 203             p -= p.lt()
 204
 205     return (qs,r)