[破碎的状态] [-21] 100633 J
题意:
求[tex]C_n^k\ mod\ m[/tex]
[tex]n<=10^{18},0<=k<=n,m<=10^6[/tex]
解法:
众所周知..
[tex]C_n^k = \frac{n!}{k! (n-k)!}[/tex]
我们根据中国剩余定理..
首先..拆成若干个[tex]p^k[/tex]的形式
对于每个[tex]p^k[/tex],我们只要对于n,k和n-k知道(以下均以n为例):
1,[tex] n! = p^{a} * q[/tex]中a的值
2,[tex]q\ mod\ p^k[/tex]的值
这时候我们只要把这个a和q'和剩下两个值进行下操作,实际上是
[tex]a_{sum} = a_1 - a_2 - a_3[/tex]
[tex]q_{sum} = \frac {q_1}{q_2 q_3} = q_1 * q_2^{\phi(p^k)-1} * q_3^{\phi(p^k)-1}[/tex]
那么..
a我们可以通过不断除的办法来求..而q..
我们可以发现....
对于每个数x,[/tex]gcd(x,p)=1[/tex],则直接乘以[tex]x%p[/tex]
对于每个数[tex]x = p * a[/tex],我们要乘以[tex]x/p^2 % p[/tex]
对于每个数[tex]x = p^2 * a[/tex],我们要乘以[tex]x/p^2 % p[/tex]
......
我们可以把所有gcd(x,p)=1的给找出来,乘以一下,因为值太多,我们可以统计有多少个x % p = 1, x % p = 2 .....
这个数字是n/p或者n/p+1...直接乘以下就好了..
那么把剩下的数字/p,就转化回了原问题
考虑到我不会中国剩余定理的快速求的姿势..我暴力了一发.....写的不太好看
综上..
代码如下:
#include<set> #include<map> #include<list> #include<queue> #include<stack> #include<string> #include<time.h> #include<math.h> #include<memory> #include<vector> #include<bitset> #include<fstream> #include<stdio.h> #include<utility> #include<sstream> #include<string.h> #include<iostream> #include<stdlib.h> #include<algorithm> using namespace std; int modo; int phi_modo; long long power(long long x,long long y) { if (y==0) return 1ll; long long t=power(x,y/2); t=t*t%modo; if (y%2==1) t=t*x%modo; return t; } int facts[1000005]; struct result { long long modo_ans; long long cnt; result () { modo_ans=1; cnt=0; } friend result operator * (const result &a,const result &b) { result c; c.modo_ans=(a.modo_ans*b.modo_ans)%modo; c.cnt=a.cnt+b.cnt; return c; } friend result operator / (const result &a,const result &b) { result c; c.modo_ans=(a.modo_ans*power(b.modo_ans,phi_modo-1))%modo; c.cnt=a.cnt-b.cnt; return c; } }; result fact(long long x,int p) { if (x==0) { return result(); } if (x<modo) { result a; a.modo_ans=1; a.modo_ans=facts[x]; a.cnt=x/p; return a*fact(x/p,p); } else { result a; a.modo_ans=1; int t=x%modo; a.modo_ans=facts[t]; long long temp=facts[modo]; a.cnt=x/p; a.modo_ans=a.modo_ans*power(temp,x/modo)%modo; return a*fact(x/p,p); } } int main() { #ifdef absi2011 freopen("input.txt","r",stdin); freopen("output.txt","w",stdout); #endif long long n,m,k; cin>>n>>m>>k; int i; int ans=1; int modo_ans=0; for (i=2;i<=k;i++) { if (k%i==0) { modo=1; int sum=0; for (;k%i==0;) { sum++; k/=i; modo*=i; } phi_modo=modo/i*(i-1); int j; facts[0]=1; for (j=1;j<=modo;j++) { facts[j]=facts[j-1]; if (j%i!=0) facts[j]=facts[j-1]*(long long)j%modo; } result a=fact(n,i)/fact(m,i)/fact(n-m,i); if (a.cnt>=sum) { a.cnt=0; a.modo_ans=0; } else { a.modo_ans=a.modo_ans*power(i,a.cnt)%modo; } for (j=0;j<modo;j++) { if ((ans*(long long)j+modo_ans)%modo==a.modo_ans) break; } modo_ans=(ans*(long long)j+modo_ans)%(ans*modo); ans*=modo; } } printf("%d\n",modo_ans); return 0; }