扩展欧几里得算法详解：从原理到C++实现与实战应用

1. 项目概述：为什么我们需要扩展欧几里得算法？

如果你写过C++，大概率用过std::gcd来求最大公约数，这背后的欧几里得算法优雅而高效。但很多朋友在遇到“求乘法逆元”或者“求解线性丢番图方程”这类问题时，就有点懵了，其实这正是扩展欧几里得算法（Extended Euclidean Algorithm, EEA）的用武之地。简单说，EEA不仅能算出两个整数的最大公约数（GCD），还能顺藤摸瓜，找出一对整数系数x和y，满足a*x + b*y = gcd(a, b)这个关键等式。

这个等式的威力超乎想象。在密码学里，它是RSA算法中计算私钥的核心；在竞赛编程中，它是求解模线性方程、求逆元的基石；甚至在你需要处理一些整数关系的数学模型时，它都能提供一种系统性的解法。网上很多EEA的代码实现，要么过于学术化让人难以理解，要么缺少关键细节导致调试时一头雾水。今天，我就结合自己十多年摸爬滚打的经验，带你从零开始，用C++彻底吃透扩展欧几里得算法。我们不只满足于写出能跑的代码，更要搞清楚每一行代码背后的数学原理，以及在实际编码中会遇到哪些坑，如何优雅地避开它们。

2. 算法核心原理与递归实现拆解

2.1 从欧几里得到扩展：关键等式的递推关系

普通的欧几里得算法基于一个核心定理：gcd(a, b) = gcd(b, a % b)。它通过递归将问题规模不断缩小，直到b == 0，此时a就是最大公约数。扩展欧几里得算法在此基础上，要求我们不仅算出gcd，还要记录下系数x和y。

算法的突破口在于观察递归的相邻两层状态。假设我们正在处理gcd(a, b)，并且我们已经通过递归求出了下一层gcd(b, a % b)的解x1和y1，即：b * x1 + (a % b) * y1 = gcd

现在，我们需要用x1和y1来表示本层的x和y，使得a * x + b * y = gcd成立。这里的关键是利用取模运算的定义：a % b = a - (a / b) * b，其中/是整数除法（向下取整）。

将a % b的表达式代入上面的等式：gcd = b * x1 + (a - (a / b) * b) * y1展开并重新分组：gcd = b * x1 + a * y1 - (a / b) * b * y1gcd = a * y1 + b * (x1 - (a / b) * y1)

将这个结果与目标形式a * x + b * y = gcd对比，我们立刻可以得到本层系数的递推关系：x = y1y = x1 - (a / b) * y1

这就是整个算法的灵魂。递归的终点（基线条件）是当b == 0时，此时gcd(a, 0) = a。对应的系数非常直观：我们需要一组x,y满足a * x + 0 * y = a。显然，x = 1,y = 0（实际上y可以是任意整数，但取0最方便）就是一组特解。

2.2 递归版C++实现与逐行解析

基于上面的推导，我们可以写出最经典的递归实现。这个版本结构清晰，直接反映了数学推导过程。

#include <iostream> using namespace std; // 扩展欧几里得算法 (递归版本) // 输入：整数 a, b // 输出：三元组 (g, x, y)，其中 g = gcd(a, b)，且满足 a*x + b*y = g int extendedGcd(int a, int b, int &x, int &y) { // 基线条件：当 b 为 0 时，gcd(a, 0) = a if (b == 0) { x = 1; // 系数 x 设为 1 y = 0; // 系数 y 设为 0 (可以是任意值，0最方便) return a; // 返回最大公约数 a } // 递归调用，处理子问题 gcd(b, a % b) // x1, y1 是子问题的解，满足：b*x1 + (a%b)*y1 = g int x1, y1; int g = extendedGcd(b, a % b, x1, y1); // 利用递推关系更新当前层的解 // x = y1 // y = x1 - (a / b) * y1 x = y1; y = x1 - (a / b) * y1; return g; // 返回计算出的最大公约数 } int main() { int a = 56, b = 15; int x, y, g; g = extendedGcd(a, b, x, y); cout << "gcd(" << a << ", " << b << ") = " << g << endl; cout << "系数 x, y 满足 " << a << "*" << x << " + " << b << "*" << y << " = " << g << endl; // 验证：56 * (-4) + 15 * 15 = -224 + 225 = 1 return 0; }

