zcmimi's blog
查看原题

点击跳转

p=p_1^{k_1}p_2^{k_2}\dots p_n^{k_n}

求出每个{n\choose m} \equiv a_i \pmod {p_i^{k_i}}

得到同余方程组

\begin{cases} {n\choose m} \equiv a_1 \pmod {p_1^{k_1}} \\ {n\choose m} \equiv a_2 \pmod {p_2^{k_2}} \\ \vdots \\ {n\choose m} \equiv a_n \pmod {p_n^{k_n}} \end{cases}

使用crt即可求出n\choose m的值

p^t不是质数,那么如何求{n\choose m} \mod p^t呢?

可以计算{n\choose m}=\frac{n!}{m!(n-m)!},分别计算出n!,m!,(n-m)!在模p^t意义下的值

假设p=3,t=2,n=19,

n!=1\cdot 2\cdot 3\cdots 19 \\ =(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8\cdot 10\cdot 11\cdot 13\cdot 14\cdot 16\cdot 17\cdot 19)\cdot (3\cdot 6\cdot 9\cdot 12\cdot 15\cdot 18) \\ =(1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8\cdot 10\cdot 11\cdot 13\cdot 14\cdot 16\cdot 17\cdot 19)\cdot 3^6\cdot(1\cdot2\cdot3\cdot4\cdot5\cdot6)

可以发现后面的部分相当于\left\lfloor \frac np \right\rfloor!,递归计算即可

前半部分以p^t为周期,即

1\cdot 2\cdot 4\cdot 5\cdot 7\cdot 8 \equiv 10\cdot 11\cdot 13\cdot 14\cdot 16\cdot 17 \pmod {3^2}

只需要计算不满足周期的数19就可以了

#include<bits/stdc++.h>
typedef long long ll;
int lim;
ll P;
ll pw(ll x,ll b,ll p){
    ll res=1;
    while(b){
        if(b&1)res=res*x%p;
        b>>=1,x=x*x%p;
    }
    return res;
}
void exgcd(ll a,ll b,ll &x,ll &y){
    if(!b){x=1,y=0;return;}
    exgcd(b,a%b,x,y);
    ll t=x;x=y,y=t-a/b*y;
}
ll inv(ll v,ll p){
    ll x,y;
    exgcd(v,p,x,y);
    return (x+=p)>p?x-p:x;
}
ll crt(ll v,ll p){return v*inv(P/p,p)%P*(P/p)%P;}
ll fac(ll n,ll p,ll pk){
    if(!n)return 1;
    ll res=1;
    for(int i=2;i<=pk;++i)
        if(i%p)res=res*i%pk;
    res=pw(res,n/pk,pk);
    for(int i=2;i<=n%pk;++i)
        if(i%p)res=res*i%pk;
    return res*fac(n/p,p,pk)%pk;
}
ll C(ll n,ll m,ll p,ll pk){
    ll cnt=0;
    for(ll i=n;i;i/=p)cnt+=i/p;
    for(ll i=m;i;i/=p)cnt-=i/p;
    for(ll i=n-m;i;i/=p)cnt-=i/p;
    return pw(p,cnt,pk)*fac(n,p,pk)%pk*inv(fac(m,p,pk),pk)%pk*inv(fac(n-m,p,pk),pk)%pk;
}
ll exlucas(ll n,ll m){
    ll res=0,t=P,pk;
    for(int i=2;i<=lim;++i)
    if(t%i==0){
        pk=1;while(t%i==0)pk*=i,t/=i;
        res=(res+crt(C(n,m,i,pk),pk))%P;
    }
    if(t>1)res=(res+crt(C(n,m,t,t),t))%P;
    return res;
}
int main(){
    ll n,m;
    scanf("%lld%lld%lld",&n,&m,&P);lim=sqrt(P)+1;
    printf("%lld\n",exlucas(n,m));
}
LG 4720 【模板】扩展卢卡斯
comment评论
Search
search