卢卡斯（lucas）定理

对于质数 \(p\)，有

\[{ \begin{aligned} & \binom{n}{m} \equiv \binom{\left \lfloor n/p \right \rfloor }{\left \lfloor m/p \right \rfloor } \binom{n\mod{p}}{m\mod p} \pmod{p} \end{aligned} } \]

引理1

\[{ \begin{aligned} & \binom{p}{n}\mod{p}=[n=0\vee n=p] \end{aligned} } \]

引理2

\[{ \begin{aligned} & (a+b)^p \equiv a^p+b^p \pmod{p} \end{aligned} } \]

P3807 【模板】卢卡斯定理/Lucas 定理

直接应用上面的。\(n、m\)较小时就可以直接暴力求。由于每个 \(p\) 不同，因此不能预处理阶乘。
又因为 \(p\) 可能比 \(n、m\) 小，因此不能递推求逆元。

点击查看代码

#include<bits/stdc++.h>
using namespace std;
#define int long long
int p,cnt;
int ksm(int x,int k)
{
	int res=1;
	while(k)
	{
		if(k&1) res=res*x%p;
		k>>=1;x=x*x%p;
	}
	return res;
}
int C(int n,int m)
{
	int jie1=1,jie2=1;
	for(int i=n;i>=n-m+1;i--) jie1=jie1*i%p;
	for(int i=1;i<=m;i++) jie2=jie2*i%p;
	jie2=ksm(jie2,p-2);return jie1*jie2%p;
}
int lucas(int n,int m)
{
	if(m==0) return 1;
	return (C(n%p,m%p)*lucas(n/p,m/p))%p;
}
void solve()
{
	int n,m;cin>>n>>m>>p;n=n+m;cnt=0;
	cout<<lucas(n,m)<<'\n';
}
signed main()
{
	ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);
	int T;cin>>T;while(T--) solve();
}

P4720 【模板】扩展卢卡斯定理/exLucas

求 \({n\choose m }\mod p\)，其中 \(p\) 不一定是质数。

先将 \(p\) 拆开成为 \(\prod_{p_i} p_i^{k_i}\)，然后将答案用普通 CRT 合并即可。

考虑如何做 \({n\choose m}\mod p_{i}^{k_i}\) 。一般而言我们是将其拆成 \(\frac{n!}{m!(n-m)!}\) 的形式。但是可能没有逆元，考虑将所有 \(p\) 的因子拿出来切换成

\[\frac{\frac{n!}{p_i^x}}{\frac{m!}{p_i^y}\frac{(n-m)!}{p_i^z}}\times p_i^{x-y-z} \mod p_{i}^{k_i} \]

的形式。这样一定互质就可以 exgcd 求逆元了。（为了方便，下面将 \(p_i\) 和 \(k_i\) 写作 \(p\) 和 \(k\)）

于是问题变成一些形如 \(\frac{n!}{p^x}\) 的东西。同时我们还有求出 \(x\) 的值。

考虑将有 \(p\) 的因子的所有数单独拿出来。有

\[p^{\left \lfloor \frac n p \right \rfloor }\times \left \lfloor \frac n p \right \rfloor !\times \prod_{i=1,i\not\equiv 0(\bmod p)}^n i \pmod{p^k} \]

考虑最后那一项很巨大。发现我们可以按 \(p^k\) 分类。考虑

\[\prod_{i=1,i\not\equiv 0(\bmod p)}^n i =\left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{p^k} i \right)^{\left \lfloor \frac n {p^k} \right \rfloor}\times \left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{n\bmod p^k} i \right)\pmod{p^k} \]

发现这两个可以 \(O(p^k)\) 预处理。
于是原来的式子变成了

\[p^{\left \lfloor \frac n p \right \rfloor }\times \left \lfloor \frac n p \right \rfloor !\times \left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{p^k} i \right)^{\left \lfloor \frac n {p^k} \right \rfloor}\times \left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{n\bmod p^k} i \right) \pmod{p^k} \]

发现有一项类似的形式即 \(\left \lfloor \frac n p \right \rfloor !\)。我们可以递归下去求解。
设

\[F(n)=\frac{n!}{p^x}= F\left(\left \lfloor \frac n p \right \rfloor \right)\times \left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{p^k} i \right)^{\left \lfloor \frac n {p^k} \right \rfloor}\times \left (\prod_{i=1,i\not\equiv 0(\bmod p)}^{n\bmod p^k} i \right) \pmod{p^k} \]

最前面那一项没了是因为是 \(p\) 的倍数会被除掉。于是一些阶乘的答案就可以这样求出来。

然后直接求逆元然后 CRT 就做完了。对于 \(p^{x-y-z}\) 那一项，我们发现每一次递归的 \(p\) 除掉的 \(p\) 的次数就是 $\left \lfloor \frac n p \right \rfloor $，于是每次递归的时候累加一下即可。

预处理复杂度上界是 \(O(p)\) 的，因为因数一定更小。单次询问的复杂度是 \(O(\sum\log_{p_i}n)\)。

code

封装了一下。注意先调用 init 初始化再询问。同时注意剪枝。具体的，判断 \(p_i^{x-y-z}\pmod p\) 是否是 0。因为有很大概率是 0，因此剪枝效果非常明显，甚至快了接近 10 倍。

点击查看代码

namespace exlucas{
    #define int long long
    void exgcd(int a,int b,int &x,int &y,int p){
        // if(a<b)swap(a,b),swap(x,y);
        if(b==0){x=1,y=0;return;}
        int x1,y1;exgcd(b,a%b,x1,y1,p);
        x=y1,y=(x1-(a/b)*y1+p)%p;
    }
    vector <int> q[250];
    void init(int p){
        int P=p,cnt=0;
        for(int i=2;i<=P;i++){
            if((P==1)||P%i)continue;
            cnt++;vector <int> ().swap(q[cnt]);int mod=1;while(P>1&&P%i==0)mod*=i,P/=i;
            q[cnt].resize(mod+5);q[cnt][0]=1;for(int j=1;j<=mod;j++){q[cnt][j]=q[cnt][j-1];if(j%i)q[cnt][j]=q[cnt][j]*j%mod;}
        }
    }
    int ksm(int x,int k,int p){
        int res=1;x%=p;while(k){
            if(k&1)res=res*x%p;
            x=x*x%p;k>>=1;
        }return res;
    }
    int pre=1,ci,id;
    int getF(int n,int P,int p){
        if(n==0)return 1;
        int res=getF(n/p,P,p)%P,tmp=ksm(pre,n/P,P);
        if(n%P>0)tmp=tmp*q[id][n%P]%P;
        return res*tmp%P;
    }
    int getG(int n,int p){if(n<p)return 0;return (n/p)+getG(n/p,p);}
    int C(int n,int m,int p){
        int P=p,ans=0,cnt=0;
        for(int i=2;i<=P;i++){
            if((P==1)||P%i)continue;
            int mod=1;cnt++,id=cnt;while(P>1&&P%i==0)mod*=i,P/=i;
            pre=q[cnt][mod];
            ci=ksm(i,getG(n,i)-getG(m,i)-getG(n-m,i),mod);
            if(!ci){continue;}int a=getF(n,mod,i),e;int b=getF(m,mod,i),c=getF(n-m,mod,i);
            exgcd(b,mod,b,e,mod),exgcd(c,mod,c,e,mod);a=a*b%mod*c%mod*ci%p; //在算出真正的 a 前都是 \mod mod 
            e=p/mod;a=a*e%p;exgcd(e,mod,e,b,mod);a=a*e%p;
            ans=(ans+a+p)%p;
        }
        return ans;
    }
    #undef int
}