HDU 5446 Unknown Treasure(lucas + 中国剩余定理 + 模拟乘法)

题目链接： http://acm.hdu.edu.cn/showproblem.php?pid=5446
题目大意：求C(n, m) % M， 其中M为不同素数的乘积,即M=p1*p2*...*pk, 1≤k≤10。1≤m≤n≤10^18。

分析： 如果M是素数，则可以直接用lucas定理来做，但是M不是素数，而是素数的连乘积。令C(n, m)为 X ,则可以利用lucas定理分别计算出 X%p1，X%p2, ... , X % pk的值，然后用中国剩余定理来组合得到所求结果。

比较坑的地方是，中间结果会超long long 范围，要用类似于快速幂的方法来模拟。

参考代码：

#include <iostream>
#include <stdio.h>
#include <string.h>
using namespace std;
long long a[110001];
long long Mod;
long long c[110001];
long long w[1000];
long long Egcd(long long a,long long b,long long &x,long long &y){
long long d;
if(b==0){
x=1;y=0;
return a;
}
d=Egcd(b,a%b,y,x);
y-=a/b*x;
return d;
}
long long cal(long long a, long long b, long long mod) {//模拟a*b%mod
long long ret = 0;
while(b) {
if(b & 1) ret = (ret + a) % mod;
a = (a + a) % mod;
b >>= 1;
}
return ret;
}

long long c_r(int len){//中国剩余定理
int i;
long long d,x,y,n,m,ret;
ret=0;
n=1;
for(i=0;i<len;i++)
n*=c[i];
for(i=0;i<len;i++){
m=n/c[i];
d=Egcd(m,c[i],y,x);
y=(c[i]+y%c[i])%c[i];
if( i&1) y -= c[i];
//if(y>c[i]-y) y=y-c[i];
ret=(ret+cal(y*m, w[i], n))%n;
}
return (n+ret%n)%n;
}

void init(){
int i;
memset(a, 0, sizeof(a));
a[0]=1;
for(i=1;i<Mod;i++)
a[i]=a[i-1]*i%Mod;
}

long long gcd(long long a,long long b){
if(b==0) return a;
return gcd(b,a%b);
}

long long choose(long long n,long long m){
if(m>n)return 0;
else if(n==m) return 1;
long long nn=a
,mm=a[m]*a[n-m]%Mod;
long long d=gcd(nn,mm);
nn/=d;
mm/=d;
long long x,y;
Egcd(mm,Mod,x,y);
x=(x+Mod)%Mod;
return (x*nn)%Mod;
}

long long work(long long n,long long m){//lucas
long long ret=1;
while(n&&m){
ret*=choose(n%Mod,m%Mod);
ret%=Mod;
n/=Mod,m/=Mod;
}
return ret;
}

int main(){
int T;
long long n,m;
int k;
scanf("%d",&T);
while(T--){
cin >> n >> m >> k;
if(m > n-m) m = n-m;
for(int i=0;i<k;i++){
cin >> c[i];
Mod=c[i];
init();
w[i]=work(n,m);
}
cout << c_r(k) << endl;
}
}