写在前面

刷题过程中，为了避免高精度的整数运算，常常会要求答案对某个数取模。但这衍生了一个新问题，如果运算过程中涉及了除法，如 $\frac{a}{b}$ ，如果先对a和b取模再做除法，会导致答案错误。比如: $\frac{20}{5}\ mod\ 10 = 4$ 但是， $\frac{20\ mod\ 10}{5\ mod\ 10} = 0$

为了解决这个问题，我们需要学习新知识：模运算下的乘法逆元。

乘法逆元

什么是逆元呢？百度百科告诉我们：逆元是指一个可以取消另一给定元素运算的元素。

比如，在实数的加法中，任意实数 a 和它的相反数 -a，互为加法逆元。因为对于任意实数 x，加上a之后，再加上-a，仍然等于 x，换言之，-a取消了x和a的加法运算。

再比如，在实数的乘法中，任意不为0的实数b和它的倒数 $\frac{1}{b}$ ，互为乘法逆元。因为对于任意实数x，乘上b之后，在乘上 $\frac{1}{b}$ ，仍然等于x，换言之， $\frac{1}{b}$ 取消了x 和 b 的乘法运算。

接下来看看模运算下的乘法逆元定义：如果两个整数a 和 b 满足 $(a*b-1)\ mod\ p = 0$ ，即 $a*b≡1(mod\ p)$ ，则称，a 和 b 互为模 p 的乘法逆元。

换言之，如果整数a, b, p 满足 $a*b≡1(mod\ p)$ ，那么 b 可以取消 a 和任意整数x在模p时的乘法，即 a 和 b 互为模 p 的乘法逆元。即 $x*a*b≡x (mod\ p)$ ，证明过程如下：

\begin{aligned} x\ mod\ p &= x(a*b-1)+x \ mod\ p \\ &= (x*a*b-x)+x \ mod \ p \\ &= x*a*b \ mod \ p \end{aligned}

举个实际的例子：比如 10 和 12 互为模7时的乘法逆元。然后再找个整数，比如 18，此时有：

18\ mod\ 7 = 4 \\

另有，

\begin{aligned} 18 * 10 * 12\ mod\ 7 &= 2160\ mod\ 7 \\ &= 4 \end{aligned}

何时存在乘法逆元

结论

先说结论：a 与 p 互质 是 存在a关于模p的乘法逆元b 的充分必要条件。

存在逆元 → 互质

证明 a 与 p 互质 是 存在乘法逆元b 的必要条件。不妨先假设 a 与 b 互为模p时的乘法逆元，且a与p的最大公约数为d，则有：

$(a*b-1)\ mod\ p = 0$

显然，a*b-1 为 p 的 n 倍，n 为整数，则上式转化为：

$(a*b-1)=np$

将等式中的a与p用d代替，则等式转化为：

\begin{aligned} (x*d*b-1)&=n*y*d \\ (x*d*b - n*y*d) &= 1 \\ (x*b - n*y)*d &= 1 \\ x*b - n*y &= \frac{1}{d} \\ \end{aligned}

显然，x*b-n*y是一个整数，所以仅当d为1时，即a与p互质时，上述等式才可能成立。

互质 → 存在逆元

证明 a 与 p 互质 是 存在乘法逆元b 的充分条件。

在证明之前先来回忆一下辗转相除法：

\begin{aligned} gcd(a,p) &= gcd(p, a\ mod\ p) \\ \end{aligned}

不妨用r来表示模运算的结果：

\begin{aligned} gcd(a,p) &= gcd(p, r_1) \\ &= gcd(r_1,r_2) \\ &= ... \\ &= gcd(r_{m-1}, r_m) \\ &= 1\ \ //因为a与p互质 \end{aligned}

将上式用除法表示，则有：

\begin{aligned} a &= n_1p + r_1\\ p &= n_2r_1 + r_2 \\ r_1 & = n_3r_2 + r_3 \\ ... &= ... \\ r_{m-3} &= n_{m-1}r_{m-2} + r_{m-1} \\ r_{m-2} &= n_{m}r_{m-1} \\ r_{m-1} &= 1\\ r_{m} &= 0 \\ \end{aligned}

由上述的m+2个等式移项可得：

\begin{aligned} r_1 &= a - n_1p \\ r_2 &= p - n_2r_1 \\ ... &= ... \\ r_{i} &= r_{i-2} - n_{i}*r_{i-1} \\ r_{m-1} &= r_{m-3} - n_{m-1}r_{m-2} \\ 特别的，1 &= r_{m-1} \\ \end{aligned}

将 r_m-1，r_m-2 ... 以及 r₁ 依次代入 $1 = r_{m-1}$

可以得到用 a，p 以及 n_i表示的等式： $ax+py=1$

其中，x 与 y 为 n₁，n₂，n₃ ... 的加减乘的运算结果，显然为整数。

得证，当 a 与 p 互质时，存在整数x与y 满足

\begin{aligned} ax+py &= 1 \\ ax-1 &= py \\ (ax-1)\ mod\ p &= 0 \\ ax &≡ 1\ (mod\ p) \end{aligned}

即，当 a 与 p 互质时，存在整数x 与 a互为模p时的乘法逆元。

如何计算乘法逆元

扩展欧几里得算法

回忆一下欧几里得的递归过程：

\begin{aligned} gcd(a,b) &= gcd(b, a\ mod\ b) \\ &= gcd(b, a - b\lfloor \frac{a}{b} \rfloor) \end{aligned}

当 a 与 b 互质时，有 x，y，x'，y' 满足：

