拓展欧几里得算法详解

记得刚上初中的时候，数学课上老师讲过一种叫“二元一次不定方程”（形如ax + by = c （a≠0，b≠0，a、b、c为常数））的东西。当时，老师说：“二元一次不定方程有无数组解，对于任意实数x，都能找到一个y值使等式成立。”当时我们研究的范围是实数，如果说：当x，y以及a，b，c都是整数的时候，对于所有的整数x还都能找到一个整数y使得等式成立吗？

思考：ax + by = c
（a≠0，b≠0，a、b、c为常数）对于所有的整数x，是否都能找到一个整数y使得等式成立。（这个问题很简单，自己想一想吧！）

显然不能，就比如说x+2y=5当x=2时，y必须等于1.5等式才能成立，但是1.5不是整数。当数的范围被缩小到整数时，方程的解就可以被看成是直线“ax + by = c”穿过的“整数点”（横纵坐标都是整数的点），一般来说，解的个数还是无数个，就比如：2x+3y=1穿过(-1,1)，(2,-1)，(5,-3)，(8,-5)……但也有一些是无解的，就比如：2x+4y=1。

思考：二元一次不定方程在整数范围内有解的条件是什么？（一定要认真思考之后，再往下看。）

还是看刚刚我们说的那个式子：2x+4y=1。无论x、y是多少，2x+4y一定是偶数，而1却是奇数，一个奇数无论如何都不会与一个偶数相等。再比如：3x+6y=2，3x+6y一定是3的倍数，而2不是三的倍数，一个3的倍数不可能与一个不是3的倍数的数相等。因此：这个方程想要有解，c一定要是gcd(a, b)的倍数。

思考：如何用数学方法求出ax + by = c的一组整数解（原方程保证有解）。（这个问题有点难）

如果方程有解，那么c一定要是gcd(a, b)的倍数。不妨令gcd(a, b)=p，令c=np（n∈Z），原式为：ax + by = np即a (x/n) + b (y/n) = p。我们可以先去解这样一个方程：ax’+ by’=p，原方程的x就等于nx’，y就等于ny’。现在我们的当务之急就是解方程：ax + by = gcd(a, b)，这个问题比较不好想。
思考：如何解方程ax + by = gcd(a, b)。（提示：考虑辗转相除法的递归过程）

递归求解这个问题：辗转相除法告诉我们gcd(a, b)=gcd(b, a%b)，我们可以尝试从这个思路解决问题。假设我们已经知道了方程：bx+(a%b)y=gcd(b, a%b)=gcd(a, b)的解为(x0, y0)，如何去求方程ax + by =gcd(a, b)的解。尝试着把“模运算”展开，m%n=m-int(m/n)n （int表示取下整）。bx0 + (a – int (a/b)*b)*y0=gcd(a, b)，继续整理：y0*a +( x0-y0 int (a/b))b=gcd(a, b)。又因为ax + by = gcd(a, b)，所以：x=y0，y=x0-y0 int (a/b)。这样我们就可以递归下去求解了，退出条件是很显然的：当递归到gcd(a, b)*x+0*y=gcd(a, b)时，x=1 y=任意值，不妨令y=0（也可以令y等于其它值，这样会解出另外的一组的解，但也是成立的，从这个“y=任意值”，也可以分析出：该二元一次不定方程的整数解有无数组）。而这个求解ax + by = gcd(a, b)的过程就是所谓拓展欧几里得算法。

现在，试着写一下拓展欧几里得的代码吧！

struct PAIR//表示一个有序整数对(x,y) {    int x,y;    PAIR(int X=0,int Y=0)//初始化函数     {        x=X;y=Y;    }};PAIR EXGCD(int a,int b)//求解ax+by=gcd(a,b) 返回(x,y){    if(b==0)//边界条件         return PAIR(1,0);    PAIR LastAns=EXGCD(b,a%b);//递归求解     return PAIR(LastAns.y,LastAns.x-LastAns.y*int(a/b));}

这就是拓展欧几里得算法的代码的主要部分，但是还有一些值得注意的细节。就比如这个EXGCD函数会默认a≥b，否则程序将无法正确运行，所以我们可以再写一个这样的函数：

PAIR exgcd(int a,int b){    if(a>=b)//a≥b直接计算         return EXGCD(a,b);    else{        PAIR Ans=EXGCD(b,a);//交换一下         return PAIR(Ans.y,Ans.x);//这里也要反过来     }}

有了拓展欧几里得算法，我们就可以写一个求二元一次不定方程整数解的程序了，自己写写看吧！写完再看下面的样例代码。

int GCD(int a,int b){    if(b==0)        return a;    return GCD(b,a%b);}int gcd(int a,int b)//辗转相除法求最大公因数{    if(a<b)        swap(a,b);    return GCD(a,b);}PAIR SoveInt(int a,int b,int c){    int p=gcd(a,b);    if(c%p!=0)//无解情况    {        cout<<"Oops!No solution!"<<endl;        return PAIR(-1,-1);}    int n=c/p;    PAIR exAns=exgcd(a,b);    return PAIR(exAns.x*n,exAns.y*n);}

试着自己写一个求二元一次方程组整数解的程序吧！

赶稿匆忙，如有谬误，请谅解。

返回文章目录：Premier Bob的博客文章目录