合成数次元NTRUに対するGentry攻撃の実装

武蔵野アドベンドカレンダー 2017の19日目です.

概要

• 復習編
  • 多項式環とその剰余環, 巡回行列
• NTRU暗号編
  • 耐量子計算機暗号の一角であるNTRU暗号を紹介する
  • sagemath8.0上でNTRU暗号を実装する
• 攻撃編
  • $n$が合成数のNTRU暗号に対するGentry攻撃 (EUROCRYPT2001) のアイデアを伝える
  • sagemath8.0上でGentry攻撃を実装する

復習編

復習編 in cocalcでjupyter notebookを共有しています. cocalcにアカウントを作って, ノートをコピーして遊んでみてください.

整数

まずはsagemathでの演算に慣れましょう.

sagemathについては 木村先生の「Sagemathと数学」が参考になります.

f = 7
g = 2
print f + g
print f * g

実行結果

9
14

整数の足し算と掛け算が出来ます.

次にこれを法10で ($\bmod{10}$の世界で) 考えましょう.
$n$で割った余りだけ考える剰余環を$\mathbb{Z}_n := \mathbb{Z}/(n\mathbb{Z})$とします.

Z_10 = Integers(10)
print ZZ
print Z_10

f = Z_10(7)
g = Z_10(2)
print f + g
print f * g

実行結果

Integer Ring
Ring of integers modulo 10
9
4

$4 = 7 * 2 \bmod{10}$が出力されていることがわかります.

多項式

次に多項式の演算をしてみましょう. $f = 1 + x^2 + x^3$と$g = 2 + x^1 + x^6$の足し算と掛け算をやってみます.

f = 1 + x^2 + x^3
g = 2 + x^1 + x^6
print f + g
print f * g

実行結果

x^6 + x^3 + x^2 + x + 3
(x^6 + x + 2)*(x^3 + x^2 + 1)

簡単に演算が出来ました.

多項式でも整数の場合と同様に剰余環を考えることができます. $x^5-1$で割った余りだけ考える剰余環を$R5 := \mathbb{Z}[x]/(x^5-1)$と書くことにします.

R.<x> = ZZ[]
R5 = R.quotient(x^5-1)
print R
print R5

f = R5(1 + x^2 + x^3)
g = R5(2 + x^1 + x^6)
print f + g
print f * g

実行結果

Univariate Polynomial Ring in x over Integer Ring
Univariate Quotient Polynomial Ring in xbar over Integer Ring with modulus x^5 - 1
xbar^3 + xbar^2 + 2*xbar + 3
2*xbar^4 + 4*xbar^3 + 2*xbar^2 + 2*xbar + 2

ここで変数名がxbarになっているのは, R5の定義で変数を指定していないからです. R5.<y> = R.quotient(x^5-1)とすると, xbarがすべてyに代わります.

今, 法$x^5-1$を考えているので, $x^5 - 1 \equiv 0 \pmod{x^5-1}$です. これは$x^5 \equiv 1 \pmod{x^5-1}$ということです. よって, この剰余環では$x^5$は$1$になってしまいます. 掛け算した結果, $x^{5+i}$という項が出てきたら, すべて$x^i$に置き換えることになります.

確認してみましょう.

\begin{align}
f * g 
&\equiv x^9 + x^8 + x^6 + x^4 + 3 x^3 + 2 x^2 + x + 2 \\
&\equiv x^4 + x^3 + x^1 + x^4 + 3 x^3 + 2 x^2 + x + 2 \\
&\equiv 2 x^4 + 4 x^3 + 2 x^2 + 2 x + 2 \pmod{x^5-1}
\end{align}

あってますね.

さらに整数環も剰余環にして, $R5_{10} := \mathbb{Z}_10[x]/(x^5-1)$を考えて演算してみましょう.

R_10.<x>  = Z_10[]
R5_10 = R_10.quotient(x^5-1)
print R_10
print R5_10

f = R5_10(1 + x^2 + x^3)
g = R5_10(2 + x^1 + x^6)
print f + g
print f * g

実行結果