代码关键点解析与注意事项：

1. 参数传递：x和y通过引用（int &）传递。这是必须的，因为函数需要修改调用者传入的变量来返回计算结果。如果传值，结果就无法带出函数。
2. 递归顺序：一定要先进行递归调用extendedGcd(b, a % b, x1, y1)，获得子问题的解x1,y1和g，然后再根据递推公式计算当前层的x和y。这个顺序不能颠倒。
3. 整数除法：(a / b)是C++中的整数除法，当a和b一正一负时，它的行为是“向零取整”（truncate toward zero）。这与数学上常用的“向下取整”（floor）不同。在绝大多数求逆元或解方程的应用场景中，a和b都是正整数，所以没有问题。但如果需要考虑负数，这里就是一个潜在的坑，需要特别处理。
4. 解的多样性：这个函数返回的(x, y)只是满足等式a*x + b*y = g的无穷多组整数解中的一组“特解”。通解形式为：x' = x + (b/g) * ty' = y - (a/g) * t，其中t是任意整数。

注意：递归版本虽然直观，但在处理极大整数（比如在密码学中常见的1024位大数）时，可能存在函数调用栈溢出的风险。对于生产环境或性能敏感的场景，迭代版本是更稳妥的选择。

3. 迭代版实现、原理与性能对比

3.1 迭代算法的状态转换推导

递归的本质是函数调用栈，我们可以手动模拟这个栈，用循环来实现，这就是迭代法。迭代法的核心是维护一组状态变量，并在循环中不断更新它们，直到达到终止条件。

我们定义两组状态变量：(old_r, r)表示余数序列，(old_s, s)表示系数x的序列，(old_t, t)表示系数y的序列。初始化如下：

• old_r = a, r = b
• old_s = 1, s = 0（对应基线条件x=1）
• old_t = 0, t = 1（对应基线条件y=0）

在每一步迭代中，我们计算商quotient = old_r / r和新的余数new_r = old_r - quotient * r。这与欧几里得算法中的old_r % r等价。关键是如何更新系数s和t？我们需要让以下关系在每一步都成立：a * old_s + b * old_t = old_ra * s + b * t = r

当计算完新的余数new_r后，我们需要更新系数，使得a * new_s + b * new_t = new_r成立。通过推导（过程类似于递归版的递推），我们可以得到更新公式：new_s = old_s - quotient * snew_t = old_t - quotient * t

然后，我们将状态向前滚动：old_r = r, r = new_r，同时old_s = s, s = new_s，old_t = t, t = new_t。循环继续，直到r == 0。此时，old_r就是gcd(a, b)，而old_s和old_t就是对应的系数x和y。

3.2 迭代版C++代码实现

#include <iostream> #include <tuple> // 用于返回多个值 using namespace std; // 扩展欧几里得算法 (迭代版本) // 返回一个元组 (g, x, y) tuple<int, int, int> extendedGcdIterative(int a, int b) { // 初始化状态变量 // r: 当前余数， old_r: 上一次的余数 int old_r = a, r = b; // s: 当前x系数， old_s: 上一次的x系数 int old_s = 1, s = 0; // t: 当前y系数， old_t: 上一次的y系数 int old_t = 0, t = 1; while (r != 0) { int quotient = old_r / r; // 计算商 // 更新余数：new_r = old_r - quotient * r // 等价于 new_r = old_r % r int new_r = old_r - quotient * r; old_r = r; r = new_r; // 更新系数 s (对应 x) int new_s = old_s - quotient * s; old_s = s; s = new_s; // 更新系数 t (对应 y) int new_t = old_t - quotient * t; old_t = t; t = new_t; } // 循环结束时，old_r 是 gcd(a, b) // old_s 和 old_t 是系数 x 和 y // 注意：当 b=0 时，循环不会执行，直接返回 (a, 1, 0) return make_tuple(old_r, old_s, old_t); } int main() { int a = 56, b = 15; auto [g, x, y] = extendedGcdIterative(a, b); // C++17 结构化绑定 cout << "gcd(" << a << ", " << b << ") = " << g << endl; cout << "系数 x, y 满足 " << a << "*" << x << " + " << b << "*" << y << " = " << g << endl; return 0; }

