fork download
copy 
  1. /**
  2. Input:
  3. Dòng 1: Số testcase
  4. Dòng 2: Gồm hai số n, m
  5. Dòng 3: Gồm n+1 hệ số của đa thức A bậc n (bao gồm cả hệ số của hạng tử x^0)
  6. Dòng 4: Gồm m+1 hệ số của đa thức B bậc m (bao gồm cả hệ số của hạng tử x^0)
  7.  
  8. Output:
  9. Với mỗi testcase in ra 1 dòng duy nhất gồm n+m+1 hệ số của đa thức tích C = A x B
  10.  
  11. Giới hạn:
  12. 1 <= n, m <= (1e6-1)
  13. -1e6 <= Hệ số hạng tử <= 1e6
  14. Tổng n+m trên tất cả các testcase không vượt quá 2e6
  15.  
  16. Ví dụ:
  17.  
  18. Input:
  19. 2
  20. 2 4
  21. 1 2 3
  22. -4 1 2 5 -2
  23. 3 2
  24. 1 4 -4 -1
  25. 1 3 5
  26.  
  27. Output:
  28. -4 -7 -8 12 14 11 -6
  29. 1 7 13 7 -23 -5
  30. **/
  31.  
  32. #include <bits/stdc++.h>
  33. #define up(i,a,b) for (int i = (int)a; i <= (int)b; i++)
  34. using namespace std;
  35. using i64 = int64_t;
  36. using u64 = uint64_t;
  37. using u128 = unsigned __int128;
  38.  
  39. constexpr u64 MOD = 2524775926340780033ULL; // p-1 = c*2^24, primitive root 3
  40. constexpr u64 G   = 3;
  41. constexpr int MAXN = (1 << 21); // du cho n+m+1 <= 2e6
  42.  
  43. u64 power_plain(u64 a, u64 e, u64 mod){
  44.     a %= mod;
  45.     u64 res = 1;
  46.     while (e){
  47.         if (e & 1) res = (u128)res * a % mod;
  48.         a = (u128)a * a % mod;
  49.         e >>= 1;
  50.     }
  51.     return res;
  52. }
  53.  
  54. struct Montgomery64 {
  55.     u64 mod, inv_mod, r2;
  56.     void set_mod(u64 m){
  57.         mod = m;
  58.         inv_mod = 1;
  59.         for (int i = 0; i < 6; i++) inv_mod *= 2 - mod * inv_mod;
  60.         u64 r = (u64)(((u128)1 << 64) % mod);
  61.         r2 = (u64)((u128)r * r % mod);
  62.     }
  63.     u64 reduce(u128 x) const {
  64.         u64 q = (u64)x * inv_mod;
  65.         u128 m = (u128)q * mod;
  66.         u64 y = (u64)((x - m) >> 64);
  67.         return (y >> 63) ? y + mod : y;
  68.     }
  69.     u64 to_mont(u64 a)   const { return reduce((u128)a * r2); }
  70.     u64 from_mont(u64 a) const { return reduce((u128)a); }
  71.     u64 mul(u64 a, u64 b) const { return reduce((u128)a * b); }
  72.     u64 add(u64 a, u64 b) const { return a + b >= mod ? a + b - mod : a + b; }
  73.     u64 sub(u64 a, u64 b) const { return a >= b ? a - b : a + mod - b; }
  74. };
  75.  
  76. Montgomery64 mt;
  77. vector<u64> root, root_inv; // bang precompute, luu dang Montgomery
  78.  
  79. // build mot lan cho n = MAXN, dung chung cho moi bound <= MAXN
  80. void precompute_root(int n){
  81.     root.resize(n); root_inv.resize(n);
  82.     root[1] = mt.to_mont(1);
  83.     root_inv[1] = mt.to_mont(1);
  84.     u64 g_inv = power_plain(G, MOD - 2, MOD);
  85.     for (int k = 2; k * 2 <= n; k <<= 1){
  86.         u64 w         = power_plain(G,     (MOD - 1) / (2 * k), MOD);
  87.         u64 w_inv     = power_plain(g_inv, (MOD - 1) / (2 * k), MOD);
  88.         u64 w_mont     = mt.to_mont(w);
  89.         u64 w_inv_mont = mt.to_mont(w_inv);
  90.         for (int j = k / 2; j < k; j++){
  91.             root[j * 2]         = root[j];
  92.             root[j * 2 + 1]     = mt.mul(root[j], w_mont);
  93.             root_inv[j * 2]     = root_inv[j];
  94.             root_inv[j * 2 + 1] = mt.mul(root_inv[j], w_inv_mont);
  95.         }
  96.     }
  97. }
  98.  
  99. void ntt_forward(vector<u64>& a){
  100.     int n = a.size();
  101.     for (int len = n; len >= 2; len >>= 1){
  102.         int half = len / 2;
  103.         for (int i = 0; i < n; i += len){
  104.             for (int j = 0; j < half; j++){
  105.                 u64 w = root[half + j];
  106.                 u64 u = a[i + j];
  107.                 u64 v = a[i + j + half];
  108.                 a[i + j]        = mt.add(u, v);
  109.                 a[i + j + half] = mt.mul(mt.sub(u, v), w);
  110.             }
  111.         }
  112.     }
  113. }
  114.  
  115. void ntt_inverse(vector<u64>& a){
  116.     int n = a.size();
  117.     for (int len = 2; len <= n; len <<= 1){
  118.         int half = len / 2;
  119.         for (int i = 0; i < n; i += len){
  120.             for (int j = 0; j < half; j++){
  121.                 u64 w = root_inv[half + j];
  122.                 u64 u = a[i + j];
  123.                 u64 v = mt.mul(a[i + j + half], w);
  124.                 a[i + j]        = mt.add(u, v);
  125.                 a[i + j + half] = mt.sub(u, v);
  126.             }
  127.         }
  128.     }
  129.     u64 n_inv = mt.to_mont(power_plain((u64)n, MOD - 2, MOD));
  130.     for (auto& x : a) x = mt.mul(x, n_inv);
  131. }
  132.  
  133. vector<u64> multiply_mod(const vector<i64>& a, const vector<i64>& b, int bound){
  134.     vector<u64> A(bound), B(bound);
  135.     up(i, 0, a.size()-1) A[i] = mt.to_mont(a[i] < 0 ? (u64)(a[i] + (i64)MOD) : (u64)a[i]);
  136.     up(i, 0, b.size()-1) B[i] = mt.to_mont(b[i] < 0 ? (u64)(b[i] + (i64)MOD) : (u64)b[i]);
  137.  
  138.     ntt_forward(A);
  139.     ntt_forward(B);
  140.     up(i, 0, bound-1) A[i] = mt.mul(A[i], B[i]);
  141.     ntt_inverse(A);
  142.  
  143.     up(i, 0, bound-1) A[i] = mt.from_mont(A[i]);
  144.     return A;
  145. }
  146.  
  147. // |result| <= 1e18 < MOD/2
  148. i64 to_signed(u64 x){
  149.     return (x > MOD / 2) ? (i64)x - (i64)MOD : (i64)x;
  150. }
  151.  
  152. void solve(){
  153.     int n, m; cin >> n >> m;
  154.     vector<i64> a(n + 1), b(m + 1);
  155.     for (auto& x : a) cin >> x;
  156.     for (auto& x : b) cin >> x;
  157.  
  158.     int need = (int)a.size() + (int)b.size() - 1;
  159.     int bound = 1;
  160.     while (bound < need) bound <<= 1;
  161.  
  162.     vector<u64> raw = multiply_mod(a, b, bound);
  163.     up(i, 0, need-1) cout << to_signed(raw[i]) << " \n"[i == need - 1];
  164. }
  165.  
  166. signed main(){
  167.     ios_base::sync_with_stdio(false);
  168.     cin.tie(0);
  169.     #define Task "A"
  170.     if (fopen(Task".inp", "r")){
  171.         freopen(Task".inp", "r", stdin);
  172.         freopen(Task".out", "w", stdout);
  173.     }
  174.  
  175.     mt.set_mod(MOD);
  176.     precompute_root(MAXN); // build 1 lan, dung cho moi bound <= MAXN
  177.  
  178.     int tt; cin >> tt;
  179.     while (tt--) solve();
  180. }
  181.  
Success #stdin #stdout 0.02s 36040KB
comments (?)
stdin
copy 
2
2 4
1 2 3
-4 1 2 5 -2
3 2
1 4 -4 -1
1 3 5
stdout
copy 
-4 -7 -8 12 14 11 -6
1 7 13 7 -23 -5
https://ideone.com/uRLldC
ntt with 1 modulos (montgomery64)
language:
C++14 (gcc 8.3)
created:
1 day ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!