#1197581 (C++23) No.3793 Multishift

提出ソース
結果

問題	No.3793 Multishift
コンテスト
ユーザー	kwm_t
提出日時	2026-11-07 14:19:40
言語	C++23 (gcc 15.3.0 + boost 1.92.0 + ACL) コンパイル: `g++-15 -O2 -fplugin=/usr/lib/yukicoder/noavx512-gcc15.so -lm -std=c++23 -Wuninitialized -DONLINE_JUDGE -o a.out _filename_` 実行: `./a.out`
結果	AC
実行時間	966 ms / 7,000 ms
	+ 629µs
コード長	22,366 bytes
記録記録タグの例: 初AC ショートコード純ショートコード純主流ショートコード最速実行時間
コンパイル時間	6,110 ms
コンパイル使用メモリ	405,076 KB
実行使用メモリ	74,004 KB
最終ジャッジ日時	2026-11-07 14:20:03
合計ジャッジ時間	16,274 ms
ジャッジサーバーID （参考情報）	judge3_0 / judge7_0
このコードへのチャレンジ
（要ログイン）
ファイルパターン	結果
sample	AC * 1
other	AC * 13
権限があれば一括ダウンロードができます
ソースコード

raw source code
#include <bits/stdc++.h>
#include <atcoder/all>
using namespace std;
using namespace atcoder;
// using mint = modint1000000007;
// const int mod = 1000000007;
using mint = modint998244353;
const int mod = 998244353;
// const int INF = 1e9;
// const long long LINF = 1e18;
#define rep(i, n) for (int i = 0; i < (n); ++i)
#define rep2(i, l, r) for (int i = (l); i < (r); ++i)
#define rrep(i, n) for (int i = (n)-1; i >= 0; --i)
#define rrep2(i, l, r) for (int i = (r)-1; i >= (l); --i)
#define all(x) (x).begin(), (x).end()
#define allR(x) (x).rbegin(), (x).rend()
//#define P pair<int, int>
template<typename A, typename B> inline bool chmax(A& a, const B& b) { if (a < b) { a = b; return true; } return false; }
template<typename A, typename B> inline bool chmin(A& a, const B& b) { if (a > b) { a = b; return true; } return false; }

/**
 * @file factorial.hpp
 * @brief 形式的冪級数の階乗関連処理
 * @author kwm_t
 * @see https://github.com/kwm-t/cp-library/blob/main/math/fps/factorial.hpp
 */
#ifndef KWM_T_MATH_FPS_FACTORIAL_HPP
#define KWM_T_MATH_FPS_FACTORIAL_HPP

#include <vector>

namespace kwm_t::math::fps {

/**
 * @brief 階乗・逆階乗テーブル
 *
 * fac[i]  = i!
 * finv[i] = (i!)^{-1}
 *
 * @tparam mint modint型
 */
template<typename mint>
struct Factorial {
	int n;
	std::vector<mint> fac, finv;

	explicit Factorial(int n = 0) : n(n), fac(n + 1, 1), finv(n + 1, 1) {
		if (n == 0) return;
		for (int i = 2; i <= n; ++i) fac[i] = fac[i - 1] * i;
		finv[n] = fac[n].inv();
		for (int i = n; i >= 1; --i) finv[i - 1] = finv[i] * i;
	}

	// nCk
	mint comb(int n, int k) const {
		if (k < 0 || n < k) return 0;
		return fac[n] * finv[k] * finv[n - k];
	}
	mint binom(int n, int k) const {
		return comb(n, k);
	}

	// nPk
	mint perm(int n, int k) const {
		if (k < 0 || n < k) return 0;
		return fac[n] * finv[n - k];
	}
};

} // namespace kwm_t::math::fps

#endif // KWM_T_MATH_FPS_FACTORIAL_HPP

/**
 * @file fps.hpp
 * @brief 形式的冪級数（FPS）
 * @author kwm_t
 * @see https://github.com/kwm-t/cp-library/blob/main/math/fps/fps.hpp
 */
#ifndef KWM_T_MATH_FPS_FPS_HPP
#define KWM_T_MATH_FPS_FPS_HPP

#include <vector>
#include <algorithm>
 // #include "factorial.hpp"
#include "atcoder/convolution"

namespace kwm_t::math::fps {


template<class T>
struct FormalPowerSeries : std::vector<T> {
	FormalPowerSeries(const std::vector<T>& vec) : std::vector<T>(vec) {}

