## "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
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)
```