十分简明易懂的FFT
浅谈 FFT
FFT(快速傅里叶变换)概述
具体FFT的原理,我就不解释了
网上大佬讲得都比我好
说白了FFT直接作用就是计算两个多项式f(x)*g(x)的结果的系数
暴力做法很容易想O(n*n) , FFT用了一些数学上的方法把它优化成了 O(nlogn)
把这个式子按照奇偶性拆开:
for(int i=0;i<n;i+=2)l[i>>1]=x[i],r[i>>1]=x[i+1];
递归处理:
fft(l,n>>1,type);fft(r,n>>1,type);
构造单位根,合并式子:
wn(cos(2pi/n),sin(type2*pi/n))
for(int i=0;i>1;i++,w=wn) t=wr[i],x[i]=l[i]+t,x[i+(n>>1)]=l[i]-t;
递归版代码:
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <algorithm>
#include <vector>
#include <cstring>
#include <queue>
#include <cmath>
#include <complex>
using namespace std;
const int N = 3e6+5;
const double pi = acos(-1.0);
typedef complex<double> dob;
int n,m;
dob a[N],b[N];
inline int gi(){
int x=0,res=1;char ch=getchar();
while(ch>'9'||ch<'0'){if(ch=='-')res*=-1;ch=getchar();}
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return x*res;
}
void FFT(dob *A,int len,int f){
if(len==1)return;
dob wn(cos(2.0*pi/len),sin(f*2.0*pi/len)),w(1,0),t;
dob A0[len>>1],A1[len>>1];
for(int i=0;i<(len>>1);i++)A0[i]=A[i<<1],A1[i]=A[i<<1|1];
FFT(A0,len>>1,f);
FFT(A1,len>>1,f);
for(int i=0;i<(len>>1);i++,w*=wn){
t=w*A1[i];
A[i]=A0[i]+t;
A[i+(len>>1)]=A0[i]-t;
}
}
signed main(){
n=gi(),m=gi();
for(int i=0;i<=n;i++)a[i]=gi();
for(int i=0;i<=m;i++)b[i]=gi();
m+=n;
for(n=1;n<=m;n<<=1);
FFT(a,n,1);
FFT(b,n,1);
for(int i=0;i<=n;i++)a[i]*=b[i];
FFT(a,n,-1);
for(int i=0;i<=m;i++)
printf("%d ",int(a[i].real()/n+0.5));
}
把n在二进制下拆开考虑
就有了非递归版,常熟也小了很多
非递归版:
#include<bits/stdc++.h>
using namespace std;
const int N=3e6+5;
const double pi=acos(-1);
typedef complex<double> dob;
int n,m;
dob a[N],b[N];
inline int gi(){
int x=0,res=1;char ch=getchar();
while(ch>'9'||ch<'0'){if(ch=='-')res*=-1;ch=getchar();}
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return x*res;
}
int R[N],L;
void FFT(dob *A,int f){
for(int i=0;i<n;i++)if(i<R[i])swap(A[i],A[R[i]]);
for(int i=1;i<n;i<<=1){
dob wn(cos(pi/i),sin(f*pi/i));
for(int p=i<<1,j=0;j<n;j+=p){
dob w(1,0);
for(int k=0;k<i;k++,w*=wn){
dob x=A[j+k],y=w*A[j+k+i];
A[j+k]=x+y,A[j+k+i]=x-y;
}
}
}
}
signed main(){
n=gi(),m=gi();
for(int i=0;i<=n;i++)a[i]=gi();
for(int i=0;i<=m;i++)b[i]=gi();
m+=n; for(n=1;n<=m;n<<=1)L++;
for(int i=0;i<n;i++)R[i]=(R[i>>1]>>1)|((i&1)<<(L-1));
FFT(a,1),FFT(b,1);
for(int i=0;i<=n;i++)a[i]*=b[i];
FFT(a,-1);
for(int i=0;i<=m;i++)printf("%d ",(int)(a[i].real()/n+0.5));
}