迭代版优势与细节剖析：

1. 无栈溢出风险：这是迭代法最大的优点，特别适合处理大整数或深度递归可能成问题的环境。
2. 性能：通常迭代版本比递归版本有稍好的性能，因为避免了函数调用的开销。但在现代编译器的优化下，对于尾递归（本算法是尾递归），递归和迭代的性能差距可能很小。
3. 状态维护：代码中同时维护了(old_r, r)、(old_s, s)、(old_t, t)三组状态。更新时必须先计算所有新值，再统一进行旧值替换。如果更新顺序错了，比如先更新了old_s再用它去算new_s，结果就会完全错误。
4. 负数处理：和递归版一样，quotient = old_r / r的向零取整特性在涉及负数时需要注意。在纯数学推导或某些数论应用中，我们可能需要的是欧几里得除法（余数非负），这时就需要自己实现一个特殊的除法运算。

3.3 递归与迭代的选择建议

对于学习和理解算法，递归版本因其与数学推导的紧密对应而更胜一筹。它清晰地展示了“分治”和“状态回溯”的思想。对于实际项目开发，尤其是需要处理不确定规模输入或对稳定性要求极高的场景（如密码学库），我强烈推荐使用迭代版本。它的行为更可预测，没有递归深度的限制，代码虽然稍长，但逻辑流是线性的，更容易进行严格的正确性证明和测试。

4. 核心应用场景：乘法逆元与线性同余方程求解

掌握了EEA的实现，我们来看看它最经典的两个应用。这才是体现它价值的战场。

4.1 求解模意义下的乘法逆元

乘法逆元是模运算中的一个核心概念。如果存在整数x，使得(a * x) % m == 1，那么x就是a在模m下的乘法逆元，记作a^(-1) mod m。逆元存在的充要条件是gcd(a, m) == 1，即a和m互质。

如何用EEA求逆元？观察等式a*x + m*y = 1。这正是EEA的标准形式！当gcd(a, m) == 1时，EEA求出的x就是满足a*x ≡ 1 (mod m)的一个特解。

但EEA直接给出的x可能不在0到m-1的范围内。我们需要将它调整到最小非负整数解。方法是对x取模m：x = (x % m + m) % m。这个操作确保了结果在[0, m-1]之间。

// 使用扩展欧几里得算法求乘法逆元 // 返回值：a 在模 m 下的逆元。如果逆元不存在（即 gcd(a, m) != 1），返回 -1。 int modInverse(int a, int m) { int x, y; int g = extendedGcd(a, m, x, y); // 使用之前定义的递归或迭代函数 // 逆元存在的条件 if (g != 1) { // cout << "逆元不存在，因为 gcd(" << a << ", " << m << ") = " << g << endl; return -1; // 用-1表示不存在，也可用异常或bool返回值 } // 将 x 调整到 [0, m-1] 范围内 int inv = (x % m + m) % m; return inv; } // 示例：求 3 在模 11 下的逆元 // 因为 3*4 = 12 ≡ 1 (mod 11)，所以逆元是4。 int main() { int a = 3, m = 11; int inv = modInverse(a, m); if (inv != -1) { cout << a << " 在模 " << m << " 下的逆元是: " << inv << endl; cout << "验证: (" << a << " * " << inv << ") % " << m << " = " << (a * inv) % m << endl; } return 0; }

实操心得：在竞赛编程中，当模数m是质数时，求逆元有更快的费马小定理方法（用快速幂计算a^(m-2) mod m）。但当m不是质数，或者需要同时求gcd和系数时，EEA是唯一通用且高效的方法。务必记住检查gcd(a, m) == 1，否则后续计算都是无意义的。

4.2 求解线性同余方程 a*x ≡ b (mod m)

线性同余方程是形如a * x ≡ b (mod m)的方程。我们可以将其转化为EEA的形式。

方程等价于存在整数y，使得a*x - m*y = b。设d = gcd(a, m)。方程有解的充要条件是d能整除b。

求解步骤：

1. 令g = gcd(a, m)。如果b % g != 0，则方程无解。
2. 用EEA求解a*x' + m*y' = g，得到一组特解(x', y')。
3. 原方程的一个特解为x0 = x' * (b / g)。
4. 方程的通解为x = x0 + (m / g) * t，其中t为任意整数。
5. 通常我们关心最小非负整数解，可以通过取模得到：x = (x0 % (m/g) + (m/g)) % (m/g)。

// 求解线性同余方程 a*x ≡ b (mod m) // 返回一个 pair<bool, int>， first 表示是否有解，second 表示最小非负整数解（如果有解） pair<bool, int> solveLinearCongruence(int a, int b, int m) { int x, y; int g = extendedGcd(a, m, x, y); if (b % g != 0) { return {false, 0}; // 无解 } // 调整系数 a /= g; b /= g; m /= g; // 此时 gcd(a, m) = 1，可以直接求逆元 // x0 = (b * a^{-1}) mod m // 先求 a 在模 m 下的逆元 int inv_a = modInverse(a, m); // 复用之前的函数，此时逆元一定存在 int x0 = (b * inv_a) % m; // 通解是 x0 + m*t，最小非负解就是 x0 本身（因为已经模过m） // 但为了确保非负，再做一次调整 x0 = (x0 % m + m) % m; return {true, x0}; } int main() { int a = 14, b = 30, m = 100; auto [hasSolution, sol] = solveLinearCongruence(a, b, m); if (hasSolution) { cout << "方程 " << a << "*x ≡ " << b << " (mod " << m << ") 的最小非负解是: " << sol << endl; cout << "验证: (" << a << " * " << sol << ") % " << m << " = " << (a * sol) % m << endl; } else { cout << "方程无解。" << endl; } return 0; }

5. 边界处理、常见陷阱与调试技巧

即使理解了原理，实现EEA时仍然会遇到一些隐蔽的坑。下面是我在多年实践中总结出的关键问题和解决方案。

5.1 整数溢出问题

这是EEA实现中最常见也最危险的问题。在递推公式y = x1 - (a / b) * y1中，(a / b)是整数除法，但(a / b) * y1可能发生溢出，尤其是当a,b,y1都很大时（例如在RSA中处理大素数）。在递归版本中，中间计算过程的溢出会导致最终结果完全错误。

解决方案：

1. 使用更大范围的整数类型：在C++中，如果int可能不够，使用long long。对于密码学应用，则需要专门的任意精度整数库（如GMP）。
2. 注意运算顺序：有时调整计算顺序可以延缓溢出的发生，但根本之道还是使用足够宽的类型。
3. 迭代版本的潜在优势：迭代版本中，我们计算new_s = old_s - quotient * s，同样存在溢出风险。但迭代版本的状态变量更直观，便于在每一步插入溢出检查。

// 使用 long long 防止溢出的递归版本示例 long long extendedGcdLL(long long a, long long b, long long &x, long long &y) { if (b == 0) { x = 1; y = 0; return a; } long long x1, y1; long long g = extendedGcdLL(b, a % b, x1, y1); x = y1; y = x1 - (a / b) * y1; // 这里 (a/b)*y1 仍可能溢出 long long！ return g; } // 对于极端大的数，需要考虑使用 __int128（如果编译器支持）或大数库。

5.2 负数输入的处理

标准的EEA通常描述为非负整数。当输入a或b为负数时，C++的%运算符和/运算符的行为可能不符合数论中的欧几里得除法定义（要求余数非负）。这会导致递归或迭代过程中出现意外的结果。

解决方案：在算法开始前，先将负数转换为正数处理，并在最后调整系数的符号。

int extendedGcdWithNegative(int a, int b, int &x, int &y) { // 处理负数，确保递归在非负参数下进行 bool swapFlag = false; if (a < 0 && b < 0) { a = -a; b = -b; } else if (a < 0) { a = -a; swapFlag = true; } else if (b < 0) { b = -b; swapFlag = true; } // 调用标准非负版本 int g = extendedGcd(a, b, x, y); // 根据原始符号调整系数 if (swapFlag) { // 如果只有一个原始参数为负，需要调整一个系数的符号 // 具体规则：保证 a*x + b*y = g 成立，其中 a, b 是原始值（可能为负） // 一种简单方法：如果a为负，则x变号；如果b为负，则y变号。 // 但更通用的方法是：在递归基线条件中，用原始a,b的符号逻辑。 // 这里提供一个简化思路：在调用前记录a,b符号，调用后修正。 // 更健壮的做法是修改算法，使其能正确处理带符号的除法和取模。 } return g; }

实际上，更干净的做法是使用一个自定义的取模函数，使其始终返回非负余数。

// 欧几里得取模，结果总是非负 int euclideanMod(int a, int b) { int r = a % b; if (r < 0) { r += (b > 0) ? b : -b; // 加上除数的绝对值 } return r; } // 欧几里得除法，商向负无穷取整 int euclideanDiv(int a, int b) { if (b == 0) throw runtime_error("Division by zero"); if (b < 0) { a = -a; b = -b; } if (a >= 0) { return a / b; } else { // 对于负数a，确保商向负无穷取整 return -((-a + b - 1) / b); } }

在EEA实现中，将所有的a % b替换为euclideanMod(a, b)，将所有的a / b替换为euclideanDiv(a, b)，就能得到一个对负数输入鲁棒的版本。不过这会增加一些计算开销。

5.3 零值输入的处理

如果输入a和b都是零，那么最大公约数gcd(0, 0)在数学上是未定义的（通常约定为0或认为不存在）。在EEA中，我们试图找到x,y使得0*x + 0*y = g，这除了g=0外无解。大多数实现将gcd(0,0)返回0，并任意指定一组系数（如x=1, y=0）。如果a非零而b为零，那么gcd(a,0)=|a|，解为(x = sign(a), y = 0)。

建议：在函数入口处添加对a和b均为零的判断，根据应用需求决定是返回一个特定值（如(0,1,0)）还是抛出异常。

5.4 调试与验证技巧

1. 基础验证：对于任何实现，首先用几组小数字测试，并手动验证等式a*x + b*y == g是否成立。例如测试(56, 15),(48, 18),(17, 5)等。
2. 随机测试：生成大量随机整数对（包括正数、负数和零），用你的EEA实现和另一个可靠的计算库（如Python的math.gcd和pow(a, -1, m)求逆元）进行对比测试。
3. 逆元测试：对于互质的随机数对(a, m)，计算逆元inv，验证(a * inv) % m == 1。
4. 方程解测试：对于随机生成的(a, b, m)，用你的solveLinearCongruence求解，并将解代入原方程验证。
5. 边界测试：专门测试(0, n),(n, 0),(0,0),(1, 大数),(大数, 1),(负数, 正数)等情况。
6. 打印中间状态：如果结果不对，在递归或迭代函数中添加打印语句，输出每一步的a,b,x,y值，与手动计算的过程对比，这是定位逻辑错误最有效的方法。

6. 项目实战：构建一个简单的模运算工具库

将EEA及其应用封装成一个易用的工具库，是巩固理解并提升代码复用性的好方法。下面是一个简单的头文件库示例。

// File: modular_arithmetic.h #ifndef MODULAR_ARITHMETIC_H #define MODULAR_ARITHMETIC_H #include <tuple> #include <stdexcept> namespace ModularArithmetic { // 使用欧几里得取模和除法（处理负数） int euclideanMod(int a, int b); int euclideanDiv(int a, int b); // 扩展欧几里得算法 (迭代版，处理负数) // 返回 (gcd, x, y) 满足 a*x + b*y = gcd std::tuple<int, int, int> extendedGcd(int a, int b); // 求乘法逆元，模 m // 如果逆元不存在（gcd(a, m) != 1），抛出 std::invalid_argument 异常 int modInverse(int a, int m); // 求解线性同余方程 a*x ≡ b (mod m) // 返回一个解 x (在 [0, m/gcd(a,m)-1] 范围内) // 如果无解，抛出 std::invalid_argument 异常 int solveCongruence(int a, int b, int m); } // namespace ModularArithmetic #endif // MODULAR_ARITHMETIC_H

// File: modular_arithmetic.cpp #include "modular_arithmetic.h" namespace ModularArithmetic { int euclideanMod(int a, int b) { if (b == 0) { throw std::invalid_argument("Modulus cannot be zero"); } int r = a % b; if (r < 0) { r += (b > 0) ? b : -b; } return r; } int euclideanDiv(int a, int b) { if (b == 0) { throw std::invalid_argument("Division by zero"); } if (b < 0) { a = -a; b = -b; } if (a >= 0) { return a / b; } else { // 向负无穷取整 return -((-a + b - 1) / b); } } std::tuple<int, int, int> extendedGcd(int a, int b) { int old_r = a, r = b; int old_s = 1, s = 0; int old_t = 0, t = 1; while (r != 0) { int quotient = euclideanDiv(old_r, r); int new_r = old_r - quotient * r; old_r = r; r = new_r; int new_s = old_s - quotient * s; old_s = s; s = new_s; int new_t = old_t - quotient * t; old_t = t; t = new_t; } // 确保 gcd 非负 if (old_r < 0) { old_r = -old_r; old_s = -old_s; old_t = -old_t; } return std::make_tuple(old_r, old_s, old_t); } int modInverse(int a, int m) { if (m <= 0) { throw std::invalid_argument("Modulus must be positive"); } auto [g, x, y] = extendedGcd(a, m); if (g != 1) { throw std::invalid_argument("Inverse does not exist (a and m are not coprime)"); } // 将 x 调整到 [0, m-1] 范围 int inv = euclideanMod(x, m); return inv; } int solveCongruence(int a, int b, int m) { if (m <= 0) { throw std::invalid_argument("Modulus must be positive"); } auto [g, x, y] = extendedGcd(a, m); if (b % g != 0) { throw std::invalid_argument("No solution to the congruence"); } // 化简方程 a /= g; b /= g; m /= g; // 此时 a 和 m 互质 int inv_a = modInverse(a, m); // 这里g=1，逆元一定存在 int x0 = (b * inv_a) % m; // 结果已经在 [0, m-1] 范围内，因为 m 现在是正数，且取模结果非负 // 但为了使用我们自己的 euclideanMod 确保非负： x0 = euclideanMod(x0, m); return x0; } } // namespace ModularArithmetic

这个工具库提供了健壮的、带错误处理的EEA及相关函数。它使用了处理负数的欧几里得除法和取模，使得函数在更广泛的输入范围内行为一致。通过封装成命名空间和抛出异常，调用方可以更容易地集成和处理错误。

在实际使用中，你可以根据需要扩展这个库，比如添加对long long类型的支持，或者增加中国剩余定理（CRT）的求解函数，后者正是基于EEA求解一系列线性同余方程组。