2194: 快速傅立叶之二

最新推荐文章于 2019-12-16 21:16:07 发布

转载最新推荐文章于 2019-12-16 21:16:07 发布 · 115 阅读

0 ·

CC 4.0 BY-SA版权

原文链接：http://www.cnblogs.com/htwx/articles/5585480.html

本文深入解析了快速傅立叶变换（FFT）算法在计算特定卷积问题中的应用，通过一个具体的编程实例展示了如何使用FFT算法来解决复杂的数学计算问题，特别是针对大规模数据集的高效处理。

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

卷积是要求两个下标的和是定值，这道题的差是定值，这样就可以将其中一个数组翻转，然后fft。。。

 1 #include<iostream>
 2 #include<cstdlib>
 3 #include<cmath>
 4 #include<cstring>
 5 #include<cstdio>
 6 #include<algorithm>
 7 #include<string>
 8 #include<map>
 9 #include<queue>
10 #include<vector>
11 #include<set>
12 #define inf 1000000000
13 #define maxn 266333
14 #define maxm 10000+5
15 #define eps 1e-10
16 #define ll long long
17 #define for0(i,n) for(int i=0;i<=(n);i++)
18 #define for1(i,n) for(int i=1;i<=(n);i++)
19 #define for2(i,x,y) for(int i=(x);i<=(y);i++)
20 #define for3(i,x,y) for(int i=(x);i>=(y);i--)
21 #define for4(i,x) for(int i=head[x],y=e[i].go;i;i=e[i].next,y=e[i].go)
22 using namespace std;
23 int read(){
24     int x=0,f=1;char ch=getchar();
25     while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
26     while(ch>='0'&&ch<='9'){x=10*x+ch-'0';ch=getchar();}
27     return x*f;
28 }
29 struct cp{
30     double x,y;
31     cp operator +(cp b){return (cp){x+b.x,y+b.y};}
32     cp operator -(cp b){return (cp){x-b.x,y-b.y};}
33     cp operator *(cp b){return (cp){x*b.x-y*b.y,x*b.y+y*b.x};}
34 };
35 const double PI=acos(-1.0);
36 cp a[maxn],b[maxn],c[maxn],y[maxn];
37 int n,m,len,rev[maxn],ans[maxn];
38 char s[maxn];
39 void fft(cp *x,int n,int flag){
40     for0(i,n-1)y[rev[i]]=x[i];
41     for0(i,n-1)x[i]=y[i];
42     for(int m=2;m<=n;m<<=1){
43         cp wn=(cp){cos(2.0*PI/m*flag),sin(2.0*PI/m*flag)};
44         for(int i=0;i<n;i+=m){
45             cp w=(cp){1.0,0};int mid=m>>1;
46             for0(j,mid-1){
47                 cp u=x[i+j],v=x[i+j+mid]*w;
48                 x[i+j]=u+v;x[i+j+mid]=u-v;
49                 w=w*wn;
50             }
51         }
52     }
53     if(flag==-1)for0(i,n-1)x[i].x/=n;
54 }
55 int main(){
56     //freopen("input.txt","r",stdin);
57     //freopen("output.txt","w",stdout);
58     n=read();
59     for0(i,n-1)a[n-1-i].x=read(),b[i].x=read();
60     n=2*n-1;m=1;
61     while(m<n)m<<=1,len++;
62     swap(n,m);
63     for0(i,n-1){
64         int x=i,y=0;
65         for1(j,len)y<<=1,y|=(x&1),x>>=1;
66         rev[i]=y;
67     }
68     fft(a,n,1);fft(b,n,1);
69     for0(i,n-1)c[i]=a[i]*b[i];
70     fft(c,n,-1);
71     for3(i,(m-1)>>1,0)
72         printf("%lld\n",(ll)(c[i].x+0.5));
73     return 0;
74 }

View Code

转载于:https://www.cnblogs.com/htwx/articles/5585480.html