【专题】数论 数论知识总结——史诗大作(这是一个flag) 【SRM20】数学场 乘法逆元的应用 欧拉函数
推荐:数论知识总结——史诗大作(这是一个flag)
---下面都是学习的笔记,还没有整理,比较凌乱,有需自取吧。---
【素数测试】Miller-Rabin算法
引用自:数论部分第一节:素数与素性测试by Matrix67
当p为素数时,有
费马小定理:a^(p-1)=1(%p)
Miller-Rabin测试:对于x^2=1(%p),有x=n-1或x=1(x^2-1=(x+1)(x-1))。
探测过程:对于一个底数(可随机)x=d*2^k,从d到x探测时候均符合d*2^k==1&&(d*2^(k-1)==1||d*2^(k-1)==n-1),则为素数探测成功。
底数2,7,61即可探测出所有int范围内的素数。
const int px[]={2,7,61}; int pow(int x,int k,int m){ int ans=1; while(k){ if(k&1)ans=1ll*ans*x%m; x=1ll*x*x%m; k>>=1; } return ans; } bool isprime(int x) { for(int i=0;i<=2;i++){ if(px[i]==x)return 1; int tmp=x-1; while(!(tmp&1))tmp>>=1; int last=pow(px[i],tmp,x),now; while(tmp<x-1){ now=1ll*last*last%x; if(now==1&&last!=1&&last!=x-1)return 0; last=now; tmp<<=1; } if(now!=1)return 0; } return 1; }View Code
【欧几里德算法】辗转相除法
唯一分解定理:任何一个大于1的整数都可以表示为若干素数的乘积。
欧几里德算法的简单解释:假设g=gcd(a,b),那么a和b都可以表示为若干个g的和,而我们的目的是求出g。
不妨a≥b(即使小于一次操作就可以自然反过来了),那么a%b可以求出g组成的最小的数字c,用b和c再做b%c直到c为0则此时得出g的最小表示b。
实际上是每次把除数和余数作为新的被除数和除数,核心是gcd(a,b)=gcd(b,a%b)。
特别的,gcd(a,0)=a。
gcd不怕0,不怕负数。
---
更相减损术:gcd中,一个数可以视为另一个数字和它们的差值,即有gcd(a,b)=gcd(a,a-b)。
感性地理解,就是一个数的对gcd的贡献可以转化为差值的贡献。
根据更相减损术的推广,若干数字的gcd可以转化为这些数的相减树+任意一个原数字,如gcd(a,b,c,d,e)=gcd(a-b,b-c,c-d,d-e,e)。
其实也很好理解,只要有一个原数就可以通过相减树得到其它所有数字。
辗转相除亦有类似推广,虽然不能还原出原数,但是信息全部继承了。
eg.已知a,a%b,b%c,则gcd(a,a%b)=gcd(a,b),则b参与贡献完毕,gcd(b,b%c)=gcd(b,c),因为b已参与贡献,只需要b%c参与就完成了c的贡献。
---
gcd可以预处理。
ll gcd(ll a,ll b) {return b==?a:gcd(b,a%b);} ll lcm(ll a,ll b) {return 1ll*a/gcd(a,b)*b;}View Code
gcd问题求解技巧:
1.设b[x]表示数列中有多少数字是x的倍数,求解b[x]可以O(n ln n)枚举每个数字的倍数或者O(n log n)线性筛后对每个数分解素因数。
序列中gcd为x的子集个数是2^b[x]-1。
例题:【SRM20】数学场
2.gcd(a,b,c,d,e)=gcd(a-b,b-c,c-d,d-e,e)
3.d|gcd(a,b)<=>d|a&&d|b。
【扩展欧几里德算法】(不定方程,模线性方程)
引用自:jumping_frog(扩欧)
刘汝佳《算法竞赛入门经典 第二版》
Matrix67: The Aha Moments(同余)
【基本算法】对于不完全为 0 的非负整数 a,b,gcd(a,b)表示 a,b 的最大公约数,必然存在整数对 x,y ,使得gcd(a,b)=ax+by。
证明:设 a>b。
1,显然当 b=0,gcd(a,b)=a。此时 x=1,y=0;
2,ab!=0 时
设 ax1+by1=gcd(a,b);
bx2+(a mod b)y2=gcd(b,a mod b);
根据朴素的欧几里德原理有 gcd(a,b)=gcd(b,a mod b);
则:ax1+by1=bx2+(a mod b)y2;
即:ax1+by1=bx2+(a-(a/b)*b)y2=ay2+bx2-(a/b)*by2;
根据恒等定理得:x1=y2; y1=x2-(a/b)*y2;
这样我们就得到了求解 x1,y1 的方法:x1,y1 的值基于 x2,y2.
上面的思想是以递归定义的,因为 gcd 不断的递归求解一定会有个时候 b=0,所以递归可以结束。
void gcd(ll a,ll b,ll& d,ll& x,ll& y) { if(!b){d=a;x=;y=;} else{gcd(b,a%b,d,y,x);y-=x*(a/b);} }View Code
上面求出ax+by=gcd(a,b)的一组解(x1,y1),任取另外一组解(x2,y2),令g=gcd(a,b),则有:
ax1+by1=ax2+by2
a(x1-x2)=b(y2-y1)两边同除以g
a'(x1-x2)=b'(y2-y1)此时a'与b'互素,因此(x1-x2)一定是b‘的整数倍
x1-x2=kb'代入得y2-y1=ka'
注意上述推导并没有用到ax+by的右边是什么,也就是说只要求出一组解就可以写成任意解,一切都只需要求出一组解。
结论:设a,b,c为任意整数,若方程ax+by=c的一组整数解为(x0,y0),则它的任意整数解都可以写成(x0+kb',y0-ka'),其中a'=a/gcd(a,b),b'=b/gcd(a,b),k取任意整数。
所以一定要先求出特解(利用右边),然后纯粹利用左边求出通解,在下面不定方程的求解中这个顺序尤为重要。
注意:再次强调,扩欧只能用于求解ax+by=gcd(a,b),特殊地,用于求解a'x+b'y=1(a'与b'互质)。
【应用一:求解不定方程】
(丢番图方程)ax+by=c求解步骤:
1.求g=gcd(a,b)。
2.若c%g≠0,方程无解。(扩欧只能求解ax+by=gcd(a,b)的整数解,如果整体乘上c/gcd(a,b)后不是整数了就无解,开始先同除g方便操作)
3.方程两边同时除以g,得到a'x+b'y=c'。(此方程完全等价,因为解没有变化)
4.通过扩展欧几里德算法求解a'x+b'y=1的一组解(x0,y0)
5.得到原不定方程的一组解(X=x0*c',Y=y0*c')
6.得到原不定方程通解(X+kb',Y-ka')注意一正一负
7.得到最小非负整数解:x=((X%b‘)+b’)%b‘ , y=((Y%a’)+a‘)%a’(+a或b转为正数)
PS:x、y可能为0,若要求正解要再次特判+a或+b。
PPS:除以g是为了方便计算,对解没有影响。除以g后扩欧只能计算a'x+b'y=1的方程,这时候解就要乘c/g得到原解。
PPPS:一个问题
{
求解不定方程时,我们通过ax+by=gcd(a,b)整体扩大c/gcd(a,b)倍来得到方程ax+by=c的解。
对于ax+by=gcd(a,b),通解是(x0+ka/g,y0-kb/g),所以我之前以为ax+by=c的通解应该由ax+by=gcd(a,b)的每一个解扩大c/gcd(a,b)得到。
换句话说,我以为ax+by=c的通解是(c/g*(x0+ka/g),c/g*(y0-kb/g)),后来发现我错了。
紫书中有句不起眼的话:“注意,上面的推导过程并没有用到“ax+by”的右边是什么”。
所以,ax+by=c的通解同样是(x0+ka/g,y0-kb/g),而(c/g*(x0+ka/g),c/g*(y0-kb/g))会让通解的数量减少。
那为什么之前会有错误的认识呢?其实(c/g*(x0+ka/g),c/g*(y0-kb/g))没有包含的那些通解就是除以c/g后不是整数的解。
换句话说,ax+by=gcd(a,b)中的有一些非整数解,扩大c/gcd(a,b)倍后就变成了ax+by=c的整数解。
综上,ax+by=c的通解求法是先求特解(X=x0*c/gcd(a,b),Y=y0*c/gcd(a,b)),后求通解(X+kb'/g,Y-ka'/g)。
}
//???紫书中似乎提到一个扩欧解x0,y0的特性:“x可能为负数,但因为|x|+|y|最小,所以加上n以后一定非负”???
【应用二:模线性方程组】
---
题外话:
减法取模:(a-b)%n=((a%n)-(b%n)+n)%n注意加n
大整数取模:自左向右逐个加入,每次乘10后取模。
幂取模:快速幂运算
---
目标:给定a,b,n,解方程ax≡b(%n)。
a≡b(%n)的含义是“a和b除以n的余数相同“,其充要条件是”a-b是n的整数倍“。
eg. a mod 3 = 1等价于a≡1(mod 3)等价于a-1=k*3。
---
同余性质证明模板:
性质:如果a≡b(mod m),x≡y(mod m),则a+x≡b+y(mod m)。
证明:条件告诉我们,可以找到p和q使得a-mp = b-mq,也存在r和s使得x-mr = y-ms。于是a-mp + x-mr = b-mq + y-ms,即a+x-m(p+r) = b+y-m(q+s),这就告诉我们a+x和b+y除以m的余数相同。
---
同余性质:
1.支持同余数的两边同加减乘幂,对于同余运算来说重要的只有余数部分,无关整的部分。
2.除法:若ac≡bc(%m),c≠0则a≡b(%m/gcd(c,m))其中gcd(c,m)表示c,m的最大公约数(简证:g=gcd(c,m),ac'g-bc'g=k*m'g,约去g,则可得到特殊情况。)
所以,同余符号左右和模数可以同时除以公约数。
特殊地 ,gcd(c,m)=1则a≡b(%m)
{证明:g=gcd(c,m)
ac'g-bc'g=k*m'g,约去g,得到ac'-bc'=k*m'那么得到特殊情况ac≡bc(%m),gcd(c,m)=1。
ac-bc=k*m即(a-b)c=k*m,由于c与m互质,则(a-b)=k*m即a≡b(%m),得证。
}
3.若a ≡ b (mod m),n|m,则 a ≡ b (mod n)
若a ≡ b (mod mi) (i=1,2...n) 则 a ≡ b (mod [m1,m2,...mn]) 其中[m1,m2,...mn]表示m1,m2,...mn的最小公倍数
---
ax≡b(%n)
ax-b=ny即ax-ny=b。令g=gcd(a,n),方程有解当且仅当g|b。
同除以g,得a'x-n'y=b'。
利用扩欧解a'x-n'y=1,本质上是在解ax≡g(%n),此时的解距离x0还差b/g倍。
x0=(x*b/g)%n,xi=x0+k*(n/g)。
1.证明方程有一解:当g|b时,x为ax≡g(%n)的解,则x0=(x*b/g)%n为ax0≡b(%n)的一解。
2.证明方程有g个解:一解相当于一个同余等价类,x0等价于所有x%n=x0的解x,这里有g个解是指有g个同余等价类。(我觉得紫书中也是这个意思,但非要用y来说明,令人浮想联翩……)
首先证明xi=x0+k*(n/g)也是方程的解,前已证得x0为ax0≡b(%n)的一解,则:
axi=a(x0+k*n/g)(%n)=ax0+ak*n/g(%n)由于g|a即a'=a/g,所以axi=ax0(%n),得证。
然后显然xi至多有g个,k的取值为0~g-1。当k≥g时,k=(g*[k/g]+k%g),[k/g]*(g*n/g)部分可以约去,则实际上与0~g-1等价。
3.解的间隔显然是n/g。
特别注意:大多数题目要求的解是非负数或正数,而扩欧解出来不一定是!!!
【乘法逆元(a-1)】
逆元详解by ACdreamer
乘法逆元小结by Tuesday..
乘法逆元:方程ax≡1(%n)的解称为a关于模n的逆(记得模n)。当gcd(a,n)=1时,该方程有唯一解(同余等价类);否则,该方程无解。
可以理解为倒数,a*x模n下等于1。在模n意义下,除以a等价于乘a的逆元,即a^(-1)≡a^(φ(n)-1) (%n)。
求逆元的方法:O(log(n))
前提:gcd(a,n)=1即a,n互素。
1.若n为素数(n=p),则由费马小定理a^(p-1)≡1(%p)得x=a^(p-2)%p。
2.若n不为素数,则由欧拉定理a^φ(n)≡1(%n)的x=a^(φ(n)-1)%n。
3.前两种方法需要快速幂运算,为了方便也可以直接用扩欧解ax-ny=1得到的解x即为逆元。
5.对于x=a!/b!(%p),b!|a!,若gcd(b!,p)>1,可以使用中国剩余定理求解(不会爆数据),过程见下面“中国剩余定理”。
ll inv(ll a,ll n) { ll d,x,y; gcd(a,n,d,x,y); //return d==1?(x+n)%n:-1; return (x%n+n)%n; }View Code
<线性求1~n的乘法逆元及其应用>
引用自:乘法逆元的应用bylzw4896s
乘法逆元函数inv(x)是完全积性函数,即inv(xy)=inv(x)*inv(y)或(xy)^(-1)=x^(-1)*y^(-1) (%p)。
过程:
1.过程中所有数都mod p,实际上在模运算中只要不涉及除法都可以模。
2.首先求出1!,2!,...,n!的值并记录,然后求出n!的逆元。
由与逆元是积性函数,有:
(n!)^(-1)=(n-1)!^(-1)*n^(-1),移项得
( (n!)^(-1) )/( n^(-1) ) = (n-1)! ^ (-1)
根据除以一个数等价于乘以这个数的逆元
(n!)^(-1)*n=(n-1)^(-1)
从而逆向求出1!,2!,...,n!的逆元。
3.接下来利用阶乘的逆求出每个数的逆。
( (n-1)! ^ (-1) ) * n ^ (-1) = (n!) ^ (-1)移项得
n ^ (-1) = ((n!) ^ (-1) ) /( (n-1)! ^ (-1) )
根据除以一个数等价于乘以这个数的逆元
n^(-1)=(n!)^(-1)*(n-1)!
注意,两次变换都是根据(n!)^(-1)=(n-1)!^(-1)*n^(-1)【(n!)-1=((n-1)!)-1*(n)-1】进行的(不同的数字除过去),这基于乘法逆元是积性函数。
void gcd(int a,int b,int &x,int &y) { if(b==){x=;y=;} else{gcd(b,a%b,y,x);y-=x*(a/b);} } void pre_inv() { fac[]=fac[]=; for(int i=;i<p;i++)fac[i]=(fac[i-]*i)%p; int xx,yy; gcd(fac[p-],p,xx,yy); inv[p-]=((xx%p)+p)%p;//扩欧解不一定是最小非负解! for(int i=p-;i>=;i--)inv[i]=(inv[i+]*(i+))%p; }View Code
【欧拉函数(phi)】
引用自:欧拉函数by handsomecui
欧拉函数求法与应用bysentimental_dog
1.在数论中,对于正整数N,少于或等于N ([1,N]),且与N互质的正整数(包括1,1与所有正整数互质)的个数,记作φ(n)。
φ(x)=x(1-1/p(1))(1-1/p(2))(1-1/p(3))(1-1/p(4))…..(1-1/p(n))其中p(1),p(2)…p(n)为x的所有质因数;
x是正整数; φ(1)=1(唯一和1互质的数,且小于等于1)。注意:每种质因数只有一个,与质因数个数无关!
φ(1)=1φ(2)=1φ(3)=2 φ(4)=2 φ(5)=4 φ(6)=2 大于一时互质的数字不含自身!
公式的简单理解:对于每个素因子pi,x范围内pi的倍数有x/pi个,即不含pi的有x*(1-1/pi),连乘起来即可。
2.p^k型欧拉函数:
对于质数p,φ(p)=p-1。(含1,不含本身)
对于质数p的k次幂,φ(p^k)=p^k-p^(k-1)=(p-1)p^(k-1)。
3.mn型欧拉函数:欧拉函数是不完全积性函数。
若m,n互质,φ(mn)=φ(m)φ(n)。
推论:当n为奇数时,φ(2n)=φ(n)。
//4.除了N=2,φ(N)都是偶数。(素数除了2减1后都是偶数,合数可以分解为素数之积)
//5.∑d|nφ(d)=n
//6.当n>1时,1~n中与n互质的整数和是为n*φ(n)/2。
int get_phi(int n) { //phi(n)=n*(1-1/p1)*(1-1/p2)...*(1-1/pk) int m=(int)sqrt(n+0.5); int ans=n; for(int i=;i<=m;i++) if(n%i==) { ans=ans/i*(i-);//在允许的情况下,先除后乘避免爆范围 while(n%i==)n/=i; } if(n>)ans=ans/n*(n-); return ans; }View Code
【欧拉定理&费马小定理】
欧拉定理:a,n为正整数且互素,a^φ(n)≡1(%n)
费马小定理:p为素数,a^(p-1)≡1(%p)
★重要应用:简化幂的模运算 a^b=a^(b%φ(p)) mod p
例如: 计算 7^222的个位数,实际上是求7^222被10除的余数。
且7与10互质,φ(10)=1,由欧拉定理知7^4= 1mod 10
故7^222=(7^4)^55*(7^2)=>(1^55)*(7^2)=>49=>9 mod 10
【欧拉线性筛】
引用自:贾志鹏线性筛
百度——贾志鹏线性筛。
<素数筛>
本质:每个数字只被最小的素数筛去,把每个数字分解为[素数*倍数]即[prime[j]*i],其中prime[j]为最小素数。
枚举i=1~n,每次选择比i小的素数组合为合数筛选。
if(i%prime[j]==0)break;如果当前素数已经能整除i,则这个倍数i中有质因数prime[j],之后筛的prime[j+1]*i这个数字的最小素数是prime[j],所以后面一定会被prime[j]乘一个倍数筛掉,此时筛就没有必要。
#include<cstdio> #include<algorithm> #include<cstring> const int maxsize=; int n,prime[maxsize],tot; bool mark[maxsize]; int main() { scanf("%d",&n); memset(mark,,sizeof(mark)); tot=; for(int i=;i<=n;i++) { if(!mark[i])prime[++tot]=i; for(int j=;j<=tot;j++) { if(i*prime[j]>n)break; mark[i*prime[j]]=; if(i%prime[j]==)break; } } printf("tot=%d\n",tot); for(int i=;i<tot;i++)printf("%d ",prime[i]); printf("%d\n",prime[tot]); return ; }View Code
找一个数字的所有质因数:记录每个数字的最小素数,然后对目标不断除最小素数就可以得到所有质因数。
线性筛求每个数的每个质因数可以用上面的方法,也可以枚举素数的方数的倍数。下面的方法复杂度也差不多。
遍历每个数字的每个质因数:c[i]记录数字i剩余数字,对每个数字i找其所有倍数j的c[j]除以c[i](c[i]=1就不执行),这样会自然的使用所有素数方数筛选,只不过多一些不影响复杂度的枚举,大概可以跑1千万。
<积性函数>
积性函数:对于所有互质的正整数a和b有f(ab)=f(a)f(b)的数论函数。
完全积性函数:对于所有正整数a和b有f(ab)=f(a)f(b)的数论函数
f(1)=1
1.性质一:对于n=∏pi^ki,有f(n)=∏f(pi^ki)
2.性质二:对于完全积性函数,还有f(n)=∏f(pi)^ki以及f(n^k)=f(n)^k
3.常见积性函数:欧拉函数,莫比乌斯函数,乘法逆元函数等。
坑:浅谈一类积性函数的前缀和
<欧拉函数>
线性筛欧拉函数本质上是因为欧拉函数是积性函数(但不是完全积性函数)。
欧拉函数是积性函数的原因可以从公式上显然看出。
欧拉函数:φ(x)=x(1-1/p(1))(1-1/p(2))(1-1/p(3))(1-1/p(4))…..(1-1/p(n))其中p(1),p(2)…p(n)为x的所有质因数
筛选过程和素数筛十分相似,因为本质上都是根据唯一分解定理,一个数由其最小的素数筛掉,只是过程中顺便求phi。
1.若i为素数,phi[i]=i-1。
2.若i%prime[j]==0,则i中有prime[j],那么phi[i*prime[j]]中有的素数phi[i]都有,差的只是最前面乘的数字多个prime[j],即phi[i*prime[j]]=phi[i]*prime[j]。然后和素数筛一样退出。
3.若i%prime!=0,则i中没有prime[j],那么phi[i*prime[j]]中缺少prime[j],最前面乘的数字多个prime[j],那么就是乘上(prime[j]-1)/prime[j]和prime[j],约分得phi[i*prime[j]]=phi[i]*(prime[j]-1)。
void pre_phi(int x) { phi[]=;int tot=; memset(mark,,sizeof(mark)); for(int i=;i<=x;i++) { if(!mark[i]) { prime[++tot]=i; phi[i]=i-; } for(int j=;j<=tot;j++) { if(i*prime[j]>x)break; mark[i*prime[j]]=;//筛 if(i%prime[j]==)phi[i*prime[j]]=phi[i]*prime[j]; else phi[i*prime[j]]=phi[i]*(prime[j]-); } } }View Code
<莫比乌斯函数>
1.μ(1)=1。
2.x含有奇数个素因子,μ(x)=-1。
3.x含有重复素因子,μ(x)=0。
4.x含有偶数个素因子,μ(x)=1。
显然,μ函数是积性函数,那么根据μ(1)=1和μ(p)=-1(p为素数),可以筛出所有莫比乌斯函数。
另有埃式筛法,μ(x)是其所有因子的μ和取反。
例题:【CodeForces】585 E. Present for Vitalik the Philatelist
【中国剩余定理(CRT)】
引用自:孙子定理_百度百科
由于公式较多,就直接贴图片了。
详细过程看百科,只要顺序看下来就能理解了,下面的过程是从那里简化的。
中国剩余定理实际上实在求解一元同余方程组(S)
<中国剩余定理的通解>
当m1...mn两两互质时有解,通解可以采用构造法获得。
注意此时逆元必定只针对模mi意义下成立(mi之间两两互质)
方程组(S)的通解形式为:
★在模M意义下只有唯一解:
<证明>
由以上两条就很容易看出构造法的神奇了,在%mi意义下,ti*Mi=1(gcd(Mi,mi)=1),而其它ajtjMj由于模数互质无法变1,又因为Mj中含mi而约去,最终剩下ai。
由于mi两两互质,显然下一个解间隔为M,故得证。
★重要应用:处理非互素逆元
【BZOJ】2142 礼物CRT