模取幂运算

**模冪**（英語：modular exponentiation）是一種對模進行的冪運算，在計算機科學，尤其是公開密鑰加密方面有一定用途。

模冪運算是指求整數 $b$ 的 $e$ 次方 $b^{e}$ 被正整數 $m$ 所除得到的餘數 $c$ 的過程，可用數學符號表示為 $c=b^{e}{\bmod {m}}$ 。由 $c$ 的定義可得 $0\leq c<m$ 。

例如，給定 $b=5$ ， $e=3$ 和 $m=13$ ， $5^{3}=125$ 被 13 除得的餘數 $c=8$ 。

指數 $e$ 為負數時可使用擴展歐幾里得算法找到 $b$ 模除 $m$ 的模逆元 $d$ 來執行模冪運算，即：

　
: $c=b^{e}{\bmod {m}}=d^{-e}{\bmod {m}}$ ，其中 $e<0$ 且 $b\cdot d\equiv 1{\pmod {m}}$ 。

即使在整數很大的情況下，上述模冪運算其實也是易於執行的。然而，計算模的離散對數（即在已知 $b$ ， $c$ 和 $m$ 時求出指數 $e$ ）則較為困難。這種類似單向函數的表現使模冪運算可用於加密算法。

## 直接算法

計算模冪的最直接方法是直接算出 $b^{e}$ ，再把結果模除 $m$ 。假設已知 $b=4$ ， $e=13$ ，以及 $m=497$ ，要求 $c$ ：

　
: $c\equiv 4^{13}{\pmod {497}}$

可用計算器算得 413 結果為 67,108,864，模除 497，可得 _c_ 等於 445。

注意到 $b$ 只有一位， $e$ 也只有兩位，但 $b^{e}$ 的值卻有 8 位。

強加密時 $b$ 通常至少有 1024 位。考慮 $b=5\times 10^{76}$ 和 $e=17$ 的情況， $b$ 的長度為 77 位， $e$ 的長度為 2 位，但是 $b^{e}$ 的值以十進位表示卻已經有 1304 位。現代計算機雖然可以進行這種計算，但計算速度卻大大降低。

用這種算法求模冪所需的時間取決於操作環境和處理器，時間複雜度為 $\mathrm {O} (e)$ 。

## 內存優化

這種方法比第一種所需要的步驟更多，但所需內存和時間均大為減少，其原理為：
給定兩個整數 $a$ 和 $b$ ，以下兩個等式是等價的：

　
: $c{\bmod {m}}=(a\cdot b){\bmod {m}}$

　
: $c{\bmod {m}}=[(a{\bmod {m}})\cdot (b{\bmod {m}})]{\bmod {m}}$

算法如下：

1. 令 $c=1$ ， $e'=0$ 。
2. $e'$ 自增 1。
3. 令 $c=(b\cdot c){\bmod {m}}$ .
4. 若 $e'<e$ ，則返回第二步；否則， $c$ 即為 $c\equiv b^{e}{\pmod {m}}$ 。

再以 $b=4$ ， $e=13$ ， $m=497$ 為例說明，算法第三步需要執行 13 次：

 * $e'=1\cdot c=(1\cdot 4){\bmod {4}}97=4{\bmod {4}}97=4$
* $e'=2\cdot c=(4\cdot 4){\bmod {4}}97=16{\bmod {4}}97=16$
* $e'=3\cdot c=(16\cdot 4){\bmod {4}}97=64{\bmod {4}}97=64$
* $e'=4\cdot c=(64\cdot 4){\bmod {4}}97=256{\bmod {4}}97=256$
* $e'=5\cdot c=(256\cdot 4){\bmod {4}}97=1024{\bmod {4}}97=30$
* $e'=6\cdot c=(30\cdot 4){\bmod {4}}97=120{\bmod {4}}97=120$
* $e'=7\cdot c=(120\cdot 4){\bmod {4}}97=480{\bmod {4}}97=480$
* $e'=8\cdot c=(480\cdot 4){\bmod {4}}97=1920{\bmod {4}}97=429$
* $e'=9\cdot c=(429\cdot 4){\bmod {4}}97=1716{\bmod {4}}97=225$
* $e'=10\cdot c=(225\cdot 4){\bmod {4}}97=900{\bmod {4}}97=403$
* $e'=11\cdot c=(403\cdot 4){\bmod {4}}97=1612{\bmod {4}}97=121$
* $e'=12\cdot c=(121\cdot 4){\bmod {4}}97=484{\bmod {4}}97=484$
* $e'=13\cdot c=(484\cdot 4){\bmod {4}}97=1936{\bmod {4}}97=445$

