因为 \(x^{\overline p} \equiv x^p - x \pmod p\)，所以设 \(n = pq + r\)，其中 \(r \in [0, p - 1]\)，则有：
\[\begin{aligned} x^{\overline n} & = (\prod\limits_{i = 0}^{q - 1} (x + ip)^{\overline p}) (x + pq)^{\overline r} \\ & = (x^p - x)^q (x + pq)^{\overline r} \\ & = (x^p - x)^q x^{\overline r} \\ & = x^q (x^{p - 1} - 1)^q x^{\overline r} \end{aligned} \]设 \(k - q = a(p - 1) + b\)，其中 \(b \in [0, p - 2]\)，那么：
\[\begin{aligned} \begin{bmatrix} n \\ k \end{bmatrix} & = [x^k] x^{\overline n} \\ & = [x^k] x^q (x^{p - 1} - 1)^q x^{\overline r} \\ & = [x^{k - q}] (x^{p - 1} - 1)^q x^{\overline r} \\ & = (-1)^{q - a} \binom{q}{a} \begin{bmatrix} r \\ b \end{bmatrix} \end{aligned} \]设 \(m - q = k(p - 1) + t\)，其中 \(t \in [0, p - 2]\)，那么：
\[\sum\limits_{i = 0}^m \begin{bmatrix} n \\ i \end{bmatrix} = \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} \sum\limits_{j = 0}^k (-1)^{q - j} \binom{q}{j} + \sum\limits_{i = t + 1}^r \begin{bmatrix} r \\ i \end{bmatrix} \sum\limits_{j = 0}^{k - 1} (-1)^{q - j} \binom{q}{j} \]又因为：
\[\sum\limits_{i = 0}^m (-1)^i \binom{n}{i} = (-1)^m \binom{n - 1}{m} \]这个式子由 \(\binom{n}{m} = \binom{n - 1}{m} + \binom{n - 1}{m - 1}\) 易证。
那么：
\[\begin{aligned} \sum\limits_{i = 0}^m \begin{bmatrix} n \\ i \end{bmatrix} & = \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} (-1)^{q - k} \binom{q - 1}{k} + \sum\limits_{i = t + 1}^r \begin{bmatrix} r \\ i \end{bmatrix} (-1)^{q - k + 1} \binom{q - 1}{k - 1} \\ & = \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} (-1)^{q - k} \binom{q - 1}{k} + \sum\limits_{i = 0}^r \begin{bmatrix} r \\ i \end{bmatrix} (-1)^{q - k + 1} \binom{q - 1}{k - 1} + \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} (-1)^{q - k} \binom{q - 1}{k - 1} \\ & = (-1)^{q - k} (\binom{q}{k} \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} - \binom{q - 1}{k - 1} \sum\limits_{i = 0}^r \begin{bmatrix} r \\ i \end{bmatrix}) \\ & = (-1)^{q - k} (\binom{q}{k} \sum\limits_{i = 0}^t \begin{bmatrix} r \\ i \end{bmatrix} - \binom{q - 1}{k - 1} r!) \end{aligned} \]那么剩下的问题是给定 \(r, t\)，求 \(x^{\overline r}\) 的前 \(t\) 次项系数和。
考虑一个多项式 \(f(x)\)，设模数的原根为 \(g\)，若我们求出 \(f(g^0), f(g^1), \ldots, f(g^{p - 2})\)，那么我们可以线性地求出 \(f(x)\) 的前 \(t\) 次项系数和。
具体地，考虑 \(\sum\limits_{i = 0}^{p - 2} g^{ki}\)，只有当 \(k = 0\) 时这个式子等于 \(p - 1 \equiv -1 \pmod p\)，当 \(k \ne 0\) 时由等比数列求和公式可知它等于 \(0\)。
那么考虑若要求 \(f(x)\) 的第 \(j\) 次项，\(-\sum\limits_{i = 0}^{p - 2} g^{-ij} f(g^i)\) 即为答案。
所以：
\[\begin{aligned} \sum\limits_{i = 0}^t [x^i] x^{\overline r} & = \sum\limits_{i = 0}^t \sum\limits_{j = 0}^{p - 2} g^{-ij} f(g^j) \\ & = \sum\limits_{j = 0}^{p - 2} f(g^j) (\sum\limits_{i = 0}^t g^{-ij}) \end{aligned} \]求 \(f(g^i)\) 可以预处理阶乘及其逆元然后 \(O(1)\) 单次计算，后面那个 \(\sum\limits_{i = 0}^t g^{-ij}\) 可以用等比数列求和公式然后预处理原根的次方做到 \(O(1)\) 单次计算。
总时间复杂度 \(O(p)\)。
code
#include <bits/stdc++.h>
#define pb emplace_back
#define fst first
#define scd second
#define mkp make_pair
#define mems(a, x) memset((a), (x), sizeof(a))
using namespace std;
typedef long long ll;
typedef double db;
typedef unsigned long long ull;
typedef long double ldb;
typedef pair<ll, ll> pii;
const int maxn = 1000100;
ll n, L, R, P, fac[maxn], ifac[maxn], N, G, q, r, inv[maxn], pw[maxn];
inline ll qpow(ll b, ll p) {
ll res = 1;
while (p) {
if (p & 1) {
res = res * b % P;
}
b = b * b % P;
p >>= 1;
}
return res;
}
inline ll C(ll n, ll m) {
if (n < m || n < 0 || m < 0) {
return 0;
}
if (n >= P) {
return C(n / P, m / P) * C(n % P, m % P) % P;
} else {
return fac[n] * ifac[m] % P * ifac[n - m] % P;
}
}
inline bool check(ll x) {
ll t = P - 1;
vector<ll> vc;
for (ll i = 2; i * i <= t; ++i) {
if (t % i == 0) {
vc.pb(i);
while (t % i == 0) {
t /= i;
}
}
}
if (t > 1) {
vc.pb(t);
}
for (ll y : vc) {
if (qpow(x, (P - 1) / y) == 1) {
return 0;
}
}
return 1;
}
inline ll g(ll n, ll m) {
ll p = 1, ip = 1, iG = qpow(G, P - 2), res = 0;
for (int i = 0; i <= P - 2; ++i, p = p * G % P, ip = ip * iG % P) {
if (p + n - 1 >= P) {
continue;
}
if (ip == 1) {
res = (res + (m + 1) * fac[p + n - 1] % P * ifac[p - 1]) % P;
} else {
res = (res + (pw[((-(m + 1) * i) % (P - 1) + P - 1) % (P - 1)] + P - 1) % P * inv[(ip + P - 1) % P] % P * fac[p + n - 1] % P * ifac[p - 1]) % P;
}
}
return (P - res) % P;
}
inline ll h(ll m) {
if (m < q) {
return 0;
}
ll k = (m - q) / (P - 1);
ll t = (m - q) % (P - 1);
ll ans = C(q, k) * g(r, t) % P;
ans = (ans - C(q - 1, k - 1) * fac[r] % P + P) % P;
if ((q - k) & 1) {
ans = (P - ans) % P;
}
return ans;
}
void solve() {
scanf("%lld%lld%lld%lld", &n, &L, &R, &P);
N = P - 1;
fac[0] = 1;
for (int i = 1; i <= N; ++i) {
fac[i] = fac[i - 1] * i % P;
}
ifac[N] = qpow(fac[N], P - 2);
for (int i = N - 1; ~i; --i) {
ifac[i] = ifac[i + 1] * (i + 1) % P;
}
inv[1] = 1;
for (int i = 2; i <= N; ++i) {
inv[i] = (P - P / i) * inv[P % i] % P;
}
for (G = 1;; ++G) {
if (check(G)) {
break;
}
}
pw[0] = 1;
for (int i = 1; i <= P - 2; ++i) {
pw[i] = pw[i - 1] * G % P;
}
q = n / P;
r = n - P * q;
printf("%lld\n", (h(R) - h(L - 1) + P) % P);
}
int main() {
int T = 1;
// scanf("%d", &T);
while (T--) {
solve();
}
return 0;
}