其实是复习笔记了…… Tab. 第 233 次复习
Part0. FFT 用来干什么
众所周知,用朴素算法实现 多项式乘法 $f(x)\cdot g(x)$:分别枚举 $f(x)$ 和 $g(x)$ 中的一个项,系数相乘、指数相加,得到答案中的一个项(的一部分)。
不妨设 $f(x),g(x)$ 的次数都是 $n$,那么上述的朴素算法就是 $O(n^2)$ 的。这样的时间复杂度是多项式乘法的一个瓶颈……
而 FFT 就是用来解决多项式乘法的问题的,可以把它的时间复杂度优化到 $O(n \log n)$。
实际上,FFT 并不是直接计算多项式乘法,而是把原来的多项式 $f(x),g(x)$ 在 $O(n\log n)$ 的复杂度内转换为它的 点值表示(后面会讲),而点值表示的多项式相乘的时间复杂度是 $O(n)$ 的。最后再用 $O(n\log n)$ 的时间复杂度把所得多项式的点值表示转化为一般形式。
Part1. 点值表示
如果完全没有了解过 FFT 的 reader 们看完 Part0 都有一点懵。那么现在就来科普一下什么是点值表示。
Part1/1. 什么是点值表示
对于一个 $n$ 次的多项式 $f(x)$,我们将 $n+1$ 个不同的 $x_i$ 代入这个 $f(x)$,可以得到 $n+1$ 个多项式的值 $y_i$。根据这 $n+1$ 个数对 $(x_i,y_i)$ 可以唯一确定一个 $n$ 次的多项式。
也就是说一个多项式与一个点值表示是一一对应的。
那么 FFT 完成的操作就是:
- 把已知的一个多项式转化成对应的点值表示;
- 把已知的点值表示转换成对应的多项式。
复杂度都是 $O(n\log n)$。
Part1/2. 点值表示的乘法
假设我们用 $\{x_0,x_1,x_2,…,x_n\}$ 来表示 $f(x)$ 为 $\{(x_0,f(x_0)),(x_1,f(x_1)),…,(x_n,f(x_n))\}$、表示 $g(x)$ 为 $\{(x_0,g(x_0)),(x_1,g(x_1)),…,(x_n,g(x_n))\}$。
那么 $f(x)\cdot g(x)$ 的点值表达式就是 $\{(x_0,f(x_0)\cdot g(x_0)),(x_1,f(x_1)\cdot g(x_1)),…,(x_n,f(x_n)\cdot g(x_n))\}$ ,其实就是对应项相乘,可见它的时间复杂度是 $O(n)$ 的。
Part2. 算法思想
Part2/1. 小科普:复数
是复数(complex)不是负数(negative)
首先了解到虚数 $i=\sqrt{-1}$。
然后一个复数 $\mathbf a$ 可以表示为实部 $x$ 和虚部 $y\times i$ 之和 —— $\mathbf a=x+y\times i$。对于一个实数,我们可以把它放在“一维数轴”上;那么对于一个复数,我们需要把它放在“二维数轴”,也就是直角坐标系上。
对于复数的乘法,比如 $\mathbf a=x+yi,\mathbf b=x’+y’i$ ,计算 $\mathbf {ab}$:
- 可以展开直接相乘:$\mathbf{ab}=(x+yi)(x’+y’i)=(xx’-yy’)+(xy’+x’y)i$ ;
- 可以从几何角度理解:
对于一个“向量”表示的复数,有如下定义(实际上不满足向量的运算法则):
那么两个复数相乘等于“幅角相加、模长相乘”
Part2/2. 单位复根
Hint. 复数的模长 $|\mathbf a|$(这个不是绝对值),表示把它表示成向量后的长度,即 $|\mathbf a|=\sqrt{x^2+y^2}$,
定义 n 次单位复根 $\mathbf \omega_n^i$ 为满足 $\mathbf z^n=1$ 的复数 $\mathbf z$。比如 $n=3$ 时,单位复根表现在直角坐标系上为:
可以发现 $n$ 次单位复根满足以下性质($\mathbf \omega_n^i$ 表示第 $i$ 个单位复根):
- $|\mathbf \omega_n^i|=1$(都在单位圆上);
- $\mathbf \omega_n^i$ 的幅角为 $\frac {2\pi i}{n}$;Tab. 这里的 $i$ 是下标不是虚数,角度为弧度制($360^\circ$ 弧度角为 $2\pi$)
- 总共有 $n$ 个单位复根,记为 $\omega_n^0,\omega_n^1,…,\omega_n^{n-1}$。
另外,还有一些比较明显的性质:
- $\omega_n^i\cdot\omega_n^j=\omega_n^{i+j}$;
- $\omega_n^{i+n}=\omega_n^i$;
- $\omega_n^i=\omega_{kn}^{ki}$;
- $\omega_n^i=-\omega_n^{i+n/2}$ (这里的“负”表示大小相等、方向相反)
以上这些结论都可以通过在单位圆上画出单位根来证明。
单位复根有什么用呢?因为 $n$ 次单位复根恰好有 $n$ 个,就可以把这 $n$ 个单位复根分别代入 $n-1$ 次多项式 $f(x)$ —— 得到点值表达 $\{(\omega_n^0,f(\omega_n^0)),(\omega_n^1,f(\omega_n^1)),…,(\omega_n^{n-1},f(\omega_n^{n-1}))\}$
Part2/3. 傅里叶正变换(一般形式 to 点值表达)
Tab. 以下思路来自 Tiw_Air_OAO 的 cnblogs 博客,但不尽相同:https://www.cnblogs.com/Tiw-Air-OAO/p/10162034.html
FFT 的实现方法是 分治:记多项式 $f(x)=a_0+a_1x+a_2x^2+…+a_{n-1}x^{n-1}$。
接下来把 $f(x)$ 的偶数次项、奇数次项分开:
那么可以得到 $f(x)=x\cdot f_o(x^2)+f_e(x^2)$ (直接展开这个式子就可以证明了)。
令 $p=n/2$,然后有以下结论:
所以,如果已经知道 $f_o(\omega_p^i),f_e(\omega_p^i)$,就可以 $O(1)$ 求出 $f(\omega_n^i)$ 和 $f(\omega_n^{i+p})$。实现方法就是递归:每次从求解 $f(\omega_n^i)$ 的单位复根到求解 $f_o(\omega_{n/2}^i)$ 和 $f_e(\omega_{n/2}^i)$ ;最后当多项式只剩两项(一次项和常数项)时直接代入计算。
由于 $n$ 每次会除以二,而要计算 $n$ 次,所以时间复杂度为 $O(n\log n)$。
Part2/4. 傅里叶逆变换(点值表达 to 一般形式)
其实傅里叶正变换是作的这样一个矩阵乘法:
(矩阵乘法是第一个矩阵的第 i 行依次乘第二个矩阵的第 j 列,得到结果的第 i 行、第 j 列)
记上面的第一个矩阵(系数矩阵)为 $V$。我们再定义一个矩阵 $D$:
于是计算 $D\times V$ 的第 $(i,j)$ 项:
① 当 $i=j$ 时:
② 当 $i\not=j$ 时:
Hint. 用等比数列……
因为 $\omega_n^{(j-i)n}=\omega_n^0=1$,所以 $(D\times V)_{i,j}=0$,于是发现:
所以我们把点值表达再乘上这个矩阵 $D$(就像正变换一样),最后再给每个值除以 $n$。
Part3. 离散傅里叶变换实现
虽然之前讲的算法思想是递归求解,然而实测的效率很低……于是改成 迭代(递推)。
先来研究一下递归的时候系数是被怎么分的,比如 n=8 时(下面的数字代表的是多项式 $f(x)$ 的系数,如 “2” 代表的是 $f(x)$ 的第 2 项的系数 $a_2$):
现在还看不出规律,换成二进制再看看:
原来的系数编号的顺序是 {000,001,010,011,100,101,110,111},最后得到的编号顺序为 {000,100,010,110,001,101,011,111}……惊奇地发现每个位置上的数字的编号都反过来了(比如原来的第 4 个系数 011 变成了 110)。所以只需要做一个二进制翻转即可。
Part3. 源代码
版子题:UOJ #34
1 | /*Lucky_Glass*/ |
The End
Thanks for reading~
Email: lucky_glass@foxmail.com,欢迎提问~