用欧几里得算法计算最大公约数

给定两个非负整数 a 和 b，我们需要求它们的最大公约数（GCD），即同时整除这两个数的最大整数。通常记为 gcd(a, b)，其数学定义是：

gcd(a, b) = max { k > 0 : k ∣ a 且 k ∣ b }

符号 ∣ 表示整除，例如 k ∣ a 表示 k 整除 a。

如果一个数为零而另一个数非零，按照定义，它们的最大公约数就是那个非零的数。如果两个数都为零，最大公约数没有定义，因为公约数可以任意大；不过，为了保留 gcd 的结合律，将 gcd(0, 0) 定义为零很方便。于是得到一条简单规则：如果其中一个数为零，最大公约数就是另一个数。

下面的欧几里得算法能以 O(log min(a, b)) 的时间求出两个数的最大公约数。由于 gcd 满足结合律，求两个以上整数的最大公约数时，可以使用 gcd(a, b, c) = gcd(a, gcd(b, c))，依此类推。

这一算法最早见于欧几里得的《几何原本》（约公元前300年），但它可能起源于更早的时期。

算法

欧几里得算法最初采用这样的表述：不断用较大的数减去较小的数，直到其中一个数为零。如果 g 整除 a 和 b，它也整除 a − b。反过来，如果 g 整除 a − b 和 b，它也整除 a = b + (a − b)。因此，{a, b} 与 {b, a − b} 的公约数集合相同。

注意，在从 a 中减去 b 至少 ⌊a / b⌋ 次之前，a 仍是较大的数。为加快计算，可以将 a − b 替换成 a − ⌊a / b⌋b = a mod b。这样，算法就可以写成非常简单的递推关系：

gcd(a, b) = a，当 b = 0；否则 gcd(a, b) = gcd(b, a mod b)。

实现

int gcd (int a, int b) {
    if (b == 0)
        return a;
    else
        return gcd (b, a % b);
}

使用 C++ 的三元运算符，可以将函数体写成一行：

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

最后，下面是非递归实现：

int gcd (int a, int b) {
    while (b) {
        a %= b;
        swap(a, b);
    }
    return a;
}

从 C++17 开始，C++ 已经提供了求最大公约数的标准函数。

时间复杂度

算法的运行时间可以用拉梅定理估计。这个定理揭示了欧几里得算法与斐波那契数列之间令人意外的联系：

若 a > b ≥ 1，并且对于某个 n 有 b < F_n，则欧几里得算法最多进行 n − 2 次递归调用。

还可以证明，这一定理给出的上界是最优的。当 a = F_n、b = F_n−1 时，gcd(a, b) 恰好进行 n − 2 次递归调用。换句话说，相邻的斐波那契数是欧几里得算法的最坏情况输入。

斐波那契数呈指数增长，因此欧几里得算法的时间复杂度为 O(log min(a, b))。

另一种估计方法是注意到：在 a ≥ b 时，a mod b 至多为 a 的一半，因此算法每次迭代都会将较大的数至少减半。将这一推理用于一组数 a₁, …, a_n ≤ C 的最大公约数计算，还能将总运行时间估计为 O(n + log C)，而非 O(n log C)。这是因为每次非平凡迭代都会将当前最大公约数候选值至少减半。

最小公倍数

计算最小公倍数（LCM）可以通过以下简单公式归约为计算最大公约数：

lcm(a, b) = (a × b) / gcd(a, b)

因此，可以利用欧几里得算法，以相同的时间复杂度计算最小公倍数。一种实现如下：先将 a 除以最大公约数，以避免先计算 a × b 带来的中间结果溢出。

int lcm (int a, int b) {
    return a / gcd(a, b) * b;
}

二进制最大公约数算法

二进制最大公约数算法是对常规欧几里得算法的一种优化。

常规算法较慢的部分是取模运算。虽然我们将取模视为 O(1) 操作，但它比加法、减法或位运算这样的简单操作慢得多。因此，最好能避开取模。

确实可以设计一种不使用取模的快速 GCD 算法。它基于以下性质：

若两个数均为偶数，可以各提出一个因子2，再求剩余数的最大公约数：gcd(2a, 2b) = 2 gcd(a, b)。
若一个数为偶数、另一个数为奇数，可以从偶数中去除因子2：当 b 为奇数时，gcd(2a, b) = gcd(a, b)。
若两个数均为奇数，用其中一个数减去另一个数不会改变最大公约数：gcd(a, b) = gcd(b, a − b)。

仅利用这些性质，以及 GCC 提供的一些快速位操作函数，就可以实现一个快速版本：

int gcd(int a, int b) {
    if (!a || !b)
        return a | b;
    unsigned shift = __builtin_ctz(a | b);
    a >>= __builtin_ctz(a);
    do {
        b >>= __builtin_ctz(b);
        if (a > b)
            swap(a, b);
        b -= a;
    } while (b);
    return a << shift;
}

注意，这种优化通常没有必要。大多数编程语言的标准库已经提供 GCD 函数。例如，C++17 的 numeric 头文件中就有 std::gcd。

练习题

来源与许可

来源：CP-Algorithms：Euclidean algorithm for computing the greatest common divisor，由 CP-Algorithms 贡献者维护，原文标记为从 e-maxx 翻译而来。本中文版本翻译正文并调整排版，保留示例源码。原文与本译文均按 CC BY-SA 4.0 提供；材料按现状提供，不附保证。原文实现以非负整数为前提，先除后乘不能保证最终 LCM 不溢出；两个输入均为零时，所示 LCM 实现会出现除零，需要调用方另行处理。

文章版权归作者所有，未经允许请勿转载。

THE END