LogoCSP Wiki By Yundou
0 MathBasics

乘法逆元

定义

如果一个线性同余方程 ax1(modb)ax \equiv 1 \pmod b,则 xx 称为 amodba \bmod b 的逆元,记作 a1a^{-1}

如何求逆元

扩展欧几里得法

void exgcd(int a, int b, int& x, int& y) {
  if (b == 0) {
    x = 1, y = 0;
    return;
  }
  exgcd(b, a % b, y, x);
  y -= a / b * x;
}

扩展欧几里得法和求解 线性同余方程 是一个原理,在这里不展开解释。

快速幂法

证明

因为 ax1(modb)ax \equiv 1 \pmod b

所以 axab1(modb)ax \equiv a^{b-1} \pmod b(根据 费马小定理);

所以 xab2(modb)x \equiv a^{b-2} \pmod b

然后我们就可以用快速幂来求了。

int qpow(long long a, int b) {
  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;
}

注意:快速幂法使用了费马小定理,要求 bb 是一个素数;而扩展欧几里得法只要求 gcd(a,b)=1\gcd(a, b) = 1

线性求逆元

求出 1,2,,n1,2,\dots,n 中每个数关于质数 pp 的逆元。

如果对于每个数进行单次求解,以上两种方法就显得慢了,很有可能超时,所以下面来讲一下如何线性(O(n)O(n))求逆元。

首先,很显然的 111(modp)1^{-1} \equiv 1 \pmod p

???+ note "证明" 对于 pZ\forall p \in \mathbf{Z},有 1×11(modp)1 \times 1 \equiv 1 \pmod p 恒成立,故在 pp11 的逆元是 11,而这是推算出其他情况的基础。

其次对于递归情况 i1i^{-1},我们令 k=pik = \lfloor \frac{p}{i} \rfloorj=pmodij = p \bmod i,有 p=ki+jp = ki + j。再放到 modp\mod p 意义下就会得到:ki+j0(modp)ki+j \equiv 0 \pmod p

pp 为质数时,可知 jj 在模 pp 意义下的乘法逆元必定存在。这时在上式 ki+j0(modp)ki+j \equiv 0 \pmod p 的两边同时乘 i1×j1i^{-1} \times j^{-1}

kj1+i10(modp)kj^{-1}+i^{-1} \equiv 0 \pmod p

i1kj1(modp)i^{-1} \equiv -kj^{-1} \pmod p

再带入 j=pmodij = p \bmod i,有 p=ki+jp = ki + j,有:

i1pi(pmodi)1(modp)i^{-1} \equiv -\lfloor\frac{p}{i}\rfloor (p \bmod i)^{-1} \pmod p

我们注意到 pmodi<ip \bmod i < i,而在迭代中我们完全可以假设我们已经知道了所有的模 pp 下的逆元 j1,j<ij^{-1}, j < i

故我们就可以推出逆元,利用递归的形式,而使用迭代实现:

i1{1,if i=1,pi(pmodi)1,otherwise.(modp)i^{-1} \equiv \begin{cases} 1, & \text{if } i = 1, \\ -\lfloor\frac{p}{i}\rfloor (p \bmod i)^{-1}, & \text{otherwise}. \end{cases} \pmod p
inv[1] = 1;
for (int i = 2; i <= n; ++i) {
  inv[i] = (long long)(p - p / i) * inv[p % i] % p;
}

使用 ppip-\lfloor \dfrac{p}{i} \rfloor 来防止出现负数。

当 p 不为质数

线性同余方程 中指出,如果 iipp 不互素时不存在相应的逆元。如果 pp 不为质数,则至少有一个 ii 的逆元不存在。此时这个建立在递推式上的方法就不能保证结果的正确性。例如当 p=8p = 8i=3i = 3 时,根据递推式,逆元需要从 inv[p % i]inv[2] 计算。而 2 在模 8 意义下不存在逆元,inv[2] 的值是未定义的(准确来讲,递推时 inv[2] 依赖于 inv[0] 的初始值)。因此 3 在模 8 意义下的乘法逆元便无法正确求出。

另外,根据线性求逆元方法的式子:i1kj1(modp)i^{-1} \equiv -kj^{-1} \pmod p

递归求解 j1j^{-1}, 直到 j=1j=1 返回 11

中间优化可以加入一个记忆化来避免多次递归导致的重复,这样求 1,2,,n1,2,\dots,n 中所有数的逆元的时间复杂度仍是 O(n)O(n)

注意:如果用以上给出的式子递归进行单个数的逆元求解,目前已知的时间复杂度的上界为 O(n13)O(n^{\frac 1 3}),具体请看 知乎讨论。算法竞赛中更好地求单个数的逆元的方法有扩展欧几里得法和快速幂法。

线性求任意 n 个数的逆元

上面的方法只能求 11nn 的逆元,如果需要求任意给定 nn 个数(1ai<p1 \le a_i < p)的逆元,就需要下面的方法:

首先计算 nn 个数的前缀积,记为 sis_i,然后使用快速幂或扩展欧几里得法计算 sns_n 的逆元,记为 svnsv_n

因为 svnsv_nnn 个数的积的逆元,所以当我们把它乘上 ana_n 时,就会和 ana_n 的逆元抵消,于是就得到了 a1a_1an1a_{n-1} 的积逆元,记为 svn1sv_{n-1}

同理我们可以依次计算出所有的 svisv_i,于是 ai1a_i^{-1} 就可以用 si1×svis_{i-1} \times sv_i 求得。

所以我们就在 O(n+logp)O(n + \log p) 的时间内计算出了 nn 个数的逆元。

s[0] = 1;
for (int i = 1; i <= n; ++i) s[i] = s[i - 1] * a[i] % p;
sv[n] = qpow(s[n], p - 2);
// 当然这里也可以用 exgcd 来求逆元,视个人喜好而定.
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;

例题:

On this page