繰り返し二乗法・ダブリング

繰り返し二乗法

べき乗を高速に求めることができるアルゴリズム。

ある数を \(N\) 乗するとき、計算量 \(\mathcal{O}(\log{N})\) でべき乗が求まる。

アルゴリズム

例として、 \(3^{45} \pmod{17}\) を求めることを考える。

指数を2のべき乗の和で表すと

\[\begin{split}\begin{eqnarray} 3^{53} &=& 3^{(1 + 4 + 16 + 32)} \\ &=& 3^1 \cdot 3^4 \cdot 3^{16} \cdot 3^{32} \end{eqnarray}\end{split}\]

のように変形できる。 \(3^k\) は \(3^{k-1}\) を2乗すれば求められるから、\(3^1, 3^2, 3^4, 3^8, \ldots\) をこの順に求め、上式に含まれる値のみ掛け算した値を解とすればよい。

実装

繰り返し二乗法の利用例として、 \(x^y \pmod z\) を求める関数の実装を示す。

なお、Pythonの組み込み関数 pow も同様に振舞う。

def modpow(x, y, z):
    L = y.bit_length()
    a = 1
    for i in range(L):
        if y & 1 == 1:
            a = a * x % z
        x = x * x % z
        y >>= 1
    return a

print(modpow(3, 53, 17), pow(3, 53, 17))

5 5

ダブリング

ある状態に対して同じ操作を \(N\) 回適用した後の状態について、操作を2のべき乗の和に分解することで高速に求めることができるアルゴリズム。

とりうる状態の数を \(M\) とすると、計算量は次の通り。

• 前処理: \(\mathcal{O}(M\log{N})\)

• クエリ; \(\mathcal{O}(\log{N})\)

また、先の \(x^y \pmod z\) はダブリングで求められる。

前処理では、法 \(z\) のもとでとりうるすべての値 \(0\leq j < z\) について繰り返し二乗法を行う。

配列 \(Z[i][j] := j \cdot 3^{2^i}\) を定義し、この配列に繰り返し二乗法の結果を格納する。

実際に \(x^y \pmod z\) を求めるときは、剰余 \(j\) とべき乗の指数 \(i\) をキーとして、配列から順に取り出せばよい。

\[\begin{split}\begin{eqnarray} 1 \cdot 3^{1} = Z[0][1] &\equiv& 3 \pmod{17} \\ 3 \cdot 3^{4} = Z[2][3] &\equiv& 5 \\ 5 \cdot 3^{16} = Z[4][5] &\equiv& 5 \\ 5 \cdot 3^{32} = Z[5][5] &\equiv& 5 \\ \therefore 3^{53} &\equiv& 5\\ \end{eqnarray}\end{split}\]

実装

\(x^y \pmod z\) のダブリングでの実装例を示す。

def modpow_doubling(x, y, z):
    L = y.bit_length()
    Z = [[0 for j in range(z)] for i in range(L)]
    for j in range(z):
        Z[0][j] = j * x % z

    for i in range(L-1):
        for j in range(z):
            Z[i+1][j] = Z[i][Z[i][j]]

    a = 1
    for i in range(L):
        if y & 1 == 1:
            a = Z[i][a]
        y >>= 1
    return a

print(modpow_doubling(3, 53, 17), pow(3, 53, 17))

5 5