之前在b站上看到了一些介绍FFT的视频
《快速傅里叶变换(FFT)——有史以来最巧妙的算法?》
《这个算法改变了世界》
于是打算写一篇记录一下qwq(本博客中的截图基本上来源于第一个视频)
Fast Fourier Transform 是一种能在O(nlogn)O(nlogn)O(nlogn)的时间内将一个用系数表示的多项式转化成点值表示。
P(x)=p0+p1x+p2x2+...+pdxdP(x)=p_0+p_1x+p_2x^2+...+p_dx^dP(x)=p0+p1x+p2x2+...+pdxd
系数表示法:[p0,p1,p2,...,pd][p_0,p_1,p_2,...,p_d][p0,p1,p2,...,pd]
点值表示法:(x0,P(x0),(x1,P(x1),...,(xd,P(xd)){(x_0,P(x_0),(x_1,P(x_1),...,(x_d,P(x_d))}(x0,P(x0),(x1,P(x1),...,(xd,P(xd))(n+1个不同点可以确定一个n次的多项式:可以用矩阵列式求解)
当我们试图算两个多项式相乘的结果:
A(x)=a0+a1x+...+adxdA(x)=a_0+a_1x+...+a_dx^dA(x)=a0+a1x+...+adxd
B(x)=b0+b1x+...+bdxdB(x)=b_0+b_1x+...+b_dx^dB(x)=b0+b1x+...+bdxd
传统的做法需要O(n2)O(n^2)O(n2),而如果我们先求出若干个点的坐标(多少个点根据最后算出的会是几次多项式而定,比如A是三次多项式,B是五次多项式,算出来的C会是八次多项式,那么就需要找九个点),最后再根据这些点把多项式还原成系数表示,但是对于每个点,知道一个横坐标,需要(d+1)次计算才能知道这个点的纵坐标,总复杂度还是会达到O(n2)O(n^2)O(n2),达咩。
可不可以通过选取一些巧妙的点来减少我们的计算次数?
我们可以每次将多项式分成奇项和偶项,通过这样的划分可以便捷地算出横坐标互为相反数的一对点对的坐标,此时我们需要求出x2x^2x2分别代入Pe(x2)P_e(x^2)Pe(x2)和Po(x2)P_o(x^2)Po(x2)算出来的值,原问题变成了两个子问题,每个子问题的次数是原问题的1/21/21/2。我们想着能不能一直这么递归下去,直到分到底层(次数为1),再自底向上,层层归并求出答案?
实际上在递归步骤,我们假设了每个多项式我们都使用相反数对来求值,然而对两个子问题而言,每个求值点都是平方数,没办法继续了(悲.jpg)于是我们希望把新的求值点也弄成相反数对QAQ
(P.S.为了方便讨论,下文的n均为2的整数次幂,如果题目的n不是2的整数次幂,我们可以加一波操作~)
于是——复数大法好!
按照上面的逻辑推演,再向下递归一层我们需要两个相反点对,即要求出x4=1x^4=1x4=1的四个解;
再向下递归一层我们需要四个相反点对,即要求出x8=1x^8=1x8=1的八个解;
现在推广到d阶多项式,我们要先取n>d个点,(并且n等于2的整数次幂),我们为求解多项式乘积所选取的点就是1的n个n次方根。
单位根的性质:
性质一:ω2n2k=ωnkω^{2k}_{2n}=ω^{k}_nω2n2k=ωnk
性质二:ωnk+n/2=−ωnkω^{k+n/2}_n=−ω^k_nωnk+n/2=−ωnk(对应的点关于原点对称)
快速傅里叶变换的数学证明
设A(x)=a0+a1x+a2x2+...+an−1xn−1A(x) = a_0 + a_1x + a_2x^2 +...+a_{n-1}x^{n-1}A(x)=a0+a1x+a2x2+...+an−1xn−1,为求离散傅里叶变换,要把一个x=ωnkx = \omega_n^kx=ωnk代入。
考虑将A(x)A(x)A(x)的每一项按照下标的奇偶分成两部分:
A(x)=(a0+a2x2+...+an−2xn−2)+(a1x+a3x3+...+an−1xn−1)A(x) = (a_0 + a_2x^2 + ... + a_{n - 2}x^{n - 2}) + (a_1x + a_3x^3 + ... + a_{n-1}x^{n-1})A(x)=(a0+a2x2+...+an−2xn−2)+(a1x+a3x3+...+an−1xn−1)
设两个多项式:
A1(x)=a0+a2x+...+an−2xn2−1A_1(x) = a_0 + a_2x + ... + a_{n - 2}x^{\frac{n}{2} - 1}A1(x)=a0+a2x+...+an−2x2n−1
A2(x)=a1+a3x+...+an−1xn2−1A_2(x) = a_1 + a_3x + ... + a_{n - 1}x^{\frac{n}{2} - 1}A2(x)=a1+a3x+...+an−1x2n−1
则:A(x)=A1(x2)+xA2(x2)A(x) = A_1(x^2) + xA_2(x^2)A(x)=A1(x2)+xA2(x2)
假设k A(ωnk)=A1(ωn2k)+ωnkA2(ωn2k)=A1(ωn2k)+ωnkA2(ωn2k)\begin{align*} A(\omega_n^k) &= A_1(\omega_n^{2k}) + \omega_n^kA_2(\omega_n^{2k}) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) + \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{align*}A(ωnk)=A1(ωn2k)+ωnkA2(ωn2k)=A1(ω2nk)+ωnkA2(ω2nk) 所以,如果我们知道两个多项式A1(x)A_1(x)A1(x)和A2(x)A_2(x)A2(x)分别在(ωn20,ωn21,ωn22,...,ωn2n2−1)(\omega_{\frac{n}{2}}^{0}, \omega_{\frac{n}{2}}^{1}, \omega_{\frac{n}{2}}^{2}, ... , \omega_{\frac{n}{2}}^{\frac{n}{2} - 1})(ω2n0,ω2n1,ω2n2,...,ω2n2n−1)的值,就可以求出A(x)A(x)A(x)在ωn0,ωn1,ωn2,...,ωnn−1\omega_n^0, \omega_n^1, \omega_n^2, ..., \omega_n^{n-1}ωn0,ωn1,ωn2,...,ωnn−1处的值了,A1(x)A_1(x)A1(x)和A2(x)A_2(x)A2(x)是规模缩小一半的子问题,分治边界n=1n=1n=1. 结论:把多项式A(x)A(x)A(x)的离散傅里叶变换结果作为另一个多项式B(x)B(x)B(x)的系数,取单位根的倒数即ωn0,ωn−1,ωn−2,...,ωn−(n−1)\omega_{n}^{0}, \omega_{n}^{-1}, \omega_{n}^{-2}, ..., \omega_{n}^{-(n - 1)}ωn0,ωn−1,ωn−2,...,ωn−(n−1)作为xxx代入B(x)B(x)B(x),得到的每个数再除以nnn,得到的是A(x)A(x)A(x)的各项系数。 参考资料: 小学生都能看懂的FFT!!! 代码:咕咕,回头补
那么对于A(ωnk+n2)A(\omega_n^{k + \frac{n}{2}})A(ωnk+2n):
A(ωnk+n2)=A1(ωn2k+n)+ωnk+n2A2(ωn2k+n)=A1(ωn2k×ωnn)+ωnk+n2A2(ωn2k×ωnn)=A1(ωn2k)−ωnkA2(ωn2k)\begin{align*} A(\omega_n^{k + \frac{n}{2}}) &= A_1(\omega_n^{2k + n}) + \omega_n^{k + \frac{n}{2}}A_2(\omega_n^{2k + n}) \\ &= A_1(\omega_{\frac{n}{2}}^{k} \times \omega_n^n) + \omega_n^{k + \frac{n}{2}} A_2(\omega_{\frac{n}{2}}^{k} \times \omega_n^n) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) - \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{align*}A(ωnk+2n)=A1(ωn2k+n)+ωnk+2nA2(ωn2k+n)=A1(ω2nk×ωnn)+ωnk+2nA2(ω2nk×ωnn)=A1(ω2nk)−ωnkA2(ω2nk)
现在我们有了这若干个点的纵坐标,问题来到如何根据这些点的坐标还原答案多项式的系数,这需要用到IFFT(Inverse Fast Fourier Transform)
《快速傅里叶变换(FFT)——有史以来最巧妙的算法?》