FFT 快速傅里叶变换 学习笔记&

十分简明易懂的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));
}
原文地址:https://www.cnblogs.com/naruto-mzx/p/12116448.html