Submit Info #62010

Problem Lang User Status Time Memory
Polynomial Interpolation cpp-acl cai_lw AC 1078 ms 33.75 MiB

ケース詳細
Name Status Time Memory
example_00 AC 1 ms 0.45 MiB
example_01 AC 1 ms 0.45 MiB
max_random_00 AC 1078 ms 33.71 MiB
max_random_01 AC 1078 ms 33.75 MiB
random_00 AC 1019 ms 33.14 MiB
random_01 AC 715 ms 17.93 MiB
random_02 AC 437 ms 13.98 MiB

#include <bits/stdc++.h> #include <atcoder/modint> #include <atcoder/convolution> using namespace std; using mint = atcoder::modint998244353; using atcoder::convolution; /* F is a field. Must support operator +, +=, -, -=, *, /. F(0) must return additive identity. F(1) must return multiplicative identity. */ template<typename F, vector<F> (*Conv)(const vector<F>&, const vector<F>&) = convolution<F>> struct Polynomial { using Self = Polynomial<F, Conv>; vector<F> coef; Polynomial() {} Polynomial(const vector<F> &_coef) : coef(_coef) {} Polynomial(vector<F> &&_coef) : coef(_coef) {} size_t size() const { return coef.size(); } void trim_leading_zeros() { while (!coef.empty() && coef.back() == F(0)) coef.pop_back(); } void operator+=(const Polynomial &other) { for (size_t i = 0; i < other.size(); i++) if (i < coef.size()) coef[i] += other.coef[i]; else coef.push_back(other.coef[i]); } Polynomial operator+(const Polynomial &other) const { Polynomial ret = *this; ret += other; return ret; } void operator-=(const Polynomial &other) { for (size_t i = 0; i < other.size(); i++) if (i < coef.size()) coef[i] -= other.coef[i]; else coef.push_back(-other.coef[i]); } Polynomial operator-(const Polynomial &other) const { Polynomial<F> ret = *this; ret -= other; return ret; } Polynomial operator*(const Polynomial &other) const { return Polynomial(Conv(this->coef, other.coef)); } void operator*=(const Polynomial &other) { *this = *this * other; } // f(x)g(x) = 1 (mod x^N) Polynomial recip_mod(size_t deg = 0) { size_t n = deg; if (n == 0) n = size(); if (n == 1) return Polynomial({ F(1) / coef[0] }); size_t n_low = (n + 1) / 2; Polynomial low_recip; if (size() <= n_low) low_recip = recip_mod(n_low); else { Polynomial low(vector<F>(coef.begin(), coef.begin() + n_low)); low_recip = low.recip_mod(); } Polynomial ret = low_recip * low_recip; if (size() <= n) ret *= *this; else { Polynomial trunc(vector<F>(coef.begin(), coef.begin() + n)); ret *= trunc; } ret.coef.resize(n, F(0)); for (size_t i = n_low; i < n; i++) ret.coef[i] = -ret.coef[i]; return ret; } Polynomial rev() const { Polynomial ret = *this; reverse(ret.coef.begin(), ret.coef.end()); return ret; } // f(x) = g(x)q(x) + r(x) deg(r) < deg(g) pair<Polynomial, Polynomial> div_mod(const Polynomial &other) const { if (size() < other.size()) return {Polynomial(), *this}; size_t quot_size = size() - other.size() + 1; Polynomial quot; if (quot_size > 1) { quot = this->rev() * other.rev().recip_mod(quot_size - 1); quot.coef.resize(quot_size, F(0)); reverse(quot.coef.begin(), quot.coef.end()); quot.coef[0] = F(0); } else { quot.coef.emplace_back(0); } Polynomial rem = *this - quot * other; quot.coef[0] = rem.coef[other.size() - 1] / other.coef.back(); for (size_t i = 0; i < other.size(); i++) rem.coef[i] -= quot.coef[0] * other.coef[i]; return {quot, rem}; } Polynomial operator/(const Polynomial &other) const { return div_mod(other).first; } void operator/=(const Polynomial &other) { *this = *this / other; } Polynomial operator%(const Polynomial &other) const { return div_mod(other).second; } void operator%=(const Polynomial &other) { *this = *this % other; } Polynomial derivative() const { if (size() <= 1) return Polynomial(); Polynomial deriv; deriv.coef.reserve(size() - 1); F multiplier(0); for (size_t i = 1; i < size(); i++) { multiplier += F(1); deriv.coef.push_back(coef[i] * multiplier); } return deriv; } vector<F> multipoint_evaluate(const vector<F> &points) const { size_t n = points.size(); unique_ptr<MultipointEvaluationTree> tree = make_unique<MultipointEvaluationTree>(); tree->build(points, 0, n); vector<F> ret; ret.reserve(n); tree->eval(*this, ret); return ret; } static Polynomial from_interpolate(const vector<F> &x, const vector<F> &y) { size_t n = x.size(); unique_ptr<MultipointEvaluationTree> tree = make_unique<MultipointEvaluationTree>(); tree->build(x, 0, n); vector<F> y_norm; y_norm.reserve(n); tree->eval(tree->p.derivative(), y_norm); for (size_t i = 0; i < n; i++) y_norm[i] = y[i] / y_norm[i]; return tree->interpolate(x, y_norm, 0, n); } private: struct MultipointEvaluationTree { unique_ptr<MultipointEvaluationTree> left, right; Polynomial p; void build(const vector<F> &points, size_t l, size_t r) { if (r - l == 1) { p = Polynomial({-points[l], 1}); return; } size_t m = l + (r - l) / 2; left = make_unique<MultipointEvaluationTree>(); left->build(points, l, m); right = make_unique<MultipointEvaluationTree>(); right->build(points, m, r); p = left->p * right-> p; } void eval(const Polynomial &f, vector<F> &out) { Polynomial fp = f % p; fp.trim_leading_zeros(); if (!left && !right) { out.push_back(fp.size() == 0 ? F(0) : fp.coef[0]); return; } if (left) left->eval(fp, out); if (right) right->eval(fp, out); } Polynomial interpolate(const vector<F> &x, const vector<F> &y_norm, size_t l, size_t r) { if (r - l == 1) return Polynomial({y_norm[l]}); size_t m = l + (r - l) / 2; Polynomial fl = left->interpolate(x, y_norm, l, m); Polynomial fr = right->interpolate(x, y_norm, m, r); return fl * right->p + fr * left->p; } }; }; // https://judge.yosupo.jp/problem/inv_of_formal_power_series void library_checker_inv_of_formal_power_series() { int n; cin >> n; vector<mint> a; a.reserve(n); for (int i = 0; i < n; i++) { int x; cin >> x; a.emplace_back(x); } Polynomial<mint> poly(move(a)); auto ans = poly.recip_mod(); for (int i = 0; i < n; i++) { cout << ans.coef[i].val(); cout << (i == n - 1 ? '\n' : ' '); } } // https://judge.yosupo.jp/problem/division_of_polynomials void library_checker_division_of_polynomials() { int n, m; cin >> n >> m; vector<mint> a, b; a.reserve(n); b.reserve(m); for (int i = 0; i < n; i++) { int x; cin >> x; a.emplace_back(x); } for (int i = 0; i < m; i++) { int x; cin >> x; b.emplace_back(x); } Polynomial<mint> f(move(a)), g(move(b)); auto [q, r] = f.div_mod(g); q.trim_leading_zeros(); r.trim_leading_zeros(); cout << q.size() << ' ' << r.size() << '\n'; for (int i = 0; i < q.size(); i++) { cout << q.coef[i].val(); if (i + 1 < q.size()) cout << ' '; } cout << '\n'; for (int i = 0; i < r.size(); i++) { cout << r.coef[i].val(); if (i + 1 < r.size()) cout << ' '; } cout << '\n'; } // https://judge.yosupo.jp/problem/multipoint_evaluation void library_checker_multipoint_evaluation() { int n, m; cin >> n >> m; vector<mint> a, p; a.reserve(n); p.reserve(m); for (int i = 0; i < n; i++) { int x; cin >> x; a.emplace_back(x); } for (int i = 0; i < m; i++) { int x; cin >> x; p.emplace_back(x); } Polynomial<mint> f(move(a)); vector<mint> v = f.multipoint_evaluate(p); for (int i = 0; i < m; i++) { cout << v[i].val(); cout << (i + 1 == m ? '\n' : ' '); } } // https://judge.yosupo.jp/problem/polynomial_interpolation void library_checker_polynomial_interpolation() { int n, m; cin >> n; vector<mint> a, b; a.reserve(n); b.reserve(n); for (int i = 0; i < n; i++) { int x; cin >> x; a.emplace_back(x); } for (int i = 0; i < n; i++) { int x; cin >> x; b.emplace_back(x); } auto f = Polynomial<mint>::from_interpolate(a, b); for (int i = 0; i < n; i++) { cout << f.coef[i].val(); cout << (i + 1 == n ? '\n' : ' '); } } int main() { ios::sync_with_stdio(false); cin.tie(nullptr); library_checker_polynomial_interpolation(); }