	using std::vector<T>::vector;
	using std::vector<T>::operator=;
	using F = FormalPowerSeries;

	F operator-() const {
		F res(*this);
		for (auto& e : res) e = -e;
		return res;
	}
	F& operator*=(const T& g) {
		for (auto& e : *this) e *= g;
		return *this;
	}
	F& operator/=(const T& g) {
		assert(g != T(0));
		*this *= g.inv();
		return *this;
	}
	F& operator+=(const F& g) {
		int n = this->size(), m = g.size();
		for (int i = 0; i < std::min(n, m); ++i) (*this)[i] += g[i];
		return *this;
	}
	F& operator-=(const F& g) {
		int n = this->size(), m = g.size();
		for (int i = 0; i < std::min(n, m); ++i)(*this)[i] -= g[i];
		return *this;
	}
	F& operator<<=(const int d) {
		int n = this->size();
		if (d >= n) *this = F(n);
		this->insert(this->begin(), d, 0);
		this->resize(n);
		return *this;
	}
	F& operator>>=(const int d) {
		int n = this->size();
		this->erase(this->begin(), this->begin() + min(n, d));
		this->resize(n);
		return *this;
	}

	// O(n log n)
	F& inv_inplace(int d = -1) {
		int n = this->size();
		assert(n != 0 && (*this)[0] != 0);
		if (d == -1) d = n;
		assert(d >= 0);

		F res{ (*this)[0].inv() };
		for (int m = 1; m < d; m *= 2) {
			F f(this->begin(), this->begin() + std::min(n, 2 * m));
			F g(res);
			f.resize(2 * m), atcoder::internal::butterfly(f);
			g.resize(2 * m), atcoder::internal::butterfly(g);
			for (int i = 0; i < 2 * m; ++i) f[i] *= g[i];
			atcoder::internal::butterfly_inv(f);
			f.erase(f.begin(), f.begin() + m);
			f.resize(2 * m), atcoder::internal::butterfly(f);
			for (int i = 0; i < 2 * m; ++i) f[i] *= g[i];
			atcoder::internal::butterfly_inv(f);
			T iz = T(2 * m).inv();
			iz *= -iz;
			for (int i = 0; i < m; ++i) f[i] *= iz;
			res.insert(res.end(), f.begin(), f.begin() + m);
		}
		res.resize(d);
		*this = std::move(res);
		return *this;
	}

	F inv(int d = -1) const { return F(*this).inv_inplace(d); }

	F& add_inplace(const F& g, int d = -1) {
		int n = this->size();
		if (d == -1) d = n;
		assert(d >= 0);
		this->resize(d);
		for (int i = 0; i < std::min(d, (int)g.size()); ++i)
			(*this)[i] += g[i];
		return *this;
	}
	F& add_resize_inplace(const F& g) { return add_inplace(g, std::max((int)this->size(), (int)g.size())); }
	F add(const F& g, int d = -1) const { return F(*this).add_inplace(g, d); }
	F add_resize(const F& g) const { return F(*this).add_resize_inplace(g); }

	F& sub_inplace(const F& g, int d = -1) {
		int n = this->size();
		if (d == -1) d = n;
		assert(d >= 0);
		this->resize(d);
		for (int i = 0; i < std::min(d, (int)g.size()); ++i)
			(*this)[i] -= g[i];
		return *this;
	}
	F& sub_resize_inplace(const F& g) { return sub_inplace(g, std::max((int)this->size(), (int)g.size())); }
	F sub(const F& g, int d = -1) const { return F(*this).sub_inplace(g, d); }
	F sub_resize(const F& g) const { return F(*this).sub_resize_inplace(g); }

	// fast: FMT-friendly modulus only
	// O(n log n)
	F& multiply_inplace(const F& g, int d = -1) {
		int n = this->size();
		if (d == -1) d = n;
		assert(d >= 0);
		*this = convolution(move(*this), g);
		this->resize(d);
		return *this;
	}
	F& multiply_resize_inplace(const F& g) { return multiply_inplace(g, (int)this->size() + (int)g.size() - 1); }
	F multiply(const F& g, const int d = -1) const { return F(*this).multiply_inplace(g, d); }
	F multiply_resize(const F& g) const { return F(*this).multiply_resize_inplace(g); }

