結果

問題 No.3182 recurrence relation’s intersection sum
ユーザー kwm_t
提出日時 2025-06-14 14:10:52
言語 C++23
(gcc 13.3.0 + boost 1.87.0)
結果
AC  
実行時間 8 ms / 2,000 ms
コード長 4,837 bytes
コンパイル時間 6,488 ms
コンパイル使用メモリ 334,260 KB
実行使用メモリ 7,844 KB
最終ジャッジ日時 2025-06-14 14:11:01
合計ジャッジ時間 8,353 ms
ジャッジサーバーID
（参考情報） 		judge2 / judge5
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
sample AC * 3
other AC * 40
権限があれば一括ダウンロードができます

ソースコード

diff # 

#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; }
#include <vector>
template<class T>
struct Formal_power_series : std::vector<T> {
	using std::vector<T>::vector;
	using std::vector<T>::operator=;
	using F = Formal_power_series;
	F operator-() const {
		F res(*this);
		for (auto &e : res) e = -e;
		return res;
	}
	F inv(int d = -1) const {
		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() + min(n, 2 * m));
			F g(res);
			f.resize(2 * m), internal::butterfly(f);
			g.resize(2 * m), internal::butterfly(g);
			for (int i = 0; i < 2 * m; ++i) f[i] *= g[i];
			internal::butterfly_inv(f);
			f.erase(f.begin(), f.begin() + m);
			f.resize(2 * m), internal::butterfly(f);
			for (int i = 0; i < 2 * m; ++i) f[i] *= g[i];
			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);
		return res;
	}
	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(const F &g, const int d = -1) const { return F(*this).divide_inplace(g, d); }
};

template<class T>
struct Bostan_mori {
	std::vector<T> p, q;
	Bostan_mori(std::vector<T> &_p, std::vector<T> &_q) : p(_p), q(_q) {}
	void rev(std::vector<T> &f) const {
		int d = f.size();
		for (int i = 0; i < d; ++i) if (i & 1) f[i] = -f[i];
	}
	void even(std::vector<T> &f) const {
		int d = (f.size() + 1) >> 1;
		for (int i = 0; i < d; ++i) f[i] = f[i << 1];
		f.resize(d);
	}
	void odd(std::vector<T> &f) const {
		int d = f.size() >> 1;
		for (int i = 0; i < d; ++i) f[i] = f[i << 1 | 1];
		f.resize(d);
	}
	T operator[] (long long n) const {
		std::vector<T> _p(p), _q(q), _q_rev(q);
		rev(_q_rev);
		for (; n; n >>= 1) {
			_p = convolution(move(_p), _q_rev);
			if (n & 1) odd(_p);
			else     even(_p);
			_q = convolution(move(_q), move(_q_rev));
			even(_q);
			_q_rev = _q; rev(_q_rev);
		}
		return _p[0] / _q[0];
	}
};
std::vector<mint> BarlekampMassey(const std::vector<mint> &A) { //線形漸化式の初めのn項の特性方程式を返す d次ならO(d^2).
	//d次->2*d項があれば十分 ΣCix^i=0のCiを返す C0=1.
	int M = A.size();
	mint y = 1; //1つ前のxの値.
	std::vector<mint> C, B;//Cはret Bは1つ前のC.
	C.reserve(M + 1), B.reserve(M + 1);
	C.emplace_back(mint(1)), B.emplace_back(mint(1));

	for (int i = 1; i <= M; i++) {
		int lc = C.size(), lb = B.size();
		mint x = 0;
		for (int k = 0; k < lc; k++) x += C.at(k)*A.at(i - lc + k); //revしている?.
		B.emplace_back(mint(0)); lb++;
		if (x == 0) continue;
		mint Multi = x / y;
		if (lc < lb) { //C<-C-x/y*x^m*B B.emplace_back(0)でx^mを担当?.
			auto memo = C;
			if (lb > lc)C.insert(C.begin(), (lb - lc), 0);
			for (int k = 0; k < lb; k++) C.at(k) -= Multi * B.at(k); //O(d^2)の「掛け算」が重い モンゴメれば2倍速.
			B = memo; y = x;
		}
		else for (int k = 0; k < lb; k++) C.at(lc - 1 - k) -= Multi * B.at(lb - 1 - k);
	}
	reverse(C.begin(), C.end());//revで考えていたので元に戻す.
	return C; //何やってるか分からないシリーズ.
}
// 線形漸化式の初めのn項からN項目を求める d次漸化式ならO(d^2+dlogdlogN).
mint bmbm(long long N, const std::vector<mint> &A) {
	if (N < 0) return 0;
	auto C = BarlekampMassey(A);
	auto CA = A;
	CA.resize(C.size() - 1);
	auto B = atcoder::convolution(CA, C);
	B.resize(C.size() - 1);

	Bostan_mori<mint> bm(B, C);
	return bm[N];
}
int main() {
	ios::sync_with_stdio(false);
	cin.tie(nullptr);

	long long K, L, R; cin >> K >> L >> R;
	int n = 5000;
	std::vector<mint> A(n);
	A[0] = 1;
	for (int i = 1; i < n; i++) A.at(i) = A.at(i - 1)*K + mint(i - 1).pow(K) + mint(K).pow(i - 1);
	for (int i = 1; i < n; ++i) A.at(i) += A.at(i - 1);
	cout << (bmbm(R, A) - bmbm(L - 1, A)).val() << endl;
	return 0;
	return 0;
}
0