[cf2015ICLFinalsDiv1J]Ceizenpok’s formula

题意:$C_n^m\% k$

解题关键:扩展lucas+中国剩余定理裸题

 1 #include<algorithm>
 2 #include<iostream>
 3 #include<cstring>
 4 #include<cstdio>
 5 #include<cmath>
 6 #include<cstdlib>
 7 typedef long long ll;
 8 using namespace std;
 9 ll mod,n,m,x,y,module[10002],piset[10002],r[10002];
10 
11 ll mod_pow(ll x,ll n,ll p){
12     ll res=1;
13     while(n){
14         if(n&1) res=res*x%p;
15         x=x*x%p;
16         n>>=1;
17     }
18     return res;
19 }
20 
21 ll extgcd(ll a,ll b,ll &x,ll &y){
22     ll d=a;
23     if(b)  d=extgcd(b,a%b,y,x),y-=a/b*x;
24     else x=1,y=0;
25     return d;
26 }
27 
28 ll inv(ll t,ll mod){ extgcd(t,mod,x,y);return (x+mod)%mod;}
29 
30 ll multi(ll n,ll pi,ll pk){//求非互质的部分 
31     if (!n) return 1;
32     ll ans=1;
33     for (ll i=2;i<=pk;i++) if(i%pi) ans=ans*i%pk;
34     ans=mod_pow(ans,n/pk,pk);
35     for (ll i=2;i<=n%pk;i++) if(i%pi) ans=ans*i%pk;
36     return ans*multi(n/pi,pi,pk)%pk;
37 }
38 
39 
40 ll exlucas(ll n,ll m,ll pi,ll pk){//组合数 c(n,m)mod pk=pi^k 
41     if(m>n) return 0;
42     ll a=multi(n,pi,pk),b=multi(m,pi,pk),c=multi(n-m,pi,pk);
43     ll k=0;
44     for(ll i=n;i;i/=pi) k+=i/pi;
45     for(ll i=m;i;i/=pi) k-=i/pi;
46     for(ll i=n-m;i;i/=pi) k-=i/pi;
47     return a*inv(b,pk)%pk*inv(c,pk)%pk*mod_pow(pi,k,pk)%pk;//组合数求解完毕 
48 }
49 
50 ll crt(int n,ll *r,ll *m){
51     ll M=1,ret=0;
52     for(int i=0;i<n;i++) M*=m[i];
53     for(int i=0;i<n;i++){
54         ll w=M/m[i];
55         ret+=w*inv(w,m[i])*r[i];
56         ret%=M;
57     }
58     return (ret+M)%M;
59 }
60 
61 ll fz(ll n,ll *m,ll *piset){//分解质因子 
62     ll num=0;
63     for (ll i=2;i*i<=n;i++){
64         if(n%i==0){
65             ll pk=1;
66             while(n%i==0) pk*=i,n/=i;
67             m[num]=pk;
68             piset[num]=i;
69             num++;
70         }
71     }
72     if(n>1) m[num]=n,piset[num]=n,num++;
73     return num;
74 }
75 
76 ll excomb(ll n,ll m,ll p){
77     ll num=fz(p,module,piset);
78     for(int i=0;i<num;i++){
79         r[i]=exlucas(n,m,piset[i],module[i]);
80     }
81     return crt(num,r,module);
82 }
83 
84 int main(){
85     cin>>n>>m>>mod;
86     printf("%d",excomb(n,m,mod));
87     return 0;
88 }

 

posted @ 2017-10-02 17:05  Elpsywk  阅读(246)  评论(0编辑  收藏  举报