洛谷P1919 【模板】A*B Problem升级版（FFT）

最新推荐文章于 2021-09-11 09:19:23 发布

weixin_34130269

最新推荐文章于 2021-09-11 09:19:23 发布

阅读量104

点赞数

CC 4.0 BY-SA版权

原文链接：http://www.cnblogs.com/bztMinamoto/p/9742468.html

本文介绍如何使用快速傅立叶变换(FFT)算法进行大数乘法，通过将数字转换为多项式，利用FFT加速乘法过程，再通过逆FFT得到结果。文章包含详细的C++代码实现。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

传送门

话说FFT该不会真的只能用来做这种板子吧……

我们把两个数字的每一位都看作多项式的系数

然后这就是一个多项式乘法

上FFT就好了

然后去掉前导零

（然而连FFT的板子都背不来orz，而且空间又开小了……）

 1 //minamoto
 2 #include<iostream>
 3 #include<cstdio>
 4 #include<cmath>
 5 using namespace std;
 6 char sr[1<<21],z[20];int C=-1,Z;
 7 inline void Ot(){fwrite(sr,1,C+1,stdout),C=-1;}
 8 inline void print(int x){
 9     if(C>1<<20)Ot();if(x<0)sr[++C]=45,x=-x;
10     while(z[++Z]=x%10+48,x/=10);
11     while(sr[++C]=z[Z],--Z);
12 }
13 const int N=2e5+5;const double Pi=acos(-1.0);
14 int r[N],l=0,limit=1,c[N],n;char sa[N],sb[N];
15 struct complex{
16     double x,y;
17     complex(double xx=0,double yy=0){x=xx,y=yy;}
18     inline complex operator +(complex b){return complex(x+b.x,y+b.y);}
19     inline complex operator -(complex b){return complex(x-b.x,y-b.y);}
20     inline complex operator *(complex b){return complex(x*b.x-y*b.y,x*b.y+y*b.x);}
21 }a[N],b[N];
22 void FFT(complex *a,int type){
23     for(int i=0;i<limit;++i)
24     if(i<r[i]) swap(a[i],a[r[i]]);
25     for(int mid=1;mid<limit;mid<<=1){
26         complex Wn(cos(Pi/mid),type*sin(Pi/mid));
27         for(int R=mid<<1,j=0;j<limit;j+=R){
28             complex w(1,0);
29             for(int k=0;k<mid;++k,w=w*Wn){
30                 complex x=a[j+k],y=w*a[j+k+mid];
31                 a[j+k]=x+y,a[j+k+mid]=x-y;
32             }
33         }
34     }
35 }
36 int main(){
37 //    freopen("testdata.in","r",stdin);
38     scanf("%d",&n),--n;
39     scanf("%s%s",sa,sb);
40     for(int i=0;i<=n;++i) a[i].x=sa[n-i]-'0',b[i].x=sb[n-i]-'0';
41     while(limit<=n*2) limit<<=1,++l;
42     for(int i=0;i<=limit;++i) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
43     FFT(a,1),FFT(b,1);
44     for(int i=0;i<=limit;++i) a[i]=a[i]*b[i];
45     FFT(a,-1);
46     for(int i=0;i<=limit;++i) c[i]=(int)(a[i].x/limit+0.5);
47     for(int i=0;i<=limit;++i)
48     if(c[i]>9){
49         c[i+1]+=c[i]/10,c[i]%=10;
50         if(i+1>limit) ++limit;
51     }
52     for(int i=limit;i>=0;--i)
53     if(c[i]==0) --limit;
54     else break;
55     for(int i=limit;i>=0;--i) print(c[i]);
56     Ot();
57     return 0;
58 }