FFT 入门
推荐博客 : https://oi.men.ci/fft-notes/
卷积的理解 : https://www.zhihu.com/question/22298352?rf=21686447
题目链接 :http://uoj.ac/problem/34
这是一道模板题。
给你两个多项式,请输出乘起来后的多项式。
输入格式
第一行两个整数 nn 和 mm,分别表示两个多项式的次数。
第二行 n+1n+1 个整数,表示第一个多项式的 00 到 nn 次项系数。
第三行 m+1m+1 个整数,表示第二个多项式的 00 到 mm 次项系数。
输出格式
一行 n+m+1n+m+1 个整数,表示乘起来后的多项式的 00 到 n+mn+m 次项系数。
样例一
input
1 2
1 2
1 2 1
output
1 4 5 2
explanation
(1+2x)?(1+2x+x2)=1+4x+5x2+2x3(1+2x)?(1+2x+x2)=1+4x+5x2+2x3。
限制与约定
0≤n,m≤105,保证输入中的系数大于等于 0 且小于等于 9。
时间限制:1s1s
空间限制:256MB
题意 : 给你两个多项式的系数,从 0 到 n 给出,求这两个多项式相乘后的系数,从小到大输出
思路分析 : 裸的 FFT ,参考kuangbin 的板子
就是要注意以下数组的大小,main中的 len 是 2^k , 因此当m+n = 2e5 左右时,此时 2^k = 260000+ , 因此要注意数组的大小
代码示例:
#include<bits/stdc++.h> using namespace std; #define ll long long const int maxn = 2e5+63000; const double pi = acos(-1.0); int n, m; struct Complex{ double x, y; Complex (double _x=0, double _y=0):x(_x), y(_y){} Complex operator -(const Complex &b)const{ return Complex(x-b.x, y-b.y); } Complex operator +(const Complex &b)const{ return Complex(x+b.x, y+b.y); } Complex operator *(const Complex &b)const{ return Complex(x*b.x-y*b.y, x*b.y+y*b.x); } }; Complex x1[maxn], x2[maxn]; int sum[maxn]; void change(Complex y[], int len){ for(int i = 1, j = len/2; i < len-1; i++){ if (i < j) swap(y[i], y[j]); int k = len/2; while(j >= k){ j -= k; k /= 2; } if (j < k) j += k; } } void fft(Complex y[], int len, int on){ change(y, len); for(int h = 2; h <= len; h <<= 1){ Complex wn(cos(-on*2*pi/h), sin(-on*2*pi/h)); for(int j = 0; j < len; j += h){ Complex w(1, 0); for(int k = j; k < j+h/2; k++){ Complex u = y[k]; Complex t = w*y[k+h/2]; y[k] = u+t; y[k+h/2] = u-t; w = w*wn; } } } if (on == -1){ for(int i = 0; i < len; i++) y[i].x /= len; } } int main () { cin >> n >> m; int len = 1; while(len <= (n+m)) len <<= 1; for(int i = 0; i <= n; i++) scanf("%lf", &x1[i].x); fft(x1, len, 1); for(int i = 0; i <= m; i++) scanf("%lf", &x2[i].x); fft(x2, len, 1); for(int i = 0; i < len; i++) x1[i] = x1[i]*x2[i]; fft(x1, len, -1); for(int i = 0; i <= n+m; i++){ sum[i] = (int)(x1[i].x+0.5); // sum[] 是最后的答案 printf("%d%c", sum[i], i ==n+m?'\n':' '); } return 0; }
____________________________________________________________________________
int rev[maxl]; void get_rev(int bit)//bit表示二进制位数,计算一个数在二进制翻转之后形成的新数 { for(int i=0;i<(1<<bit);i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1)); } void fft(cd *a,int n,int dft)//n表示我的多项式位数 { for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i],a[rev[i]]); //中间的那个if保证了每个数做多只被交换了1次 //如果不写那么会有一些数被交换两次,导致最终的位置没有变 for(int step=1;step<n;step<<=1)//模拟一个合并的过程 { cd wn=exp(cd(0,dft*PI/step));//计算当前单位复根 for(int j=0;j<n;j+=step<<1) { cd wnk(1,0);//计算当前单位复根 for(int k=j;k<j+step;k++) {//蝴蝶操作 cd x=a[k]; cd y=wnk*a[k+step]; a[k]=x+y;//这就是上文中F(x)=G(x)+ωH(x)的体现 a[k+step]=x-y; //后半个“step”中的ω一定和“前半个”中的成相反数 //“红圈”上的点转一整圈“转回来”,转半圈正好转成相反数 //一个数相反数的平方与这个数自身的平方相等.. wnk*=wn; } } } if(dft==-1) for(int i=0;i<n;i++) a[i]/=n; //考虑到如果是IDFT操作,整个矩阵中的内容还要乘上1/n }
东北日出西边雨 道是无情却有情