Bzoj2142礼物:组合数取模

题目链接：2142:礼物

组合数学推公式是很简单的，关键是取模

注意p并不是质数，而且p可能很大，所以lucas在这道题上并没有什么卵用

存在这样一个事实:设P分解质因数后有一项为pi^ai，那么C(x,y)%p%(pi^ai)=x%(pi^ai)，设为xi

而题目中已经给出pi^ai<=100000，所以我们只要计算出所有的xi，然后CRT一下就可以算出C(x,y)了

现在问题转化为怎么求xi

考虑把n!拆开得到1*2*3*4*5*6*...*n，显然和pi^ai不互质，在%(pi^ai)的意义下并没有逆元

所以我们把这n项全%(pi^ai)，得到1*2*3*...*pi^ai-1*1*2*3*...*pi^ai-1...*(pi^ai)^?*(1+2+3+....)

前面一段是循环的，中间有一段是循环节的一部分，后面(pi^ai)^(balabala)在组合数公式中是可以上下相减来消去的

所以前面一段快速幂，递归处理即可

#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<iostream>
#include<algorithm>
#define LL long long
using namespace std;
const int maxn=500010;
LL n,m,mod,cnt=0,num=0,sum=0;
struct Yinz{LL v,s;}yz[maxn];
LL w[maxn],jc[maxn];

void get_devide(){
LL tmp=mod;
for (int i=2;i<=sqrt(tmp);++i)
if (tmp%i==0){
yz[++cnt].v=i; yz[cnt].s=1;
while(tmp%i==0)yz[cnt].s*=i,tmp/=i;
}
if (tmp>1) yz[++cnt].v=tmp,yz[cnt].s=tmp;
}

LL powe(LL x,LL y,LL mod){
LL ret=1;
while (y){
if (y&1) ret=ret*x%mod;
x=x*x%mod; y>>=1;
}return ret;
}

LL get_val(LL N,LL x,LL y){
if (N<x) return jc
; num+=N/x;
return jc[N%y]*powe(jc[y-1],N/y,y)%y*get_val(N/x,x,y)%y;
}

void exgcd(LL a,LL b,LL &x,LL &y){
int tmp;
if(!b){x=1;y=0;return;}
exgcd(b,a%b,x,y);
tmp=x; x=y; y=tmp-(a/b)*y;
}

LL get_inv(LL a,LL b){
LL x,y; exgcd(a,b,x,y);
return (x%b+b)%b;
}

LL solve(LL x,LL y){
jc[0]=1; for(int i=1;i<y;++i) jc[i]=jc[i-1]*(i%x?i:1)%y;
num=0; LL tmp1=get_val(n,x,y),tmp2=1; LL tmp=num; num=0;
for (int i=1;i<y;++i) jc[i]=jc[i-1]*(i%x?get_inv(i,y):1)%y;
for (int i=1;i<=m;++i) tmp2=tmp2*get_val(w[i],x,y)%y;
return tmp1*tmp2%y*powe(x,tmp-num,y)%y;
}

int main(){
scanf("%lld%lld%lld",&mod,&n,&m);
LL sum=0;
for (int i=1;i<=m;++i){
scanf("%lld",&w[i]);
sum+=w[i];
if (sum>n) {printf("Impossible\n");return 0;}
}
if (sum<n) w[++m]=n-sum;
get_devide();//mod=pi^ai;
LL ans=0;
for (int i=1;i<=cnt;++i){
int now=mod/yz[i].s;
LL tmp1=solve(yz[i].v,yz[i].s);
LL tmp2=get_inv(now,yz[i].s);
ans=(ans+tmp1*tmp2%mod*now%mod)%mod;
}
printf("%lld",ans);
}