Univariate Polynomial Ring in x over Ring of integers modulo 10
Univariate Quotient Polynomial Ring in xbar over Ring of integers modulo 10 with modulus x^5 + 9
xbar^3 + xbar^2 + 2*xbar + 3
2*xbar^4 + 4*xbar^3 + 2*xbar^2 + 2*xbar + 2

はい, できました.

R_10やR5_10の要素をRの要素に__自然に__戻すには以下のようにします.

f = R_10(1 + x^2 + x^3)
print f.change_ring(ZZ)
print f.change_ring(ZZ).parent()

g = R5_10(2 + x^1 + x^6)
# print g.change_ring(ZZ) -> Error
print g.lift()
print g.lift().parent()
print g.lift().change_ring(ZZ)
print g.lift().change_ring(ZZ).parent()

実行結果

x^3 + x^2 + 1
Univariate Polynomial Ring in x over Integer Ring
2*x + 2
Univariate Polynomial Ring in x over Ring of integers modulo 10
2*x + 2
Univariate Polynomial Ring in x over Integer Ring

多項式とベクトルと巡回行列

多項式環は加法だけ考えて群とみなせば, （無限次元）ベクトル空間となります. 適当な次数以上の項を無視してベクトルに表してみましょう.

# 多項式をベクトルに直す関数
# f_0 + f_1 x + ... + f_{n-1} x^{n-1} -> (f_0,f_1,...,f_{n-1})
# R(x^i)としておかないと, i=0のときにエラー
def vectorize(f,n):
    return vector(ZZ,[f.monomial_coefficient(R(x^i)) for i in range(n)])

f = 1 + x^2 + x^3
g = 2 + x^1 + x^6
print f + g

fv= vectorize(f,10)
gv= vectorize(g,10)
print fv + gv

実行結果

x^6 + x^3 + x^2 + x + 3
(3, 1, 1, 1, 0, 0, 1, 0, 0, 0)

法$x^n-1$でも同様です. 都合上, R_10やR5_10の要素はいったんRの要素とみなしてから, ベクトルに直します.

f = R10(1 + x^2 + x^3)
g = R5_10(2 + x^1 + x^6)
# vectorize(f,5) -> TypeError

fv = vectorize(f.change_ring(ZZ),5)
gv = vectorize(g.lift().change_ring(ZZ),5)
print fv
print gv

# ベクトルを多項式 in Rに直す関数
def vec_to_poly(z):
    tmp = R(0)
    for i in range(len(z)):
        tmp += R(z[i]*x**i)
    return tmp

print vec_to_poly(fv)
print vec_to_poly(gv)

実行結果

(1, 0, 1, 1, 0)
(2, 2, 0, 0, 0)
x^3 + x^2 + 1
2*x + 2

次に掛け算を考えましょう. ベクトルとベクトルの間には掛け算という演算がないので, うまく定義する必要があります.

多項式の掛け算をうまく分解してやると,

\begin{align}
f * g
&\equiv \sum_{i=0}^{4} f_i x^i \cdot g \\
&\equiv \sum_{i=0}^{4} f_i \cdot (x^i \cdot g) \pmod{x^5-1}.
\end{align}

となります. これを踏まえて,

(f_0,f_1,f_2,f_3,f_4)
 \cdot 
\begin{pmatrix}
g \\
x g \bmod{x^5-1} \\
x^2 g \bmod{x^5-1} \\
x^3 g \bmod{x^5-1} \\
x^4 g \bmod{x^5-1} \\
\end{pmatrix}
= \sum_{i=0}^{4} f_i \cdot (x^i \cdot g \bmod{x^5-1})
= f * g

とベクトルと行列の積として考えることが出来ます.

今回は法が$x^5-1$なので, $g$に対応する行列は__巡回行列__になります.

f = R5_10(1 + x^2 + x^3)
g = R5_10(2 + x^1 + x^6)
print f * g

# 多項式 in Rを巡回行列に直す関数
def circulant_matrix(f,n):
    return matrix(ZZ,n,n,lambda i,j: f.monomial_coefficient(R(x^((j-i) % n))))

