B. 2194: 快速傅立叶之二解题报告

$$\begin{eqnarray}&c[k] = \sum_{i}^{n}a[i]b[i-k] \\&c[k] = \sum_{i}^{n}a[n-i]b[i-k] (倒序保存a) \\&c[n-k]= \sum_{i}^{n}a[n-i]b[i-k] (倒序保存c) \\&通过卷积 o (nlog(n))得到c\end{eqnarray}$$

#include<bits/stdc++.h>
using namespace std;
const int N=135000;
const double Pi=acos(-1.0);
int n,m;
inline int read()
{
    int s=0,w=1;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-')w=-1;ch=getchar();}
    while(isdigit(ch)){s=s*10+ch-'0';ch=getchar();}
    return s*w;
}
inline int _min(int a,int b){return a<b?a:b;}
inline int _max(int a,int b){return a>b?a:b;}
struct CP
{
    double x,y;
    CP (double xx=0,double yy=0)
    {
        x=xx;y=yy;
    }
    CP operator + (const CP &B)    const
    {
        return CP(x+B.x,y+B.y);
    }
    CP operator - (const CP &B) const
    {
        return CP(x-B.x,y-B.y);
    }
    CP operator * (const CP &B) const
    {
        return CP(x*B.x-y*B.y,x*B.y+y*B.x);
    }
}f[N<<1];
int tr[N<<1];
void FFT(CP *f,bool flag)
{
    for(int i=0;i<n;i++)
    {
        if(i<tr[i])
            swap(f[i],f[tr[i]]);
    }
    for(int p=2;p<=n;p<<=1)
    {
        int len=p>>1;
        CP tG(cos(2*Pi/p),sin(2*Pi/p));
        if(flag==0)
        {
            tG.y*=-1;
        }
        for(int k=0;k<n;k+=p)
        {
            CP buf(1,0);
            for(int l=k;l<k+len;l++)
            {
                CP tt=buf*f[len+l];
                f[len+l]=f[l]-tt;
                f[l]=f[l]+tt;
                buf=buf*tG;
            }
        }
    }
}
int main()
{
    n=read();m=n;
    for(int i=0;i<n;i++)
    {
        f[n-i].x=read();
        f[i].y=read();
    }
    int tmp=n;
    for(n=1;n<=tmp*2;n<<=1);
    for(int i=0;i<n;i++)
    {
        tr[i]=(tr[i>>1]>>1)|(i&1?(n>>1):0);
    }
    FFT(f,1);
    for(int i=0;i<n;i++)
    {
        f[i]=f[i]*f[i];
    }
    FFT(f,0);
    //存的是c[n-k],0<=k<=n-1,1<=n-k<=n,所以输出1~n 
    for(int i=m;i>=1;i--)
    {
        printf("%d\n",(int)(f[i].y/n/2+0.49));
    }
    return 0;
}
原文地址:https://www.cnblogs.com/HKHbest/p/14448001.html