扩展欧几里得算法

费马定理（Fermat’s theorem）使我们能够通过二进制指数法（binary exponentiation）在 $O(\log n)$ 的操作内计算模数除法的乘法逆元（modular multiplicative inverses），但这只适用于对质数取模。它的一个推广是欧拉定理（Euler’s theorem），此定理表明，如果 m 和 a 互质（即他们之间除了1，没有其他公约数），那么

$a^{\phi(m)} ≡ 1(\mod m)$

其中， $\phi(m)$ 是欧拉函数，定义为小于 $m$ 的正整数 $x$ 与 $m$ 互质的个数。当 $m$ 是一个质数时，所有小于 $m$ 的 $m-1$ 个正整数与 $m$ 都是互质的，所以 $\phi(m)=m-1$ ，这就得出了费马定理。

这允许我们通过 $\phi (m)$ 来计算 $a$ 的模除乘法逆元，记作 $a^{\phi (m)−1}$ ，但反过来，计算 $\phi (m)$ 并不那么快速：通常你需要获得 $m$ 的因数分解才能做到这一点。有一种更通用的方法，它是通过修改欧几里得算法（Euclidean algorthm）实现的。

算法

扩展欧几里得算法（Extended Euclidean algorithm），除了用来寻找最大公约数 $g = \gcd (a,b)$ ，也可以用来寻找整数 $x$ 和 $y$ 满足：

$a \cdot x + b \cdot y = g$

这可以通过将 $b$ 替换为 $m$ 和 $g$ 替换为 $1$ 来解决了寻找模运算逆元的问题：

$a^{−1}⋅a+k⋅m=1$

如果 $a$ 和 $m$ 不是互质的，那么就没有解，因为 $a$ 和 $m$ 的任何整数组合都不能得到一个不是他们最大公约数的倍数的数。

算法同样是递归的：它先计算出 $\gcd (b,a \mod b)$ 的系数 $x'$ 和 $y'$ ，然后恢复初始数字对的解。如果我们有一对 $(b,a \mod b)$ 的解 $(x',y')$ ：

$b⋅x′+(a \mod b)⋅y′=g$

那么，为了得到初始输入的解，我们可以将表达式 $(a \mod b)$ 重写为 $(a−⌊\frac{b}{a}⌋⋅b)$ 并将它代入前述方程:

$b⋅x′+(a−⌊\frac{b}{a}⌋⋅b)⋅y′=g$

现在我们按 $a$ 和 $b$ 进行分组，重新排列表达式可以得到：

$a \cdot \underbrace{y'}_x + b \cdot \underbrace{(x' - ⌊\frac{b}{a}⌋ \cdot y')}_y = g$

将上式与初始表达式相比较，可以得出 $a$ 和 $b$ 的系数即为初始的 $x$ 和 $y$ 。

实现

我们将算法实现为递归函数。由于它的输出不是一个而是三个整数，我们通过引用将系数传递给它：

int gcd(int a, int b, int &x, int &y) {
    if (a == 0) {
        x = 0;
        y = 1;
        return b;
    }
    int x1, y1;
    int d = gcd(b % a, a, x1, y1);
    x = y1 - (b / a) * x1;
    y = x1;
    return d;
}

为了计算模运算的乘法逆元，我们只需传递 $a$ 和 $m$ ，然后返回算法找到的 $x$ 系数。由于我们传递了两个正数，其中一个系数将是正的，另一个系数将是负的（具体哪一个取决于迭代次数是奇数还是偶数），所以我们需要检查 $x$ 是否为负，并加上 $m$ 以获得正确的余数。

int inverse(int a) {
    int x, y;
    gcd(a, M, x, y);
    if (x < 0)
        x += M;
    return x;
}

这段代码的执行时间约为160ns，比用二进制指数法求逆元快10ns。为了进一步优化它，我们同样可以用迭代方式实现，这样的执行时间为135ns：

int inverse(int a) {
    int b = M, x = 1, y = 0;
    while (a != 1) {
        y -= b / a * x;
        b %= a;
        swap(a, b);
        swap(x, y);
    }
    return x < 0 ? x + M : x;
}