自动搬运

查看原文

来自洛谷，原作者为

	Fading AFO

搬运于2025-08-24 22:02:57，当前版本为作者最后更新于2019-07-05 10:52:08，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

扩展 Lucas 定理 exLucas

其他题解讲的都很不清楚，而且很奇怪。。。

所以我来发一篇良心题解，求大家点赞~

问题

求

$$C_n^m \bmod{p}$$

其中$p\leq 10^6$

算法

设

$$p=p_1^{\alpha_1}p_2^{\alpha_2}...$$

求出

$$\left\{\begin{aligned}C_n^m \bmod\ {p_1^{\alpha_1}} \\C_n^m \bmod\ {p_2^{\alpha_2}}\\ ......\end{aligned} \right. $$

最后用中国剩余定理合并即可。

假设我们现在要求

$$C_n^m\bmod{P^k}([P\in {\texttt{prime}}])$$ $$\therefore\ =\frac {n!}{m!(n-m)!}\bmod{P^k}$$

我们无法求得$m!,(n-m)!$的逆元，理由很简单，不一定有解。

怎么办呢？发现$a$对$p$有逆元的充要条件为

$$\gcd(a,p)=1$$

所以我们换一个式子：

$$=\frac {\frac {n!}{P^{x}}}{\frac {m!}{P^{y}}\frac {(n-m)!}{P^{z}}}P^{x-y-z}\bmod{P^k} $$

其中

$x$为$n!$中$P$因子的个数（包含多少$P$因子）。

其他都是同理的。

所以如果我们对于一个$n$，可以求出$\frac {n!}{P^x}$，那不就完事了吗？（这样我们就可以求逆元了）

问题等价于求

$$\frac {n!}{P^x}\bmod{P^k}$$

我们对$n!$进行变形：

$$n!=1\cdot 2\cdot 3\cdot ...\cdot n$$ $$=(P\cdot 2P\cdot 3P...)(1\cdot 2\cdot ...)$$

左边是$1\sim n$中所有$\leq n$的$P$的倍数。

右边反之。

因为$1\sim n$中有$\lfloor \frac nP\rfloor$个$P$的倍数，

$$\therefore\ =P^{\lfloor \frac nP\rfloor}(1 \cdot 2 \cdot 3...)(1\cdot 2\cdot ...) $$$$=P^{\lfloor \frac nP\rfloor}(\lfloor \frac nP\rfloor)!\prod_{i=1,i\not\equiv 0\pmod {P}}^ni $$

显然后面这个$\bmod\ P$是有循环节的。

$$=P^{\lfloor \frac nP\rfloor}(\lfloor \frac nP\rfloor)!(\prod_{i=1,i\not\equiv 0\pmod {P}}^{P^k}i)^{\lfloor \frac n{P^k}\rfloor}(\prod_{i=P^k\lfloor \frac n{P^k}\rfloor,i\not\equiv 0\pmod {P}}^ni) $$

其中

$$(\prod_{i=1,i\not\equiv 0\pmod {P}}^{P^k}i)^{\lfloor \frac n{P^k}\rfloor} $$

是循环节$1\sim P$中所有无$P$因子数的乘积。

$$(\prod_{i=P^k\lfloor \frac n{P^k}\rfloor,i\not\equiv 0\pmod {P}}^ni) $$

是余项的乘积。

比如

$$22!\equiv(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8)(10\cdot 11\cdot 13\cdot 14\cdot 16\cdot 17)(19\cdot 20\cdot 22)(3\cdot 6\cdot 9\cdot 12\cdot 15\cdot 18\cdot 21)\pmod {3^2} $$$$\equiv(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8)^2(19\cdot 20\cdot 22)3^7(1\cdot 2\cdot 3\cdot 4\cdot 5\cdot 6\cdot 7)\pmod {3^2} $$$$\equiv3^77!(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8)^2(19\cdot 20\cdot 22)\pmod {3^2} $$

发现正好是一样的。$\lfloor \frac {22}{3^2}\rfloor=2$

理解了吧...

$$=P^{\lfloor \frac nP\rfloor}(\lfloor \frac nP\rfloor)!(\prod_{i=1,i\not\equiv 0\pmod {P}}^{P^k}i)^{\lfloor \frac n{P^k}\rfloor}(\prod_{i=P^k\lfloor \frac n{P^k}\rfloor,i\not\equiv 0\pmod {P}}^ni) $$

发现前面这个$P^{\lfloor \frac nP\rfloor}$是要除掉的。

然而$(\lfloor \frac nP\rfloor)!$里可能还包含$P$。

所以，我们定义函数：

$$f(n)=\frac {n!}{P^x}$$ $$\therefore\ f(n)=f(\lfloor \frac nP\rfloor)(\prod_{i=1,i\not\equiv 0\pmod {P}}^{P^k}i)^{\lfloor \frac n{P^k}\rfloor}(\prod_{i=P^k\lfloor \frac n{P^k}\rfloor,i\not\equiv 0\pmod {P}}^ni) $$

这样就好了。时间复杂度是$O(\log_{P}n)$。

边界$f(0)=1$

看看原式

$$=\frac {\frac {n!}{P^{x}}}{\frac {m!}{P^{y}}\frac {(n-m)!}{P^{z}}}P^{x-y-z}\bmod{P^k} $$$$=\frac {f(n)}{f(m)f(n-m)}P^{x-y-z}\bmod{P^k}$$

由于$f(m)$一定与$P^k$互质，所以我们可以直接用$exgcd$求逆元啦。

最后我们还有一个问题，$P^{x-y-z}$咋求?

其实很好求。

比如说，我们要求$f(n)=\frac {n!}{P^x}$中的$x$

设$g(n)=x$(上述的$x$)

再看看阶乘的式子

$$n!=P^{\lfloor \frac nP\rfloor}(\lfloor \frac nP\rfloor)!(\prod_{i=1,i\not\equiv 0\pmod {P}}^{P^k}i)^{\lfloor \frac n{P^k}\rfloor}(\prod_{i=P^k\lfloor \frac n{P^k}\rfloor,i\not\equiv 0\pmod {P}}^ni) $$

很显然我们要的就是$P^{\lfloor \frac nP\rfloor}$

而且$(\lfloor \frac nP\rfloor)!$可能还有$P$因子。

所以我们可以得到以下递推式

$$g(n)=\lfloor \frac nP\rfloor+g(\lfloor \frac nP\rfloor) $$

时间复杂度是$O(\log_{P}n)$。

边界$g(n)=0(n<P)$

所以答案就是

$$=\frac {\frac {n!}{P^{x}}}{\frac {m!}{P^{y}}\frac {(n-m)!}{P^{z}}}P^{g(n)-g(m)-g(n-m)}\bmod{P^k} $$

最后用中国剩余定理合并答案即可得到$C_n^m$。

代码:

#include<bits/stdc++.h>
#define ll long long
using namespace std;
#ifndef Fading
inline char gc(){
    static char now[1<<16],*S,*T;
    if (T==S){T=(S=now)+fread(now,1,1<<16,stdin);if (T==S) return EOF;}
    return *S++;
}
#endif
#ifdef Fading
#define gc getchar
#endif
void exgcd(ll a,ll b,ll &x,ll &y){
    if (!b) return (void)(x=1,y=0);
    exgcd(b,a%b,x,y);
    ll tmp=x;x=y;y=tmp-a/b*y;
}
ll gcd(ll a,ll b){
    if (b==0) return a;
    return gcd(b,a%b); 
}
inline ll INV(ll a,ll p){
    ll x,y;
    exgcd(a,p,x,y);
    return (x+p)%p;
}
inline ll lcm(ll a,ll b){
    return a/gcd(a,b)*b;
}
inline ll mabs(ll x){
    return (x>0?x:-x);
}
inline ll fast_mul(ll a,ll b,ll p){
    ll t=0;a%=p;b%=p;
    while (b){
        if (b&1LL) t=(t+a)%p;
        b>>=1LL;a=(a+a)%p;
    }
    return t;
}
inline ll fast_pow(ll a,ll b,ll p){
    ll t=1;a%=p;
    while (b){
        if (b&1LL) t=(t*a)%p;
        b>>=1LL;a=(a*a)%p;
    }
    return t;
}
inline ll read(){
    ll x=0,f=1;char ch=gc();
    while (!isdigit(ch)) {if (ch=='-') f=-1;ch=gc();}
    while (isdigit(ch)) x=x*10+ch-'0',ch=gc();
    return x*f;
}
inline ll F(ll n,ll P,ll PK){
    if (n==0) return 1;
    ll rou=1;//循环节
    ll rem=1;//余项 
    for (ll i=1;i<=PK;i++){
        if (i%P) rou=rou*i%PK;
    }
    rou=fast_pow(rou,n/PK,PK);
    for (ll i=PK*(n/PK);i<=n;i++){
        if (i%P) rem=rem*(i%PK)%PK;
    }
    return F(n/P,P,PK)*rou%PK*rem%PK;
}
inline ll G(ll n,ll P){
    if (n<P) return 0;
    return G(n/P,P)+(n/P);
}
inline ll C_PK(ll n,ll m,ll P,ll PK){
    ll fz=F(n,P,PK),fm1=INV(F(m,P,PK),PK),fm2=INV(F(n-m,P,PK),PK);
    ll mi=fast_pow(P,G(n,P)-G(m,P)-G(n-m,P),PK);
    return fz*fm1%PK*fm2%PK*mi%PK;
}
ll A[1001],B[1001];
//x=B(mod A)
inline ll exLucas(ll n,ll m,ll P){
    ll ljc=P,tot=0;
    for (ll tmp=2;tmp*tmp<=P;tmp++){
        if (!(ljc%tmp)){
            ll PK=1;
            while (!(ljc%tmp)){
                PK*=tmp;ljc/=tmp;
            }
            A[++tot]=PK;B[tot]=C_PK(n,m,tmp,PK);
        }
    }
    if (ljc!=1){
        A[++tot]=ljc;B[tot]=C_PK(n,m,ljc,ljc);
    }
    ll ans=0;
    for (ll i=1;i<=tot;i++){
        ll M=P/A[i],T=INV(M,A[i]);
        ans=(ans+B[i]*M%P*T%P)%P;
    }
    return ans;
}
signed main(){
    ll n=read(),m=read(),P=read();
    printf("%lld\n",exLucas(n,m,P));
    return 0;
}