中国剩余定理与扩展 Lucas定理与扩展 学习笔记
2017-01-16 09:29
375 查看
中国剩余定理
问题
求同余方程组⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪x≡c1(modm1)x≡c2(modm2)x≡c3(modm3)...x≡ck(modmk){x≡c1(modm1)x≡c2(modm2)x≡c3(modm3)...x≡ck(modmk)
其中满足(mi,mj)=1,1<=i≠j<=k(mi,mj)=1,1<=i≠j<=k
x的最小正(非负)整数解
结论
令M=m1∗m2∗m3∗...∗mkM=m1∗m2∗m3∗...∗mk则x=(∑i=1kci∗Mmi∗inv(Mmi,mi))%Mx=(∑i=1kci∗Mmi∗inv(Mmi,mi))%M
证明
a.在模M意义下x只有唯一解 (有多解那还了得)b.令nini满足ni%Mmi=0ni%Mmi=0且ni%mi=1ni%mi=1,则N=∑i=1kci∗niN=∑i=1kci∗ni为原问题的一个解
c.根据上面的式子容易得出ni=Mmix,ni=miy+1ni=Mmix,ni=miy+1,则Mmix=miy+1Mmix=miy+1,即Mmix−miy=1Mmix−miy=1
d.由于mimi两两互质,所以MmiMmi与mimi也互质,令a=Mmi,b=mia=Mmi,b=mi,则ax−by=1ax−by=1
e.可以发现我们已经将其化简成扩展欧几里得的基本形式ax+by=(a,b)ax+by=(a,b),其中也可以认为x为a在模b意义下的逆元
代码
codevs 3990#include<algorithm> #include<iostream> #include<cstring> #include<cstdio> #include<cmath> using namespace std; #define LL long long LL k,l,r,n,M,x,y,Min,ans,m[15],c[15]; void exgcd(LL a,LL b,LL &x,LL &y) { if (!b) x=1,y=0; else exgcd(b,a%b,y,x),y-=a/b*x; } int main() { scanf("%d%lld%lld",&k,&l,&r); M=1LL; for (int i=1;i<=k;++i) { scanf("%lld%lld",&m[i],&c[i]); M*=m[i]; } for (int i=1;i<=k;++i) { LL a=M/m[i],b=m[i]; exgcd(a,b,x,y); x=(x%b+b)%b; if (!x) x+=b; n+=c[i]*a*x; } n%=M; if (!n) n+=M; if (r>=n) ans=(r-n)/M+1; if (l>=n) ans=ans-((l-n)/M+1); if ((l-n)%M==0) ++ans; if (ans) { if (l<=n) Min=n; else Min=n+((l-n)/M+1)*M; } printf("%lld\n%lld\n",ans,Min); }
扩展中国剩余定理
问题
求同余方程组⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪x≡c1(modm1)x≡c2(modm2)x≡c3(modm3)...x≡ck(modmk){x≡c1(modm1)x≡c2(modm2)x≡c3(modm3)...x≡ck(modmk)
x的最小正(非负)整数解
结论
对于两个方程x≡c1(modm1)x≡c1(modm1)
x≡c2(modm2)x≡c2(modm2)
将其合并成一个方程,有解条件为(m1,m2)|(c2−c1)(m1,m2)|(c2−c1)
m=m1m2(m1,m2)m=m1m2(m1,m2)
c=(inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2))%m2(m1,m2)∗m1+c1c=(inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2))%m2(m1,m2)∗m1+c1
最终得出一个式子x≡c(modm)x≡c(modm)
x=c%mx=c%m即为原问题的一个解
证明
a.将两个方程写成x=m1k1+c1,x=m2k2+c2x=m1k1+c1,x=m2k2+c2,两式联立得m1k1+c1=m2k2+c2m1k1+c1=m2k2+c2,移项m1k1=c2−c1+m2k2m1k1=c2−c1+m2k2,b.根据贝祖定理,以上等式有解充要条件为(m1,m2)|(c2−c1)(m1,m2)|(c2−c1)
c.将等式两边同除(m1,m2)(m1,m2)得m1(m1,m2)∗k1=(c2−c1)(m1,m2)+m2(m1,m2)∗k2m1(m1,m2)∗k1=(c2−c1)(m1,m2)+m2(m1,m2)∗k2,即m1(m1,m2)k1≡(c2−c1)(m1,m2)(modm2(m1,m2))m1(m1,m2)k1≡(c2−c1)(m1,m2)(modm2(m1,m2)),进一步化简k1≡inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)(modm2(m1,m2))k1≡inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)(modm2(m1,m2)),k1=inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)+ym2(m1,m2)k1=inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)+ym2(m1,m2)
d.将其回代x=m1k1+c1x=m1k1+c1,得x=inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)∗m1+ym1m2(m1,m2)+c1x=inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)∗m1+ym1m2(m1,m2)+c1,即x≡inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)∗m1+c1(modm1m2(m1,m2))x≡inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2)∗m1+c1(modm1m2(m1,m2))
e.至此我们又将其化简成了x≡c(modm)x≡c(modm)的形式,其中m=m1m2(m1,m2)m=m1m2(m1,m2),c=(inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2))%m2(m1,m2)∗m1+c1c=(inv(m1(m1,m2),m2(m1,m2))∗(c2−c1)(m1,m2))%m2(m1,m2)∗m1+c1
代码
pku 2891#include<algorithm> #include<iostream> #include<cstring> #include<cstdio> #include<cmath> using namespace std; #define LL long long #define N 1005 int k; LL c ,m ,c1,c2,m1,m2,t; bool flag; LL gcd(LL a,LL b) { if (!b) return a; else return gcd(b,a%b); } void exgcd(LL a,LL b,LL &x,LL &y) { if (!b) x=1LL,y=0LL; else exgcd(b,a%b,y,x),y-=a/b*x; } LL inv(LL a,LL b) { LL x=0LL,y=0LL; exgcd(a,b,x,y); x=(x%b+b)%b; if (!x) x+=b; return x; } int main() { while (~scanf("%d",&k)) { flag=true; for (int i=1;i<=k;++i) scanf("%I64d%I64d",&m[i],&c[i]); for (int i=2;i<=k;++i) { m1=m[i-1],m2=m[i],c1=c[i-1],c2=c[i]; t=gcd(m1,m2); if ((c2-c1)%t!=0) {flag=false;break;} m[i]=m1*m2/t; c[i]=inv(m1/t,m2/t)*((c2-c1)/t)%(m2/t)*m1+c1; c[i]=(c[i]%m[i]+m[i])%m[i]; } if (!flag) puts("-1"); else printf("%I64d\n",c[k]); } }
Lucas定理
问题
求Cmn%pCnm%p,其中p为质数结论
令n=nk∗pk+nk−1∗pk−1+...+n2∗p2+n1∗p+n0n=nk∗pk+nk−1∗pk−1+...+n2∗p2+n1∗p+n0
m=mk∗pk+mk−1∗pk−1+...+m2∗p2+m1∗p+m0m=mk∗pk+mk−1∗pk−1+...+m2∗p2+m1∗p+m0
则Cmn=∏i=0kCminiCnm=∏i=0kCnimi
证明
去问卢卡斯听说要用二项式定理什么的
我怎么可能会这种东西
记住就行了
代码
zoj 3557#include<algorithm> #include<iostream> #include<cstring> #include<cstdio> #include<cmath> using namespace std; #define LL long long LL n,m,Mod; LL fast_pow(LL a,LL p) { LL ans=1LL; for (;p;p>>=1,a=a*a%Mod) if (p&1) ans=ans*a%Mod; return ans; } LL inv(LL x) { return fast_pow(x,Mod-2); } LL C(LL n,LL m) { if (m>n) return 0LL; LL up=1LL,down=1LL; for (LL i=n-m+1;i<=n;++i) up=up*i%Mod; for (LL i=1;i<=m;++i) down=down*i%Mod; return up*inv(down)%Mod; } LL lucas(LL n,LL m) { if (m>n) return 0LL; LL ans=1; for (;m;n/=Mod,m/=Mod) ans=ans*C(n%Mod,m%Mod)%Mod; return ans; } int main() { while (~scanf("%lld%lld%lld",&n,&m,&Mod)) printf("%lld\n",lucas(n-m+1,m)); }
扩展Lucas定理
问题
求Cmn%pCnm%p结论
令p=p1k1∗p2k2∗...∗pqkqp=p1k1∗p2k2∗...∗pqkq
列出同余方程组
⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪ans≡c1(modp1k1)ans≡c2(modp2k2)...ans≡cq(modpqkq){ans≡c1(modp1k1)ans≡c2(modp2k2)...ans≡cq(modpqkq)
其中c1...cqc1...cq是对于每一个Cmn%pikiCnm%piki求出的答案
然后根据中国剩余定理合并
可见pikipiki并不是质数,而是一个质数的幂的形式
对于如何求Cmn%pikiCnm%piki将在证明中给出
证明
a.由于同于方程组在模p=p1k1∗p2k2∗...∗pqkqp=p1k1∗p2k2∗...∗pqkq意义下有唯一解,可以证明上面做法的正确性b.由于p=p1k1∗p2k2∗...∗pqkqp=p1k1∗p2k2∗...∗pqkq实质上是将p进行质因数分解,所以pikipiki满足两两互质,可以直接用中国剩余定理合并
c.对于如何求Cmn%pikiCnm%piki
根据Cmn=n!m!(n−m)!Cnm=n!m!(n−m)!,只要分别求出n!%piki,m!%piki,(n−m)!%pikin!%piki,m!%piki,(n−m)!%piki就可以通过逆元求出Cmn%pikiCnm%piki
d.对于如何求n!%pikin!%piki
我们以n=19,pi=3,ki=2n=19,pi=3,ki=2为例
n!=1∗2∗3∗4∗5∗6∗7∗8∗9∗10∗11∗12∗13∗14∗15∗16∗17∗18∗19n!=1∗2∗3∗4∗5∗6∗7∗8∗9∗10∗11∗12∗13∗14∗15∗16∗17∗18∗19
=(1∗2∗4∗5∗7∗8∗10∗11∗13∗14∗16∗17∗19)∗36∗(1∗2∗3∗4∗5∗6)=(1∗2∗4∗5∗7∗8∗10∗11∗13∗14∗16∗17∗19)∗36∗(1∗2∗3∗4∗5∗6)
根据这个例子发现,求解n!可以分为3部分:第一部分是pipi的幂的部分,也就是3636即pi⌊npi⌋pi⌊npi⌋,可以直接求解;第二部分是一个新的阶乘,也就是6!即⌊npi⌋!⌊npi⌋!,可以递归下去求解;第三部分是除前两部分之外剩下的数
考虑第三部分如何求解
发现第三部分在模pikipiki意义下是以pikipiki为周期的,即(1∗2∗4∗5∗7∗8)≡(10∗11∗13∗14∗16∗17)(modpiki)(1∗2∗4∗5∗7∗8)≡(10∗11∗13∗14∗16∗17)(modpiki),所以只求pikipiki长度的即可;但是还剩下一个孤立的19,可以发现剩下孤立的数长度不会超过pikipiki,只需要暴力求解即可
e.最后一个问题是对于求出的m!%pikim!%piki和(n−m)!%piki(n−m)!%piki有可能与pikipiki不互质,无法求逆元
所以要将m!%pikim!%piki和(n−m)!%piki(n−m)!%piki中质因子pipi先全部除去,求出逆元后再全部乘回去
计算n!中质因子p的个数x的公式为x=⌊np⌋+⌊np2⌋+⌊np3⌋+...x=⌊np⌋+⌊np2⌋+⌊np3⌋+...
递推式也可以写为f(n)=f(⌊np⌋)+⌊np⌋f(n)=f(⌊np⌋)+⌊np⌋
代码
codeforces2015ICL,Finals,Div.1#J#include<algorithm> #include<iostream> #include<cstring> #include<cstdio> #include<cmath> using namespace std; #define LL long long LL n,m,MOD,ans; LL fast_pow(LL a,LL p,LL Mod) { LL ans=1LL; for (;p;p>>=1,a=a*a%Mod) if (p&1) ans=ans*a%Mod; return ans; } void exgcd(LL a,LL b,LL &x,LL &y) { if (!b) x=1LL,y=0LL; else exgcd(b,a%b,y,x),y-=a/b*x; } LL inv(LL A,LL Mod) { if (!A) return 0LL; LL a=A,b=Mod,x=0LL,y=0LL; exgcd(a,b,x,y); x=((x%b)+b)%b; if (!x) x+=b; return x; } LL Mul(LL n,LL pi,LL pk) { if (!n) return 1LL; LL ans=1LL; if (n/pk) { for (LL i=2;i<=pk;++i) if (i%pi) ans=ans*i%pk; ans=fast_pow(ans,n/pk,pk); } for (LL i=2;i<=n%pk;++i) if (i%pi) ans=ans*i%pk; return ans*Mul(n/pi,pi,pk)%pk; } LL C(LL n,LL m,LL Mod,LL pi,LL pk) { if (m>n) return 0LL; LL a=Mul(n,pi,pk),b=Mul(m,pi,pk),c=Mul(n-m,pi,pk); LL k=0LL,ans; 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; ans=a*inv(b,pk)%pk*inv(c,pk)%pk*fast_pow(pi,k,pk)%pk; return ans*(Mod/pk)%Mod*inv(Mod/pk,pk)%Mod; } int main() { scanf("%I64d%I64d%I64d",&n,&m,&MOD); for (LL x=MOD,i=2;i<=MOD;++i) if (x%i==0) { LL pk=1LL; while (x%i==0) pk*=i,x/=i; ans=(ans+C(n,m,MOD,i,pk))%MOD; } printf("%I64d\n",ans); }
相关文章推荐
- 中国剩余定理及其扩展——学习(复习)笔记
- 学习笔记:中国剩余定理(CRT)及扩展中国剩余定理
- [BZOJ2142]礼物(扩展Lucas定理+中国剩余定理)
- [BZOJ 1951][Sdoi2010]古代猪文:Lucas定理|中国剩余定理|费马小定理|扩展欧几里得
- 【学习笔记-中国剩余定理】POJ1006 Biorhythms
- [BZOJ 1951] 古代猪文【Lucas定理/费马小定理/中国剩余定理/扩展欧几里得】
- codeforces 100633J. Ceizenpok’s formula (扩展Lucas定理+中国剩余定理)
- 「中国剩余定理CRT」学习笔记
- [BZOJ2142]礼物-扩展lucas定理-中国剩余定理
- 【bzoj2142】礼物 扩展Lucas定理+中国剩余定理
- [BZOJ2142] 礼物 - Lucas定理及扩展 - 中国剩余定理 - 扩展欧几里得算法
- BZOJ1951 [Sdoi2010]古代猪文 【费马小定理 + Lucas定理 + 中国剩余定理 + 逆元递推 + 扩展欧几里得】
- [BZOJ2142]-扩展Lucas+中国剩余定理
- [Codeforces2015ICL,Finals,Div. 1#J]Ceizenpok's formula(扩展Lucas定理+中国剩余定理)
- [中国剩余定理]【学习笔记】
- 【BZOJ2142】礼物(扩展lucas定理,中国剩余定理合并方程)
- BZOJ_P1951&Codevs_P1830 [SDOI2010]古代猪文(Lucas定理+扩展欧几里得+中国剩余定理)
- 学习笔记——中国剩余定理(CRT)
- 中国剩余定理学习笔记
- 【学习笔记】中国剩余定理