CRT, lucas及其扩展形式
CRT, lucas及其扩展形式
exgcd
int exgcd(int a, int b, int &x, int &y) {
if (b == 0) return a, x = 1, y = 0;
int y = exgcd(b, a % b, x, y), t;
t = x, x = y, y = t - a / b * y;
}
证明:
gcd的过程中, 假设我们已经求出了\(b * x + (a~\%~b) * y = gcd(a, b)\)推导到\(a*x + b*y = gcd(a, b)\)的形式
中国剩余定理(CRT)
问题: 求最小整数解x, 其中\(P_1,P_2\cdots P_n\)互素
设\(M = \Pi^n_{i=1} P_i~~~ M_i=\frac M{P_i}~~~M_i^{i-1}\)为\(modM_i\)意义下的逆元 那么
扩展中国剩余定理(exCRT)
问题同中国剩余定理但各模数P不互质
假设已经求出前\(k-1\)个方程组成的同余方程组的一个解为\(x\), 且有 \(M=LCM_{i-1}^{k-1}P_i\), 则满足前k个方程的解为\(x + t * M\)
加入第k个方程得\(x + t * M = a_i \mod P_k\)
即 \(t*M=a_i-x \mod P_i\)
\(exgcd\)求解即可, 同时可以判是否有解
Lucas定理
解决的问题, P是质数且不大, 而n, m可以很大
证明咕咕咕, 其实没啥用
其实可以理解为转化为p进制下
ll C(int n,int m) {
if (m > n) return 0;
return jie[n] * fpw(jie[m], P-2) % P * fpw(jie[n-m], P-2) % P;
}
ll lau(ll x, ll y) {
if (!y) return 1;
return C(x % P, y % P) * lau(x / P, y / P) % P;
}
exLucas定理
这个比较有用, 不要求P是质数了
设\(P= \prod p_i^k\), 如果求出了\(\binom{n}{m}\bmod{p_i^{k_i}}\) 就可以利用CRT合并答案了
把p都提取出来, 可以先计算上面有多少p, 下面有多少p
然后考虑其他的数, 即求\(K!\mod P\), 发现在模\(p^k\)的时候有循环节, 因为\(p^k+i \equiv p^k\)
所以对于与p互质的数直接暴力扫描一个\(p^k\)然后快速幂即可
对与不与p互质的数, 去除所有p后, 仍会有其他因子, 不妨递归解决子问题, 即先提取出一个p来, p的倍数为\(1p, 2p, 3p \cdots \frac npp\) ,都去掉一个p然后相乘等于\(\frac np!\) , 递归下去就好喽
#include<iostream>
#include<cstring>
#include<cmath>
#include<cstdio>
#define ll long long
using namespace std;
ll m, n, p;
ll exgcd(ll a, ll b, ll &x, ll &y) {
if (!b) return x = 1, y = 0, a;
ll gcd = exgcd(b, a % b, x, y), t;
t = x, x = y, y = t - a / b * y;
return gcd;
}
ll inline inv(ll n, ll mod) {
ll x, y; exgcd(n, mod, x, y);
return (x % mod + mod) % mod;
}
ll fpw(ll x, ll mi, ll mod) {
ll res = 1;
while (mi) {
if (mi & 1) res = res * x % mod;
mi >>= 1;
x = x * x % mod;
}
return res;
}
ll fac(ll n, ll pi, ll pk) {
if (!n) return 1;
ll res = 1;
for (ll i = 2;i <= pk; i++)
if (i % pi) res = res * i % pk;
res = fpw(res, n / pk, pk);
for (ll i = 2;i <= n % pk; i++)
if (i % pi) res = res * i % pk;
return res * fac(n / pi, pi, pk) % pk;
}
inline ll ptime(ll x, ll pi) {
ll res = 0;
for (;x ; x /= pi) res += x / pi;
return res;
}
ll C(ll n, ll m, ll pi, ll pk) {
ll up = fac(n, pi, pk), d1 = fac(m, pi, pk), d2 = fac(n - m, pi, pk);
ll k = ptime(n, pi) - ptime(n-m, pi) - ptime(m, pi);
return up * inv(d1, pk) % pk * inv(d2, pk) % pk * fpw(pi, k, pk) % pk;
}
ll CRT(ll b, ll mod) {
return b * inv(p / mod, mod) % p * (p / mod) % p;
}
inline ll exlucus(ll n, ll m) {
ll res = 0, tmp = p, pk;
ll lim = sqrt(p);
for (int i = 2;i <= lim; i++) {
if (tmp % i == 0) {
pk = 1;
while(tmp % i == 0) pk *= i, tmp /= i;
res = (res + CRT(C(n, m, i, pk), pk)) % p;
// cout << tt << endl;
}
}
if (tmp > 1) res = (res + CRT(C(n, m, tmp, tmp), tmp)) % p;
return res;
}
int main() {
cin >> n >> m >> p;
printf ("%lld\n", exlucus(n, m));
return 0;
}