【数论】[总结]扩展欧几里德算法

NOIP2012 Day2 T1 的 mod 一题会用到扩展欧几里德，本人理解能力有限，花了两三天终于理解了，这里总结一下

 

（以下的所有符号以C/C++中的符号为准，特别是 / 在下面的情况下相当于Pascal的div）

 

欧几里德的辗转相除法

int gcd(int a,int b) { return (b)?gcd(b,a%b):a; }

 

现在再看这样一个等式ax+by=c，如果 c % gcd(a,b) =0 的话，那么方程就必定有整数解

证明很容易，ax % gcd(a,b) =0，而且 bx % gcd(a,b) =0，那么 (ax+by) % gcd(a,b) =0，又 c % gcd(a,b) =0 ，所以ax+by=c必定有整数解

 

那么 ax+by=1 也必定有整数解（因为 c%1=0 是显然成立的 ），所以如果我们能求出 ax+by=1 的一组整数解x0,y0，那么 ax+by=c 的整数解就可以求得 x=c*x0,y=c*y0

 

现在我们求得了一组整数解x,y，但是这个x和y是其中的一组，比如我们要找最小的x该怎么办呢？我们可以由原方程ax+by得到 ax+by=a(x+b)+b(y-a) ，那么我们求的x最小整数解就可以通过+b来求得了，但是时间一下子又慢了（如果求得的x是一个很小的负数），所以我们可以根据负数取余的特性（绝对值取余后再加上符号，例如 -13%3 = -(13%3) ）快速得到最小正整数解：(x%b+b)%b                              <由于只求x，所以就没有处理y>

 

现在的问题就转化成了如何求 ax+by=1 的一组整数解了

 

欧几里得的算法是 gcd(a,b)=gcd(b,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)，那么同理 bx+(a%b)y=gcd(b,a%b)  成

立。为了表示方便，我们直接把gcd(a,b)写成1来表示（这里gcd(a,b)的值就是1）

 

由数学上的知识很容易得到 a%b = a-a/b*b

 

上面方程bx+(a%b)y=1就可以写成            bx+(a-a/b*b)y=1                化简得到  ==>     bx+ay-(a/b*b)y=1

 

好了关键就在这里，后面的a和b和前面ax+by=1的a和b是相同的，但是x和y不同，那么x和y之间有什么关系呢？我们把它们整理一下

ax+by=1

ay+b(x-a/b*y)=1

也就是说在进行了一次gcd后，以前的x变成了y，而y变成了x-a/b*y          。。。。。。。。。这个结论先放着，等会再用

 

gcd进行到边界，也就是b=0的时候返回a，这个时候我们看看方程ax+by=1，就成了ax=1，我们再想想这个时候的a是多少？对了，a肯定等于1！因为这时候的a=gcd(a,b)=1 ，所以很自然的我们解得 x=1,y=0

 

现在就把刚才的结论用到这里开始进行倒推，最后求得一组x,y的整数解

 

 

这个时候就可以把扩展欧几里德的代码放上来了

int gcd(int a,int b,int &x,int &y)//一个函数，返回gcd(a,b)
{
	if(b==0)
	{
		x=1;y=0;//边界时候的赋值
		return a;
	}
	else{
		int ans=gcd(b,a%b);
		int t=x;x=y;//这一句和后一句就是刚才推到的求法
		y=t-a/b*y;
		return ans;
	}
}

 

当然，如果你的目的仅仅是为了求得一组x,y的话，那么就可以定义为void，而不用返回值了，当然定义成int也没错。。。