FFT其实没什么需要特别了解的,了解下原理,(特别推荐算法导论上面的讲解),模板理解就行了。重在运用吧。
处理过程中要特别注意精度。
先上个练习的地址吧:
http://vjudge.net/vjudge/contest/view.action?cid=53596#overview
Problem A A * B Problem Plus
A*B的大数乘法,似乎大数模拟乘法不行的,得用FFT优化到nlogn,将一个数AnAn-1....A1A0,看做An*10^n+An-1*10^n-1+....A1*10+A0*10^0,这样就可以将两个数相乘当成多项式乘法了。
代码:
1 #include <cstdio> 2 #include <iostream> 3 #include <cstring> 4 #include <string> 5 #include <cstdlib> 6 #include <algorithm> 7 #include <vector> 8 #include <cmath> 9 #include <queue> 10 #include <map> 11 #include <set> 12 #include <complex> 13 #define pb push_back 14 #define in freopen("solve_in.txt", "r", stdin); 15 #define out freopen("solve_out.txt", "w", stdout); 16 #define pi (acos(-1.0)) 17 #define bug(x) printf("line %d :>>>>>>>>>>>>>>>\n", x); 18 #define pb push_back 19 using namespace std; 20 #define esp 1e-8 21 22 typedef complex<double> CD; 23 const int maxn = 50000+100; 24 char s1[maxn], s2[maxn]; 25 26 struct Complex{ 27 double x, y; 28 Complex(){} 29 Complex(double x, double y):x(x), y(y){} 30 }; 31 Complex operator + (const Complex &a, const Complex &b){ 32 Complex c; 33 c.x = a.x+b.x; 34 c.y = a.y+b.y; 35 return c; 36 } 37 Complex operator - (const Complex &a, const Complex &b){ 38 Complex c; 39 c.x = a.x-b.x; 40 c.y = a.y-b.y; 41 return c; 42 } 43 Complex operator * (const Complex &a, const Complex &b){ 44 Complex c; 45 c.x = a.x*b.x-a.y*b.y; 46 c.y = a.x*b.y+a.y*b.x; 47 return c; 48 } 49 50 inline void FFT(vector<Complex> &a, bool inverse) 51 { 52 int n = a.size(); 53 for(int i = 0, j = 0; i < n; i++) 54 { 55 if(j > i) 56 swap(a[i], a[j]); 57 int k = n; 58 while(j & (k>>=1)) j &= ~k; 59 j |= k; 60 } 61 double PI = inverse ? -pi : pi; 62 for(int step = 2; step <= n; step <<= 1) 63 { 64 double alpha = 2*PI/step; 65 Complex wn(cos(alpha), sin(alpha)); 66 for(int k = 0; k < n; k += step) 67 { 68 Complex w(1, 0); 69 for(int Ek = k; Ek < k+step/2; Ek++) 70 { 71 int Ok = Ek + step/2; 72 Complex u = a[Ek]; 73 Complex t = a[Ok]*w; 74 a[Ok] = u-t; 75 a[Ek] = u+t; 76 w = w*wn; 77 } 78 } 79 } 80 if(inverse) 81 for(int i = 0; i < n; i++) 82 a[i].x = (a[i].x/n); 83 } 84 vector<double> operator * (const vector<double> &v1, const vector<double> &v2) 85 { 86 int S1 = v1.size(), S2 = v2.size(); 87 int S = 2; 88 while(S < S1+S2) S <<= 1; 89 vector<Complex> a(S), b(S); 90 for(int i = 0; i < S; i++) 91 a[i].x = a[i].y = b[i].x = b[i].y = 0.0; 92 for(int i = 0; i < S1; i++) 93 a[i].x = v1[i]; 94 for(int i = 0; i < S2; i++) 95 b[i].x = v2[i]; 96 FFT(a, false); 97 FFT(b, false); 98 for(int i = 0; i < S; i++) 99 a[i] = a[i] * b[i]; 100 FFT(a, true); 101 vector<double> res(S1+S2-1, 0.0); 102 for(int i = 0; i < S1+S2-1; i++) 103 res[i] = a[i].x; 104 return res; 105 } 106 int ans[maxn*2+100]; 107 vector<double > v1, v2; 108 int main() 109 { 110 111 while(scanf("%s%s", s1, s2) != -1) 112 { 113 114 v1.clear(); 115 v2.clear(); 116 int len1 = strlen(s1); 117 int len2 = strlen(s2); 118 for(int i = 0; s1[i]; i++) 119 v1.pb((double)(s1[len1-1-i]-'0')); 120 for(int i = 0; s2[i]; i++) 121 v2.pb((double)(s2[len2-1-i]-'0')); 122 v1 = v1*v2; 123 memset(ans, 0, sizeof(ans)); 124 int carry = 0, top = 0; 125 for(int i = 0; i < v1.size(); i++){ 126 carry += (int)(v1[i]+0.5); 127 ans[top++] = carry%10; 128 carry /= 10; 129 } 130 while(carry){ 131 ans[top++] = carry%10; 132 carry /= 10; 133 } 134 while(top) 135 { 136 if(ans[top]) 137 break; 138 top--; 139 } 140 for(int i = top; i >= 0; i--) 141 printf("%d", ans[i]); 142 cout<<endl; 143 } 144 return 0; 145 }