【HDU1402】【FNT版】A * B Problem Plus

  1 /*
  2 宋代朱敦儒
  3 《西江月·世事短如春梦》
  4 世事短如春梦，人情薄似秋云。不须计较苦劳心。万事原来有命。
  5 幸遇三杯酒好，况逢一朵花新。片时欢笑且相亲。明日阴晴未定。 
  6 */
  7 #include <cstdio>
  8 #include <cstring>
  9 #include <algorithm>
 10 #include <cmath>
 11 #include <queue>
 12 #include <vector>
 13 #include <iostream>
 14 #include <string>
 15 #include <ctime>
 16 #define LOCAL
 17 const int MAXN = 202000 + 10;
 18 const long long MOD = (479ll<<21ll) + 1;//费马数论变换的费马素数 
 19 long long G = 3;//元根 
 20 using namespace std;
 21 typedef long long ll;
 22 ll x1[MAXN], x2[MAXN];
 23 ll Inv[MAXN], wn[MAXN];
 24 ll Inv2[MAXN];
 25 char a[MAXN], b[MAXN];
 26 
 27 ll pow(ll a, ll b){
 28    if (b == 0) return 1 % MOD;
 29    if (b == 1) return a % MOD;
 30    ll tmp = pow(a, b / 2);
 31    if (b % 2 == 0) return (tmp * tmp) % MOD;
 32    else return ((tmp * tmp) % MOD * (a % MOD)) % MOD;
 33 }
 34 ll exgcd(ll a, ll b, ll &x, ll &y){
 35    if (b == 0){x = 1ll; y = 0; return a;}
 36    ll tmp = exgcd(b, a % b, y, x);
 37    y -= x * (a / b);
 38    return tmp;
 39 }
 40 ll inv(ll a, ll p){
 41    ll x, y;
 42    ll tmp = exgcd(a, p, x, y);
 43    return ((x % MOD) + MOD) % MOD; 
 44 }
 45 void change(ll *x, int len, int loglen){
 46      for (int i = 0; i < len; i++){
 47          int t = i, k = 0, tmp = loglen;
 48          while (tmp--) {k = (k << 1) + (t & 1); t >>= 1;}
 49          if (k < i) swap(x[i], x[k]);
 50      }
 51      return;
 52 }
 53 void FNT(ll *x, int len, int loglen, int type){
 54      if (type) change(x, len, loglen);
 55      int t;
 56      t = (type ? 1: (1 << loglen));
 57      for (int i = 0; i < loglen; i++){
 58          if (!type) t >>= 1;
 59          int l = 0, r = l + t;
 60          while (l < len){
 61                ll a, b;
 62                ll tmp = 1ll, w = wn[t] % MOD;
 63                if (!type) w = Inv[t];
 64                for (int j = l; j < l + t; j++){
 65                    if (type){
 66                       a = x[j] % MOD;
 67                       b = (x[j + t] * (tmp % MOD)) % MOD;
 68                       x[j] = (a + b) % MOD;
 69                       x[j + t] = ((a - b) % MOD + MOD) % MOD; 
 70                    }else{
 71                       a = (x[j] + x[j + t]) % MOD;
 72                       b = ((((x[j] - x[j + t]) % MOD + MOD) % MOD) * tmp) % MOD;
 73                       x[j] = a;
 74                       x[j + t] = b;
 75                    }
 76                    tmp = (tmp * w) % MOD;
 77                }
 78                l = r + t;
 79                r = l + t;
 80          } 
 81          if (type) t <<= 1;
 82      }
 83      if (!type){
 84         change(x, len, loglen);
 85         for (int i = 0; i < len; i++) x[i] = (x[i] % MOD * Inv2[len]) % MOD;
 86      }
 87 }
 88 int work(char *a, char *b){
 89      int len1, len2, loglen = 0;
 90      len1 = strlen(a);
 91      len2 = strlen(b);
 92      while ((1 << loglen) < (max(len1, len2) << 1)) loglen++;
 93      for (int i = 0; i < len1; i++) x1[i] = 1ll * (a[i] - ‘0‘);
 94      for (int i = 0; i < len2; i++) x2[i] = 1ll * (b[i] - ‘0‘);
 95      for (int i = len1; i < (1<<loglen); i++) x1[i] = 0;
 96      for (int i = len2; i < (1<<loglen); i++) x2[i] = 0;
 97      FNT(x1, (1<<loglen), loglen, 1);
 98      FNT(x2, (1<<loglen), loglen, 1);
 99      for (int i = 0; i < (1<<loglen); i++) x1[i] = (x1[i] * x2[i]) % MOD;
100      FNT(x1, (1<<loglen), loglen, 0);
101      
102      int len;
103      ll over, t;
104      over = 0;
105      len = 0;
106      for (int i = (len1 + len2) - 2; i >= 0; i--){
107          t = x1[i] + over;
108          a[len++] = t % 10;
109          over = t / 10;
110      }
111      while (over){
112            a[len++] = over % 10;
113            over /= 10;
114      }
115      return len;
116 }
117 void print(char *t, int len){
118      for (; len >= 0 && !t[len]; len--);
119      if (len < 0) printf("0");
120      else for (int i = len; i >= 0; i--) printf("%c", t[i] + ‘0‘);
121      printf("\n");
122 }
123 void prepare(){
124      wn[0] = MOD - 1;
125      for (int i = 1; i < 200010; i++){
126          if (i != 1 && (MOD - 1) / (i * 2) == (MOD - 1) / ((i - 1) * 2)){
127             wn[i] = wn[i - 1];
128             Inv[i] = Inv[i - 1];
129          }else{
130             wn[i] = pow(G, (MOD - 1) / (2 * i));
131             Inv[i] = inv(wn[i], MOD);
132          }
133          Inv2[i] = inv(i, MOD); 
134      }
135 }
136 
137 int main() {
138     memset(a, 0, sizeof(a));
139     memset(b, 0, sizeof(b));
140     prepare();
141     while (scanf("%s%s", a, b) != EOF){
142           print(a, work(a, b));
143           memset(a, 0, sizeof(a));
144           memset(b, 0, sizeof(b));
145     }
146     return 0;
147 }