Exhaustive Search法
いわゆる全探索法であり, 力任せな方法. 与えられる楕円曲線$E$とその上の点$P$, $Q$について$m = \#P$とする. この時, ECDLPの解$d$は$1 \leq d \leq m $の不等式を満たす. このため, この$d$の範囲を全て計算してみることでECDLPは解かれる. 計算量的には$O(n)$となり, 一番ナイーブな方法.
Baby-step Giant-step法
一般的に適用できる中で最も有名なアルゴリズムだと思われる. これは任意の巡回群に対して適用することができる方法で, ここではECDLPに限定して解説する. 基本的な考え方としては, $d$を求めたい離散対数とすると, $d$はある$r$について$d = ir+j$と表せるので, これを使って $$ \begin{align} dP &= Q\\ irP + jP &= Q\\ jP &= Q - irP \end{align} $$ とすることができる. ここで, $j < r$となることから次のように考えられる. まず最初に$P, \ldots, rP$を計算して格納しておく. 次に, $T = -rP$は事前に計算できるのでこれを計算する. そして$Q + iT$を$i$を増やしながら計算していくことによっていつかは$jP = Q + iT$となるような$i, j$を手に入れることができる. すると, 最初の式から$d = ir + j$なのでこれでECDLPは解けたことになる.
多くの場合, $r = \lceil\sqrt{\#P}\rceil$が用いられる. このため, 計算量は$O(\sqrt{\#P})$となる. これはExhastive Searchよりも十分良い結果となっている.
from ecpy import * # https://github.com/elliptic-shiho/ecpy/ import math def bsgs(E, P, Q, ordP): # m = sqrt(#P) (round up) m = int(math.sqrt(ordP) + 0.5) # Calculate Baby-step table baby = [] for i in xrange(m): baby += [(i, i * P)] h = -P * m t = Q # Main loop: Giant-step for j in xrange(m): b = filter(lambda x: x[1] == t, baby) if len(b) == 1: return b[0][0] + j * m t = t + h def main(): F = FiniteField(17) E = EllipticCurve(F, 3, 7) P = E(4, 10, 1) Q = E(10, 0, 1) ordP = 14 x = bsgs(E, P, Q, ordP) assert x*P == Q print x # 7 if __name__ == "__main__": main()
Pollard Rho法
ECDLPについては最も成果を挙げている手法. Baby-step Giant-stepと同様, これも任意の巡回群について適用できる方法だが, Baby-step Giant-step法とは異なり確率的(乱択)アルゴリズムである. 実際には素因数分解アルゴリズムとしてのPollard Rhoアルゴリズムのほうが有名なので, 考え方をそちらで説明する. 実際の所, 考え方は共通しているので素因数分解の方法のみ理解していればECDLPへの応用もすぐに可能である.
このアルゴリズムでは乱数列を生成しながらその衝突(素因数分解の対象の数との最大公約数が1以外の値であること)を検出し, 衝突時の値から失敗又は素因数を返す. フロイドの循環検出法を用いて, $i$番目と$2i$番目の乱数を生成するのみとすることで, $O(1)$の記憶容量が存在すれば良いように作られている. この衝突する数列の様子を図にした際にギリシャ文字の$\rho$に見えることから, Pollard Rho法と呼ばれるようになった.
from gmpy import gcd def rand(x, n): return (x*x+3) % n def pollard_rho_factor(n): x = 2 y = 2 d = 1 while d == 1: x = rand(x, n) y = rand(rand(y, n), n) d = gcd(abs(x-y), n) if d == n: return False return d def main(): print pollard_rho_factor(18446744073709551617) print pollard_rho_factor(10023859281455311421) if __name__ == "__main__": main()
これは演算を差し替えれば任意の有限巡回群に適用できることが分かる.よって, これをECDLPへ適用することを考える. ECDLPへ適用する際には衝突の定義を変えなければならない. 以降, 楕円曲線$E$上の点$P$, $Q$について$Q = dP$となる$d$を求めるアルゴリズムとして考える.
- 乱数列$a_i$, $b_i$について $R_i = a_i P + b_i Q$を定義する.
- $R_i = R_{2i}$となる$i$が存在するとする.
- $R_i = a_i P + b_i Q = (a_i + b_i d)P$なので, 以下が成り立つ.
$$ R_i = R _ {2i} \iff a_i + b_i d = a _ {2i} + b _ {2i} d \mod \#E $$
- よって, $d = (a_i - a _ {2i})\cdot(b _ {2i} - b_i)^{-1} \mod \#E$と計算できる.
最後の段階で, 逆元が存在しない場合には失敗を返す.
以上を踏まえてコードにすると以下のようになる. 関数$f$は乱数を生成するために少々ごちゃごちゃしたことをしているが, これについてはPollard's rho algorithm for logarithms - Wikipediaを参考にした. (よく見れば, 先の条件を満たす乱数に見える値を作っているだけということが読み取れる. )
from ecpy import * # https://github.com/elliptic-shiho/ecpy/ from ecpy.util import modinv def f(alpha, beta, x, a, b, order): if x.x % 3 == 0: return (beta+x, a, (b+1) % order) elif x.x % 3 == 1: return (x * 2, (a*2) % order, (b*2) % order) elif x.x % 3 == 2: return (alpha+x, (a+1) % order, b) def pollard_rho_ECDLP(alpha, beta, curve, order): a, b, x = 0, 0, curve(0, 1, 0) A, B, X = a, b, x i = 0 while i < order: x, a, b = f(alpha, beta, x, a, b, order) X, A, B = f(alpha, beta, X, A, B, order) X, A, B = f(alpha, beta, X, A, B, order) if x == X and not A == a == b == B == 0: # except x == X == O r = (B - b) % order return (modinv(r, order) * (a - A)) % order i += 1 def main(): p = 100003 r = 100379 # order of `E` a = -7 b = 104 Gx = 91802 Gy = 58807 F = FiniteField(p) E = EllipticCurve(F, a, b) x = 12345 G = E(Gx, Gy) P = G * x print pollard_rho_ECDLP(G, P, E, r) if __name__ == "__main__": main()
また, このアルゴリズムは不完全である. 位数が合成数の場合には逆元が計算できなくなってしまうため, もしちゃんとしたECDLP Solverを書くのであれば位数を素因数分解し各因数それぞれについてECDLPを解き, 中国人剰余定理を用いて元の値を計算するような工夫が必要となる.
Pollard Kangaroo法 (Lambda法)
考案者の名前が同じなことから見当が付くかと思うが, このアルゴリズムはRho法とよく似ている. このアルゴリズムも楕円曲線上の点を整数へと写す写像$f$を用いるが, 衝突の判定が微妙に異なる. ここでは, Rho法と同様に楕円曲線$E$上の点$P, Q$について, $Q = xP$を満たすような$a \lt x \lt b$を探すことを目的とする.
- $N$を決める. これはステップ数とも言える値で, 任意に定めることができる. (ここではまず100とする. )
- 数列$T$を以下の規則に従って定める. ここで, $f$は楕円曲線上の点から整数へのハッシュ関数である (適当な乱数関数のシード値にx座標を渡すだけでも良い.). $$ \begin{align} T_0 &= b P\newline T _ {i + 1} &= T_i + f(T_i) P\newline \end{align} $$
- 同様に, 数列$W$を以下の規則に従って定める. $$ \begin{align} W_0 &= Q\newline W _ {i + 1} &= W_i + f(W_i) P\newline \end{align} $$
- 関数$D_A(k)$を以下のように定める. ただし, $D_A(0) := 0$である. $$ D_A(k) := \sum _ {i = 1} ^ k f(A_i) $$
- ここで, $T_i = W_j$となるような$i$, $j$が存在したとする. すると, $T_i = (b + D_T(i)) P$, $W_j = (x + D_W(j)) P$であることから$x = b + D_T(i) - D_W(j)$となり, ECDLPが解ける. この時, $D_W(j) \gt b - a + D_T(i)$が成立している時は範囲内に離散対数が存在しないと言えるため, 失敗と判定する.
- 失敗した時には, $N$や$f$を取り替えて再度最初からやり直す.
これを$T$, $W$をカンガルーに例えた話に言い換える. 手懐けた(Tame)カンガルー$T$と野生の(Wild)カンガルー$W$が居るときに, $W$がどこから来たのかを知るために捕まえることを考える. そこで, $T$を$N$回ジャンプさせ, ジャンプした地点それぞれに罠を仕掛ける. どこかでその罠に$W$が引っかかると, $W$を捕まえることができる. この時, 元の位置からの移動距離$D_T(i)$, $D_W(j)$がわかっているので, その距離を引くことで元の位置が分かる. *1 Lambda($\lambda$)法の名前の由来は, 手懐けたカンガルーがジャンプしていった地点のどこかで野生のカンガルーが罠にかかる様子を図にすると$\lambda$の形のように見えることから来ている.
以下に実際のコードを示す.
from ecpy import * # https://github.com/elliptic-shiho/ecpy/ def f(x, n): return int((x.x ** 2 + 5) % n) def pollard_lambda(curve, P, Q, order, a=0, b=None): if b is None: b = order - 1 N = 100 INCR = 1 print "[+] Constructing xs/ys..." xs = [P * b] for i in xrange(1, N): xs += [xs[-1] + P * f(xs[-1], order)] ys = [Q] for i in xrange(N - 1): ys += [ys[-1] + P * f(ys[-1], order)] while True: D = [sum(map(lambda x: f(x, order), xs[:n])) % order for n in xrange(N)] ds = [sum(map(lambda x: f(x, order), ys[:n])) % order for n in xrange(N)] print "[*] Find Collision" for j, d in enumerate(ds): Y = ys[j] if Y in xs: xj = xs.index(Y) if d > (b-a + D[xj]): print "[-] Failed... Retry with new N" break else: m = (b + D[xj] - ds[j]) % order if m * P == Q: print "[+] Found!" return m else: print "[-] Failed... Can't tame a 'roo!" N += INCR # additional xs/ys for i in xrange(INCR): xs += [xs[-1] + P * f(xs[-1], order)] for i in xrange(INCR): ys += [ys[-1] + P * f(ys[-1], order)] def main(): p = 100003 r = 100379 # order of `E` a = -7 b = 104 Gx = 91802 Gy = 58807 F = FiniteField(p) E = EllipticCurve(F, a, b) x = 12345 G = E(Gx, Gy) P = G * x print pollard_lambda(E, G, P, r, a=10000, b=20001) if __name__ == "__main__": main()
このアルゴリズムは, ステップ数に比例した記憶領域を必要とする.
指数計算法(Index-Calculus)
離散対数問題(DLP)に対しては非常に効果的だが, そのままではECDLPに対してはそれほど効果がない方法. 複雑な方法になるためここでは実際のコードを書くことはしない.
まず, 因子基底(factor base)と呼ばれる基底を取る. これはDLPの場合は素数を適当にいくつか取ってくれば良い. これらは(基底という名前の通り)$\mathbb{Z}$係数の線形空間をなすことが必要なことに注意する. 以降, 因子基底の集合を$\mathfrak{D} = \{F_1, F_2, \ldots, F_n\}$, 因子基底のなす空間を$\mathcal{F}$と表す. 次に, ランダムに選んだ$\alpha _ i$, $\beta _ i$を用いて$\alpha_i P + \beta_i Q$を計算し, これが$\mathcal{F}$に含まれていることを確認する. もし含まれていないならば$\alpha_i$, $\beta_i$を取り直す. $\alpha_i$, $\beta_i$はそれぞれ$n$個ずつ集める必要がある. そうすると, $\mathbf{\alpha} = (\alpha_1, \ldots, \alpha _ n)$, $\mathbf{\beta} = (\beta_1, \ldots, \beta_n)$というベクトルに対して以下の等式を満たすような$\mathbf{f} _ i = (f _ {1i}, \ldots, f _ {ni})$が存在する. $$ \mathbf{\alpha} P + \mathbf{\beta} Q = \sum _ {i = 1} ^ n \mathbf{f} _ i F _ i $$ 次に, $\mathbf{f} _ i$から$n$次正方行列$ M $を作る.
$$ M := \begin{pmatrix} f _ {11} & f _ {12} & \cdots & f _ {1n}\newline f _ {21} & f _ {22} & \cdots & f _ {2n}\newline \vdots & \vdots & \ddots & \vdots\newline f _ {n1} & f _ {n2} & \cdots & f _ {nn}\newline \end{pmatrix} $$
この$ M $は $M \begin{pmatrix}F _ 1, F _ 2, \ldots, F _ n\end{pmatrix}^T = \mathbf{\alpha} P + \mathbf{\beta} Q$をを満たす. この関係式に対して行基本変形を施し, $ M $の一番下の行が全て0になる(上三角行列になる)ようにする. この変形された行列を$ M' $として, 対応する$\mathbf{\alpha}$, $\mathbf{\beta}$を$\mathbf{\alpha}^\ast$, $\mathbf{\beta}^\ast$とする.
$M'$の一番下の行は全て0であるため, $\mathbf{\alpha}^\ast _ n P + \mathbf{\beta}^\ast _ n Q = \mathbf{0}$という式ができる. この式を$Q$について解くことで, $Q = xP$を満たす$x$を求めることができる.
この指数計算法はDLPでは他と比べて早い(準指数時間)で終了するアルゴリズムであり, 非常に有用と言える. しかし, 今回考えているECDLPでは主に点を分解する段階において非常に時間がかかることが分かっており, 用いられることは少ない.
Pohlig-Hellman法
この方法は巡回群の位数が小さな素因数に分解される時に有効な手法である. 楕円曲線$E$上の点$P$, $Q$について$Q = xP$を満たす$x$を求めることを考える.
- 楕円曲線$E$の位数$\#E$を素因数分解する. ここでは$\#E = p_1^{e_1}p_2^{e_2}\cdots p_k^{e_k}$と素因数分解されるとする.
- 各素因数$p_i^{e_i}$について, $\beta_i = \#E / p_i^{e_i}$を計算する.
- 全探索やRho法等を用いて$x_i \cdot \beta_i P = \beta_i Q$を満たす$x_i$をそれぞれ求める. $\beta_i$をかけたことにより$\beta_i P$, $\beta_i Q$はそれぞれ位数$p_i^{e_i}$を持つ点になるため, 素因数が小さければ$x_i$は非常に高速に求めることができる.
- 各$x_i$はそれぞれ$\mod p_i^{e_i}$されているので, 中国人剰余定理より全ての$i$について$x \equiv x_i \mod p_i^{e_i}$を満たす$x$を求める. これにより, ECDLPが解ける.
以上を実装したものを以下に示す.
from ecpy import * # https://github.com/elliptic-shiho/ecpy/ from ecpy.util import crt, modinv from primefac import factorint # pip install primefac def factor(n): ret = [] factors = factorint(n) return factors, sorted([p**factors[p] for p in factors.keys()]) def exhaustive_search(curve, P, Q, a, b): i = a while i < b: if P * i == Q: return i i += 1 return None def pohlig_hellman(curve, P, Q, order): if tuple(P) == (0, 1, 0): # P * x == O => x == #P return order _, factors = factor(order) d = [] for x in factors: m = (order / x) d += [exhaustive_search(curve, P * m, Q * m, 0, x)] return crt(d, factors) def main(): p = 100003 F = FiniteField(p) E = EllipticCurve(F, -7, 1) order = 100080 Gx = 71233 Gy = 31138 G = E(Gx, Gy) x = 12345 P = G * x print pohlig_hellman(E, G, P, order) if __name__ == "__main__": main()
ここでは全探索関数を実装して利用しているが, ここは任意のECDLPを解くアルゴリズムを適用できる.
参考文献
- The Arithmetic of Elliptic Curves - J.H.Silverman 著
- クラウドを支えるこれからの暗号技術 - 光成滋生 著
- Modern Cryptanalysis: Techniques for Advanced Code Breaking - Christopher Swenson 著
- 暗号理論入門 原書第3版 - J.A.ブーフマン 著 林芳樹 訳
- Monte Carlo Methods for Index Computation $\pmod p$ - J.M.Pollard
- デファクト暗号の評価 - 宮地充子
- Attacks on the Elliptic Curve Discrete Logarithm Problem - Kenneth J.Giuliani
*1:多少わかりづらいかも知れないが, これは元論文にも書かれているちゃんとした解説である