fork download
copy 
  1.  
  2. #include <bits/stdc++.h>
  3. using namespace std;
  4.  
  5. template <typename T>
  6. class Modular {
  7.  public:
  8.   using Type = typename decay<decltype(T::value)>::type;
  9.  
  10.   constexpr Modular() : value() {}
  11.   template <typename U>
  12.   Modular(const U& x) {
  13.     value = normalize(x);
  14.   }
  15.  
  16.   template <typename U>
  17.   static Type normalize(const U& x) {
  18.     Type v;
  19.     if (-mod() <= x && x < mod()) v = static_cast<Type>(x);
  20.     else v = static_cast<Type>(x % mod());
  21.     if (v < 0) v += mod();
  22.     return v;
  23.   }
  24.  
  25.   const Type& operator()() const { return value; }
  26.   template <typename U>
  27.   explicit operator U() const { return static_cast<U>(value); }
  28.   constexpr static Type mod() { return T::value; }
  29.  
  30.   Modular& operator+=(const Modular& other) { if ((value += other.value) >= mod()) value -= mod(); return *this; }
  31.   Modular& operator-=(const Modular& other) { if ((value -= other.value) < 0) value += mod(); return *this; }
  32.   template <typename U> Modular& operator+=(const U& other) { return *this += Modular(other); }
  33.   template <typename U> Modular& operator-=(const U& other) { return *this -= Modular(other); }
  34.   Modular& operator++() { return *this += 1; }
  35.   Modular& operator--() { return *this -= 1; }
  36.   Modular operator++(int) { Modular result(*this); *this += 1; return result; }
  37.   Modular operator--(int) { Modular result(*this); *this -= 1; return result; }
  38.   Modular operator-() const { return Modular(-value); }
  39.  
  40.   template <typename U = T>
  41.   typename enable_if<is_same<typename Modular<U>::Type, int>::value, Modular>::type& operator*=(const Modular& rhs) {
  42. #ifdef _WIN32
  43.     uint64_t x = static_cast<int64_t>(value) * static_cast<int64_t>(rhs.value);
  44.     uint32_t xh = static_cast<uint32_t>(x >> 32), xl = static_cast<uint32_t>(x), d, m;
  45.     asm(
  46.       "divl %4; \n\t"
  47.       : "=a" (d), "=d" (m)
  48.       : "d" (xh), "a" (xl), "r" (mod())
  49.     );
  50.     value = m;
  51. #else
  52.     value = normalize(static_cast<int64_t>(value) * static_cast<int64_t>(rhs.value));
  53. #endif
  54.     return *this;
  55.   }
  56.   template <typename U = T>
  57.   typename enable_if<is_same<typename Modular<U>::Type, int64_t>::value, Modular>::type& operator*=(const Modular& rhs) {
  58.     int64_t q = static_cast<int64_t>(static_cast<long double>(value) * rhs.value / mod());
  59.     value = normalize(value * rhs.value - q * mod());
  60.     return *this;
  61.   }
  62.   template <typename U = T>
  63.   typename enable_if<!is_integral<typename Modular<U>::Type>::value, Modular>::type& operator*=(const Modular& rhs) {
  64.     value = normalize(value * rhs.value);
  65.     return *this;
  66.   }
  67.  
  68.   Modular inverse() const {
  69.     Type a = value, b = mod(), u = 0, v = 1;
  70.     while (a != 0) {
  71.       Type t = b / a;
  72.       b -= t * a; swap(a, b);
  73.       u -= t * v; swap(u, v);
  74.     }
  75.     assert(b == 1);
  76.     return Modular(u);
  77.   }
  78.   Modular& operator/=(const Modular& other) { return *this *= other.inverse(); }
  79.  
  80.   template <typename U>
  81.   friend bool operator==(const Modular<U>& lhs, const Modular<U>& rhs);
  82.  
  83.   template <typename U>
  84.   friend std::istream& operator>>(std::istream& stream, Modular<U>& number);
  85.  
  86.  private:
  87.   Type value;
  88. };
  89.  
  90. template <typename T> bool operator==(const Modular<T>& lhs, const Modular<T>& rhs) { return lhs.value == rhs.value; }
  91. template <typename T, typename U> bool operator==(const Modular<T>& lhs, U rhs) { return lhs == Modular<T>(rhs); }
  92. template <typename T, typename U> bool operator==(U lhs, const Modular<T>& rhs) { return Modular<T>(lhs) == rhs; }
  93.  
  94. template <typename T> bool operator!=(const Modular<T>& lhs, const Modular<T>& rhs) { return !(lhs == rhs); }
  95. template <typename T, typename U> bool operator!=(const Modular<T>& lhs, U rhs) { return !(lhs == rhs); }
  96. template <typename T, typename U> bool operator!=(U lhs, const Modular<T>& rhs) { return !(lhs == rhs); }
  97.  
  98. template <typename T> Modular<T> operator+(const Modular<T>& lhs, const Modular<T>& rhs) { return Modular<T>(lhs) += rhs; }
  99. template <typename T, typename U> Modular<T> operator+(const Modular<T>& lhs, U rhs) { return Modular<T>(lhs) += rhs; }
  100. template <typename T, typename U> Modular<T> operator+(U lhs, const Modular<T>& rhs) { return Modular<T>(lhs) += rhs; }
  101.  
  102. template <typename T> Modular<T> operator-(const Modular<T>& lhs, const Modular<T>& rhs) { return Modular<T>(lhs) -= rhs; }
  103. template <typename T, typename U> Modular<T> operator-(const Modular<T>& lhs, U rhs) { return Modular<T>(lhs) -= rhs; }
  104. template <typename T, typename U> Modular<T> operator-(U lhs, const Modular<T>& rhs) { return Modular<T>(lhs) -= rhs; }
  105.  
  106. template <typename T> Modular<T> operator*(const Modular<T>& lhs, const Modular<T>& rhs) { return Modular<T>(lhs) *= rhs; }
  107. template <typename T, typename U> Modular<T> operator*(const Modular<T>& lhs, U rhs) { return Modular<T>(lhs) *= rhs; }
  108. template <typename T, typename U> Modular<T> operator*(U lhs, const Modular<T>& rhs) { return Modular<T>(lhs) *= rhs; }
  109.  
  110. template <typename T> Modular<T> operator/(const Modular<T>& lhs, const Modular<T>& rhs) { return Modular<T>(lhs) /= rhs; }
  111. template <typename T, typename U> Modular<T> operator/(const Modular<T>& lhs, U rhs) { return Modular<T>(lhs) /= rhs; }
  112. template <typename T, typename U> Modular<T> operator/(U lhs, const Modular<T>& rhs) { return Modular<T>(lhs) /= rhs; }
  113.  
  114. template<typename T, typename U>
  115. Modular<T> power(const Modular<T>& a, const U& b) {
  116.   assert(b >= 0);
  117.   Modular<T> x = a, res = 1;
  118.   U p = b;
  119.   while (p > 0) {
  120.     if (p & 1) res *= x;
  121.     x *= x;
  122.     p >>= 1;
  123.   }
  124.   return res;
  125. }
  126.  
  127. template <typename T>
  128. string to_string(const Modular<T>& number) {
  129.   return to_string(number());
  130. }
  131.  
  132. template <typename T>
  133. std::ostream& operator<<(std::ostream& stream, const Modular<T>& number) {
  134.   return stream << number();
  135. }
  136.  
  137. template <typename T>
  138. std::istream& operator>>(std::istream& stream, Modular<T>& number) {
  139.   typename common_type<typename Modular<T>::Type, int64_t>::type x;
  140.   stream >> x;
  141.   number.value = Modular<T>::normalize(x);
  142. //  stream >> number.value;
  143. //  number.value = Modular<T>::normalize(number.value);
  144.   return stream;
  145. }
  146.  
  147. /*
  148. using ModType = int;
  149.  
  150. struct VarMod { static ModType value; };
  151. ModType VarMod::value;
  152. ModType& md = VarMod::value;
  153. using Mint = Modular<VarMod>;
  154. */
  155.  
  156. constexpr int md = 998244353;
  157. using Mint = Modular<std::integral_constant<decay<decltype(md)>::type, md>>;
  158.  
  159. template <typename T>
  160. class NTT {
  161.  public:
  162.   using Type = typename decay<decltype(T::value)>::type;
  163.  
  164.   static Type md;
  165.   static Modular<T> root;
  166.   static int base;
  167.   static int max_base;
  168.   static vector<Modular<T>> roots;
  169.   static vector<int> rev;
  170.  
  171.   static void clear() {
  172.     root = 0;
  173.     base = 0;
  174.     max_base = 0;
  175.     roots.clear();
  176.     rev.clear();
  177.   }
  178.  
  179.   static void init() {
  180.     md = T::value;
  181.     assert(md >= 3 && md % 2 == 1);
  182.     auto tmp = md - 1;
  183.     max_base = 0;
  184.     while (tmp % 2 == 0) {
  185.       tmp /= 2;
  186.       max_base++;
  187.     }
  188.     root = 2;
  189.     while (power(root, (md - 1) >> 1) == 1) {
  190.       root++;
  191.     }
  192.     assert(power(root, md - 1) == 1);
  193.     root = power(root, (md - 1) >> max_base);
  194.     base = 1;
  195.     rev = {0, 1};
  196.     roots = {0, 1};
  197.   }
  198.  
  199.   static void ensure_base(int nbase) {
  200.     if (md != T::value) {
  201.       clear();
  202.     }
  203.     if (roots.empty()) {
  204.       init();
  205.     }
  206.     if (nbase <= base) {
  207.       return;
  208.     }
  209.     assert(nbase <= max_base);
  210.     rev.resize(1 << nbase);
  211.     for (int i = 0; i < (1 << nbase); i++) {
  212.       rev[i] = (rev[i >> 1] >> 1) + ((i & 1) << (nbase - 1));
  213.     }
  214.     roots.resize(1 << nbase);
  215.     while (base < nbase) {
  216.       Modular<T> z = power(root, 1 << (max_base - 1 - base));
  217.       for (int i = 1 << (base - 1); i < (1 << base); i++) {
  218.         roots[i << 1] = roots[i];
  219.         roots[(i << 1) + 1] = roots[i] * z;
  220.       }
  221.       base++;
  222.     }
  223.   }
  224.  
  225.   static void fft(vector<Modular<T>> &a) {
  226.     int n = (int) a.size();
  227.     assert((n & (n - 1)) == 0);
  228.     int zeros = __builtin_ctz(n);
  229.     ensure_base(zeros);
  230.     int shift = base - zeros;
  231.     for (int i = 0; i < n; i++) {
  232.       if (i < (rev[i] >> shift)) {
  233.         swap(a[i], a[rev[i] >> shift]);
  234.       }
  235.     }
  236.     for (int k = 1; k < n; k <<= 1) {
  237.       for (int i = 0; i < n; i += 2 * k) {
  238.         for (int j = 0; j < k; j++) {
  239.           Modular<T> x = a[i + j];
  240.           Modular<T> y = a[i + j + k] * roots[j + k];
  241.           a[i + j] = x + y;
  242.           a[i + j + k] = x - y;
  243.         }
  244.       }
  245.     }
  246.   }
  247.  
  248.   static vector<Modular<T>> multiply(vector<Modular<T>> a, vector<Modular<T>> b) {
  249.     if (a.empty() || b.empty()) {
  250.       return {};
  251.     }
  252.     int eq = (a.size() == b.size() && a == b);
  253.     int need = (int) a.size() + (int) b.size() - 1;
  254.     int nbase = 0;
  255.     while ((1 << nbase) < need) nbase++;
  256.     ensure_base(nbase);
  257.     int sz = 1 << nbase;
  258.     a.resize(sz);
  259.     b.resize(sz);
  260.     fft(a);
  261.     if (eq) b = a; else fft(b);
  262.     Modular<T> inv_sz = 1 / static_cast<Modular<T>>(sz);
  263.     for (int i = 0; i < sz; i++) {
  264.       a[i] *= b[i] * inv_sz;
  265.     }
  266.     reverse(a.begin() + 1, a.end());
  267.     fft(a);
  268.     a.resize(need);
  269.     return a;
  270.   }
  271. };
  272.  
  273. template <typename T> typename NTT<T>::Type NTT<T>::md;
  274. template <typename T> Modular<T> NTT<T>::root;
  275. template <typename T> int NTT<T>::base;
  276. template <typename T> int NTT<T>::max_base;
  277. template <typename T> vector<Modular<T>> NTT<T>::roots;
  278. template <typename T> vector<int> NTT<T>::rev;
  279.  
  280. template <typename T>
  281. vector<Modular<T>> inverse(const vector<Modular<T>>& a) {
  282.   assert(!a.empty());
  283.   int n = (int) a.size();
  284.   vector<Modular<T>> b = {1 / a[0]};
  285.   while ((int) b.size() < n) {
  286.     vector<Modular<T>> x(a.begin(), a.begin() + min(a.size(), b.size() << 1));
  287.     x.resize(b.size() << 1);
  288.     b.resize(b.size() << 1);
  289.     vector<Modular<T>> c = b;
  290.     NTT<T>::fft(c);
  291.     NTT<T>::fft(x);
  292.     Modular<T> inv = 1 / static_cast<Modular<T>>((int) x.size());
  293.     for (int i = 0; i < (int) x.size(); i++) {
  294.       x[i] *= c[i] * inv;
  295.     }
  296.     reverse(x.begin() + 1, x.end());
  297.     NTT<T>::fft(x);
  298.     rotate(x.begin(), x.begin() + (x.size() >> 1), x.end());
  299.     fill(x.begin() + (x.size() >> 1), x.end(), 0);
  300.     NTT<T>::fft(x);
  301.     for (int i = 0; i < (int) x.size(); i++) {
  302.       x[i] *= c[i] * inv;
  303.     }
  304.     reverse(x.begin() + 1, x.end());
  305.     NTT<T>::fft(x);
  306.     for (int i = 0; i < ((int) x.size() >> 1); i++) {
  307.       b[i + ((int) x.size() >> 1)] = -x[i];
  308.     }
  309.   }
  310.   b.resize(n);
  311.   return b;
  312. }
  313.  
  314. template <typename T>
  315. vector<Modular<T>> inverse_old(vector<Modular<T>> a) {
  316.   assert(!a.empty());
  317.   int n = (int) a.size();
  318.   if (n == 1) {
  319.     return {1 / a[0]};
  320.   }
  321.   int m = (n + 1) >> 1;
  322.   vector<Modular<T>> b = inverse_old(vector<Modular<T>>(a.begin(), a.begin() + m));
  323.   int need = n << 1;
  324.   int nbase = 0;
  325.   while ((1 << nbase) < need) {
  326.     ++nbase;
  327.   }
  328.   NTT<T>::ensure_base(nbase);
  329.   int size = 1 << nbase;
  330.   a.resize(size);
  331.   b.resize(size);
  332.   NTT<T>::fft(a);
  333.   NTT<T>::fft(b);
  334.   Modular<T> inv = 1 / static_cast<Modular<T>>(size);
  335.   for (int i = 0; i < size; ++i) {
  336.     a[i] = (2 - a[i] * b[i]) * b[i] * inv;
  337.   }
  338.   reverse(a.begin() + 1, a.end());
  339.   NTT<T>::fft(a);
  340.   a.resize(n);
  341.   return a;
  342. }
  343.  
  344. template <typename T>
  345. vector<Modular<T>> operator*(const vector<Modular<T>>& a, const vector<Modular<T>>& b) {
  346.   if (a.empty() || b.empty()) {
  347.     return {};
  348.   }
  349.   if (min(a.size(), b.size()) < 100) {
  350.     vector<Modular<T>> c(a.size() + b.size() - 1, 0);
  351.     for (int i = 0; i < (int) a.size(); i++) {
  352.       for (int j = 0; j < (int) b.size(); j++) {
  353.         c[i + j] += a[i] * b[j];
  354.       }
  355.     }
  356.     return c;
  357.   }
  358.   return NTT<T>::multiply(a, b);
  359. }
  360.  
  361. template <typename T>
  362. vector<Modular<T>>& operator*=(vector<Modular<T>>& a, const vector<Modular<T>>& b) {
  363.   return a = a * b;
  364. }
  365.  
  366. template <typename T>
  367. vector<T>& operator+=(vector<T>& a, const vector<T>& b) {
  368.   if (a.size() < b.size()) {
  369.     a.resize(b.size());
  370.   }
  371.   for (int i = 0; i < (int) b.size(); i++) {
  372.     a[i] += b[i];
  373.   }
  374.   return a;
  375. }
  376.  
  377. template <typename T>
  378. vector<T> operator+(const vector<T>& a, const vector<T>& b) {
  379.   vector<T> c = a;
  380.   return c += b;
  381. }
  382.  
  383. template <typename T>
  384. vector<T>& operator-=(vector<T>& a, const vector<T>& b) {
  385.   if (a.size() < b.size()) {
  386.     a.resize(b.size());
  387.   }
  388.   for (int i = 0; i < (int) b.size(); i++) {
  389.     a[i] -= b[i];
  390.   }
  391.   return a;
  392. }
  393.  
  394. template <typename T>
  395. vector<T> operator-(const vector<T>& a, const vector<T>& b) {
  396.   vector<T> c = a;
  397.   return c -= b;
  398. }
  399.  
  400. template <typename T>
  401. vector<T> operator-(const vector<T>& a) {
  402.   vector<T> c = a;
  403.   for (int i = 0; i < (int) c.size(); i++) {
  404.     c[i] = -c[i];
  405.   }
  406.   return c;
  407. }
  408. const int N=1e6+10;
  409. vector<Mint>v[N],tree[4*N];
  410. void build(int idx,int l,int r){
  411.     if(l>r)return;
  412.     if(l==r){tree[idx]=v[l];return;}
  413.     int mid=(l+r)/2;
  414.     build(2*idx+1,l,mid);
  415.     build(2*idx+2,mid+1,r);
  416.     tree[idx]=tree[2*idx+1]*tree[2*idx+2];
  417. }
  418. #define debug(...) fprintf(stdout, __VA_ARGS__), fflush(stderr)    
  419. #define time__(d)for(long blockTime=0;(blockTime==0?(blockTime=clock())!=0:false); debug("%s:%.4fs",d,(double)(clock()-blockTime)/CLOCKS_PER_SEC))
  420.  
  421. int main(){
  422.   time__("time required"){
  423.   long long n=1000000;
  424.   for(int i=0;i<N;i++){
  425.     v[i].resize(2);
  426.     v[i][0]=1;
  427.     v[i][1]=i+1;
  428.    }
  429.    cout<<(n*(n+1)/2)%md<<" ";
  430.    build(0,0,n-1);
  431.    //for(auto xx:tree[0])cout<<xx<<" ";
  432.  }
  433. }
Success #stdin #stdout 1.96s 302376KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
878323500 time required:1.7918s
https://ideone.com/LbMYij
language:
C++14 (gcc 8.3)
created:
4 years 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!