欧几里得计算最大公因数数学原理及实现

发表于 | 分类于

最大公因数定义

最大公因数又名最大公约数，是指两个整数（本文约定为正整数）共同的最大约数。在实际编程中按照维基百科的说法主要有穷举法、素因数分解法、短除法、辗转相除法（欧几里得算法）。本文主要尝试分析实现穷举法和辗转相除法。具体分析实现见下小节。

穷举法

穷举法的思路非常简单，即在所有可能的解中遍历找出符合要求的解。按照遍历思路可以从小往大遍历，也可以从两个数中较小的一个数往前遍历。本文使用第二种思路实现穷举法暴力求解最大公约数，具体代码如下所示：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
#include <algorithm>
#include <<utility>>
using namespace std;

unsigned int BruteGCD(unsigned int num1, unsigned int num2)
{
        if (num1 > num2) {
                swap(num1, num2);
        }
        unsigned int result = 0;
        for (unsigned int i = num1; i >= 1; i--) {
                if (num1 % i == 0 && num2 % i == 0) {
                        result = i; break;
                }
        }
        return result;
}

欧几里得算法

欧几里得算法的主要依据是$gcd(a,b)=gcd(b, a\%b)$。该公式证明如下：

  1. 令$g$为$(a,b)$的任意一个公因数，则有$a=k_1g, b=k_2g$，假设$a \mod b=r$，则有$a = kb+r$
  2. 将前两式代入第三式得出$k_1g = kk_2g + r=>r = (k_1g-kk_2g)g$，可得出$r\mod g=0$
  3. 结合上述可知$a,b$任意的公因数也是$r$的因数，即可推知$(a,b)$的任意公因数是$(b,r)$的公因数
  4. 同理对于任意$(b,r)$的公因数$h$，$b=m_1h,r=m_2h=>a=km_1h+m_2h$，亦可推出$(b,r)$的任意公因数也是$(a,b)$的公因数
  5. 最终得出两者具有相同的公因数序列，进一步得知最大公因数也相等

根据上述数学结论如何得出最大公因数算法呢？需要注意的是，对于两个数，如果其中一个数能被另外一个数整除，那么被除数一定是两个数的最大公因数，因此我们只要重复$gcd(a,b)=gcd(b, a\%b)$过程直到出现两数可以整除为止，以上思路的实现代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
/*显而易见的递归算法*/
unsigned int EuclideanAlgorithm(unsigned int num1, unsigned int num2)
{
        if (num1 % num2 == 0) {
                return num2;
        }
        else {
                return EuclideanAlgorithm(num2, num1 % num2);
        }
}

/*非递归写法*/
unsigned int EuclideanAlgorithm(unsigned int num1, unsigned int num2)
{
        while (num1 % num2 != 0) {
                unsigned int tmp = num1 % num2;
                num1 = num2;
                num2 = tmp;
        }
        return num2;
}

以上算法的时间复杂度如何计算呢？注意根据取模运算特点，对于$M>N$，那么有$M\mod N<M/2$，因此可以大致认为欧几里得算法的时间复杂度为$O(\log M)$。

算法测试

上述算法的测试代码如下，随机生成两个正整数，然后调用每个函数进行测试，为了测试代码正确性，我们直接抄录书上的代码作为测试基准，教材实现代码如下Gcd函数所示

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
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
#include <vector>
#include <random>
#include <ctime>
#include <iostream>
#include <algorithm>
#include <chrono>
#include <functional>
using namespace std;
using namespace chrono;

unsigned int Gcd(unsigned int M, unsigned int N)
{
        unsigned int Rem;

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

int main()
{
	default_random_engine randEngine(unsigned(time(nullptr)));
	uniform_int_distribution<unsigned int> uIntDis(unsigned int(1), unsigned int(10e4));

        int testTimes = 15;

	while ((testTimes--) > 0) {
                unsigned int num1 = uIntDis(randEngine);
                unsigned int num2 = uIntDis(randEngine);


                auto startTime = system_clock::now();

                unsigned int testResult = 0;
                //testResult = BruteGCD(num1, num2); /*调用具体的算法*/

                auto endTime = system_clock::now();

                unsigned int correctResult = Gcd(num1, num2);

                cout << "num1 = " << num1 << "\tnum2 = " << num2;
                cout << "\tgcd = " << correctResult << "\t";
                cout << ((testResult == correctResult) ? "Correct;" : "Wrong;");

                auto duration = duration_cast<microseconds>(endTime - startTime);
                cout << "TimeCost:" << duration.count() << endl;
        }
	return 0;
}

参考文章

最大公因数-维基百科

辗转相除法-维基百科

最大公约数(Gcd)两种算法