阅读本文您不需要掌握的知识或阅读的书籍有
- 复变函数
- 算法导论
- 欧拉公式
阅读本文您需要具备
- 全日制小学生学历及其同等学历
- 知道复数的单位是根号负一,用 i 表示
0 概要
快速傅立叶可以在最坏时间 O(nlgn) 内计算两个多项式的乘积, 本文以
为例, 计算它们的乘积
1 最小知识
我们首先了解学习 fft 所需的最小知识 点值表示法 和 单位复根。
1.1 点值表示法
多项式 A 是由它的系数构成, 用系数来表示多项式称为系数表示法:
对于 A 来说,它是一个三次多项式,有 4 个系数,想要唯一的确定这个多项式,至少需要四个点。
所以对于有 n 个系数的多项式至少需要 n 个点才能够唯一确定。
用点
1.2 单位复根
我们不加证明的给出结论, 复数域上
例如,当 n 等于 8 时,有 8 个单位复根:
那么问题来了,我们都知道
我们可以使用小学数学的知识进行分析:
我们发现,单位复根都在复平面上以原点为圆心,半径是1的圆上。
对于单位复根
若以顺时针为正,
特别,当 n 是偶数时,
(相消引理)
2 fast fourier transform
- 计算 A 和 B 在单位复根处的点值表达式;O(nlgn)
- 点乘得到 C 的点值表达式;O(n)
- 插值,将点值表达式转为系数表达式。O(nlgn)
2.1 计算点值
C 的最高次为 6 次,所以最小有6个点才能确定。
所以至少要使用 8 次单位复根,即我们刚刚已经计算好的
直接算的话,复杂度至少是O{n^2},所以我们按照每项次数的次数的奇偶性分成两组来计算
(蝴蝶操作)
单位复根代入之
由复根的性质,转化为 4 次单位复根
没错,计算规模神奇的减少了一半
同理,也在B上使用同样的方法。
开局只有一把 Python 解释器,写个程序计算一下
def fft(para=[1,2,3,4,0,0,0,0]):
if 2 == len(para):
[a0, a1] = para
return [a0 + a1 * 1, a0 - a1 * 1]
n = len(para)
odd = para[::2]
odd_r = compute(odd)
eve = para[1::2]
eve_r = compute(eve)
res = []
for i in range(n):
xi = complex(math.cos(i*1.0/n * 2 * math.pi), math.sin(i*1.0/n * 2 * math.pi))
yi = odd_r[i%(n/2)] + xi * eve_r[i%(n/2)]
res.append(yi)
return res
解得
A = fft([1,2,3,4,0,0,0,0])
B = fft([4,6,2,7,0,0,0,0])
In [4]: A
Out[4]:
[(10+0j),
(-0.41421356237309426+7.242640687119286j),
(-2.0000000000000004-1.9999999999999996j),
(2.4142135623730954+1.2426406871192848j),
(-2+7.347880794884119e-16j),
(2.4142135623730945-1.2426406871192857j),
(-1.9999999999999991+2.0000000000000004j),
(-0.4142135623730968-7.242640687119285j)]
In [5]: B
Out[5]:
[(19+0j),
(3.292893218813454+11.192388155425117j),
(1.9999999999999991-0.9999999999999998j),
(4.707106781186548+7.1923881554251174j),
(-7+1.592040838891559e-15j),
(4.707106781186546-7.1923881554251174j),
(2.000000000000001+1.0000000000000002j),
(3.29289321881345-11.192388155425117j)]
2.2 点乘
这个没啥好说的
C_dot = []
for a, b in zip(A, B):
C_dot.append(a * b)
In [7]: C_dot
Out[7]:
[(190+0j),
(-82.42640687119285+19.213203435596448j),
(-5.999999999999998-1.9999999999999973j),
(2.4264068711928566+23.213203435596427j),
(14-8.327598234202001e-15j),
(2.4264068711928424-23.21320343559642j),
(-6.000000000000001+2.000000000000003j),
(-82.42640687119285-19.213203435596387j)]
可以得到
2.3 插值
将点值表达式转为系数表达式即插值。
$$
$$
即求转化矩阵
容易得证
所以插值和fft算法基本相同,直接上代码吧
def rfft():
if 2 == len(para):
[a0, a1] = para
return [a0 + a1 * 1, a0 - a1 * 1]
n = len(para)
odd = para[::2]
odd_r = fft(odd)
eve = para[1::2]
eve_r = fft(eve)
res = []
for i in range(n):
xi = complex(math.cos(((-i) % n)*1.0/n * 2 * math.pi), math.sin(((-i) % n)*1.0/n * 2 * math.pi))
yi = odd_r[(-i)%(n/2)] + xi * eve_r[(-i)%(n/2)]
res.append(yi/n)
return res
In [10]: rfft(C_dot)
Out[10]:
[(4+8.118390173882291e-15j),
(13.999999999999996+5.6847126901744435e-15j),
(25.999999999999996+1.3762094118649165e-17j),
(45-7.538383680317139e-15j),
(44-1.1206293398714794e-14j),
(29.000000000000004-3.197071506826809e-15j),
(28.000000000000004-3.783880944511323e-15j),
(3.552713678800501e-15+8.448827874285115e-15j)]
来手动算一下验证一下结果吧
3 reference
https://en.wikipedia.org/wiki/Discrete_Fourier_transform_%28general%29
若你觉得我的文章对你有帮助,欢迎点击上方按钮对我打赏
扫描二维码,分享此文章