因此最終結果 $c$ 為 445，與第一種方法所求結果相等。

與第一種方法相同，這種算法需要 $\mathrm {O} (e)$ 的時間才能完成。但是，由於在計算過程中處理的數字比第一種算法小得多，因此計算時間至少減少了 $\mathrm {O} (e)$ 倍。

算法偽代碼如下：

```
**function** modular_pow(b, e, m)
    **if** m = 1
        **then** **return** 0
    c := 1
    **for** e' = 0 **to** e-1
        c := (c * b) **mod** m
    **return** c
```

## 從右到左二位算法

第三種方法結合了第二種算法和平方求冪原理，使所需步驟大大減少，同時也與第二種方法一樣減少了內存占用量。

首先把 $e$ 表示成二進位，即：

　
: $e=\sum _{i=0}^{n-1}a_{i}2^{i}$

此時 $e$ 的長度為 $n$ 位。對任意 $i$ （ $0\leq i<n$ ）， $a_{i}$ 可取 0 或 1 任一值。由定義有 $a_{n-1}=1$ 。

 $b^{e}$ 的值可寫作：

　
: $b^{e}=b^{\left(\sum _{i=0}^{n-1}a_{i}2^{i}\right)}=\prod _{i=0}^{n-1}\left(b^{2^{i}}\right)^{a_{i}}$

因此答案 $c$ 即為：

　
: $c\equiv \prod _{i=0}^{n-1}\left(b^{2^{i}}\right)^{a_{i}}\ ({\mbox{mod}}\ m)$

### 偽代碼

下述偽代碼基於布魯斯・施奈爾所著《應用密碼學》。其中 base，exponent 和 modulus 分別對應上式中的 $b$ ， $e$ 和 $m$ 。

```
**function** modular_pow(base, exponent, modulus)
    **if** modulus = 1 **then** **return** 0
    Assert :: (modulus - 1) * (modulus - 1) does not overflow base
    result := 1
    base := base **mod** modulus
    **while** exponent > 0
        **if** (exponent **mod** 2 == 1):
           result := (result * base) **mod** modulus
        exponent := exponent >> 1
        base := (base * base) **mod** modulus
    **return** result
```

注意到在首次進入循環時，變量 base 等於 $b$ 。在第三行代碼中重複執行平方運算，會確保在每次循環結束時，變量 base 等於 $b^{2^{i}}{\bmod {m}}$ ，其中 $i$ 是循環執行次數。

本例中，底數 $b$ 的指數為 $e=13$ 。
指數用二進位表示為 1101，有 4 位，故循環執行 4 次。
4 位數字從右到左依次為 1，0，1，1。

首先，初始化結果 $R$ 為 1，並將 _b_ 的值保存在變量 $x$ 中：

　
: $R\leftarrow 1\,(=b^{0}){\text{且 }}x\leftarrow b$ .

　
: 第 1 步 第 1 位為 1，故令 $R\leftarrow R\cdot x{\text{ }}(=b^{1})$ ；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{2})$ 。

　
: 第 2 步 第 2 位為 0，故不給 _R_ 賦值；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{4})$ 。

　
: 第 3 步 第 3 位為 1，故令 $R\leftarrow R\cdot x{\text{ }}(=b^{5})$ ；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{8})$ 。

　
: 第 4 步 第 4 位為 1，故令 $R\leftarrow R\cdot x{\text{ }}(=b^{13})$ ；

　
: 這是最後一步，所以不需要對 _x_ 求平方。

綜上， $R$ 為 $b^{13}$ 。

以下計算 $b=4$ 的 $e=13$ 次方對 497 求模的結果。

