FFT

递归版

#include<bits/stdc++.h>
using namespace std;
#define N 4000005
int n,m;
const double Pi=acos(-1);
int tr[N];
struct CP{
    double x,y;
}f[N],g[N];
CP operator + (CP a,CP b){return (CP){a.x+b.x,a.y+b.y};}
CP operator - (CP a,CP b){return (CP){a.x-b.x,a.y-b.y};}
CP operator * (CP a,CP b){return (CP){a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};}
void FFT(CP *f,int len,bool flag){
    if(len==1) return ;
    FFT(f,len/2,flag);
    FFT(f+len/2,len/2,flag);
    CP tG=(CP){cos(2*Pi/len),sin(2*Pi/len)},buf=(CP){1,0};
    if(!flag) tG.y*=-1;
    for(int k=0;k<len/2;k++){
        CP tt=buf*f[k+len/2];
        f[k+len/2]=f[k]-tt;
        f[k]=f[k]+tt;
        buf=buf*tG;
    } 
}
int main(){
    scanf("%d%d",&n,&m);
    for(int i=0;i<N-3;i++) f[i].x=f[i].y=g[i].x=g[i].y=0;
    for(int i=0;i<=n;i++) scanf("%lf",&f[i].x);
    for(int i=0;i<=m;i++) scanf("%lf",&g[i].x);
    for(m+=n,n=1;n<=m;n<<=1);
    for(int i=0;i<n;i++) tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);
    for(int i=0;i<n;i++){
        if(i<tr[i]) swap(f[i],f[tr[i]]),swap(g[i],g[tr[i]]);
    }
    FFT(f,n,1),FFT(g,n,1);
    for(int i=0;i<n;i++) f[i]=f[i]*g[i];
    for(int i=0;i<n;i++){
        if(i<tr[i]) swap(f[i],f[tr[i]]);
    }
    FFT(f,n,0);
    for(int i=0;i<=m;i++) printf("%d ",(int)(f[i].x/n+0.49));
    return 0;
}

非递归

#include<bits/stdc++.h>
using namespace std;
#define N 4000005
int n,m;
const double Pi=acos(-1);
int tr[N];
struct CP{
    double x,y;
}f[N],g[N];
CP operator + (CP a,CP b){return (CP){a.x+b.x,a.y+b.y};}
CP operator - (CP a,CP b){return (CP){a.x-b.x,a.y-b.y};}
CP operator * (CP a,CP b){return (CP){a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};}
void FFT(CP *f,bool flag){
    for(int len=2;len<=n;len<<=1){
        CP tG=(CP){cos(2*Pi/len),sin(2*Pi/len)};
        if(!flag) tG.y*=-1;
        for(int qi=0;qi<n;qi+=len){
            CP buf=(CP){1,0};
            for(int k=qi;k<qi+len/2;k++){
                CP tt=buf*f[k+len/2];
                f[k+len/2]=f[k]-tt;
                f[k]=f[k]+tt;
                buf=buf*tG;             
            }
        }       
    }

}
int main(){
    scanf("%d%d",&n,&m);
    for(int i=0;i<N-3;i++) f[i].x=f[i].y=g[i].x=g[i].y=0;
    for(int i=0;i<=n;i++) scanf("%lf",&f[i].x);
    for(int i=0;i<=m;i++) scanf("%lf",&g[i].x);
    for(m+=n,n=1;n<=m;n<<=1);
    for(int i=0;i<n;i++) tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);
    for(int i=0;i<n;i++){
        if(i<tr[i]) swap(f[i],f[tr[i]]),swap(g[i],g[tr[i]]);
    }
    FFT(f,1),FFT(g,1);
    for(int i=0;i<n;i++) f[i]=f[i]*g[i];
    for(int i=0;i<n;i++){
        if(i<tr[i]) swap(f[i],f[tr[i]]);
    }
    FFT(f,0);
    for(int i=0;i<=m;i++) printf("%d ",(int)(f[i].x/n+0.49));
    return 0;
}
posted @ 2023-09-21 21:50  hubingshan  阅读(54)  评论(0)    收藏  举报  来源