乘法逆元的解法及证明

写在前面

刷题过程中，为了避免高精度的整数运算，常常会要求答案对某个数取模。但这衍生了一个新问题，如果运算过程中涉及了除法，如$\frac{a}{b}$，如果先对a和b取模再做除法，会导致答案错误。比如:
$$\frac{20}{5}mod10=4$$
但是，
$$\frac{20mod10}{5mod10}=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\ast b-1)modp=0$，即 $a\ast b\equiv 1(modp)$，则称，a 和 b 互为模 p 的乘法逆元。

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

$$\begin{array}{rl}xmodp& =x(a\ast b-1)+xmodp\\ & =(x\ast a\ast b-x)+xmodp\\ & =x\ast a\ast bmodp\end{array}$$

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

$$18mod7=4$$

另有，

$$\begin{array}{rl}18\ast 10\ast 12mod7& =2160mod7\\ & =4\end{array}$$

何时存在乘法逆元

结论

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

存在逆元 → 互质

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

$$(a\ast b-1)modp=0$$

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

$$(a\ast b-1)=np$$

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

$$\begin{array}{rl}(x\ast d\ast b-1)& =n\ast y\ast d\\ (x\ast d\ast b-n\ast y\ast d)& =1\\ (x\ast b-n\ast y)\ast d& =1\\ x\ast b-n\ast y& =\frac{1}{d}\end{array}$$

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

互质 → 存在逆元

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

在证明之前先来回忆一下辗转相除法：
$$\begin{array}{rl}gcd(a,p)& =gcd(p,amodp)\end{array}$$

不妨用r来表示模运算的结果：
$$\begin{array}{rl}gcd(a,p)& =gcd(p,r_1)\\ & =gcd(r_1,r_2)\\ & =...\\ & =gcd(r_{m-1},r_m)\\ & =1//因为a与p互质\end{array}$$

将上式用除法表示，则有：
$$\begin{array}{rl}a& =n_{1}p+r_1\\ p& =n_2{r}_1+r_2\\ r_1& =n_3{r}_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{array}$$

由上述的m+2个等式移项可得：
$$\begin{array}{rl}{r}_1& =a-n_{1}p\\ r_2& =p-n_2{r}_1\\ ...& =...\\ r_i& =r_{i-2}-n_i\ast r_{i-1}\\ r_{m-1}& =r_{m-3}-n_{m-1}{r}_{m-2}\\ 特别的，1& =r_{m-1}\end{array}$$

将 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{array}{rl}ax+py& =1\\ ax-1& =py\\ (ax-1)modp& =0\\ ax& \equiv 1(modp)\end{array}$$

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

如何计算乘法逆元

扩展欧几里得算法

回忆一下欧几里得的递归过程：
$$\begin{array}{rl}gcd(a,b)& =gcd(b,amodb)\\ & =gcd(b,a-b\lfloor \frac{a}{b}\rfloor )\end{array}$$

当 a 与 b 互质时，有 x，y，x'，y' 满足：
$$\{\begin{array}{r}ax+by=1\\ b{x}^{\prime }+(a-b\lfloor \frac{a}{b}\rfloor )y^{\prime }=1\end{array}$$

将上式合并移项：
$$\begin{array}{rl}b{x}^{\prime }+(a-b\lfloor \frac{a}{b}\rfloor )y^{\prime }& =ax+by\\ a(x-y^{\prime })+b(y-(x^{\prime }-\lfloor \frac{a}{b}\rfloor )y^{\prime })& =0\end{array}$$

令x，y，x'，y' 满足下述关系，则对于任意 a 和 b 上述公式都能成立。
$$\{\begin{array}{rl}x& =y^{\prime }\\ y& =(x^{\prime }-\lfloor \frac{a}{b}\rfloor )y^{\prime }\end{array}$$

😯好巧，这就是我们要寻找的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{array}{r}{a}^{p-1}\equiv 1(modp)\\ a\ast a^{p-2}\equiv 1(modp)\end{array}$$

于是，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{array}{rl}k& =\lfloor \frac{p}{i}\rfloor \\ j& =pmodi\end{array}$$
则有：
$$\{\begin{array}{rl}p& \equiv 0(modp)\\ (k\ast i+j)& \equiv 0(modp)\\ (k\ast i+j)\ast i^{-1}\ast j^{-1}& \equiv 0(modp)\\ k\ast j^{-1}+i^{-1}& \equiv 0(modp)\\ -k\ast j^{-1}& \equiv i^{-1}(modp)\\ -\lfloor \frac{p}{i}\rfloor \ast (pmodi{)}^{-1}& \equiv i^{-1}(modp)\end{array}$$

而且，在递推过程中，当我们要求$i^{-1}$时，对于任意的$(pmodi{)}^{-1}$都是已知的啦(因为 p mod i肯定小于 i)，所以可根据$(pmodi{)}^{-1}$ 通过常数次运算得到$i^{-1}$。
$$i^{-1}\equiv \{\begin{array}{rl}1& ,i=1\\ -\lfloor \frac{p}{i}\rfloor (pmodi{)}^{-1}& ,i>2\end{array}(modp)$$

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{array}{rl}\prod _{i=1}^n{a}_i\ast {a_i}^{-1}& \equiv 1modp\\ \prod _{i=1}^n{a}_i\ast \prod _{i=1}^n{a_i}^{-1}& \equiv 1modp\\ S_n\ast {S_n}^{-1}& \equiv 1modp\end{array}$$

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

又有
$$\begin{array}{rl}{S_i}^{-1}\ast {a_{i+1}}^{-1}\ast a_{i+1}& \equiv {S_i}^{-1}(modp)\\ {S_{i+1}}^{-1}\ast a_{i+1}& \equiv {S_i}^{-1}(modp)\end{array}$$

换言之，我们可以利用 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}\equiv S_{i-1}\ast {S_i}^{-1}(modp)$$

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;