fv = vectorize(f.lift().change_ring(ZZ),5)
fm = circulant_matrix(f.lift().change_ring(ZZ),5)
gm = circulant_matrix(g.lift().change_ring(ZZ),5)

print fv
print fm
print gm
print fv * gm == vectorize((f*g).lift().change_ring(ZZ),5)
print fm * gm == circulant_matrix((f*g).lift().change_ring(ZZ),5)

実行結果

2*xbar^4 + 4*xbar^3 + 2*xbar^2 + 2*xbar + 2
(1, 0, 1, 1, 0)
[1 0 1 1 0]
[0 1 0 1 1]
[1 0 1 0 1]
[1 1 0 1 0]
[0 1 1 0 1]
[2 2 0 0 0]
[0 2 2 0 0]
[0 0 2 2 0]
[0 0 0 2 2]
[2 0 0 0 2]
True
True

$f$の方も巡回行列にしてみましょう.
$f$に対応する巡回行列と$g$に対応する巡回行列の積が$f*g$に対応する巡回行列になります. 可換環$K$について, $K[x]/(x^n-1)$は, $n$次の$K$要素巡回行列がなす環$\mathrm{Cir}(K,n)$と同型ということがいえます.

NTRU暗号

この章と次の章に対応するノートブックはNTRU-CS in cocalcにあります.

NTRU暗号はJ. Hoffstein, J. Pipher, and J.H. Silvermanの3人によって1990年代後半に提案された公開鍵暗号です. 多項式環に基づいており, これまでさまざまな攻撃に耐えてきています. 最近, 耐量子計算機暗号が盛り上がっていることもあり, 注目を集めています.

パラメータ設定

正整数$n$, $q$について,

\begin{align}
R_{n} &:= \mathbb{Z}[x]/(x^n-1) \\
R_{q,n} &:= \mathbb{Z}_q[x]/(x^n-1)
\end{align}

と定義します.
また, $R_{n}$の部分集合として, ${0,1}$係数多項式のうち係数$1$が丁度$d$個ある多項式の集合を$S(n,d)$と書きます.

# parameters ====================
# NTRU256.2
n,p,q,df,dg,dr=[256,2,127,35,35,22]
# ====================
R.<x>  = ZZ[]
Rq.<x> = Integers(q)[]
Rqn = Rq.quotient(x^n-1)
Rp.<x> = Integers(p)[]
Rpn = Rp.quotient(x^n-1)
# ====================
print "n,p,q", n,p,q
print "df,dg, dr", df, dg, dr

NTRU256.2というのは後で攻撃対象になるパラメータです. 当時のNTRU論文ではNTRU263.2 n,p,q,df,dg,dr=[263,2,127,35,35,22] 等のパラメータ集合が提案されています. 約20年前の設定なので数字が小さいですね.

アルゴリズム詳細

以上のパラメータを用いた場合, NTRU論文によると, $p=2$のときの鍵の作り方は以下です.

1. $f \in S(n,d_f)$と$g \in S(n,d_g)$をランダムに選ぶ
2. $f$が$R_{p,n}$と$R_{q,n}$で逆元を持たない場合, 1に戻る
3. $h = g/f \bmod{q}$とする
4. $f$を秘密鍵とする
5. $h$を公開鍵とする

暗号化はこんな感じです

1. $m \in R_{p,n}$を平文とする
2. $r \in S(n,d_r)$をランダムに選ぶ
3. $c := p * h * r + m \bmod{q}$とする
4. $c$を暗号文とする

復号はこう.

1. $d := f * c \bmod{q}$を計算し, $d$を$R_n$の多項式だと思う
2. $m' := d/f \bmod{p}$を計算する

$d \equiv f * c \equiv p * g * r + m * f \pmod{q}$
なので, $d = pgr + mf$が整数上でも成立するなら,
$m' \equiv d/f \equiv (pgr + mf)/f \equiv mf/f \equiv m \pmod{p}$
として復号できます.
$g, r, m, f$はすべて$0,1$係数の多項式なので, 掛けても足しても小さいということで, $d = prg + mf$が整数上でも成立します.

実装

実装していきます.

