(qqqqc)自用模板
目录
dsu
组合数
取模
线段树
lazy线段树
2-sat
树状数组
FFT
SAM
Hash
int128
NTT
dsu
struct DSU {std::vector f, siz;DSU(int n) : f(n), siz(n, 1) { std::iota(f.begin(), f.end(), 0); }int leader(int x) {while (x != f[x]) x = f[x] = f[f[x]];return x;}bool same(int x, int y) { return leader(x) == leader(y); }bool merge(int x, int y) {x = leader(x);y = leader(y);if (x == y) return false;siz[x] += siz[y];f[y] = x;return true;}int size(int x) { return siz[leader(x)]; }
};
组合数
const int N = 2e5 + 7;
const int mod = 998244353;
int fac[N], inv[N];
void init(int n) {fac[0] = 1;for (int i = 1; i <= n; ++i) fac[i] = fac[i - 1] * i % mod;inv[n] = qpow(fac[n], mod - 2);for (int i = n - 1; ~i; --i) inv[i] = inv[i + 1] * (i + 1) % mod;
}
int C(int n, int m) {return n < m ? 0 : fac[n] * inv[m] % mod * inv[n - m] % mod;
}
取模
constexpr int P = 1000000007;
// assume -P <= x < 2P
int norm(int x) {if (x < 0) {x += P;}if (x >= P) {x -= P;}return x;
}
template
T power(T a, int b) {T res = 1;for (; b; b /= 2, a *= a) {if (b % 2) {res *= a;}}return res;
}
struct Z {int x;Z(int x = 0) : x(norm(x% P)) {}int val() const {return x;}Z operator-() const {return Z(norm(P - x));}Z inv() const {assert(x != 0);return power(*this, P - 2);}Z& operator*=(const Z& rhs) {x = x * rhs.x % P;return *this;}Z& operator+=(const Z& rhs) {x = norm(x + rhs.x);return *this;}Z& operator-=(const Z& rhs) {x = norm(x - rhs.x);return *this;}Z& operator/=(const Z& rhs) {return *this *= rhs.inv();}friend Z operator*(const Z& lhs, const Z& rhs) {Z res = lhs;res *= rhs;return res;}friend Z operator+(const Z& lhs, const Z& rhs) {Z res = lhs;res += rhs;return res;}friend Z operator-(const Z& lhs, const Z& rhs) {Z res = lhs;res -= rhs;return res;}friend Z operator/(const Z& lhs, const Z& rhs) {Z res = lhs;res /= rhs;return res;}friend std::istream& operator>>(std::istream& is, Z& a) {int v;is >> v;a = Z(v);return is;}friend std::ostream& operator<<(std::ostream& os, const Z& a) {return os << a.val();}
};std::vector fac(N + 1), invfac(N + 1);
fac[0] = 1;
for (int i = 1; i <= N; i++) {fac[i] = fac[i - 1] * i;
}
invfac[N] = fac[N].inv();
for (int i = N; i; i--) {invfac[i - 1] = invfac[i] * i;
}
线段树
struct Max {int x;Max(int x = 0) : x(x) {}
};Max operator+(const Max& a, const Max& b) {return std::max(a.x, b.x);
}template>struct SegmentTree {const int n;const Merge merge;std::vector info;SegmentTree(int n) : n(n), merge(Merge()), info(4 << std::__lg(n)) {}SegmentTree(std::vector init) : SegmentTree(init.size()) {std::function build = [&](int p, int l, int r) {if (r - l == 1) {info[p] = init[l];return;}int m = (l + r) / 2;build(2 * p, l, m);build(2 * p + 1, m, r);pull(p);};build(1, 0, n);}void pull(int p) {info[p] = merge(info[2 * p], info[2 * p + 1]);}void modify(int p, int l, int r, int x, const Info& v) {if (r - l == 1) {info[p] = v;return;}int m = (l + r) / 2;if (x < m) {modify(2 * p, l, m, x, v);} else {modify(2 * p + 1, m, r, x, v);}pull(p);}void modify(int p, const Info& v) {modify(1, 0, n, p, v);}Info rangeQuery(int p, int l, int r, int x, int y) {if (l >= y || r <= x) {return Info();}if (l >= x && r <= y) {return info[p];}int m = (l + r) / 2;return merge(rangeQuery(2 * p, l, m, x, y), rangeQuery(2 * p + 1, m, r, x, y));}Info rangeQuery(int l, int r) {return rangeQuery(1, 0, n, l, r);}
};
lazy线段树
struct Max {int x;Max(int x = -inf) : x(x) {}
};Max operator+(const Max& a, const Max& b) {return std::max(a.x, b.x);
}void apply(Max& a, const Max& b) {a.x = std::max(a.x, b.x);
}template>struct LazySegmentTree {const int n;const Merge merge;std::vector info;std::vector tag;LazySegmentTree(int n) : n(n), merge(Merge()), info(4 << std::__lg(n)), tag(4 << std::__lg(n)) {}LazySegmentTree(std::vector init) : LazySegmentTree(init.size()) {std::function build = [&](int p, int l, int r) {if (r - l == 1) {info[p] = init[l];return;}int m = (l + r) / 2;build(2 * p, l, m);build(2 * p + 1, m, r);pull(p);};build(1, 0, n);}void pull(int p) {info[p] = merge(info[2 * p], info[2 * p + 1]);}void apply(int p, const Tag& v) {::apply(info[p], v);::apply(tag[p], v);}void push(int p) {apply(2 * p, tag[p]);apply(2 * p + 1, tag[p]);tag[p] = Tag();}void modify(int p, int l, int r, int x, const Info& v) {if (r - l == 1) {info[p] = v;return;}int m = (l + r) / 2;push(p);if (x < m) {modify(2 * p, l, m, x, v);} else {modify(2 * p + 1, m, r, x, v);}pull(p);}void modify(int p, const Info& v) {modify(1, 0, n, p, v);}Info rangeQuery(int p, int l, int r, int x, int y) {if (l >= y || r <= x) {return Info();}if (l >= x && r <= y) {return info[p];}int m = (l + r) / 2;push(p);return merge(rangeQuery(2 * p, l, m, x, y), rangeQuery(2 * p + 1, m, r, x, y));}Info rangeQuery(int l, int r) {return rangeQuery(1, 0, n, l, r);}void rangeApply(int p, int l, int r, int x, int y, const Tag& v) {if (l >= y || r <= x) {return;}if (l >= x && r <= y) {apply(p, v);return;}int m = (l + r) / 2;push(p);rangeApply(2 * p, l, m, x, y, v);rangeApply(2 * p + 1, m, r, x, y, v);pull(p);}void rangeApply(int l, int r, const Tag& v) {return rangeApply(1, 0, n, l, r, v);}
};
2-sat
struct TwoSat {int n;std::vector> e;std::vector ans;TwoSat(int n) : n(n), e(2 * n), ans(n) {}void addClause(int u, bool f, int v, bool g) {e[2 * u + !f].push_back(2 * v + g);e[2 * v + !g].push_back(2 * u + f);}bool satisfiable() {std::vector id(2 * n, -1), dfn(2 * n, -1), low(2 * n, -1);std::vector stk;int now = 0, cnt = 0;std::function tarjan = [&](int u) {stk.push_back(u);dfn[u] = low[u] = now++;for (auto v : e[u]) {if (dfn[v] == -1) {tarjan(v);low[u] = std::min(low[u], low[v]);} else if (id[v] == -1) {low[u] = std::min(low[u], dfn[v]);}}if (dfn[u] == low[u]) {int v;do {v = stk.back();stk.pop_back();id[v] = cnt;} while (v != u);++cnt;}};for (int i = 0; i < 2 * n; ++i) if (dfn[i] == -1) tarjan(i);for (int i = 0; i < n; ++i) {if (id[2 * i] == id[2 * i + 1]) return false;ans[i] = id[2 * i] > id[2 * i + 1];}return true;}std::vector answer() { return ans; }
};
树状数组
template
struct Fenwick {int n;std::vector a;Fenwick(int n) : n(n), a(n + 1) {}void add(int x, T v) {for (int i = x; i <= n; i += i & -i) {a[i] += v;}}T sum(int x) {T ans = 0;for (int i = x; i > 0; i -= i & -i) {ans += a[i];}return ans;}T rangeSum(int l, int r) {return sum(r) - sum(l);}T find_kth(int k) {T ans = 0, cnt = 0;for (int i = 20; i >= 0; i--) {ans += (1 << i);if (ans > n || cnt + a[ans] >= k) {ans -= (1 << i);} else {cnt += a[ans];}}return ++ans;}
};
FFT
constexpr int P = 998244353;
// assume -P <= x < 2P
int norm(int x) {if (x < 0) {x += P;}if (x >= P) {x -= P;}return x;
}
template
T power(T a, int b) {T res = 1;for (; b; b /= 2, a *= a) {if (b % 2) {res *= a;}}return res;
}
struct Z {int x;Z(int x = 0) : x(norm(x)) {}int val() const {return x;}Z operator-() const {return Z(norm(P - x));}Z inv() const {assert(x != 0);return power(*this, P - 2);}Z& operator*=(const Z& rhs) {x = x * rhs.x % P;return *this;}Z& operator+=(const Z& rhs) {x = norm(x + rhs.x);return *this;}Z& operator-=(const Z& rhs) {x = norm(x - rhs.x);return *this;}Z& operator/=(const Z& rhs) {return *this *= rhs.inv();}friend Z operator*(const Z& lhs, const Z& rhs) {Z res = lhs;res *= rhs;return res;}friend Z operator+(const Z& lhs, const Z& rhs) {Z res = lhs;res += rhs;return res;}friend Z operator-(const Z& lhs, const Z& rhs) {Z res = lhs;res -= rhs;return res;}friend Z operator/(const Z& lhs, const Z& rhs) {Z res = lhs;res /= rhs;return res;}friend std::istream& operator>>(std::istream& is, Z& a) {int v;is >> v;a = Z(v);return is;}friend std::ostream& operator<<(std::ostream& os, const Z& a) {return os << a.val();}
};std::vector rev;
std::vector roots{ 0, 1 };
void dft(std::vector& a) {int n = a.size();if ((int) rev.size() != n) {int k = __builtin_ctz(n) - 1;rev.resize(n);for (int i = 0; i < n; i++) {rev[i] = rev[i >> 1] >> 1 | (i & 1) << k;}}for (int i = 0; i < n; i++) {if (rev[i] < i) {std::swap(a[i], a[rev[i]]);}}if ((int) roots.size() < n) {int k = __builtin_ctz(roots.size());roots.resize(n);while ((1 << k) < n) {Z e = power(Z(3), (P - 1) >> (k + 1));for (int i = 1 << (k - 1); i < (1 << k); i++) {roots[2 * i] = roots[i];roots[2 * i + 1] = roots[i] * e;}k++;}}for (int k = 1; k < n; k *= 2) {for (int i = 0; i < n; i += 2 * k) {for (int j = 0; j < k; j++) {Z u = a[i + j];Z v = a[i + j + k] * roots[k + j];a[i + j] = u + v;a[i + j + k] = u - v;}}}
}
void idft(std::vector& a) {int n = a.size();std::reverse(a.begin() + 1, a.end());dft(a);Z inv = (1 - P) / n;for (int i = 0; i < n; i++) {a[i] *= inv;}
}
struct Poly {std::vector a;Poly() {}Poly(const std::vector& a) : a(a) {}Poly(const std::initializer_list& a) : a(a) {}int size() const {return a.size();}void resize(int n) {a.resize(n);}Z operator[](int idx) const {if (idx < size()) {return a[idx];} else {return 0;}}Z& operator[](int idx) {return a[idx];}Poly mulxk(int k) const {auto b = a;b.insert(b.begin(), k, 0);return Poly(b);}Poly modxk(int k) const {k = std::min(k, size());return Poly(std::vector(a.begin(), a.begin() + k));}Poly divxk(int k) const {if (size() <= k) {return Poly();}return Poly(std::vector(a.begin() + k, a.end()));}friend Poly operator+(const Poly& a, const Poly& b) {std::vector res(std::max(a.size(), b.size()));for (int i = 0; i < (int) res.size(); i++) {res[i] = a[i] + b[i];}return Poly(res);}friend Poly operator-(const Poly& a, const Poly& b) {std::vector res(std::max(a.size(), b.size()));for (int i = 0; i < (int) res.size(); i++) {res[i] = a[i] - b[i];}return Poly(res);}friend Poly operator*(Poly a, Poly b) {if (a.size() == 0 || b.size() == 0) {return Poly();}int sz = 1, tot = a.size() + b.size() - 1;while (sz < tot) {sz *= 2;}a.a.resize(sz);b.a.resize(sz);dft(a.a);dft(b.a);for (int i = 0; i < sz; ++i) {a.a[i] = a[i] * b[i];}idft(a.a);a.resize(tot);return a;}friend Poly operator*(Z a, Poly b) {for (int i = 0; i < (int) b.size(); i++) {b[i] *= a;}return b;}friend Poly operator*(Poly a, Z b) {for (int i = 0; i < (int) a.size(); i++) {a[i] *= b;}return a;}Poly& operator+=(Poly b) {return (*this) = (*this) + b;}Poly& operator-=(Poly b) {return (*this) = (*this) - b;}Poly& operator*=(Poly b) {return (*this) = (*this) * b;}Poly deriv() const {if (a.empty()) {return Poly();}std::vector res(size() - 1);for (int i = 0; i < size() - 1; ++i) {res[i] = (i + 1) * a[i + 1];}return Poly(res);}Poly integr() const {std::vector res(size() + 1);for (int i = 0; i < size(); ++i) {res[i + 1] = a[i] / (i + 1);}return Poly(res);}Poly inv(int m) const {Poly x{ a[0].inv() };int k = 1;while (k < m) {k *= 2;x = (x * (Poly{ 2 } - modxk(k) * x)).modxk(k);}return x.modxk(m);}Poly log(int m) const {return (deriv() * inv(m)).integr().modxk(m);}Poly exp(int m) const {Poly x{ 1 };int k = 1;while (k < m) {k *= 2;x = (x * (Poly{ 1 } - x.log(k) + modxk(k))).modxk(k);}return x.modxk(m);}Poly pow(int k, int m) const {int i = 0;while (i < size() && a[i].val() == 0) {i++;}if (i == size() || 1LL * i * k >= m) {return Poly(std::vector(m));}Z v = a[i];auto f = divxk(i) * v.inv();return (f.log(m - i * k) * k).exp(m - i * k).mulxk(i * k) * power(v, k);}Poly sqrt(int m) const {Poly x{ 1 };int k = 1;while (k < m) {k *= 2;x = (x + (modxk(k) * x.inv(k)).modxk(k)) * ((P + 1) / 2);}return x.modxk(m);}Poly mulT(Poly b) const {if (b.size() == 0) {return Poly();}int n = b.size();std::reverse(b.a.begin(), b.a.end());return ((*this) * b).divxk(n - 1);}std::vector eval(std::vector x) const {if (size() == 0) {return std::vector(x.size(), 0);}const int n = std::max((int) x.size(), size());std::vector q(4 * n);std::vector ans(x.size());x.resize(n);std::function build = [&](int p, int l, int r) {if (r - l == 1) {q[p] = Poly{ 1, -x[l] };} else {int m = (l + r) / 2;build(2 * p, l, m);build(2 * p + 1, m, r);q[p] = q[2 * p] * q[2 * p + 1];}};build(1, 0, n);std::function work = [&](int p, int l, int r, const Poly& num) {if (r - l == 1) {if (l < (int) ans.size()) {ans[l] = num[0];}} else {int m = (l + r) / 2;work(2 * p, l, m, num.mulT(q[2 * p + 1]).modxk(m - l));work(2 * p + 1, m, r, num.mulT(q[2 * p]).modxk(r - m));}};work(1, 0, n, mulT(q[1].inv(n)));return ans;}
};//cdq + fft
//vector f, c, p, invp, w, sumw, invsumw;
// void cdq(int l, int r) {
// if (l == r) {
// f[l] += (f[l - 1] + c[l - 1]) * invp[l - 1];
// return;
// }// int mid = l + r >> 1;
// cdq(l, mid);// vector a(r - l + 2), b(r - l + 2);
// for (int i = l; i <= r; i++) {
// a[i - l] = w[i - l];
// b[i - l + 1] = (i <= mid) ? f[i] : 0;
// }// auto res = Poly(a) * Poly(b);// for (int i = mid + 1; i <= r; i++) {
// Z tmp = invp[i - 1] * invsumw[i - 1];
// tmp *= (1 - p[i - 1] + P);
// tmp *= res[i - l];
// f[i] -= tmp;
// }// cdq(mid + 1, r);
// }std::vector fac(N + 1), invfac(N + 1);
fac[0] = 1;
for (int i = 1; i <= N; i++) {fac[i] = fac[i - 1] * i;
}
invfac[N] = fac[N].inv();
for (int i = N; i; i--) {invfac[i - 1] = invfac[i] * i;
}
SAM
struct SuffixAutomaton {static constexpr int ALPHABET_SIZE = 26, N = 1e5;struct Node {int len;int link;int next[ALPHABET_SIZE];Node() : len(0), link(0), next{} {}} t[2 * N];int cntNodes;SuffixAutomaton() {cntNodes = 1;std::fill(t[0].next, t[0].next + ALPHABET_SIZE, 1);t[0].len = -1;}int extend(int p, int c) {if (t[p].next[c]) {int q = t[p].next[c];if (t[q].len == t[p].len + 1)return q;int r = ++cntNodes;t[r].len = t[p].len + 1;t[r].link = t[q].link;std::copy(t[q].next, t[q].next + ALPHABET_SIZE, t[r].next);t[q].link = r;while (t[p].next[c] == q) {t[p].next[c] = r;p = t[p].link;}return r;}int cur = ++cntNodes;t[cur].len = t[p].len + 1;while (!t[p].next[c]) {t[p].next[c] = cur;p = t[p].link;}t[cur].link = extend(p, c);return cur;}
};
Hash
const int L = 1e6 + 5;
const int HASH_CNT = 2;int hashBase[HASH_CNT] = { 29, 31 };
int hashMod[HASH_CNT] = { (int) 1e9 + 9, 998244353 };struct StringWithHash {char s[L];int ls;int hsh[HASH_CNT][L];int pwMod[HASH_CNT][L];void init() {ls = 0;for (int i = 0; i < HASH_CNT; ++i) {hsh[i][0] = 0;pwMod[i][0] = 1;}}StringWithHash() { init(); }void extend(char c) {s[++ls] = c;for (int i = 0; i < HASH_CNT; ++i) {pwMod[i][ls] =1ll * pwMod[i][ls - 1] * hashBase[i] % hashMod[i];hsh[i][ls] = (1ll * hsh[i][ls - 1] * hashBase[i] + c) % hashMod[i];}}vector getHash(int l, int r) {vector res(HASH_CNT, 0);for (int i = 0; i < HASH_CNT; ++i) {int t =(hsh[i][r] - 1ll * hsh[i][l - 1] * pwMod[i][r - l + 1]) % hashMod[i];t = (t + hashMod[i]) % hashMod[i];res[i] = t;}return res;}
};bool equal(const vector& h1, const vector& h2) {assert(h1.size() == h2.size());for (unsigned i = 0; i < h1.size(); i++)if (h1[i] != h2[i]) return false;return true;
}
int128
using i128 = __int128;std::ostream &operator<<(std::ostream &o, i128 n) {if (!n) {return o << 0;} else {std::string s;while (n) {s += n % 10 + '0';n /= 10;}std::reverse(s.begin(), s.end());return o << s;}
}i128 gcd(i128 a, i128 b) {while (b) {a %= b;std::swap(a, b);}return a;
}
NTT
const int MOD = 998244353;struct mod_int {int val;mod_int(long long v = 0) {if (v < 0)v = v % MOD + MOD;if (v >= MOD)v %= MOD;val = v;}static int mod_inv(int a, int m = MOD) {// https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Exampleint g = m, r = a, x = 0, y = 1;while (r != 0) {int q = g / r;g %= r; swap(g, r);x -= q * y; swap(x, y);}return x < 0 ? x + m : x;}explicit operator int() const {return val;}mod_int& operator+=(const mod_int& other) {val += other.val;if (val >= MOD) val -= MOD;return *this;}mod_int& operator-=(const mod_int& other) {val -= other.val;if (val < 0) val += MOD;return *this;}static unsigned fast_mod(uint64_t x, unsigned m = MOD) {/*#if !defined(_WIN32) || defined(_WIN64)return x % m;#endif// Optimized mod for Codeforces 32-bit machines.// x must be less than 2^32 * m for this to work, so that x / m fits in a 32-bit integer.unsigned x_high = x >> 32, x_low = (unsigned) x;unsigned quot, rem;asm("divl %4\n": "=a" (quot), "=d" (rem): "d" (x_high), "a" (x_low), "r" (m));return rem;*/return x % m;}mod_int& operator*=(const mod_int& other) {val = fast_mod((uint64_t) val * other.val);return *this;}mod_int& operator/=(const mod_int& other) {return *this *= other.inv();}friend mod_int operator+(const mod_int& a, const mod_int& b) { return mod_int(a) += b; }friend mod_int operator-(const mod_int& a, const mod_int& b) { return mod_int(a) -= b; }friend mod_int operator*(const mod_int& a, const mod_int& b) { return mod_int(a) *= b; }friend mod_int operator/(const mod_int& a, const mod_int& b) { return mod_int(a) /= b; }mod_int& operator++() {val = val == MOD - 1 ? 0 : val + 1;return *this;}mod_int& operator--() {val = val == 0 ? MOD - 1 : val - 1;return *this;}mod_int operator++(int) { mod_int before = *this; ++* this; return before; }mod_int operator--(int) { mod_int before = *this; --* this; return before; }mod_int operator-() const {return val == 0 ? 0 : MOD - val;}bool operator==(const mod_int& other) const { return val == other.val; }bool operator!=(const mod_int& other) const { return val != other.val; }mod_int inv() const {return mod_inv(val);}mod_int pow(long long p) const {assert(p >= 0);mod_int a = *this, result = 1;while (p > 0) {if (p & 1)result *= a;a *= a;p >>= 1;}return result;}friend ostream& operator<<(ostream& stream, const mod_int& m) {return stream << m.val;}
};
void debug(vector a) {for (auto i : a)cerr << i << ' ';cerr << endl;
}
void r_debug(vector a) {for (auto i : a)cerr << 1 / i << ' ';cerr << endl;
}
vector inv, factorial, inv_factorial;void prepare_factorials(int maximum) {inv.assign(maximum + 1, 1);// Make sure MOD is prime, which is necessary for the inverse algorithm below.for (int p = 2; p * p <= MOD; p++)assert(MOD % p != 0);for (int i = 2; i <= maximum; i++)inv[i] = inv[MOD % i] * (MOD - MOD / i);factorial.resize(maximum + 1);inv_factorial.resize(maximum + 1);factorial[0] = inv_factorial[0] = 1;for (int i = 1; i <= maximum; i++) {factorial[i] = i * factorial[i - 1];inv_factorial[i] = inv[i] * inv_factorial[i - 1];}
}
inline mod_int C(int n, int m) {assert(n >= m); assert(n < factorial.size());return factorial[n] * inv_factorial[m] * inv_factorial[n - m];
}
namespace NTT {vector roots = { 0, 1 };vector bit_reverse;int max_size = -1;mod_int root;bool is_power_of_two(int n) {return (n & (n - 1)) == 0;}int round_up_power_two(int n) {assert(n > 0);while (n & (n - 1))n = (n | (n - 1)) + 1;return n;}// Given n (a power of two), finds k such that n == 1 << k.int get_length(int n) {assert(is_power_of_two(n));return __builtin_ctz(n);}// Rearranges the indices to be sorted by lowest bit first, then second lowest, etc., rather than highest bit first.// This makes even-odd div-conquer much easier.void bit_reorder(int n, vector& values) {if ((int) bit_reverse.size() != n) {bit_reverse.assign(n, 0);int length = get_length(n);for (int i = 0; i < n; i++)bit_reverse[i] = (bit_reverse[i >> 1] >> 1) + ((i & 1) << (length - 1));}for (int i = 0; i < n; i++)if (i < bit_reverse[i])swap(values[i], values[bit_reverse[i]]);}void find_root() {int order = MOD - 1;max_size = 1;while (order % 2 == 0) {order /= 2;max_size *= 2;}root = 2;// Find a max_size-th primitive root of MOD.while (!(root.pow(max_size) == 1 && root.pow(max_size / 2) != 1))root++;}void prepare_roots(int n) {if (max_size < 0)find_root();assert(n <= max_size);if ((int) roots.size() >= n)return;int length = get_length(roots.size());roots.resize(n);// The roots array is set up such that for a given power of two n >= 2, roots[n / 2] through roots[n - 1] are// the first half of the n-th primitive roots of MOD.while (1 << length < n) {// z is a 2^(length + 1)-th primitive root of MOD.mod_int z = root.pow(max_size >> (length + 1));for (int i = 1 << (length - 1); i < 1 << length; i++) {roots[2 * i] = roots[i];roots[2 * i + 1] = roots[i] * z;}length++;}}void fft_iterative(int N, vector& values) {assert(is_power_of_two(N));prepare_roots(N);bit_reorder(N, values);for (int n = 1; n < N; n *= 2)for (int start = 0; start < N; start += 2 * n)for (int i = 0; i < n; i++) {mod_int even = values[start + i];mod_int odd = values[start + n + i] * roots[n + i];values[start + n + i] = even - odd;values[start + i] = even + odd;}}const int FFT_CUTOFF = 150;vector mod_multiply(vector left, vector right) {int n = left.size();int m = right.size();// Brute force when either n or m is small enough.if (min(n, m) < FFT_CUTOFF) {const uint64_t ULL_BOUND = numeric_limits::max() - (uint64_t) MOD * MOD;vector result(n + m - 1);for (int i = 0; i < n; i++)for (int j = 0; j < m; j++) {result[i + j] += (uint64_t) ((int) left[i]) * ((int) right[j]);if (result[i + j] > ULL_BOUND)result[i + j] %= MOD;}for (uint64_t& x : result)if (x >= MOD)x %= MOD;return vector(result.begin(), result.end());}int N = round_up_power_two(n + m - 1);left.resize(N);right.resize(N);bool equal = left == right;fft_iterative(N, left);if (equal)right = left;elsefft_iterative(N, right);mod_int inv_N = mod_int(N).inv();for (int i = 0; i < N; i++)left[i] *= right[i] * inv_N;reverse(left.begin() + 1, left.end());fft_iterative(N, left);left.resize(n + m - 1);return left;}vector mod_multiply_all(const vector>& polynomials) {if (polynomials.empty())return { 1 };struct compare_size {bool operator()(const vector& x, const vector& y) {return x.size() > y.size();}};priority_queue, vector>, compare_size> pq(polynomials.begin(), polynomials.end());while (pq.size() > 1) {vector a = pq.top(); pq.pop();vector b = pq.top(); pq.pop();pq.push(mod_multiply(a, b));}return pq.top();}
}
本文来自互联网用户投稿,文章观点仅代表作者本人,不代表本站立场,不承担相关法律责任。如若转载,请注明出处。 如若内容造成侵权/违法违规/事实不符,请点击【内容举报】进行投诉反馈!