\begin{aligned} ax+by=1 \\ bx'+ (a - b\lfloor \frac{a}{b} \rfloor)y' = 1\\ \end{aligned} \right.

将上式合并移项：

\begin{aligned} bx'+ (a - b\lfloor \frac{a}{b} \rfloor)y'&=ax+by \\ a(x-y') + b(y-(x'-\lfloor \frac{a}{b} \rfloor)y')&=0 \end{aligned}

令x，y，x'，y' 满足下述关系，则对于任意 a 和 b 上述公式都能成立。

\begin{aligned} x &= y' \\ y &= (x'-\lfloor \frac{a}{b} \rfloor)y' \\ \end{aligned} \right.

😯好巧，这就是我们要寻找的x，y和x'，y'的计算公式。另外，当 b = 0，即到达欧几里得算法的递归边界时，令 x = 1, y = 0 即可满足ax+by=1(此时的a显然为1)。

于是我们得到了下述代码，最终得到的 x 即为 a 在模 b 时的乘法逆元，当然前提是函数的返回值为1，即 a 与 b 互质。

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 t = y;
    y = x-(a/b)*y;
    x = t;
    return r;
}

费马小定理 & 快速幂

费马小定理：若 p 为素数，且 gcd(a, p) = 1，则满足

\begin{aligned} a^{p-1} ≡ 1\ (mod\ p) \\ a*a^{p-2} ≡ 1\ (mod\ p) \end{aligned}

于是，a 在模 p 时的乘法逆元为 $a^{p-2}$ ，可用快速幂求得。

inline int quick_power(int a, int b, int p) {
  int ans = 1;
  a = (a % p + p) % p;
  for (; b; b >>= 1) {
    if (b & 1) ans = (a * ans) % p;
    a = (a * a) % p;
  }
  return ans;
}
// quick_power(a, p-2, p) 即为 a 的乘法逆元。

值得注意的是，该算法仅在 p 为素数时有效，而扩欧并没有该限制。

线性求1到n中每个数的逆元

现在要求1到n中每个数i在模p时的乘法逆元。

特别的，i = 1时，其乘法逆元i^-1为1。

当 i > 1 时，不妨先设：

\begin{aligned} k &= \lfloor\frac{p}{i}\rfloor \\ j &= p\ mod\ i \\ \end{aligned} \right.

则有：

\begin{aligned} p &≡ 0\ (mod\ p) \\ (k*i +j) &≡ 0\ (mod\ p) \\ (k*i +j)*i^{-1}*j^{-1} &≡ 0\ (mod\ p) \\ k*j^{-1}+i^{-1} &≡ 0\ (mod\ p) \\ -k*j^{-1} &≡i^{-1}\ (mod\ p) \\ -\lfloor\frac{p}{i}\rfloor *(p\ mod\ i)^{-1} &≡i^{-1} \ (mod\ p) \\ \end{aligned} \right.

而且，在递推过程中，当我们要求 $i^{-1}$ 时，对于任意的 $(p\ mod\ i)^{-1}$ 都是已知的啦(因为 p mod i肯定小于 i)，所以可根据 $(p\ mod\ i)^{-1}$ 通过常数次运算得到 $i^{-1}$ 。

\begin{aligned} 1&, i = 1 \\ -\lfloor\frac{p}{i}\rfloor(p\ mod\ i)^{-1}&, i > 2 \\ \end{aligned} \right. \ (mod\ p)

inv[1] = 1;
for (int i = 2; i <= n; ++i) {
  inv[i] = (p-p/i)*inv[p%i]%p;
}

求任意n个整数的乘法逆元

设有 n 个整数，分别为 a₁，a₂，a₃...，对应的逆元为a₁^-1，a₂^-1，a₃^-1...。

设S_i 为 a₁ 到 a_i i 个数的乘积，S_i^-1 为 a₁^-1 到 a_i^-1 i 个数的乘积。那么 S_i 与 S_i^-1 互为逆元。证明过程如下：

\begin{aligned} \prod_{i = 1}^{n}{a_i*{a_i}^{-1}} &≡ 1\ mod\ p \\ \prod_{i = 1}^{n}a_i*\prod_{i = 1}^{n}{a_i}^{-1} &≡ 1\ mod\ p \\ S_n*{S_n}^{-1} &≡ 1\ mod\ p \\ \end{aligned}

接下来，我们可以先求出 S₁ 到 S_n，然后通过扩展欧几里得或者快速幂求出 S_n 的乘法逆元 S_n^-1。

又有

\begin{aligned} {S_i}^{-1}*{a_{i+1}}^{-1}*a_{i+1} &≡ {S_i}^{-1}\ (mod\ p) \\ {S_{i+1}}^{-1}*a_{i+1} &≡ {S_i}^{-1}\ (mod\ p) \end{aligned}

换言之，我们可以利用 a_i可以取消a_i^-1乘法运算的性质，由S_i^-1 计算出 S_i-1^-1。

这样就可以O(n)的计算出所有的 S_i^-1。

最后，由S_i^-1以及S_i 计算出所有的 a_i^-1：

{a_i}^{-1} ≡ S_{i-1}*{S_{i}}^{-1}\ (mod\ p)

s[0] = 1;
for (int i = 1; i <= n; ++i) s[i] = s[i - 1] * a[i] % p;
sv[n] = quick_power(s[n], p - 2);
for (int i = n; i >= 1; --i) sv[i - 1] = sv[i] * a[i] % p;
for (int i = 1; i <= n; ++i) inv[i] = sv[i] * s[i - 1] % p;

四种方法求乘法逆元！四种！