sample(range(n),d)で 集合${0,\dots,n-1}$から$d$個重複無しで選んだリストが生成出来るので, これを使って, $S(n,d)$からのランダム元生成を実現できます. (random_tpoly(d))

# ランダム多項式の生成
def random_tpoly(d):
    f = R(0)
    for i in sample(range(n),d):
        f += R(x**i)
    return f

# 鍵生成
def skgen():
    f,g = [0,0]
    while not Rpn(f).is_unit() or not Rqn(f).is_unit():
        f = random_tpoly(df)
        g = random_tpoly(dg)
    return f, g

def pkgen(f,g):
    h = Rqn(g)/Rqn(f)
    return h

# 暗号化
def enc(h,m):
    r = random_tpoly(dr)
    c = p * h * Rqn(r) + Rqn(m)
    return c

# 復号
def dec(f,c):
    d = (c * Rqn(f)).lift()
    m = Rpn(d)/Rpn(f)
    return R(m.lift())

# テスト ====================
# 擬似乱数のシードを0に固定
set_random_seed(0)

f,g = skgen()
h = pkgen(f,g)

# ランダムな平文を作って, 暗号化して復号する
m = random_tpoly(int(n/2))
c = enc(h,m)
d = dec(f,c)

# 表示用にh.lift()とRの要素にする
print "pk", h.lift()
print "m", m
print "d", d
print "flag,", m == d

これを実行すると,

実行結果

pk :=  62*x^255 + 59*x^254 + (snip) + 74*x + 51
m :=  x^248 + x^247 + (snip) + x^4 + x^3
d :=  x^248 + x^247 + (snip) + x^4 + x^3
flag, True

m == dとなり, 無事に復号できました.

Coppersmith-Shamir攻撃

NTRUはCRYPTO 1996のランプセッションで公開され, 論文自体はANTS1998で出版されています. CoppersmithとShamirによる攻撃論文が, EUROCRYPT __1997__に採択されています (Lattice Attacks on NTRU | SpringerLink).

格子

Lattice Attacksというからには Lattice (格子) を使っています.
基底ベクトル$b_1,\dots,b_n$を適当な整数係数で足したベクトル全体の集合を格子と呼びます.
すなわち,

\Lambda = \{ \sum_i w_i b_i \mid w_i \in \mathbb{Z}\}

が格子です. あ, この稿ではベクトルはすべて列ベクトルです.

鍵の関係式

NTRU暗号の鍵をみると,

h = g/f \bmod{q}

という関係があります. CoppersmithとShamirは$H$を$h$に対応する巡回行列として,

L = 
\begin{pmatrix}
I_n & H \\
O_n & q I_n 
\end{pmatrix}
\in \mathbb{Z}^{2n \times 2n}

という行列の列ベクトルを基底とする格子を考えたときに, $(f,g)$というベクトルがこの格子に含まれることに気づきました. 証明は以下です.

1. 鍵の関係式から, $hf \equiv g \pmod{q}$が成立
2. なので, ある$k \in R_n$が存在して, $hf = g + qk$
3. よって, $hf - qk = g$
4. これと基底$L$を考えると, $(f,-k) L = (f, fh) + (0, -k q) = (f,g)$より, $(f,g)$は$L$が張る格子に含まれる.

LLLやBKZといった格子基底縮小アルゴリズムと呼ばれるアルゴリズムを使うと, 格子の基底を与えられたときに（そこそこ）短いベクトルを探せることが多いです.

$L$を格子基底縮小アルゴリズムに掛けて短いベクトルを探してやると, それが$(f,g)$なり$(fx^{i},gx^{n-i})$であることが期待できます.
これをプログラムにしてみましょう.

def vectorize(f,n):
    return vector(ZZ,[f.monomial_coefficient(R(x^i)) for i in range(n)])

def vec_to_poly(z):
    tmp = R(0)
    for i in range(len(z)):
        tmp += R(z[i]*x**i)
    return tmp

def circulant_matrix(f,n):
    return matrix(ZZ,n,n,lambda i,j: f.monomial_coefficient(R(x^((j-i) % n))))