	// O(n log n)
	F& divide_inplace(const F& g, int d = -1) {
		int n = this->size();
		if (d == -1) d = n;
		assert(d >= 0);
		*this = convolution(move(*this), g.inv(d));
		this->resize(d);
		return *this;
	}
	F& divide_resize_inplace(const F& g) { return divide_inplace(g, (int)this->size() + (int)g.size() - 1); }
	F divide(const F& g, const int d = -1) const { return F(*this).divide_inplace(g, d); }
	F divide_resize(const F& g) const { return F(*this).divide_resize_inplace(g); }

	// sparse
	// O(nk)
	F& multiply_inplace(std::vector<std::pair<int, T>> g) {
		int n = this->size();
		auto [d, c] = g.front();
		if (d == 0) g.erase(g.begin());
		else c = 0;

		for (int i = (n - 1); i >= 0; --i) {
			(*this)[i] *= c;
			for (auto& [j, b] : g) {
				if (j > i) break;
				(*this)[i] += (*this)[i - j] * b;
			}
		}
		return *this;
	}
	F multiply(const std::vector<std::pair<int, T>>& g) const { return F(*this).multiply_inplace(g); }
	// O(nk)
	F& divide_inplace(std::vector<std::pair<int, T>> g) {
		int n = this->size();
		auto [d, c] = g.front();
		assert(d == 0 && c != T(0));
		T ic = c.inv();
		g.erase(g.begin());
		for (int i = 0; i < n; ++i) {
			for (auto& [j, b] : g) {
				if (j > i) break;
				(*this)[i] -= (*this)[i - j] * b;
			}
			(*this)[i] *= ic;
		}
		return *this;
	}
	F divide(const std::vector<std::pair<int, T>>& g) const { return F(*this).divide_inplace(g); }

	// divide by x^shift * h(x), where h(0) != 0
	// return the removed shift amount
	// O(nk)
	int divide_shifted_inplace(std::vector<std::pair<int, T>> g) {
		int n = this->size();
		assert(!g.empty());

		int shift = g.front().first;

		// remove x^shift
		*this >>= shift;

		for (auto& [d, c] : g) d -= shift;

		assert(g.front().first == 0);
		assert(g.front().second != T(0));

		divide_inplace(std::move(g));

		return shift;
	}

	F divide_shifted(const std::vector<std::pair<int, T>>& g) const { return F(*this).divide_shifted_inplace(g); }

	// multiply and divide (1 + cz^d)
	// O(n)
	void multiply_inplace(const int d, const T c) {
		int n = this->size();
		if (c == T(1)) for (int i = (n - d - 1); i >= 0; --i) (*this)[i + d] += (*this)[i];
		else if (c == T(-1)) for (int i = (n - d - 1); i >= 0; --i) (*this)[i + d] -= (*this)[i];
		else for (int i = (n - d - 1); i >= 0; --i) (*this)[i + d] += (*this)[i] * c;
	}
	// O(n)
	void divide_inplace(const int d, const T c) {
		int n = this->size();
		if (c == T(1)) for (int i = 0; i < n - d; ++i) (*this)[i + d] -= (*this)[i];
		else if (c == T(-1)) for (int i = 0; i < n - d; ++i) (*this)[i + d] += (*this)[i];
		else for (int i = 0; i < n - d; ++i) (*this)[i + d] -= (*this)[i] * c;
	}

	// O(n)
	T eval(const T& a) const {
		T x(1), res(0);
		for (auto e : *this) res += e * x, x *= a;
		return res;
	}

	// O(n)
	F& integ_inplace() {
		int n = this->size();
		assert(n > 0);
		if (n == 1) return *this = F{ 0 };
		this->insert(this->begin(), 0);
		this->pop_back();
		std::vector<T> inv(n);
		inv[1] = 1;
		int p = T::mod();
		for (int i = 2; i < n; ++i) inv[i] = -inv[p % i] * (p / i);
		for (int i = 2; i < n; ++i) (*this)[i] *= inv[i];
		return *this;
	}
	F integ() const { return F(*this).integ_inplace(); }

	// O(n)
	F& deriv_inplace() {
		int n = this->size();
		assert(n > 0);
		for (int i = 2; i < n; ++i) (*this)[i] *= i;
		this->erase(this->begin());
		this->push_back(0);
		return *this;
	}
	F deriv() const { return F(*this).deriv_inplace(); }

