CTFのためのPARI/GP

PARI/GPは整数論的計算を行うことに特化した計算環境である。

PARI

計算ライブラリ

GP

スクリプト言語

言語についてちゃんと学ぶには、PARI/GP documentationの"User's Guide to PARI/GP"と"GP Tutorial"を参照。 関数については"Online User's Guide"を参照。 GP Reference Cardは期待するほどには役に立たない。

Syntax of GP

ヘルプ

? funcname

コメント

\\ コメント

なかなか珍しいバックスラッシュによるコメントアウト。

基本演算

• 加算 +

• 減算 -

• 乗算 *

• 冪算 ^

• 除算 /

  • 整数に対しては分数になる

• 少数切り捨て除算 \

• 剰余 %

• 四捨五入 \/

• ビットシフト <<, >>

• 比較 <, >, <=, >=, !=(<>), ==

  • 注意: [] == 0

• 厳密な比較 ===

• 否定 前置!

• 論理演算 &(&&), |(||)

  • ビット演算はbit*関数群を使用

ベクトル(aka タプル・配列)

• 行ベクトル: [1, 2, 3]

• 列ベクトル: [1, 2, 3]~

• 参照: vector[i] (1-origin)

• 要素数: #vector

? [2, 4, 6][1]
%1 = 2
? #[2, 4, 6]
%2 = 3

行列

• 行列: [1, 2; 3, 4]

• 参照: matrix[i,j]

  • 行参照: matrix[i,]

  • 列参照: matrix[,j]

? [1, 2; 3, 4]
%1 = 
[1 2]

[3 4]

? [1, 2; 3, 4][1, 1]
%2 = 1
? [1, 2; 3, 4][1,]
%3 = [1, 2]
? [1, 2; 3, 4][,1]
%4 = [1, 3]~

初等整数論的演算 ($\mathbb{Z}/n\mathbb{Z}$)

素因数分解

• factor(x): xを素因数分解

• factor(x, lim): lim未満の素数を使用

  • 事前にlim以下の素数を計算しておく必要がある: default(primelimit, 5000000)

結果は行列として返される。形式: [素因数1, 素因数1の個数; 素因数2, 素因数2の個数; ...]

? factor(2^17+1)
%1 = 
[3 1]

[43691 1]

? default(primelimit, 10^6)
? factor(3*nextprime(10^6)*nextprime(10^7), 10^6)
%2 = 
[3 1]

[10000049000057 1]

制限時間付き素因数分解

Programming in GP: other specific functionsに次のコードが掲載されている。

default(factor_add_primes,1);
default(primelimit,16777216);
timefact(N,sec)=
{
  trap(alarmer,factor(N,0),alarm(sec);my(F=factor(N));alarm(0);F);
}

? timefact(2*nextprime(10^30)*nextprime(10^31), 1)
%2 = 
[2 1]

[10000000000000000000000000000603000000000000000000000000001881 1]

乱数生成

random(n)

nが整数であれば[0, n)の範囲の整数をランダムに返す。

大きい素数を生成

nextprime(n)を用いればnとほぼ同じ大きさの素数が得られる。

\\ 10^10の次の素数
? nextprime(10^10)
%1 = 10000000019
\\ 128bitの素数をランダムに生成
? nextprime(2^127 + random(2^127))
%2 = 267183150806405707118987831126661135943

$\mathbb{Z}/n\mathbb{Z}$上の点

$n$を法とした$x$: Mod(x, n)

? Mod(3, 7) * 3
%1 = Mod(2, 7)

注意: 単に剰余を計算するなら%演算子がある。

Mod(x, n)からxとnを取り出す

g = Mod(x, n)に対して

• xの取り出し lift(g)

• nの取り出し g.mod

? lift(Mod(2, 7))
%1 = 2
? Mod(2, 7).mod
%2 = 7

逆元

割り算するだけで良い。指数に-1を渡しても良い。

