<a href="http://zh.wikipedia.org/wiki/File:Euklid.jpg" target="_blank"></a>
算法描述
两个数a,b的最大公约数记为GCD(a,b)。a,b的最大公约数是两个数的公共素因子的乘积。如462可以分解成2 × 3 × 7 × 11;1071可以分解成3 × 3 × 7 × 17。462和1071的最大公约数等于它们共有的素因数的乘积3 × 7 = 21。如果两数没有公共的素因数,那么它们的最大公约数是1,也即这两个数互素,即GCD(a,b)=1。另g=GCD(a,b),则a=mg, b=ng,其中m,n均为正整数。由上述分析可知,m,n互素。因为m,n没有公共素因子,GCD(m,n)=1。
辗转相除法是一种递归算法。设k表示步骤数(从0开始计数),算法的计算过程如下。每一步的输入是都是前两次计算的余数rk−1和rk−2。因为每一步计算出的余数都在不断减小,所以,rk−1小于rk−2。在第k步中,算法计算出满足以下等式的商qk和余数 rk:
<dl><dd>rk−2 = qk rk−1 + rk</dd></dl>
其中rk < rk−1。也就是rk−2要不断减去rk−1直到比rk−1小。
在第一步计算时(k = 0),设r−2和r−1分别等于a和b,第2步(此时k = 1)时计算r−1(即b)和r0(第一步计算产生的余数)相除产生的商和余数,以此类推。整个算法可以用如下等式表示:
<dl></dl>
<dd>a = q0 b + r0</dd>
<dd>b = q1 r0 + r1</dd>
<dd>r0 = q2 r1 + r2</dd>
<dd>r1 = q3 r2 + r3</dd>
<dd>…</dd>
如果输入参数a小于b,则第一步计算的结果是交换两个变量的值:因为a < b,所以a和b相除得到的商q0等于0,余数r0等于a。所以在运算的每一步中得出的余数一定小于上一步计算的余数(rk一定小于rk−1)。由于每一步的余数都在减小并且不为负数,必然存在第N步时rN等于0,使算法终止,rN−1 就是a和b的最大公约数。其中N不可能无穷大,因为在r0和0之间只有有限个自然数。
辗转相除法的正确性可以用两步来证明。首先,算法的最终结果rN−1同时整除a和b。因为它是一个公约数,所以必然小于或者等于最大公约数g。然后,任何a和b的公约数都能整除rN−1。所以g一定小于或等于rN−1。两个不等式只在rN−1 = g是同时成立。
算法实现
递归版本
<a></a>
1 function gcd(a, b)
2 if a<b
3 swap(a,b);
4 if b==0
5 then return a;
6 else
7 return gcd(b, a mod b);
8 end
循环版本
1 function gcd(a,b)
3 then swap(a,b);
4 while(b!=0)
5 {
6 c = a mod b;
7 a = b;
8 b = c;
9 }
10 return a;
11 end
扩展欧几里德算法
<dd>sk = sk−2 − qk−1sk−1</dd>
<dd>tk = tk−2 − qk−1tk−1</dd>
算法开始时:
<dd>s−2 = 1, t−2 = 0</dd>
<dd>s−1 = 0, t−1 = 1</dd>
加上这个递归式后,当算法终止于rN = 0,贝祖等式的整数s和t分别由sN和tN给出。
这个算法的正确性可以用数学归纳法来证明。假设递归至第k−1步是正确的,也就是假设:
<dl><dd>rj = sj a + tj b</dd></dl>
对所有j小于k成立。则第k步运算得出以下等式:
<dl><dd>rk = rk−2 − qk−1rk−1</dd></dl>
因为rk−2和rk−1被假定是正确的,所以可以用s和t表示:
<dl><dd>rk = (sk−2 a + tk−2 b) − qk−1(sk−1 a + tk−1 b)</dd></dl>
整理后得到第k步的结果,和我们期望得到的结果一致:
<dl><dd>rk = sk a + tk b = (sk−2 − qk−1sk−1) a + (tk−2 − qk−1tk−1) b</dd></dl>
算法分析
加百利·拉梅于1884年指出,一个算法的效率可以用计算所需步数乘以每步计算的开销表示。用辗转相除法计算两个数的最大公约数所需的步数不会超过其中较小数的位数h的5倍。因为每一步的计算开销通常也是h数量级的,所以辗转相除法的复杂度是h2。
计算两个自然数a和b的最大公约数所需的步数可以表示为T(a, b)。如果a和b的最大公约数是g,m和n是两个互素整数,使a = mg,b = ng,那么:
<dl><dd>T(a, b) = T(m, n)</dd></dl>
这可以通过在辗转相除法的计算过程中的每一步都除以g来证明。同样,当a和b同时乘以w时,计算步数不变:T(a, b) = T(wa, wb)。所以,对于数值上相近的数,如T(a, b)和T(a, b + 1),计算步数可能相差很大。
根据辗转相除法的递归性质可以得出另一个公式:
<dl><dd>T(a, b) = 1 + T(b, r0) = 2 + T(r0, r1) = … = N + T(rN−2, rN−1) = N + 1</dd></dl>
并且定义T(x, 0) = 0。
辗转相除法的平均步骤数有三种不同的定义。第一种定义是计算已知自然数a和从0到a − 1范围内随机选取的自然数b的最大公约数所需的时间T(a):
<dl><dd></dd></dl>
但是因为T(a, b)在连续整数间变化非常剧烈,所以T(a)的平均值也会显得很杂乱。
为了解决这个问题,第二种定义规定τ(a)只要取遍其中所有和a互素的数即可:
因为第一种定义可以通过用第二种定义的求和来完成:
所以也可以由以下公式近似:
第三种定义Y(n)定义为从1到n间随机选取a和b时计算他们的最大公约数所需的平均步骤数:
将T(a)的近似公式代入,得到Y(n)的近似值:
在辗转相除法的每一步中,商qk和余数rk都通过rk−2和rk−1求出:
<dl><dd>rk−2 = qk rk−1 + rk.</dd></dl>
所以每一步的计算开销主要与计算商qk的算法有关,因为余数rk可以很迅速地从rk−2、rk−1和qk计算出来:
<dl><dd>rk = rk−2 − qk rk−1.</dd></dl>
综合考虑算法需要的步数和每一步的计算开销,辗转相除法的随两个数字a和b的平均位数成平方级的速度增长(h2)。设h0、h1、…、hN−1表示计算过程中的余数r0、r1、…、rN−1的位数,因为算法的步数N随h线性增长,所以算法的运算时间为:
本文转自NewPanderKing51CTO博客,原文链接: http://www.cnblogs.com/newpanderking/archive/2011/07/25/2116323.html,如需转载请自行联系原作者