Chapter III. 漸入佳境
基礎數學
基礎數論
作者WiwiHo
先備知識常用數學演算法

本文包含許多比較細瑣的證明,實際上不會那些證明也不太會怎麼樣,有興趣的讀者再參考即可,對證明沒有興趣的讀者可以不用點開預設是收合的證明框框。

數論(number theory)是數學中研究整數性質的領域,整數作為一種在電腦中容易儲存和運算的東西,在整個資訊領域有很廣泛的應用,在程式競賽中也很重要。

常用數學演算法中,我們提過的模運算與質因數都屬於數論的範圍。之前我們提到了「同餘」是一個類似於等於的概念,兩個同餘的整數即便實際上不一樣,我們在模運算之下也可以將它們視為相同的東西,拿去做加、減、乘運算,不過我們也特別說了除法並沒有這種好性質,像是 $6$ 和 $12$ 在模 $6$ 下同餘,但除以 $3$ 之後就變成 $2$ 和 $4$,它們 $\bmod 6$ 的結果是不一樣的。這件事情其實不難理解,畢竟我們熟知的除法,是會把兩個整數的運算結果變成非整數的,而我們的模運算是整數之間的運算,要求除法之後還要保持同餘性質當然不太合理。

不過,同餘之下可以加也可以減,但可以乘卻不能除,那感覺就缺少了一塊。本文會以「想辦法在同餘之下做除法」為主軸,同時讓讀者認識程式競賽中常用的各種數論知識。

費馬小定理


要做到除法,我們應該先來嘗試暸解一下除法的好朋友:乘法。模 $m$ 下的運算其實是在「讓世界上的數字只剩下 $0,1,\dots,m-1$」,所以要是我們不斷做同一種運算,總有一天肯定會在運算過程中出現重複的數字,舉一個很顯然的例子,一開始我們讓 $x=0$,不斷把這個數字 $+1$,加了 $m$ 次後就會回到 $0$ 了,接下來就會不斷循環。

把 $+1$ 換成加一個固定的數字 $a$ 也差不多,總有一天一定會回到 $0$,因為 $m \times a \equiv 0 \pmod{m}$ 嘛,所以最慢就是 $m$ 步以後一定會回到 $0$(搞不好會比 $m$ 步更早,像是每次 $+5$、$m=10$,只要 2 步就回到 $0$ 了)。

那麼乘法呢?既然我們是要乘法,應該要從 $x=1$ 開始,每次乘一個固定的數字 $a$,看看什麼時候會回到 $1$。舉些例子:

跟加法不一樣的是,有時候會陷入永遠回不到 $1$ 的迴圈之中,而且 $a^m$ 好像可能會是 $1$ 也有可能不是……,看似沒什麼規律,我們先來看比較特殊的例子:

Theorem 1
費馬小定理

對於質數 $p$ 和不是它的倍數的 $a$,
\[ a^{p-1} \equiv 1 \pmod{p} \]

像我們剛剛 $a=3,\ m=5$ 的例子,我們發現了 $a^4 \equiv 1 \pmod 5$,而 $4$ 正好就是 $5-1$,$a=4,\ m=5$ 同樣也是。

$p-1$ 這個數字其實很有道理,畢竟一直乘一個不是 $0$(不同餘)的數字,總不可能跑出 $0$ 嘛,所以可以跑出來的數字頂多就 $p-1$ 種,只是不一定每種都會出現就是了,但大致上可以這樣理解。因為質數有這樣的特殊性質,所以很多題目都會使用質數作為模數,像 $10^9+7$ 和 $998244353$ 都是常見的質數。至於這個性質實際上有什麼用,馬上就會知道了。

費馬小定理有一個好用的應用,$a^{p-1} \equiv 1 \pmod p$ 的左右可以同時乘 $a$,所以 $a^n \equiv a^{n \bmod (p-1)} \pmod{p}$,如果想要算的 $a^n \bmod p$ 的 $n$ 超級大,那就可以先把 $n$ 模 $p-1$ 再用快速冪。

模運算下的除法


除法在很多問題,例如要數數量時,是難以避免的存在。舉個例子,像是我們只是想要算 $P^n_k=\frac{n!}{(n-k)!} \bmod 11$,$n,k \leq C$,有很多筆詢問,如果我們可以輕鬆在 $O(1)$ 時間進行大數運算,那只要先預處理好 $1!,2!,\dots,C!$,就直接算完除法後 $\bmod$,就能在 $O(C)$ 的預處理之後,用 $O(1)$ 的時間回答每一個詢問了。然而世界沒這麼美好,我們要大數運算的話既不能 $O(1)$ 也不會輕鬆,可是如果算階乘的時候先 $\bmod$,又不能直接用除法,不然像是 $5! \bmod 11 = 10$、$3! \bmod 11 = 6$,$\frac{5!}{3!} \bmod 11 = 9$,但 $10/6$ 不知道是什麼東西。當然也有個辦法是硬算 $(n-k+1) \times (n-k+2) \times \dots \times n \bmod 11$,就可以每乘完一次就模一次,但這樣一個詢問就得花上 $O(C)$ 的時間。