# g = h * f mod qのL1ﾉﾙﾑがdgに一致するかどうかチェック
def check_g(h,f):
    return vectorize((h*Rqn(f)).lift().change_ring(ZZ),n).norm(1) == dg

# 公開鍵から格子の基底を作る関数 (Coppersmith and Shamirが提案)
def CS(h,n):
    Id = identity_matrix(n)
    Ze = zero_matrix(n)
    HH = circulant_matrix(h.lift().change_ring(ZZ),n)
    return block_matrix(ZZ,2,2,[Id,HH,Ze,q*Id])

def find_fv(h):
    # 基底を得て, それを格子基底簡約アルゴリズム (今回はBKZ) で簡約する
    L = CS(h,n)
    L = L.BKZ()

    # Lの横ベクトルを(f,g)と思って, fに対応するベクトルのL1ノルムがdfに等しいなら採用
    fv_candidates = []
    for i in range(2*n):
        fv = L[i][0:n]
        if fv.norm(1) == df:
            if min(fv.list()) >= 0:
                fv_candidates.append(fv)
            if max(fv.list()) <= 0:
                fv_candidates.append(-fv)
    if len(fv_candidates) == 0:
        raise ValueError('Oops, I cannot find any fv_candidates')
    else:
        return fv_candidates
    
def test_CS(seed,debug=true):
    set_random_seed(seed)
    f,g = skgen()
    h = pkgen(f,g)
    print h
    
    # Find fv
    if debug:
        print "start find_fv"
    try:
        fv_candidates = find_fv(h)
    except ValueError as err:
        print("Value Error: {0}".format(err))
    if debug:
        print "end find_fv"

    # f0s is a list of polynomials f0's in R = ZZ[x]
    f0s = [vec_to_poly(fv_candidates[i]) for i in range(len(fv_candidates))]
    # check_gでフィルター
    f0s = [f0 for f0 in f0s if check_g(h,f0)]
    return f,g,h,f0s

# 今回のパラメータ設定ではCS攻撃はうまくいかない
# print "!!! Start !!!"
# %time f,g,h,f0s = testCS(0,true)
# for i in range(len(f0s)):
#    print "i:", i
#    print "fi:", f0s[i]
#    print "gi:", R((h*Rqn(f0s[i])).lift())
# print "!!! End !!!"

Gentry攻撃

最後に dblp: Craig Gentryの2001年の論文であるKey Recovery and Message Attacks on NTRU-Composite (pdf)の攻撃方法を実装してみます.

この章に対応するノートブックはNTRU-Composite in cocalcにあります.

前置き

SilvermanはNTRU Tech. Rep. #11 の中で, nを2の冪にするとFFTが使えるので, 暗号化・復号が高速化出来るだろう, と述べました.

それに対してGentryが新しい攻撃を提案し, （当時の技術で）実装し, NTRU256.2というパラメータセットでも秘密鍵回復攻撃が成功することを示しました. 以下ではこの秘密鍵回復攻撃を追試してみます.

アイデア

$n$が$256$ということは, 非自明な約数$d$を持ちます. 例えば$d = 64$としましょう. すると, $x^n-1$は$x^d-1$で割り切れます.

$\theta \colon R_{n} \to R_{d}$を

\theta(\sum_{i=0}^{n-1} f_i x^i)
 := \sum_{i=0}^{n-1} f_i x^i \bmod{x^d - 1}
 = \sum_{i=0}^{d-1} \sum_{j=0}^{n/d-1} f_{jd+i} x^i

と定義すると, これは環準同型です.

環準同型性のチエック

$n = 10$, $d = 5$で試してみます.

Rq.<x>  = Integers(127)[]
Rqn.<xx> = Rq.quotient(x^10-1)
Rqd.<xx> = Rq.quotient(x^5-1)

f = Rqn(1 + x^2 + x^3)
g = Rqn(2 + x^1 + x^6)

print "Addition"
print f + g
print Rqd(f + g)
print Rqd(f) + Rqd(g)

print "Multiplication"
print f * g
print Rqd(f * g)
print Rqd(f) * Rqd(g)

実行結果

