#!/usr/bin/python2
# -*- coding: utf-8 -*-
# †
from math import cos, sin, pi
def myfft(a, inv):
n = len(a)
if n == 1:
return a
m = n / 2
a0 = myfft(a[0::2], inv)
a1 = myfft(a[1::2], inv)
theta = 2 * pi / n
if inv:
theta *= -1
z = complex(cos(theta), sin(theta))
pz = 1
r = [None] * n
for i in xrange(n):
r[i] = a0[i%m] + pz * a1[i%m]
pz *= z
return r
def fft(a):
return myfft(a, False)
def ifft(a):
n = len(a)
b = myfft(a, True)
for i in xrange(n):
b[i] /= n
return b
def convol(a, b):
n = 1
while n < len(a) + len(b):
n <<= 1
a += [0] * (n - len(a))
b += [0] * (n - len(b))
a = fft(a)
b = fft(b)
c = [ai * bi for ai, bi in zip(a, b)]
c = ifft(c)
return c
|
|
#!/usr/bin/python2
# -*- coding: utf-8 -*-
# †
O = 3 # 原始根 \omega
M = (5 << 37) * 211 + 1
def myfmt(a, inv):
n = len(a)
if n == 1:
return a
m = n / 2
a0 = myfmt(a[0::2], inv)
a1 = myfmt(a[1::2], inv)
z = pow(O, (M-1)/n, M)
if inv:
z = pow(z, M-2, M) # 乗法逆元を求めている
pz = 1
r = [None] * n
for i in xrange(n):
r[i] = a0[i%m] + pz * a1[i%m]
r[i] %= M
pz = (pz * z) % M
return r
def fmt(a):
return myfmt(a, False)
def ifmt(a):
n = len(a)
b = myfmt(a, True)
inv = pow(n, M-2, M) # 乗法逆元を求めている
for i in xrange(n):
b[i] = b[i] * inv % M
return b
def convol(a, b):
n = 1
while n < len(a) + len(b):
n <<= 1
a += [0] * (n - len(a))
b += [0] * (n - len(b))
a = fmt(a)
b = fmt(b)
c = [ai * bi % M for ai, bi in zip(a, b)]
c = ifmt(c)
return c
|