算法学习笔记(2): 逆元及其应用

逆元

定义

逆元素，是指一个可以取消另一给定元素运算的元素

具体来说，对于实际的一些应用，如：

当我们想要求(11 / 3) % 10时

明显可以看出，是没有办法直接算的，这时就需要引入逆元

(a) 在模(p)意义下的逆元记作 (a^{-1})，也可以用inv(a)表示

应当满足

[a * a^{-1} equiv 1 pmod p ]

则此时，(11 / 3) % 10就可以写成(11 * inv(3)) % 10

可以求出，inv(3)在模10意义下= 7

[begin{align} 3 times inv(3) &= 21 \ 21 &equiv 1 pmod p end{align} ]

故(11 / 3) % 10 = (11 * 7) % 10 = ((11 % 10) * (7 % 10)) = (1 * 7) % 10 = 7

为什么我要多此一举在第三步再变换一次？

在实际应用中，a * b可能会很大以至于溢出，导致错误，所以分开来乘以减小数据规模

---

如何求？

费马小定理

依据费马小定理(需要注意先决条件，(a)与(p)互质且(p)是质数)

费马小定理可以通过欧拉定理求解，详见后文欧拉定理

[gcd(a, p) == 1 ; and ; text{p is prime} implies a^p equiv a pmod p ]

故此时可以有

[a^{-1} = a^{p-2} ]

扩展欧几里得算法

如果不满足先决条件呢？

这是相对来说的通发，但是总会有数据可以卡

根据观察

[a^{-1},a equiv 1 pmod p ]

令(i = a^{-1})换成等式可以知道

[ia + rp = 1 ]

由于已知(a, p)，则此时可以通过扩展欧几里得算法求解 (i) 的值

扩展欧几里得算法可以参考这篇文章：扩展欧几里得算法。

是我认为写的非常好的一篇文章。

---

欧拉定理

再推广一下？若 (p) 不为质数呢？

那么就要有欧拉定理来了

[gcd(k, p) == 1 implies k^{varphi(p)} equiv 1 pmod p ]

(varphi{(p)})指 ([1, p]) 中与(p)互质的数的个数。特别的，(1)也算。

举个例子：

• (varphi(7) = 6) ，因为7是质数（所以在(p)为质数的时候就退化成费马小定理了）

• (varphi(6) = 2)，因为只有1, 5和它互质

但是如何求(varphi(p))呢？

1. 将(p)分解质因数，于是有 (p = a_1^{c_1} , a_2^{c_2} , a_3 ^{c_3} ldots a_n^{c_n})

2. 此时(varphi(p) = p prodlimits_{i=1}^{n}frac {a_i -1}{a_i})

---

欧拉定理证明

令集合(A)为 ([1, p]) 中所有与(p)互质的数，即

[A_1 = {a_1, a_2, a_3, ldots, a_{varphi(p)}} ]

将(A)中每一个元素在模(p)意义下乘(k)，由于(A)中元素与(p)互质，且(k)也与(p)互质，可知

[A_2 = {ka_1 % p, ka_2 % p, ka_3 % p, ldots, ka_{varphi(p)} %p} ]

也满足为 ([1, p]) 中所有与p互质的数，故可知 (A_1 = A_2)

于是

[prodlimits_{i=1}^{varphi(p)} {a_i} equiv prodlimits_{i=1}^{varphi(p)} k{a_i}pmod p ]

即是

[prodlimits_{i=1}^{varphi(p)} {a_i} equiv k^{varphi(p)} prodlimits_{i=1}^{varphi(p)} {a_i}pmod p ]

左右相减，变形即可知 (k^{varphi(p)} equiv 1 pmod p)

扩展欧拉定理

[a^k equiv a^{k bmod varphi(p) + varphi(p)} pmod p ]

想必证明很简单，这里就不展开叙述了

---

补充：快速幂

可以看出，如果要利用欧拉定理，需要求(a^k)，当(k)非常大的时候，就需要快速幂的帮助了

推荐阅读：快速幂

这里给出一种参考代码

// (a**x) % p int quickPow(int a, int x, int p) {     int r = 1;     while (x) {         // no need to use quickMul when p*p can be smaller than int64.max !!!         if (x & 1) r = (r * a) % p;         a = (a * a) % p, x >>= 1;     }     return r; } 

至于其中的那一行注释，主要是考虑到当(a), (p)都很大（如：a = 1e15, p = 1e17 + 1时，a * a一定会溢出，所以需要“快速”乘来辅助）

实际上“快速”乘特别慢，是O(logn)的复杂度……所以叫龟速乘也不为过

推荐阅读：快速乘总结 - 一只不咕鸟，里面有更详细的阐述

这里给出快速乘的一种参考代码

// a*b % p O(log b) int quickMul(int a, int b, int p) {     // let b < a, to reduce a little time to process.     if (a < b) std::swap(a, b);      int r = 0;     while (b) {         if (b & 1) r = (r + a) % p;         a = (a<<1) % p, b >>= 1;     }     return r; } 

---

notice: 适当的使用long long

线性求逆元

不妨设我们需要求(i)在模(p)意义下的逆元

很容易知道，1的逆元为1，所以边界条件就有了

令 (p = k i + r), 放在模 (p) 意义下则有 (ki + r equiv 0 pmod p)

两边同时乘以 (i^{-1}r^{-1}) 可以得到 (kr^{-1} + i^{-1} equiv 0 pmod p)

变换一下

[begin{aligned} i^{-1} &equiv -kr ^{-1} pmod p \ i^{-1} &equiv -lfloor frac pi rfloor (p mod i)^{-1} pmod p \ inv(i) &equiv (p - lfloor frac pi rfloor)inv(p % i) pmod p end{aligned} ]

所以，有了递推式

inv[i] = (p - p/i) * inv[p % i] % p; 

线性求阶乘逆元

这个东西一般用于求组合数

我们先预处理出阶乘

fac[0] = 1; for (int i = 1; i <= n; ++i)     fac[i] = (fac[i - 1] * i) % p; 

根据逆元定义(i frac 1i equiv 1 pmod p)

所以 (inv(i!) equiv frac 1 {i!} pmod p)

稍微变换一下

[frac 1 {i!} equiv frac 1 {(i + 1)!}(i + 1) pmod p ]

所以有了递推式

ifac[i] = ifac[i + 1] * (i + 1) % p 

我们逆着推，假设最大需要到(n)

ifac[n] = quickPow(fac[n], p - 2); for (int i = n; i; i--)     ifac[i - 1] = ifac[i] * i % p; 

同时求逆元与阶乘逆元

还是逆元的本质是求倒数

[inv(i) equiv frac 1i pmod p ]

稍微变换一下

[inv(i) equiv frac 1 {i!} (i - 1)! equiv inv(i!) (i - 1)! pmod p ]

所以

inv[i] = ifac[i] * fac[i - 1] % p 

合起来就是

for (int i = n; i; i--) {     inv[i] = ifac[i] * fac[i - 1] % p;     ifac[i - 1] = ifac[i] * i % p; } 

就可以在较少的常数下同时求得两者了