Addition
xbar^6 + xbar^3 + xbar^2 + xbar + 3
xbar^3 + xbar^2 + 2*xbar + 3
xbar^3 + xbar^2 + 2*xbar + 3
Multiplication
xbar^9 + xbar^8 + xbar^6 + xbar^4 + 3*xbar^3 + 2*xbar^2 + xbar + 2
2*xbar^4 + 4*xbar^3 + 2*xbar^2 + 2*xbar + 2
2*xbar^4 + 4*xbar^3 + 2*xbar^2 + 2*xbar + 2

環準同型性が確認できました.

鍵の関係式

Coppersmith-Shamir攻撃法では 鍵の関係式

h f + q k = g

から格子の基底を作成していました. この関係式に準同型写像$\theta$を適用すると,

\theta(h) \cdot \theta(f) + q \theta(k) = \theta(g)

となります.
自然なアイデアとして,

L_d = 
\begin{pmatrix}
I_d & H_d \\
O_d & q I_d 
\end{pmatrix}
\in \mathbb{Z}^{2d \times 2d}

を考えます. $(\theta(f),\theta(k)) L_d = ... = (\theta(f),\theta(g))$となり, $(\theta(f),\theta(g))$も小さくなった格子に含まれます.

さて, この環準同型は$f = f_0 + f_1 x^d + f_2 x^{2d} + f_3 x^{3d}$という$f \in R_{n}$を$\theta(f) = f_0 + f_1 + f_2 + f_3 \in R_d$に写しています.
秘密鍵$f$は$0,1$係数なので, $\lVert \theta(f)\rVert_{\infty}$は高々$4$ですし, $\lVert \theta(f) \rVert_{1} = d_f$です.

まとめると, $2d$次元に小さくなった基底$L_d$が張る格子の中に, まだまだ短いベクトル$(\theta(f),\theta(g))$が含まれています. これなら格子基底縮小アルゴリズムで見つけられそうです.

こうして小さい格子から秘密鍵の情報を求め, それを大きい格子から秘密鍵を求める際のヒントとして用いるというのがGentry攻撃のアイデアです.

実装

$d1 = 128$と$d2 = 64$と非自明な約数$d$を二つ使うので,
それぞれに環を定義していきます.

# ====================
# parameters
# ====================
debug = false
n,p,q,df,dg,dr=[256,2,127,35,35,22]
d2 = int(n/4)
d1 = int(n/2)
R.<x>  = ZZ[]
Rq.<x> = Integers(q)[]
Rqn = Rq.quotient(x^n-1)
Rqd2 = Rq.quotient(x^d2-1)
Rqd1 = Rq.quotient(x^d1-1)
Rp.<w> = Integers(p)[]
Rpn = Rp.quotient(w^n-1)
# ====================
# NTRU 
# ====================
def random_tpoly(p):
    f = R(0)
    for i in sample(range(n),df):
        f += R(x**i)
    return f

def skgen():
    f,g = [0,0]
    while not Rpn(f).is_unit() or not Rqn(f).is_unit():
        f = random_tpoly(p)
        g = random_tpoly(p)
    return f,g

def pkgen(f,g):
    h = Rqn(g)/Rqn(f)
    return h

def enc(h,m):
    r = random_tpoly(dr)
    c = p * h * Rqn(r) + Rqn(m)
    return c

def dec(f,c):
    d = (c * Rqn(f)).lift()
    m = Rpn(d)/Rpn(f)
    return R(m.lift())

# ====================
# Attack
# ====================
def vectorize(b,n):
    return vector(ZZ,[b.monomial_coefficient(R(x^i)) for i in range(n)])

def circulant_matrix(b,n):
    return matrix(ZZ,n,n,lambda i,j: b.monomial_coefficient(R(x^((j-i) % n))))

def vec_to_poly(z):
    tmp = R(0)
    for i in range(len(z)):
        tmp += R(z[i]*w**i)
    return tmp

def check_g(h,f):
    return vectorize((h*Rqn(f)).lift().change_ring(ZZ),n).norm(1) == dg

# ====================
# Gentry's attack
# ====================

