Ed25519

# はじめに（Introduction）

Ed25519とはエドワーズ曲線デジタル署名アルゴリズム（Edwards-curve Digital Signature Algorithm）を実装したものです。
楕円曲線暗号の一つで、他のものより**高速**で**コンパクト**である事が特徴です。

ここでは、仕様をみながら基本的な実装していきます。

# Ed25519

## EdDSAパラメータ(EdDSA parameters)

EdDSAには7つのパラメーターがあります。
Ed25519のパラメータも共に記載します。

|EdDSA|Ed25519|
|:--|:--|
|整数 $b \geq 10$ |$b = 256$|
|$2b$ビット出力する暗号化ハッシュ関数 $\text{H}$ | $\text{H}$ は **SHA-512** |
|$1 \pmod{4}$ となる素数 $q$| $q$ は素数 $2^{255} − 19$ |
|有限体 $F_q$ 要素の $(b-1)$ ビットエンコーディング|エンコーディングはリトルエンディアン |
|$F_q$ の非平方要素 $d$ |$d = −\frac{121665}{121666} \in F_q$|
|**※**を満たす $2^{b-4}$ と $2^{b-3}$ の間の素数 $\ell$|$\ell$は素数 $2^{252} + 27742317777372353535851937790883648493$|
|集合の要素 $B\ne (0,1)$|$x$ が正の一意の点 $(x, \frac{4}{5}) \in E$|

以下の楕円曲線を使います。

```math
E = \{ (x, y) \in F_q \times F_q : −x^2 + y^2 = 1 + dx^2y^2 \}
```

楕円曲線上の点に対して、以下の加算則を定めると、**中立元**（neutral element）$0=(0,1)$ で**群**をなします。

```math
(x_1, y_1) + (x_2, y_2) = \left( \frac{x_1y_2 + x_2y_1}{1 + dx_1x_2y_1y_2}\ ,\ \frac{y_1y_2 + x_1x_2}{1 − dx_1x_2y_1y_2} \right)
```

パラメータの**※**は、$\ell B=0$ を満たす事です。（$nB$ は群における $B$ の $n$ 倍を意味します。）

要素 $(x,y)\in E$ を、$256$ビット（$b$ビット）でエンコードします。
$y$ は $255$ビット（$(b-1)$ビット）でエンコードされ、最後のビットは符号ビットです。
$x$ が負の場合、符号ビットは $1$ です。

```python
import hashlib


class Ed25519P():
    p = 2**255 - 19
    L = 2**252 + 27742317777372353535851937790883648493
    d = (-121665 * pow(121666, p-2, p)) % p

    def __init__(self, x, y):
        if x is None or y is None:
            raise ValueError("x or y are None")
        self.x = x % self.p
        self.y = y % self.p
        # : −x^2 + y^2 = 1 + dx^2y^2
        if (-1 * pow(self.x, 2, self.p) + pow(self.y, 2, self.p)) % self.p != (1 + self.d * pow(self.x, 2, self.p) * pow(self.y, 2, self.p)) % self.p:
            raise ValueError("x and y are not on the curve")

    def __add__(self, other):
        x1 = self.x
        y1 = self.y
        x2 = other.x
        y2 = other.y
        #        x1y2 + x2y1
        # x3 = ---------------
        #       1 + dx1x2y1y2
        x3 = ((x1*y2 + x2*y1) * pow((1 + self.d*x1*x2*y1*y2) %
                                    self.p, self.p-2, self.p)) % self.p
        #        y1y2 + x1x2
        # y3 = ---------------
        #       1 − dx1x2y1y2
        y3 = ((y1*y2 + x1*x2) * pow((1 - self.d*x1*x2*y1*y2) %
                                    self.p, self.p-2, self.p)) % self.p
        return self.__class__(x3, y3)

    def __sub__(self, other):
        return self + other.__class__(self.p - self.x, self.y)

    def __rmul__(self, coefficient):
        coef = coefficient % self.L
        current = self.__class__(self.x, self.y)
        result = self.__class__(0, 1)
        while coef:
            if coef & 1:
                result += current
            current += current
            coef >>= 1
        return result

    def __eq__(self, other):
        return (self.x == other.x) and (self.y == other.y)

    def __str__(self):
        return '({:064x},{:064x})'.format(self.x, self.y)

    def encode(self):
        return int.to_bytes(self.y | ((self.x & 1) << 255), 32, 'little')
```

> コードについて、いくつか解説します。
> 除算に関しては、乗法逆元の乗算としています。（参照：付録（Appendix）- 乗法逆元（multiplicative inverse））
> 減算（`__sub__`）の部分は、点を反転（正負を逆転）して加算しています。
> 点の反転は、点 $(x_1,y_1)$ の $x_1$ を反転したもの 点$(-x_1,y_1)$ です。（参照：付録（Appendix）- 符号反転（negation））
> スカラー倍算（`__rmul__`）の部分は二進展開を使います。（参照：付録（Appendix）- 二進展開（binary expansion））

