卢卡斯定理&扩展卢卡斯

@

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;
}

例题

luogu2183 [国家集训队] 礼物

基本裸的,不难看出组合数的方案,重点在于取模

#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;                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                 
}
posted @ 2019-07-27 23:12  lcyfrog  阅读(312)  评论(0编辑  收藏  举报