扩展欧几里得算法 (exgcd)

扩展欧几里得算法 (exgcd)

问题 1

求不定方程 \(ax + by = \gcd(a,b)\) 的一组或特定范围内的解

特解

前提 / 特解

当 \(b = 0\) 时
原式 = \(ax + by = a\)，易得 \(x = 1, y = 0\)

引理

由欧几里得算法：\(\gcd(a,b) = \gcd(b,a \% b)\)

推式子

由裴蜀定理：

\[\gcd(a , b) = ax + by \]

用欧几里得算法拆开

\[\gcd(b,a\%b) = (a \% b)x_1 + by_1 %这时方程变化也就有了新的解 \]

在给 \(a\%b\) 换种表示方式

\[bx_1 + (a - \lfloor \frac{a}{b} \rfloor y_1) \]

去括号，并把 \(b\) 提出来得

\[a y_1 + b(x_1 - \frac{a}{b}y_1) \]

根据原方程各个未知数的系数
\(x=y_1 , y=x_1 - \frac{a}{b}y_1\)

实现

考虑递归，先求出下一轮的解 \(x_{1},y_{1}\)
递归终止条件：\(b = 0\) 时，我们可以求一个固定解 (\(x=1,y=0\))
然后回代得到特解 \((x_0,y_0)\)

普通写法

int exgcd(int a, int b, int &x, int &y) {
	if (b == 0) {
		x = 1 , y = 0;
		return a;
	}
	int r = exgcd(b, a % b, x, y);
	int tmp = x;
    x = y , y = temp - a / b * y
	return r;
}

优美的写法

区别

• exgcd(b , a%b , y , x)

• y = y - a / b * x

int exgcd(int a, int b, int &x, int &y) {
	if (b == 0) {
		x = 1 , y = 0;
		return a;
	}
	int r = exgcd(b, a % b, y, x);
	y = y - a / b * x;
	return r;
}

:::info[why]
递归传参和 y -= x*(a/b) 的逻辑，为什么能对应理论上的 \(x=y_1、y=x_1 - \lfloor \frac{a}{b} \rfloor y_1\) 推导公式，核心是引用传参的特性 + 理论推导的变量映射

1. 先回顾核心理论推导（关键前提）

我们的目标是求 \(a \cdot x + b \cdot y = \gcd(a,b)\) 的解 \((x,y)\)，根据欧几里得算法，先求解子问题：\(\gcd(b, a\%b) = b \cdot x_1 + (a\%b) \cdot y_1\) 其中 \(a\%b = a - \lfloor \frac{a}{b} \rfloor \cdot b\)，代入子问题公式并整理：

\[\begin{align*} \gcd(a,b) &= b \cdot x_1 + (a - \lfloor \frac{a}{b} \rfloor \cdot b) \cdot y_1 \\ &= a \cdot y_1 + b \cdot (x_1 - \lfloor \frac{a}{b} \rfloor \cdot y_1) \end{align*} \]

对比原问题 \(\gcd(a,b) = a \cdot x + b \cdot y\)，可得理论核心等式：

\[\begin{cases} x = y_1 \quad \text{（原问题x = 子问题的y₁）} \\ y = x_1 - \lfloor \frac{a}{b} \rfloor \cdot y_1 \quad \text{（原问题y = 子问题的x₁ - (a/b) * y）} \end{cases} \]

2. 拆解代码的执行逻辑（对应理论）

代码中用了引用传参（&x、&y），这是关键 —— 递归调用会直接修改传入的变量值

代码逻辑描述

函数定义为 void exgcd(LL a, LL b, LL &x, LL &y)，核心逻辑如下：

• 终止条件：当 b=0 时，方程变为 \(a \cdot x + 0 \cdot y = a\)，直接赋值 x=1, y=0；
• 递归步骤：当 b≠0 时，先调用 exgcd(b, a % b, y, x) 求解子问题，再执行 y -= x * (a / b) 更新 y 值。

步骤 1：递归调用 exgcd(b, a%b, y, x) 的含义

递归调用的参数是 (b, a%b, y, x)，而非 (b, a%b, x, y)，这是为了让子问题的解直接映射到理论中的 \(x_1、y_1\)：

• 子问题的目标：求 \(b \cdot x_1 + (a\%b) \cdot y_1 = \gcd(b, a\%b)\) 的解 \((x_1, y_1)\)。
• 代码中把当前的 y 传给子问题的 x 参数（即子问题的 \(x_1\) 会赋值给当前的 y）；
• 把当前的 x 传给子问题的 y 参数（即子问题的 \(y_1\) 会赋值给当前的 \(x\)）。

递归返回后：

• 当前的 x = 子问题的 (y_1)（正好对应理论公式 (x = y_1)，无需额外修改）；
• 当前的 y = 子问题的 (x_1)（这是中间值，还需要按公式修正）。

步骤 2：执行 y -= x * (a / b) 的含义

此时：

• 当前 x 是子问题的 $(y_1)；
• 当前 y 是子问题的 (x_1)；
• a/b 是 (\lfloor \frac {a}{b} \rfloor)（整数除法）。