	// O(n log n)
	F& log_inplace(int d = -1) {
		int n = this->size();
		assert(n > 0 && (*this)[0] == 1);
		if (d == -1) d = n;
		assert(d >= 0);
		this->resize(d);
		F f_inv = this->inv();
		this->deriv_inplace();
		this->multiply_inplace(f_inv);
		this->integ_inplace();
		return *this;
	}
	F log(const int d = -1) const { return F(*this).log_inplace(d); }

	// O(n log n)
	// https://arxiv.org/abs/1301.5804 (Figure 1, right)
	F& exp_inplace(int d = -1) {
		int n = this->size();
		assert(n > 0 && (*this)[0] == 0);
		if (d == -1) d = n;
		assert(d >= 0);
		F g{ 1 }, g_fft{ 1, 1 };
		(*this)[0] = 1;
		this->resize(d);
		F h_drv(this->deriv());
		for (int m = 2; m < d; m *= 2) {
			// prepare
			F f_fft(this->begin(), this->begin() + m);
			f_fft.resize(2 * m), atcoder::internal::butterfly(f_fft);

			// Step 2.a'
			// {
			F _g(m);
			for (int i = 0; i < m; ++i) _g[i] = f_fft[i] * g_fft[i];
			atcoder::internal::butterfly_inv(_g);
			_g.erase(_g.begin(), _g.begin() + m / 2);
			_g.resize(m), atcoder::internal::butterfly(_g);
			for (int i = 0; i < m; ++i) _g[i] *= g_fft[i];
			atcoder::internal::butterfly_inv(_g);
			_g.resize(m / 2);
			_g /= T(-m) * m;
			g.insert(g.end(), _g.begin(), _g.begin() + m / 2);
			// }

			// Step 2.b'--d'
			F t(this->begin(), this->begin() + m);
			t.deriv_inplace();
			// {
			  // Step 2.b'
			F r{ h_drv.begin(), h_drv.begin() + m - 1 };
			// Step 2.c'
			r.resize(m); atcoder::internal::butterfly(r);
			for (int i = 0; i < m; ++i) r[i] *= f_fft[i];
			atcoder::internal::butterfly_inv(r);
			r /= -m;
			// Step 2.d'
			t += r;
			t.insert(t.begin(), t.back()); t.pop_back();
			// }

			// Step 2.e'
			if (2 * m < d) {
				t.resize(2 * m); atcoder::internal::butterfly(t);
				g_fft = g; g_fft.resize(2 * m); atcoder::internal::butterfly(g_fft);
				for (int i = 0; i < 2 * m; ++i) t[i] *= g_fft[i];
				atcoder::internal::butterfly_inv(t);
				t.resize(m);
				t /= 2 * m;
			}
			else { // この場合分けをしても数パーセントしか速くならない
				F g1(g.begin() + m / 2, g.end());
				F s1(t.begin() + m / 2, t.end());
				t.resize(m / 2);
				g1.resize(m), atcoder::internal::butterfly(g1);
				t.resize(m), atcoder::internal::butterfly(t);
				s1.resize(m), atcoder::internal::butterfly(s1);
				for (int i = 0; i < m; ++i) s1[i] = g_fft[i] * s1[i] + g1[i] * t[i];
				for (int i = 0; i < m; ++i) t[i] *= g_fft[i];
				atcoder::internal::butterfly_inv(t);
				atcoder::internal::butterfly_inv(s1);
				for (int i = 0; i < m / 2; ++i) t[i + m / 2] += s1[i];
				t /= m;
			}

			// Step 2.f'
			F v(this->begin() + m, this->begin() + std::min<int>(d, 2 * m)); v.resize(m);
			t.insert(t.begin(), m - 1, 0); t.push_back(0);
			t.integ_inplace();
			for (int i = 0; i < m; ++i) v[i] -= t[m + i];

			// Step 2.g'
			v.resize(2 * m); atcoder::internal::butterfly(v);
			for (int i = 0; i < 2 * m; ++i) v[i] *= f_fft[i];
			atcoder::internal::butterfly_inv(v);
			v.resize(m);
			v /= 2 * m;

			// Step 2.h'
			for (int i = 0; i < std::min(d - m, m); ++i)(*this)[m + i] = v[i];
		}
		return *this;
	}
	F exp(const int d = -1) const { return F(*this).exp_inplace(d); }

