fork(4) download
copy 
  1. #include <vector>
  2. #include <cassert>
  3.  
  4. // ------------ Modulo Class ------------ //
  5. template<unsigned mod>
  6. class modulo {
  7. private:
  8. 	unsigned x;
  9. public:
  10. 	modulo() : x(0) {};
  11. 	modulo(unsigned x_) : x(x_) {};
  12. 	operator unsigned() { return x; }
  13. 	modulo operator==(const modulo& m) const { return x == m.x; }
  14. 	modulo operator!=(const modulo& m) const { return x != m.x; }
  15. 	modulo& operator+=(const modulo& m) { x = (x + m.x >= mod ? x + m.x - mod : x + m.x); return *this; }
  16. 	modulo& operator-=(const modulo& m) { x = (x < m.x ? x - m.x + mod : x - m.x); return *this; }
  17. 	modulo& operator*=(const modulo& m) { x = 1ULL * x * m.x % mod; return *this; }
  18. 	modulo operator+(const modulo& m) const { return modulo(*this) += m; }
  19. 	modulo operator-(const modulo& m) const { return modulo(*this) -= m; }
  20. 	modulo operator*(const modulo& m) const { return modulo(*this) *= m; }
  21. };
  22.  
  23. // ------------ Matrix Functions ------------ //
  24. typedef std::vector<modulo<998244353> > matrix_base;
  25. typedef std::vector<matrix_base> matrix;
  26. matrix mul(const matrix& a, const matrix& b) {
  27. 	assert(a[0].size() == b.size());
  28. 	matrix ret(a.size(), matrix_base(b[0].size(), 0));
  29. 	for (int i = 0; i < a.size(); i++) {
  30. 		for (int j = 0; j < b[0].size(); j++) {
  31. 			for (int k = 0; k < b.size(); k++) ret[i][j] += a[i][k] * b[k][j];
  32. 		}
  33. 	}
  34. 	return ret;
  35. }
  36. matrix unit(int n) {
  37. 	matrix ret(n, matrix_base(n, 0));
  38. 	for (int i = 0; i < n; i++) ret[i][i] = 1;
  39. 	return ret;
  40. }
  41. matrix power(const matrix& a, long long b) {
  42. 	assert(a.size() == a[0].size());
  43. 	matrix f = a, ret = unit(a.size());
  44. 	while (b) {
  45. 		if (b & 1) ret = mul(ret, f);
  46. 		f = mul(f, f);
  47. 		b >>= 1;
  48. 	}
  49. 	return ret;
  50. }
  51.  
  52. // ------------ Modpower Algorithm ------------ //
  53. inline int modpow(int a, int b, int m) {
  54. 	int ret = 1;
  55. 	while (b) {
  56. 		if (b & 1) ret = 1LL * ret * a % m;
  57. 		a = 1LL * a * a % m;
  58. 		b >>= 1;
  59. 	}
  60. 	return ret;
  61. }
  62.  
  63. // ------------ Number Theoretic Transform ------------ //
  64. inline static std::vector<int> FastModuloTransform(std::vector<int> v, int base, int root) {
  65. 	int n = v.size();
  66. 	for (int i = 0, j = 1; j < n - 1; j++) {
  67. 		for (int k = n >> 1; k >(i ^= k); k >>= 1);
  68. 		if (i < j) std::swap(v[i], v[j]);
  69. 	}
  70. 	for (int b = 1; b <= n / 2; b *= 2) {
  71. 		int x = modpow(root, (base - 1) / (b << 1), base);
  72. 		for (int i = 0; i < n; i += (b << 1)) {
  73. 			int p = 1;
  74. 			for (int j = i; j < i + b; j++) {
  75. 				int t1 = v[j], t2 = 1LL * v[j + b] * p % base;
  76. 				v[j] = t1 + t2; v[j] = (v[j] < base ? v[j] : v[j] - base);
  77. 				v[j + b] = t1 - t2 + base; v[j + b] = (v[j + b] < base ? v[j + b] : v[j + b] - base);
  78. 				p = 1LL * p * x % base;
  79. 			}
  80. 		}
  81. 	}
  82. 	return v;
  83. }
  84. inline static std::vector<int> FastConvolutionMod(std::vector<int> v1, std::vector<int> v2, int mod, int tr) {
  85. 	int n = v1.size() * 2; // v1 and v2 must be the same size!!
  86. 	v1.resize(n);
  87. 	v2.resize(n);
  88. 	v1 = FastModuloTransform(v1, mod, tr);
  89. 	v2 = FastModuloTransform(v2, mod, tr);
  90. 	for (int i = 0; i < n; i++) v1[i] = 1LL * v1[i] * v2[i] % mod;
  91. 	v1 = FastModuloTransform(v1, mod, modpow(tr, mod - 2, mod));
  92. 	int t = modpow(n, mod - 2, mod);
  93. 	for (int i = 0; i < n; i++) v1[i] = 1LL * v1[i] * t % mod;
  94. 	return v1;
  95. }
  96.  
  97. // ------------ Lagrange Interpolation ------------ //
  98. std::vector<int> lagrange_interpolation(std::vector<int> &v, int m) {
  99. 	int n = v.size() - 1;
  100. 	std::vector<int> inv(n + 2); inv[1] = 1;
  101. 	for (int i = 2; i <= n; i++) inv[i] = 1LL * inv[m % i] * (m - m / i) % m;
  102. 	std::vector<int> ret(n + 1);
  103. 	int q = 1;
  104. 	for (int i = 1; i <= n; i++) q = 1LL * q * inv[i] % m;
  105. 	if (n % 2 == 1) q = (m - q) % m;
  106. 	for (int i = 0; i <= n; i++) {
  107. 		ret[i] = 1LL * v[i] * q % m;
  108. 		q = 1LL * q * (m - n + i) % m * inv[i + 1] % m;
  109. 	}
  110. 	return ret;
  111. }
  112. int lagrange_function(int x, std::vector<int> &v, int m) {
  113. 	int n = v.size() - 1;
  114. 	int mul = 1;
  115. 	for (int i = 0; i <= n; i++) mul = 1LL * mul * (x - i + m) % m;
  116. 	int ret = 0;
  117. 	for (int i = 0; i <= n; i++) ret = (ret + 1LL * v[i] * modpow(x - i + m, m - 2, m)) % m;
  118. 	return 1LL * ret * mul % m;
  119. }
  120.  
  121. // ------------ Fibonacci Number ------------ //
  122. int nth_fibonacci(long long x) {
  123. 	matrix e(2, matrix_base(2));
  124. 	e[0][0] = e[0][1] = e[1][0] = 1;
  125. 	matrix p(2, matrix_base(1));
  126. 	p[0][0] = 1; p[1][0] = 0;
  127. 	return mul(power(e, x), p)[1][0];
  128. }
  129.  
  130. #include <iostream>
  131. using namespace std;
  132. const int mod = 998244353;
  133. int fact[1000009], inv[1000009], factinv[1000009];
  134. int combination(int a, int b) {
  135. 	if (a < 0 || b < 0 || a < b) return 0;
  136. 	return 1LL * fact[a] * factinv[b] % mod * factinv[a - b] % mod;
  137. }
  138. int n; long long m;
  139. int main() {
  140. 	cin >> n >> m; n--; m--;
  141. 	if (n == 0) {
  142. 		cout << nth_fibonacci(m + 1) << endl;
  143. 	}
  144. 	else {
  145. 		fact[0] = 1;
  146. 		for (int i = 1; i < 2 * n; i++) fact[i] = 1LL * fact[i - 1] * i % mod;
  147. 		inv[1] = 1;
  148. 		for (int i = 2; i < 2 * n; i++) inv[i] = 1LL * inv[mod % i] * (mod - mod / i) % mod;
  149. 		factinv[0] = 1;
  150. 		for (int i = 1; i < 2 * n; i++) factinv[i] = 1LL * factinv[i - 1] * inv[i] % mod;
  151. 		vector<int> a(3 * n); a[0] = a[1] = 1;
  152. 		for (int i = 2; i < 3 * n; i++) a[i] = (a[i - 1] + a[i - 2]) % mod;
  153. 		int e = 1; while (e < n) e <<= 1;
  154. 		vector<int> b(e), bt(a.begin(), a.begin() + e);
  155. 		for (int i = 0; i < n; i++) b[i] = combination((n - 1) + i, i);
  156. 		vector<int> res = FastConvolutionMod(b, bt, mod, 3);
  157. 		if (m < n) {
  158. 			cout << res[m] << endl;
  159. 		}
  160. 		else {
  161. 			vector<int> diff(n);
  162. 			for (int i = 0; i < n; i++) diff[i] = (a[2 * n + i] - res[i] + mod) % mod;
  163. 			vector<int> poly = lagrange_interpolation(diff, mod);
  164. 			int f = nth_fibonacci(m + 2 * n + 1);
  165. 			int ret = lagrange_function(m % mod, poly, mod);
  166. 			cout << (f - ret + mod) % mod << endl;
  167. 		}
  168. 	}
  169. 	return 0;
  170. }
Success #stdin #stdout 0s 15200KB
comments (?)
stdin
copy 
3 20
stdout
copy 
46345
https://ideone.com/NLbvPn
language:
C++14 (gcc 8.3)
created:
8 years ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!