所以 \(y -= x*(a/b)\) 等价于：\(y = x_1 - \lfloor \frac{a}{b} \rfloor \cdot y_1\) 这对应理论推导的第二个公式。

总结

1. 引用传参的巧思：递归时交换 y 和 x 传入，让子问题的 (y_1) 直接赋值给原问题的 x，满足 (x = y_1)；
2. 公式的直接映射：递归返回后，y 初始是子问题的 (x_1)，执行 y -= x*(a/b) 等价于 ( y = x_1 - (a / b) y_1 )；
3. 核心本质：代码通过参数顺序交换 + 简单的减法操作，把理论推导的两个公式用极简的方式实现，没有额外的临时变量，是扩展欧几里得算法的经典优化写法。
   :::

构造通解

如何构造呢？我们要让和不变，即始终等于 \(\gcd(a,b)\)。
即构造 \(ax + by = 0\)
不废话了，通解为

\[x = x_0 + \frac{b}{\gcd(a,b)}*k \]

\[y = y_0 - \frac{a}{\gcd(a,b)}*k \]

(\(k \ni \mathbb{Z}\))
为什么这么做 \(ax + by\) 不会变呢？

通解和特解不同的地方

\(\frac{b}{\gcd(a,b)}\) 和 \(\frac{a}{\gcd(a,b)}\)
注意到两式在分别乘上 \(a，b\) 后二式相等。

:::info[ 例子]

扩展欧几里得算法实例：求解 ( 8x + 6y = \gcd (8,6) )

1. 基础信息

• 目标方程：\(( 8x + 6y = \gcd(8,6) )\)
• 最大公约数计算：\(( \gcd(8,6) = 2 )\)，因此实际求解方程为 ( \(8x + 6y = 2\))

2. 递归求解特解（核心步骤）

终止条件：当 \(( b = 0 )\) 时，方程 \(( ax + 0·y = a )\) 的解为 \(( x=1, y=0 )\)；
回代公式：若子问题 \(( bx_1 + (a\%b)y_1 = \gcd(b,a\%b) )\) 的解为 ( (x_1,y_1) )，则原问题的解为：\([ x = y_1, \quad y = x_1 - \lfloor \frac{a}{b} \rfloor · y_1 ]\)

分步递归 & 回代

第一层（初始问题）：求解 (8x + 6y = 2)

• (a = 8, b = 6)，(b≠0)，计算 (a% b = 8%6 = 2)
• 递归求解子问题：(6x_1 + 2y_1 = 2)

第二层（子问题）：求解 (6x_1 + 2y_1 = 2)

• (a = 6, b = 2)，(b≠0)，计算 (a% b = 6%2 = 0)
• 递归求解子问题：(2x_2 + 0y_2 = 2)

第三层（终止条件）：求解 (2x_2 + 0y_2 = 2)

• (b = 0)，触发终止条件，解为：(x_2 = 1, y_2 = 0)（验证：(2×1 + 0×0 = 2)）

回代到第二层

• 子问题解：(x_2 = 1, y_2 = 0)
• 当前层 (a = 6, b = 2)，(\lfloor \frac {6}{2} \rfloor = 3)
• 回代公式计算：(x_1 = y_2 = 0)，(y_1 = x_2 - 3・y_2 = 1 - 0 = 1)
• 验证：(6×0 + 2×1 = 2)，正确。

回代到第一层（初始问题）

• 子问题解：(x_1 = 0, y_1 = 1)
• 当前层 (a = 8, b = 6)，(\lfloor \frac {8}{6} \rfloor = 1)
• 回代公式计算：(x = y_1 = 1)，(y = x_1 - 1・y_1 = 0 - 1 = -1)
• 验证：( 8×1 + 6×(-1) = 2 )，正确。

最终特解

方程 (8x + 6y = 2) 的一组特解为：( (x_0, y_0) = (1, -1) )

3. 构造通解

通解公式

设 ( d = \gcd (a,b) )，通解需满足 (ax + by = d)，因此：
[
\begin {cases}
x = x_0 + \frac {b}{d}・k \
y = y_0 - \frac {a}{d}・k
\end {cases}
]
（注：若 (y = y_0 + \frac {a}{d}・k) 会导致和改变，正确符号为减号）

代入实例计算

• ( d=2, \frac{b}{d} = \frac{6}{2}=3, \frac{a}{d} = \frac{8}{2}=4 )
• 通解：
  [
  \begin {cases}
  x = 1 + 3k \
  y = -1 - 4k
  \end {cases}
  ]
  （(k) 为任意整数）

通解验证（选不同 (k) 值）

(k) 值	( x = 1 + 3k )	( y = -1 - 4k )	验证 (8x + 6y)
-1	( 1 - 3 = -2 )	( -1 + 4 = 3 )	( 8×(-2) + 6×3 = -16 + 18 = 2 )
0	( 1 + 0 = 1 )	( -1 - 0 = -1 )	( 8×1 + 6×(-1) = 8 - 6 = 2 )
1	( 1 + 3 = 4 )	( -1 - 4 = -5 )	( 8×4 + 6×(-5) = 32 - 30 = 2 )
:::

