欧几里德与扩展欧几里德算法----数论
转载自https://www.cnblogs.com/hadilo/p/5914302.html
一、欧几里得算法(重点是证明,对后续知识有用)
欧几里得算法,也叫辗转相除,简称 gcd,用于计算两个整数的最大公约数
定义 gcd(a,b) 为整数 a 与 b 的最大公约数
引理:gcd(a,b)=gcd(b,a%b)
先把这个引理的结论解释一下:
- 假设a=mc,b=nc,如果a,b的最大公约数是c,则m,n一定互质;
- 假设r=a%b=a-pb=mc-pnc=(m-pn)c,由1可知,r,b的最大公约数是c,且m-pn,n互质,
所以gcd(a,b)==gcd(b,a%b)
证明:
设 r=a%b , c=gcd(a,b)
则 a=xc , b=yc , 其中x , y互质
r=a%b=a-pb=xc-pyc=(x-py)c
而b=yc
可知:y 与 x-py 互质
证明:
假设 y 与 x-py 不互质
设 y=nk , x-py=mk , 且 k>1 (因为互质)
将 y 带入可得
x-pnk=mk
x=(pn+m)k
则 a=xc=(pn+m)kc , b=yc=nkc
那么此时 a 与 b 的最大公约数为 kc 不为 k
与原命题矛盾,则 y 与 x-py 互质
因为 y 与 x-py 互质,所以 r 与 b 的最大公约数为 c
即 gcd(b,r)=c=gcd(a,b)
得证
当a%b=0时,gcd(a,b)=b
这样我们可以写成递归形式
int gcd(int a, int b)//最大公约数 { return b == 0 ? a : gcd(b, a % b); } int lcm(int a, int b)//最小公倍数 { return a / gcd(a, b) * b; }
二、扩展欧几里得算法
扩展欧几里得算法,简称 exgcd,一般用来求解不定方程,求解线性同余方程,求解模的逆元等
引理:存在 x , y 使得 gcd(a,b)=ax+by
证明:
当 b=0 时,gcd(a,b)=a,此时 x=1 , y=0
当 b!=0 时,
设 ax1+by1=gcd(a,b)=gcd(b,a%b)=bx2+(a%b)y2
又因 a%b=a-a/b*b (参见文章开头的两个假设,这里的a/b取的是整数(假设里面的p),a-a/b*b!=0)
则 ax1+by1=bx2+(a-a/b*b)y2
ax1+by1=bx2+ay2-a/b*by2
ax1+by1=ay2+bx2-b*a/b*y2
ax1+by1=ay2+b(x2-a/b*y2)
解得 x1=y2 , y1=x2-a/b*y2
因为当 b=0 时存在 x , y 为最后一组解,gcd(a,b)=a,此时 x=1 , y=0
而每一组的解可根据后一组得到
所以第一组的解 x , y 必然存在
得证
根据上面的证明,在实现的时候采用递归做法
先递归进入下一层,等到到达最后一层即 b=0 时就返回x=1 , y=0
再根据 x=y’ , y=x’-a/b*y’ ( x’ 与 y’ 为下一层的 x 与 y ) 得到当层的解
不断算出当层的解并返回,最终返回至第一层,得到原解
这里的exgcd求得的x是 ax+by=1的解 (即a,b互质,gcd(a,b)==c)
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 niYuan(LL a, LL b) //求a对b取模的逆元 { LL x, y; exgcd(a, b, x, y); return (x + b) % b; }
模板题:https://www.cnblogs.com/-citywall123/p/10693036.html
三、exgcd 解不定方程
---------(使用不将a与b转为互质的方法)-----(gcd(a,b)!=c)
对于 ax+by=c 的不定方程,设 r=gcd(a,b)
当 c%r!=0 时无整数解
当 c%r=0 时,将方程右边 *r/c 后转换为 ax+by=r 的形式
可以根据扩展欧几里得算法求得一组整数解 x0 , y0
而这只是转换后的方程的解,原方程的一组解应再 *c/r 转变回去
(如 2x+4y=4 转换为 2x+4y=2 后应再将解得的 x , y 乘上2)
则原方程解为 x1=x0*c/r , y1=y0*(c/r);
通解 x=x1+b/r*t , y=y1-a/r*t ,其中 t 为整数
证明:
将 x , y 带入方程得
ax+ab/r*t+by-ab/r*t=c
ax+by=c
此等式恒成立
得证
这里 b/r 与 a/r 为最小的系数,所以求得的解是最多最全面的
证明:
为了推出证明中的 ax+by=c ,且想达到更小的系数,只能将 b/r 与 a/r 同除以一个数 s
而 b/r 与 a/r 互质,且 s 为整数,则 s=1 ,不影响通解
那么 b/r 与 a/r 就为最小的系数
得证
ll x,y,r,s; 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 gcd(ll a,ll b) { return b==0?a:gcd(b,a%b); } void BD(ll a,ll b,ll c,ll r) { exgcd(a,b,x,y); x=x*c/r;//得到原方程的解x和y y=y*c/r; }
模板题:https://www.cnblogs.com/-citywall123/p/10698205.html
模板题:http://www.cnblogs.com/hadilo/p/5917173.html
四、exgcd 解线性同余方程
线性同余方程有解的充分必要条件是当且仅当 c 能够被 a 与 b的最大公约数整除 c%gcd(a,b)==0
关于 x 的模方程 ax%b=c 的解---------要保证x的系数a为正,如果a小于0,等号两边乘以-1 eg:https://www.cnblogs.com/-citywall123/p/10811824.html
方程转换为 ax+by=c 其中 y 一般为非正整数
(PS怎么转换:ax%b=c%b
假设一个整数y; 因为ax%b=(ax+by)%b 所以(ax+by)%b=c%b --> ax+by=c)
则问题变为用 exgcd 解不定方程
解得 x1=x0*c/r
通解为 x=x1+b/r*t
设 s=b/r (已证明 b/r 为通解的最小间隔,r=gcd(a,b))
则 x 的最小正整数解为 (x1%s+s)%s
同理可得y的最小正整数解为(y1%ss+ss)%ss;------ss=a/r;
证明:
若 x1>0,则 (x1%s+s)%s=x1%s%s+s%s=x1%s=x1-ts (t∈N)
若 x1<0,因在 C++ 里 a%b=-(-a%b)<0 (a<0 , b>0) 如 -10%4=-2
则 (x1%s+s)%s=(-(-x1%s)+s)%s=(-(ts-x1)+s)%s=ts-x1 (t∈N)
即为 x1 通过加或减上若干个 s 后得到的最小正整数解
得证
亦可伪证 x1<0 的情况:设 x1=-5 , s=2
则 (x1%s+s)%s=(-5%2+2)%2=(-1+2)%2=3%2=1
即为 x1 加上 3 个 s 后的到的最小正整数解
ll x,y,r,s; 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 gcd(ll a,ll b) { return b==0?a:gcd(b,a%b); } void T(ll a,ll b,ll c) { r=gcd(a,b); s=b/r; exgcd(a,b,x,y);//得到x0 x=x*c/r; //得到x1 x=(x%s+s)%s; //得到最小正整数解 }
模板题:https://www.cnblogs.com/-citywall123/p/10694160.html
模板题:https://www.cnblogs.com/-citywall123/p/10680585.html
模板题:http://www.cnblogs.com/hadilo/p/5951091.html