	// O(n log n)
	F& pow_inplace(const long long k, int d = -1) {
		int n = this->size();
		if (d == -1) d = n;
		assert(d >= 0 && k >= 0);
		if (k == 0) {
			*this = F(d);
			if (d > 0) (*this)[0] = 1;
			return *this;
		}
		int l = 0;
		while (l < n && (*this)[l] == 0) ++l;
		if (l > (d - 1) / k || l == n) return *this = F(d);
		T c = (*this)[l];
		this->erase(this->begin(), this->begin() + l);
		*this /= c;
		this->log_inplace(d - l * k);
		*this *= k;
		this->exp_inplace();
		*this *= c.pow(k);
		this->insert(this->begin(), l * k, 0);
		return *this;
	}
	F pow(const long long k, const int d = -1) const { return F(*this).pow_inplace(k, d); }

	// O(n log n)
	F& shift_inplace(const T c) {
		int n = this->size();
		auto fc = Factorial<T>(n);
		for (int i = 0; i < n; ++i) (*this)[i] *= fc.fac[i];
		reverse(this->begin(), this->end());
		F g(n);
		T cp = 1;
		for (int i = 0; i < n; ++i) g[i] = cp * fc.finv[i], cp *= c;
		this->multiply_inplace(g, n);
		reverse(this->begin(), this->end());
		for (int i = 0; i < n; ++i) (*this)[i] *= fc.finv[i];
		return *this;
	}
	F shift(const T c) const { return F(*this).shift_inplace(c); }

	F operator*(const T& g) const { return F(*this) *= g; }
	F operator/(const T& g) const { return F(*this) /= g; }
	F operator+(const F& g) const { return F(*this) += g; }
	F operator-(const F& g) const { return F(*this) -= g; }
	F operator<<(const int d) const { return F(*this) <<= d; }
	F operator>>(const int d) const { return F(*this) >>= d; }
	F operator*(std::vector<std::pair<int, T>> g) const { return F(*this) *= g; }
	F operator/(std::vector<std::pair<int, T>> g) const { return F(*this) /= g; }
};


} // namespace kwm_t::math::fps

#endif // KWM_T_MATH_FPS_FPS_HPP

/**
 * @file builtin.hpp
 * @brief コンパイラ組み込み機能の補助
 * @author kwm_t
 * @see https://github.com/kwm-t/cp-library/blob/main/base/builtin.hpp
 */
#ifndef KWM_T_BASE_BUILTIN_HPP
#define KWM_T_BASE_BUILTIN_HPP

 /**
  * @brief ローカル環境用(MSC)builtin
  */
#ifdef _MSC_VER
#include <intrin.h>
#define __builtin_popcount (int)__popcnt
#define __builtin_popcountll (int)__popcnt64
inline int __builtin_ctz(unsigned int n) { unsigned long i; _BitScanForward(&i, n); return i; }
inline int __builtin_ctzll(unsigned long long n) { unsigned long i; _BitScanForward64(&i, n); return i; }
inline int __builtin_clz(unsigned int n) { unsigned long i; if (_BitScanReverse(&i, n)) { return 31 - i; }return 32; }
inline int __builtin_clzll(unsigned long long n) { unsigned long i; if (_BitScanReverse64(&i, n)) { return 63 - i; }return 64; }
#endif

#endif // KWM_T_BASE_BUILTIN_HPP

/**
 * @file multipoint_evaluation.hpp
 * @brief 多点評価
 * @author kwm_t
 * @see https://github.com/kwm-t/cp-library/blob/main/math/fps/multipoint_evaluation.hpp
 */
#ifndef KWM_T_MATH_FPS_MULTIPOINT_EVALUATION_HPP
#define KWM_T_MATH_FPS_MULTIPOINT_EVALUATION_HPP

#include <vector>
#include <algorithm>
#include <atcoder/convolution>

 // #include "base/builtin.hpp"

