"Fossies" - the Fresh Open Source Software Archive

Member "NZMATH-1.2.0/nzmath/poly/multivar.py" (19 Nov 2012, 13770 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 "multivar.py" see the Fossies "Dox" file reference documentation.

    1 """
    2 Class definitions of multivariate polynomials.
    3 
    4 The polynomials are immutable data types, under the public API.  If
    5 one tries to manipulate its underlying data attributes, immutability
    6 will be able to be broken.
    7 """
    8 
    9 from __future__ import division
   10 import logging
   11 import warnings
   12 import nzmath.ring as _ring
   13 import nzmath.poly.formalsum as formalsum
   14 
   15 
   16 _log = logging.getLogger('nzmath.poly.multivar')
   17 
   18 
   19 class TermIndices(object):
   20     """
   21     Indices of terms of multivariate polynomials.
   22     """
   23     def __init__(self, indices):
   24         """
   25         TermIndices(indices)
   26 
   27         The constructor does not check the validity of indices,
   28         such as integerness, nonnegativity, etc.
   29         """
   30         self._tuple = tuple(indices)
   31 
   32     def __hash__(self):
   33         """
   34         hash(indices)
   35         """
   36         return hash(self.__class__.__name__) ^ hash(self._tuple)
   37 
   38     def __len__(self):
   39         """
   40         len(indices)
   41         """
   42         return len(self._tuple)
   43 
   44     def __iter__(self):
   45         """
   46         iter(indices)
   47 
   48         Return iterator yielding successive elements.
   49         """
   50         return iter(self._tuple)
   51 
   52     def __repr__(self): # debug
   53         return repr(self._tuple)
   54 
   55     def __eq__(self, other):
   56         """
   57         self == other
   58         """
   59         if hash(self) != hash(other):
   60             return False
   61         return self._tuple == other._tuple
   62 
   63     def __ne__(self, other):
   64         """
   65         self != other
   66         """
   67         return not self.__eq__(other)
   68 
   69     def __le__(self, other):
   70         """
   71         self <= other
   72         """
   73         return self._tuple <= other._tuple
   74 
   75     def __lt__(self, other):
   76         """
   77         self < other
   78         """
   79         return self._tuple < other._tuple
   80 
   81     def __gt__(self, other):
   82         """
   83         self > other
   84         """
   85         return self._tuple > other._tuple
   86 
   87     def __ge__(self, other):
   88         """
   89         self >= other
   90         """
   91         return self._tuple >= other._tuple
   92 
   93     def __add__(self, other):
   94         """
   95         (i1, ..., in) + (j1, ..., jn) = (i1 + j1, ..., in + jn)
   96         """
   97         if len(self) != len(other):
   98             raise TypeError("different length indices")
   99         return self.__class__([i + j for (i, j) in zip(self, other)])
  100 
  101     def __sub__(self, other):
  102         """
  103         (i1, ..., in) - (j1, ..., jn) = (i1 - j1, ..., in - jn)
  104         """
  105         if len(self) != len(other):
  106             raise TypeError("different length indices")
  107         return self.__class__([i - j for (i, j) in zip(self, other)])
  108 
  109     def __mul__(self, scale):
  110         """
  111         (i1, ..., in) * s = (i1 * s, ..., in * s)
  112         """
  113         return self.__class__([i * scale for i in self])
  114 
  115     def __getitem__(self, index):
  116         """
  117         (i1, ..., in)[k] = ik
  118         """
  119         return self._tuple[index]
  120 
  121     def pop(self, pos):
  122         """
  123         Return the index at 'pos' and a new TermIndices object as the
  124         omitting-the-'pos' indices.
  125         """
  126         index = self.__class__(self._tuple[:pos] + self._tuple[pos + 1:])
  127         return self._tuple[pos], index
  128 
  129     def gcd(self, other):
  130         """
  131         Return the greatest common divisor.
  132 
  133         gcd((i1, ..., in), (j1, ..., jn)) = (min(i1, j1), ..., min(in, jn))
  134         """
  135         if len(self) != len(other):
  136             raise TypeError("different length indices")
  137         return self.__class__([min(i, j) for (i, j) in zip(self, other)])
  138 
  139     def lcm(self, other):
  140         """
  141         Return the least common multiple.
  142 
  143         lcm((i1, ..., in), (j1, ..., jn)) = (max(i1, j1), ..., max(in, jn))
  144         """
  145         if len(self) != len(other):
  146             raise TypeError("different length indices")
  147         return self.__class__([max(i, j) for (i, j) in zip(self, other)])
  148 
  149     def total_degree(self):
  150         """
  151         Return the total degree of indices, i.e. the sum of indices.
  152         """
  153         return sum(self)
  154 
  155 
  156 class PolynomialInterface(formalsum.FormalSumContainerInterface):
  157     """
  158     Base class for all multivariate polynomials.
  159     """
  160     def total_degree(self):
  161         """
  162         Return the maximum total degree of terms.
  163         """
  164         return max([b.total_degree() for b in self.iterbases()])
  165 
  166     def __neg__(self):
  167         """
  168         -self
  169         """
  170         return self.construct_with_default([(d, -c) for (d, c) in self])
  171 
  172     def __pos__(self):
  173         """
  174         +self
  175         """
  176         return self.construct_with_default(self._coefficients)
  177 
  178 
  179 class BasicPolynomial(PolynomialInterface):
  180     """
  181     The class for basic multivariate polynomials.
  182     """
  183 
  184     def __init__(self, coefficients, **kwds):
  185         """
  186         BasicPolynomial(coefficients [, keyword_arguments...])
  187 
  188         coefficients can be any dict initial values.
  189         """
  190         PolynomialInterface.__init__(self)
  191         self._coefficients = dict([(TermIndices(i), c) for (i, c) in dict(coefficients).iteritems()])
  192         if "number_of_variables" in kwds:
  193             self.number_of_variables = kwds.pop("number_of_variables")
  194         else:
  195             self.number_of_variables = 0
  196             for i in self.iterbases():
  197                 if len(i) > self.number_of_variables:
  198                     self.number_of_variables = len(i)
  199         self._init_kwds = kwds
  200 
  201     def __mul__(self, other):
  202         """
  203         self * other
  204 
  205         If type of other is Polynomial, do multiplication in ring.
  206         Otherwise, do scalar multiplication.
  207         """
  208         if isinstance(other, PolynomialInterface):
  209             return self.ring_mul(other)
  210         else:
  211             return self.scalar_mul(other)
  212 
  213     def __rmul__(self, other):
  214         """
  215         other * self
  216         """
  217         return self.scalar_mul(other)
  218 
  219     def ring_mul(self, other):
  220         """
  221         self * other
  222 
  223         Both self and other must have the same length tuples of
  224         indices for every term.
  225         """
  226         mul_coeff = {}
  227         for ds, cs in self:
  228             if not cs:
  229                 continue
  230             for do, co in other:
  231                 if not co:
  232                     continue
  233                 assert len(ds) == len(do)
  234                 indices = ds + do
  235                 if indices in mul_coeff:
  236                     mul_coeff[indices] += cs*co
  237                 else:
  238                     mul_coeff[indices] = cs*co
  239         return self.construct_with_default([(d, c) for (d, c) in mul_coeff.iteritems() if c])
  240 
  241     def scalar_mul(self, scale):
  242         """
  243         Return the result of scalar multiplication.
  244         """
  245         return self.construct_with_default([(i, c * scale) for (i, c) in self if c])
  246 
  247     def term_mul(self, term):
  248         """
  249         Return the result of multiplication with the given term.
  250         The term can be given as a tuple (degree indices, coeff)
  251         or as a Polynomial instance.
  252         """
  253         if isinstance(term, PolynomialInterface):
  254             degrees, coeff = iter(term).next()
  255         else:
  256             degrees, coeff = term
  257         return self.construct_with_default([(d + degrees, c * coeff) for (d, c) in self])
  258 
  259     def __pow__(self, index):
  260         """
  261         self ** index
  262 
  263         pow with three arguments is not supported by default.
  264         """
  265         if index < 0:
  266             raise ValueError("negative index is not allowed.")
  267         elif index == 0:
  268             indices = (0,)
  269             for i, c in self:
  270                 if c:
  271                     indices = (0,) * len(i)
  272                     one = _ring.getRing(c).one
  273                     break
  274             else:
  275                 one = 1
  276             return self.construct_with_default({indices: one})
  277         elif index == 1:
  278             return self
  279         elif index == 2:
  280             return self.square()
  281         # special polynomials
  282         if not self:
  283             return self
  284         elif len(self) == 1:
  285             return self.construct_with_default([(d*index, c**index) for (d, c) in self])
  286         # general
  287         for i, c in self:
  288             if c:
  289                 zero = (0,) * len(i)
  290                 one = _ring.getRing(c).one
  291                 break
  292         power_product = self.construct_with_default({zero: one})
  293         power_of_2 = self
  294         while index:
  295             if index & 1:
  296                 power_product *= power_of_2
  297             index //= 2
  298             if index:
  299                 power_of_2 = power_of_2.square()
  300         return power_product
  301 
  302     def square(self):
  303         """
  304         Return squared polynomial.
  305         """
  306         # zero
  307         if not self:
  308             return self
  309 
  310         polynomial = self.construct_with_default
  311         data_length = len(self)
  312         # monomial
  313         if data_length == 1:
  314             for i, c in self:
  315                 result = polynomial([(i * 2, c ** 2)])
  316         # binomial
  317         elif data_length == 2:
  318             (i1, c1), (i2, c2) = [(i, c) for (i, c) in self]
  319             result = polynomial({i1 * 2: c1**2, i1 + i2: c1*c2*2, i2 * 2: c2**2})
  320         # general (recursive)
  321         else:
  322             half = data_length >> 1
  323             coefficients = [(i, c) for (i, c) in self]
  324             left, right = polynomial(coefficients[:half], **self._init_kwds), polynomial(coefficients[half:])
  325             result = left.square() + left * right * 2 + right.square()
  326         return result
  327 
  328     def __hash__(self):
  329         """
  330         hash(self)
  331 
  332         Return the hash satisfying hash(self) == hash(other)
  333         if self == other.
  334         """
  335         hashvalue = hash(self.__class__.__name__)
  336         for i, c in self:
  337             hashvalue = hashvalue ^ (hash(i) | hash(c))
  338         return hashvalue & 0x7fffffff
  339 
  340     def __call__(self, target, value):
  341         """
  342         Substitute 'value' to 'target' index variable.
  343         If 'target' is a tuple of indices, it has to be sorted and
  344         'value' also has to be a tuple of the same length.
  345 
  346         Note that the result will not be a univar polynomial nor a
  347         scalar.
  348         """
  349         result = {}
  350         if isinstance(target, (int, long)):
  351             for i, c in self:
  352                 deg, index = i[target], i[:target] + (0,) + i[target + 1:]
  353                 if index in result:
  354                     result[index] += c * value ** deg
  355                 else:
  356                     result[index] = c * value ** deg
  357         elif len(target) == len(value):
  358             substituted = self
  359             for var, val in zip(target[::-1], value[::-1]):
  360                 substituted = substituted(var, val)
  361             return substituted
  362         else:
  363             raise TypeError("argument lengths mismatsch")
  364         return self.__class__([(i, c) for (i, c) in result.iteritems() if c], **self._init_kwds)
  365 
  366     def __len__(self):
  367         """
  368         Return the number of data entries.
  369         """
  370         return len(self._coefficients)
  371 
  372     def __getitem__(self, index):
  373         """
  374         Return the coefficient of specified index (tuple of degrees).
  375         If there is no term of index, return 0.
  376         """
  377         return self._coefficients.get(index, 0)
  378 
  379     def iterterms(self):
  380         """
  381         iterator for (degree, coefficient) pairs.
  382         """
  383         return self._coefficients.iteritems()
  384 
  385     def itercoefficients(self):
  386         """
  387         iterator for coefficients.
  388         """
  389         return self._coefficients.itervalues()
  390 
  391     def iterbases(self):
  392         """
  393         iterator for degrees.
  394         """
  395         return self._coefficients.iterkeys()
  396 
  397     def partial_differentiate(self, target):
  398         """
  399         Return the polynomial obtained by partial differentiation with
  400         the 'target' index variable.
  401         """
  402         partial = {}
  403         for i, c in self:
  404             index_diffed = list(i)
  405             target_degree = index_diffed[target]
  406             index_diffed[target] -= 1
  407             index_diffed = tuple(i)
  408             partial[index_diffed] = target_degree * c
  409         return self.construct_with_default([(i, c) for (i, c) in partial.iteritems() if c])
  410 
  411     def erase_variable(self, target=0):
  412         """
  413         Erase a variable from the polynomial.  The target variable is
  414         specified by the position in indices.
  415 
  416         The method takes no care about resulting polynomial type, i.e.
  417         the result remains as the same type even if their indices have
  418         length less than 2.
  419         """
  420         result = dict()
  421         for term, coeff in self:
  422             term = term[:target] + term[target + 1:]
  423             if term in result:
  424                 result[term] += coeff
  425             else:
  426                 result[term] = coeff
  427 
  428         return self.__class__([(d, c) for (d, c) in result.iteritems() if c],
  429                               number_of_variables=(self.number_of_variables - 1),
  430                               **self._init_kwds)
  431 
  432     def combine_similar_terms(self, target):
  433         """
  434         Combine similar terms and return the resulting univariate
  435         polynomial with polynomial coefficients in the form of list of
  436         (degree, coefficient) pairs.  The target variable is specified
  437         by the position in indices.
  438         """
  439         zero = self.construct_with_default(())
  440         result = {}
  441         for i, c in self:
  442             result[i[target]] = result.get(i[target], zero) + self.__class__([(i, c)], **self._init_kwds)
  443         for i, c in result.iteritems():
  444             result[i] = c.erase_variable(target)
  445         return result.items()
  446 
  447     def __repr__(self): # debug use
  448         return "BasicPolynomial(%s)" % repr(self._coefficients)
  449 
  450     def construct_with_default(self, maindata):
  451         """
  452         Create a new multivar polynomial of the same class with self,
  453         with given only the maindata and use copy of self's data if
  454         necessary.
  455         """
  456         return self.__class__(maindata,
  457                               **self._init_kwds)