This documentation is automatically generated by online-judge-tools/verification-helper
View the Project on GitHub ebi-fly13/Library
#include "math/sum_of_powers_iota.hpp"
$\sum_{i = 0}^{n-1} i ^ k$ を $O(k)$ で求める。
#pragma once #include <cassert> #include "../math/lagrange_interpolation.hpp" #include "../math/linear_sieve.hpp" #include "../modint/base.hpp" namespace ebi { template <Modint mint> mint sum_of_powers_iota(long long n, int k) { assert(n > 0 && k >= 0); linear_sieve sieve(k + 1); auto pow_table = sieve.pow_table<mint>(k + 1, k); for (int i = 0; i < k + 1; i++) { pow_table[i + 1] += pow_table[i]; } return lagrange_interpolation(pow_table, n - 1); } } // namespace ebi
#line 2 "math/sum_of_powers_iota.hpp" #include <cassert> #line 2 "math/lagrange_interpolation.hpp" #include <vector> /* reference: https://atcoder.jp/contests/abc208/editorial/2195 verify: https://atcoder.jp/contests/abc208/tasks/abc208_f */ namespace ebi { template <class mint> mint lagrange_interpolation(const std::vector<mint> &f, long long n) { const int d = int(f.size()) - 1; // Nのd次以下の多項式 mint fact = 1; std::vector<mint> inv_fact(d + 1); for (int i = 1; i < d + 1; ++i) { fact *= i; } inv_fact[d] = fact.inv(); for (int i = d; i > 0; i--) { inv_fact[i - 1] = inv_fact[i] * i; } std::vector<mint> l(d + 1), r(d + 1); l[0] = 1; for (int i = 0; i < d; ++i) { l[i + 1] = l[i] * (n - i); } r[d] = 1; for (int i = d; i > 0; --i) { r[i - 1] = r[i] * (n - i); } mint res = 0; for (int i = 0; i < d + 1; ++i) { res += mint((d - i) % 2 == 1 ? -1 : 1) * f[i] * l[i] * r[i] * inv_fact[i] * inv_fact[d - i]; } return res; } } // namespace ebi #line 2 "math/linear_sieve.hpp" #line 2 "modint/base.hpp" #include <concepts> #include <iostream> #include <utility> namespace ebi { template <class T> concept Modint = requires(T a, T b) { a + b; a - b; a * b; a / b; a.inv(); a.val(); a.pow(std::declval<long long>()); T::mod(); }; template <Modint mint> std::istream &operator>>(std::istream &os, mint &a) { long long x; os >> x; a = x; return os; } template <Modint mint> std::ostream &operator<<(std::ostream &os, const mint &a) { return os << a.val(); } } // namespace ebi #line 2 "template/int_alias.hpp" #include <cstdint> namespace ebi { using ld = long double; using std::size_t; using i8 = std::int8_t; using u8 = std::uint8_t; using i16 = std::int16_t; using u16 = std::uint16_t; using i32 = std::int32_t; using u32 = std::uint32_t; using i64 = std::int64_t; using u64 = std::uint64_t; using i128 = __int128_t; using u128 = __uint128_t; } // namespace ebi #line 5 "math/linear_sieve.hpp" /* reference: https://37zigen.com/linear-sieve/ verify: https://atcoder.jp/contests/abc162/submissions/25095562 */ #line 13 "math/linear_sieve.hpp" namespace ebi { struct linear_sieve { private: using u64 = std::uint64_t; int n; std::vector<int> sieve; std::vector<int> prime; public: linear_sieve(int _n) : n(_n), sieve(std::vector<int>(_n + 1, -1)) { for (int i = 2; i <= n; i++) { if (sieve[i] < 0) { sieve[i] = i; prime.emplace_back(i); } for (auto p : prime) { if (u64(p) * u64(i) > u64(n) || p > sieve[i]) break; sieve[p * i] = p; } } } std::vector<int> prime_table() const { return prime; } std::vector<std::pair<int, int>> prime_power_table(int m) const { assert(m <= n); std::vector<std::pair<int, int>> table(m + 1, {1, 1}); for (int i = 2; i <= m; i++) { int p = sieve[i]; table[i] = {p, p}; if (sieve[i / p] == p) { table[i] = table[i / p]; table[i].second *= p; } } return table; } std::vector<std::pair<int, int>> factorize(int x) { assert(x <= n); std::vector<std::pair<int, int>> res; while (x > 1) { int p = sieve[x]; int exp = 0; if (p < 0) { res.emplace_back(x, 1); break; } while (sieve[x] == p) { x /= p; exp++; } res.emplace_back(p, exp); } return res; } std::vector<int> divisors(int x) { assert(x <= n); std::vector<int> res; res.emplace_back(1); auto pf = factorize(x); for (auto p : pf) { int sz = (int)res.size(); for (int i = 0; i < sz; i++) { int ret = 1; for (int j = 0; j < p.second; j++) { ret *= p.first; res.emplace_back(res[i] * ret); } } } return res; } template <class T> std::vector<T> fast_zeta(const std::vector<T> &f) { std::vector<T> F = f; int sz = f.size(); assert(sz <= n + 1); for (int i = 2; i < sz; i++) { if (sieve[i] != i) continue; for (int j = (sz - 1) / i; j >= 1; j--) { F[j] += F[j * i]; } } return F; } template <class T> std::vector<T> fast_mobius(const std::vector<T> &F) { std::vector<T> f = F; int sz = F.size(); assert(sz <= n + 1); for (int i = 2; i < sz; i++) { if (sieve[i] != i) continue; for (int j = 1; j * i < sz; j++) { f[j] -= f[j * i]; } } return f; } template <Modint mint> std::vector<mint> pow_table(int m, int k) { assert(m <= n && k >= 0); std::vector<mint> table(m + 1, 1); table[0] = (k == 0); for (int i = 2; i <= m; i++) { if (sieve[i] == i) { table[i] = mint(i).pow(k); continue; } table[i] = table[sieve[i]] * table[i / sieve[i]]; } return table; } template <Modint mint> std::vector<mint> inv_table() { return pow_table(mint::mod() - 2); } }; } // namespace ebi #line 8 "math/sum_of_powers_iota.hpp" namespace ebi { template <Modint mint> mint sum_of_powers_iota(long long n, int k) { assert(n > 0 && k >= 0); linear_sieve sieve(k + 1); auto pow_table = sieve.pow_table<mint>(k + 1, k); for (int i = 0; i < k + 1; i++) { pow_table[i + 1] += pow_table[i]; } return lagrange_interpolation(pow_table, n - 1); } } // namespace ebi