FFT
Posted __AiR_H
tags:
篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了FFT相关的知识,希望对你有一定的参考价值。
1 #include <iostream> 2 #include <cstdio> 3 #include <cstring> 4 #include <cstdlib> 5 #include <algorithm> 6 #include <queue> 7 #include <vector> 8 #include <stack> 9 #include <map> 10 #include <set> 11 #include <cmath> 12 #include <cctype> 13 #include <ctime> 14 #include <ext/rope> 15 16 using namespace std; 17 using namespace __gnu_cxx; 18 19 #define REP(i, n) for (int i = 0; i < (n); ++i) 20 #define eps 1e-9 21 #define SZ(x) ((int)x.size()) 22 23 typedef long long ll; 24 typedef pair<int, int> pii; 25 26 #define PI acos(-1.0) 27 const int maxn = 1e5 + 10; 28 struct Complex { 29 double x, y; 30 Complex(double _x = 0.0, double _y = 0.0) { x = _x; y = _y; } 31 Complex operator - (const Complex &b) const { return Complex(x - b.x, y - b.y); } 32 Complex operator + (const Complex &b) const { return Complex(x + b.x, y + b.y); } 33 Complex operator * (const Complex &b) const { return Complex(x * b.x - y * b.y, x * b.y + y * b.x); } 34 }; 35 char str_A[maxn], str_B[maxn]; 36 Complex A[maxn * 4], B[maxn * 4], ans[maxn * 4]; 37 int ans_final[maxn * 4], len_final = 0; 38 int N, M, n = 1; 39 40 void fft(Complex a[], int len, int oper) { 41 for (int i = 1, j = 0; i < len - 1; ++i) { 42 for (int s = len; ~j & s;) { j ^= s >>= 1; } 43 if (i < j) { swap(a[i], a[j]); } 44 } 45 Complex omega_m, omega, t, u; 46 int m, m2; 47 for (int d = 0; (1 << d) < len; ++d) { 48 m = (1 << d); m2 = m << 1; 49 omega_m = Complex(cos(PI / m * oper), sin(PI / m * oper)); 50 for (int i = 0; i < len; i += m2) { 51 omega = Complex(1, 0); 52 for (int j = 0; j < m; ++j) { 53 u = a[i + j]; t = omega * a[i + j + m]; 54 a[i + j] = u + t; a[i + j + m] = u - t; 55 omega = omega * omega_m; 56 } 57 } 58 } 59 if (oper == -1) { REP(i, N + M - 1) { a[i].x /= len; } } 60 } 61 62 63 int main() { 64 #ifdef __AiR_H 65 freopen("in.txt", "r", stdin); 66 // freopen("out.txt", "w", stdout); 67 #endif // __AiR_H 68 while (scanf("%s %s", str_A, str_B) != EOF) { 69 memset(ans, 0, sizeof(ans)); memset(ans_final, 0, sizeof(ans_final)); 70 N = strlen(str_A); M = strlen(str_B); n = 1; 71 while (n < N + M - 1) { n <<= 1; } 72 REP(i, N) { A[i] = Complex(str_A[N - 1 - i] - ‘0‘, 0); } 73 for (int i = N; i < n; ++i) { A[i] = Complex(0, 0); } 74 REP(i, M) { B[i] = Complex(str_B[M - 1 - i] - ‘0‘, 0); } 75 for (int i = M; i < n; ++i) { B[i] = Complex(0, 0); } 76 fft(A, n, 1); fft(B, n, 1); 77 REP(i, n) { ans[i] = A[i] * B[i]; } 78 fft(ans, n, -1); 79 REP(i, N + M) { ans_final[i] = ans[i].x + 0.5; } 80 REP(i, N + M - 1) { ans_final[i + 1] += ans_final[i] / 10; ans_final[i] %= 10; } 81 len_final = N + M - 1; 82 if (ans_final[len_final]) { ++len_final; } 83 while (ans_final[len_final - 1] == 0 && len_final > 1) { --len_final; } 84 for (int i = len_final - 1; i > 0; --i) { printf("%d", ans_final[i]); } printf("%d\n", ans_final[0]); 85 } 86 return 0; 87 }
以上是关于FFT的主要内容,如果未能解决你的问题,请参考以下文章