要怎麼做到「在同餘下的除法」呢?從剛剛的 $P^n_k$ 當出發點,我們最初步的基本要求是,有兩個整數 $a,b$,我們知道 $\frac{a}{b}$ 是整數,那我們要有辦法在只知道 $a \bmod m$ 和 $b \bmod m$ 的情況下,算出 $\frac{a}{b} \bmod m$。

我們要先回來探究一下除法的本質。在我們熟悉的一般算數世界,除法是一種「乘法的逆向操作」,一個數字 $\times a$ 再 $\div a$ 之後,就會變回本來的數字,甚至除法根本就是一種乘法,$\div a$ 的意思就是 $\times \frac{1}{a}$,也就是「乘以 $a$」和「乘以 $a$ 的倒數」互為逆向操作。

回到模運算底下的世界,既然乘法會保持同餘性質,那我們可不可以幫數字找到一個「類似於倒數的東西」,乘上它就類似於除法呢?在實數的世界裡有 $a \times \frac{1}{a} = 1$,感覺倒數就應該要是跟本來的數乘起來等於 $1$,所以我們在同餘世界如法炮製看看:一個整數 $a$ 在模 $m$ 下的「倒數」是 $b$,必須要滿足 $ab \bmod m = 1$。

舉例來說,$a=6,\ m=11$,因為我們的模數 $11$ 是個質數,所以我們馬上就可以用費馬小定理找到答案!我們知道 $6^{10} \equiv 1 \pmod{11}$,所以 $6 \times 6^{9} \bmod 11 = 1$,$6^9 \bmod 11$ 就是 $2$。要是我們直接把 $2$ 當成 $6$ 的「倒數」,用 $5! \equiv 10 \pmod{11}$ 和 $\frac{1}{3!} \stackrel{好像}{\equiv} 2 \pmod{11}$ 去算 $\frac{5!}{3!}$ 的話,答案就是 $20 \equiv 9 \pmod{11}$,算出正確答案了!

這個「$a$ 的倒數」$b$ 我們把它稱作模逆元模反元素(modular multiplicative inverse,或簡稱為 modular inverse),直接記作 $a^{-1}$,不過注意這裡的意思和實數運算之下的 $a^{-1}$ 是不一樣的。

費馬小定理告訴了我們,當模數是個質數 $p$ 時,$0 < a < p$ 的模逆元 $a^{-1}$ 就是 $a^{p-2} \bmod p$,只要用快速冪就可以在 $O(\log p)$ 的時間內找到一個數字的模逆元了。

對於某一種二元運算,如果跟 $a$ 做運算,再跟 $b$ 做運算之後,等同於沒有做運算,那我們就說 $b$ 是 $a$ 在這種運算底下的「反元素」(inverse)。舉例來說,$5$ 和 $-5$ 互為實數加法底下的反元素,$4$ 和 $0.25$ 互為實數乘法底下的反元素,剛剛舉例的 $6$ 和 $2$ 則互為「模 $11$ 下的乘法的反元素」,所以說得精確一點其實 $2$ 是 $6$ 的「模 $11$ 下的乘法反元素」,只是通常提到模運算底下的反元素就是指乘法,所以經常直接簡稱叫模反元素,「逆元」是 inverse 的另一個翻譯。

歐拉函數與歐拉定理


我們學會了怎麼在模數是質數的時候找到模逆元,那不是質數的時候怎麼辦?我們先來認識一個函數:

Definition 1
歐拉函數

