BZOJ 2142 礼物 数论

这道题是求组合数终极版.

C(n,m) mod P

n>=1e9 m>=1e9 P>=1e9且为合数且piqi<=1e5

拓展lucas定理.

实际上就是一点数论小知识的应用.

这篇文章对于CRT和lucas定理的学习非常不错.

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define FILE "dealing"
#define up(i,j,n) for(int i=j;i<=n;i++)
#define db double
#define pii pair<ll,ll>
#define pb push_back
template<class T> inline bool cmin(T& a,T b){return a>b?a=b,true:false;}
template<class T> inline bool cmax(T& a,T b){return a<b?a=b,true:false;}
template<class T> inline T squ(T a){return a*a;}
const int maxn=101000+10,inf=1e9+10;
ll read(){
ll x=0,f=1,ch=getchar();
while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+ch-'0',ch=getchar();
return x*f;
}
ll n,m,mod,w[maxn];
struct node{
ll p,a,num;
}a[maxn];int tot=0;
void divide(ll n){
for(ll i=2;i*i<=n;i++){
if(n%i==0){
n/=i;
a[++tot].p=i,a[tot].a=1,a[tot].num=i;
while(n%i==0){
a[tot].a++;
a[tot].num*=i;
n/=i;
}
}
if(n==1)break;
}
if(n!=1)a[++tot].p=n,a[tot].num=n,a[tot].a=1;
}
ll qpow(ll a,ll b,ll mod){
ll ans=1;
for(;b;b>>=1,a=a*a%mod)
if(b&1)ans=ans*a%mod;
return ans;
}
pii deal(ll n,int pos){
pii t;t.first=1,t.second=0;
if(n==0)return t;
ll mod=a[pos].num,p=a[pos].p;
for(int i=1;i<mod;i++)
if(i%p)t.first=t.first*i%mod;
t.first=qpow(t.first,n/mod,mod);
for(int i=1;i<=n%mod;i++)if(i%p)t.first=t.first*i%mod;
t.second+=n/p;
pii w=deal(n/p,pos);
t.first=t.first*w.first%mod;
t.second+=w.second;
return t;
}
ll A[maxn],M[maxn],Ans=1;
void exgcd(ll a,ll b,ll& d,ll& x,ll& y){
if(b==0){d=a;x=1,y=0;return;}
exgcd(b,a%b,d,x,y);
ll t=x;
x=y;
y=t-a/b*x;
}
void CRT(){
ll AA=A[1],MM=M[1];
for(int i=2;i<=tot;i++){
ll d,k,x,y;
exgcd(MM,M[i],d,x,y);
x*=(A[i]-AA);
x=(x%M[i]+M[i])%M[i];
AA=(AA+x*MM)%(MM*M[i]);
MM=MM*M[i];
}
printf("%lld\n",AA);
}
int main(){
//freopen(FILE".in","r",stdin);
//freopen(FILE".out","w",stdout);
mod=read(),n=read(),m=read();
ll sum=0;
up(i,1,m)w[i]=read(),sum+=w[i];
if(sum<n)w[++m]=n-sum;
else if(sum>n){
printf("Impossible\n");
return 0;
}
divide(mod);
for(int i=1;i<=tot;i++){
ll mod=a[i].num,p=a[i].p;
pii t=deal(n,i);
for(int j=1;j<=m;j++){
pii k=deal(w[j],i);
t.second-=k.second;
ll x,y,d;
exgcd(k.first,mod,d,x,y);
x=(x%mod+mod)%mod;
t.first=t.first*x%mod;
}
A[i]=t.first*qpow(p,t.second,mod)%mod,M[i]=mod;
}
CRT();
return 0;
}

View Code