poly.py

@@ -1,115 +1,115 @@
-from fractions import gcd
+from fractions import gcd
-class Polynomial():
+class Polynomial():
-    def __init__(self, lst, mod):
+    def __init__(self, lst, mod):
-        self.poly = list(map(lambda x: x % mod, lst))    
+        self.poly = list(map(lambda x: x % mod, lst))    
-        self.mod = mod
+        self.mod = mod
-        self.normalize()
+        self.normalize()
-    def normalize(self):
+    def normalize(self):
-        while self.poly and self.poly[-1] == 0:
+        while self.poly and self.poly[-1] == 0:
-            self.poly.pop() 
+            self.poly.pop() 
-
+
-    #zwraca jednomian stopnia n 
+    #zwraca jednomian stopnia n 
-    @staticmethod
+    @staticmethod
-    def Monomial(n, c, mod):
+    def Monomial(n, c, mod):
-        zeros = [0]*n
+        zeros = [0]*n
-
+
-        zeros.append(c)
+        zeros.append(c)
-        return Polynomial(zeros, mod)
+        return Polynomial(zeros, mod)
-
+
-    def __add__(self, p2):
+    def __add__(self, p2):
-        p1 = self
+        p1 = self
-        len_p1, len_p2= len(p1.poly), len(p2.poly)
+        len_p1, len_p2= len(p1.poly), len(p2.poly)
-        res = [0] * max(len_p1, len_p2)
+        res = [0] * max(len_p1, len_p2)
-        if len_p1 > len_p2:
+        if len_p1 > len_p2:
-            for _ in range(len_p1-len_p2):
+            for _ in range(len_p1-len_p2):
-                p2.poly.append(0)
+                p2.poly.append(0)
-        else:
+        else:
-            for _ in range(len_p2-len_p1):
+            for _ in range(len_p2-len_p1):
-                p1.poly.append(0)
+                p1.poly.append(0)
-
+
-        for i in range(len(res)):
+        for i in range(len(res)):
-            res[i] = (p1.poly[i] + p2.poly[i]) % self.mod
+            res[i] = (p1.poly[i] + p2.poly[i]) % self.mod
-        return Polynomial(res, self.mod)
+        return Polynomial(res, self.mod)
-
+
-    def __sub__(self, p2):
+    def __sub__(self, p2):
-        p1 = self
+        p1 = self
-        res = []
+        res = []
-        len_p2 = len(p2.poly)
+        len_p2 = len(p2.poly)
-        for i in range(len(p1.poly)):
+        for i in range(len(p1.poly)):
-            if i < len_p2:
+            if i < len_p2:
-                res.append(p1.poly[i] - p2.poly[i] % self.mod)
+                res.append(p1.poly[i] - p2.poly[i] % self.mod)
-            else:
+            else:
-                res.append(p1.poly[i])
+                res.append(p1.poly[i])
-        return Polynomial(res, self.mod)
+        return Polynomial(res, self.mod)
-
+
-    def __mul__(self, p2):
+    def __mul__(self, p2):
-        res = [0]*(len(self.poly)+len(p2.poly)-1)
+        res = [0]*(len(self.poly)+len(p2.poly)-1)
-        for i, x1 in enumerate(self.poly):
+        for i, x1 in enumerate(self.poly):
-            for j, x2 in enumerate(p2.poly):
+            for j, x2 in enumerate(p2.poly):
-                res[i+j] += x1 * x2 % self.mod
+                res[i+j] += x1 * x2 % self.mod
-        return Polynomial(res, self.mod)
+        return Polynomial(res, self.mod)
-
+
-    def __eq__(self, p2):
+    def __eq__(self, p2):
-        p1 = self
+        p1 = self
-        return p1.poly == p2.poly and p1.mod == p2.mod
+        return p1.poly == p2.poly and p1.mod == p2.mod
-
+
-    def __pow__(self, n):
+    def __pow__(self, n):
-        p1 = self
+        p1 = self
-        for i in range(n):
+        for i in range(n):
-            p1 = p1 * p1
+            p1 = p1 * p1
-        return p1
+        return p1
-    
+    
-    def __truediv__(self, p2):
+    def __mod__(self, p2):
-        p1 = self
+        p1 = self
-        m = self.mod
+        m = self.mod
-        
+        
-        if len(p1.poly) < len(p2.poly):
+        if len(p1.poly) < len(p2.poly):
-            return p1
+            return p1
-        
+        
-        if len(p2.poly) == 0:
+        if len(p2.poly) == 0:
-            raise ZeroDivisionError
+            raise ZeroDivisionError
-            
+            
-        divisor_coeff = p2.poly[-1]
+        divisor_coeff = p2.poly[-1]
-        divisor_exp = len(p2.poly) - 1 
+        divisor_exp = len(p2.poly) - 1 
-        while len(p1.poly) >= len(p2.poly):   
+        while len(p1.poly) >= len(p2.poly):   
-            max_coeff_p1 = p1.poly[-1] #wspolczynnik przy najwyzszej potedze
+            max_coeff_p1 = p1.poly[-1] #wspolczynnik przy najwyzszej potedze
-            try:
+            try:
-              tmp_coeff = modDiv(max_coeff_p1, divisor_coeff, m)
+              tmp_coeff = modDiv(max_coeff_p1, divisor_coeff, m)
-            except ZeroDivisionError as e:
+            except ZeroDivisionError as e:
-                raise e
+                raise e
-            
+            
-            tmp_exp = len(p1.poly)-1 - divisor_exp
+            tmp_exp = len(p1.poly)-1 - divisor_exp
-            tmp = [0] * tmp_exp
+            tmp = [0] * tmp_exp
-            tmp.append(tmp_coeff)
+            tmp.append(tmp_coeff)
-            sub = Polynomial(tmp, m) * p2
+            sub = Polynomial(tmp, m) * p2
-            p1 = p1 - sub
+            p1 = p1 - sub
-        p1.normalize()
+        p1.normalize()
-        return Polynomial(p1.poly, m)
+        return Polynomial(p1.poly, m)
-
+
-    def poly_gcd(self, p2):
+    def poly_gcd(self, p2):
-        p1 = self
+        p1 = self
-        try:
+        try:
-            divisible = p2
+            divisible = p2
-        except ZeroDivisionError as e:
+        except ZeroDivisionError as e:
-            raise e
+            raise e
-        if p2.poly == []:
+        if p2.poly == []:
-            return p1
+            return p1
-        
+        
-        return p2.poly_gcd(p1 / p2)
+        return p2.poly_gcd(p1 % p2)
-
+
-def modDiv(a, b, m): # a*b^-1 (mod m)
+def modDiv(a, b, m): # a*b^-1 (mod m)
-    if gcd(b, m) != 1:
+    if gcd(b, m) != 1:
-        raise ZeroDivisionError
+        raise ZeroDivisionError
-    else:
+    else:
-        return (a * modinv(b, m)) % m
+        return (a * modinv(b, m)) % m
-#rozszerzony algorytm euklidesa
+#rozszerzony algorytm euklidesa
-def egcd(a, b):
+def egcd(a, b):
-    if a == 0:
+    if a == 0:
-        return (b, 0, 1)
+        return (b, 0, 1)
-    else:
+    else:
-        g, y, x = egcd(b % a, a)
+        g, y, x = egcd(b % a, a)
-        return (g, x - (b // a) * y, y)
+        return (g, x - (b // a) * y, y)
-
+
-def modinv(a, m):
+def modinv(a, m):
-    g, x, y = egcd(a, m)
+    g, x, y = egcd(a, m)
     return x % m