attachment:f5.pyx of days10/CodingSprint

   1 """
   2 Faugere's F_5 algorithm
   3 
   4 This implementation is along the lines of John Perry's pseudocode 
   5 and Singular implementation. It was inspired by Martin Albrecht's
   6 implementation and discussions with him and Ludovic Perret.
   7 
   8 See http://www.math.usm.edu/perry/Research/ for details.
   9 
  10 AUTHOR:
  11     -- 20081022 Simon King
  12 EXAMPLE:
  13     sage: R.<x,y,z> = PolynomialRing(GF(29))
  14     sage: I =  R* [3*x^4*y + 18*x*y^4 + 4*x^3*y*z + 20*x*y^3*z + 3*x^2*z^3, \
  15                    3*x^3*y^2 + 7*x^2*y^3 + 24*y^2*z^3, 
  16                    12*x*y^4 + 17*x^4*z + 27*y^4*z + 11*x^3*z^2]
  17     sage: J = I.homogenize()
  18 
  19     sage: f5 = F5() # original F5
  20     sage: gb = f5(J)
  21     Generator 2/3
  22     1 critical pairs in degree 7
  23     1 critical pairs in degree 9
  24     1 critical pairs in degree 11
  25     Generator 3/3
  26     1 critical pairs in degree 6
  27     2 critical pairs in degree 7
  28     5 critical pairs in degree 8
  29     7 critical pairs in degree 9
  30     11 critical pairs in degree 10
  31     10 critical pairs in degree 11
  32     WARNING: top-reduction to zero
  33     WARNING: reduction to zero
  34     WARNING: reduction to zero
  35     3 critical pairs in degree 12
  36     5 critical pairs in degree 13
  37     2 critical pairs in degree 14
  38     1 critical pairs in degree 16
  39     Note that we encountered 3 S-polynomials reducing to zero
  40     sage: f5.zeroes(), len(gb)
  41     (3, 18)
  42 """
  43 
  44 import sage
  45 import sage.all
  46 from sage.all import copy
  47 import sys
  48 from sage.structure.sequence import Sequence
  49 from sage.rings.ideal import Ideal
  50 #from sage.rings.polynomial.multi_polynomial_libsingular cimport MPolynomial_libsingular
  51 
  52 class CounterClass:
  53     def __init__(self):
  54         self.Gen = 0
  55         self.Del = 0
  56     def __call__(self,i):
  57         if i>0:
  58             self.Gen += 1
  59         else:
  60             self.Del += 1
  61 Counter = CounterClass()
  62 CounterB = CounterClass()
  63 
  64 def dehomogenize(I,R):
  65     "dehomogenize(I,R): dehomogenize an ideal I that was obtained by homogenizing an ideal of ring R"
  66     h = I[0].parent()('h')
  67     return ([R(X(h=1)) for X in I]*R).reduced_basis()
  68 
  69 cdef class DecoratedPolynomial:
  70     """
  71     A decorated polynomial is a polynomial p in an ideal I = <f_1,...,f_m> of a polynomial
  72     ring R, together with a 'signature' (u,i), where u is a monomial in R and 1 <= i <= m
  73     is an integer. The signature has the following meaning:
  74     * One can express p as an ideal combination $\sum_{k=1}^i t_k\cdot f_k$, and u is the leading
  75       term of $t_i$ of f_1,...,f_i.
  76     The signature allows for avoiding some useless critical pairs in the Buchberger algorithm.
  77     """
  78     cdef int nu
  79     cdef object mu  # the multiplier
  80     cdef tuple Emu  # exponents of the multiplier
  81     #cdef MPolynomial_libsingular p
  82     cdef object p
  83 
  84     def __init__(self,mu, int nu,  p):
  85         self.mu = mu
  86         self.Emu = mu.degrees()
  87         self.nu = nu
  88         self.p = p
  89         #self.rule = -1
  90         #self.hist = []
  91     #def __init__(self,mu, int nu, object p):
  92     #    pass
  93     def sig(DecoratedPolynomial self):
  94         return (self.mu,self.nu)
  95     def poly(DecoratedPolynomial self):
  96         return self.p
  97     #def info(self):
  98     #    return self.hist
  99 
 100     def __richcmp__(DecoratedPolynomial self, DecoratedPolynomial other, op):
 101         # < 0, <= 1, == 2, != 3, > 4, >= 5
 102         # Idea for applications: When considering a pair of decorated monomials and
 103         # reducing them, we always reduce the "bigger" one, 
 104         if op == 0:
 105             if self.nu<other.nu:
 106                 return True
 107             if (self.nu==other.nu) and (self.mu<other.mu):
 108                 return True
 109             return False
 110         if op == 1:
 111             if self.nu<other.nu:
 112                 return True
 113             if (self.nu==other.nu) and (self.mu<=other.mu):
 114                 return True
 115             return False
 116         if op == 2:
 117             return (self.nu==other.nu) and (self.Emu==other.Emu) and (self.p==other.p)
 118         if op == 3:
 119             return not ((self.nu==other.nu) and (self.Emu==other.Emu) and (self.p==other.p))
 120         if op == 4:
 121             if self.nu>other.nu:
 122                 return True
 123             if (self.nu==other.nu) and (self.mu>other.mu):
 124                 return True
 125             return False
 126         if op == 5:
 127             if self.nu>other.nu:
 128                 return True
 129             if (self.nu==other.nu) and (self.mu>=other.mu):
 130                 return True
 131             return False
 132 
 133 
 134 cdef class F5:
 135     cdef dict Rules       # rewriting rules (indexed by leading monomial, yielding a monomial ideal)
 136     cdef object Rgb       # reduced GB of the first few input generators
 137     cdef object RgbMon    # their leading ideal
 138     cdef list Inputgen    # the ideal generators we got, interreduced and increasingly sorted
 139     cdef list Basis       # list of decorated polynomials that yields self.Rgb
 140     cdef object Ring      # the polynomial ring we started with
 141     cdef object HRing     # homogenized ring
 142     cdef object One       # "One" in the ring we work with (may be homogenisation of self.Ring
 143     cdef int nvars        # number of variables of HRing
 144     cdef int homog        # 1, if the input was homogeneous; 0 otherwise (then we homogenize/dehomogenize)
 145     cdef int Zeroes       # Number of S-polynomials reducing to zero, to test the efficiency of 
 146                           # our implementation of the F5 criterion.
 147                           # We are not keeping track of the corresponding critical pairs.
 148     def __init__(F5 self):
 149         "F = F5() set up the machinery to compute Groebner bases with the F5 algorithm"
 150         self.Rules={}
 151         self.Rgb = None
 152         self.RgbMon = None
 153         self.Inputgen = []
 154         self.Basis = []
 155         self.Ring = None
 156         self.HRing = None
 157         self.One = None
 158         self.nvars = 0
 159         self.homog = 0
 160         self.Zeroes = 0
 161 
 162     #def __init__(self):
 163     #    pass
 164 
 165     #def __dealloc__(self):
 166     #    self.Rules={}
 167     #    self.Rgb = None
 168     #    self.RgbMon = None
 169     #    self.Inputgen = []
 170     #    self.Basis = []
 171     #    self.Ring = None
 172     #    self.HRing = None
 173     #    self.One = None
 174     #    self.nvars = 0
 175     #    self.homog = 0
 176     #    self.Zeroes = 0
 177     #    print "Deleted"
 178 
 179     def __call__(F5 self, I):
 180         """
 181         F(I): Apply the F5 algorithm to an ideal I and return its reduced Groebner basis
 182         """
 183         if not isinstance(I, sage.rings.polynomial.multi_polynomial_ideal.MPolynomialIdeal):
 184             I = Ideal(I)
 185         self.Ring = I.ring()
 186         self.Zeroes = 0 # S-polynomials reducing to zero.
 187                          # We do not keep track of the corresponding critical pair.
 188         self.Rgb            = Sequence([I.ring()(0)],I.ring()) # (reduced) GB of the first few input generators
 189         self.RgbMon = self.Rgb
 190         # TODO: is homogenization really needed?
 191         # Certainly it is needed in the F5 matrix version. But here??
 192         if not I.is_homogeneous():
 193             print "homogenizing the input"
 194             J = I.homogenize()
 195             self.homog = 0
 196         else:
 197             J = I
 198             self.homog = 1
 199         if J.gens()==(0,):
 200             self.Inputgen = [J.ring()(0)]     # reduced_basis would yield a segfault
 201         else:
 202             self.Inputgen       = list(J.reduced_basis())   # the ideal generators we got, sorted increasingly
 203         self.HRing = J.ring()
 204         self.Inputgen.sort()
 205         self.Rules          = {} # rewriting rules
 206         self.Basis          = [] # List of decorated polynomials that will eventually form a
 207                                  # a Groebner basis of J
 208         self.One            = J.ring()(1) # self.One.parent() =self.HRing != self.Ring, if there was a homogenisation
 209         self.nvars          = len(J.ring().gens())
 210         return self.basis()
 211 
 212     ## Some methods that reveal the C-attributes of self.
 213     def rules(self):
 214         return self.Rules
 215     def zeroes(self):
 216         return self.Zeroes
 217 
 218     #############
 219     ## Main procedure
 220     ## Return the previously computed basis, or do the computation when first called
 221     def basis(slf, maxgen=None):
 222         """
 223         F.basis(): Compute a reduced Groebner basis using the F5 algorithm.
 224         """
 225         cdef F5 self = slf
 226         if self.Basis:
 227             if self.homog:
 228                 return self.Rgb # or [X.poly for X in self.Basis]*self.HRing, currently without reduction
 229             else:
 230                 return dehomogenize(self.Rgb,self.Ring)
 231         if not self.homog:
 232             h = self.HRing('h')
 233         if self.Inputgen[0] == 0:
 234             return Sequence([self.HRing(0)], self.HRing)
 235         cdef DecoratedPolynomial newDP = DecoratedPolynomial(self.One,0,self.Inputgen[0]/self.Inputgen[0].lc())
 236         
 237         self.Rules[0]=[]
 238         self.Basis.append(newDP)
 239         cdef int i
 240         cdef int UpTo
 241         if maxgen is None:
 242             UpTo = len(self.Inputgen)
 243         else:
 244             UpTo = maxgen
 245         cdef int laststop = 0
 246         cdef DecoratedPolynomial X,Y
 247         for i in range(1,UpTo):
 248             print "Generator %d/%d"%(i+1,UpTo)
 249             # by induction, we have a reduced GB for the generators 1,...,i-1.
 250             self.Rgb = ([Y.p for Y in self.Basis]*self.HRing).reduced_basis()
 251             self.RgbMon = ([Y.p.lm() for Y in self.Basis]*self.HRing).reduced_basis()
 252             self.insert(i)
 253             for X in self.Basis[laststop:]:
 254                 # it it contains a constant polynomial, the ideal is trivial
 255                 if self.homog:
 256                     if X.p.degree()==0:
 257                         self.Rgb = Sequence([self.One], self.HRing)
 258                         return Sequence([self.Ring(1)], self.Ring)
 259                 else:
 260                     if X.p(h=1).degree()==0:
 261                         self.Rgb = Sequence([self.One], self.HRing)
 262                         return Sequence([self.Ring(1)], self.Ring)
 263             laststop = len(self.Basis)
 264         # Now, there only remains to toss out the reduced GB
 265         if self.Zeroes:
 266             print "Note that we encountered %d S-polynomials reducing to zero"%(self.Zeroes)
 267         self.Rgb = ([Y.p for Y in self.Basis]*self.HRing).reduced_basis()
 268         if self.homog:
 269             return self.Rgb
 270         else:
 271             return dehomogenize(self.Rgb,self.Ring)
 272 
 273     cdef insert(F5 self,int i):
 274         """F.insert(i): Compute a Groebner basis for the ideal spanned by generators #0,...,i,
 275            when a Groebner basis for generators #0,...,i-1 is already known.
 276         """
 277         # We may reduce the i-th generator f_i by everything that is in the "old" ideal (f_0,...,f_{i-1})
 278         f = self.Inputgen[i].reduce(self.Rgb)
 279         self.Rules[i]=[]
 280         if f.lc()==0:
 281             return
 282         cdef DecoratedPolynomial newDP = DecoratedPolynomial(self.One, i, f/f.lc())
 283 
 284         cdef int newindices = len(self.Basis)  # self.Basis[:newindices] was our starting point
 285         self.Basis.append(newDP)
 286         cdef int j,d
 287         cdef list P
 288         cdef list Pd,tmpL,X
 289         cdef DecoratedPolynomial Y,Z
 290         # generate new S-polynomials
 291         # The argument i is used for the F5 criterion (it appears in the signature of 
 292         # polynomials we want to reduce with)
 293         P=self.crit_pairs(newDP,newindices,i,newindices)  # only crit pairs involving newDP and old GB elements
 294         if not self.homog:
 295             h = self.HRing('h')
 296         while P:
 297             #print i
 298             #print "len P =",len(P)
 299             d = min([X[0].degree() for X in P])
 300             tmpL,Pd = [],[]
 301             for X in P:
 302                 if X[0].degree() == d:
 303                     Pd.append(X)
 304                 else:
 305                     tmpL.append(X)
 306             print "%d critical pairs in degree %d"%(len(Pd),d)
 307             P = tmpL
 308             S = self.SPolys(Pd)
 309             R = self.reduction(S, i, newindices)
 310             for Y in R:
 311                 self.Basis.append(Y)
 312                 if self.homog:
 313                     if Y.p.degree()==0:
 314                         return
 315                 else:
 316                     if Y.p(h=1).degree()==0:
 317                         return
 318                 P.extend(self.crit_pairs(Y,len(self.Basis),i,newindices))   # crit pairs for all of self.Basis (including the new ones)
 319 
 320 
 321     cdef list crit_pair(self, DecoratedPolynomial X, DecoratedPolynomial Y, int i, int newindices):
 322         """
 323         F.crit_pair(X,Y,i,newindices) returns (t,u_X,X,u_Y,Y) corresponding to a critical pair
 324             X,Y when necessary for the computation of a Groebner basis of (f_1,...,f_i). Here,
 325             the F5 criterion for label i is used. By switching X and Y, if necessary, we will
 326             have X<Y in the output. 
 327             Assumption: F.Basis[:newindices] (a private attribute) corresponds to a GB of 
 328                         (f_1,...,f_{i-1}).
 329         """
 330         R = self.HRing
 331         t_X = X.p.lm()
 332         t_Y = Y.p.lm()
 333         t = t_X.lcm(t_Y)
 334         u_X = R.monomial_quotient(t,t_X)
 335         u_Y = R.monomial_quotient(t,t_Y)
 336         umuX = X.mu*u_X
 337         if X.nu==i and self.topreducible(umuX):
 338             #X.hist.append(('topred',umuX,i))
 339             return [] # since the row with signature (umuX,nu) is kicked out by F5, we don't need the S-poly
 340         umuY = Y.mu*u_Y
 341         if (X.nu==Y.nu) and (umuX==umuY):  # found in Perry's f5_library.lib
 342             # I guess this is part of "rewrite criterion".
 343             # But I confess I don't understand why this is correct
 344             return []
 345         if Y.nu==i and self.topreducible(umuY):
 346             #Y.hist.append(('topred',umuY,i))
 347             return []
 348         # adopt minor optimization (check is_rewritable)
 349         if self.is_rewritable(u_X,X) or self.is_rewritable(u_Y,Y):
 350             return []
 351         # Eventually we try to reduce from "top" to "bottom".
 352         # We return the pair X,Y, if u_X*X is on top of u_Y*Y.
 353         # Hence, later we will replace u_Y*Y by the S-poly of X and Y
 354         if (X.nu>Y.nu) or ((X.nu==Y.nu) and (umuX>umuY)): 
 355             return [[t,Y,u_Y,X,u_X]]
 356         return [[t,X,u_X,Y,u_Y]]
 357 
 358     cdef list crit_pairs(self, DecoratedPolynomial X, int ub, int i, int newindices):
 359         "compare with crit_pair, but here Y runs on self.Basis[:ub]"
 360         cdef list CP = []
 361         cdef DecoratedPolynomial Y
 362         for Y in self.Basis[:ub]:
 363             CP.extend(self.crit_pair(X,Y,i,newindices))
 364         return CP
 365 
 366     cdef list SPolys(self,list P):
 367         """
 368         F.SPolys(P), P a list of critical pairs, given by a 5-tuple of 
 369           1) a monomial t (the least common multiple of the leading monomials of X.poly() and Y.poly()), 
 370           2) the cofactor u_X = t/lm(X.poly()) for X, 
 371           3) a decorated polynomial X,
 372           4) the cofactor u_Y = t/lm(X.poly()) for Y,
 373           5) a decorated polynomial Y.
 374         Moreover, (u_X*X)<(u_Y*Y).
 375         The S-polynomials are computed and tested for rewritability. The non-rewritable S-polynomials
 376         are added to the rewrite rules and are returned in a list.
 377         """
 378         # P is a list of tuples (lcm, u_X,X, u_Y,Y), with u_X*lm(X) = u_Y*lm(Y) = lcm
 379         P.sort() # which means increasing by lcm
 380         cdef list S = []
 381         cdef DecoratedPolynomial X,Y,Z
 382         for (t,X,u_X,Y,u_Y) in P:
 383             if (not self.is_rewritable(u_X,X)) and (not self.is_rewritable(u_Y,Y)):
 384                 s = (u_X*X.p*Y.p.lc() - u_Y*Y.p*X.p.lc()) # it will be reduced soon, in a different method...
 385                 if s!=0:
 386                     # We know that u_X*X is on top of u_Y*Y.
 387                     # Hence, thinking abouz the F5 matrix version, we would 
 388                     # create rows labeled (X.mu*u_X,X.nu) and (Y.mu*u_Y,Y.nu),
 389                     # and a reduction from top to bottom would yield
 390                     # the S-polynomial at (Y.mu*u_Y,Y.nu). 
 391                     # TODO: Is this thinking correct?
 392                     Z = DecoratedPolynomial(Y.mu*u_Y,Y.nu, s/s.lc())
 393                     #Z.hist.append(('S-poly',X,X.sig(),Y,Y.sig()))
 394                     self.add_rule(Z)
 395                     S.append(Z)
 396         S.sort() # by increasing signature, hence, from top to bottom in the F5 matrix version
 397         return S
 398 
 399     cdef list reduction(self,list S, int i, int newindices):
 400         """
 401         F.reduction(S,i,newindices): S a list of decorated polynomials, and F.Basis[:newindices] (a private
 402         attribute) provides a Groebner basis for the ideal generatd by the generators #0,...,i-1. 
 403         It is assumed that the list is sorted by increasing signature. 
 404         We will reduce each element of S by the GB and the elements of S of lower signature.
 405         Returns the list of completely reduced decorated polynomials.
 406         """
 407         cdef DecoratedPolynomial X,Y
 408         cdef list truncGB = list(self.Rgb)
 409         # We use self.Rgb, since it is a reduced GB for the ideal of the first i-1 indices.
 410         # Later on, we add to it new reduced S-polynomials.
 411         cdef list completed  = []
 412         cdef list newly_completed = []
 413         cdef list redo = []
 414         #sys.stdout.write('r')
 415         #sys.stdout.flush()
 416         while S:
 417             #sys.stdout.write('.')
 418             #sys.stdout.flush()
 419             X = S.pop(0) # this is the one with the smallest signature
 420             # By construction, the signature of all our input polynomials has index i:
 421             #sys.stdout.write('r')
 422             #sys.stdout.flush()
 423             X.p = X.p.reduce(truncGB)
 424             #sys.stdout.write('.')
 425             #sys.stdout.flush()
 426             #X.hist.append('reduction')
 427             if X.p.lc()!=0:
 428                 X.p = X.p/X.p.lc()
 429             newly_completed, redo = self.top_reduction(X,newindices,completed)
 430             completed.extend(newly_completed)
 431             truncGB.extend([Y.p for Y in newly_completed]) # hence basis is "old basis" plus "completed"
 432             S.extend(redo)
 433             S.sort() # by increasing signature, since reductions must be done from top to bottom.
 434             #print "S",len(S)
 435             #print "rule",len(self.Rules[i])
 436             #print "Basis", len(self.Basis)
 437         #sys.stdout.write('\n')
 438         #sys.stdout.flush()
 439         return completed
 440 
 441     cdef tuple top_reduction(self, DecoratedPolynomial X, int newindices,list completed):
 442         if X.p.lc()==0:
 443             print "WARNING: reduction to zero"
 444             self.Zeroes = self.Zeroes+1
 445             return ([],[])
 446         #X.hist.append(('test',X.p))
 447         #sys.stdout.write('t')
 448         #sys.stdout.flush()
 449         J = self.find_reductor(X,newindices,completed)
 450         #sys.stdout.write('.')
 451         #sys.stdout.flush()
 452         if not J:
 453             return ([X],[])
 454         cdef DecoratedPolynomial Y = J[1]
 455         
 456         u = J[0]
 457         p = X.p-u*Y.p
 458         #sys.stdout.write('\n')
 459         #sys.stdout.flush()
 460         if p.lc()!=0:
 461             p = p/p.lc()
 462         else:
 463             print "WARNING: top-reduction to zero"
 464             self.Zeroes = self.Zeroes+1 
 465             return ([],[])
 466         newmono = Y.mu*u
 467         # TODO: is it newmono<X.mu or newmono<=X.mu?
 468         # At least it seems there is an infinite loop if it is "<" only
 469         if (Y.nu<X.nu) or (Y.nu==X.nu and newmono<=X.mu):
 470             # i.e., (signature of Y)*u is not bigger, hence, it is on top of sig(X)
 471             # hence, X ought to be replaced
 472             X.p = p
 473             #X.hist.append(('top-red by',Y,'with',u,X.p.lm(),Y.p.lm(),X.p.lm()))
 474             #X.hist.append('redo1')
 475             return ([],[X])
 476         # Otherwise, (signature of Y)*u is below sig(X). Hence, 
 477         # the corresponding row (that is now being created) must be reduced
 478         # (which is done right on the spot).
 479         cdef DecoratedPolynomial Z = DecoratedPolynomial(newmono,Y.nu,p)
 480         
 481         #Z.hist.append(('top-red obtained from',X,Y))
 482         #X.hist.append(('top-red 2',Y,Z))
 483         self.add_rule(Z)
 484         return ([],[X,Z])
 485 
 486     cdef list find_reductor(self, DecoratedPolynomial X, int newindices, list completed):
 487         """
 488         F.find_reductor(X,i,C): returns [], or [[u,Y]] if there is a decorated polynomial Y 
 489         in self.Basis+C (private attribute) and X is *safely* top-reducible by Y with multiplier u. 
 490         That's to say, u*Y has the same leading monomial as X, Y is neither detected by the F5 
 491         criterion nor the rewritten criterion, and u*Y has smaller signature than X.
 492         """
 493         cdef DecoratedPolynomial Y
 494         t = X.p.lm()
 495         R = self.HRing
 496         for Y in self.Basis:
 497             if Y<X:
 498                 tt = Y.p.lm()
 499                 if R.monomial_divides(tt,t):
 500                     u = R.monomial_quotient(t,tt)
 501                     umu = Y.mu*u
 502                     # Y is a reductor of X if u*Y has different signature from X,
 503                     # (u,Y) is not rewritable, and umu is not topreducible.
 504                     if ((Y.nu!=X.nu) or (umu!=X.mu)) and (not self.topreducible(umu)) and (not self.is_rewritable(u,Y)):
 505                         return [u,Y]
 506         for Y in completed: 
 507             if Y<X:
 508                 tt = Y.p.lm()
 509                 if R.monomial_divides(tt,t):
 510                     u = R.monomial_quotient(t,tt)
 511                     umu = Y.mu*u
 512                     if ((Y.nu!=X.nu) or (umu!=X.mu)) and (not self.topreducible(umu)) and (not self.is_rewritable(u,Y)):
 513                         return [u,Y]
 514         return []
 515 
 516     def topreducible(self, u):
 517         """
 518         F.topreducible(u): Tests whether the F5 criterion would apply to decorated polynomials 
 519         with signature (u,i), where u is a monomial and self.Rgb (private attribute) provides
 520         a Groebner basis for the ideal spanned by the generators #0,...,i-1.
 521         """
 522         return u.reduce(self.RgbMon)==0
 523 
 524     cdef add_rule(self, DecoratedPolynomial X):
 525         "F.add_rule(X) adds the rewriting rule corresponding to the decorated polynomial X"
 526         #print "Add", X.nu, X.mu
 527         self.Rules[X.nu].append(X)
 528 
 529     def is_rewritable(self, u, DecoratedPolynomial X):
 530         """
 531         F.is_rewritable(u,X) tests if the is a rewriting (corresponding to Gaussian 
 532         elimination "from top") for u*X, where u is a monomial and X a decorated polynomial.
 533         """
 534         cdef DecoratedPolynomial Y
 535         cdef list R = self.Rules[X.nu]
 536         cdef int l = len(R)
 537         HR = self.HRing
 538         umu = X.mu*u
 539         while (1):
 540             l-=1
 541             if l<0:
 542                 return False
 543             Y = R[l]
 544             if Y is X:
 545                 return False
 546             if HR.monomial_divides(Y.mu,umu):
 547                 return True
 548         #return (X.mu*u).reduce([Y.mu for Y in self.Rules[X.nu][X.rule+1:]])==0