@(卢卡斯定理&扩展卢卡斯)
Lucas
卢卡斯定理:
\((^n_m)=(^{\frac{n}{p}}_{\frac{m}{p}}) * (^{n \pmod p}_{m \pmod {p}})\pmod p\)
后面部分可以递归,就是n和m在p进制下每一位的对应的组合数
证明要二项式定理
条件是P是质数
代码:
ll C(ll s,ll t){
if(s>t) return 0;
if(s==t) return 1;
return (f[t]*g[s])%p*g[t-s]%p;
}
ll lucas(ll s,ll t){
if(s==0) return 1;
pp=s%p;qq=t%p;
//printf("%lld %lld\n",s,t);
return C(pp,qq)*lucas(s/p,t/p)%p; //C处理出0~p-1的阶乘与其逆元就好
}
复杂度大概是\(O \log p\)?
EXlucas
处理p不为质数,考虑把p分解为
\[ p=\prod^{k}_{i-1}p_i^{e_i} \]
然后对于每个质数的乘积求出\((^n_m)\pmod {p_i^{e_i}}\),用CRT合并
每个质因子乘积的\((^n_m)\pmod {p_i^{e_i}}\)为
\[ \frac{n!}{m!*(n-m)!} \pmod {p_i^{e_i}} \]
考虑对于\(x! \pmod {p_i^{e_i}}\),有
\[ x! = {\prod^{p^k}_{i=1,i\pmod p !=0} i\pmod p}^{\frac{x}{p}}*\prod_{i=\frac{x}{p}+1}^x i\pmod p *{p^{\frac{n}{p}}} *\frac{x}{p} ! \]
举个例子,
当\(n=19,p=3,k=2\)时:
\(n!=1*2*3*\cdots*19\)
\(=(1*2*4*5*7*8*10*11*13*14*16*17*19)*3^{6}*6!\)
\(\equiv(1*2*4*5*7*8)^2*19*3^6*6!\)
然后发现x必须先要和p互质,
于是上面\(p^{\frac{n}{p}}\)就不算了
然后最后在外面补一个
复杂度PlogP的
看代码:
#include<bits/stdc++.h>
#define rep(i,a,b) for(register int i=(a);i<=(b);i++)
#define ll long long
ll test[15]={0,2,61,23,37};
inline ll read(){
ll x=0;int pos=1;char ch=getchar();
for(;!isdigit(ch);ch=getchar()) if(ch=='-') pos=0 ;
for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';
return pos?x:-x;
}
ll n,m,p;
void exgcd(ll a,ll b,ll &x,ll &y){
if(!b) x=1,y=0;
else{
exgcd(b,a%b,y,x);
y-=x*(a/b);
}
}
ll ksm(ll aa,ll bb,ll mod){
ll nans=1;
while(bb){
if(bb&1) nans=(nans*aa)%mod;
aa=(aa*aa)%mod;
bb>>=1;
}
return nans;
}
ll inv(ll xx,ll mod){
ll x,y;exgcd(xx,mod,x,y);
return (x%mod+mod)%mod;
}
ll fac(ll xx,ll pi,ll pk){
if(!xx) return 1;
ll res=1;
for(ll i=2;i<=pk;++i){
if(i%pi!=0) res=(res*i)%pk;
}
res=ksm(res,xx/pk,pk);
for(ll i=2;i<=xx%pk;++i){
if(i%pi!=0) res=(res*i)%pk;
}
return res*fac(xx/pi,pi,pk)%pk;
}
ll C(ll n,ll m,ll pi,ll pk){
ll up=fac(n,pi,pk),d1=fac(m,pi,pk),d2=fac(n-m,pi,pk);
ll resp=0;
for(ll i=n;i;i/=pi) resp+=i/pi;
for(ll i=m;i;i/=pi) resp-=i/pi;
for(ll i=n-m;i;i/=pi) resp-=i/pi;
return up*inv(d1,pk)%pk*inv(d2,pk)*ksm(pi,resp,pk)%pk;
}
ll CRT(ll b,ll mod){
return b*inv(p/mod,mod)%p*(p/mod)%p;
}
ll out(ll aa){
if(aa>9) out(aa/10);
putchar(aa%10+'0');
}
ll work(ll n,ll m){
ll res=0,tmp=p,pk=1;
int len=sqrt(p);
for(register int i=2;i<=len;++i){
if(tmp%i==0){
pk=1;while(tmp%i==0) pk*=i,tmp/=i;
int pp=i;
res=(res+CRT(C(n,m,pp,pk),pk))%p;
}
}
if(tmp>1){
int pp=tmp;pk=tmp;
res=(res+CRT(C(n,m,pp,pk),pk))%p;
}
return res;
}
int main(){
n=read();m=read();p=read();
out(work(n,m));
return 0;
}
例题
基本裸的,不难看出组合数的方案,重点在于取模
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<queue>
#include<cmath>
#define rep(i,a,b) for(long long i=(a);i<=(b);++i)
#define per(i,a,b) for(long long i=(a);i>=(b);--i)
#define ll long long
long long read(){
long long x=0,pos=1;char ch=getchar();
for(;!isdigit(ch);ch=getchar()) if(ch=='-') pos=0;
for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';
return pos?x:-x;
}
using namespace std;
ll n,m,mod,w[10];
ll ksm(ll a,ll b,ll mo){
ll na=1;
while(b){
if(b&1){
na=(na*a)%mo;
}
a=(a*a)%mo;
b>>=1;
}
return na;
}
ll fac(ll n,ll pi,ll pk){
if(!n) return 1;
ll res=1;
rep(i,2,pk) if(i%pi) res=(res*i)%pk;
res=ksm(res,n/pk,pk);
rep(i,2,n%pk) if(i%pi) res=(res*i)%pk;
return res*fac(n/pi,pi,pk)%pk;
}
void exgcd(ll a,ll b,ll &x,ll &y){
if(b==0) x=1,y=0;
else exgcd(b,a%b,y,x),y-=x*(a/b);
}
ll inv(ll now,ll mo){
ll x,y;
exgcd(now,mo,x,y);
return (x%mo+mo)%mo;
}
ll C(ll n,ll m,ll pi,ll pk){
ll up=fac(n,pi,pk),d1=fac(m,pi,pk),d2=fac(n-m,pi,pk);
ll resp=0;
for(int i=n;i;i/=pi) resp+=i/pi;
for(int i=m;i;i/=pi) resp-=i/pi;
for(int i=n-m;i;i/=pi) resp-=i/pi;
return up*inv(d1,pk)%pk*inv(d2,pk)%pk*ksm(pi,resp,pk)%pk;
}
ll CRT(ll a,ll b){
return a*inv(mod/b,b)%mod*(mod/b)%mod;
}
ll exlucas(ll n,ll m,ll p){
ll pk=1,res=0,tmp=p;
int len=sqrt(p);
rep(i,2,len){
if(tmp%i==0){
pk=1;
while(tmp%i==0) tmp/=i,pk*=i;
ll Ci=C(n,m,i,pk);
res=(res+CRT(Ci,pk))%p;
}
}
if(tmp>1){
ll Ci=C(n,m,tmp,tmp);
res=(res+CRT(Ci,tmp))%p;
}
return res;
}
int main(){
mod=read();
n=read();m=read();
rep(i,1,m){
w[i]=read();
}
ll ans=1,tot=0;
rep(i,1,m){
if(n<tot+w[i]){
printf("Impossible");
return 0;
}
ans=ans*(exlucas(n-tot,w[i],mod))%mod;
tot+=w[i];
}
printf("%lld",ans);
return 0;
}