最大公约数算法

两个(或更多个)整数的最大公约数(Greatest Common Divisor, GCD)是指能同时整除它们的最大正整数,求解最大公约数最常用的算法是辗转相除法(又称欧几里得算法)及Stein算法,至于直接穷举的暴力算法,因为过于低效,实际中是不会使用的。下面来总结下这两种方法。

问题描述

给定两个正整数M和N,求其最大公约数。函数原型:

1
unsigned int GCD(unsigned int M, unsigned int N);

调用此函数时保证M≥N>0。

辗转相除法

辗转相除法的介绍可见维基百科,其理论基础就是下面这个公式:

GCD(M, N) = GCD(N, M%N)

此公式的证明在此不做介绍,直接给出算法实现:

1
2
3
4
5
6
7
8
9
10
unsigned int GCD(unsigned int M, unsigned int N) {
unsigned int Rem;

while (N > 0) {
Rem = M % N;
M = N;
N = Rem;
}
return M;
}

上述算法可以进一步优化,算法中使用了一个中间变量Rem,它用于交换两数,其实这是没有必要的。注意到:

GCD(M, N) = GCD(N, M % N) = GCD(M % N, N % (M % N))

在第一步计算出M % N后,保留N的值其实是没有必要的,可以直接跳过这一步,使用最后一个式子进行计算即可:

1
2
M = M % N;  // 这一步后,新的M的值即为M%N
N = N % M; // 相当于公式中的N%(M%N)

因为这相当于计算了2次,截止条件变为MN中的任意一个为0,否则就会出现%0这样的错误。返回结果是MN中不为0那个,不过由于另一个数为0,可以用M + N替代。这样就可以写出代码:

1
2
3
4
unsigned int GCD(unsigned int M, unsigned int N) {
while ((M %= N) && (N %= M));
return M + N;
}

如此简单而优雅!

Stein算法

此算法又被称为Binary GCD algorithm,它是针对欧几里得算法的一些缺陷提出的。

欧几里德算法是计算两个数最大公约数的传统算法,无论从理论还是从实际效率上都是很好的。但是却有一个致命的缺陷,这个缺陷在素数比较小的时候一般是感觉不到的,只有在大素数时才会显现出来。
一般实际应用中的整数很少会超过64位(当然现在已经允许128位了),对于这样的整数,计算两个数之间的模是很简单的。对于字长为32位的平台,计算两个不超过32位的整数的模,只需要一个指令周期,而计算64位以下的整数模,也不过几个周期而已。但是对于更大的素数,这样的计算过程就不得不由用户来设计,为了计算两个超过64位的整数的模,用户也许不得不采用类似于多位数除法手算过程中的试商法,这个过程不但复杂,而且消耗了很多CPU时间。对于现代密码算法,要求计算128位以上的素数的情况比比皆是,设计这样的程序迫切希望能够抛弃除法和取模。
由J. Stein 1961年提出的Stein算法很好的解决了欧几里德算法中的这个缺陷,Stein算法只有整数的移位和加减法。

算法的理论和实现方法在此从略,此处仅摘录一下实现代码以供参考,引用自WiKi

递归实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
unsigned int gcd(unsigned int u, unsigned int v)
{
// simple cases (termination)
if (u == v)
return u;

if (u == 0)
return v;

if (v == 0)
return u;

// look for factors of 2
if (~u & 1) // u is even
{
if (v & 1) // v is odd
return gcd(u >> 1, v);
else // both u and v are even
return gcd(u >> 1, v >> 1) << 1;
}

if (~v & 1) // u is odd, v is even
return gcd(u, v >> 1);

// reduce larger argument
if (u > v)
return gcd((u - v) >> 1, v);

return gcd((v - u) >> 1, u);
}

非递归实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
unsigned int gcd(unsigned int u, unsigned int v)
{
int shift;

/* GCD(0,v) == v; GCD(u,0) == u, GCD(0,0) == 0 */
if (u == 0) return v;
if (v == 0) return u;

/* Let shift := lg K, where K is the greatest power of 2
dividing both u and v. */
for (shift = 0; ((u | v) & 1) == 0; ++shift) {
u >>= 1;
v >>= 1;
}

while ((u & 1) == 0)
u >>= 1;

/* From here on, u is always odd. */
do {
/* remove all factors of 2 in v -- they are not common */
/* note: v is not zero, so while will terminate */
while ((v & 1) == 0) /* Loop X */
v >>= 1;

/* Now u and v are both odd. Swap if necessary so u <= v,
then set v = v - u (which is even). For bignums, the
swapping is just pointer movement, and the subtraction
can be done in-place. */
if (u > v) {
unsigned int t = v; v = u; u = t;} // Swap u and v.
v = v - u; // Here v >= u.
} while (v != 0);

/* restore common factors of 2 */
return u << shift;
}

扩展问题

如果需要求多个数的最大公约数,可以先求出任意两个数的最大公约数,然后再求此数与其它数间的最大公约数,即有以下等式成立:

GCD(a, b, c) = GCD(a, GCD(b, c)) = GCD(c, GCD(a,b))

由此,易写出求N个数最大公约数的算法:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#include <stdarg.h>

unsigned int GCDN(unsigned int parmN, ...) {
unsigned int M;
va_list ap;
va_start(ap, parmN);

M = va_arg(ap, unsigned int);

while (--parmN) {
M = GCD(M, va_arg(ap, unsigned int));
}

va_end(ap);
return M;
}

求出最大公约数的同时也就得到了最小公倍数(Lowest Common Multiple, LCM),二者存在以下关系:

GCD(M, N) * LCM(M, N) = M * N

故易得:

1
2
3
unsigned int LCM(unsigned int M, unsigned int N) {
return M * N / GCD(M, N);
}
文章目录
  1. 1. 问题描述
  2. 2. 辗转相除法
  3. 3. Stein算法
  4. 4. 扩展问题