多項式の最大公約数

というか、最大公約式?
整係数の範囲でLCDを求めるプログラムを作った。ヒマがあれば解説も後ほど。



from poly import Poly

def LCD(f, g):
if isinstance(f, Poly):
if isinstance(g, Poly):
return LCDPoly(f, g)
else:
return LCDPolyNum(f, g)
else:
if isinstance(g, Poly):
return LCDPolyNum(g, f)
else:
return LCDNum(f, g)

def LCDNum(a, b):
a = abs(a)
b = abs(b)
while 1:
mod = a % b
if mod == 0:
return b
a = b
b = mod

def LCDList(L):
n = len(L)
a = L[0]
for i in range(1, n):
if a == 0:
a = L[i]
elif L[i] != 0:
a = LCDNum(a, L[i])
if a == 1:
return 1
return a

def LCDPoly(f, g):
a = LCDNum(LCDList(f), LCDList(g))
if a > 1:
f /= a
g /= a

if g > f:
tmp = f
f = g
g = tmp
while 1:
g /= LCDList(g)
mod = f % g
if mod == 0:
return g * a
elif f == mod:
return a
f = g
g = mod

def LCDPolyNum(f, a):
return LCDNum(LCDList(f), a)

def LCDPolyList(L):
d = L[0]
for i in range(1, len(L)):
if d == 0:
d = L[i]
elif L[i] != 0:
d = LCD(d, L[i])
return d

f = Poly([8, -8, 2])
g = Poly([0, 4, -2])
h = Poly([-2, -2, 4])

print LCD(f, g)
print LCD(-12, -4)
print LCDList(h)
print LCD(f, 2)


# poly.py
class Poly(list):
def __add__(self, a):
res = Poly([])
m = len(self)
if isinstance(a, Poly):
n = len(a)
if m < n:
res[0:n] = a[0:n]
for i in range(m):
res[i] += self[i]
else:
res[0:m] = self[0:m]
for i in range(n):
res[i] += a[i]
else:
res[0:m] = self[0:m]
if len(res) == 0:
res.append(a)
else:
res[0] += a
res.normalize()
return res

def __radd__(self, a):
return self + a

def __sub__(self, a):
res = Poly([])
m = len(self)
if isinstance(a, Poly):
n = len(a)
if m < n:
res[0:n] = [ -x for x in a[0:n] ]
for i in range(m):
res[i] -= self[i]
else:
res[0:m] = self[0:m]
for i in range(n):
res[i] -= a[i]
else:
res[0:m] = self[0:m]
if len(res) == 0:
res.append(-a)
else:
res[0] -= a
res.normalize()
return res

def __rsub__(self, a):
return -self + a

def __mul__(self, a):
if isinstance(a, Poly):
m = len(self)
n = len(a)
res = Poly([ 0 for i in range(m + n - 1) ])
for i in range(m):
for j in range(n):
res[i+j] += self[i] * a[j]
return res
else:
return Poly([ x * a for x in self ])

def __rmul__(self, a):
return self * a

def __div__(self, a):
if isinstance(a, Poly):
m = len(self)
n = len(a)
res = Poly([ 0 for i in range(m - n + 1) ])
b = Poly([ x for x in self ])
for i in range(m - n, -1, -1):
d = b[i+n-1] / a[n-1]
for j in range(n):
b[i+j] -= d * a[j]
res[i] = d
else:
res = Poly([ x / a for x in self])

res.normalize()
return res

def __mod__(self, a):
m = len(self)
n = len(a)
res = Poly([ 0 for i in range(m - n + 1) ])
b = Poly([ x for x in self ])
for i in range(m - n, -1, -1):
d = b[i+n-1] / a[n-1]
for j in range(n):
b[i+j] -= d * a[j]
res[i] = d

b.normalize()
return b

def __neg__(self):
return Poly([ -x for x in self ])

def __eq__(self, a):
if isinstance(a, Poly):
if len(self) != len(a):
return False
for i in range(0, len(self)):
if self[i] != a[i]:
return False
return True
else:
if len(self) == 0:
return a == 0
elif len(self) == 1:
return self[0] == a
else:
return False

def __ne__(self, a):
return not self.__eq__(a)

def __gt__(self, a):
if isinstance(a, Poly):
if len(self) == len(a):
return abs(self[len(self)-1]) > abs(a[len(a)-1])
else:
return len(self) > len(a)
else:
return True

def normalize(self):
for i in range(len(self) - 1, -1, -1):
if self[i] == 0:
self[i:] = [ ]
else:
return

def __str__(self):
s = ""
n = len(self)
for i in range(0, n):
if self[i] == 0:
continue
if i == 0:
s = str(self[i])
else:
tmp = ""
if self[i] == 1:
pass
elif self[i] == -1:
tmp = "-"
else:
tmp = str(self[i])
if s != "" and self[i] > 0:
tmp = "+" + tmp
tmp += "x"
if i > 1:
tmp += "^" + str(i)
s += tmp
if s == "":
s = "0"

return s