问题 2

求不定方程 \(ax + by = c\) 的一组或特定范围内的解，即 P5656。

无解情况

由裴蜀定理，当 \(c \% \gcd(a,b) \neq 0\) 时，原方程无解。

有解时的转化

因为方程有解，所以我们可以把 \(\gcd(a,b)\) 转化成 \(c\) 即乘 \(c / \gcd(a,b)\)
显然，等号左边也要乘 \(c / \gcd(a,b)\)

求通解的各种问题

有正整数解

• 正整数解的数量
• 所有正整数解中 \(x\) 的最小值，所有正整数解中 \(y\) 的最小值。
• 所有正整数解中 \(x\) 的最大值，以及所有正整数解中 \(y\) 的最大值。

无正整数解

• 所有整数解中 \(x\) 的最小正整数值，\(y\) 的最小正整数值。

通解构造方式

令 \(d=\gcd(a,b)\), \(k \% [c / \gcd(a,b)] = 0\) 即乘了 \(\gcd(a,b)\), 已经变成 \(ax + by = c\) 的解，且 \(k \ni \mathbb{Z}\)。

\[x = x_0 + \frac{b}{d}*k \]

\[y = y_0 - \frac{a}{d}*k \]

然后就可以愉快地推式子啦！（解不等式大家应该都会把）
当 \(x > 0\) 时

\[x_0 + \frac {b}{d}*k > 0 \]

移项得

\[k > \frac {-x_0 * d}{b} \]

又由于 \(k \ni \mathbb{Z}\)
\(x\) 为正整数时，有

\[\lceil {\frac {-x_0 * d+1}{b}} \rceil \le k \]

同理得到当 \(y\) 为正整数时，有

\[k \le \lfloor {\frac {-y_0 * d-1}{a}} \rfloor \]

综合二式得到

\[\lceil {\frac {-x_0 * d + 1}{b}} \rceil \le k \le \lfloor {\frac {-y_0 * d - 1}{a}} \rfloor \]

若该不等式成立，则方程有正整数解，数量为 n 的个数（右边界 - 左边界 + 1）
若该不等式不成立，则方程有整数解但无正整数解
有了这两个，别的就很好推算了（见代码）。

 exgcd(a, b, x, y);
        x = x * c / d, y = y * c / d;
        ll minn = ceil((-x * d + 1) * 1.0 / b);
        ll maxx = floor((y * d - 1) * 1.0 / a);
        
        if (minn <= maxx) {
            ll cnt = maxx - minn + 1;
            ll x_max = x + b / d * maxx;
            ll x_min = x + b / d * minn;
            ll y_max = y - a / d * minn;
            ll y_min = y - a / d * maxx;
            cout << cnt << " " << x_min << " " << y_min << " " << x_max << " " << y_max << '\n';
        } else {
            ll x_min = x + b / d * minn;
            ll y_min = y - a / d * maxx;
            cout << x_min << " " << y_min << '\n';
        }

同余方程

题目描述

求关于 \(x\) 的同余方程 \(a x \equiv 1 \pmod {b}\) 的最小正整数解。

转化

设商是 \(y\)，则我们有 \(ax = by + 1\)。
移项得：

\[ax - by = 1 \]

这个式子又可以变成

\[ax + by = 1 \]

因为 \(y\) 是变量

最小正整数解

(x % mod + mod) % mod 即可
其余就按上面的做法做就行了

---

【模板】有理数取余

题目描述

求 \(bx\equiv a\pmod{19260817}\) 的解。

说明/提示

对于所有数据，保证 \(0\leq a \leq 10^{10001}\)，\(1 \leq b \leq 10^{10001}\)，且 \(a, b\) 不同时是 \(19260817\) 的倍数。

转化

由于 \(19260817\) 为质数，所以用 \(m\) 表示

\[my + a = bx \]

移项得

\[my - bx = a \]

无解判断

• 当 $ a \neq 0 $ 时：
  $ 0 < a < m $，而 $ m $ 的最小正整数倍是$ m $ 本身（$ 1 \times m = m $），显然 $ a < m $，因此 $ m $ 不可能整除 $ a $，不满足有解条件，同余方程无解。这是绝大多数场景，此时代码 if(b % mod == 0) 输出 angry（无解）是完全正确的。

• 当 $a = 0 $ 时：
  此时 $ m \mid 0 $（0 能被任意非零整数整除），理论上同余方程 $ bx \equiv 0 \pmod{m} $ 有解（且任意 $ x $ 都是解，因为 $ b = k \times m $，则 $ bx = k \times m \times x $，必然是 $ m $ 的倍数）。

补充代码判断的合理性总结

代码中的 if(b % mod == 0) 本质是一种快速的无解判定优化，其合理性来源于：

1. $ b % mod == 0 \implies \gcd(b, mod) = mod $；
2. 代码中 $ a $ 的取值范围 $ 0 \leq a < mod $，使得 $ mod \mid a $ 几乎无法满足（仅当 $ a=0 $ 时满足，无实际求解意义）；