bzoj 2142: 礼物(组合数取模终极版) 组合数学+中国剩余定理+exgcd

题意

一年一度的圣诞节快要来到了。每年的圣诞节小E都会收到许多礼物，当然他也会送出许多礼物。不同的人物在小E心目中的重要性不同，在小E心中分量越重的人，收到的礼物会越多。小E从商店中购买了n件礼物，打算送给m个人，其中送给第i个人礼物数量为wi。请你帮忙计算出送礼物的方案数（两个方案被认为是不同的，当且仅当存在某个人在这两种方案中收到的礼物不同）。由于方案数可能会很大，你只需要输出模P后的结果。

m<=5,n,p<=10^9，p不一定为素数。

分析

显然答案就等于若干个组合数搞起来，这里就不推了。

因为p不是素数，所以显然不能用lucas定理来搞，那么我们就考虑把p分解质因数，分开求组合数然后用天朝剩余定理合并即可。

那么对于一个组合数Cmnmodpa=n!m!(n−m)!

考虑如何求n!

复制一波别人的题解：

c(n,m)我们把他看做三个阶乘，分别为n!,m!,(n-m)!

但是如果直接用逆元的话并不可以，因为有可能在这些阶乘中出现取模的数，这样的话就直接变为0了。所以我们要把这些数提取出来。

引用一名大牛的方法

n!=t1∗piu1

m!=t2∗piu2

(n−m)!=t3∗piu3

所以t1,t2,t3我们怎么求呢？

假设19!,pi=3（popoqqq大爷举的栗子）

19!=1*2…*19

19!=1∗2∗4∗5∗7∗8∗10∗11∗13∗14∗16∗17∗19∗36∗(1∗2∗3∗4∗5∗6)

所以后面的123456又是个子问题可以递归求解，而这个前面的乱七八糟一堆东西我们只需要把它乘到一起就好了！并且我们提取出来了6个pi，现在继续递归下去的话，乱七八糟的东西再次乘到一起，pi继续提取，然后继续递归就搞定n！了！

现在我们只需要t1∗niyuan(t2)∗niyuan(t3)∗piu1−u2−u3

这样我们就可以搞定所有的c(n,m)最终乘到一起就好了！

代码

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;

typedef long long LL;

LL n,m,p,a[20];
int tot;

LL ksm(LL x,LL y,LL mo)
{
LL ans=1;
while (y)
{
if (y&1) ans=ans*x%mo;
x=x*x%mo;y>>=1;
}
return ans;
}

LL work(LL p,LL pa,LL n,LL op)
{
if (n<=1) return 1;
tot+=op*(n/p);
LL x=1;
for (int i=1;i<=pa;i++)
if (i%p>0) x=x*i%pa;
x=ksm(x,n/pa,pa);
for (int i=1;i<=n%pa;i++)
if (i%p>0) x=x*i%pa;
return x*work(p,pa,n/p,op)%pa;
}

LL gcd(LL x,LL y)
{
if (!y) return x;
else return gcd(y,x%y);
}

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

void merge(LL &a1,LL &m1,LL a2,LL m2)
{
LL a=m1,b=m2,c=a2-a1,x,y;
exgcd(a,b,x,y);
LL d=gcd(a,b);c/=d;
x*=c;y*=c;
a1=(x*m1%(m1*m2)+a1+m1*m2)%(m1*m2);m1*=m2;
}

LL get_ny(LL x,LL p,LL mo)
{
return ksm(x,mo/p*(p-1)-1,mo);
}

LL solve(LL n,LL m)
{
LL p1=sqrt(p),pp=p,a=-1,mo=-1;
for (int i=2;i<=p1;i++)
if (pp%i==0)
{
LL w=1;
while (pp%i==0) pp/=i,w*=i;
tot=0;
LL x=work(i,w,n,1),y=work(i,w,m,-1),z=work(i,w,n-m,-1);
x=x*get_ny(y,i,w)%w*get_ny(z,i,w)%w*ksm(i,tot,w)%w;
if (a==-1) a=x,mo=w;
else merge(a,mo,x,w);
}
if (pp>1)
{
LL w=pp;
tot=0;
LL x=work(pp,w,n,1),y=work(pp,w,m,-1),z=work(pp,w,n-m,-1);
x=x*get_ny(y,pp,w)%w*get_ny(z,pp,w)%w*ksm(pp,tot,w)%w;
if (a==-1) a=x,mo=w;
else merge(a,mo,x,w);
}
return a;
}

int main()
{
scanf("%lld",&p);
scanf("%lld%lld",&n,&m);
LL sum=0;
for (int i=1;i<=m;i++)
{
scanf("%lld",&a[i]);
sum+=a[i];
}
if (sum>n)
{
puts("Impossible");
return 0;
}
LL ans=1;
for (int i=1;i<=m;i++)
{
ans=ans*solve(n,a[i])%p;
n-=a[i];
}
printf("%lld",ans);
return 0;
}