library
This documentation is automatically generated by online-judge-tools/verification-helper
fps/poly_find_root.hpp
- View this file on GitHub
- Last update: 2025-04-05 12:46:42+09:00
- Include:
#include "fps/poly_find_root.hpp"
Depends on
Verified with
Code
#ifndef KK2_FPS_POLY_FIND_ROOT_HPP
#define KK2_FPS_POLY_FIND_ROOT_HPP 1
#include <vector>
#include "../random/gen.hpp"
#include "mod_pow.hpp"
#include "poly_gcd.hpp"
namespace kk2 {
template <class FPS, class mint = typename FPS::value_type>
std::vector<mint> find_root(const FPS &f) {
long long p = mint::getmod();
std::vector<mint> res;
if (p == 2) {
for (int i = 0; i < 2; i++) {
if (f.eval(mint(i)) == mint(0)) { res.push_back(mint(i)); }
}
return res;
}
std::vector<FPS> fs;
fs.push_back(poly_gcd(mod_pow(p, FPS{0, 1}, f) - FPS{0, 1}, f));
while (!fs.empty()) {
auto g = fs.back();
fs.pop_back();
if (g.size() == 2) res.push_back(-g[0]);
if (g.size() <= 2) continue;
FPS s = FPS{random::rng(0, p), 1};
FPS t = poly_gcd(mod_pow((p - 1) / 2, s, g) - FPS{1}, g);
fs.push_back(t);
if (g.size() != t.size()) fs.push_back(g / t);
}
return res;
}
} // namespace kk2
#endif // KK2_FPS_POLY_FIND_ROOT_HPP#line 1 "fps/poly_find_root.hpp"
#include <vector>
#line 1 "random/gen.hpp"
#include <algorithm>
#include <cassert>
#include <numeric>
#include <random>
#include <unordered_set>
#line 10 "random/gen.hpp"
#line 1 "random/seed.hpp"
#include <chrono>
namespace kk2 {
namespace random {
using u64 = unsigned long long;
inline u64 non_deterministic_seed() {
u64 seed = std::chrono::duration_cast<std::chrono::nanoseconds>(
std::chrono::high_resolution_clock::now().time_since_epoch())
.count();
seed ^= reinterpret_cast<u64>(&seed);
seed ^= seed << 5;
seed ^= seed >> 41;
seed ^= seed << 20;
return seed;
}
inline u64 deterministic_seed() { return 5801799128519729247ull; }
inline u64 seed() {
#if defined(KK2_RANDOM_DETERMINISTIC)
return deterministic_seed();
#else
return non_deterministic_seed();
#endif
}
} // namespace random
} // namespace kk2
#line 12 "random/gen.hpp"
namespace kk2 {
namespace random {
using i64 = long long;
using u64 = unsigned long long;
inline u64 rng() {
static std::mt19937_64 mt(kk2::random::seed());
return mt();
}
// [l, r)
inline i64 rng(i64 l, i64 r) {
assert(l < r);
return l + rng() % (r - l);
}
// [l, r)
template <class T> std::vector<T> random_vector(int n, T l, T r) {
std::vector<T> res(n);
for (int i = 0; i < n; i++) res[i] = rng(l, r);
return res;
}
// [l, r)
std::vector<i64> distinct_rng(i64 l, i64 r, i64 n) {
assert(l < r and n <= r - l);
std::unordered_set<i64> st;
for (i64 i = n; i; --i) {
i64 m = rng(l, r + 1 - i);
if (st.find(m) != st.end()) m = r - i;
st.insert(m);
}
std::vector<i64> res(st.begin(), st.end());
std::sort(res.begin(), res.end());
return res;
}
template <class Iter> void shuffle(Iter first, Iter last) {
if (first == last) return;
int len = 1;
for (auto it = first + 1; it != last; ++it) {
len++;
int j = rng(0, len);
if (j != len - 1) std::iter_swap(first + j, it);
}
}
template <class T> std::vector<T> perm(int n) {
std::vector<T> res(n);
std::iota(res.begin(), res.end(), T(0));
shuffle(res.begin(), res.end());
return res;
}
template <class T> std::vector<T> choices(int l, int r, int k) {
assert(l < r and k <= r - l);
std::vector<T> res(r - l);
std::iota(res.begin(), res.end(), T(l));
shuffle(res.begin(), res.end());
res.resize(k);
return res;
}
} // namespace random
} // namespace kk2
#line 1 "fps/mod_pow.hpp"
#line 5 "fps/mod_pow.hpp"
namespace kk2 {
// return f ^ k mod g
template <class FPS, class mint = typename FPS::value_type, class T>
FPS mod_pow(T k, const FPS &f, const FPS &g) {
// assert(!is_signed_v<T> || k >= 0);
assert(!g.empty());
auto inv = g.rev().inv();
auto quo = [&](const FPS &poly) {
if (poly.size() < g.size()) return FPS{};
int n = poly.size() - g.size() + 1;
return (poly.rev().pre(n) * inv.pre(n)).pre(n).rev();
};
FPS res{1}, b(f);
while (k) {
if (k & 1) {
res *= b;
res -= quo(res) * g;
res.shrink();
}
b *= b;
b -= quo(b) * g;
b.shrink();
k >>= 1;
}
return res;
}
} // namespace kk2
#line 1 "fps/poly_gcd.hpp"
#line 5 "fps/poly_gcd.hpp"
#include <array>
#include <utility>
namespace kk2 {
namespace poly_gcd_impl {
template <class FPS> using Vec = std::array<FPS, 2>;
template <class FPS> struct mat_poly {
FPS a00, a01, a10, a11;
mat_poly() = default;
mat_poly(FPS a00_, FPS a01_, FPS a10_, FPS a11_) : a00(a00_), a01(a01_), a10(a10_), a11(a11_) {}
mat_poly &operator*=(const mat_poly &r) {
FPS A00 = a00 * r.a00 + a01 * r.a10;
FPS A01 = a00 * r.a01 + a01 * r.a11;
FPS A10 = a10 * r.a00 + a11 * r.a10;
FPS A11 = a10 * r.a01 + a11 * r.a11;
A00.shrink();
A01.shrink();
A10.shrink();
A11.shrink();
std::swap(a00, A00);
std::swap(a01, A01);
std::swap(a10, A10);
std::swap(a11, A11);
return *this;
}
static mat_poly identity() { return mat_poly(FPS{1}, FPS(), FPS(), FPS{1}); }
mat_poly operator*(const mat_poly &r) const { return mat_poly(*this) *= r; }
};
template <class FPS> Vec<FPS> operator*(const mat_poly<FPS> &a, const Vec<FPS> &b) {
FPS x0 = a.a00 * b[0] + a.a01 * b[1];
FPS x1 = a.a10 * b[0] + a.a11 * b[1];
x0.shrink();
x1.shrink();
return {x0, x1};
};
template <class FPS> void inner_naive_gcd(mat_poly<FPS> &a, Vec<FPS> &b) {
FPS quo = b[0] / b[1];
FPS rem = b[0] - quo * b[1];
FPS x10 = a.a00 - quo * a.a10;
FPS x11 = a.a01 - quo * a.a11;
rem.shrink();
x10.shrink();
x11.shrink();
std::swap(x10, a.a10);
std::swap(x11, a.a11);
std::swap(x10, a.a00);
std::swap(x11, a.a01);
b = {b[1], rem};
}
template <class FPS> mat_poly<FPS> inner_half_gcd(Vec<FPS> b) {
int n = (int)b[0].size(), m = (int)b[1].size();
int k = (n + 1) >> 1;
if (m <= k) return mat_poly<FPS>::identity();
mat_poly<FPS> m1 = inner_half_gcd(Vec<FPS>{b[0] >> k, b[1] >> k});
b = m1 * b;
if ((int)b[1].size() <= k) return m1;
inner_naive_gcd(m1, b);
if ((int)b[1].size() <= k) return m1;
int l = (int)b[0].size() - 1;
int j = 2 * k - l;
b[0] = b[0] >> j;
b[1] = b[1] >> j;
return inner_half_gcd(b) * m1;
}
template <class FPS> mat_poly<FPS> inner_poly_gcd(const FPS &a, const FPS &b) {
Vec<FPS> c{a, b};
c[0].shrink();
c[1].shrink();
int n = (int)c[0].size(), m = (int)c[1].size();
if (n < m) {
mat_poly<FPS> ret = inner_poly_gcd(c[1], c[0]);
std::swap(ret.a00, ret.a01);
std::swap(ret.a10, ret.a11);
return ret;
}
mat_poly<FPS> res = mat_poly<FPS>::identity();
while (1) {
mat_poly<FPS> m1 = inner_half_gcd(c);
c = m1 * c;
if (c[1].empty()) return m1 * res;
inner_naive_gcd(m1, c);
if (c[1].empty()) return m1 * res;
res = m1 * res;
}
}
template <class FPS> FPS poly_gcd(FPS a, FPS b) {
Vec<FPS> c{a, b};
mat_poly<FPS> m = inner_poly_gcd(a, b);
c = m * c;
if (!c[0].empty()) {
auto coeff = c[0].back().inv();
for (auto &x : c[0]) x *= coeff;
}
return c[0];
}
// f ^ {-1} mod g
template <class FPS> std::pair<bool, FPS> poly_inv(const FPS &f, const FPS &g) {
Vec<FPS> c{f, g};
mat_poly<FPS> m = inner_poly_gcd(f, g);
FPS gcd_ = (m * c)[0];
if (gcd_.size() != 1) return {0, FPS()};
Vec<FPS> x{FPS{1}, g};
return {1, ((m * x)[0] % g) * gcd_[0].inv()};
}
} // namespace poly_gcd_impl
using poly_gcd_impl::poly_gcd;
using poly_gcd_impl::poly_inv;
} // namespace kk2
#line 9 "fps/poly_find_root.hpp"
namespace kk2 {
template <class FPS, class mint = typename FPS::value_type>
std::vector<mint> find_root(const FPS &f) {
long long p = mint::getmod();
std::vector<mint> res;
if (p == 2) {
for (int i = 0; i < 2; i++) {
if (f.eval(mint(i)) == mint(0)) { res.push_back(mint(i)); }
}
return res;
}
std::vector<FPS> fs;
fs.push_back(poly_gcd(mod_pow(p, FPS{0, 1}, f) - FPS{0, 1}, f));
while (!fs.empty()) {
auto g = fs.back();
fs.pop_back();
if (g.size() == 2) res.push_back(-g[0]);
if (g.size() <= 2) continue;
FPS s = FPS{random::rng(0, p), 1};
FPS t = poly_gcd(mod_pow((p - 1) / 2, s, g) - FPS{1}, g);
fs.push_back(t);
if (g.size() != t.size()) fs.push_back(g / t);
}
return res;
}
} // namespace kk2