bzoj2194 快速傅立叶之二 FFT

2194: 快速傅立叶之二

Time Limit: 10 Sec  Memory Limit: 259 MB
Submit: 1697  Solved: 1004
[Submit][Status][Discuss]

Description

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

 

Input

第一行一个整数N,接下来N行，第i+2..i+N-1行，每行两个数，依次表示a[i],b[i] (0 < = i < N)。

Output

输出N行，每行一个整数，第i行输出C[i-1]。

Sample Input

5
3 1
2 4
1 1
2 4
1 4

Sample Output

24
12
10
6
1

HINT

卷积的形式

这样的形式化出来才可以fft

我们把b数组翻转，i-->-i
c[i+j]+=a[i]*b[j]

 

#pragma GCC optimize(2)
#pragma G++ optimize(2)
#include<cstdio>
#include<cmath>
#include<iostream>
#include<algorithm>
#include<cstring>

#define N 300007
#define pi acos(-1)
using namespace std;
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}
    while(isdigit(ch)){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
    return x*f;
}

int n,m,L;
int rev[N];
struct comp
{
    double r,v;
    inline comp operator+(const comp &a){return (comp){r+a.r,v+a.v};}
    inline comp operator-(const comp &a){return (comp){r-a.r,v-a.v};}
    inline comp operator*(const comp &a){return (comp){r*a.r-v*a.v,r*a.v+v*a.r};}
}a[N],b[N];

void FFT(comp *a,int f)
{
    for (int i=0;i<n;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    for (int i=1;i<n;i<<=1)
    {
        comp wn=(comp){cos(pi/i),f*sin(pi/i)};
        for (int j=0;j<n;j+=(i<<1))
        {
            comp w=(comp){1,0};
            for (int k=0;k<i;k++,w=w*wn)
            {
                comp x=a[j+k],y=w*a[j+k+i];
                a[j+k]=x+y,a[j+k+i]=x-y;
            }
        }
    }
    if(f==-1)for (int i=0;i<n;i++)a[i].r/=n;
}
int main()
{
    n=read()-1,m=n*2;
    for(int i=0;i<=n;i++)
        a[i].r=(double)read(),b[n-i].r=(double)read();
    for(n=1;n<=m;n<<=1,L++);
    if(L)L--;
    for(int i=0;i<n;i++)
        rev[i]=(rev[i>>1]>>1)|((i&1)<<L);
    FFT(a,1),FFT(b,1);
    for (int i=0;i<n;i++)a[i]=a[i]*b[i];
    FFT(a,-1);
    for (int i=m/2;i<=m;i++)
        printf("%d\n",(int)(a[i].r+0.5));
}