namespace kwm_t::math::fps {

namespace detail {

// 複数FPSの積（divide & conquer）
template <class T>
std::vector<T> product_sequence(const std::vector<std::vector<T>>& polys) {
	if (polys.empty()) return { T(1) };

	auto dfs = [&](auto&& self, int l, int r) -> std::vector<T> {
		if (l + 1 == r) return polys[l];
		int m = (l + r) >> 1;
		return atcoder::convolution(
			self(self, l, m),
			self(self, m, r)
		);
	};

	return dfs(dfs, 0, (int)polys.size());
}


// 区間畳み込み（必要な範囲だけ取る）
template <class T>
std::vector<T> convolution_range(std::vector<T> a,
	std::vector<T> b,
	int l, int r) {
	if (l >= r) return {};

	if (l < 0) l = 0;

	while ((int)a.size() > r) a.pop_back();
	while ((int)b.size() > r) b.pop_back();

	int need = (int)a.size() + (int)b.size() - 1;
	need -= std::min(need - r, l);

	int z = 1;
	while (z < need) z <<= 1;

	a.resize(z, T(0));
	b.resize(z, T(0));

	atcoder::internal::butterfly(a);
	atcoder::internal::butterfly(b);

	for (int i = 0; i < z; i++) a[i] *= b[i];

	atcoder::internal::butterfly_inv(a);

	static std::vector<T> inv = { T(1), T(1) / T(2) };
	while ((int)inv.size() <= __builtin_ctz(z)) {
		inv.push_back(inv.back() * inv[1]);
	}

	std::vector<T> res(r - l);
	T coef = inv[__builtin_ctz(z)];
	for (int i = l; i < r; i++) {
		res[i - l] = a[i] * coef;
	}

	return res;
}

} // namespace detail


/**
 * @brief 多点評価
 *
 * f(x) を p[i] で評価する
 *
 * @tparam T modint
 * @param f 多項式
 * @param p 評価点
 * @return std::vector<T> f(p[i])
 *
 * 計算量:
 *   O((n + m) log^2 n)
 *
 * verified
 *   https://atcoder.jp/contests/fps-24/submissions/74606689
 */
template <class T>
std::vector<T> multipoint_evaluation(std::vector<T> f,
	const std::vector<T>& p) {
	int N = (int)p.size();
	int M = (int)f.size();

	// reverse
	std::reverse(f.begin(), f.end());

	// g[i] = (x - p[i])
	std::vector<std::vector<T>> leaf(N);
	for (int i = 0; i < N; i++) {
		leaf[i] = { T(1), -p[i] };
	}

	// 全体積
	auto prod = detail::product_sequence(leaf);

	// 逆元
	auto inv = (FormalPowerSeries<T>(prod)).inv(M);

	// f / prod
	auto tmp = atcoder::convolution(inv, f);
	tmp.resize(M);

	// セグ木構築
	int size = 1;
	while (size < N) size <<= 1;

	std::vector<std::vector<T>> seg(2 * size, { T(1) });
	for (int i = 0; i < N; i++) {
		seg[size + i] = leaf[i];
	}
	for (int i = size - 1; i > 0; i--) {
		seg[i] = atcoder::convolution(seg[i << 1], seg[i << 1 | 1]);
	}

	std::vector<T> res(N);

	// remainder tree
	auto dfs = [&](auto&& self, int node, std::vector<T> cur) -> void {
		if ((int)seg[node].size() == 1) return;

		if (node >= size) {
			res[node - size] = cur.back();
			return;
		}

		const auto& L = seg[node << 1];
		const auto& R = seg[node << 1 | 1];

		self(self, node << 1,
			detail::convolution_range(
				cur, R,
				(int)cur.size() - (int)L.size() + 1,
				(int)cur.size()
			));

		self(self, node << 1 | 1,
			detail::convolution_range(
				cur, L,
				(int)cur.size() - (int)R.size() + 1,
				(int)cur.size()
			));
	};

	dfs(dfs, 1, tmp);

	return res;
}

} // namespace kwm_t::math::fps

#endif // KWM_T_MATH_FPS_MULTIPOINT_EVALUATION_HPP

/**
 * @file lagrange.hpp
 * @brief Lagrange 補間
 * @author kwm_t
 * @see https://github.com/kwm-t/cp-library/blob/main/math/modint/lagrange.hpp
 */
#ifndef KMW_T_MATH_MODINT_LAGRANGE_HPP
#define KMW_T_MATH_MODINT_LAGRANGE_HPP

#include <vector>

 /**
  * @brief ラグランジュ補間（0,1,...,Nでの値から f(x) を求める）
  *
  * 典型用途:
  *   f(0), f(1), ..., f(N) が分かっているときに f(x) を求める
  *
  * 計算量:
  *   O(N)
  *
  * @tparam Mint modint型
  *
  * 制約 / 注意:
  *   - mod は素数
  *   - x は long long
  *   - y.size() >= 1
  * verified
  *   https://atcoder.jp/contests/arc230/submissions/79418709
  */

