卢卡斯定理
常用于求解 ( n m ) ( m o d P ) \binom{n}{m} (mod \ \ P) (mn)(mod P) , 且 n n n, m m m值域特别大, 而 P P P 的值域比较小
P P P为素数
直接记定理就好了 :
(
n
m
)
(
m
o
d
P
)
=
(
n
/
P
m
/
P
)
⋅
(
n
%
P
m
%
P
)
(
m
o
d
P
)
\binom{n}{m} (mod \ \ P)= \binom{n / P}{m / P} \cdot \binom {n \% P}{m \% P} (mod\ \ P)
(mn)(mod P)=(m/Pn/P)⋅(m%Pn%P)(mod P)
然后我们预处理出
1
∼
P
1 \sim P
1∼P的阶乘和逆元, 那么求这个组合数的复杂度就是
O
(
P
+
l
o
g
P
n
)
O(P + log_P^{n})
O(P+logPn)
代码 :
#include <bits/stdc++.h>
using namespace std;
const int maxn = 1e6 + 5;
int n, m, fac[maxn], P;
inline int read(){
char ch = getchar(); int u = 0, f = 1;
while(!isdigit(ch)){if(ch == '-')f = -1; ch = getchar();}
while(isdigit(ch)){u = u * 10 + ch - 48; ch = getchar();} return u * f;
}
inline int ksm(int x, int k){
int cnt = 1;
while(k){
if(k & 1)cnt = 1ll * cnt * x % P;
x = 1ll * x * x % P; k >>= 1;
} return cnt;
}
inline int C(int n, int m, int P){
if(n < m)return 0;
return (int)(1ll * fac[n] * ksm(fac[n - m], P - 2) % P * ksm(fac[m], P - 2) % P);
}
inline Lucas(int n, int m, int P){
if(n == 0 || m == 0)return 1;
return (int)(1ll * Lucas(n / P, m / P, P) * C(n % P, m % P, P) % P);
}
int main(){
n = read(); m = read(); P = read();
fac[0] = 1;
for(register int i = 1; i <= P; i++)
fac[i] = 1ll * fac[i - 1] * i % P;
cout << Lucas(n, m, P);
return 0;
}
P P P 不为素数
其实这个东西和卢卡斯定理并没有什么关系
因为
P
P
P不为素数, 所以 :
P
=
∏
P
i
a
i
P = \prod P_i^{a_i}
P=∏Piai
我们考虑求出
(
n
m
)
(
m
o
d
P
i
a
i
)
\binom {n}{m} (mod \ \ P_i^{a_i})
(mn)(mod Piai)
然后直接用
C
R
T
CRT
CRT合并就是答案了
那么现在问题转化为了如何快速求出 ( n m ) ( m o d P i a i ) \binom{n}{m}(mod\ \ P_i^{a_i}) (mn)(mod Piai)了
因为实在膜 P i a i P_i^{a_i} Piai 的意义下, 所以阶乘其实是有循环节的, 我们可以考虑直接暴力处理循环节和不完整的循环
难道这样就好了吗, 事实上并没有这么简单
因为
P
i
a
i
P_i^{a_i}
Piai是
P
i
P_i
Pi的倍数, 所以很有可能再膜意义下就直接变为零了, 但是事实上并不是这样的, 我们观察组合数的公式 :
(
n
m
)
=
n
!
(
n
−
m
)
!
m
!
\binom {n}{m} = \frac{n!}{(n - m)!m!}
(mn)=(n−m)!m!n!
所以分母和分子的
P
i
P_i
Pi是可已被约掉的, 所以我们不能直接取膜, 我们再考虑把是
P
i
P_i
Pi的倍数的数给单独提出来搞一下, 我们设
P
k
=
P
i
a
i
P_k = P_i ^ {a_i}
Pk=Piai, 那么则有
f
a
c
(
n
)
=
(
f
a
c
(
⌊
n
P
i
⌋
)
⋅
P
i
⌊
n
P
i
⌋
)
⋅
(
∏
i
=
1
P
k
i
)
⌊
n
P
k
⌋
⋅
(
∏
i
=
1
n
m
o
d
P
k
i
)
(
m
o
d
P
i
a
i
)
fac(n) = \left(fac(\lfloor \frac{n}{P_i} \rfloor) \cdot P_i^{\lfloor \frac{n}{P_i} \rfloor}\right) \cdot \left(\prod_{i = 1}^{P_k} i\right)^{\lfloor\frac{n}{P_k}\rfloor} \cdot \left(\prod_{i = 1}^{n ~mod~ P_k} i\right) \pmod {P_i^{a_i}}
fac(n)=(fac(⌊Pin⌋)⋅Pi⌊Pin⌋)⋅(i=1∏Pki)⌊Pkn⌋⋅(i=1∏n mod Pki)(modPiai)
我们发现这个是一个递归式, 那么直接递归就好了, 特别的如果
n
=
0
n = 0
n=0就返回
1
1
1
需要注意的是 P i P_i Pi这个东西不能在递归的时候直接乘上, 我们考虑, 在递归完之后, 或者递归之前直接处理出来就好了
代码 :
#include <bits/stdc++.h>
typedef long long LL;
using namespace std;
LL n, m, P;
inline LL read(){
char ch = getchar(); LL u = 0, f = 1;
while(!isdigit(ch)){if(ch == '-')f = -1; ch = getchar();}
while(isdigit(ch)){u = u * 10 + ch - 48; ch = getchar();} return u * f;
}
inline LL ksm(LL x, LL k, LL P){
LL cnt = 1;
while(k){
if(k & 1)cnt = cnt * x % P;
x = x * x % P; k >>= 1;
} return cnt;
}
inline LL Exgcd(LL a, LL b, LL &x, LL &y){
if(b == 0){x = 1; y = 0; return a;}
LL x0 = 0, y0 = 0;
LL d = Exgcd(b, a % b, y0, x0);
x = x0; y = y0 - ((a / b) * x0);
return d;
}
inline LL fac(LL n, LL Pi, LL Pk){
if(n == 0)return 1;
LL rnt = 1;
for(register int i = 2; i <= Pk; i++)
if(i % Pi)rnt = rnt * i % Pk;
rnt = ksm(rnt, n / Pk, Pk);
for(register int i = 2; i <= n % Pk; i++)
if(i % Pi)rnt = rnt * i % Pk;
return rnt * fac(n / Pi, Pi, Pk) % Pk;
}
inline LL inv(LL x, LL P){
LL inv, tmp;
Exgcd(x, P, inv, tmp);
return (inv + P) % P;
}
inline LL ExLucas(LL n, LL m, LL Pi, LL Pk){
LL C1 = fac(n, Pi, Pk), C2 = fac(n - m, Pi, Pk), C3 = fac(m, Pi, Pk);
LL tim = 0;
for(register LL i = n; i; i /= Pi)tim += i / Pi;
for(register LL i = n - m; i; i /= Pi)tim -= i / Pi;
for(register LL i = m; i; i /= Pi)tim -= i / Pi;
return C1 * inv(C2, Pk) % Pk * inv(C3, Pk) % Pk * ksm(Pi, tim, Pk) % Pk;
}
inline LL China(LL n, LL m, LL P){
LL A = 0, M = P, Pk = 1;
for(register LL i = 2; i <= P; i++)
if(P % i == 0){
Pk = 1;
while(P % i == 0)Pk *= i, P /= i;
A = (A + ExLucas(n, m, i, Pk) * inv(M / Pk, Pk) % M * (M / Pk) % M) % M;
}
return (A + M) % M;
}
int main(){
n = read(); m = read(); P = read();
cout << China(n, m, P);
return 0;
}