歐拉函數(Euler's totient function)$\varphi(n)$ 定義為 $[1,n-1]$ 之中和 $n$ 互質的整數數量。

歐拉函數直接有個公式:

\[ \varphi(n)=n \prod_{p \mid n} \frac{p-1}{p} \]

其中 $p$ 是 $n$ 的質因數。這個公式看起來很神奇,但我們先暫時接受它,至於證明我們以後再說。

舉例來說,$\phi(20)=20 \times \frac{1}{2} \times \frac{4}{5} = 8$,跟 $20$ 互質的數有 $1,3,7,9,11,13,17,19$,的確是 8 個、對任何質數 $p$,所有 $1$ 到 $p-1$ 的整數都跟它互質,所以 $\varphi(p)=p \times \frac{p-1}{p}=p-1$。

歐拉函數是數論中很重要的函數,有很多很酷的性質,其中一個性質便是:

Theorem 2
歐拉定理

如果 $a,m$ 互質,則

\[ a^{\varphi(m)} \equiv 1 \pmod{m} \]

我們前面舉過 $a=5,\ m=6$ 的例子,$\varphi(6)=2$,所以 $5^2 \equiv 1 \pmod 6$。

費馬小定理其實是歐拉定理的特別版,在 $m$ 是質數的時候,$\varphi(m)=m-1$,所以 $a^{m-1} \equiv 1 \pmod m$,就跟費馬小定理說的一樣。相似的,歐拉定理告訴了我們,對於和 $m$ 互質的 $a$,$a^{-1} \equiv a^{\varphi(m)-1} \pmod{m}$,像在模 $6$ 下,$5$ 的模逆元就是 $5$ 自己。

什麼時候會有模逆元


我們學會在 $a,m$ 互質的時候找 $a$ 模 $m$ 下的模逆元了,那其他狀況怎麼辦?先直接說結論:沒有別的狀況了,$a,m$ 不互質就是沒有模逆元,像剛剛提過的 $a=4,m=6$ 的例子,枚舉一下 $1,2,3,4,5$ 就會發現,沒有一個乘上 $4$ 會同餘 $1$ 的。

要知道為什麼,我們就要來研究一下模逆元到底是個什麼東西了,前面我們都是直接湊出一個滿足 $ab \equiv 1 \pmod{m}$ 的 $b$ 當作 $a^{-1}$,但滿足這個條件的 $b$ 到底會長怎麼樣?寫出同餘的定義,我們要的是 $ab - 1$ 是 $m$ 的倍數,也就是說,存在一個整數 $k$,滿足 $ab-1=mk$,移項一下後就變成 $ab-mk=1$,只要這個式子存在 $b,k$ 都是整數的解,那這個 $b$ 就可以是模逆元。

什麼條件下會有整數的解呢?

Lemma 1
貝祖定理(Bézout's identity)的其中一個方向

對於整數 $a,b,c$,如果 $g=\gcd(a,b)$ 不是 $c$ 的因數,那麼 $ax+by=c$ 就沒有整數解。

證明:

所以 $ab-mk=1$ 要有整數解,就得要滿足 $\gcd(a,m)$ 是 $1$ 的因數,也就是 $a,m$ 互質,因此不互質就沒救了。

雖然我們前面學到的方法已經可以幫我們在任何有模逆元的狀況下找出模逆元,但是如果 $m$ 不是質數的話,我們就得用歐拉定理,算歐拉函數得要質因數分解,要是我們需要算很多不同超大模數之下的模逆元,那肯定不夠快。有沒有什麼快速的方法能幫我們找到模逆元呢?

擴展歐幾里得演算法


剛剛我們遇到的問題,講得一般化一點,就是:

例題
二元一次方程式求整數解
Source:經典題

給三個整數 $a,b,c$,求 $ax+by=c$ 的任意一組 $x,y$ 都是整數的解,或判斷沒有這種解。

貝祖定理告訴我們,要是 $g=\gcd(a,b)$ 不是 $c$ 的因數,那就沒有整數解,那是不是 $g$ 是 $c$ 的因數時就肯定有解?答案是對,至於要怎麼做,就是擴展歐幾里得演算法(extended Euclidean algorithm)要解決的問題了。擴展歐幾里得演算法要解決的問題是如下的形式:對於整數 $a,b$,求出 $ax+by=\gcd(a,b)$ 的整數解,在本來的問題中,$c$ 是 $g=\gcd(a,b)$ 的倍數時,就只要先找到 $ax+by=g$ 的整數解,然後把 $x,y$ 都乘上 $c/g$,就可以得到原問題的整數解了,因此我們只要會解 $ax+by=g$ 就好。

既然它叫作擴展「歐幾里得演算法」,問題中又出現了 $\gcd$,那就要想到我們在基礎數學 / 常用數學演算法學過的歐幾里得演算法(輾轉相除法)了!輾轉相除法有一種很好寫的實作方式是

int gcd(int a, int b) {
    if (b == 0) return a;
    return gcd(b, a % b);
}

擴展歐幾里得演算法就是在嘗試用 $\gcd(a,b)=\gcd(b, a \bmod b)$ 這個遞迴式,來找出 $ax+by=g$ 的解。既然要遞迴,那我們就嘗試把「$ax+by=\gcd(a,b)$」裡的 $a$ 直接換成 $b$、$b$ 直接換成 $a \bmod b$,也就是說我們想辦法用
\[ bx' + (a \bmod b)y' = \gcd(b, a \bmod b) \]
的整數解來湊出 $ax+by=\gcd(a, b)$ 的解。注意 $\gcd(b, a \bmod b)$ 和 $\gcd(a,b)$ 是一樣的,所以要是我們知道一個整數 $(x',y')$ 解的話,只要把它們變成「$a$ 乘一坨東西加 $b$ 乘一坨東西」的形式,那兩坨東西就分別是我們想要的 $x,y$。有個 $a \bmod b$ 在那邊不太好做,要是 $a,b \geq 0$ 的話那 $a \bmod b$ 就是 $a - \lfloor \frac ab \rfloor b$,這樣整理一下就可以得到
\[ bx' + \left(a - \left\lfloor\frac ab\right\rfloor b \right)y' = ay' - b\left( x' - \left\lfloor \frac ab \right\rfloor y' \right) \]
所以我們就知道
\[ x = y' \qquad y = x' - \left\lfloor \frac ab \right\rfloor y' \]
是一組 $ax+by=g$ 的整數解。別忘記遞迴要有終止條件,上面的輾轉相除法終止條件是 $b=0$,在 $b=0$ 的時候,$ax+by=\gcd(a,b)=a$ 的一個整數解就是 $(1,0)$,這樣我們就可以遞迴地算出一個 $ax+by=g$ 的整數解:

pair<int, int> exgcd(int a, int b) { // 回傳 {x, y}
    if (b == 0) return {1, 0};
    auto [x, y] = exgcd(b, a % b); // a % b == a - a / b * b
    return {y, x - a / b * y};
}

遞迴的過程和本來的輾轉相除法是一樣的,所以複雜度同樣是 $O(\log a + \log b)$。

既然擴展歐幾里得演算法總是求得出解,那我們就知道:

Theorem 3
貝祖定理

對於整數 $a,b,c$,$ax+by=c$ 有整數解若且唯若 $c$ 是 $g=\gcd(a,b)$ 的倍數。

負的係數

剛才我們偷偷假設 $a,b \geq 0$,如果 $a,b$ 是負的怎麼辦?負整數的最大公因數是什麼?STL 規定 std::gcd(a, b) 的回傳值是 $\lvert a \rvert$ 與 $\lvert b \rvert$ 的最大公因數,而 $a=b=0$ 時回傳 $0$,所以理當總要是個非負整數吧,但如果讀者自己試試看執行上面的程式碼,呼叫 gcd(10, -12),會得到 $-2$,而呼叫 exgcd(10, -12),會得到 $x=y=1$,代入 $10x-12y$ 後也是得到 $-2$,似乎跟我們想像的最大公因數定義不太一樣。其實會這樣只是因為負號不影響整除,所以輾轉相除法根本不在乎,結果有可能會多個負號,但取絕對值後是正確的。解決方法也很簡單,發現得到的結果多一個負號,就把 $x,y$ 也都變相反數就好了。

注意到 a / b * y 中的 a / b 其實跟我們剛剛說好的 $\lfloor \frac ab \rfloor$ 不一樣,可以這樣寫是因為 C++ 裡面的 a % b 無論正負,就是 a - a / b * b,所以把上面式子的 $a \bmod b$ 都當成 a % b,$\lfloor \frac ab \rfloor$ 都當成 a / b,就可以發現這個程式碼是會給出正確結果的。

所有整數解

整數解很可能不只一個吧!要怎麼知道其他的整數解呢?把 $ax+by=c$ 畫成一條二維平面上的直線的話,它的斜率是 $-b/a$,所以把直線上的一個點右移 $b$ 單位,再下移 $a$ 單位,就會得到另一個也在直線上的點,而且這個位移量都是整數,我們還可以幫它們除以 $g=\gcd(a,b)$,右移 $b/g$ 再下移 $a/g$ 也可以。

Observation 1
整數解的一般解

如果 $ax+by=c$ 的一個整數解是 $(x_0,y_0)$,令 $g=\gcd(a,b),\ a'=a/g,\ b'=b/g$,那麼所有整數解都可以被表示成

\[ (x_0 + kb', y_0 - ka') \]

其中 $k$ 是任意整數。所有這種形式的 $(x,y)$ 也都是整數解。

證明:

運算過程與答案的範圍

在純算 $\gcd(a,b)$ 的輾轉相除法中,我們知道 $a,b$ 只會一直往 $0$ 靠近而已,過程中不會出現絕對值比本來還大的數字,但擴展歐幾里得演算法好像就不太好說了,畢竟推答案的時候出現了 * y,還出現了各種不同運算,有沒有可能最後算出來的解超級大,或者是雖然最後會得到的解很小,但運算過程中用到了很大的數字呢?

Lemma 2
擴展歐幾里得演算法解的範圍

用擴展歐幾里得演算法求到的 $ax+by=g=\gcd(a,b)$ 的解 $(x,y)$,$a,b \neq 0$ 的話,滿足 $\lvert x \rvert \leq \lvert b/g \rvert,\ \lvert y \rvert \leq \lvert a/g \rvert$。

而 $a$ 或 $b$ 為 $0$ 時答案就是 $(0,1)$ 或 $(1,0)$,因此,我們就可以放心地相信答案和過程中出現的數字不會超過輸入本身的值域。

擴展歐幾里得與模逆元


我們學會了快速找 $ax+by=c$ 的解,所以我們只要把擴展歐幾里得演算法拿去解 $ab-mk=1$,得到一個 $b$ 後,把它 $\bmod m$ 一下弄進 $[0,m-1]$ 的範圍內,就可以在 $O(\log m)$ 的時間找模逆元了!

剛剛學到的事情也多告訴了我們一些關於模逆元的性質,像是:

Lemma 3
模逆元的唯一性

對於一個正整數 $m$ 和非負整數 $0 \leq a < m$,要是它在模 $m$ 下有模逆元的話,那模逆元在 $[0,m-1]$ 內是唯一的。

證明:

模逆元:小結


如果是第一次接觸到這種概念的讀者,可能還是心理會覺得怪怪的,剛剛我們舉了 $\frac{5!}{3!} \bmod 11$ 的例子,但用個莫名其妙的方法算對了感覺就像巧合一般。其實只要理解「乘上 $a^{-1}$ 和乘 $a$ 互為逆向操作」這件事,就比較可以接受一點了:在一般運算的世界,$\frac{5!}{3!}$ 乘上它的分母 $3!$ 後就是 $5!$,而在模 $11$ 運算下的 $5! \times (3!)^{-1}$ 再乘上 $3!$ 的話,因為 $(3!)^{-1} \times 3! \equiv 1 \pmod{11}$,$5! \times (3!)^{-1} \times 3!$ 就也是 $5!$,很合理地 $5 \times (3!)^{-1} \bmod m$ 就應該要和 $\frac{5!}{3!} \bmod m$ 是一樣的。

總而言之,模逆元就是在模運算之下的除法,不僅可以用在要計算本來就是整數的結果,有些輸出本來可能不是整數的題目會說「可以證明,答案一定可以被表示成 $\frac{P}{Q}$ 的形式,$P,Q$ 互質,請輸出 $P \times Q^{-1} \bmod m$ 的結果」,例如 $\frac{2}{5} \bmod 11$ 就是 $7$,在這種題目就是放心地模運算,要用除法就模逆元,最後算出來的答案就會自然而然的是 $P \times Q^{-1} \bmod m$,也不用管什麼 $P,Q$ 互質,畢竟模逆元有某種除法效果,分子分母的公因數會自己被消掉。

在要求答案要模的題目,通常為了確保 $[1,m-1]$ 都有模逆元,$m$ 會是一個質數,像常見的 $10^9+7$ 和 $998244353$ 都是質數。不過如果題目要模的數字不是質數,那有可能是在暗示你不要用除法,或者得用一些技巧來湊出答案,這點要特別注意。

習題


習題
Line

給 $[-2 \times 10^9, 2 \times 10^9]$ 內的三個整數 $A,B,C$,求一組 $Ax+By+C=0$ 的 $x,y$ 都是 $[-5 \times 10^{18}, 5 \times 10^{18}]$ 內整數的解。不存在的話輸出 $-1$。

習題
Euler Totient Function
Source:SPOJ ETF

有 $T$ 筆詢問,每個詢問給一個整數 $n$,求 $\varphi(n)$。

條件限制
  • $T \leq 20000$
  • $1 \leq n \leq 10^6$
習題
Exponentiation II
Source:CSES 1712

有 $n$ 筆測資,每筆測資有三個整數 $a,b,c$,求 $a^{b^c} \bmod 10^9+7$。

條件限制
  • $1 \leq n \leq 10^5$
  • $0 \leq a,b,c \leq 10^9$
習題
乘法逆元
Source:LOJ 110

給正整數 $n$ 和質數 $p$,求 $1$ 到 $n$ 的整數模 $p$ 下的模逆元。

條件限制
  • $1 \leq n \leq 3 \times 10^6$
  • $n < p < 20\,000\,528$