# d2次元に圧縮されたfを探す
def CS2(h):
    Id = identity_matrix(d2)
    Ze = zero_matrix(d2)
    HH = circulant_matrix(Rqd2(h).lift().change_ring(ZZ),d2)
    return block_matrix(ZZ,2,2,[Id,HH,Ze,q*Id])

def find_fv2(h,debug=false):
    L2 = CS2(h)
    L2 = L2.BKZ()

    fv2_candidates = []
    for i in range(d1):
        fv2 = L2[i][0:d2]
        if fv2.norm(1) == df:
            if min(fv2.list()) >= 0:
                fv2_candidates.append(fv2)
            if max(fv2.list()) <= 0:
                fv2_candidates.append(-fv2)
    if len(fv2_candidates) == 0:
        raise ValueError('Oops, I cannot find any fv2_candidates')
    else:
        return fv2_candidates

# d2次元に圧縮されたヒントからd1次元に復元する

def CS2to1(h,fd):
    HH = circulant_matrix(Rqd1(h).lift().change_ring(ZZ),d1)
    t = vector(ZZ,list(zero_vector(d2))+list(fd)) * HH
    Hnew = block_matrix(ZZ,[[identity_matrix(d2),-identity_matrix(d2)]]) * HH
    CS1  = block_matrix([\
                  [identity_matrix(d2),Hnew],\
                  [zero_vector(d2).row(),t.row()],\
                  [zero_matrix(d1,d2),q*identity_matrix(d1)]\
                 ])
    return CS1[0:(d1+1),0:d1]

def find_fv1(h,fv2,debug=false):
    L1 = CS2to1(h,fv2)
    L1 = L1.BKZ()

    fv1_candidates = []
    # range(1,d1+1) because L1[0] = 0
    for i in range(1,L1.nrows()):
        fv1 = zero_vector(d1)
        fv1_left = L1[i][0:d2]
        if max(fv1_left.list()) <= 0:
            fv1_left = -fv1_left
        if min(fv1_left.list()) >= 0:
            fv1 = vector(ZZ,list(fv1_left)+list(fv2-fv1_left))
        if fv1.norm(1) == df:
            fv1_candidates.append(fv1)
    if len(fv1_candidates) == 0:
        raise ValueError('Oops, I cannot find any fv1_candidates')
    else:
        return fv1_candidates

# d1次元に圧縮されたヒントからn次元に復元する

def CS1to0(h,fd):
    HH = circulant_matrix(Rqn(h).lift().change_ring(ZZ),n)
    t = vector(ZZ,list(zero_vector(d1))+list(fd))*HH;
    Hnew = block_matrix(ZZ,[[identity_matrix(d1),-identity_matrix(d1)]]) * HH
    CS0  = block_matrix([\
                  [identity_matrix(d1),Hnew],\
                  [zero_vector(d1).row(),t.row()],\
                  [zero_matrix(n,d1),q*identity_matrix(n)]\
                 ])
    return CS0[0:(n+1),0:n]

def get_supports_of_0(z):
    s0 = []
    for i in range(len(z)):
        if z[i] == 0:
            s0.append(i)
    return s0

def find_fv0(h,fv1,debug=false):
    supports_of_0 = get_supports_of_0(fv1)
    # Gentryの勧めに従って, fv1の係数が0の場合, その列ベクトルは使わないので, 消す
    L0 = CS1to0(h,fv1).delete_rows(supports_of_0)
    # デフォルトはblock_size=10だが遅くなるので, 5に変更. 精度が必要らしいのでlong doubleを使うよう設定
    L0 = L0.BKZ(block_size=5,fp='ld')

    fv0_candidates = []
    for i in range(1,L0.nrows()):
        fv0 = zero_vector(n)
        fv0_left = L0[i][0:d1]
        if max(fv0_left.list()) <= 0:
            fv0_left = -fv0_left
        if min(fv0_left.list()) >= 0:
            fv0 = vector(ZZ,list(fv0_left)+list(fv1-fv0_left))
            if fv0.norm(1) == df:
                fv0_candidates.append(fv0)
    if len(fv0_candidates) == 0:
        raise ValueError('Oops, I cannot find any fv0_candidates')
    else:
        return fv0_candidates

