FFT目害学笔记

我的YY理解

FFT可以处理形如:

\[h(x)=g(x)f(x)->h(x)=\sum\limits_{i=0}^{n}\sum\limits_{j=0}^{m}g_if_jx^{i+j}=\sum\limits_{i=0}^{n+m}\sum\limits_{j=0}^{i}g_jf_{i-j}x^i \]

的卷积形式。暴力是\(O(n^2)\)而FFT可以做到\(O(nlogn)\)级别.
具体用到了复数在二维平面上的表示,形如\(cos\theta+isin\theta\)的形式加减乘,实现点值表示下的\(O(n)\)快速卷积。
以单位根为x代入到多项式中,由于其特殊性质就可以求得系数。但是只能求长度(最大次数)为\(2^n-1\)的数,这也没有关系,只取需要的项即可。
为了计算需要将下标按照奇偶性分到\(1\)~\(n/2\),\(n/2+1\)~\(n\),
这样会发现每个数被分到的位置就是二进制下的数倒转过来表示的数。
例如在\(n=3\)时,\(6(110)\)->\(3(011)\),\(n=4\)时\(5(0101)\)->\(10(1010)\)等。
于是可以预处理出来每个数要分到哪里去,FFT时直接交换即可。

多项式乘法

是FFT板子题。
放个码:

#include<cstring>
#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
namespace EMT{
	typedef long long ll;typedef double db;//(double)clock() / (double)CLOCKS_PER_SEC;
	#define pf printf
	#define F(i,a,b) for(register int i=a;i<=b;i++)
	#define D(i,a,b) for(register int i=a;i>=b;i--)
	inline int read(){int x=0,f=1;char ch=getchar();while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}while(ch>='0'&&ch<='9')x=x*10+ch-'0',ch=getchar();return x*f;}
	inline void file(){freopen("in.in","r",stdin);freopen("my.out","w",stdout);}
	inline int max(int a,int b){return a>b?a:b;}inline int min(int a,int b){return a<b?a:b;}
	inline void pi(int x){pf("%d ",x);}inline void pn(){pf("\n");}
	const int N=1e7+10;
	const db PI=acos(-1.0);
	struct comp{
		db x,y;
		friend comp operator +(comp a,comp b){return (comp){a.x+b.x,a.y+b.y};}
		friend comp operator -(comp a,comp b){return (comp){a.x-b.x,a.y-b.y};}
		friend comp operator *(comp a,comp b){return (comp){a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};}
	}a[N],b[N];
	int lim=1,n,m,ans[N],l,r[N];
	inline void fft(comp *a,int tp){
		F(i,0,lim-1)if(i<r[i])std::swap(a[i],a[r[i]]);
		for(int mid=1;mid<lim;mid<<=1){
			comp wn=(comp){cos(PI/mid),tp*sin(PI/mid)};
			for(int R=mid<<1,j=0;j<lim;j+=R){
				comp w=(comp){1,0};
				for(int k=0;k<mid;k++,w=w*wn){
					comp x=a[j+k],y=w*a[j+mid+k];
					a[j+k]=x+y;
					a[j+mid+k]=x-y;
				}
			}
		}
	}
	inline short main(){
		n=read(),m=read();
		while(lim<=n+m)lim<<=1,l++;
		F(i,0,n)a[i].x=read();
		F(i,0,m)b[i].x=read();
		F(i,0,lim-1)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
		fft(a,1),fft(b,1);
		F(i,0,lim)a[i]=a[i]*b[i];
		fft(a,-1);
		F(i,0,n+m)pi((int)(a[i].x/lim+0.5));
		return 0;
	}
}
signed main(){return EMT::main();}

快速傅立叶之二

将b数组倒置,发现

\[C_k=\sum\limits_{i=k}^{n}a_i\times b_{n+k-i} \]

因为\(b_{n+1}\)~\(b_{n+k}\)都是0
\(a_{n+1}\)~\(a_{n+k}\)都是0
所以C还可以写成

\[C_k=\sum\limits_{i=0}^{n+k}a_i \times b_{n+k-i} \]

\[D_k=\sum\limits_{i=0}^{k}a_i \times b_{k-i} \]

\[C_k=D_{n+k} \]

其中D满足卷积形式。于是可以用fft求出D下标+n输出即可。

上一篇:CF717A Festival Organization


下一篇:连续时间的马尔可夫链