拓展卢卡斯定理模板完整注释

#define _CRT_SECURE_NO_WARNINGS
#include<bits/stdc++.h>
using namespace std;
//#define ll long long
typedef long long ll;
//拓展gcd,求解ax+by=gcd(a,b)
void ex_gcd(ll a, ll b, ll & x, ll &y) {
    if (!b) {
        x = 1; y = 0;
    }
    else {
        ex_gcd(b, a%b, y, x);
        y -= x * (a / b);
    }
}
//求解a关于p的逆元，其中gcd(a,p)==1
ll inverse(ll a, ll p) {
    if (!a)return 0ll;
    ll x, y;
    ex_gcd(a, p, x, y);
    x = (x%p + p) % p;
    return x;
}
//快速幂a^b
ll qpow(ll a, ll b, ll p) {
    ll res = 1;
    while (b) {
        if (b & 1)res = res * a%p;
        a = a * a%p;
        b >>= 1;
    }
    return res;
}
//中国剩余定理，求解x=a[i]%m[i]方程组
ll CRT(int n, ll *a, ll *m) {
    ll M = 1, res = 0;
    for (int i = 0; i < n; i++)M = M * m[i];
    for (int i = 0; i < n; i++) {
        ll x, y;
        ll tm = M / m[i];
        ex_gcd(tm, m[i], x, y);
        res = (res + tm * x*a[i] % M) % M;
    }
    return (res + M) % M;
}
//x表示质数，P表示分解的x^t
//calc计算n！%P，P不一定质数,复杂度P*log(x)(n)
ll calc(ll n, ll x, ll P) {
    if (!n)return 1;
    ll s = 1;
    //线性时间内计算1~P中不是x的倍数的乘积
    //由于P不是质数，不能用威尔逊定理
    //如果多组，这里可以提前打表优化，不过好像也不太方便，P总数太多
    for (ll i = 2; i <= P; i++) {
        if (i%x)s = s * i%P;
    }
    s = qpow(s, n / P, P);//前面一部分的数关于P下面的循环
    //这里是计算多余的，也就是n%P剩下的
    for (ll i =2; i <=n%P; i++)
        if (i%x)s = s*i%P;
    return s * calc(n / x, x, P) % P;
}
//x表示质数，P表示分解的x^t
//计算C(n,m)%P,P为一个质数的幂
ll multilucas(ll n, ll m, ll x, ll P) {
    ll cnt = 0;
    //下面分别用log的复杂度三个阶乘中含x的指数数目
    //这一块不懂得话可以了解n！%p的相关性质
    for (register ll i = n; i>0; i /= x)
        cnt += i / x;
    for ( ll i = m; i>0; i /= x)cnt -= i / x;
    for (ll i = n - m; i>0; i /= x)cnt -= i / x;
    return qpow(x, cnt, P) % P*calc(n, x, P) % P*inverse(calc(m, x, P), P) % P*inverse(calc(n - m, x, P), P) % P;
}
//求解C(n,m)%P的最终答案
ll ex_lucas(ll n, ll m, ll P) {
    ll cnt = 0;
    //a数组保存C(n,m)对每一种分解之后质因数求余的结果
    //由于每一种质因数的t次方是互质的，所以用CRT即可求解
    //p数组保存每一类质因数分解之后的答案
    ll p[20], a[20];//最多质因数分解20个
    for (ll i = 2; i*i <= P; i++) {
        //提前打素数表可以优化一点
        if (P%i == 0) {
            p[cnt] = 1;
            while (P%i == 0) { p[cnt] = p[cnt] * i; P /= i; }
            a[cnt] = multilucas(n, m, i, p[cnt]);
            cnt++;
        }
    }
    //最后是一个质数
    if (P > 1)p[cnt] = P, a[cnt++] = multilucas(n, m, P, P);
    return CRT(cnt, a, p);
}
int main()
{
    ll n, m, p;
    scanf("%lld%lld%lld", &n, &m, &p);
    printf("%lld\n", ex_lucas(n, m, p));
    cin >> n;
    return 0;
}