您的位置:首页 > 其它

HDU5446 Unknown Treasure Lucas+中国剩余定理

2015-09-16 17:11 204 查看
题目链接:http://acm.hdu.edu.cn/showproblem.php?pid=5446

题目大意:求C(n,m) mod M的值,其中M=p1*p2*...*pk<=10^18且pi为素数。

分析:首先我们由Lucas定理可以解决大组合数取模的问题,即,我们先求出Ai=C(n,m) mod pi,然后得到k组线性同余方程组,接下来中国剩余定理解出方程组的解即可。需要注意的是,由于n和m都是10^18级的,long
long乘法会溢出,要用到快速乘。

实现代码如下:

#include <cstdio>
#include <iostream>
using namespace std;
typedef long long ll;
const int maxn=15;
ll prime[maxn]; //M=∑p
ll A[maxn]; //ans[i]=C(n,m)%pi

ll multi(ll a,ll b,ll m)
{//快速乘
ll ans=0;
while(b)
{
if(b&1) ans=(ans+a)%m;
b>>=1;
a=(a+a)%m;
}
return ans;
}
ll quick_mod(ll a,ll b,ll m)
{
ll ans=1;
a%=m;
while(b)
{
if(b&1) ans=ans*a%m;
b>>=1;
a=a*a%m;
}
return ans;
}
ll C(ll n,ll m,int cnt)
{
ll p=prime[cnt];
if(m>n) return 0;
if(m>n-m) m=n-m;
ll ans=1;
for(ll i=1;i<=m;i++)
{
ll a=(n+i-m)%p;
ll b=i%p;
ans=multi(ans, multi(a,quick_mod(b,p-2,p),p) ,p);
}
return ans%p;
}
ll lucas(ll n,ll m,int cnt)
{
int p=prime[cnt];
if(m==0) return 1;
return multi( C(n%p,m%p,cnt) , lucas(n/p,m/p,cnt) , p );
}

//求n组线性同余方程x=r(mod a)的解
ll a,b,c,d,a1,r1,a2,r2,x0,y0;
void exgcd(ll a,ll b,ll &d,ll &x,ll &y)
{
if(!b)
{
x=1,y=0,d=a;
return ;
}
else
{
exgcd(b,a%b,d,x,y);
ll temp=x;
x=y;
y=temp-(a/b)*y;
}
}
ll solve(int n)
{
ll mod=1,x=0;
for(int i=0;i<n;i++)
mod*=prime[i];
for(int i=0;i<n;i++)
{
ll tmp=mod/prime[i];
exgcd(prime[i],tmp,d,x0,y0);x0=d;
x=(x+multi( multi(y0,tmp,mod),A[i],mod) )%mod;

}
return (x+mod)%mod;
}
int main()
{
int t,k;
ll m,n;
scanf("%d",&t);
while(t--)
{
scanf("%I64d%I64d%d",&n,&m,&k);
for(int i=0;i<k;i++)
{
scanf("%I64d",&prime[i]);
A[i]=lucas(n,m,i)%prime[i];
}
printf("%I64d\n",solve(k));
}
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: