[学习笔记]扩展LUCAS定理

可以先做这个题[SDOI2010]古代猪文

此算法和LUCAS定理没有半毛钱关系。

【模板】扩展卢卡斯

不保证P是质数。

$C_n^m=\frac{n!}{m!(n-m)!}$

麻烦的是分母。

如果互质就有逆元了。

所以可以考虑把分子分母不互质的数单独提出来处理。

然鹅P太一般,直接处理要考虑的东西太多。

我们不妨令$p=p_1^{q_1}*p_2^{q_2}*...*p_k^{q_k}$

对每一个$p_i^{q_i}$分别求解(不妨叫这个数为$pk$)(这样会容易很多)

即求ai满足:$\frac{n!}{m!(n-m)!} = a_i\space mod \space pk$

然后可以$CRT$合并

(CRT可以合并的原因是,我们可以求出满足这些同余方程的通解。发现这些解mod lcm都是同一个数x。$C_n^m$一定是这些个解之一,不管是哪一个,$mod\space lcm$即mod p都是x。我们也就求出了答案)

$\frac{n!}{m!(n-m)!} = a_i\space mod \space pk$

现在不互质的就是pi的倍数

首先我们可以把分子分母的所有$pi$质因子都提出来,然后上下次数相消。

对于$n!$中p的质因个数,就是不断除以p^i下取整。

剩下的都是和pk互质的了。存在逆元

以求$19!\space mod\space 3^2$为例

$19!=1*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19$

提完质数3之后,变成:

$=(1*2*3*4*5*6)*3^2*(1*2*4*5*7*8*10*11*13*14*16*17*19)$

前面直接递归下去处理。后面的每一项,都可以%pk处理。

(不论分母还是分子位置,如果是分母位置,因为存在逆元,

10*inv =1 mod 9

1*inv = 1 mod 9

这两个inv显然是同一个inv.所以把所有项%pk没有问题

然后变成:

$=(1*2*3*4*5*6)*3^2*(1*2*4*5*7*8)^2*1)$

后面那个1是多出来的19

其实pk长度的循环节有n/pk个。直接算。剩下(例如这里的19),个数少于pk,直接算。

(所以,扩展LUCAS的重要适用条件是,$p_i^{q_i}$不能太大(1e5左右))

递归算出来即可。

对于分母位置的两个阶乘,算出来结果之后,再处理inv

(这里可以先乘完之后再找inv,不用一边找inv一边乘。)

(注意inv处理要用exgcd,不保证质数,不能用费马)

(CRT可以不用保存结果,Mi=p/pk 可以一次到位)

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define int long long
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
char ch;x=;bool fl=false;
while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
for(x=numb;isdigit(ch=getchar());x=x*+numb);
(fl==true)&&(x=-x);
}
namespace Miracle{
ll p;
ll qm(ll x,ll y,ll pk){
x%=pk;
ll ret=;
while(y){
if(y&) ret=(ret*x)%pk;
x=(x*x)%pk;
y>>=;
}
return ret;
}
ll calc(ll n,ll pi,ll pk){//计算阶乘部分 (质因子已经提前处理这里不予考虑)
if(!n) return ;
ll res=;
for(reg i=;i<pk;++i)//每个循环节
if(i%pi) res=(res*i)%pk;
res=qm(res,n/pk,pk);
for(reg i=;i<=n%pk;++i)
if(i%pi) res=(res*i)%pk;
return res*calc(n/pi,pi,pk)%pk;
}
void exgcd(ll a,ll b,ll &x,ll &y){//exgcd
if(!b){
x=,y=;return;
}
exgcd(b,a%b,y,x);
y-=(a/b)*x;
}
ll inv(ll n,ll pk){//逆元
ll x,y;exgcd(n,pk,x,y);
x=(x%pk+pk)%pk;
return x;
}
ll C(ll n,ll m,ll pi,ll pk){//计算C(n,m)mod pi^k
ll up=calc(n,pi,pk),d1=calc(m,pi,pk),d2=calc(n-m,pi,pk);
ll k=;
for(reg i=n;i;i/=pi) k+=i/pi;//处理质因子个数
for(reg i=m;i;i/=pi) k-=i/pi;
for(reg i=n-m;i;i/=pi) k-=i/pi;
return up*inv(d1,pk)%pk*inv(d2,pk)%pk*qm(pi,k,pk)%pk;
} ll CRT(ll b,ll mod){//CRT每步算出来了之后直接合并
return (b*inv(p/mod,mod)%p*(p/mod))%p;
}
ll EXLUCAS(ll n,ll m){//质因数分解+开始处理C
ll ret=;
ll tmp=p;
for(reg i=;(ll)i*i<=tmp;++i){
if(tmp%i==){
ll pi=i,pk=;
while(tmp%i==) pk*=i,tmp/=i;
(ret+=CRT(C(n,m,pi,pk),pk))%=p;
}
}
if(tmp>) (ret+=CRT(C(n,m,tmp,tmp),tmp))%=p;
return ret;
}
int main(){
ll n,m;
scanf("%lld%lld%lld",&n,&m,&p);
printf("%lld",EXLUCAS(n,m));
return ;
} }
signed main(){
Miracle::main();
return ;
} /*
Author: *Miracle*
Date: 2018/12/1 8:44:42
*/

这个算法的核心思路是:

1.不互质的要提出来单独处理

2.直接处理P,不互质的太多了

3.分成质因子处理,CRT合并

4.对于阶乘上下提出不互质的部分(质因子),转化成互质存在逆元的情况

5.观察剩余部分,后面可以对pk取模。

发现一部分还是阶乘,递归处理。

另一部分发现有循环节,利用循环节加速处理。

剩下的边角考虑一下。

(5本质上就是对每个数提取pi质因子,剩下的再乘起来。不过用递归和循环节加速了一下)

还有一个无聊的题:

[国家集训队]礼物

简单的组合数学,非要考你扩展LUCAS。。。。

[学习笔记]扩展LUCAS定理的相关教程结束。

相关推荐:

如何在 NumPy 中将二维数组沿新轴扩展为三维数组

本文介绍多种高效、内存友好的方法,将二维NumPy数组沿新增的第一维(axis=0)重复扩展为三维数组,涵盖np.expand_dims+np.repeat、np.tile、广播机制等主流方案,并对比其性能与适用场景。 本文介绍多种高效、内存友好的方法,将二维numpy数组沿新增的第一维(axis=...

生成每个元素可正、可负、可不选的扩展幂集(Extended Powerset)

本文介绍一种高效生成“带符号幂集”的方法:对输入序列中每个元素,独立选择保留原值、取负值或完全排除,从而生成所有合法子集,避免后续过滤,时间复杂度最优为o(3ⁿ)。 本文介绍一种高效生成“带符号幂集”的方法:对输入序列中每个元素,独立选择保留原值、取负值或完全排除,从而生成所有合法子集,避免后续过滤...

生成每个元素可正、可负、可不选的扩展幂集(三元选择)

本文介绍一种高效生成“扩展幂集”的方法:对输入可迭代对象中的每个元素,允许其以正值、负值或完全不出现于子集中——共三种独立选择,避免生成冗余组合后再过滤。 本文介绍一种高效生成“扩展幂集”的方法:对输入可迭代对象中的每个元素,允许其以正值、负值或完全不出现于子集中——共三种独立选择,避免生成冗余组合...