初始化：

　
: $R\leftarrow 1\,(=b^{0})$ 且 $x\leftarrow b=4$ 。

　
: 第 1 步 第 1 位為 1，故令 $R\leftarrow R\cdot 4{\pmod {497}}\equiv 4{\pmod {497}}$ ；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{2})\equiv 4^{2}\equiv 16{\pmod {497}}$ 。

　
: 第 2 步 第 2 位為 0，故不給 _R_ 賦值；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{4})\equiv 16^{2}{\pmod {497}}\equiv 256{\pmod {497}}$ 。

　
: 第 3 步 第 3 位為 1，故令 $R\leftarrow R\cdot x{\text{ }}(=b^{5})\equiv 4\cdot 256{\pmod {497}}\equiv 30{\pmod {497}}$ ；

　
: 令 $x\leftarrow x^{2}{\text{ }}(=b^{8})\equiv 256^{2}{\pmod {497}}\equiv 429{\pmod {497}}$ 。

　
: 第 4 步 第 4 位為 1，故令 $R\leftarrow R\cdot x{\text{ }}(=b^{13})\equiv 30\cdot 429{\pmod {497}}\equiv 445{\pmod {497}}$ ；

綜上， $R$ 為 $4^{13}{\pmod {497}}\equiv 445{\pmod {497}}$ ，與先前算法中所得結果相同。

該算法時間複雜度為 O(log exponent)。指數 exponent 值較大時，這種算法與前兩種 O(exponent) 算法相比具有明顯的速度優勢。例如，如果指數為 220 = 1048576，此算法只需執行 20 步，而非 1,048,576 步。

### Lua 實現

```
**function** modPow(b,e,m)
        **if** m == 1 **then**
                **return** 0
        **else**
                **local** r = 1
                b = b % m
                **while** e > 0 **do**
                        **if** e % 2 == 1 **then**
                                r = (r*b) % m
                        **end**
                        e = e >> 1     --Lua 5.2或更早版本使用e = math.floor(e / 2)
                        b = (b^2) % m
                **end**
                **return** r
        **end**
**end**
```

## 軟體實現

鑑於模冪運算是計算機科學中的重要操作，並且已有高效算法，所以許多程式語言和高精度整數庫都有執行模冪運算的函數：

* Python 內置的 `pow()`（求冪）函數 [1] （頁面存檔備份，存於網際網路檔案館）
* .NET 框架的 `BigInteger` 類的 `ModPow()` 方法
* Java 的 `java.math.BigInteger` 類的 `modPow()` 方法
* Perl 的 `Math::BigInt` 模塊的 `bmodpow()` 方法 [2] （頁面存檔備份，存於網際網路檔案館）
* Raku 內置的 `expmod` 例程
* Go 的 `big.Int` 類的 `Exp()`（求冪）方法 [3] （頁面存檔備份，存於網際網路檔案館）
* PHP 的 BC Math 庫的 `bcpowmod()` 函數 [4] （頁面存檔備份，存於網際網路檔案館）
* GNU 多重精度運算庫（GMP）的 `mpz_powm()` 函數 [5] （頁面存檔備份，存於網際網路檔案館）
* FileMaker Pro 的 `@PowerMod()` 函數（以 1024 位 RSA 加密為例）
* Ruby 的 `openssl` 包的 `OpenSSL::BN#mod_exp` 方法 [6] （頁面存檔備份，存於網際網路檔案館）
* HP Prime 計算器的 CAS.powmod () 函數 [7]

## 另見

* 蒙哥馬利算法，用於計算模很大時的餘數。

## 參考資料

## 外部連結

* Schneier, Bruce. Applied Cryptography: Protocols, Algorithms, and Source Code in C, Second Edition 2nd. Wiley. 1996. ISBN 978-0-471-11709-4.
* Paul Garrett, Fast Modular Exponentiation Java Applet （頁面存檔備份，存於網際網路檔案館）
* Gordon, Daniel M. A Survey of Fast Exponentiation Methods (PDF). Journal of Algorithms (Elsevier BV). 1998, **27** (1): 129–146 [2019-11-30]. ISSN 0196-6774. doi:10.1006/jagm.1997.0913. （原始內容存檔 (PDF) 於 2019-08-27）.