扩展Lucas定理
原问题
求:
其中, (n) , (m) 很大,不能够直接求阶乘。(p) 大小可接受但是不保证是质数。
学过 (Lucas) 定理后,我们知道其限制是 (p) 必须为质数。
我们可以分解问题:
-
既然去了这个条件,我们可以将 (p) 进行质因数分解,模每一个 (p_i^{k_i}),然后就相当于得到一组同余方程组。
用 (CRT) 解出即可。
本人在这里因为 (CRT) 学的不太好产生了一个错误思路,导致纠结好久,总结在此 。
若有人同样误入歧途,希望有些许帮助。(如果没有纠结什么奇怪想法最好先不看
-
那么现在问题是解:
[C^{n}_{m}impliesfrac{m!}{(n)!(m-n)!} quad(mod;p_i^{k_i}) ]问题有二: 一,(m,n) 都很大,不能直接求出阶乘。
二,(n!,(m-n)!) 内可能有 (p_i) 这个因子,无法求出逆元。
很难解决,但是我们发现,阶乘取模一个数后呈现一个以长 (p_i^{k_i}) 循环节不断循环的模式,似乎可以只处理一个循环节,再处理余项,两者长度都不超过 (p_i^{k_i}) ,可以接受。
另外,从其中我们可以提出所有的 (p_i) 使式子变成这样:
[frac{frac{m!}{p_i^{a_1}}}{frac{n!}{p_i^{a_2}} imesfrac{(m-n)!}{p_i^{a_3}}} imes p_i^{a1-a2-a3} quad(mod;p_i^{k_i}) ]其中 (a_1-a_2-a_3) 一定是大于 (0) 的,因为 (a) 就是从 (1) 开始一定长的自然数中 (p) 这个因子的个数。数越大的区域,(p) 因子的密度越大(因为可能有 (p) 的几次方项),那么两个小数段的 (a) 值之和,一定比等长的连续数段的 (a) 值小。
这样就避开了逆元问题。
那么我们来具体对一个 (n!) 展开进行量化:
[n!equiv left. egin{align} [1 imes 2 imescdots (p_i-1) imes (p_i+1) imescdots imes (p_i^{k_i}-1)] \ imes[(1+p_i^{k_i}) imes(2+p_i^{k_i}) imescdots imes(p_i^{k_i}-1+p_i^{k_i})] end{align} ight}f 循环节 m \ imescdots imes[(1+lfloorfrac{n}{p_i^{k_i}} floor) imes p_i^{k_i}) imescdots imes n]quadf余项 m\ imes(1 imes2 imesdots imeslfloorfrac{n}{p_i} floor) imes p_i^{lfloorfrac{n}{p_i} floor}quadf 阶乘项 m\ (mod;p_i^{k_i}) ]出现了三部分,一部分是将有 (p) 因子的项都提出一个,形成一个新的阶乘。
第二部分是循环的,在模 (p_i^{k_i}) 意义下,所有加在后面的 (i imes p_i^{k_i}) 都可以不要。
式子就变成这样:
[n!equiv p_i^{lfloor frac{n}{p_i} floor}*lfloor frac{n}{p_i} floor!*(prod_{i=1;&;p_i mid i}^{p_i^{k_i}}i)^{lfloor frac{n}{p_i^{k_i}} floor}*(prod_{i=1;&;p_i mid i}^{n\%p_i^{k_i}}i)quad(mod;p_i^{k_i}) ]可以在函数中处理掉余项和循环节,阶乘项中, (p_i) 的幂次是要被提出去的不用管,而 (lfloor frac{n}{p_i} floor!) 项可以再递归下去,像 (Lucas) 那样,边界是输入为 (0) 。
这样就完成了对阶乘的处理。
从这里就能看出,(Exlucas) 精髓是一种 快速对模任意数意义下的阶乘进行处理 的方法。
这样处理出三个阶乘项,并且求出下面两个的逆元,乘起来即可。
-
还没完,我们还没处理后面的:
[p_i^{a1-a2-a3} ]可以见到,上面处理阶乘的递归函数,每一次递归可以提取出(p) 的 ({lfloorfrac{n}{p_i} floor}) 次方。
那么我们再建立一个递归函数 (G) :
[G(x)=lfloorfrac{n}{p_i} floor+G(lfloorfrac{n}{p_i} floor) ]递归边界是 (n<p_i) ,这时就没有 (p_i) 可以提取了,返回 (0) 。
之后快速幂求出,乘在 (2.) 的那一坨后面。
-
最后,把每个 (p_i) 的解用 (CRT) 组合起来。
-
完毕
(Exlucas) 的重点是一种 快速对模任意数意义下的阶乘进行处理 的方法,以及 用 (CRT) 进行质因数分解的组合 这一思路。
MOD
板子题,直接贴代码了。
(frak code)
#include<cstdio>
#include<algorithm>
using namespace std;
typedef int int_;
#define int long long
int n,m,p;
int gcd(int a,int b){
return b==0?a:gcd(b,a%b);
}
void exgcd(int a,int b,int &x,int &y){
if(b==0){
x=1,y=0;
}
else{
exgcd(b,a%b,y,x);
y-=(a/b)*x;
}
}
int get_inv(int x,int mod){
int d=gcd(x,mod);
if(d!=1) printf("err inv %lld %lld
",x,mod),exit(0);
int e,f;
exgcd(x,mod,e,f);
e=((e%mod)+mod)%mod;
return e;
}
int ksm(int x,int q,int mod){
int ret=1;
x%=mod;
while(q>0){
if(q&1) (ret *= x ) %= mod ;
(x *= x) %= mod ;
q>>=1;
}
return ret;
}
int fac(int x,int pi,int pik){
if(x == 0) return 1ll;
int mul=1;
for(int i=1;i < pik;i++){
if(i%pi != 0) (mul *= i) %= pik;
}
mul = ksm(mul,x/pik,pik);
for(int i=(x/pik)*pik+1;i <= x;i++){
if(i%pi != 0) (mul *= (i%pik)) %= pik;
}
(mul *= fac(x/pi,pi,pik)) %= pik;
return mul;
}
int g(int x,int pi){
if(x < pi) return 0;
else return x/pi+g(x/pi,pi);
}
int exlucas(int pi,int pik){
return ( ( ( fac(m,pi,pik) * get_inv(fac(n,pi,pik),pik) ) %pik * get_inv(fac(m-n,pi,pik),pik) ) %pik * ksm(pi,g(m,pi)-g(n,pi)-g(m-n,pi),pik) ) %pik;
}
int Crt(){
int ans=0;
int mod=p;
int t,ai,mi,e,f;
for(int i=2;i*i <= mod;i++){
if(mod%i != 0) continue;
t=1;
while(mod%i == 0){
mod/=i;
t*=i;
}
ai=exlucas(i,t);
mi=p/t;
exgcd(mi,t,e,f);
e=((e%t)+t)%t;
(e *= ai) %= p;
(ans += (e*mi)%p) %= p ;
}
if(mod != 1){
t=mod;
ai=exlucas(mod,t);
mi=p/t;
exgcd(mi,t,e,f);
e=((e%t)+t)%t;
(e *= ai) %= p;
(ans += (e*mi)%p) %= p ;
}
return ans;
}
int_ main()
{
//freopen("data.in","r",stdin);
//freopen("my.out","w",stdout);
scanf("%lld %lld %lld",&m,&n,&p);
printf("%lld",Crt());
return 0;
}
内容采用“知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议”进行许可。请您在转载时注明来源及链接。
(frak by;thorn\_)