摘自：https://www.cnblogs.com/RabbitHu/p/FFT.html
快速傅里叶变换(FFT)： 一种能在O(nlogn)的时间内将一个多项式转换成它的点值表示的算法。
点值表示：设A(x)是一个n−1次多项式，那么把n个不同的x代入，会得到n个y。这n对(x,y)唯一确定了该多项。由多项式可以求出其点值表示，而由点值表示也可以求出多项式。
设有两个n−1次多项式A(x)和B(x))，我们的目标是——把它们乘起来。普通的多项式乘法是O(n^2)，但有趣的是，两个用点值表示的多项式相乘，复杂度是O(n)的！具体方法：C(xi)=A(xi)×B(xi)。
高精度乘法：等价于x=10的A(x)*B(x)
1 #include <cstdio> 2 #include <cmath> 3 #include <cstring> 4 #include <algorithm> 5 #include <complex> 6 #define space putchar(' ') 7 #define enter putchar(' ') 8 using namespace std; 9 typedef long long ll; 10 template <class T> 11 void read(T &x){ 12 char c; 13 bool op = 0; 14 while(c = getchar(), c < '0' || c > '9') 15 if(c == '-') op = 1; 16 x = c - '0'; 17 while(c = getchar(), c >= '0' && c <= '9') 18 x = x * 10 + c - '0'; 19 if(op) x = -x; 20 } 21 template <class T> 22 void write(T x){ 23 if(x < 0) putchar('-'), x = -x; 24 if(x >= 10) write(x / 10); 25 putchar('0' + x % 10); 26 } 27 const int N = 1000005; 28 const double PI = acos(-1); 29 typedef complex <double> cp; 30 char sa[N], sb[N]; 31 int n = 1, lena, lenb, res[N]; 32 cp a[N], b[N], omg[N], inv[N]; 33 void init(){ 34 for(int i = 0; i < n; i++){ 35 omg[i] = cp(cos(2 * PI * i / n), sin(2 * PI * i / n)); 36 inv[i] = conj(omg[i]); 37 } 38 } 39 void fft(cp *a, cp *omg){ 40 int lim = 0; 41 while((1 << lim) < n) lim++; 42 for(int i = 0; i < n; i++){ 43 int t = 0; 44 for(int j = 0; j < lim; j++) 45 if((i >> j) & 1) t |= (1 << (lim - j - 1)); 46 if(i < t) swap(a[i], a[t]); // i < t 的限制使得每对点只被交换一次（否则交换两次相当于没交换） 47 } 48 for(int l = 2; l <= n; l *= 2){ 49 int m = l / 2; 50 for(cp *p = a; p != a + n; p += l) 51 for(int i = 0; i < m; i++){ 52 cp t = omg[n / l * i] * p[i + m]; 53 p[i + m] = p[i] - t; 54 p[i] += t; 55 } 56 } 57 } 58 int main(){ 59 scanf("%s%s", sa, sb); 60 lena = strlen(sa), lenb = strlen(sb); 61 while(n < lena + lenb) n *= 2; 62 for(int i = 0; i < lena; i++) 63 a[i].real(sa[lena - 1 - i] - '0'); 64 for(int i = 0; i < lenb; i++) 65 b[i].real(sb[lenb - 1 - i] - '0'); 66 init(); 67 fft(a, omg); 68 fft(b, omg); 69 for(int i = 0; i < n; i++) 70 a[i] *= b[i]; 71 fft(a, inv); 72 for(int i = 0; i < n; i++){ 73 res[i] += floor(a[i].real() / n + 0.5); 74 res[i + 1] += res[i] / 10; 75 res[i] %= 10; 76 } 77 for(int i = res[lena + lenb - 1] ? lena + lenb - 1: lena + lenb - 2; i >= 0; i--) 78 putchar('0' + res[i]); 79 enter; 80 return 0; 81 } 82 83 //可以预处理ωkn和ω−kn，分别存在omg和inv数组中。调用fft时，如果无需取倒数，则传入omg；如果需要取倒数，则传入inv。