简单介绍
扩展欧几里得算法是用于求解类似于这样的式子的,据说这个式子叫做贝祖等式。
a∗x+b∗y=gcd(a,b)
其中,a,b,x,y 都为整数。不过一般为了方便,都会将它化成这样
a∗x+b∗y=1
其中,a,b 互质(其实就是最上面那个式子中a,b都除以一个 gcd(a,b) 就行了)
欧几里得算法(辗转相除法)
- 要讲扩展欧几里得算法,那么就先讲讲欧几里得算法,之后再将它扩展出来。
- 所谓欧几里得算法,就是用来求两个整数的最大公因数的方法,相信大家都会,不过这里还是给出简单的证明吧。
- 首先先证明辗转相减法的正确性,即
gcd(a,b)=gcd(b,a−b)
- 不妨设 a=p∗m, b=p∗n, gcd(m,n)=1
- 那么就有 a−b=(m−n)∗p ,接下来只需证明 gcd(m−n,n)=1 就行了
- 利用反证法, 若 m−n=q∗k1, n=q∗k2, (q≠1),则 m=(m−n)+n=q∗(k1+k2),这样一来gcd(n,m) 就变成 q了, 矛盾,固 gcd(m,n)=1
- 所以,辗转相减法是正确的
- 之后,辗转相除法就很简单了,它的内容是这样的
gcd(a,b)=gcd(b,a%b)
其实就相当于做多次辗转相减法嘛~
扩展欧几里得算法
大致思路
扩欧的思路大概就是对a,b做一次辗转相除,得到一堆式子,再利用这些式子后往前推出二元一次不定方程的一个特解。
方法&证明
对于要求解的式子,为了方便,不妨令
G=gcd(a,b)
a∗x+b∗y=G
我们对 a,b 做一次辗转相除法,可以得到多个这样的式子:
{a∗x1+b∗y1b∗x0+(a%b)∗y0=G=G
由于是从下到上,即 x0 和 y0 是已知的,那么经过 剧烈的运算,可以解出 x1 和 x2
{x1y1==y0x0−(a/b)∗y0
注意,这里的 a 和 b 可能和原式子中的 a、b 不同,它们是辗转相除中的 a, b。
这样,当b==0时,a肯定是 G 了,此时返回x1=1,x2=0,其余情况就用上面解出的公式求,回到最开始的原式时,得到的x1、x2就是所求的特解了。
关于通解
对于一次不定方程
a∗x+b∗y=gcd(a,b)
只要我们求出一组特解 (x′,y′) ,那么 我们可以得出一个通解
{xy=x′+t∗b=y′−t∗a
其中, t 是一个整数。
证明
首先,对于所有 (x,y) 符合上面的要求,那么它们肯定是一组解
a∗(x′+t∗b)+b∗(y′−t∗a)=a∗x′+b∗y′+a∗t∗b−b∗t∗a=a∗x′+b∗y′然后再证明一下数对 (x,y) 为此不定方程的解当且仅当符合上面的要求。用反证法可以很容易地证明出来。
程序
#include <cstdio>
int a,b,x,y;
void gcd(int &x1,int &y1,int aa,int bb){ //为了不和主程序中的a和b混淆,就用aa和bb
if (bb==0){x1=1; y1=0; return;} //关键的一步,相当于辗转相除法求出最大公约数的那个式子,其中aa就是最大公约数了
int x0,y0;
gcd(x0,y0,bb,aa%bb); //求出x0和y0
x1=y0; y1=x0-aa/bb*y0; //根据上面的公式求出x1和y1
if (x1<0) x1+=bb,y1-=aa; //为了防止范围溢出,利用通解性质适当调整一下
if (x1>bb) x1-=bb,y1+=aa; //同上
}
int main(){
scanf("%d%d",&a,&b); //题目中的a和b
gcd(x,y,a,b);
printf("%d %d",x,y); //一组特解x1,y1就出来了
}
例子模拟
下面对于
19∗x+15∗y=1
来模拟一次程序实现的过程
(我是递归到bb==0时才开始按回溯的顺序记录的,你们可以自己做一次辗转相除,看看除下来对应的aa、bb是不是相同)
aa | bb | x1 | y1 | x0 | y0 |
---|---|---|---|---|---|
1 | 0 | 1 | 0 | ||
3 | 1 | 0 | 1 | 1 | 0 |
4 | 3 | 1 | -1 | 0 | 1 |
15 | 4 | -1 | 4 | 1 | -1 |
19 | 15 | 4 | -5 | -1 | 4 |
所以,最终得出的特解就是
{xy==−54
验算一下, 4∗19+(−5)∗15=1 ,确实没错!实际上,在上面每一步中,你算一算
x1∗aa+y1∗bb
的值,它们都是等于1(即gcd(a,b))的。
除此之外每一行的 x0 和 y0 其实就是上一行的 x1 和 y1 。
应用
扩欧应用广泛,下面举一些例子
求解同余方程
a∗x≡gcd(a,n) (mod n)
下面是更常见的一个形式,不过要求 a 和 n 互质
a∗x≡1 (mod n)可以将其看成 a∗x+n∗y=gcd(a,n) ,就和上面讲的东西一样了
求乘法逆元
http://blog.csdn.net/jackypigpig/article/details/56846964求解同余方程组
⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪xxx...x≡a1modb1≡a2modb2≡a3modb3≡anmodbn
值得注意
- 我们可以发现,当 a 和 b 是斐波那契数列里相邻的两个元素时,时间复杂度是相对更高的。
- 在实际程序实现时,可能会出现溢出的现象,所以可以在适当的地方“模”一下 b。
小小的说明
对于本篇博文中的一些为了书写方便而不规范的符号说明一下
- 有些地方取模运算(mod)用了
%
代替,因为 mod... 太长了。 - 文中的
/
符号实际上是“整除”的意思,例如在“方法&证明”中,a/b实际上是 ⌊ab⌋ ,由于它实在太高了,所以就用矮一点的/
代替了。