$y$ と 符号ビットから、以下の式で $x$ を求める事ができます。
$x=\pm\sqrt{(y^2-1)/(dy^2+1)}$

```python
def DecodePoint(encode_point):
    if len(encode_point) != 32:
        raise ValueError("invalid length of encode point")
    y = int.from_bytes(encode_point, 'little')
    sign = y >> 255
    y &= (1 << 255) - 1
    # x^2 = (y^2−1)/(dy^2+1)
    x2 = ((y*y-1)
          * pow(Ed25519P.d*y*y+1, Ed25519P.p-2, Ed25519P.p)) % Ed25519P.p
    if x2 == 0:
        if sign:
            raise ValueError("invalid sign of encode point")
    # x = ±√(y^2−1)/(dy^2+1)
    x = pow(x2, (Ed25519P.p+3) // 8, Ed25519P.p)
    if (x*x - x2) % Ed25519P.p != 0:
        x = (x * pow(2, (Ed25519P.p-1) // 4, Ed25519P.p)) % Ed25519P.p
        if (x*x - x2) % Ed25519P.p != 0:
            raise ValueError("invalid square of encode point")
    if (x & 1) != sign:
        x = Ed25519P.p - x
    return Ed25519P(x, y)

By = (4*pow(5, Ed25519P.p-2, Ed25519P.p)) % Ed25519P.p
B = DecodePoint(By.to_bytes(32, 'little'))
```


> コードについて、いくつか解説します。
> パラメータの長さ $32$ バイト（$256$ビット）をチェックします。
> リトルエンディアンで数値化した後、$y$ と符号（偶数or奇数）を取得します。
> $y$ から、$x^2$ を計算し、平方根を求めます。（参照：付録（Appendix）- 平方根（square root））
> 符号（偶数or奇数）から、値を確定します。

## EdDSAの鍵と署名検証(EdDSA keys and signatures)

### 秘密鍵と公開鍵（secret key and public key）

EdDSA __秘密鍵__は $32$ バイト（$256$ビット）のデータ __$k$__ です。
これをハッシュした値 $\text{H}(k)=(h_0,h_1,\ldots ,h_{511})$ から数値を求めます。

```math
a = 2^{254}+\sum_{3 \leq i \leq 253}2^ih_i\in\{ 2^{254},2^{254}+8,\ldots ,2^{255}-8 \}
```

> 捕捉します。
> $254$ビット目を $1$ 、最初の $3$ビットを $0$ としたものをリトルエンディアンで数値にしたものです。
>
|0|1|2|3|4|$\cdots$|253|254|255|
|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|
|$0$|$0$|$0$|$h_3$|$h_4$|$\cdots$|$h_{253}$|$1$|$0$|
>
> したがって、最小値は $2^{254}$ 、下位 $3$ビットは $0$ なので、$8$ 刻みの数値となり、最大値は $2^{255}-8$ となります。

$a$ から、$A=aB$ を計算します。
EdDSA __公開鍵__は $A$ です。

### 署名（signature）

__秘密鍵 $k$__によるメッセージ $M$ の署名は、次のように定義されます。
秘密鍵 $k$ から求めたハッシュ値の後半部分とメッセージ $M$ を使って以下を求めます。
$r = \text{H}(h_{256},\ldots ,h_{511},M) \in \\{0,1,\ldots , 2^{512}-1 \\}$ 
ハッシュ値をリトルエンディアンで数値にして、$\\{0,1,\ldots , 2^{512}-1 \\}$ の整数とます。
$r$ から、$R=rB$ を計算します。
$a,A,r,R,M$ を使って以下を求めます。
$S = (r + \text{H}(R, A, M)a) \mod \ell$

秘密鍵 $k$ によるメッセージ $M$ の署名は、$64$ バイト（$512$ビット）の $(R,S)$ となります。
$S$ は、$32$ バイト（$256$ビット）リトルエンディアンエンコーディングです。

```python
def Sign(secret, m):
    if len(secret) != 32:
        raise ValueError("invalid length of secret")
    h = hashlib.sha512(secret).digest()
    b = (h[0] & 0xf8).to_bytes(1, byteorder='little') + \
        h[1:31] + \
        ((h[31] & 0x7f) | 0x40).to_bytes(1, byteorder='little')
    s = int.from_bytes(b, 'little')
    A = s*B
    prefix = h[32:]
    r = int.from_bytes(hashlib.sha512(prefix + m).digest(), 'little')
    R = r*B
    k = int.from_bytes(hashlib.sha512(
        R.encode() + A.encode() + m).digest(), 'little')
    S = (r + k * s) % Ed25519P.L
    return R.encode() + S.to_bytes(32, 'little')
```

