はてなブログ謎オラクルによるブログのスパム認定が解除されたので、yoshi-campの参加記続きを書いていきまーす。
Day 0とDay 1はこっち:
Day 2
[mitsu] ネコちゃん式 ~安全な曲線の生成~
yoshi-camp一番の鬼講師ことmitsu先生が今回も楕円曲線についてお話してくれました。
タイトル通り安全な楕円曲線を作る方法を考えるのですが、これは楕円曲線の位数を計算すると確認できます。 つまり、この講義では楕円曲線の位数を計算するためのSchoofのアルゴリズムの原理と高速化についてお勉強しました。
Frobenius写像とHasseの定理
Hasseの定理は楕円曲線の本には必ずといっていいほど書いてある定理です。
ここで、左辺の絶対値の中身はFrobenius写像という写像のトレースになっています。 Frobenius写像も有名ですが、次のように定義される写像のことです。
楕円曲線の位数がFrobenius写像のトレースから計算できるというのは、pwnerの私でも知っている話でした。 しかし、具体的にトレースをどう計算するかは一切知らなかったので、ここから新しい内容に入ります。
等分多項式
Schoofのアルゴリズム本体に入る前に、アルゴリズムで必要となる楕円曲線の等分多項式を計算します。 つまり、n-ねじれ点のx座標を根に持つ多項式を求めます。
多項式自体は簡単に求められます。 単純に楕円曲線の加算公式を使い、 を計算して方程式が得られるからです。
簡単でも面倒なので、mitsuくんからの天啓を使ってみましょう。 楕円曲線に付随する多項式を計算して出来上がったものがこちらになります。
これ自体は等分多項式ではなく、等分多項式は、次のようにを消去した関数として定義されています。
の定義に照らしあわせると、は帰納的に計算できます。実装してみましょう。
from functools import cache @cache def division_polynomial(E, m): f = lambda m: division_polynomial(E, m) a, b = E.a4(), E.a6() PR = PolynomialRing(E.base(), 'x') x = PR.gen() if m == 0: return PR(0) elif m == 1 or m == 2: return PR(1) elif m == 3: return 3*x^4 + 6*a*x^2 + 12*b*x - a^2 elif m == 4: return 2*(x^6 + 5*a*x^4 + 20*b*x^3 - 5*a^2*x^2 - 4*a*b*x - 8*b^2 - a^3) elif m % 2 == 0: n = m // 2 return (f(n+2)*f(n-1)^2 - f(n-2)*f(n+1)^2) * f(n) else: n = (m-1) // 2 F = 4*(x^3 + a*x + b) if n % 2 == 1: return f(n+2)*f(n)^3 - F^2*f(n-1)*f(n+1)^3 else: return F^2*f(n+2)*f(n)^3 - f(n-1)*f(n+1)^3
この関数で求めた方程式の根が、ねじれ点のx座標になっていることを確認します。
E = EllipticCurve(GF(13), [2, 3]) m = 9 f = division_polynomial(E, m) print(f) for x, _ in f.roots(): try: P = E.lift_x(x) print(P, m*P) except ValueError: pass
結果:
9*x^40 + 8*x^38 + x^37 + 3*x^36 + 10*x^35 + 12*x^33 + 10*x^32 + 11*x^31 + 3*x^30 + 5*x^29 + 9*x^28 + 2*x^27 + x^26 + x^23 + 4*x^21 + 9*x^20 + 6*x^19 + 6*x^18 + 2*x^17 + 12*x^16 + 2*x^15 + 6*x^14 + 5*x^13 + 3*x^12 + 9*x^11 + 7*x^10 + 9*x^9 + 2*x^8 + 10*x^7 + 3*x^6 + 3*x^5 + x^4 + 4*x^3 + 11*x^2 + 6*x (0 : 9 : 1) (0 : 1 : 0) (11 : 2 : 1) (0 : 1 : 0) (9 : 3 : 1) (0 : 1 : 0) (3 : 7 : 1) (0 : 1 : 0)
できました!
Schoofのアルゴリズム
Schoofのアルゴリズムの本質は中国人剰余定理です。 Frobenius写像のトレースは、未満の小さな素数を使って
となるようなを計算し、これらからを中国人剰余定理で復元します。 中国人剰余定理とHasseの定理から、必要な素数は
を満たすまで取ってくれば良いです。 (Hasse boundに入れるために、実際にはこれで求まったからを引く必要があります。)
もう1つの重要情報として、Frobenius写像の特性多項式は
で表されます。綺麗ですね。 受験生のお供、Cayley-Hamiltonの定理より
ですので、ここに-ねじれ点Pを代入すれば、
つまり、とおくと
であり、とおくと、
が成り立ちます。 この式で計算できないのはだけですが、この値は十分に小さいので全探索します。 (実際にはもっと賢い方法があり、後で登場します。)
実装の上での問題
ここまでの理論を実装しようと思うと、いくつかの問題が発生します。
1つ目の問題が、多項式の計算方法です。 さきほどのを求めるための方程式は具体的な点に対して計算するのではなく、一般的なねじれ点に対して計算するので、多項式環の上で計算する必要があります。 これは、「等分多項式を満たす」「楕円曲線上の演算である」という制約を加えれば良いので、剰余環の上で計算することで解決します。
2つ目の問題が、点の計算です。 一般には大きいので、といった計算は愚直にやると終わりません。
まず、ですが、今考えている上では
なので、の部分を楕円曲線の式で置き換えれば、y座標もほぼの計算で終わります。
同様にも
と計算できます。
点の加算についても高速化を考えてみましょう。 2つの点を(ただし、)とおき、を考えます。
まず、のとき、
とおくと、
となります。ではの項にが登場しますが、2乗されているので[tex: y2]となり、ここが楕円曲線のに関する式に置換できます。
(です。今後もはこの定義で書きます。) 一方、ではは消えませんが、中のと外側のをまとめて
と書けます。
のときも、
と表せるので同様にの形で書けます。
3つ目がで剰余が取れない問題です。 の計算で分数があるので逆元を計算する必要がありますが、これがを法として存在しないことがあります。 どうすればいいでしょうか?
逆元を取りたい多項式をとおきます。 具体的には、またはです。
今、は法の下で逆元を持たないと仮定しているので、です。 のとき、の定義よりです。 のとき、まずとの関係についてですが、の根は、つまり位数が2の点のx座標になります。 しかし、のは3以上の素数なので位数は奇数であり、とは互いに素です。 そのためはの因子となるわけですが、の定義よりです。 したがって、いずれのケースでもです。
つまり、の逆元が存在しない場合は
上で計算すれば良いです。(の元であれば何でもいい。)
最後に、講義中は特に言及されませんでしたが、実装している上で一番詰まった点として、の扱いがあります。
私は上での楕円曲線演算を、SymbolicEllipticCurveOverQuotientRing
というクラスを作って実装したのですが、この際sageの機能を使って
self.PR = PolynomialRing(E.base(), 'xx') self.xx = self.PR.gen() self.QR = self.PR.quotient(f, 'x') self.x = self.QR.gen()
と定義しました。はこれで良いのですが、楕円曲線の加算の式を見返すとが登場します。困ったなぁ。 これには、優秀な先人のyoshi-camp参加記を参考にを1にする方法で解決しました。
というのも、今回の計算においてはどれだけ加算や乗算を繰り返しても、y座標は必ずという形になることを知っています。 そして、Schoofのアルゴリズムにおいては点の比較のみが重要となるため、の部分を1に置き換えても問題ありません。
Schoofのアルゴリズムの実装
長かったよ......
class SymbolicEllipticCurvePointOverQuotientRing(object): def __init__(self, E, x, y): self.E = E self.x, self.y = x, y def __eq__(self, other): return self.E == other.E \ and self.x == other.x \ and self.y == other.y def __add__(self, other): P, Q, O = self, other, self.E(0, 0) if P == O: return Q elif Q == O: return P a1, b1 = P.x, P.y a2, b2 = Q.x, Q.y if a1 != a2: n = b1 - b2 d = a1 - a2 else: n = 3 * a1^2 + self.E.a d = 2 * b1 * self.E.F dd = self.E.f.parent()(d.list()) g = self.E.f.gcd(dd) if g == 1: r = n / d a3 = r^2 * self.E.F - a1 - a2 b3 = r * (a1 - a3) - b1 return self.E(a3, b3) else: EE = SymbolicEllipticCurveOverQuotientRing(self.E.base, g) PP = EE(P.x.list(), P.y.list()) QQ = EE(Q.x.list(), Q.y.list()) RR = PP + QQ return self.E(RR.x.list(), RR.y.list()) def __mul__(self, n): P = self Q = self.E(0, 0) while n > 0: if n % 2 == 1: Q = P + Q P = P + P n >>= 1 return Q def __rmul__(self, n): return self * n class SymbolicEllipticCurveOverQuotientRing(object): def __init__(self, E, f): self.base = E self.PR = PolynomialRing(E.base(), 'xx') self.xx = self.PR.gen() self.QR = self.PR.quotient(f, 'x') self.x = self.QR.gen() self.f = f self.a, self.b = self.PR(E.a4()), self.PR(E.a6()) self.F = self.x^3 + self.a*self.x + self.b def __call__(self, x, y): return SymbolicEllipticCurvePointOverQuotientRing( self, self.QR(x), self.QR(y) ) def frobenius_trace(E): a, b, p = E.a4(), E.a6(), E.base().order() l_prod, li = 1, 2 t_list, l_list = [], [] while l_prod < 4 * sqrt(p): li = next_prime(li) fl = division_polynomial(E, li) EE = SymbolicEllipticCurveOverQuotientRing(E, fl) pl = p % li L = EE(EE.x^(p^2), EE.F^((p^2-1)//2)) + pl*EE(EE.x, 1) for tli in range(1, li): R = tli * EE(EE.x^p, EE.F^((p-1)//2)) if L == R: t_list.append(tli) l_list.append(li) l_prod *= li break t = CRT(t_list, l_list) if t >= 2 * sqrt(p): t -= l_prod return t def schoof(E): return -frobenius_trace(E) + E.base().order() + 1
mitsuくんは真面目にPolynomialRingのPolynomialRingを作り、剰余環を真面目に実装したそうです。sage力が欲しい。
実行してみましょう。
p = random_prime(1<<31, 1<<32) a = randrange(0, p) b = randrange(0, p) F = GF(p) E = EllipticCurve(F, [a, b]) print("sage's order:", timeit('E.order()')) print("ptr-schoof:", timeit('schoof(E)')) print("Order:", E.order()) print("Order:", schoof(E))
やはり、tを総当りしていることもあり、sageの実装には多くの場合勝てません。
sage's order: 625 loops, best of 3: 146 ns per loop ptr-schoof: 5 loops, best of 3: 350 ms per loop Order: 435017356 Order: 435017356
素数が64-bitくらいなら数十秒以内に求まります。
Schoof-Elkie-Atkinsという高速化されたアルゴリズムもお勉強しました。 これのwriteupも途中まで書いていましたが、書けば書くほど無限に内容が増えるので全部コメントアウトして消しました。
モジュラー多項式と楕円曲線の関係を、複素平面上の格子を経由して学ぶことでアハ体験できる講義内容になっていました。 興味のある方はyoshi-campに参加しましょう。
[Xornet] Bilinyar Map in CTF
ここから人間味あふれる講義に戻ります。 Xornetさんの講義はタイトルがよくわかりませんが、楕円曲線でおなじみのペアリングについてお勉強できます。
ペアリングの重要な性質として、双線型性というものがあります。 双線型性は2変数関数に対する性質で、次のようなものです。
これより、スカラー倍に対して次の性質があります。
わかりやすいですね。 楕円曲線上のペアリングでは、次のようになります。
楽しいですね。 ペアリングというのは楕円曲線に対して1つ定まるようなものではなく、こういう性質さえ満たせば良いのでいろいろと定義できます。
Weil Pairing
親のペアリングより見たペアリングですね。 名前の読み方が初見殺しのペアリングランキングでは常にNo. 1です。
n-ねじれ点を2つ渡すと、有限体上のn乗根が返ってきます。
n-ねじれ点どうしでしかペアリングは計算できず、という性質があるため、線形従属な2点では1が返ってきてしまいます。
Tate Pairing
名前がかわいいペアリングランキングでいつも上位にいる常連さんです。
という定義です。怖いですね。
といった性質があります。
DDH問題
ペアリングを使って解ける代表的な問題としてDDH問題があります。 これは、が既知のときに、与えられた点がかどうかを判別する問題です。
Tate Pairingを使うと
なので、のときに両者の値が一致します。
確認してみましょう。
p = 0x3a05ce0b044dade60c9a52fb6a3035fc9117b307ca21ae1b6577fef7acd651c1f1c9c06a644fd82955694af6cd4e88f540010f2e8fdf037c769135dbe29bf16a154b62e614bb441f318a82ccd1e493ffa565e5ffd5a708251a50d145f3159a5 k = 2 E = EllipticCurve(GF(p), [1, 0]) assert (p^k - 1) % E.order() == 0 P = E.random_point() aP = 314159265 * P bP = 2718281828 * P Q1 = 314159265 * 2718281828 * P Q2 = randrange(0, E.order()) * P print(aP.tate_pairing(bP, n=E.order(), k=k)) print(P.tate_pairing(Q1, n=E.order(), k=k)) print(P.tate_pairing(Q2, n=E.order(), k=k))
埋め込み次数*3を知っている必要があるので、講義で用意されていた埋め込み次数の小さい楕円曲線を使っています。
15128028... 15128028... 85119939...
と、DDH問題が解けていることがわかります。
MOV Attack
楕円曲線上の離散対数問題について、ベースポイントの埋め込み次数が小さいときに使える攻撃です。 古のyoshi-campで「猫でもわかる楕円曲線暗号」を発表したときに理論無視実装コードだけやりました。
の位数をとすると、ペアリングが定義できます。 適当な点に対して
が成り立つので、拡大体上での離散対数問題として解くことができ、楕円曲線上でやるよりはマシに計算できます。*4
以下、パラメータは「Cyber Apocalypse CTF 2022: Intergalactic Chase / MOVs Like Jagger」という問題から取っているそうです。
a = -35 b = 98 p = 434252269029337012720086440207 E = EllipticCurve(GF(p), [a, b]) ec_order = 434252269029337012720086440208 G = E(16378704336066569231287640165, 377857010369614774097663166640) x = randint(1, ec_order) xG = x * G k = 1 while (p^k - 1) % ec_order: k += 1 print("Embedding degree:", k) Ee = EllipticCurve(GF(p^k), [a, b]) Ge, xGe = Ee(G), Ee(xG) T = Ee.random_point() Q = (T.order() // gcd(T.order(), G.order())) * T assert G.order() / Q.order() in ZZ a = Ge.weil_pairing(Q, G.order()) b = xGe.weil_pairing(Q, G.order()) print("e(G, Q) =", a) print("e(xG, Q) =", b) y = b.log(a) print(x) print(y)
1分くらいで解けました。
Embedding degree: 2 e(G, Q) = 20183620868354029403390095407*z2 + 340156830587302922110931927752 e(xG, Q) = 382029521819594423509410822768*z2 + 81650503901539005236114419028 25571984854458949809049655853 25571984854458949809049655853
BLS
BLS(Boneh-Lynn-Shacham)は署名方式の名前で、CDH問題を解くのが難しいがDDH問題は簡単に解けるような曲線を利用した署名です。 検証者はDDH問題をペアリングで解くことで、署名の正当性を検証します。
まず、鍵生成で署名者は鍵ペアを生成します。
次に、この鍵で署名します。 署名したいメッセージをハッシュ関数を通して楕円曲線上の点に変換し、署名を作ります。
- 署名:
最後に、検証者は公開鍵と署名、メッセージから署名の正当性を調べます。
- 検証:を検証する。
正しい署名であれば、
より両者の値は一致するはずです。
実装してみましょう。
import hashlib def h(m: bytes): h = hashlib.sha256(m).digest() x = Integer(int.from_bytes(h, "little")) while True: rhs = K(x**3 + a*x + b) if rhs.is_square(): return E.lift_x(x) x += 1 q = 0x3a05ce0b044dade60c9a52fb6a3035fc9117b307ca21ae1b6577fef7acd651c1f1c9c06a644fd82955694af6cd4e88f540010f2e8fdf037c769135dbe29bf16a154b62e614bb441f318a82ccd1e493ffa565e5ffd5a708251a50d145f3159a5 a, b = 1, 0 o = 21992493417575896428286087521674334179336251497851906051131955410904158485314789427947788692030188502157019527331790513011401920585195969087140918256569620608732530453375717414098148438918130733211117668960801178110820764957628836 r = 2344807743588794553832391292516279194397209456764712786969868894104465782493871625440983981162219279755855675661203 k = 2 assert (q^2 - 1) % r == 0 assert o == 2**2 * r**2 K = GF(q) E = EllipticCurve(K, [a, b]) G = E.gens()[0] msg = b"neko to innu" H = h(msg) x = 5963 P = x * G S = x * H wp1 = S.tate_pairing(G, n=o, k=k) wp2 = H.tate_pairing(P, n=o, k=k) wp3 = E.random_point().tate_pairing(G, n=o, k=k) print(wp1 == wp2) print(wp3 == wp2)
実行すると、署名が正しいときのみTrueになります。
BLS12-381
BLS12-381はBLSで使える楕円曲線の一種です。*5 埋め込み次数12、381-bitの素数qで定義されており、上の楕円曲線でペアリングします。
12次の拡大体についてですが、既約多項式で逐次的に拡大して構成するのが効率的だそうで、は次のような多段の拡大になっています。
埋め込み次数に対してには位数の部分群が複数存在し、これにはペアリングもニッコリ。
しかし、は計算が重いという問題があり、当然その上での楕円曲線の演算も重くなってしまいます。
実際、381-bitの素数でGF(q^12)
を作ろうとするとPCがキュイィィィインと泣きます。
なぜ計算量を増やしてこんな面倒なことをするかというと、計算機上で実装する際の効率が良いみたいな恩恵もあるらしいです。
でも計算できなければ使い物になりません。 BLSで登場する演算は、ペアリングを除いて点のスカラー倍だけです。 そこでBLS12-381では、次の2つの曲線が使われます。
ここでは-1の平方根なので、次のようにも書けます。
さて、1つ目の楕円曲線上の点は上の点とみなしても良いですが、2つ目の楕円曲線は曲線が違うのでダメそうです。 そこで、次のような写像を考えます。ただし、です。
をの式に代入すると
となり、に移ります。*6 (これでにした理由がわかりましたね。)
実装パート
Circle City Con 2021というCTFのRandom is Also Validという問題を解きます。 この問題では謎の方法でBLS12-381が実装されており、その上でDDH問題を解かなければなりません。 問題の詳細はXornetさんのブログを読みましょう。
こんな感じ。
import json out = json.load(open("out.txt", "r")) q = 0x1a0111ea397fe69a4b1ba7b6434bacd764774b84f38512bf6730d2a0f6b0f6241eabfffeb153ffffb9feffffffffaaab r = 0x73eda753299d7d483339d80809a1d80553bda402fffe5bfeffffffff00000001 k = 12 assert (q^k - 1) % r == 0 G1_x = 4 G1_y = 0x0a989badd40d6212b33cffc3f3763e9bc760f988c9926b26da9dd85e928483446346b8ed00e1de5d5ea93e354abe706c G1_mul = 0x396c8c005555e1568c00aaab0000aaab G2_x = [2] G2_y = [0x013a59858b6809fca4d9a3b6539246a70051a3c88899964a42bc9a69cf9acdd9dd387cfa9086b894185b9a46a402be73, 0x02d27e0ec3356299a346a09ad7dc4ef68a483c3aed53f9139d2f929a3eecebf72082e5e58c6da24ee32e03040c406d4f] G2_mul = 0x5d543a95414e7f1091d50792876a202cd91de4547085abaa68a205b2e5a7ddfa628f1cb4d9e82ef21537e293a6691ae1616ec6e786f0c70cf1c38e31c7238e5 F12 = GF(q^12, 'w', modulus=x^12-2*x^6+2) w = F12.gen() E = EllipticCurve(F12, [0, 4]) def twist(Px, Py): z = w^-1 if isinstance(Px, list): Px = sum(map(lambda v: v[1]*(w^6-1)^v[0], enumerate(Px))) if isinstance(Py, list): Py = sum(map(lambda v: v[1]*(w^6-1)^v[0], enumerate(Py))) return z^2 * F12(Px), z^3 * F12(Py) G1 = G1_mul * E(G1_x, G1_y) G2 = G2_mul * E(twist(G2_x, G2_y)) o = E.order() dec = lambda v: list(map(lambda x: int(x, 16), v)) m = 0 for i, bitinfo in enumerate(out): PA = E(dec(bitinfo['P_A_x']), dec(bitinfo['P_A_y'])) PB = E(twist(dec(bitinfo['P_B_x']), dec(bitinfo['P_B_y']))) wp2 = PA.weil_pairing(PB, n=r) if len(bitinfo['P_C_x']) == 1 and len(bitinfo['P_C_y']) == 1: PC = E(dec(bitinfo['P_C_x']), dec(bitinfo['P_C_y'])) wp1 = PC.weil_pairing(G2, n=r) else: PC = E(twist(dec(bitinfo['P_C_x']), dec(bitinfo['P_C_y']))) wp1 = G1.weil_pairing(PC, n=r) m = (m << 1) | int(wp1 == wp2) if i % 8 == 7: print(int.to_bytes(int(m), i//8+1, 'little')[::-1])
詰まりポイントとしては、PC.weil_pairing(G)
なのかG.weil_pairing(PC)
なのかで結果が変わるので注意しましょう。(1敗)
CCC{Diffie_0r_Hellman_wh4t_d0_y0U_d3c1d3?}
ツイストできたので満足ですが、実は宿題のzipfelを解かなければなりません。 これはそのうちやりますやります。
[ネコチャン] Tornado Cats
マネーロンダリング Tornado Cashというミキサープロトコルとゼロ知識証明についての勉強です。
いろいろとスクリプトを組んでいたのですが、PCを入れ替えて全部消えたのでここで公開できる成果物がありません!
しかし、すべての講義の中でたぶん一番資料が充実しており、なんか10,000円くらいの鈍器本として売ってそうな内容でした。 時間がなくてTornado Cash本体に到達していませんが、いつかkitakitsune*7ですべてが解明されると思います。
Day 3
帰りは基本バスの中ですやすやしてました。 箱根駅に着いたら、大通りから離れた場所にめちゃくちゃおいしい焼肉屋があり、もぐもぐしてました。
あとはなんか源泉かけ流し謎温泉に入ったり、1個50円くらいのお菓子を買ったりしていました。
写真集
おわりに
楽しかったです。 また春か夏にやりましょう。 @everyone
あと、まだ宿泊費を振り込んでない人は、これを見たら振り込みましょう。宿から焼肉からすべて負担させられると破産しちゃいます。