? 1/Mod(2, 7)
%1 = Mod(4, 7)
? Mod(2, 7)^-1
%2 = Mod(4, 7)

$\mathbb{Z}/p\mathbb{Z}$上の離散n乗根

sqrtn(x, n, &z)

存在すればzにxのn乗根の生成元が入る。存在しなければzは0になる。関数自体は主値(最小の元？)を返す。

? sqrtn(Mod(1, 7), 3, &z)
%1 = Mod(1, 7)
? z
%2 = Mod(2, 7)
? [z^0, z^1, z^2]
%3 = [Mod(1, 7), Mod(2, 7), Mod(4, 7)]

User's Guide to PARI/GPに全ての根を得るための方法として次のコードが掲載されている。

sqrtnall(x,n)=
{ my(V,r,z,r2);
  r = sqrtn(x,n, &z);
  if (!z, error("Impossible case in sqrtn"));
  if (type(x) == "t_INTMOD" || type(x)=="t_PADIC" ,
  r2 = r*z; n = 1;
  while (r2!=r, r2*=z;n++));
  V = vector(n); V[1] = r;
  for(i=2, n, V[i] = V[i-1]*z);
  V
}
addhelp(sqrtnall,"sqrtnall(x,n):compute the vector of nth-roots of x");

? sqrtnall(Mod(1, 7), 3)
%3 = [Mod(1, 7), Mod(2, 7), Mod(4, 7)]

$\mathbb{Z}/pq\mathbb{Z}$上の非自明な平方根

