[SDOI2013]方程(exlucas+容斥)
题目
给一个方程和一些限制,要求一些数必须小于\(a_i\),一些数必须大于\(b_i\),求不同的正整数解的数量
思路
扩展卢卡斯定理模板题
对于第二类限制\(x_i\geq a_i\),将\(m\)减去\(a_i-1\),即\(x_i\geq 1\)
对于第一类限制\(x_i\leq a_i\),由于限制数很少可以容斥;不满足条件即\(x_i\geq a_i+1\),当成第二类限制做即可
用插板法求解为\(C_{m-1}^{n-1}\% p\),由于\(p\)不是质数,需要用扩展卢卡斯定理,码量瞬间多了几十排qwq
Code
#include<bits/stdc++.h>
#define Min(x,y) ((x)<(y)?(x):(y))
#define Max(x,y) ((x)>(y)?(x):(y))
using namespace std;
typedef long long ll;
int T,n1,n2;
ll n,m,p,a[20],f[1000005];
template <class T>
void read(T &x)
{
char c; int sign=1;
while((c=getchar())>'9'||c<'0') if(c=='-') sign=-1; x=c-48;
while((c=getchar())>='0'&&c<='9') x=(x<<1)+(x<<3)+c-48; x*=sign;
}
ll quickpow(ll a,ll b,ll mod)
{
ll ret=1;
while(b)
{
if(b&1) ret=ret*a%mod;
a=a*a%mod;
b>>=1;
}
return ret;
}
ll exgcd(ll a,ll b,ll &x,ll &y)
{
if(!b) { x=1; y=0; return a; }
ll d=exgcd(b,a%b,x,y);
ll t=x; x=y; y=t-a/b*y;
return d;
}
ll inv(ll a,ll p)
{
ll x,y;
exgcd(a,p,x,y);
return (x%p+p)%p;
}
ll fac(ll n,ll pi,ll pk)
{
if(!n) return 1LL;
if(n<pi) return f[n];
return quickpow(f[pk-1],n/pk,pk)*f[n%pk]%pk*fac(n/pi,pi,pk)%pk;
}
ll C(ll n,ll m,ll pi,ll pk)
{
if(n<m) return 0;
f[0]=1;
for(int i=1;i<=pk;i++)
if(i%pi!=0) f[i]=f[i-1]*i%pk;
else f[i]=f[i-1];
ll jn=fac(n,pi,pk),jm=fac(m,pi,pk),jnm=fac(n-m,pi,pk);
int k=0;
for(ll i=n;i;i/=pi) k+=i/pi;
for(ll i=m;i;i/=pi) k-=i/pi;
for(ll i=n-m;i;i/=pi) k-=i/pi;
return jn * inv(jm,pk)%pk * inv(jnm,pk)%pk * quickpow(pi,k,pk)%pk;
}
ll crt(ll a,ll pk)
{
ll x=p/pk;
return a*x%p*inv(x,pk)%p;//关于pk的逆元
}
ll solve(ll n,ll m,ll pi,ll pk)
{
ll ret=0;
for(int i=0,t=(1<<n1);i<t;++i)
{
int opt=1;
ll nown=n;
for(int j=0;j<n1;++j)
if(i>>j&1)
opt=-opt,nown-=a[j+1];
ret=(ret+opt*C(nown,m,pi,pk))%pk;
}
return ret;
}
ll exlucas(ll n,ll m,ll P)
{
if(n<m) return 0;
ll ret=0;
for(ll i=2;i*i<=P;++i)
{
if(P%i==0)
{
ll pk=1;
while(P%i==0)
{
pk*=i;
P/=i;
}
ret=(ret+crt(solve(n,m,i,pk),pk))%p;
}
}
if(P!=1) ret=(ret+crt(solve(n,m,P,P),P))%p;
return (ret%p+p)%p;
}
int main()
{
read(T);read(p);
while(T--)
{
read(n);read(n1);read(n2);read(m);
for(int i=1;i<=n1+n2;++i) read(a[i]);
for(int i=n1+1;i<=n1+n2;++i) m-=(a[i]-1);
printf("%lld\n",exlucas(m-1,n-1,p));
}
return 0;
}