首页 > 代码库 > [Zjoi2014]力 FFT

[Zjoi2014]力 FFT

#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
using namespace std;
typedef double dd;
const dd pi=acos(0.0)*2;
#define N 400005
struct P{
    dd x,y;
    P(dd A=0.0,dd B=0.0){
        x=A;
        y=B;
    }
    P operator+(P a){
        return P(x+a.x,y+a.y);
    }
    P operator-(P a){
        return P(x-a.x,y-a.y);
    }
    P operator*(P a){
        return P(x*a.x-y*a.y,x*a.y+y*a.x);
    }
}a[N],b[N],w[2][N];
int n,na,nb,rev[N];
void FFT(P *a,int f){
    int i,j,k,t,l;
    P x,y;
    for(i=0;i<n;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
    for(i=1;i<n;i<<=1)
        for(j=0,t=n/(i<<1);j<n;j+=i<<1)
        for(k=l=0;k<i;k++,l+=t){
            x=w[f][l]*a[j+k+i];    
            y=a[j+k];
            a[j+k]=y+x;
            a[j+k+i]=y-x;
        }
        if(f)for(i=0;i<n;i++)a[i].x/=n;
}
void prepare(){
    int i,j,x,y;
    for(i=0;i<n;i++){
        x=i,y=0;
        for(j=1;j<n;x>>=1,j<<=1)(y<<=1)|=x&1;
        rev[i]=y;
    }
    for(i=0;i<n;i++){
        w[0][i]=P(cos(2*pi*i/n),sin(2*pi*i/n));
        w[1][i]=P(cos(2*pi*i/n),-sin(2*pi*i/n));
    }
}
void mult(){
    int i;
    prepare();
    FFT(a,0);
    FFT(b,0);
    for(i=0;i<n;i++)a[i]=a[i]*b[i];
    FFT(a,1);    
}
dd x[N],ans[N];
int nn;
int main(){
    scanf("%d",&nn);
    int i;
    for(i=0;i<nn;i++)scanf("%lf",&x[i]);
    for(n=1;n<nn;n<<=1);
    n<<=1;
    for(i=1;i<nn;i++)b[i]=P(1.0/i/i,0);
    for(i=0;i<nn;i++)a[i]=P(x[i],0);
    mult();
    for(i=0;i<nn;i++)ans[i]+=a[i].x;
    for(i=0;i<n;i++)a[i]=b[i]=P(0,0);
    for(i=1;i<nn;i++)b[i]=P(1.0/i/i,0);
    for(i=0;i<nn;i++)a[i]=P(x[nn-i-1],0);
    mult();
    for(i=0;i<nn;i++)ans[i]-=a[nn-i-1].x;
    for(i=0;i<nn;i++)printf("%.3lf\n",ans[i]);
}