pqを法としたxの平方根rは $$ \left\{\begin{array}{l} r\equiv\pm\sqrt{x} \pmod{p} \\ r\equiv\pm\sqrt{x} \pmod{q} \end{array} \right. $$ より、中国人剰余定理を用いて求めることができる。

? chinese(sqrt(Mod(4, 3)), sqrt(Mod(4, 5)))
%1 = Mod(7, 15)
? chinese(sqrt(Mod(4, 3)), -sqrt(Mod(4, 5)))
%2 = Mod(13, 15)
? chinese(-sqrt(Mod(4, 3)), sqrt(Mod(4, 5)))
%3 = Mod(2, 15)
? chinese(-sqrt(Mod(4, 3)), -sqrt(Mod(4, 5)))
%4 = Mod(8, 15)

オイラーのφ関数

eulerphi(n)

n以下の正整数のうち、nと互いに素なものの個数を返す。

カーマイケル関数

znstar(n)[2][1]

$ k ^ x \equiv 1 \pmod{n}$が任意のk(!=0)に成立するような最小のxを返す。

乗法群の位数

znorder(g)

gを生成元とする$\mathbb{Z}/n\mathbb{Z}$の部分巡回群の位数を返す。

離散対数

znlog(x, g, znorder(g))

gを底としたxの離散対数を求める。 gは$(\mathbb{Z}/n\mathbb{Z})^*$の元である必要がある。

第三引数は$g$が$(\mathbb{Z}/n\mathbb{Z})^*$の生成元である場合に限り省略できる。 そうでない時に省略した場合、エラーとなることがある。

\\ 2^4 = 3 (mod 13)
? znlog(3, Mod(2, 13))
%1 = 4
\\ 4^2 = 3 (mod 13)
? g = Mod(4, 13); znlog(3, g, znorder(g))
%2 = 2

生成元

znprimroot(n)

Z/nZの生成元を返す。

? znprimroot(13)
%1 = Mod(2, 13)

最大公約数

• gcd(x, y)

• gcd([x1, x2, ...])

2引数で渡しても、配列で渡しても良い。

? gcd(4, 6)
%1 = 2
? gcd(gcd(12, 18), 30)
%2 = 6
? gcd([12, 18, 30])
%3 = 6

拡張ユークリッドの互除法

bezout(x, y)

$xu + yv = d$ なる[u, v, d]を返す。

\\ gcd(4, 6) = 2
\\ 4*-1 + 6*1 = 2
? bezout(4, 6)
%1 = [-1, 1, 2]

中国人剰余定理

• chinese(x, y): x, yは$\mathbb{Z}/n\mathbb{Z}$の元。

• chinese([x1, x2, ...]): x1, x2, ...は$\mathbb{Z}/n\mathbb{Z}$の元。

2引数で渡しても、配列で渡しても良い。

? chinese(chinese(Mod(2, 3), Mod(3, 5)), Mod(2, 7))
%1 = Mod(23, 105)
? chinese([Mod(2, 3), Mod(3, 5), Mod(2, 7)])
%2 = Mod(23, 105)

RSA

正当な計算例

\\ 鍵生成
p = 13
q = 17
e = 7
n = p*q
phi = (p-1)*(q-1)
d = lift(1/Mod(e, phi))
pk = [n, e]
sk = [n, d]

\\ 平文
m = Mod(8, n)
\\ 暗号化
c = m^e
\\ 復号
m = c^d

攻撃例: 離散対数問題を解ける場合の既知平文攻撃(KPA)

d = znlog(m, c, znorder(c))

ElGamal暗号

正当な計算例

\\ 鍵生成
\\ p=2*q+1 is prime
\\ #<g> = q
q = 29
p = 2*q+1
g = znprimroot(p)^2
x = random(q)
h = g^x
pk = [p, q, g, h]
sk = [p, x]
\\ 平文
\\ m \in <g> に注意
m = g ^ 2
\\ 暗号化
r = random(q)
c1 = g^r
c2 = m * h^r
c = [c1, c2]
\\ 復号
m = c2 / c1^x

q, pのランダムな生成

q = nextprime(2^127 + random(2^127))
p = 1; until(isprime(p), q = nextprime(q+2); p = 2*q+1)

攻撃例: 離散対数問題を解ける場合の唯鍵攻撃(KOA)

x = znlog(h, g, znorder(g))

攻撃例: 離散対数問題を解ける場合の暗号文単独攻撃(COA)

KOAできるので意味なし

r = znlog(c1, g, znorder(g))
m = c2 / h^r

楕円曲線

楕円曲線の生成

ellinit([a1, a2, a3, a4, a6])

楕円曲線 $y^2+a_1xy+a_3y=x^3+a_2x^2+a_4x+a_6$ を表す構造を生成する。 すなわちワイエルシュトラスの標準形のパラメータである。

そのままでは有理数上の楕円曲線となる。

? E = ellinit([0,0,0,-1,1])
%1 = [0, 0, 0, -1, 1, 0, -2, 4, -1, 48, -864, -368, -6912/23, [-1.3247179572447460259609088544780973407, 0.66235897862237301298045442723904867037 - 0.56227951206230124389918214490937306150*I, 0.66235897862237301298045442723904867037 + 0.56227951206230124389918214490937306150*I]~, 4.7070877612301855618837521175707435935, 2.3535438806150927809418760587853717967 - 1.0982915250610051220258220794831228014*I, 2.4199001261583493071188336079560183889 - 1.1754943508222875080 E-38*I, 1.2099500630791746535594168039780091945 + 0.77020648244262125847219661990494373329*I, 5.1697545958774928400543891205186285388]

素体$\mathbf{F}_p$上の楕円曲線の生成

ellinit([a1, a2, a3, a4, a6] * Mod(1, p))

$\mathbb{Z}/p\mathbb{Z}$ 上の楕円曲線 $y^2+a_1xy+a_3y=x^3+a_2x^2+a_4x+a_6$ を表す構造を生成する。

? E = ellinit([0,0,0,-1,1] * Mod(1, 7))
%1 = [Mod(0, 7), Mod(0, 7), Mod(0, 7), Mod(6, 7), Mod(1, 7), Mod(0, 7), Mod(5, 7), Mod(4, 7), Mod(6, 7), Mod(6, 7), Mod(4, 7), Mod(3, 7), Mod(2, 7), 0, 0, 0, 0, 0, 0]

ガロア体$\mathbf{F}_{p^f}$上の楕円曲線の生成

(version >= 2.7)

ellinit([a1, a2, a3, a4, a6], t)

tの属するガロア体上の楕円曲線 $y^2+a_1xy+a_3y=x^3+a_2x^2+a_4x+a_6$ を表す構造を生成する。

? t = ffgen(2^4);
? ellinit([1, 0, 0, 0, 1], t)
%1 = [1, 0, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, Vecsmall([4]), [x, [Vecsmall([0]), Vecsmall([0, 1]), [Vecsmall([0, 1]), Vecsmall([0]), Vecsmall([0]), Vecsmall([0])]]], [0, 0, 0, 0]]

点の表現

• 通常の点: [x, y]

• 無限遠点: [0]

群上の基本演算

• 加算 $P_1+P_2$ : elladd(E, P1, P2)

• 減算 $P_1-P_2$ : ellsub(E, P1, P2)

• $n$倍算(冪算) $[n]P$ : ellpow(E, P, n)

x座標として$x$を持つ点

ellordinate(E, x)

楕円曲線$E$上の、x座標として$x$を持つ点のy座標の一覧を返す。

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? ellordinate(E, 1)
%2 = [Mod(1, 7), Mod(6, 7)]
? ellordinate(E, 2)
%3 = [Mod(0, 7)]
? ellordinate(E, 4)
%4 = []

楕円曲線上の点かどうかの判定

ellisoncurve(E, z)

点$z$が楕円曲線$E$上の点なら1、そうでないなら0を返す。

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? ellisoncurve(E, [0, 1])
%2 = 1
? ellisoncurve(E, [0, 0])
%3 = 0

ランダムな点

random(E)

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? random(E)
%2 = [Mod(5, 7), Mod(4, 7)]

群の位数

ellorder(E, P)

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? ellorder(E, [0])
%2 = 1
? ellorder(E, [1, 1])
%3 = 12

群上の離散対数

elllog(E, P, G, ellorder(E, G))

$G$を底とした$P$の離散対数を求める。

第四引数は$G$が$E$の全ての点の生成元である場合に限り省略できる。 そうでない時に省略した場合、エラーとなることがある。

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? G = [6, 1]; P = ellpow(E, G, 3);
? elllog(E, P, G, ellorder(E, G))
%73 = 3

濃度($\mathbb{Z}/p\mathbb{Z}$上の楕円曲線上の点の数)

p+1-ellap(E, p)

ellap(E, p)は$E$のFrobenius写像のトレースを計算する。

? E = ellinit([0,0,0,-1,1] * Mod(1, 7));
? 7+1-ellap(E, 7)
%2 = 12

楕円曲線上のElGamal暗号

正当な計算例

\\ 鍵生成
p = 7
E = ellinit([0,0,0,-1,1] * Mod(1, p))
g = [1, 1]
ord = ellorder(E, g)

x = random(ord)
h = ellpow(E, g, x)
pk = [p, g, h]
sk = [p, x]
\\ 平文
\\ m \in E(F_p) に注意
m = [3, ellordinate(E, 3)[1]]
\\ 暗号化
r = random(ord)
c1 = ellpow(E, g, r)
c2 = elladd(E, m, ellpow(E, h, r))
c = [c1, c2]
\\ 復号
m = ellsub(E, c2, ellpow(E, c1, x))

攻撃例: 離散対数問題を解ける場合の唯鍵攻撃(KOA)

x = elllog(E, h, g, ellorder(E, g))

多項式

Pari/GPにおいては文字式はt_POL型となり多項式として扱われる。 また、t_POLMOD型として多項式を法とする多項式環も扱うことができる。

多項式環上の多項式

Mod(x, y)

多項式$y$を法とした多項式$x$

? Mod(x, x^2 + 1) * x
%1 = Mod(-1, x^2 + 1)

注意: 単に剰余を計算するなら%演算子がある。

Mod(x, y)からxとyを取り出す

普通のModとほぼ同じ。

z = Mod(x, y)に対して

• xの取り出し z.pol または lift(z)

• yの取り出し z.mod

? Mod(x, x^2 + 1).pol
%1 = x
? lift(Mod(x, x^2 + 1))
%2 = x
? Mod(x, x^2 + 1).mod
%3 = x^2 + 1

多項式の特定の係数を取り出す

polcoeff(x, n, v)

多項式$x$の$v^n$の係数。 vを省略した場合にはvariable(x)で得られる多項式xの主要な変数が選択される。

? polcoeff(3*x^2, 2)
%1 = 3
? polcoeff(3*x^2, 1)
%2 = 0
? polcoeff(x*y, 1)
%3 = y
? polcoeff(x*y, 1, y)
%4 = x

see also: 多項式の係数をベクトルで取り出す

多項式の係数をベクトルで取り出す

Vec(pol)

多項式polの係数をベクトルで次数が高い順に返す。

? Vec(2*x^2+1)
%1 = [2, 0, 1]

なおVecrev(pol)とすると次数が低い順になる。

与えた係数を持つ多項式

Pol(vec)

ベクトルvecの内容を次数が高い順に係数に持つような多項式を返す。

? Pol([2, 0, 1])
%1 = 2*x^2 + 1

多項式の次数

poldegree(x, v)

多項式xの次数を返す。 vを省略した時はx主要な変数が選択される。

? poldegree(x^2+1)
%1 = 2
? poldegree(x^2+1, x)
%2 = 2

変数の置き換える

subst(x, y, z)

多項式xの変数yをzに置き換える。

? subst(x^2 + x, x, y)
%1 = y^2 + y
? subst(x*y, x, y)
%2 = y^2
? subst(x^2 + x, x, 2)
%3 = 6

式の評価

eval(x)

x中の変数を定義済み変数の内容で置き換えて評価する。

? x = y^2; y = 2; eval(x)
%1 = 4

xは文字列でもよい。

? x = y^2; y = 2;
? x
%1 = y^2
? eval("x^2")
%2 = y^4
? eval(eval("x^2"))
%3 = 16

注意: substの利用を検討すべき

ガロア体(有限体)

多項式環を用いて$\mathbf{F}_p$の拡大体$\mathbf{F}_{p^f}$を定義する。

既約多項式の生成

ffinit(p, n, v)

$\mathbf{F}_p[v]$上の$n$次の既約多項式を生成する。vを省略した場合はxが用いられる。

? ffinit(2, 4)
%1 = Mod(1, 2)*x^4 + Mod(1, 2)*x^3 + Mod(1, 2)*x^2 + Mod(1, 2)*x + Mod(1, 2)
? ffinit(2, 4, 'y)
%2 = Mod(1, 2)*y^4 + Mod(1, 2)*y^3 + Mod(1, 2)*y^2 + Mod(1, 2)*y + Mod(1, 2)

注意: ガロア体を利用するにあたって生成多項式を指定したい場合は用いない。

位数を指定したガロア体の生成

ffgen(q, v)

位数qのガロア体の生成元を返す。 生成元を表す変数名としてvに指定した変数名が用いられる。vを省略した場合はxが用いられる。 戻り値はt_FFENT型の値となる。

? g=ffgen(2^3, 't)
%1 = t
? g^3
%2 = t^2 + 1

注意1: 生成多項式を指定したい場合は次項を参照。

注意2: ここで言う生成元とは乗法的生成元ではない。つまり、ffgen(q)^i (i=1,...,q-1)によりガロア体上の0を除く全ての元が生成できるとは限らない。

生成多項式を指定したガロア体の生成

ffgen(x, v)

有限体$\mathbf{F}_p$上の多項式xを生成多項式としたガロア体の生成元を返す。 生成元を表す変数名としてvに指定した変数名が用いられる。vを省略した場合はxが用いられる。 戻り値はt_FFENT型の値となる。

? g=ffgen((x^3+x^2+1)*Mod(1,2), 't)
%1 = t
? g^3
%2 = t^2 + 1

注意: ここで言う生成元とは乗法的生成元ではない。つまり、ffgen(x)^i (i=0,1,...)によりガロア体上の0を除く全ての元が生成できるとは限らない。

ガロア体のパラメータの取り出し

xをガロア体上の元として

• xを表す多項式(t_POL型): x.pol

• 生成多項式: x.mod
• 生成多項式の$\mathbf{F}_p$の$p$: x.p

• 拡大次数($\mathbf{F}_{p^f}$の$f$): x.f

? g = ffgen(2^3);
? x = g^3
%1 = x^2 + 1
? x.pol
%2 = x^2 + 1
? type(x)
%3 = "t_FFELT"
? type(x.pol)
%4 = "t_POL"
? x.mod
%5 = x^3 + x^2 + 1
? x.p
%6 = 2
? x.f
%7 = 3

ガロア体上の離散対数

fflog(x, g, fforder(g))

gを底としたxの離散対数を求める。

第三引数を省略するとgが原始根であると仮定する。解が存在しない場合の動作は未定義である。

? g = ffgen(2^4)
%1 = x
? fflog(g^4, g, fforder(g))
%2 = 4

乗法的位数

fforder(g)

乗法的位数、つまり$g^e=1$なる正で最小の整数を求める。

? g = ffgen(2^4)
%1 = x
? fforder(g)
%2 = 5
? fforder(g+1)
%3 = 5

原子根(乗法的生成元)

ffprimroot(x, &o)

xの属するガロア体の乗法的生成元を得る。 同時にoを指定すると、乗法的位数に関する情報をoに格納する。 oはfforderに代えて用いることが出来る。

? t = ffgen(2^4);
? g = ffprimroot(t, &o)
%1 = x^2 + x + 1
? fforder(g, o)
%2 = 15

多項式と数値の一般的な対応について

ガロア体$\mathbf{F}_{p^f}$上の多項式$\sum_i a_i x^i$は、 一般に$x$に$p$を代入した整数と同一視される。

次のような関数を定義しておくと良い。

int2ff(n, g) = { subst(Pol(digits(n, g.p)), 'x, g) };
ff2int(x) = { subst(x.pol, variable(x.pol), x.p) };

? t = ffgen(2^4);
? int2ff(10, t)
%1 = x^3 + x
? ff2int(t^3 + t)
%2 = 10

ガロア体上の楕円曲線上のECDSA

正当な計算例

int2ff(n, g) = { subst(Pol(digits(n, g.p)), 'x, g) };
ff2int(x) = { subst(x.pol, variable(x.pol), x.p) };

\\ 楕円曲線 K283
t = ffgen((t^283 + t^12 + t^7 + t^5 + 1) * Mod(1, 2), 't)
E = ellinit([1,0,0,0,1], t)
G_x = 9737095673315832344313391497449387731784428326114441977662399932694280557468376967222;
G_y = 3497201781826516614681192670485202061196189998012192335594744939847890291586353668697;
G = [int2ff(G_x, t), int2ff(G_y, t)]
n = ellorder(E, G)

\\ 鍵生成
d_A = random(n)
Q_A = ellpow(E, G, d_A)
\\ メッセージのハッシュ
z = 1234567
\\ 署名
k = random(n)
kG = ellpow(E, G, k)
r = Mod(ff2int(kG[1]), n)
s = Mod(k, n)^-1 * (z + r * d_A)
sig = lift([r, s])
\\ 検証
w = s^-1
u1 = z*w
u2 = r*w
r == ff2int(elladd(E, ellpow(E, G, lift(u1)), ellpow(E, Q_A, lift(u2)))[1])

攻撃例: 同じkを使いまわした場合の鍵回復

k = (z1 - z2) / (s1 - s2)
d_A = (s1*k - z1) / r

格子

ペアリング