結果
| 問題 | No.3620 Compositional Power with Schröder Coordinate 2 |
| コンテスト | |
| ユーザー |
 askr58
|
| 提出日時 | 2026-08-11 02:40:03 |
| 言語 | C++23 (gcc 15.2.0 + boost 1.90.0) |
| 結果 |
AC
|
| 実行時間 | 3,089 ms / 10,000 ms |
| + 294µs | |
| コード長 | 16,308 bytes |
| 記録 | |
| コンパイル時間 | 2,899 ms |
| コンパイル使用メモリ | 235,992 KB |
| 実行使用メモリ | 42,280 KB |
| 最終ジャッジ日時 | 2026-08-11 02:40:25 |
| 合計ジャッジ時間 | 20,140 ms |
|
ジャッジサーバーID （参考情報） |
judge1_0 / judge2_0 |
（要ログイン）
| ファイルパターン | 結果 |
|---|---|
| sample | AC * 2 |
| other | AC * 7 |
ソースコード
#include <iostream>
#include <ranges>
#include <algorithm>
#include <vector>
using namespace std;
using ll=long long;
//geneerated by Codex outside of contest time.
#include <algorithm>
#include <cassert>
#include <cstddef>
#include <vector>
#include <atcoder/convolution>
#include <atcoder/modint>
template <class Field = atcoder::modint998244353>
class FPS {
public:
using value_type = Field;
using size_type = std::size_t;
explicit FPS(size_type n = 0) : coefficients_(n) {}
// 長さ n の零形式的冪級数を作る。
template <class U>
explicit FPS(size_type n = 0, U v = 0) : coefficients_(n,v) {}
// Field に変換できる型 U の係数列から作る。
template <class U>
explicit FPS(const std::vector<U>& coefficients) {
coefficients_.reserve(coefficients.size());
for (const auto& coefficient : coefficients) {
coefficients_.emplace_back(coefficient);
}
}
size_type size() const noexcept {
return coefficients_.size();
}
void resize(size_type n) {
coefficients_.resize(n);
}
Field& operator[](size_type i) {
return coefficients_[i];
}
const Field& operator[](size_type i) const {
return coefficients_[i];
}
FPS& operator+=(const FPS& rhs) {
const size_type old_size = size();
resize(std::max(size(), rhs.size()));
for (size_type i = 0; i < rhs.size(); ++i) {
if (i < old_size) {
coefficients_[i] += rhs[i];
} else {
coefficients_[i] = rhs[i];
}
}
return *this;
}
FPS& operator-=(const FPS& rhs) {
const size_type old_size = size();
resize(std::max(size(), rhs.size()));
for (size_type i = 0; i < rhs.size(); ++i) {
if (i < old_size) {
coefficients_[i] -= rhs[i];
} else {
coefficients_[i] = -rhs[i];
}
}
return *this;
}
FPS& operator*=(const FPS& rhs) {
return *this = *this * rhs;
}
// 形式的微分を求める。
FPS derivative() const {
if (size() == 0) {
return FPS();
}
FPS result(size() - 1);
for (size_type i = 1; i < size(); ++i) {
result[i - 1] = coefficients_[i] * Field(i);
}
return result;
}
// 競技プログラミングで一般的な名前も用意する。
FPS diff() const {
return derivative();
}
// 定数項を 0 として形式的積分を求める。
FPS integral() const {
FPS result(size() + 1);
for (size_type i = 0; i < size(); ++i) {
result[i + 1] = coefficients_[i] / Field(i + 1);
}
return result;
}
// this * result == 1 (mod x^size()) となる乗法逆元を求める。
FPS inv() const {
assert(size() > 0 && coefficients_[0] != Field(0));
FPS result(std::vector<Field>{coefficients_[0].inv()});
while (result.size() < size()) {
const size_type next_size = std::min(size(), result.size() * 2);
FPS prefix(next_size);
for (size_type i = 0; i < next_size; ++i) {
prefix[i] = coefficients_[i];
}
FPS correction = prefix * result;
correction.resize(next_size);
for (auto& coefficient : correction.coefficients_) {
coefficient = -coefficient;
}
correction[0] += Field(2);
result *= correction;
result.resize(next_size);
}
return result;
}
// log(this) (mod x^size()) を O(N log N) で求める。
FPS log() const {
assert(size() > 0 && coefficients_[0] == Field(1));
FPS result = derivative() * inv();
result.resize(size() - 1);
result = result.integral();
result.resize(size());
return result;
}
// exp(this) (mod x^size()) を Newton 法により O(N log N) で求める。
FPS exp() const {
if (size() == 0) {
return FPS();
}
assert(coefficients_[0] == Field(0));
FPS result(std::vector<Field>{Field(1)});
while (result.size() < size()) {
const size_type next_size = std::min(size(), result.size() * 2);
FPS prefix(next_size);
for (size_type i = 0; i < next_size; ++i) {
prefix[i] = coefficients_[i];
}
result.resize(next_size);
FPS correction = prefix - result.log();
correction[0] += Field(1);
result *= correction;
result.resize(next_size);
}
return result;
}
FPS& operator/=(const FPS& rhs) {
return *this = *this / rhs;
}
friend FPS operator+(FPS lhs, const FPS& rhs) {
lhs += rhs;
return lhs;
}
friend FPS operator-(FPS lhs, const FPS& rhs) {
lhs -= rhs;
return lhs;
}
// atcoder::convolution への依存は積の定義だけに閉じ込める。
friend FPS operator*(const FPS& lhs, const FPS& rhs) {
return FPS(atcoder::convolution(lhs.coefficients_, rhs.coefficients_));
}
friend FPS operator/(const FPS& lhs, FPS rhs) {
if (lhs.size() == 0) {
return FPS();
}
rhs.resize(std::max(rhs.size(), lhs.size()));
FPS quotient = lhs * rhs.inv();
quotient.resize(std::max(lhs.size(), rhs.size()));
return quotient;
}
private:
std::vector<Field> coefficients_;
};
// generated by Codex outside of contest time.
#include <algorithm>
#include <cassert>
#include <cstddef>
#include <limits>
#include <vector>
//#include "fps.cpp"
// i = 0, 1, ..., m - 1 に対する [x^k] g(x)f(x)^i を列挙する。
//
// P(x,y) / Q(x,y) = g(x) / (1-yf(x)) に Graeffe 変換を適用する。
// 二変数多項式は配列にはせず、x のブロック幅 stride を用いて
// [x^a y^b] を 1 本の FPS の [x^(a+stride*b)] に格納する。
template <class Field>
std::vector<Field> power_projection(const FPS<Field>& f,
const FPS<Field>& g,
std::size_t m,
std::size_t k) {
if (m == 0) {
return {};
}
// 目標次数を n-1 (n は 2 冪) にそろえる。
std::size_t n = 1;
while (n <= k) {
assert(n <= std::numeric_limits<std::size_t>::max() / 2);
n *= 2;
}
const std::size_t shift = n - 1 - k;
// 各 y ブロックの間隔を 2n にすると、積における x 次数の和が
// 隣の y ブロックへ桁上がりしない。
std::size_t stride = 2 * n;
std::size_t p_y_size = 1;
std::size_t q_y_size = std::min<std::size_t>(m, 2);
FPS<Field> p(n);
for (std::size_t i = 0; i < g.size() && i + shift < n; ++i) {
p[i + shift] = g[i];
}
FPS<Field> q(q_y_size == 1 ? n : stride + n);
q[0] = Field(1);
if (q_y_size == 2) {
for (std::size_t i = 0; i < f.size() && i < n; ++i) {
q[stride + i] = -f[i];
}
}
while (n > 1) {
// q(-x,y)。符号はブロック内の x 次数だけで決まる。
FPS<Field> q_minus = q;
for (std::size_t y = 0; y < q_y_size; ++y) {
const std::size_t base = y * stride;
for (std::size_t x = 1; x < n; x += 2) {
q_minus[base + x] = -q_minus[base + x];
}
}
FPS<Field> pq = p * q_minus;
FPS<Field> qq = q * q_minus;
const std::size_t next_n = n / 2;
const std::size_t next_stride = n;
const std::size_t next_p_y_size =
std::min(m, p_y_size + q_y_size - 1);
const std::size_t next_q_y_size =
std::min(m, q_y_size + q_y_size - 1);
FPS<Field> next_p((next_p_y_size - 1) * next_stride + next_n);
FPS<Field> next_q((next_q_y_size - 1) * next_stride + next_n);
// P(x,y)Q(-x,y) の x 奇数次部分と、
// Q(x,y)Q(-x,y) の x 偶数次部分を抜き出す。
for (std::size_t y = 0; y < next_p_y_size; ++y) {
const std::size_t old_base = y * stride;
const std::size_t new_base = y * next_stride;
for (std::size_t x = 1; x < n; x += 2) {
if (old_base + x < pq.size()) {
next_p[new_base + x / 2] = pq[old_base + x];
}
}
}
for (std::size_t y = 0; y < next_q_y_size; ++y) {
const std::size_t old_base = y * stride;
const std::size_t new_base = y * next_stride;
for (std::size_t x = 0; x < n; x += 2) {
if (old_base + x < qq.size()) {
next_q[new_base + x / 2] = qq[old_base + x];
}
}
}
p = next_p;
q = next_q;
p_y_size = next_p_y_size;
q_y_size = next_q_y_size;
n = next_n;
stride = next_stride;
}
// x の幅が 1 になったので、ブロック先頭を集めて通常の FPS に戻す。
FPS<Field> numerator(m);
FPS<Field> denominator(m);
for (std::size_t y = 0; y < p_y_size; ++y) {
numerator[y] = p[y * stride];
}
for (std::size_t y = 0; y < q_y_size; ++y) {
denominator[y] = q[y * stride];
}
assert(denominator[0] != Field(0));
FPS<Field> answer = numerator * denominator.inv();
answer.resize(m);
std::vector<Field> result(m);
for (std::size_t i = 0; i < m; ++i) {
result[i] = answer[i];
}
return result;
}
// generated by Codex outside of contest time.
#include <algorithm>
#include <cassert>
#include <cstddef>
#include <vector>
//#include "fps.cpp"
namespace multipoint_evaluation_internal {
template <class Field>
FPS<Field> trim(FPS<Field> f) {
while (f.size() > 0 && f[f.size() - 1] == Field(0)) {
f.resize(f.size() - 1);
}
return f;
}
// a mod b を、反転した FPS の除算を用いて求める。
// このファイル内で用いる b はすべてモニックである。
template <class Field>
FPS<Field> polynomial_remainder(FPS<Field> a, const FPS<Field>& b) {
a = trim(a);
if (a.size() < b.size()) {
return a;
}
const std::size_t quotient_size = a.size() - b.size() + 1;
FPS<Field> reversed_a(quotient_size);
FPS<Field> reversed_b(quotient_size);
for (std::size_t i = 0; i < quotient_size; ++i) {
reversed_a[i] = a[a.size() - 1 - i];
if (i < b.size()) {
reversed_b[i] = b[b.size() - 1 - i];
}
}
FPS<Field> reversed_quotient = reversed_a * reversed_b.inv();
reversed_quotient.resize(quotient_size);
FPS<Field> quotient(quotient_size);
for (std::size_t i = 0; i < quotient_size; ++i) {
quotient[i] = reversed_quotient[quotient_size - 1 - i];
}
FPS<Field> remainder = a - quotient * b;
if (remainder.size() >= b.size()) {
remainder.resize(b.size() - 1);
}
return trim(remainder);
}
} // namespace multipoint_evaluation_internal
// 任意の相異なるとは限らない点 x[i] で f を評価する。
// 積木と剰余計算により O(M(N) log N) = O(N log^2 N)。
template <class Field, class T>
std::vector<Field> multipoint_evaluation(const FPS<Field>& f,
const std::vector<T>& x) {
const std::size_t point_count = x.size();
if (point_count == 0) {
return {};
}
std::size_t leaf_count = 1;
while (leaf_count < point_count) {
leaf_count *= 2;
}
std::vector<FPS<Field>> product_tree(2 * leaf_count);
for (std::size_t i = 0; i < leaf_count; ++i) {
if (i < point_count) {
product_tree[leaf_count + i] =
FPS<Field>(std::vector<Field>{-Field(x[i]), Field(1)});
} else {
product_tree[leaf_count + i] =
FPS<Field>(std::vector<Field>{Field(1)});
}
}
for (std::size_t i = leaf_count; i-- > 1;) {
product_tree[i] = product_tree[2 * i] * product_tree[2 * i + 1];
}
std::vector<FPS<Field>> remainder_tree(2 * leaf_count);
remainder_tree[1] =
multipoint_evaluation_internal::polynomial_remainder(f, product_tree[1]);
for (std::size_t i = 1; i < leaf_count; ++i) {
remainder_tree[2 * i] = multipoint_evaluation_internal::polynomial_remainder(
remainder_tree[i], product_tree[2 * i]);
remainder_tree[2 * i + 1] =
multipoint_evaluation_internal::polynomial_remainder(
remainder_tree[i], product_tree[2 * i + 1]);
}
std::vector<Field> result(point_count, Field(0));
for (std::size_t i = 0; i < point_count; ++i) {
const FPS<Field>& value = remainder_tree[leaf_count + i];
if (value.size() > 0) {
result[i] = value[0];
}
}
return result;
}
// x[i] = a r^i であることを仮定して f(x[i]) を列挙する。
// r^(ij) = r^C(i+j,2) / (r^C(i,2) r^C(j,2)) を用いる chirp 変換で、
// 1 回の FPS の積に帰着する。計算量は O(M(N)) = O(N log N)。
template <class Field, class T>
std::vector<Field> multipoint_evaluation_geometric(
const FPS<Field>& f, const std::vector<T>& x) {
const std::size_t point_count = x.size();
std::vector<Field> result(point_count, Field(0));
if (point_count == 0) {
return result;
}
const Field constant = f.size() == 0 ? Field(0) : f[0];
const Field first = Field(x[0]);
if (first == Field(0)) {
// 等比数列で初項が 0 なら全評価点が 0。
std::fill(result.begin(), result.end(), constant);
return result;
}
if (point_count == 1) {
Field power = Field(1);
for (std::size_t j = 0; j < f.size(); ++j) {
result[0] += f[j] * power;
power *= first;
}
return result;
}
const Field ratio = Field(x[1]) / first;
if (ratio == Field(0)) {
// a, 0, 0, ... の場合。逆元を使う chirp 変換から外して扱う。
Field power = Field(1);
for (std::size_t j = 0; j < f.size(); ++j) {
result[0] += f[j] * power;
power *= first;
}
std::fill(result.begin() + 1, result.end(), constant);
return result;
}
const std::size_t degree_size = f.size();
if (degree_size == 0) {
return result;
}
// chirp[t] = r^{t(t-1)/2} を乗算だけで列挙する。
const std::size_t chirp_size = point_count + degree_size - 1;
std::vector<Field> chirp(chirp_size, Field(1));
Field ratio_power = Field(1);
for (std::size_t t = 1; t < chirp_size; ++t) {
chirp[t] = chirp[t - 1] * ratio_power;
ratio_power *= ratio;
}
FPS<Field> left(degree_size);
Field first_power = Field(1);
Field inverse_chirp = Field(1);
const Field inverse_ratio = ratio.inv();
Field inverse_ratio_power = Field(1);
for (std::size_t j = 0; j < degree_size; ++j) {
left[degree_size - 1 - j] = f[j] * first_power * inverse_chirp;
first_power *= first;
inverse_chirp *= inverse_ratio_power;
inverse_ratio_power *= inverse_ratio;
}
FPS<Field> product = left * FPS<Field>(chirp);
inverse_chirp = Field(1);
inverse_ratio_power = Field(1);
for (std::size_t i = 0; i < point_count; ++i) {
result[i] = product[degree_size - 1 + i] * inverse_chirp;
inverse_chirp *= inverse_ratio_power;
inverse_ratio_power *= inverse_ratio;
}
return result;
}
int main(){
cin.tie(nullptr);
ios::sync_with_stdio(false);
int n,m;
cin>>n>>m;
FPS f(n),g(n),h(n);
for(int i=0;i<n;i++){
int t;
cin>>t;
f[i]=t;
}
for(int i=0;i<n;i++){
int t;
cin>>t;
g[i]=t;
}
for(int i=0;i<n;i++){
int t;
cin>>t;
h[i]=t;
}
auto v=power_projection(g,FPS(1,1),n,n-1);
FPS r(n);
using mint=atcoder::modint998244353;
for(int i=0;i<n;i++)r[i]=h[i]*v[i];
vector<mint> x(m,1);
for(int i=1;i<m;i++)x[i]=x[i-1]*f[1];
auto ans=multipoint_evaluation_geometric(r,x);
for(int i=0;i<m;i++)cout<<ans[i].val()<<(i+1==m?"\n":" ");
}