"Fossies" - the Fresh Open Source Software Archive

Member "NZMATH-1.2.0/nzmath/poly/termorder.py" (19 Nov 2012, 13940 Bytes) of package /linux/misc/old/NZMATH-1.2.0.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) Python source code syntax highlighting (style: standard) with prefixed line numbers. Alternatively you can here view or download the uninterpreted source code file. For more information about "termorder.py" see the Fossies "Dox" file reference documentation.

    1 """
    2 Term Order for polynomials.
    3 """
    4 
    5 import re
    6 import nzmath.compatibility
    7 import nzmath.ring as ring
    8 
    9 
   10 _INTERFACE_MSG = "%s is interface"
   11 
   12 class TermOrderInterface (object):
   13     """
   14     (abstract term order)
   15 
   16     A term order is primalily a function, which determines precedence
   17     between two terms (or monomials).  By the precedence, all terms
   18     are ordered.
   19 
   20     More precisely in terms of Python, a term order accepts two tuples
   21     of integers, each of which represents power indices of the term,
   22     and returns 0, 1 or -1 just like cmp built-in function.
   23 
   24     A TermOrder object provides not only the precedence function, but
   25     also a method to format a string for a polynomial, to tell degree,
   26     leading coefficients, etc.
   27     """
   28 
   29     _PLUS_MINUS = re.compile(r"(^-|\+ -)")
   30 
   31     def __init__(self, comparator):
   32         """
   33         'comparator' accepts two tuples of integers, each of which
   34         represents power indices of the term, and returns 0, 1 or -1
   35         just like cmp built-in function.
   36         """
   37         if type(self) is TermOrderInterface:
   38             raise NotImplementedError(_INTERFACE_MSG % self.__class__.__name__)
   39         self.comparator = comparator
   40 
   41     def cmp(self, left, right):
   42         """
   43         Compare two indices left and right and determine precedence by
   44         self.comparator.
   45         """
   46         return self.comparator(left, right)
   47 
   48     def bisect(self, array, elem, lo=0, hi=None):
   49         """
   50         Return the index where to insert item `elem' in a list `array',
   51         assuming array is sorted with the order.  (In the following,
   52         a refers `array' and x refers `elem')
   53 
   54         The return value i is such that all e in a[:i] have e <= x, and
   55         all e in a[i:] have e > x.  So if x already appears in the list,
   56         a.insert(x) will insert just after the rightmost x already there.
   57 
   58         Optional args lo (default 0) and hi (default len(a)) bound the
   59         slice of a to be searched.
   60 
   61         This function is based on the bisect.bisect_right of the Python
   62         standard library.
   63         """
   64         if hi is None:
   65             hi = len(array)
   66         while lo < hi:
   67             mid = (lo + hi) >> 1
   68             if self.cmp(elem, array[mid]) < 0:
   69                 hi = mid
   70             else: lo = mid + 1
   71         return lo
   72 
   73     def format(self, polynom, **kwds):
   74         """
   75         Return the formatted string of the polynomial.
   76         """
   77         return self.comparator(left, right)
   78 
   79     def bisect(self, array, elem, lo=0, hi=None):
   80         """
   81         Return the index where to insert item `elem' in a list `array',
   82         assuming array is sorted with the order.  (In the following,
   83         a refers `array' and x refers `elem')
   84 
   85         The return value i is such that all e in a[:i] have e <= x, and
   86         all e in a[i:] have e > x.  So if x already appears in the list,
   87         a.insert(x) will insert just after the rightmost x already there.
   88 
   89         Optional args lo (default 0) and hi (default len(a)) bound the
   90         slice of a to be searched.
   91 
   92         This function is based on the bisect.bisect_right of the Python
   93         standard library.
   94         """
   95         if hi is None:
   96             hi = len(array)
   97         while lo < hi:
   98             mid = (lo + hi) >> 1
   99             if self.cmp(elem, array[mid]) < 0:
  100                 hi = mid
  101             else: lo = mid + 1
  102         return lo
  103 
  104     def leading_coefficient(self, polynom):
  105         """
  106         Return the leading coefficient of polynomial 'polynom' with
  107         respect to the term order.
  108         """
  109         raise NotImplementedError(_INTERFACE_MSG % self.__class__.__name__)
  110 
  111     def leading_term(self, polynom):
  112         """
  113         Return the leading term of polynomial 'polynom' as tuple of
  114         (degree index, coefficient) with respect to the term order.
  115         """
  116         raise NotImplementedError(_INTERFACE_MSG % self.__class__.__name__)
  117 
  118 
  119 class UnivarTermOrder (TermOrderInterface):
  120     """
  121     term order for univariate polynomials.
  122 
  123     One thing special to univariate case is that powers are not tuples
  124     but integers.
  125     """
  126     def __init__(self, comparator):
  127         """
  128         UnivarTermOrder(comparator)
  129 
  130         'comparator' can be any callable that accepts two integers and
  131         returns 0, 1 or -1 just like cmp, i.e. if they are equal it
  132         returns 0, first one is greater 1, and otherwise -1.
  133         Theoretically acceptable comparator is only the cmp function.
  134         """
  135         TermOrderInterface.__init__(self, comparator)
  136 
  137     def format(self, polynom, varname="X", reverse=False):
  138         """
  139         Return the formatted string of the polynomial.
  140 
  141         - 'polynom' must be a univariate polynomial.
  142         - 'varname' can be set to the name of the variable (default to
  143           'X').
  144         - 'reverse' can be either True or False. If it's True, terms
  145           appear in reverse (descending) order.
  146         """
  147         degrees = _sort_with_cmp(
  148             [base for base in polynom.iterbases()],
  149             self.comparator)
  150         if reverse:
  151             degrees.reverse()
  152         str_terms = [("%s * %s ** %d" % (polynom[d], varname, d)) for d in degrees if polynom[d]]
  153         # constant
  154         if 0 in degrees and polynom[0]:
  155             const_term = str(polynom[0])
  156             if (hasattr(polynom, "getCoefficientRing") and
  157                 polynom[0] == polynom.getCoefficientRing().one):
  158                 const_term = "1"
  159             str_terms[str_terms.index("%s * %s ** 0" % (polynom[0], varname))] = const_term
  160         # degree 1
  161         if 1 in degrees and polynom[1]:
  162             str_terms[str_terms.index("%s * %s ** 1" % (polynom[1], varname))] = "%s * %s" % (polynom[1], varname)
  163         result = " + ".join(str_terms)
  164         # minus terms
  165         result = self._PLUS_MINUS.sub("- ", result)
  166         # coefficient is 1 (or -1)
  167         if hasattr(polynom, "getCoefficientRing"):
  168             one_times_x = re.compile(r"(^| )%s \* %s" % (polynom.getCoefficientRing().one, varname))
  169         else:
  170             one_times_x = re.compile(r"(^| )1 \* %s" % varname)
  171         result = one_times_x.sub(" " + varname, result)
  172         result = result.lstrip()
  173         return result
  174 
  175     def degree(self, polynom):
  176         """
  177         Return the degree of the polynomial 'polynom'.
  178         """
  179         if hasattr(polynom, "degree"):
  180             return polynom.degree()
  181         degree = -1
  182         for d in polynom.iterbases():
  183             if self.comparator(degree, d) < 0:
  184                 degree = d
  185         return degree
  186 
  187     def leading_coefficient(self, polynom):
  188         """
  189         Return the leading coefficient of polynomial 'polynom' with
  190         respect to the term order.
  191         """
  192         if hasattr(polynom, 'leading_coefficient'):
  193             return polynom.leading_coefficient()
  194         degree, lc = -1, 0
  195         for d, c in polynom:
  196             if self.comparator(degree, d) < 0:
  197                 degree, lc = d, c
  198         return lc
  199 
  200     def leading_term(self, polynom):
  201         """
  202         Return the leading term of polynomial 'polynom' as tuple of
  203         (degree, coefficient) with respect to the term order.
  204         """
  205         if hasattr(polynom, 'leading_term'):
  206             return polynom.leading_term()
  207         degree, lc = -1, 0
  208         for d, c in polynom:
  209             if self.comparator(degree, d) < 0:
  210                 degree, lc = d, c
  211         return degree, lc
  212 
  213     def tail_degree(self, polynom):
  214         """
  215         Return the least degree among all terms of the polynomial
  216         'polynom'.
  217 
  218         This method is EXPERIMENTAL.
  219         """
  220         if hasattr(polynom, "tail_degree"):
  221             return polynom.tail_degree()
  222         degree = -1
  223         for d in polynom.iterbases():
  224             if degree == -1 or self.comparator(degree, d) > 0:
  225                 degree = d
  226                 if degree == 0:
  227                     break
  228         return degree
  229 
  230 
  231 ascending_order = UnivarTermOrder(cmp)
  232 
  233 
  234 class MultivarTermOrder (TermOrderInterface):
  235     """
  236     A class of term orders for multivariate polynomials.
  237     """
  238     def __init__(self, comparator):
  239         """
  240         'comparator' accepts two tuples of integers, each of which
  241         represents power indices of the term, and returns 0, 1 or -1
  242         just like cmp built-in function.
  243         """
  244         self.comparator = comparator
  245 
  246     def format(self, polynom, varnames=None, reverse=False, **kwds):
  247         """
  248         Return the formatted string of the polynomial.
  249 
  250         An additional keyword argument 'varnames' is required to name
  251         variables.
  252         """
  253         if varnames is None:
  254             raise TypeError("keyword argument 'varnames' is required")
  255 
  256         bases = _sort_with_cmp(polynom.bases(), self.comparator)
  257         if reverse:
  258             bases.reverse()
  259 
  260         result = " + ".join([self._format_term((base, polynom[base]), varnames) for base in bases if polynom[base]])
  261         # minus terms
  262         result = self._PLUS_MINUS.sub("- ", result)
  263         # coefficient is 1 (or -1)
  264         if hasattr(polynom, "getCoefficientRing"):
  265             one_times = re.compile(r"(^| )%s \* " % polynom.getCoefficientRing().one)
  266         else:
  267             one_times = re.compile(r"(^| )1 \* ")
  268         result = one_times.sub(" ", result)
  269 
  270         result = result.lstrip()
  271         if not result:
  272             result = "0"
  273         return result
  274 
  275     def _format_term(self, term, varnames):
  276         """
  277         Return formatted term string.
  278 
  279         'term' is a tuple of indices and coefficient.
  280         """
  281         if not term[1]:
  282             return ""
  283         if term[1] == ring.getRing(term[1]).one:
  284             powlist = []
  285         else:
  286             powlist = [str(term[1])]
  287         for v, d in zip(varnames, term[0]):
  288             if d > 1:
  289                 powlist.append("%s ** %d" % (v, d))
  290             elif d == 1:
  291                 powlist.append(v)
  292 
  293         if not powlist:
  294             # coefficient is plus/minus one and every variable has degree 0
  295             return str(term[1])
  296         return " * ".join(powlist)
  297 
  298     def leading_coefficient(self, polynom):
  299         """
  300         Return the leading coefficient of polynomial 'polynom' with
  301         respect to the term order.
  302         """
  303         if hasattr(polynom, 'leading_coefficient'):
  304             return polynom.leading_coefficient()
  305         return polynom[self._max(polynom.bases())]
  306 
  307     def leading_term(self, polynom):
  308         """
  309         Return the leading term of polynomial 'polynom' as tuple of
  310         (degree index, coefficient) with respect to the term order.
  311         """
  312         if hasattr(polynom, 'leading_term'):
  313             return polynom.leading_term()
  314         max_indices = self._max(polynom.bases())
  315         return max_indices, polynom[max_indices]
  316 
  317     def _max(self, indices_list):
  318         """
  319         Return the maximum indices with respect to the comparator.
  320         """
  321         if not indices_list:
  322             raise ValueError("max() arg is an empty sequence")
  323         it = iter(indices_list)
  324         maxi = it.next()
  325         for indices in it:
  326             if self.comparator(maxi, indices) < 0:
  327                 maxi = indices
  328         return maxi
  329 
  330 
  331 def weight_order(weight, tie_breaker=None):
  332     """
  333     Return a comparator of weight ordering.
  334 
  335     The weight ordering is defined for arguments x and y that x < y
  336     if w.x < w.y (dot products) or w.x == w.y and tie breaker tells
  337     x < y.
  338 
  339     The option 'tie_breaker' is another comparator that will be used
  340     if dot products with the weight vector leaves arguments tie.  If
  341     the option is None (default) and a tie breaker is indeed necessary
  342     to order given arguments, a TypeError is raised.
  343     """
  344     def _order(mono, mial):
  345         if mono == mial:
  346             return 0
  347         wx = sum(w * x for (w, x) in zip(weight, mono))
  348         wy = sum(w * y for (w, y) in zip(weight, mial))
  349         if wx != wy:
  350             return cmp(wx, wy)
  351         return tie_breaker(mono, mial)
  352     return _order
  353 
  354 
  355 def _total_degree_lexicographic(left, right):
  356     """
  357     Total degree lexicographic (or graded lexicographic) term order :
  358       L < R iff
  359       (1) sum(li) < sum(ri) or
  360       (2) sum(li) = sum(ri) and
  361           there exists i s.t. l0 == r0, ..., l(i-1) == r(i-1), li < ri.
  362     """
  363     sum_left, sum_right = sum(left), sum(right)
  364     if sum_left != sum_right:
  365         return cmp(sum_left, sum_right)
  366     return cmp(left, right)
  367 
  368 def _total_degree_reverse_lexicographic(left, right):
  369     """
  370     Total degree reverse lexicographic (or graded reverse
  371     lexicographic) term order :
  372       L < R iff
  373       (1) sum(li) < sum(ri) or
  374       (2) sum(li) = sum(ri) and
  375           there exists i s.t. ln == rn, ..., l(i+1) == r(i+1), li > ri.
  376     """
  377     sum_left, sum_right = sum(left), sum(right)
  378     if sum_left != sum_right:
  379         return cmp(sum_left, sum_right)
  380     return cmp(right[::-1], left[::-1])
  381 
  382 
  383 lexicographic_order = MultivarTermOrder(cmp)
  384 total_degree_lexicographic_order = MultivarTermOrder(_total_degree_lexicographic)
  385 total_degree_reverse_lexicographic_order = MultivarTermOrder(_total_degree_reverse_lexicographic)
  386 
  387 
  388 ### for compatibility function
  389 def _sort_with_cmp(seq, mycmp):
  390     """
  391     Return the sorted seq by using the comparator function mycmp.
  392     """
  393     
  394     try:
  395         return sorted(seq, cmp=mycmp)
  396     except TypeError: # cmp is no longer available
  397         class internal_cmp:
  398             def __init__(self, obj):
  399                 self.obj = obj
  400             def __gt__(self, other):
  401                 return mycmp(self.obj, other.obj) > 0
  402             def __ge__(self, other):
  403                 return mycmp(self.obj, other.obj) >= 0
  404             def __lt__(self, other):
  405                 return mycmp(self.obj, other.obj) < 0
  406             def __le__(self, other):
  407                 return mycmp(self.obj, other.obj) <= 0
  408             def __cmp__(self, other):
  409                 return mycmp(self.obj, other.obj)
  410         return sorted(seq, key=internal_cmp)
  411