快速傅里叶变换 FFT 【bzoj2194】 快速傅立叶之二

题目大意:
请计算C[k]=sigma(a[i]*b[i-k]) 其中 k < = i < n ，并且有 n < = 10 ^ 5。 a,b中的元素均为小于等于100的非负整数。

题目分析：（看题目猜题意：快速傅里叶变换FFT）
我们通过傅里叶变换，实际上是对两个高次函数求了卷积。
以下是来自度娘的解释。

我们可以发x和h的下标的和是相同的。
所以说我们使用做傅里叶变换其实就做了这样一个事er：
把所有下标和为0的两项的积求和；
把所有下标和为1的两项的积求和;
把所有下标和为2的两项的积求和；
……

现在再来看C[k]=sigma(a[i]*b[i-k]) ，我们发现它们的下标和并不相等，但是它们的下标差相等啊。那我们可以把a[i]变成a[-i]，可是我们不能用负的下标，所以可以把a[i]变成a[n-i-1]（因为下标是从0开始的）,也就是把a数组翻转过来啦。

这样我们可以发现所求的东西的两项的下标和就为定值了。
（n-i-1）+（i-k） = n-k-1
这个定值就是 n-k-1（k < n ），所以要把我们把这两个数组用傅里叶变换求卷积，
然后得到的数组的前n项倒序输出。

代码如下：

#include<cstdio>
#include<algorithm>
#include<cmath>
#define N 1<<18
using namespace std;
const double pi=acos(-1);
const double DFT=2.0,IDFT=-2.0;

int pos[N];
void initialization(int len)
{
    for(int i=0;i<len;i++)
    {
        pos[i]=pos[i>>1]>>1;
        if(i&1) pos[i]|=len>>1;
    }
    return;
}


struct complex{
    double a,b;
    complex(const double &a=0.0,const double &b=0.0):a(a),b(b){}
    complex operator + (const complex &c) {return complex(a+c.a,b+c.b); }
    complex operator - (const complex &c) {return complex(a-c.a,b-c.b); }
    complex operator * (const complex &c) {return complex(a*c.a-b*c.b,a*c.b+b*c.a); }
}factor1[N],factor2[N],product[N];

void Fast_Fourier_Transform(complex x[],int len,double mode)
{
    for(int i=0;i<len;i++)
    if(i<pos[i]) swap(x[i],x[pos[i]]);

    for(int i=2;i<=len;i<<=1)
    {
        int step=i>>1;
        complex wm(cos(2.0*pi/i),sin(mode*pi/i));
        for(int j=0;j<len;j+=i)
        {
            int lim=j+step;
            complex w(1,0);
            for(int k=j;k<lim;k++)
            {
                complex a=x[k],b=w*x[k+step];
                x[k]=a+b;x[k+step]=a-b;
                w=w*wm;
            }
        }
    }

    if(mode==IDFT)
    for(int i=0;i<len;i++)
    x[i].a/=len;
    return;
}
#define FFT Fast_Fourier_Transform

int n,len;
int ans[N];
int main()
{
    scanf("%d",&n);
    len=1;
    while(len<(n<<1)) len<<=1;
    initialization(len);

    for(int i=0;i<n;i++) scanf("%lf%lf",&factor1[n-i-1].a,&factor2[i]);
    FFT(factor1,len,DFT);
    FFT(factor2,len,DFT);

    for(int i=0;i<len;i++) product[i]=factor1[i]*factor2[i];
    FFT(product,len,IDFT);
    for(int i=0;i<n;i++) ans[i]=floor(product[i].a+0.5);
    for(int i=n-1;~i;i--) printf("%d\n",ans[i]);

    return 0;
}