Code[VS] 3123 高精度练习之超大整数乘法

FFT 做高精度乘法
  1 #include <bits/stdc++.h>
  2 
  3 const double pi = acos(-1);
  4 
  5 struct complex
  6 {
  7     double a, b;
  8     
  9     inline complex(
 10         double _a = 0,
 11         double _b = 0) 
 12     {
 13         a = _a;
 14         b = _b;
 15     }
 16     
 17     inline friend complex operator +
 18     (const complex &a, const complex &b)
 19     {
 20         return complex(a.a + b.a, a.b + b.b);
 21     }
 22     
 23     inline friend complex operator - 
 24     (const complex &a, const complex &b)
 25     {
 26         return complex(a.a - b.a, a.b - b.b);
 27     }
 28     
 29     inline friend complex operator *
 30     (const complex &a, const complex &b)
 31     {
 32         return complex(a.a*b.a - a.b*b.b, a.a*b.b + a.b*b.a);
 33     }
 34     
 35     inline friend complex & operator +=
 36     (complex &a, const complex &b)
 37     {
 38         return a = a+b;
 39     }
 40     
 41     inline friend complex & operator -=
 42     (complex &a, const complex &b)
 43     {
 44         return a = a-b;
 45     }
 46     
 47     inline friend complex & operator *= 
 48     (complex &a, const complex &b)
 49     {
 50         return a = a*b;
 51     }
 52 };
 53 
 54 inline complex alpha(double a)
 55 {
 56     return complex(cos(a), sin(a));
 57 }
 58 
 59 typedef std::vector<complex> vec;
 60 
 61 vec FFT(const vec &a)
 62 {
 63     int n = a.size();
 64     
 65     if (n == 1)return a;
 66     
 67     complex w_k(1, 0);
 68     complex w_n = alpha(pi*2/n);
 69     
 70     vec ar[2], yr[2], y(n);
 71     
 72     for (int i = 0; i < n; ++i)
 73         ar[i & 1].push_back(a[i]);
 74         
 75     for (int i = 0; i < 2; ++i)
 76         yr[i] = FFT(ar[i]);
 77         
 78     for (int i = 0; i < n/2; ++i, w_k *= w_n)
 79     {
 80         y[i] = yr[0][i] + w_k * yr[1][i];
 81         y[i + n/2] = yr[0][i] - w_k * yr[1][i];
 82     }
 83         
 84     return y;
 85 }
 86 
 87 vec IFFT(const vec &a)
 88 {
 89     int n = a.size();
 90     
 91     if (n == 1)return a;
 92     
 93     complex w_k(1, 0);
 94     complex w_n = alpha(-pi*2/n);
 95     
 96     vec ar[2], yr[2], y(n);
 97     
 98     for (int i = 0; i < n; ++i)
 99         ar[i & 1].push_back(a[i]);
100         
101     for (int i = 0; i < 2; ++i)
102         yr[i] = IFFT(ar[i]);
103         
104     for (int i = 0; i < n/2; ++i, w_k *= w_n)
105     {
106         y[i] = yr[0][i] + w_k * yr[1][i];
107         y[i + n/2] = yr[0][i] - w_k * yr[1][i];
108     }
109         
110     return y;
111 }
112 
113 char s1[100005]; int len1;
114 char s2[100005]; int len2;
115 
116 vec v1, v2, p1, p2, mul, ans;
117 
118 signed main(void)
119 {
120     scanf("%s", s1); len1 = strlen(s1);
121     scanf("%s", s2); len2 = strlen(s2);
122     
123     int len = len1 + len2;
124     
125     while (len != (len&-len))++len;
126     
127     for (int i = len1 - 1; ~i; --i)v1.push_back(complex(s1[i] - '0', 0));
128     for (int i = len2 - 1; ~i; --i)v2.push_back(complex(s2[i] - '0', 0));
129     
130     while ((int)v1.size() < len)v1.push_back(complex());
131     while ((int)v2.size() < len)v2.push_back(complex());
132     
133     p1 = FFT(v1);
134     p2 = FFT(v2);
135     
136     for (int i = 0; i < len; ++i)
137         mul.push_back(p1[i] * p2[i]);
138         
139     ans = IFFT(mul);
140     
141     std::vector<int> ret;
142     
143     for (int i = 0; i < len; ++i)
144         ret.push_back((int)round(ans[i].a / len));
145         
146     for (int i = 0; i < len; ++i)
147         if (ret[i] >= 10)
148         {
149             ret[i + 1] += ret[i] / 10;
150             ret[i] %= 10;
151         }
152         
153     while (ret.size() != 1 && !ret[ret.size() - 1])
154         ret.pop_back();
155         
156     for (int i = ret.size() - 1; i >= 0; --i)
157         putchar('0' + ret[i]);
158     putchar('
');
159 }
@Author: YouSiki