之前在洛谷看某位神犇的NOI省选题解的时候呀,
注意到了FFT可以很有效地优化我们在做e.g. 多项式乘法相关的程序;
各位喜欢计算机科学/应用数学的客官们一定要学会它嗷!
(不太喜欢的话也康康呗www
快速傅里叶变换:即利用计算机计算离散傅里叶变换(DFT)的高效、快速计算方法的统称,简称FFT。快速傅里叶变换是1965年由J.W.库利和T.W.图基提出的。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少,特别是被变换的抽样点数N越多,FFT算法计算量的节省就越显著。
From Moe度娘
快速傅里叶变换(FFT)在计算机科学里常被用来加速多项式乘法;
对于此类问题,朴素高精度算法的时间复杂度为O(n^2);而使用FFT则可以达到O(n log n);
可能一个带有傅里叶老爷爷的名字会让很多同志们望而止步,但是实际上我们只要拥有一点点高中数学的知识就可以完全掌握呢!
需要的Pre-requisite有:
(聪明的你在高/初中时可以吊打的)三角函数
一点点复数的基础知识
odk那大家快上车啦!
# (I) 多项式的两种表达方法
大家耳熟能详的系数表达法
b. 容易理解的点值表示法
把一个(n-1)次的多项式放到平面直角坐标系里面,看成一个函数
把n个不同的x代入,会得出n个不同的y,在坐标系内就是n个不同的点
那么这n个点可以确定唯一的一个f(x),因为可以列出n个方程来求解f(x)的n个系数嘛,并且都有唯一的解
# (II) 在朴素算(bao)法(suan)下两种表达的乘法
系数表达法
对于两个(n-1)次也就是n项的多项式来说,两者的乘法相当于枚举每一个系数并相乘,可以得到最终的乘积;
aka 时间复杂度为O(n^2)
点值表示法
在取同样的n个x值的条件下,两个多项式的乘法运算只需要把相对应的 f(x)值相乘;
aka 时间复杂度为O(n)
瞬间对点值表达法肃然起敬有木有呀qwq
BUT! 注意到...咱们一般处理的多项式都是用系数表达法来表示的嘛,所以想要用点值表达法来计算的话必然涉及到两种表达法的转换问题:
DFT: *朴素*系数转点值的算法
IDFT: 也就是点值转系数的算法啦
Note: DFT - Discrete Fourier Transform离散傅里叶变换
IDFT - Inverse Discrete Fourier Transform
然而以上的两种算法的时间复杂度都是O(n^2)
也就是说我们十分费劲地把系数转换成了点值之后想要优化算法的行为使我们看起来像一个憨批??
这时聪明的你可能就要问了:
难道高精度乘法就只能O(n^2)爆算了么?
# (III) 一个hin简短的DFT介绍
对于系数转点值的操作,我们当然可以随便取任意n个x值代入计算
但是,正如上文提到过的,计算这些x的从0到(n-1)次方当然是O(n^2)的
然而,傅里叶姥爷说:
咱可不可以代入一组特定的x值,使得我们不用做这么多个次方运算呢?
By 莫名被cue的大神Fourier
也就是说,存不存在一些数字,它们的无论多少次方都是同一个(或有限的几个)数呢?
大家应该已经猜到了叭?
这些数字一定包括正负1;
了解过复数的同志们会发现:正负i也是满足条件的;
或者进一步来说,
在Re-Im这个复数平面里,以原点为圆心的单位圆上的所有点都是满足条件的!
原因就是,这些复数(或者说在复平面上对应的position vector)的modulus模长都为1,而不知道根据什么什么定理 这些复数的argument都代表了由实数1在复平面上的旋转角度;i次方的操作会使复数的模长变成原来的i次方,而1的i次方永远是1,也会使argument翻i倍,所以就是说在经历了次方的操作后的复数x,仍然是在单位圆上的;
对于一个n项(n-1)次的多项式,想要找到合适的x值,也就相当于把这个单位圆n等分,之后就会相应地得到n个对应的在复平面上的单位根x;
把得到的这些x代入进要转换为点值表示法的多项式即可
# (IV) FFT
虽然DFT算法想到了要选择一些很厉害的x,但是我们最终还是需要把单位根暴力代入并计算,时间复杂度仍是O(n^2)
但是呢,我们可以利用单位根的一些性质,使用Divide & Conquer分治来优化DFT,也就是我们要介绍的FFT算法
根据一大波根据(函数奇偶性和很多代换)的推导之后,可以知道:
对于一个n项(n-1)次的多项式f(x)来说,设所有x的偶数次方的和为f1(x),把所有奇数次方的和提出来一个x之后,设它为f2(x)
我们有
其中omega为代入的复数x值,k为对应的次方数,下标n为函数的项数
这时,聪明的你会发现:
这两个多项式后面一坨东西只有符号不同;
就是说,如果已知以上等式左边的两个函数值,即可同时求出等式右边的两个值;
这样,我们可以利用分治的思想,每次调用function时只扫当前前面一半的f(x)序列,即可得出后面一半序列的答案
注意到当n=1时,函数f(x)只有一个常数项,可以直接return
现在我们的时间复杂度已经变成O(n log n)啦嘿嘿
# (V) 你想要的FFT板子~
C++的FFT实现:
#include#define cp complexvoid fft(cp *a,int n,int inv) //inv => complex conjugate{ if (n==1)return; int mid=n/2; static cp b[MAXN]; fo(i,0,mid-1)b[i]=a[i*2],b[i+mid]=a[i*2+1]; fo(i,0,n-1)a[i]=b[i]; fft(a,mid,inv),fft(a+mid,mid,inv); fo(i,0,mid-1) { cp x(cos(2*pi*i/n),inv*sin(2*pi*i/n)); b[i]=a[i]+x*a[i+mid],b[i+mid]=a[i]-x*a[i+mid]; } fo(i,0,n-1)a[i]=b[i];}
# (VI) Conclusion
希望各位可以通过这篇文章学到一些新的知识w
之后也会更FFT的优化板子啦~
本人以后应该会不定时地更一些CS相关的科普或奇奇怪怪的东西!
CF还有洛谷的题解可能也会写一写x
(鞠躬