您的位置:首页 > 其它

HDU 5446——Unknown Treasure——————【CRT+lucas+exgcd+快速乘+递推求逆元】

2015-09-15 20:33 337 查看
Each test case starts with three integers n,m,k(1≤m≤n≤1018,1≤k≤10) on a line where k is the number of primes. Following on the next line are kdifferent primes p1,...,pk. It is guaranteed that M=p1⋅p2⋅⋅⋅pk≤1018 and pi≤105 for every i∈{1,...,k}.

[align=left]Output[/align]
For each test case output the correct combination on a line.

[align=left]Sample Input[/align]

1

9 5 2

3 5

[align=left]Sample Output[/align]

6

[align=left]Source[/align]
2015 ACM/ICPC Asia Regional Changchun Online

题目大意:让求C(n,m)%(∏pi) 这个式子的值。





中国剩余定理:



解题思路:首先用lucas定理将求C(a,b)%p转化成求解∏C(bi,ai),这样,我们可以得到c[i]数组。然后用中国剩余定理来求x0的值,即为答案。在求解的过程中需要用到扩展欧几里得来求解Mi的逆元,由于Mi比较大,所以在乘积的时候会爆数据范围,所以改成快速乘取模的方式代替直接乘积。

#include<bits/stdc++.h>
using namespace std;
typedef long long INT;
const int maxp=1e5+200;
INT p[15],c[15];
INT fac[maxp],inv[maxp];
INT powmod(INT a,INT n,INT mod){//快速幂取模
INT ret=1;
while(n){
if(n&1){
ret=ret*a%mod;
}
n>>=1;
a = a*a%mod;
}
return ret;
}
INT mulmod(INT a,INT b,INT mod){//快速乘取模
a = (a%mod + mod) % mod;    //用扩展欧几里得求出的值可能为负值
b = (b%mod + mod) % mod;    //用扩展欧几里得求出的值可能为负值
INT ret=0;
while(b){
if(b&1){
ret = (ret+a)%mod;
}
b >>= 1;
a = (a<<1) % mod;
}
return ret;
}
void init(INT n){   //递推出来阶乘和逆元数组
fac[0]=1;
for(int i=1;i<n;i++){
fac[i]=fac[i-1]*i % n;
}
inv[n-1]=powmod(fac[n-1],n-2,n);
for(int i=n-2;i>=0;i--){
inv[i] = inv[i+1] * (i+1) % n;
//fac
*inv[fac
]≡1%p  ==>  fac[n-1]*(n*inv[fac
])≡1%p
}
}
INT cm(INT n,INT m,INT mod){    //用逆元求组合数取模
if(n<0||m<0||m>n){
return 0;
}
return fac
*inv[n-m]%mod*inv[m]%mod;
}
INT lucas(INT n,INT m,INT mod){//lucas递归求P进制时的c
if(m==0){
return 1;
}
return lucas(n/mod,m/mod,mod) * cm(n%mod,m%mod,mod) % mod;
}
INT exgcd(INT a,INT b,INT &x,INT &y){   //求b关于模a的逆元。放在y中
if(b==0) { x = 1; y = 0; return a; }
INT d = exgcd(b, a%b , y, x);
y -= x * (a / b);
return d;
}
void CRT(INT k){//中国剩余定理求解一元线性同余方程组
INT M=1,x,y;
INT ans=0;
for(int i=1;i<=k;i++){
M *= p[i];
}
for(int i=1;i<=k;i++){
INT Mi=M/p[i];
exgcd(p[i],Mi,x,y);
ans = (ans+mulmod(mulmod(y,Mi,M),c[i],M))%M                                             ;
}
printf("%I64d\n",ans);
}
int main(){
INT n,m,k;
int t;
scanf("%d",&t);
while(t--){
scanf("%I64d%I64d%I64d",&n,&m,&k);
for(int i=1;i<=k;i++){
scanf("%I64d",&p[i]);
init(p[i]);
c[i] = lucas(n,m,p[i]);
}
CRT(k);
}
return 0;
}


  
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: