扩展欧几里得算法:从GCD到模逆元与线性同余方程求解

发布时间:2026/7/17 1:53:21
扩展欧几里得算法:从GCD到模逆元与线性同余方程求解 1. 项目概述从GCD到扩展欧几里得如果你写过C或者对算法稍有涉猎那么“最大公约数”GCD这个概念你一定不陌生。它不仅是小学数学里的老朋友更是密码学、计算机图形学、乃至解决某些特定数学问题的基石。而求解GCD最经典、最高效的算法莫过于欧几里得算法也叫辗转相除法。这个算法的思想简洁到令人惊叹两个整数的最大公约数等于其中较小的数和两数相除余数的最大公约数。用代码实现往往只需要寥寥几行递归或循环。但今天我们要聊的远不止于此。停留在基础的GCD求解就像只学会了汽车的启动和刹车却不知道它还能挂挡、转向、开空调。扩展欧几里得算法才是这个古老算法在现代编程中真正大放异彩的部分。它不仅能求出GCD还能顺带找到一组整数解 (x, y)满足著名的贝祖等式a*x b*y gcd(a, b)。这个等式的解是解锁模逆元计算、求解线性同余方程、乃至RSA加密算法等高级应用的钥匙。我见过很多初学者能默写出欧几里得算法的代码但一提到扩展形式就感到困惑更别提在实际项目中灵活运用了。这就像手里有一把瑞士军刀却只用来开瓶盖。本文将带你从最基础的欧几里得算法实现开始一步步深入到扩展欧几里得算法的原理、推导和C实现并结合具体的编程实践场景让你不仅“知其然”更“知其所以然”最终能自信地将这些知识应用到你的代码中。2. 欧几里得算法原理与基础实现2.1 算法核心思想与数学证明欧几里得算法的核心基于一个非常直观的数学观察如果一个整数d能同时整除a和b那么它也一定能整除a - b甚至能整除a除以b的余数r其中a b*q r,0 r b。更形式化地说对于任意整数a和bb ! 0有gcd(a, b) gcd(b, a % b)为什么这个结论成立我们可以这样理解设d gcd(a, b)那么d能整除a也能整除b。对于a b*q r因为d能整除a和b*q那么它必然也能整除它们的差a - b*q也就是r。所以d是b和r的一个公约数。反过来任何b和r的公约数也一定能整除b*q r a因此也是a和b的公约数。所以a和b的公约数集合与b和r的公约数集合完全相同它们的最大公约数自然也就相等。这个性质的美妙之处在于它将求解两个较大数的GCD问题转化为了求解两个更小数的GCD问题因为r b。反复应用这个过程余数会越来越小最终当余数为0时上一步的除数就是我们要找的最大公约数。2.2 C递归与迭代实现对比理解了原理实现起来就水到渠成了。最直观的方式是递归。// 递归实现 int gcd_recursive(int a, int b) { if (b 0) { return a; } return gcd_recursive(b, a % b); }递归的代码几乎就是数学定义的直接翻译如果b为0那么a就是最大公约数因为任何数和0的最大公约数就是它本身否则递归计算gcd(b, a % b)。代码简洁优雅但对于极深层的递归比如处理非常大的整数存在栈溢出的风险。因此迭代实现是更稳健、效率也丝毫不差的选择。// 迭代实现 int gcd_iterative(int a, int b) { while (b ! 0) { int temp b; b a % b; a temp; } return a; }迭代版本使用一个while循环不断用(a, b)更新为(b, a % b)直到b变为0此时的a即为结果。它避免了函数调用的开销是生产环境中更常用的形式。注意在实际编码中我们通常需要处理负数。根据定义最大公约数总是非负的。一个简单的处理方法是在计算前对输入取绝对值a abs(a); b abs(b);。或者在取模运算时C的%运算符会保留被除数的符号为了确保余数非负可以稍作调整r (a % b b) % b。但在纯粹的GCD计算中由于算法只关心余数是否为0符号通常不影响最终结果gcd(a, b) gcd(|a|, |b|)所以很多简单实现会忽略这一点。不过在扩展欧几里得算法中正确处理符号就至关重要了。2.3 边界条件与性能考量边界条件b为 0根据定义gcd(a, 0) |a|。我们的递归基线条件或循环终止条件正是基于此。a和b均为 0gcd(0, 0)在数学上通常是未定义的因为任何数都是0的约数。我们的实现无论是递归还是迭代在b0时会直接返回a如果此时a也是0则返回0。这符合一些编程库如C17的std::gcd的约定但使用时需明确其含义。负数输入如前所述建议先取绝对值或确保算法能正确处理。性能考量 欧几里得算法的时间复杂度大致是O(log(min(a, b)))。更精确地说它是关于输入数值位数的线性时间算法效率非常高。对于绝大多数应用场景其性能都不是瓶颈。然而在处理大整数如上百位的整数常见于密码学时取模运算%的成本会变得显著。虽然算法复杂度依然优秀但底层的大数运算库如GMP的效率会成为关键。一个有趣的优化是二进制GCD算法也称为Stein算法。它利用移位和减法来代替耗时的取模运算对于硬件不支持快速取模指令的环境或非常大的整数可能更有优势。其核心思想是利用以下性质gcd(a, b) 2 * gcd(a/2, b/2)如果a和b都是偶数。gcd(a, b) gcd(a/2, b)如果a是偶数b是奇数反之亦然。gcd(a, b) gcd(|a-b|, min(a, b))如果两者都是奇数。3. 扩展欧几里得算法从GCD到贝祖等式3.1 贝祖等式与算法目标基础欧几里得算法给了我们GCD但扩展欧几里得算法Extended Euclidean Algorithm能给我们更多。它的目标是找到整数x和y使得a*x b*y gcd(a, b)这个方程被称为贝祖等式Bézout‘s identity。找到的(x, y)被称为贝祖系数。为什么这个很重要试想一下如果gcd(a, b) 1即a和b互质那么贝祖等式就变成了a*x b*y 1。对等式两边同时模b我们得到a*x ≡ 1 (mod b)。看到了吗x就是a在模b意义下的乘法逆元这是RSA解密、椭圆曲线加密等密码学核心操作的基础。此外它还能直接用于求解形如a*x ≡ c (mod b)的线性同余方程。3.2 算法推导与递归关系理解扩展算法的推导紧跟着欧几里得算法的步骤。我们回顾一下辗转相除的过程a b * q1 r1 b r1 * q2 r2 r1 r2 * q3 r3 ... r_{n-2} r_{n-1} * q_n r_n (这里 r_n gcd) r_{n-1} r_n * q_{n1} 0我们知道最后一步r_n gcd(a, b)。扩展算法的思想是从最后一步r_n gcd开始反向代入将gcd表示为a和b的线性组合。更优雅的方式是采用与递归计算GCD同步进行的正向推导。我们假设在递归的每一步我们都能知道当前两个数与其GCD的线性组合关系。设我们递归计算gcd(a, b)并假设我们已经知道了对于(b, a % b)的贝祖系数(x1, y1)即b * x1 (a % b) * y1 gcd(b, a % b)由于gcd(a, b) gcd(b, a % b)所以上式右边就是我们要的gcd(a, b)。现在我们需要用a和b来表示这个等式。我们知道a % b a - b * (a / b)这里/是整数除法。代入上式b * x1 (a - b * (a / b)) * y1 gcd(a, b)整理一下a * y1 b * (x1 - (a / b) * y1) gcd(a, b)看如果我们令x y1y x1 - (a / b) * y1那么我们就得到了a*x b*y gcd(a, b)。这给了我们一个清晰的递归关系要计算(a, b)的贝祖系数(x, y)我们先递归计算(b, a % b)的系数(x1, y1)然后通过上述公式回溯计算出(x, y)。递归的基线条件是当b 0时此时gcd(a, 0) a等式a*x 0*y a显然有一组解(x1, y0)注意y可以是任意整数通常取0。3.3 C递归实现详解根据上面的推导我们可以写出清晰的递归C代码。函数将返回gcd(a, b)并通过引用参数返回贝祖系数x和y。#include iostream using namespace std; // 扩展欧几里得算法 (递归版本) // 返回 gcd(a, b)并通过 x, y 返回一组贝祖系数 int extended_gcd_recursive(int a, int b, int x, int y) { // 基线条件b 0 if (b 0) { x 1; y 0; // y可以是任意值通常取0 return a; // gcd(a, 0) a } int x1, y1; // 用于存储递归调用的结果 int gcd extended_gcd_recursive(b, a % b, x1, y1); // 回溯根据递归关系更新 x, y x y1; y x1 - (a / b) * y1; return gcd; } int main() { int a 30, b 18; int x, y, g; g extended_gcd_recursive(a, b, x, y); cout GCD( a , b ) g endl; cout Coefficients (x, y): ( x , y ) endl; cout Check: a * x b * y (a*x b*y) endl; return 0; }运行这段代码对于a30, b18你会得到输出GCD6,(x, y) (-1, 2)验证30*(-1) 18*2 6完全正确。关键点解析参数传递x和y通过引用int 传递使得函数内部对它们的修改能反映到外部。递归基线当b0我们设置x1, y0使得a*1 0*0 a成立并返回a作为GCD。回溯计算递归调用后我们得到了子问题(b, a%b)的解(x1, y1)和gcd。然后利用公式x y1,y x1 - (a/b)*y1计算出当前问题的解。解的多样性贝祖等式的解不是唯一的。如果(x, y)是一组解那么(x k*b/g, y - k*a/g)对于任意整数k也都是解其中g gcd(a, b)。我们的算法返回的是众多解中的一组。3.4 迭代实现与状态追踪递归实现虽然清晰但同样有栈深度限制。迭代实现扩展欧几里得算法需要一点技巧我们需要在循环中手动维护状态。思路是模拟递归的回溯过程但用循环正向完成。我们维护两组系数(x1, y1)对应当前的a和b(x2, y2)对应上一步的a和b。初始化时我们对应(a, b)和(b, a % b)的状态。// 扩展欧几里得算法 (迭代版本) int extended_gcd_iterative(int a, int b, int x, int y) { // 初始化状态 // 对应基线条件: s * a t * b a int s 1, t 0; // 对应 (x1, y1)初始为(1, 0) // 对应: u * a v * b b int u 0, v 1; // 对应 (x2, y2)初始为(0, 1) while (b ! 0) { int quotient a / b; // 更新 (a, b) 为 (b, a % b) int temp_a a; a b; b temp_a - quotient * b; // 即 temp_a % b // 同步更新系数 (s, t) 和 (u, v) // 新的系数应该是旧的 (u, v) 和 (s - quotient*u, t - quotient*v) int temp_s s; int temp_t t; s u; t v; u temp_s - quotient * u; v temp_t - quotient * v; } // 循环结束时a 就是 gcds 和 t 是贝祖系数 x s; y t; return a; // gcd }这个迭代版本看起来比递归复杂但它完全避免了递归开销。其核心是同时进行辗转相除和系数更新。在每一步除法a b*q r中我们不仅更新(a, b)也按照线性组合的关系更新两组系数(s, t)和(u, v)使得它们始终满足s*a_original t*b_original current_a和u*a_original v*b_original current_b。当循环结束b0时current_a就是GCD对应的系数(s, t)就是我们要的(x, y)。实操心得第一次接触迭代实现可能会觉得绕一个很好的理解方式是拿一个小例子比如a30, b18在纸上手动模拟一遍循环过程跟踪a, b, s, t, u, v每一步的变化你会发现它完美地模拟了递归的回溯逻辑。在性能要求高或递归深度可能很大的场景下迭代版本是更安全的选择。4. 核心应用场景与编程实践掌握了扩展欧几里得算法的实现我们来看看它能解决哪些实际问题。这些场景远比单纯求GCD要丰富和有趣。4.1 求解线性同余方程线性同余方程形式为a * x ≡ c (mod m)。意思是寻找整数x使得a*x除以m的余数等于c。解法步骤设g gcd(a, m)。根据贝祖等式扩展欧几里得算法可以找到s,t使得a*s m*t g。方程a*x ≡ c (mod m)有解的充要条件是g能整除c。因为等式左边a*x和右边c模m同余意味着a*x - c是m的倍数即a*x - c m*y整理得a*x m*(-y) c。这要求c必须是a和m的线性组合而a和m的所有线性组合构成的集合其最小正元素正是g。所以c必须是g的倍数。如果有解那么特解x0可以通过s * (c / g)得到注意模m。方程在模m意义下的所有解为x ≡ x0 k * (m / g) (mod m)其中k 0, 1, ..., g-1。也就是说有g个不同的解在模m的完全剩余系中。C代码示例// 求解线性同余方程 a*x ≡ c (mod m) // 返回一个解 x (在 [0, m) 范围内)如果无解则返回 -1 int solve_linear_congruence(int a, int c, int m) { int x, y; int g extended_gcd_iterative(a, m, x, y); // 使用迭代版本 if (c % g ! 0) { // 无解 return -1; } // 求出一个特解 // 注意extended_gcd 给出的 x 满足 a*x m*y g // 我们需要的是 a * (x * (c/g)) m * (y * (c/g)) c int x0 (x * (c / g)) % m; // 保证解在 [0, m) 范围内 if (x0 0) { x0 m; } // 注意这里返回的是模 m 意义下的一个特解。 // 通解为 x0 k * (m/g) k0,...,g-1 return x0; }4.2 计算模逆元模逆元是密码学和数论中的高频操作。对于整数a和模数m通常要求m 1a在模m下的乘法逆元是一个整数x满足a * x ≡ 1 (mod m)。记作a^{-1} mod m。存在条件a在模m下有乘法逆元的充要条件是gcd(a, m) 1即a与m互质。这其实就是线性同余方程a*x ≡ 1 (mod m)的特例。因此我们可以直接用扩展欧几里得算法求解。C代码示例// 计算 a 在模 m 下的乘法逆元 // 返回逆元如果不存在则返回 -1 (或者抛出异常根据设计) int mod_inverse(int a, int m) { int x, y; int g extended_gcd_iterative(a, m, x, y); if (g ! 1) { // a 和 m 不互质逆元不存在 // 在实际项目中这里可能需要更健壮的错误处理 return -1; } // 扩展欧几里得给出的 x 满足 a*x m*y 1 // 所以 a*x ≡ 1 (mod m)x 即为我们所求的逆元 // 但 x 可能为负数需要调整到 [0, m) 范围 int inv x % m; if (inv 0) { inv m; } return inv; }应用实例在RSA解密过程中私钥指数d的计算就是公钥指数e关于欧拉函数φ(n)的模逆元即d ≡ e^{-1} mod φ(n)。扩展欧几里得算法在这里扮演了核心角色。4.3 将分数转换为模意义下的整数在某些计算中比如组合数取模我们可能会遇到分数(a/b) mod m其中m是一个质数。直接进行除法在模运算中是没有定义的。我们需要将除法转换为乘以模逆元(a / b) mod m (a * b^{-1}) mod m其中b^{-1}是b模m的逆元可以用扩展欧几里得算法求出。例如计算(7 / 3) mod 5。首先求3 mod 5的逆元。因为gcd(3,5)1用扩展算法求得逆元为2因为3*26≡1 mod 5。所以(7/3) mod 5 (7 * 2) mod 5 14 mod 5 4。int fraction_mod(int numerator, int denominator, int mod) { int inv_denom mod_inverse(denominator, mod); if (inv_denom -1) { // 处理逆元不存在的情况例如分母与模数不互质 // 这可能意味着原始分数在模意义下无定义或需要其他处理 cerr Denominator has no modular inverse! endl; return -1; } return ( (numerator % mod) * inv_denom ) % mod; }4.4 实际编程中的封装与优化建议在实际的C项目中我们很少会裸写这些算法。更好的做法是将它们封装成工具函数或类并考虑以下优化点使用模板支持多种整数类型基础的int可能范围不够我们可以使用long long甚至__int128如果编译器支持或集成大数库如Boost.Multiprecision。templatetypename T T gcd_template(T a, T b) { while (b ! 0) { T t b; b a % b; a t; } return a; } templatetypename T T extended_gcd_template(T a, T b, T x, T y) { if (b 0) { x 1; y 0; return a; } T x1, y1; T g extended_gcd_template(b, a % b, x1, y1); x y1; y x1 - (a / b) * y1; return g; }处理负数与零在通用函数中必须仔细处理边界情况。对于GCD可以返回绝对值。对于扩展算法要明确约定当输入为负数或零时系数x,y的行为。迭代与递归的选择对于确定的小范围输入如竞赛编程递归的简洁性可能是首选。对于库函数或处理不可预知输入的程序迭代版本更稳健。与标准库结合C17在numeric头文件中引入了std::gcd和std::lcm。如果你的项目只需要GCD/LCM直接使用标准库是最好、最可移植的选择。但标准库没有提供扩展欧几里得算法这就需要你自己实现或使用第三方数学库如Boost。性能与可读性的权衡在绝大多数情况下这些算法的性能都不是瓶颈。代码的清晰性、正确性和健壮性更为重要。使用有意义的变量名添加必要的注释并进行充分的单元测试。5. 常见问题、调试技巧与进阶思考5.1 典型错误与排查清单即使理解了原理实现时也难免会遇到问题。下面是一些常见坑点及其解决方法问题现象可能原因排查与解决递归版本导致栈溢出输入数值过大递归深度太深。改用迭代实现。递归深度约为O(log n)对于极大的数如10^1000仍可能溢出。迭代版本无此问题。结果不正确尤其是系数1. 整数溢出。2. 递归回溯公式写错。3. 对负数取模处理不当。1. 使用更大范围的整数类型long long并在乘法时注意是否可能溢出考虑使用__int128或模运算特性。2. 仔细核对公式x y1; y x1 - (a/b)*y1;。确保a/b是整数除法。3. C/C中%运算结果符号与被除数相同。在扩展算法中我们通常需要非负余数。可以在递归基线或迭代循环中对输入a, b取绝对值并在最后根据原始符号调整系数x, y。一个简单策略是始终对a, b的绝对值进行计算最后如果原始a为负则x -x如果原始b为负则y -y。模逆元函数返回错误或无效值1. 没有检查gcd(a, m) 1。2. 返回的逆元没有调整到[0, m)范围。1.必须检查g 1。如果g ! 1逆元不存在应返回错误标识或抛出异常。2. 计算inv x % m后如果inv 0执行inv m。求解线性同余方程得到多个解如何获取全部算法只返回一个特解。记住通解公式x ≡ x0 k * (m / g) (mod m)k 0, 1, ..., g-1。如果需要所有解用一个循环生成即可。在模运算中计算(a/b) % m时结果不对直接使用了(a % m) / (b % m)这是错误的。必须将除法转换为乘以模逆元(a * mod_inverse(b, m)) % m。并且要确保b与m互质。5.2 调试技巧小数据模拟与打印中间状态对于扩展欧几里得这样逻辑稍复杂的算法调试的最佳方法是小数据模拟。手工计算选择一对小数字如a30, b18在纸上完整走一遍递归或迭代过程写下每一步的a, b, x, y或s, t, u, v。这是理解算法最有效的方式。打印日志在递归函数或迭代循环的关键步骤插入打印语句输出当前的参数和中间变量。对比你的输出和手工计算的结果不一致的地方就是bug所在。int extended_gcd_debug(int a, int b, int x, int y, int depth0) { string indent(depth*2, ); cout indent Call: gcd( a , b ) endl; if (b 0) { x 1; y 0; cout indent Base: return a , coeff (1, 0) endl; return a; } int x1, y1; int g extended_gcd_debug(b, a % b, x1, y1, depth1); x y1; y x1 - (a / b) * y1; cout indent Back: a a , b b , q (a/b) endl; cout indent (x1,y1)( x1 , y1 ) - (x,y)( x , y ) endl; return g; }单元测试编写一系列测试用例包括正数、负数、零、互质的数、倍数关系等用已知结果或暴力枚举对于小数字来验证你的函数。例如测试贝祖等式是否成立测试模逆元性质(a * inv) % m 1。5.3 从整数到多项式算法的泛化思考欧几里得算法及其扩展形式并不局限于整数。任何支持类似“带余除法”的代数结构称为欧几里得整环都可以应用此算法。一个重要的例子是多项式。对于多项式A(x)和B(x)我们也可以进行多项式除法A(x) B(x) * Q(x) R(x)其中deg(R) deg(B)。这里“大小”的比较由多项式的次数deg扮演了整数中“绝对值大小”的角色。同样地我们可以定义多项式的最大公因式GCD并使用辗转相除法求解。扩展算法也能找到多项式U(x),V(x)使得A(x)*U(x) B(x)*V(x) GCD(A(x), B(x))。这在信号处理如寻找滤波器、编码理论如Reed-Solomon码和计算机代数系统中都有应用。虽然实现细节更复杂涉及多项式运算但核心思想一脉相承。理解整数域上的算法为你学习这些更抽象领域的概念打下了坚实的基础。5.4 关于C标准库与第三方库如前所述C17提供了std::gcd和std::lcm定义在numeric头文件中。它们通常是编译器内置或高度优化的对于通用整数类型应优先使用。#include numeric #include iostream int main() { int a 30, b 18; std::cout GCD: std::gcd(a, b) std::endl; // C17 std::cout LCM: std::lcm(a, b) std::endl; // C17 return 0; }对于扩展欧几里得算法、模逆元等C标准库并未提供。如果你在项目中频繁使用这些数论工具可以考虑以下选项自己封装如本文所示编写一套健壮的模板函数。使用Boost库Boost的boost::integer模块提供了gcd、lcm以及extended_euclidean算法功能强大且经过充分测试。使用专门的数学库如GMPGNU Multiple Precision Arithmetic Library用于大数运算它自然也包含了这些数论函数。选择哪种方式取决于你的项目需求。对于学习、竞赛或轻量级应用自己实现是很好的锻炼。对于大型、要求高可靠性的项目使用成熟的第三方库是更稳妥的选择。最后再分享一个我个人的体会算法学习尤其是像扩展欧几里得这样精巧的算法切忌死记硬背代码。一定要亲手推导一遍数学公式用一个小例子在纸上或调试器中一步步跟踪直到你能够清晰地讲述出每一步“为什么”要这么做。当你真正理解了其背后的数论原理代码就自然而然地成为你表达思想的工具而不仅仅是记忆的负担。无论是解决一道编程竞赛题还是设计一个加密模块这份深刻的理解都会让你更加游刃有余。