本文介绍模意义下乘法运算的逆元,并讨论它的常见求解方法。

基本概念

非零实数 的乘法逆元就是它的倒数 。类似地,数论中也可以定义一个整数 在模 意义下的逆元 ,或简单地记作 。这就是 模逆元(modular multiplicative inverse),也称作 数论倒数

逆元

对于非零整数 ,如果存在 使得 ,就称 在模 意义下的 逆元(inverse)。

这相当于说, 是线性同余方程 的解。根据 线性同余方程 的性质可知,当且仅当 ,即 互素时,逆元 存在,且在模 的意义下是唯一的。

单个逆元的求法

利用扩展欧几里得算法或快速幂法,可以在 时间内求出单个整数的逆元。

扩展欧几里得算法

求解逆元,就相当于求解线性同余方程。因此,可以使用 扩展欧几里得算法 时间内求解逆元。同时,由于逆元对应的线性方程比较特殊,可以适当地简化相应的步骤。

参考实现

[list2tab]

  • C++

    // Extended Euclidean algorithm.
    void ex_gcd(int a, int b, int& x, int& y) {
      if (!b) {
        x = 1;
        y = 0;
      } else {
        ex_gcd(b, a % b, y, x);
        y -= a / b * x;
      }
    }
     
    // Returns the modular inverse of a modulo m.
    // Assumes that gcd(a, m) = 1, so the inverse exists.
    int inverse(int a, int m) {
      int x, y;
      ex_gcd(a, m, x, y);
      return (x % m + m) % m;
    }
  • Python

    # Extended Euclidean algorithm.
    def ex_gcd(a, b):
        if b == 0:
            return 1, 0
        else:
            x1, y1 = ex_gcd(b, a % b)
            x = y1
            y = x1 - (a // b) * y1
            return x, y
     
     
    # Returns the modular inverse of a modulo m.
    # Assumes that gcd(a, m) = 1, so the inverse exists.
    def inverse(a, m):
        x, y = ex_gcd(a, m)
        return (x % m + m) % m

这一算法适用于所有逆元存在的情形。

快速幂法

这一方法主要适用于模数是素数 的情形。此时,由 费马小定理 可知对于任意 都有

根据逆元的唯一性可知,逆元 就等于 ,因此可以直接使用 快速幂 时间内计算:

参考实现

[list2tab]

  • C++

    // Binary exponentiation.
    int pow(int a, int b, int m) {
      long long res = 1, po = a;
      for (; b; b >>= 1) {
        if (b & 1) res = res * po % m;
        po = po * po % m;
      }
      return res;
    }
     
    // Returns the modular inverse of a prime modulo p.
    int inverse(int a, int p) { return pow(a, p - 2, p); }
  • Python

    # Returns the modular inverse of a prime modulo p.
    # Use built-in pow function.
    def inverse(a, p):
        return pow(a, p - 2, p)

当然,理论上,这一方法可以利用 欧拉定理 推广到一般的模数 的情形,即利用 计算逆元。但是,单次求解 欧拉函数 并不容易,因此该算法在一般情况下效率不高。

多个逆元的求法

有些场景下,需要快速处理出多个整数 在模 意义下的逆元。此时,逐个求解逆元,总共需要 的时间。实际上,如果将它们统一处理,就可以在 的时间内求出所有整数的逆元。

考虑序列 的前缀积:

只要每个 都与 互素,它们的乘积 就与 互素。因此,可以通过前文所述算法求出 的值。因为乘积的逆元就是逆元的乘积,所以,从 出发,反向遍历序列就能求出每个 的逆元:

由此,单个 的逆元可以通过下式计算:

参考实现如下:

参考实现

[list2tab]

  • C++

    // Returns the modular inverses for each x in a modulo m.
    // Assume x mod m exists for each x in a.
    std::vector<int> batch_inverse(const std::vector<int>& a, int m) {
      int n = a.size();
      std::vector<int> prod(n);
      long long s = 1;
      for (int i = 0; i < n; ++i) {
        // prod[i] = product of a[0...i-1]; prod[0] = 1.
        prod[i] = s;
        s = s * a[i] % m;
      }
      // s = product of all elements in a.
      s = inverse(s, m);
      std::vector<int> res(n);
      for (int i = n - 1; i >= 0; --i) {
        res[i] = s * prod[i] % m;
        s = s * a[i] % m;
      }
      return res;
    }
  • Python

    # Returns the modular inverses for each x in a modulo m.
    # Assume x mod m exists for each x in a.
    def batch_inverse(a, m):
        n = len(a)
        prod = [0] * n
        s = 1
        for i in range(n):
            # prod[i] = product of a[0...i-1]; prod[0] = 1.
            prod[i] = s
            s = s * a[i] % m
        # s = product of all elements in a.
        s = inverse(s, m)
        res = [0] * n
        for i in reversed(range(n)):
            res[i] = s * prod[i] % m
            s = s * a[i] % m
        return res

算法中,只求了一次单个元素的逆元,因此总的时间复杂度是 的。

线性时间预处理逆元

如果要预处理前 个正整数在素数模 下的逆元,还可以通过本节将要讨论的递推关系在 时间内计算。这一方法常用于组合数计算中前 个正整数的阶乘的倒数的预处理。

对于 的正整数 ,考察带余除法:

将该等式对素数 取模,就得到

将等式两边同时乘以 就得到

这就是用于线性时间递推求逆元的公式。由于 ,这一公式将求解 的问题转化为规模更小的问题 。因此,从 开始,对每个 顺次应用该公式,就可以在 时间内获得前 个整数的逆元。

参考实现如下:

参考实现

[list2tab]

  • C++

    // Precomputes modular inverses of all integers from 1 to n modulo prime p.
    std::vector<int> precompute_inverses(int n, int p) {
      std::vector<int> res(n + 1);
      res[1] = 1;
      for (int i = 2; i <= n; ++i) {
        res[i] = (long long)(p - p / i) * res[p % i] % p;
      }
      return res;
    }
  • Python

    # Precomputes modular inverses of all integers from 1 to n modulo prime p.
    def precompute_inverses(n, p):
        res = [0] * (n + 1)
        res[1] = 1
        for i in range(2, n + 1):
            res[i] = (p - p // i) * res[p % i] % p
        return res

这一算法只适用于模数是素数的情形。对于模数 不是素数的情形,无法保证递推公式中得到的 仍然与 互素,因而递推所需要的 可能并不存在。一个这样的例子是 。此时,,不存在模 的逆元。

另外,得到该递推公式后,一种自然的想法是直接递归求解任意一个数 的逆元。每次递归时,都利用递推公式将它转化为更小的余数 的逆元,直到余数变为 时停止。目前尚不清楚这样做的复杂度 1,因此,推荐使用前文所述的常规方法求解。

习题

参考资料与注释

Footnotes

  1. riteme 在知乎上的回答 中指出,这样做理论上已知的复杂度的上界是 ,而在实际随机数据中的表现接近于