# 攻撃本体
def test(seed,debug=false):
    set_random_seed(seed)
    f,g = skgen()
    h = pkgen(f,g)

    # Find fv2
    if debug:
        print "start find_fv2"
    try:
        fv2_candidates = find_fv2(h)
    except ValueError as err:
        print("Value Error: {0}".format(err))
    if debug:
        print "end find_fv2"
    fv2 = fv2_candidates[0]

    # Find fv1
    if debug:
        print "start find_fv1"
    try:
        fv1_candidates = find_fv1(h,fv2)
    except ValueError as err:
        print("Value Error: {0}".format(err))
    if debug:
        print "end find_fv1"
    fv1 = fv1_candidates[0]

    # Find fv0
    if debug:
        print "start find_fv0"
    try:
        fv0_candidates = find_fv0(h,fv1)
    except ValueError as err:
        print("Value Error: {0}".format(err))
    if debug:
        print "end find_fv0"

    # f0s is a list of polynomials f0's in R = ZZ[x]
    f0s = [vec_to_poly(fv0_candidates[i]) for i in range(len(fv0_candidates))]
    f0s = [f0 for f0 in f0s if check_g(h,f0)]
    return f,g,h,f0s



print "n,p,q", n,p,q
print "df,dg,dr", df, dg, dr
print "d1,d2", d1,d2

# %time ... で時間を計測
%time f,g,h,f0s = test(0,true)

for i in range(len(f0s)):
    print "i:", i
    print "fi:", f0s[i]
    print "gi:", (h*Rqn(f0s[i])).lift().change_ring(ZZ)
print "!!!! End !!!!"

n,p,q 256 2 127
df,dg,dr 35 35 22
d1,d2 128 64
start find_fv2
end find_fv2
start find_fv1
end find_fv1
start find_fv0
end find_fv0
CPU times: user 42.3 s, sys: 187 ms, total: 42.5 s
Wall time: 43.8 s
i: 0
fi: x^251 + x^234 + x^229 + (snip) + x^20 + x^7 + 1
gi: x^254 + x^253 + x^247 + (snip) + x^25 + x^23 + x^6
i: 1
fi: x^253 + x^235 + x^233 + (snip) + x^12 + x^10 + 1
gi: x^248 + x^238 + x^234 + (snip) + x^23 + x^12 + x^7
!!!! End !!!!

無事にfi, giが復元できました.
f0s[0]とf0s[1]で復号できるかどうか確かめてみます.

# ランダムな平文を作って, 暗号化して復号する
m = random_tpoly(int(n/2))
c = enc(h,m)
print "m:", m
for i in range(len(f0s)):
    d = dec(f0s[i],c)
    print "i:", i
    print "d:", d
    print "flag,", d == m

m: x^247 + x^246 + x^228 + (snip) + x^9 + x^5 + x^4
i: 0
d: x^247 + x^246 + x^228 + (snip) + x^9 + x^5 + x^4
flag, True
i: 1
d: x^247 + x^246 + x^228 + (snip) + x^9 + x^5 + x^4
flag, True

本来の秘密鍵とはずれていますが, ちゃんと復号出来ることがわかります.

$h = g/f$なので, 任意の$i$で$h = (g x^i)/(f x^i)$が成立します. よって $f x^i$も秘密鍵になっています. 今回はこれを二つ引いてきたことがわかります.

for i in range(len(f0s)):
    print Rqn(f0s[i])/Rqn(f)

実行結果

xbar^9
xbar^137

まとめ

• 多項式環とその剰余環, 巡回行列をざっと振り返りました
• 耐量子計算機暗号の一角であるNTRU暗号をsagemath8.0上で実装しました.
• CoppersmithとShamirの格子攻撃を振り返りました.
• $n$が合成数のNTRU暗号に対するGentry攻撃 (EUROCRYPT2001) を振り返り, sagemath8.0上でGentry攻撃を実装しました.

sagemathは実際の暗号方式を攻撃する際にも非常に使えるツールです.