### 検証（verication）

__公開鍵 $A$__によるメッセージ $M$ と署名 $(R,S)$ の検証は、次のように行います。
$A$ は $A \in E$ として、$(R,S)$ は解析して $R \in E$ と $S \in \\{0,1,\ldots , \ell-1\\}$ にします。
検証者は $\text{H}(R,A,M)$ を計算し、以下をチェックします。

```math
8SB \overset{?}{=}8R + 8\text{H}(R, A, M)A
```

```python
def Verify(A, signature, m):
    if len(signature) != 64:
        raise ValueError("invalid length of signature")
    R = DecodePoint(signature[:32])
    S = int.from_bytes(signature[32:], 'little')
    H = int.from_bytes(hashlib.sha512(
        R.encode() + A.encode() + m).digest(), 'little')
    return 8*S*B == 8*R + 8*H*A
```

# おわりに（Conclusion）

単純な実装はできましたが、さらに高速に演算することができます。
また、公開鍵や秘密鍵を $32$ バイト単位にする事によって、CPUの演算的にも高速に演算可能となっています。
それにより、コンパクトにもなっています。

個人的には、署名がSchnorrなので署名の加算が可能である点、署名のロジックにナンス生成も含まれている事から、かなり使いやすい実装だと感じました。

https://github.com/tnakagawa/pyref/blob/master/ed25519/ed25519.py

# 参考（References）

- High-speed high-security signatures<br>http://ed25519.cr.yp.to/ed25519-20110926.pdf

- RFC 8032 - Edwards-Curve Digital Signature Algorithm (EdDSA)<br>https://tools.ietf.org/html/rfc8032

- Quadratic residue - Wikipedia<br>https://en.wikipedia.org/wiki/Quadratic_residue#Prime_or_prime_power_modulus

# 付録（Appendix）

## 乗法逆元（multiplicative inverse）

有限体における乗法の逆元、いわゆる $a$ に対して、$\frac{1}{a}$ となる値です。
フェルマーの小定理 $a^{p-1} \equiv 1 \pmod{p}$ より、

```math
a^{p-2} \equiv a^{-1} \pmod{p}
```

で求める事ができます。

## 二進展開（binary expansion）

ポイントの加算において、$nB$ のように $n$倍する場合、単純には $n$回 $B$ を加算すればいいのですが、
$n$ が大きな数の場合（$2^{252} \leq \ell \leq 2^{253}$）、加算に時間がかかってしまいます。
$n$ を二進数に展開し、$B$ を倍々にし $1$ のところだけを加算する事で $\log_2(n)$回ですます事ができます。

例えば、$100$を考えてみます。

```math
\begin{align*}
100_{10} &= 1100100_2 \\
&= 1000000_2 \\
&+ \ \ 100000_2 \\
&+ \ \ \ \ \ \ \ \ 100_2 \\
\end{align*}
```

通常は $100$ 回加算する必要がありますが、この方法を使えば
二倍の計算が $6$ 回、加算が $3$ 回の合計 $9$ 回で済ます事ができます。

## 符号反転（negation）

点 $(x_1,y_1)$ と 点$(-x_1,y_1)$ を加算すると、点 $(0,1)$ となるか見てみます。
点 $(x_1,y_1)$ が楕円曲線上の点である式と加算の式で求められます。

```math
\begin{align*}
−x_1^2 + y_1^2 &= 1 + dx_1^2y_1^2 \\
(x_1, y_1) + (-x_1, y_1) &= \left( \frac{x_1y_1 - x_1y_1}{1 - dx_1^2y_1^2}\ ,\ \frac{y_1^2 - x_1^2}{1 + dx_1^2y_1^2} \right) \\
&= \left( 0\ ,\ \frac{y_1^2 - x_1^2}{1 + dx_1^2y_1^2} \right) \\
&= \left( 0\ ,\ 1 \right) \\
\end{align*}
```

## 平方根（square root）

有限体 $F_p$ での $p \equiv 5 \pmod{8}$ の場合、$a$ が平方数である場合（$a=x^2$）、$x$ は以下の式で求める事ができます。

```math
x =
\begin{cases}
\pm a^{\frac{p+3}{8}}\\
\pm a^{\frac{p+3}{8}}2^{\frac{p-1}{4}}\\
\end{cases}
```

数学的な証明については割愛します。（書ききれないTT)

プログラムとしてのロジックは、まず平方数 $a$ の $\frac{p+3}{8}$乗、すなわち $a^{\frac{p+3}{8}}$ を計算します。
これを二乗して、元の平方数 $a$ との差が$0$であれば、その値となります。
$0$ でない場合、算出した値に$2^{\frac{p-1}{4}}$を乗算します。
再びこれを二乗して、元の平方数 $a$ との差が$0$であれば、その値となります。
$0$ でない場合、$a$ は平方数ではありません。