[算法笔记] 扩展欧几里得算法

1. 扩展欧几里得算法

扩展欧几里得算法用于求解这样一个问题：给定两个非零整数 \(a\) 和 \(b\) ，求一组整数解 \((x,y)\) 使得 \(ax+by = gcd(a,b)\) 成立。

易知：

\[ax+by = gcd(a,b) = gcd(b,a\%b) = ... = gcd(a',0) = a' \\ a'x+0y = gcd(a',0) = a' => x=1 \]

下面说明求解过程。

\[\begin{aligned} gcd(a,b) &= gcd(a,a\%b) \\ ax+by &= gcd(a,b) \end{aligned} \]

所以设：

\[\begin{aligned} ax_1+by_1 &= gcd(a,b) \\ bx_2+(a\%b)y_2 &= gcd(b,a\%b) \end{aligned} \]

因此有：

\[bx_2+(a\%b)y_2 = ax_1+by_1 \]

又因为（/表示整除）：

\[a\%b = a - (a/b)*b \]

所以：

\[\begin{aligned} bx_2+(a - (a/b)*b)y_2 &= ax_1+by_1 \\ ay_2+b(x_2 - (a/b)y_2) &= ax_1+by_1 \end{aligned} \]

所以：

\[\begin{aligned} x_1 &= y_2 \\ y_1 &= x_2 - (a/b)y_2 \end{aligned} \]

重复进行以上步骤，到达递归边界时必然可以得一组解。

回到原式 \(ax+by = gcd(a,b)\) ，可得到递归边界：

when b=0: ax = gcd(a,0) = a 
so: x=1, y=random

代码实现：

//solve: ax + by = gcd(a,b)
int extGcd(int a, int b, int &x, int &y)
{
    if (b == 0)
    {
        // the equation has more than one solution
        // here, y can be a random val, which decide {x,y} at the last
        x = 1, y = (int)random() % 10;
        return a;
    }
    int gcd = extGcd(b, a % b, x, y);
    int t = y;
    y = x - (a / b) * y;
    x = t;
    return gcd;
}

下面继续探讨如何得出 ax+by=gcd(a,b) 的所有解。

先说结论：

\[x' = x_0 + \frac{b}{gcd(a,b)}K \\ y' = y_0 - \frac{a}{gcd(a,b)}K \]

下面看求解过程。

设新的解为 \(x_0+s_1\) 和 \(y_0 - s_2\) :

\[\begin{aligned} a(x_0+s_1)+b(y_0-s_2) &= ax_0 + by_0 \\ as_1-bs_2 &= 0 \\ \frac{s_1}{s_2} &= \frac{b}{a} = \frac{b/gcd(a,b)}{a/gcd(a,b)} \end{aligned} \]

显然，\(\frac{b}{gcd(a,b)}\) 和 \(\frac{a}{gcd(a,b)}\) 是互质的（没有大于 1 的公约数）。

所以取：：

\[s_1 = \frac{b}{gcd(a,b)}*K \\ s_2 = \frac{a}{gcd(a,b)}*K \]

其中，\(K\) 是任意的整数。

那么问题又来了，方程中的 \(x\) 的最小非负整数解是什么呢？

从通解式 \(x' = x_0 + \frac{b}{gcd(a,b)}K\) 上看，应当是 \(x' \% \frac{b}{gcd(a,b)} = x_0 \% \frac{b}{gcd(a,b)}\) 。

但是由于在递归边界时，\(y\) 可以取任意值，所得的特解 \(x_0\) 可能为负，不能保证 \(x_0 \% \frac{b}{gcd(a,b)}\) 是非负的。

如果 \(x_0 \% \frac{b}{gcd(a,b)}\) 是负数，那么其取值范围是：

\[(-\frac{b}{gcd(a,b)},0) \]

所以，\(x\) 的最小非负整数解为：

\[\begin{aligned} x_0 \% \frac{b}{gcd(a,b)} + \frac{b}{gcd(a,b)} \quad &if \quad x_0<0 \\ x_0 \% \frac{b}{gcd(a,b)} \quad &if \quad x_0 ≥ 0 \end{aligned} \]

综合一下，也就是：

\[(x_0 \% \frac{b}{gcd(a,b)} + \frac{b}{gcd(a,b)}) \% \frac{b}{gcd(a,b)} \]

2. ax+by=c求解

\(ax+by=c\) 有解的充要条件是 \(c \% gcd(a,b) == 0\) 。

如果 \(ax+by=gcd(a,b)\) 有一组解为：\((x_0,y_0)\)，即：

\[ax_0+by_0=gcd(a,b) \]

两边同乘以 \(\frac{c}{gcd(a,b)}\)：

\[a\frac{cx_0}{gcd(a,b)} + b\frac{cy_0}{gcd(a,b)} = c \]

所以：

\[(\frac{cx_0}{gcd(a,b)},\frac{cy_0}{gcd(a,b)}) \]

是方程 \(ax+by=c\) 的一个特解。

同理可得：

\[\begin{aligned} & a(x'+s_1)+b(y'-s_2) = c \\ & ax'+by'= c \\ & \frac{s_1}{s_2} = \frac{b}{a} = \frac{b/gcd(a,b)}{a/gcd(a,b)} \end{aligned} \]

所以，通解为：

\[x = \frac{cx_0}{gcd(a,b)} + b/gcd(a,b)*K \\ y = \frac{cy_0}{gcd(a,b)} - a/gcd(a,b)*K \]

3. 同余式 ax ≡ c (MOD m) 求解

首先解释何为「同余式」，给定整数 \(a\) 和 \(b\) ，如果说二者「模 \(m\) 同余」，即满足：\(a\%m == b\%m\)，

（或者是满足 \((a-b)\%m==0\) ）。

显然，每个整数各自与 \([0,m)\) 内的唯一整数满足同余关系。

现在需要求解的问题是：同余式 \(ax = c \quad (MOD \quad m)\) 的解。

由同余式的定义，可得：

\[(ax-c)\%m=0 \]

那么我们设：

\[\begin{aligned} ax-c &= my \\ ax - my &= c \end{aligned} \]

令 \(y = -y\)，可得：

\[ax+my=c \]

使用上面第二小节的结论可得：

1. 如果 \(c \% gcd(a,m) \neq 0\) ，则同余式 \(ax=c(MOD\quad m)\) 无解。

2. 如果 \(c\%gcd(a,m)=0\)，那么设 \(ax_0+my_0 = gcd(a,m)\) ，则有：

\[\begin{aligned} x &= \frac{cx_0 + m*K}{gcd(a,m)},\quad K = 0,1,..., gcd(a,m)-1 \end{aligned} \]

4. 逆元求解以及 (b/a)%m 计算

本文的最后一个问题：设 \(a,m\) 是整数，求 \(a\) 模 \(m\) 的逆元。

首先解释一下，何为「模 \(m\) 的逆元」：如果 \(a,b\) 满足：\(a*b = 1(MOD \quad m)\)，那么我们就说 \(a,b\) 互为模 \(m\) 的逆元。

易知，逆元还有以下等价定义：

\[\begin{aligned} (a*b-1)\%m &= 0 \\ a*b\%m &= 1 \\ a &= \frac{1}{b} (MOD \quad m) \\ b &= \frac{1}{a} (MOD \quad m) \end{aligned} \]

那逆元到底有什么用呢？当我们计算 \((a*b)\%m\) 时，编程实现时为了防止溢出会计算其等价式：\(((a\%m)*(b\%m))\%m\) 。但是计算 \((b/a)\%m\) 时，该编程方法就无能为力了。

如果可以找到 \(a\) 模 \(m\) 的逆元 \(x\)，即：\(a*x=1(MOD \quad m)\)，那么：

\[(b/a)\%m = (b*x)\%m \]

要求解 \(x\)，实质上就是求解：

\[ax=1(MOD\quad m) \]

根据第三小节的结论：

1. 如果 \(1\%gcd(a,m)\neq0\) ，即 \(gcd(a,m)\neq1\)，即 \(a,m\) 互质，无解。

2. 如果 \(gcd(a,m)=1\)，那么 \(ax \equiv 1(MOD \quad m)\) 在 (0,m) 上有唯一解：

求解 \(ax \equiv 1(MOD \quad m)\)，实质上是：

\[\begin{aligned} ax - 1 = my \quad &\Rightarrow \quad ax+my=1=gcd(a,m) \\ \end{aligned} \]

使用第一小节的扩展欧几里得算法，可得 \(x_0\) ，最后最小非负数解就是逆元：

\[a^{-1} = (x_0\%m+m)\%m \]

cpp代码：

int inverse(int a, int m)
{
    int x, y;
    int gcd = extGcd(a, m, x, y);
    return (x % m + m) % m;
}

5. 费马小定理

费马小定理：

设 \(m\) 是素数 ，且 \(a\%m \neq 0\)，那么 \(a^{m-1} \equiv 1(MOD \quad m)\) 。

显然 \(a*a^{m-2} \equiv 1 (MOD \quad m)\)，\(a^{m-2}\) 就是 \(a\) 模 \(m\) 的逆元，可使用快速幂求解。