FFT(快速傅立叶变换)：HDU 1402 A * B Problem Plus

Calculate A * B.

Input

Each line will contain two integers A and B. Process to end of file.

Note: the length of each integer will not exceed 50000.

Output

For each case, output A * B in one line.

Sample Input

Sample Output

2
2000
　　
　　唉，模板题，膜的邝斌的模板。

  1 #include <iostream>
  2 #include <cstdio>
  3 #include <cstring>
  4 #include <cmath>
  5 #include <algorithm>
  6 using namespace std;
  7 
  8 const double PI = acos(-1.0);
  9 
 10 struct complex{
 11     double r,i;
 12     complex(double r_=0.0,double i_=0.0)
 13     {
 14         r=r_;i=i_;
 15     }
 16     complex operator +(const complex &b)
 17     {
 18         return complex(r+b.r,i+b.i);
 19     }
 20     complex operator -(const complex &b)
 21     {
 22         return complex(r-b.r,i-b.i);
 23     }
 24     complex operator *(const complex &b)
 25     {
 26         return complex(r*b.r-i*b.i,i*b.r+r*b.i);
 27     }
 28 };
 29  
 30 void Rader(complex *a,int len)
 31 {
 32     int k;
 33     for(int i=1,j=len/2;i<len-1;i++)
 34     {
 35         if(i<j)swap(a[i],a[j]);
 36         k=len/2;
 37         while(j>=k)
 38         {
 39             j-=k;
 40             k>>=1;
 41         }
 42         j+=k;
 43     }
 44 }
 45 
 46 void FFT(complex *a,int len,int on)
 47 {
 48     Rader(a,len);
 49     for(int h=2;h<=len;h<<=1)
 50     {
 51         complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
 52         for(int j=0;j<len;j+=h)
 53         {
 54             complex w(1,0);
 55             for(int k=j;k<j+h/2;k++)
 56             {
 57                 complex u=a[k];
 58                 complex v=a[k+h/2]*w;
 59                 a[k]=u+v;
 60                 a[k+h/2]=u-v;
 61                 w=w*wn;
 62             }
 63         } 
 64     }
 65     if(on==-1)
 66         for(int i=0;i<len;i++)
 67             a[i].r/=len;
 68 }
 69 
 70 const int maxn = 200010;
 71 complex Array1[maxn],Array2[maxn];
 72 char str1[maxn],str2[maxn];
 73 int sum[maxn],len,len1,len2;
 74 
 75 int main()
 76 {
 77     while(~scanf("%s%s",str1,str2))
 78     {
 79         len1=strlen(str1);
 80         len2=strlen(str2);    
 81         len=1;
 82         while(len<len1*2||len<len2*2)len<<=1;
 83         for(int i=0;i<len1;i++)
 84             Array1[i]=complex(str1[len1-i-1]-'0',0);
 85         for(int i=0;i<len2;i++)
 86             Array2[i]=complex(str2[len2-i-1]-'0',0);
 87         
 88         for(int i=len1;i<len;i++)
 89             Array1[i]=complex(0,0);
 90         for(int i=len2;i<len;i++)
 91             Array2[i]=complex(0,0);
 92         
 93         FFT(Array1,len,1);
 94         FFT(Array2,len,1);
 95         for(int i=0;i<len;i++)
 96             Array1[i]=Array1[i]*Array2[i];
 97         FFT(Array1,len,-1);
 98         memset(sum,0,sizeof(sum));    
 99         for(int i=0;i<len;i++){
100             sum[i]+=(int)(Array1[i].r+0.5);
101             sum[i+1]+=sum[i]/10;
102             sum[i]%=10;
103         }
104         int p=len;
105         while(!sum[p]&&p)p--;
106         for(;p!=-1;p--)
107             printf("%d",sum[p]);
108         printf("
");                    
109     }
110     return 0;
111 }

　　两个月后重打了一遍。

 1 #include <iostream>
 2 #include <cstring>
 3 #include <cstdio>
 4 #include <cmath>
 5 using namespace std;
 6 const int maxn=400010;
 7 const double PI=acos(-1.0);
 8 struct complex{
 9     double r,i;
10     complex(double r_=0.0,double i_=0.0){
11         r=r_;i=i_;
12     }
13     complex operator +(complex a){
14         return complex(r+a.r,i+a.i);
15     }
16     complex operator -(complex a){
17         return complex(r-a.r,i-a.i);
18     }
19     complex operator *(complex a){
20         return complex(r*a.r-i*a.i,r*a.i+a.r*i);
21     }
22 }A[maxn],B[maxn];
23 
24 void Rader(complex *a,int len){
25     int k;
26     for(int i=1,j=len>>1;i<len-1;i++){
27         if(i<j)swap(a[i],a[j]);
28         k=len>>1;
29         while(j>=k){
30             j-=k;
31             k>>=1;
32         }
33         j+=k;
34     }
35 }
36 
37 void FFT(complex *a,int len,int on){
38     Rader(a,len);
39     for(int h=2;h<=len;h<<=1){
40         complex wn(cos(-on*PI*2.0/h),sin(-on*PI*2.0/h));
41         for(int j=0;j<len;j+=h){
42             complex w(1,0);
43             for(int k=j;k<j+(h>>1);k++){
44                 complex x=a[k];
45                 complex y=a[k+(h>>1)]*w;
46                 a[k]=x+y;
47                 a[k+(h>>1)]=x-y;
48                 w=w*wn;
49             }    
50         }
51     }
52     if(on==-1)
53         for(int i=0;i<len;i++)
54             a[i].r/=len;
55     return;    
56 }
57 
58 char s[maxn],t[maxn];
59 int ans[maxn];
60 int main(){
61     while(~scanf("%s%s",s,t)){
62     int lens=strlen(s);
63     int lent=strlen(t),len=1;
64     while(len<=lens+lent)len<<=1;
65     memset(A,0,sizeof(A));
66     memset(B,0,sizeof(B));
67     for(int i=0;i<lens;i++)
68         A[i].r=1.0*(s[lens-i-1]-'0');
69     for(int i=0;i<lent;i++)
70         B[i].r=1.0*(t[lent-i-1]-'0');
71     FFT(A,len,1);FFT(B,len,1);
72     for(int i=0;i<len;i++)
73         A[i]=A[i]*B[i];
74     FFT(A,len,-1);
75 
76     memset(ans,0,sizeof(ans));
77     int in=0;
78     for(int i=0;i<len;i++){
79         ans[i]=in+floor(A[i].r+0.5);
80         in=ans[i]/10;
81         ans[i]%=10;
82     }
83     int i;
84     for(i=len;i&&!ans[i];i--);
85     while(~i)printf("%d",ans[i]),i--;
86     printf("
");
87     }
88     return 0;
89 }

尽最大的努力，做最好的自己！