/* * Copyright (c) 2019 Ryotaro Sato * All functions and data structures prior to the main function are released under the MIT license. * https://github.com/hitonanode/cplib-cpp/blob/master/LICENSE */ #include #include #include #include #include #include #include using ll = long long; // CUT begin // Solve ax+by=gcd(a, b) template Int extgcd(Int a, Int b, Int& x, Int& y) { Int d = a; if (b != 0) { d = extgcd(b, a % b, y, x), y -= (a / b) * x; } else { x = 1, y = 0; } return d; } // Calculate a^(-1) (MOD m) s if gcd(a, m) == 1 // Calculate x s.t. ax == gcd(a, m) MOD m template Int mod_inverse(Int a, Int m) { Int x, y; extgcd(a, m, x, y); x %= m; return x + (x < 0) * m; } // Require: 1 <= b // return: (g, x) s.t. g = gcd(a, b), xa = g MOD b, 0 <= x < b/g template /* constexpr */ std::pair inv_gcd(Int a, Int b) { a %= b; if (a < 0) a += b; if (a == 0) return { b, 0 }; Int s = b, t = a, m0 = 0, m1 = 1; while (t) { Int u = s / t; s -= t * u, m0 -= m1 * u; auto tmp = s; s = t, t = tmp, tmp = m0, m0 = m1, m1 = tmp; } if (m0 < 0) m0 += b / s; return { s, m0 }; } template /* constexpr */ std::pair crt(const std::vector& r, const std::vector& m) { assert(r.size() == m.size()); int n = int(r.size()); // Contracts: 0 <= r0 < m0 Int r0 = 0, m0 = 1; for (int i = 0; i < n; i++) { assert(1 <= m[i]); Int r1 = r[i] % m[i], m1 = m[i]; if (r1 < 0) r1 += m1; if (m0 < m1) { std::swap(r0, r1); std::swap(m0, m1); } if (m0 % m1 == 0) { if (r0 % m1 != r1) return { 0, 0 }; continue; } Int g, im; std::tie(g, im) = inv_gcd(m0, m1); Int u1 = m1 / g; if ((r1 - r0) % g) return { 0, 0 }; Int x = (r1 - r0) / g % u1 * im % u1; r0 += x * m0; m0 *= u1; if (r0 < 0) r0 += m0; } return { r0, m0 }; } // 蟻本 P.262 // 中国剰余定理を利用して,色々な素数で割った余りから元の値を復元 // 連立線形合同式 A * x = B mod M の解 // Requirement: M[i] > 0 // Output: x = first MOD second (if solution exists), (0, 0) (otherwise) template std::pair linear_congruence(const std::vector& A, const std::vector& B, const std::vector& M) { Int r = 0, m = 1; assert(A.size() == M.size()); assert(B.size() == M.size()); for (int i = 0; i < (int)A.size(); i++) { assert(M[i] > 0); const Int ai = A[i] % M[i]; Int a = ai * m, b = B[i] - ai * r, d = std::__gcd(M[i], a); if (b % d != 0) { return std::make_pair(0, 0); // 解なし } Int t = b / d * mod_inverse(a / d, M[i] / d) % (M[i] / d); r += m * t; m *= M[i] / d; } return std::make_pair((r < 0 ? r + m : r), m); } int pow_mod(int x, long long n, int md) { if (md == 1) return 0; long long ans = 1; while (n > 0) { if (n & 1) ans = ans * x % md; x = (long long)x * x % md; n >>= 1; } return ans; } #line 5 "number/combination.hpp" // nCr mod m = p^q (p: prime, q >= 1) // Can be used for n, r <= 1e18, m <= 1e7 // Complexity: O(m) (construction), O(log(n)) (per query) // https://ferin-tech.hatenablog.com/entry/2018/01/17/010829 struct combination_prime_pow { int p, q, m; std::vector fac, invfac, ppow; long long _ej(long long n) const { long long ret = 0; while (n) ret += n, n /= p; return ret; } combination_prime_pow(int p_, int q_) : p(p_), q(q_), m(1), ppow{ 1 } { for (int t = 0; t < q; ++t) m *= p, ppow.push_back(m); fac.assign(m, 1); invfac.assign(m, 1); for (int i = 1; i < m; ++i) fac[i] = (long long)fac[i - 1] * (i % p ? i : 1) % m; invfac[m - 1] = fac[m - 1]; // Same as Wilson's theorem assert(1LL * fac.back() * invfac.back() % m == 1); for (int i = m - 1; i; --i) invfac[i - 1] = (long long)invfac[i] * (i % p ? i : 1) % m; } int nCr(long long n, long long r) const { if (r < 0 or n < r) return 0; if (p == 2 and q == 1) return !((~n) & r); // Lucas long long k = n - r; long long e0 = _ej(n / p) - _ej(r / p) - _ej(k / p); if (e0 >= q) return 0; long long ret = ppow[e0]; if (q == 1) { // Lucas while (n) { ret = __int128(ret) * fac[n % p] * invfac[r % p] * invfac[k % p] % p; n /= p, r /= p, k /= p; } return (int)ret; } else { if ((p > 2 or q < 3) and (_ej(n / m) - _ej(r / m) - _ej(k / m)) & 1) ret = m - ret; while (n) { ret = __int128(ret) * fac[n % m] * invfac[r % m] * invfac[k % m] % m; n /= p, r /= p, k /= p; } return (int)ret; } } }; // nCr mod m // Complexity: O(m) space worst (construction), O(log(n) log(m)) (per query) // Input: pairs of (prime, degree), such as vector> and map // https://judge.yosupo.jp/problem/binomial_coefficient struct combination { std::vector cpps; std::vector ms; template combination(const Map& p2deg) { for (auto f : p2deg) { cpps.push_back(combination_prime_pow(f.first, f.second)); ms.push_back(cpps.back().m); } } int operator()(long long n, long long r) const { if (r < 0 or n < r) return 0; std::vector rs; for (const auto& cpp : cpps) rs.push_back(cpp.nCr(n, r)); return crt(rs, ms).first; } }; struct Sieve { std::vector min_factor; std::vector primes; Sieve(int MAXN) : min_factor(MAXN + 1) { for (int d = 2; d <= MAXN; d++) { if (!min_factor[d]) { min_factor[d] = d; primes.emplace_back(d); } for (const auto& p : primes) { if (p > min_factor[d] or d * p > MAXN) break; min_factor[d * p] = p; } } } // Prime factorization for 1 <= x <= MAXN^2 // Complexity: O(log x) (x <= MAXN) // O(MAXN / log MAXN) (MAXN < x <= MAXN^2) template std::map factorize(T x) const { std::map ret; assert(x > 0 and x <= ((long long)min_factor.size() - 1) * ((long long)min_factor.size() - 1)); for (const auto& p : primes) { if (x < T(min_factor.size())) break; while (!(x % p)) x /= p, ret[p]++; } if (x >= T(min_factor.size())) ret[x]++, x = 1; while (x > 1) ret[min_factor[x]]++, x /= min_factor[x]; return ret; } // Enumerate divisors of 1 <= x <= MAXN^2 // Be careful of highly composite numbers https://oeis.org/A002182/list // https://gist.github.com/dario2994/fb4713f252ca86c1254d#file-list-txt (n, (# of div. of n)): // 45360->100, 735134400(<1e9)->1344, 963761198400(<1e12)->6720 template std::vector divisors(T x) const { std::vector ret{ 1 }; for (const auto p : factorize(x)) { int n = ret.size(); for (int i = 0; i < n; i++) { for (T a = 1, d = 1; d <= p.second; d++) { a *= p.first; ret.push_back(ret[i] * a); } } } return ret; // NOT sorted } // Euler phi functions of divisors of given x // Verified: ABC212 G https://atcoder.jp/contests/abc212/tasks/abc212_g // Complexity: O(sqrt(x) + d(x)) template std::map euler_of_divisors(T x) const { assert(x >= 1); std::map ret; ret[1] = 1; std::vector divs{ 1 }; for (auto p : factorize(x)) { int n = ret.size(); for (int i = 0; i < n; i++) { ret[divs[i] * p.first] = ret[divs[i]] * (p.first - 1); divs.push_back(divs[i] * p.first); for (T a = divs[i] * p.first, d = 1; d < p.second; a *= p.first, d++) { ret[a * p.first] = ret[a] * p.first; divs.push_back(a * p.first); } } } return ret; } // Moebius function Table, (-1)^{# of different prime factors} for square-free x // return: [0=>0, 1=>1, 2=>-1, 3=>-1, 4=>0, 5=>-1, 6=>1, 7=>-1, 8=>0, ...] https://oeis.org/A008683 std::vector GenerateMoebiusFunctionTable() const { std::vector ret(min_factor.size()); for (unsigned i = 1; i < min_factor.size(); i++) { if (i == 1) { ret[i] = 1; } else if ((i / min_factor[i]) % min_factor[i] == 0) { ret[i] = 0; } else { ret[i] = -ret[i / min_factor[i]]; } } return ret; } // Calculate [0^K, 1^K, ..., nmax^K] in O(nmax) // Note: **0^0 == 1** template std::vector enumerate_kth_pows(long long K, int nmax) const { assert(nmax < int(min_factor.size())); assert(K >= 0); if (K == 0) return std::vector(nmax + 1, 1); std::vector ret(nmax + 1); ret[0] = 0, ret[1] = 1; for (int n = 2; n <= nmax; n++) { if (min_factor[n] == n) { ret[n] = MODINT(n).pow(K); } else { ret[n] = ret[n / min_factor[n]] * ret[min_factor[n]]; } } return ret; } }; /* * Main Code */ int main(void) { ll l, r, m; std::cin >> l >> r >> m; combination nCr(Sieve(1 << 20).factorize(m)); ll ans = 0; for (ll i = l; i <= r; i++) { ans += (ll)nCr(2 * i, i) + m - 2ll; ans %= m; } std::cout << ans << std::endl; return 0; }