演算法學習筆記(2): 逆元及其應用

jeefies發表於2023-01-13

演算法筆記

逆元

定義

逆元素，是指一個可以取消另一給定元素運算的元素

具體來說，對於實際的一些應用，如：

當我們想要求(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{aligned} 3 \times inv(3) &= 21 \\ 21 &\equiv 1 \pmod p \end{aligned} \]

故(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)\)呢？

將\(p\)分解質因數，於是有 \(p = a_1^{c_1} \, a_2^{c_2} \, a_3 ^{c_3} \ldots a_n^{c_n}\)
此時\(\varphi(p) = p \prod\limits_{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\)

於是

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

即是

\[\prod\limits_{i=1}^{\varphi(p)} {a_i} \equiv k^{\varphi(p)} \prod\limits_{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 \]