namespace kmw_t::math::modint {

template <class Mint>
Mint lagrange_interpolation(const std::vector<Mint>& y, long long x) {
	int N = (int)y.size() - 1;
	if (x <= N) return y[x];

	std::vector<Mint> invfac(N + 1, 1);
	Mint fac = 1;
	for (int i = 1; i <= N; i++) fac *= i;

	invfac[N] = fac.inv();
	for (int i = N; i > 0; i--) {
		invfac[i - 1] = invfac[i] * i;
	}

	std::vector<Mint> dp(N + 1, 1), pd(N + 1, 1);
	for (int i = 0; i < N; i++) {
		dp[i + 1] = dp[i] * (x - i);
	}
	for (int i = N; i > 0; i--) {
		pd[i - 1] = pd[i] * (x - i);
	}

	Mint res = 0;
	for (int i = 0; i <= N; i++) {
		Mint term = y[i] * dp[i] * pd[i] * invfac[i] * invfac[N - i];
		if ((N - i) & 1) res -= term;
		else res += term;
	}

	return res;
}

/**
 * @brief Lagrange 補間による多項式の係数復元
 *
 * y[i] = f(i) (0 <= i <= N) から
 * f(x) = c[0] + c[1]x + ... + c[N]x^N
 * の係数を復元する
 *
 * 計算量:
 *   O(N log^2 N)
 *
 * @tparam Mint modint型
 *
 * 制約 / 注意:
 *   - mod は素数
 *   - N < mod
 *   - y.size() >= 1
 *
 * @return c[i] = x^i の係数
 */
template <class Mint>
std::vector<Mint> lagrange_coefficients(const std::vector<Mint>& y) {
	int N = (int)y.size();
	if (N == 0) return {};

	std::vector<std::vector<Mint>> p(4 * N);

	auto build = [&](auto&& self, int node, int l, int r) -> void {
		if (r - l == 1) {
			p[node] = { -Mint(l), Mint(1) };
			return;
		}
		int mid = (l + r) / 2;
		self(self, node * 2, l, mid);
		self(self, node * 2 + 1, mid, r);
		p[node] = atcoder::convolution(p[node * 2], p[node * 2 + 1]);
	};
	build(build, 1, 0, N);

	std::vector<Mint> invfac(N, 1);
	Mint fac = 1;
	for (int i = 1; i < N; i++) fac *= i;
	invfac[N - 1] = fac.inv();
	for (int i = N - 1; i > 0; i--) {
		invfac[i - 1] = invfac[i] * i;
	}

	std::vector<Mint> res;

	auto combine = [&](auto&& self, int node, int l, int r) -> std::vector<Mint> {
		if (r - l == 1) {
			Mint c = y[l] * invfac[l] * invfac[N - 1 - l];
			if ((N - 1 - l) & 1) c = -c;
			return { c };
		}

		int mid = (l + r) / 2;
		std::vector<Mint> a = self(self, node * 2, l, mid);
		std::vector<Mint> b = self(self, node * 2 + 1, mid, r);

		std::vector<Mint> c1 = atcoder::convolution(a, p[node * 2 + 1]);
		std::vector<Mint> c2 = atcoder::convolution(b, p[node * 2]);

		std::vector<Mint> res(std::max(c1.size(), c2.size()), 0);
		for (int i = 0; i < (int)c1.size(); i++) res[i] += c1[i];
		for (int i = 0; i < (int)c2.size(); i++) res[i] += c2[i];
		return res;
	};

	res = combine(combine, 1, 0, N);
	res.resize(N);
	return res;
}
} // namespace kmw_t::math::modint

#endif // KMW_T_MATH_MODINT_LAGRANGE_HPP

int main() {
	std::ios::sync_with_stdio(false);
	std::cin.tie(nullptr);
	int n, m; cin >> n >> m;
	vector<mint>a(n + 1), b(m + 1);
	rep(i, n + 1) {
		int x; cin >> x;
		a[i] = x;
	}
	rep(i, m + 1) {
		int x; cin >> x;
		b[i] = x;
	}
	// -m,n
	vector<mint>x(n + m + 1);
	rep(i, n + m + 1)x[i] = i - m;
	auto f = kwm_t::math::fps::multipoint_evaluation(a, x);
	auto g2 = atcoder::convolution(f, b);
	vector<mint>g(n + 1);
	rep(i, n + 1) g[i] = g2[m + i];
	auto ans = kmw_t::math::modint::lagrange_coefficients(g);
	for (auto e : ans)cout << e.val() << " ";
	cout << endl;
	return 0;
}
yukicoder

結果

ソースコード
