lucas求组合数C(n,k)%p

showproblem.php?pid=3037

 1 #include<cstdio>
 2 typedef __int64 LL;
 3 const int M=100010;
 4 class LUCAS { //lucas求组合数C(n,k)%p
 5     LL F[M];
 6     LL inv(LL a,LL mod) {
 7         if(a==1) return 1;
 8         return inv(mod%a,mod)*(mod-mod/a)%mod;
 9     }
10     void init(LL p) {
11         F[0]=1;
12         for(int i=1; i<=p; i++) {
13             F[i]=F[i-1]*i%p;
14         }
15     }
16     LL Lucas(LL n,LL k,LL p) {
17         LL ans=1;
18         while(n&&k) {
19             LL a=n%p;
20             LL b=k%p;
21             if(a<b) return 0;
22             ans=ans*F[a]%p*inv(F[b]*F[a-b]%p,p)%p;
23             n/=p;
24             k/=p;
25         }
26         return ans;
27     }
28 public:
29     LL solve(LL n,LL k,LL p){
30         init(p);
31         return Lucas(n,k,p);
32     }
33 }gx;
34 int main() {
35     int t,n,m,p;
36     while(~scanf("%d",&t)) {
37         while(t--) {
38             scanf("%d%d%d",&n,&m,&p);
39             printf("%d
",(int)gx.solve(n+m,n,p));
40         }
41     }
42     return 0;
43 }

View Code

DP? http://acm.hdu.edu.cn/showproblem.php?pid=3944

 1 #include<cstdio>
 2 #include<cstring>
 3 #include<vector>
 4 #define mt(a,b) memset(a,b,sizeof(a))
 5 using namespace std;
 6 typedef __int64 LL;
 7 const int M=10010;
 8 class LUCASM { //比较适合n,k<=10^9,p<=10^4的情况
 9     int flag[M*10],prime[1300],pcnt;
10     vector<int> rev[M],fac[M];
11     int quickpow(int a,int b,int c) { //快速幂求(a^b)%c
12         int ret=1%c;
13         for(; b; a=a*a%c,b>>=1) {
14             if(b&1) {
15                 ret=ret*a%c;
16             }
17         }
18         return ret;
19     }
20 public:
21     void init() {//先初始化一次即可
22         pcnt=0;
23         mt(flag,-1);
24         for(int i=2; i<=10000; i++) {
25             if(flag[i]) {
26                 prime[pcnt++]=i;
27                 rev[i].clear();
28                 fac[i].clear();
29             }
30             for(int j=0; j<pcnt&&prime[j]<=10000/i; j++) {
31                 flag[i*prime[j]]=0;
32                 if(!(i%prime[j])) break;
33             }
34         }
35         for(int i=0; i<pcnt; i++) {
36             int tnum=1;
37             rev[prime[i]].push_back(1);
38             fac[prime[i]].push_back(1);
39             for (int j=1; j<prime[i]; j++) {
40                 tnum=(tnum*j)%prime[i];
41                 int now=quickpow(tnum,prime[i]-2,prime[i]);
42                 fac[prime[i]].push_back(tnum);
43                 rev[prime[i]].push_back(now);
44             }
45         }
46     }
47     int lucas(int n,int k,int p) {
48         int ret=1;
49         while (n && k) {
50             int num1=n%p;
51             int num2=k%p;
52             n/=p;
53             k/=p;
54             if (num1<num2) return 0;
55             int num=(fac[p][num1]*rev[p][num2])%p;//计算c(num1,num2)%p
56             num=(num*rev[p][num1-num2])%p;
57             ret=(ret*num)%p;
58         }
59         return ret;
60     }
61 } gx;
62 int main() {
63     int n,k,p,cas=1;
64     gx.init();
65     while(~scanf("%d%d%d",&n,&k,&p)) {
66         printf("Case #%d: ",cas++);
67         if(k>n/2) k=n-k;
68         int o=gx.lucas(n+1,k,p);
69         printf("%d
",(n-k+o)%p);
70     }
71     return 0;
72 }

View Code

end