Lucas定理(求组合数取模) 扩展Lucas定理(解决模数非质情况)

在比赛时 , 如果遇到$C_n^m$的n比较大 , 我们不能通过预处理阶乘和逆元来计算 , 而题目又要求对答案取一个质数模的时候 , 我们可以用Lucas定理来简化计算

Lucas 定理：

定义 : n,m是非负整数，p是素数时 , $Lucas(n,m)=C_n^m\;\% p$

公式 :

$$Lucas(n,m)=C_{n\%p}^{m\%p}*Lucas(n/p,m/p)\;\% p$$

代码 :

LL Lucas(LL n,LL m){
    if(m==0)return 1ll;
    return C(n%mod,m%mod)*Lucas(n/mod,m/mod)%mod;
}

$\\$

---

$\\$

扩展Lucas定理:

中国剩余定理 (CRT)

扩展Lucas定理用于解决模数非质数时的情况

懒得敲了,直接手写上图好了(好久没写作业不会写字了多多包涵)

也就是说 , 我们对于模数M , 在非质数的时候 , 把它分成多个两两互质的数(按照质因子分解一定保证两两互质) , 求出$C\%p_i^{q_i}所得a_i$ , 然后把多组$p_i^{q_i}和a_i$求CRT就是答案

接下来要解决的是$C\%p_i^{q_i}$ $\\$ C实际上的三个阶乘的乘除 , 所以当下解决的就变成了$n!\%p_i^{q_i}$ $\\$ 假设为$19!\%3^2$ $\\$ $19!=1*2*3*4*5*6*7*8* ……*19$ $\quad\;\;=[1*2*4*5*7*8*… *16*17*19]*(3*6*9*12*15*18)$ $\quad\;\;=[1*2*4*5*7*8*… *16*17*19]*3^6(1*2*3*4*5*6)$ $\\$ 1. 后半部分也是阶乘 , 直接重新递归即可 2. 前半部分中 , 有一个循环节 , $[1*2*4*5*7*8]\%9=[(1+9k)*(2+9k)*(4+9k)*(5+9k)*(7+9k)*(8+9k)]\%9$ , 对于后面的零头暴力就行了 3. 对于中间p的此法 , 记录一下数量 , 最后三个阶乘的数量相加减即可

模板 :

#include<bits/stdc++.h>
using namespace std;
#define i64 long long


i64 POW(i64 a,i64 b,i64 mod)
{
    i64 ans=1;
    while(b)
    {
        if(b&1) ans=ans*a%mod;
        a=a*a%mod;
        b>>=1;
    }
    return ans;
}

i64 POW(i64 a,i64 b)
{
    i64 ans=1;
    while(b)
    {
        if(b&1) ans=ans*a;
        a=a*a;
        b>>=1;
    }
    return ans;
}


i64 exGcd(i64 a,i64 b,i64 &x,i64 &y)
{
    i64 t,d;
    if(!b)
    {
        x=1;
        y=0;
        return a;
    }
    d=exGcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return d;
}

bool modular(i64 a[],i64 m[],i64 k)
{
    i64 d,t,c,x,y,i;

    for(i=2;i<=k;i++)
    {
        d=exGcd(m[1],m[i],x,y);
        c=a[i]-a[1];
        if(c%d) return false;
        t=m[i]/d;
        x=(c/d*x%t+t)%t;
        a[1]=m[1]*x+a[1];
        m[1]=m[1]*m[i]/d;
    }
    return true;
}



i64 reverse(i64 a,i64 b)
{
    i64 x,y;
    exGcd(a,b,x,y);
    return (x%b+b)%b;
}

i64 C(i64 n,i64 m,i64 mod)
{
    if(m>n) return 0;
    i64 ans=1,i,a,b;
    for(i=1;i<=m;i++)
    {
        a=(n+1-i)%mod;
        b=reverse(i%mod,mod);
        ans=ans*a%mod*b%mod;
    }
    return ans;
}

i64 C1(i64 n,i64 m,i64 mod)
{
    if(m==0) return 1;
    return C(n%mod,m%mod,mod)*C1(n/mod,m/mod,mod)%mod;
}

i64 cal(i64 n,i64 p,i64 t)
{
    if(!n) return 1;
    i64 x=POW(p,t),i,y=n/x,temp=1;
    for(i=1;i<=x;i++) if(i%p) temp=temp*i%x;
    i64 ans=POW(temp,y,x);
    for(i=y*x+1;i<=n;i++) if(i%p) ans=ans*i%x;
    return ans*cal(n/p,p,t)%x;
}

i64 C2(i64 n,i64 m,i64 p,i64 t)
{
    i64 x=POW(p,t);
    i64 a,b,c,ap=0,bp=0,cp=0,temp;
    for(temp=n;temp;temp/=p) ap+=temp/p;
    for(temp=m;temp;temp/=p) bp+=temp/p;
    for(temp=n-m;temp;temp/=p) cp+=temp/p;
    ap=ap-bp-cp;
    i64 ans=POW(p,ap,x);
    a=cal(n,p,t);
    b=cal(m,p,t);
    c=cal(n-m,p,t);
    ans=ans*a%x*reverse(b,x)%x*reverse(c,x)%x;
    return ans;
}

//计算C(n,m)%mod
i64 Lucas(i64 n,i64 m,i64 mod)
{
    i64 i,t,cnt=0;
    i64 A[205],M[205];
    for(i=2;i*i<=mod;i++) if(mod%i==0)
    {
        t=0;
        while(mod%i==0)
        {
            t++;
            mod/=i;
        }
        M[++cnt]=POW(i,t);
        if(t==1) A[cnt]=C1(n,m,i);
        else A[cnt]=C2(n,m,i,t);
    }
    if(mod>1)
    {
        M[++cnt]=mod;
        A[cnt]=C1(n,m,mod);
    }
    modular(A,M,cnt);
    return A[1];
}


int main(){
    i64 n,m,mod;
    while(scanf("%lld%lld%lld",&n,&m,&mod))
    printf("%lld\n",Lucas(n,m,mod));
}