仿射坐标下椭圆曲线算术新算法

仿射坐标下椭圆曲线算术的新快速算法

摘要

我们提出了新的算法，通过去除其公式中分子和分母的相同部分来计算 $3P$ 和 $2P+Q$，给定定义在素域和二元域上仿射坐标中的椭圆曲线上的两个点 $P$ 和 $Q$。我们的算法相比西埃、乔伊、劳特和 Montgomery 提出的方法节省了一到两次域乘法。由于 $2P+Q$ 在使用非相邻形式、二进制/三进制方法和树形方法进行标量乘法计算时，分别占所有点运算的 $1/3$ 比例、28.5% 比例和 25.8% 比例，$3P$ 在使用二进制/三进制方法和树形方法进行标量乘法计算时分别占所有点运算的 42.9% 比例和 33.4% 比例，因此利用我们提出的新的 $2P+Q$ 和 $3P$ 公式，通过非相邻形式、二进制/三进制方法和树形方法计算的标量乘法得到了改进。

关键词 ：椭圆曲线 · Double-base数制 · Elliptic曲线算术 · Scalar multiplication

1 引言

自从科布利茨[1]和米勒[2]在1985年独立提出将椭圆曲线用于密码学以来，针对加速椭圆曲线密码系统（ECC）的研究在信息安全与密码学领域得到了广泛而积极的开展。与其他公钥密码系统（如RSA）相比，椭圆曲线密码学能够以短得多的密钥长度提供同等的安全性。因此，它更适合于在智能卡或射频识别等内存和功耗受限的嵌入式系统上的密码学应用。此类密码学应用的瓶颈在于密码学方案的时间和空间消耗。因此，椭圆曲线算术的速度和内存优化一直被持续探索，成为最具吸引力的研究主题之一。

在投影坐标和雅可比坐标下的点运算得到了极大的发展[3–6]。由于投影坐标和雅可比坐标占用更多存储空间，我们专注于仿射坐标。目前已有一些高效计算标量乘法的算法。非邻接形式（NAF）[7]是在所有二进制表示中计算标量乘法最高效的表示方法，其需要计算给定椭圆曲线上的两个点 $2P$, $2P+Q/P+Q$ 的运算。$P$ 和 $Q$ 双基数系统在过去十年中被广泛应用于加速标量乘法[8–15]，其需要计算给定椭圆曲线上的两个点 $2P$, $2P+Q/P+Q$, $3P$ 的运算。$P$ 和 $Q$。

本工作聚焦于在素域和二元域上的椭圆曲线仿射坐标中加速计算 $2P+Q$, $3P$。我们的主要思想是消除其公式中分子和分母的相同部分。我们改进了定义在素域和二元域上的椭圆曲线的 $2P+Q$ 和 $3P$ 的公式。然而，该方法无法加速操作 $P+Q$, $2P$。令 $P(x_1, y_1)$, $Q(x_2, y_2)$, $d = (x_2 - x_1)^2 \cdot (2x_1 + x_2) - (y_2 - y_1)^2$，在 [16] 中计算 $2P+Q$，其中求逆为 $\frac{1}{d \cdot (x_2 - x_1)}$。考虑到 $(x_3, y_3) = 2P+Q$，求逆仅为 $\frac{1}{d}$，从而在素域上节省2次域乘法。将此技巧应用于二元域上的 $2P+Q$、素域上的 $3P$ 以及二元域上的 $3P$，这些公式分别节省1次域乘法、1次域平方，以及2次域乘法加2次域平方。

$2P+Q$ 占 NAF 中所有点运算的 $1/3$ [7] 和 $3P$ 在二进制/三进制中占 33.3%，[16] 在树形方法中占 26.6%[10]。总之，$2P+Q$ 和 $3P$ 的改进带来了标量乘法的加速。

2 预备知识

在椭圆曲线密码学中，主要使用两类曲线。其中一类是定义在二元域上的非超奇异椭圆曲线：
$$
E_2: y^2 + xy = x^3 + ax^2 + b, \tag{1}
$$
其中 $a, b \in \mathbb{F}_{2^n}, b \neq 0$。另一个是定义在大特征域上的椭圆曲线（在本工作的其余部分中，指特征大于三的域，记为 $\mathbb{F}_p$）：
$$
E_p: y^2 = x^3 + ax + b, \tag{2}
$$
其中 $a, b \in \mathbb{F}_p, \Delta = 4a^3 + 27b^2 \neq 0$。

众所周知，椭圆曲线上的所有有理点在由弦切规则定义的加法法则下构成一个阿贝尔群。对于椭圆曲线(1)和(2)上的两个点 $P=(x_1, y_1)$, $Q=(x_2, y_2)$ ($x_1 \neq x_2$)，计算 $2P = P + P$ 和 $P + Q$ 的基本公式如下（参见 [6,17,18]）。这些运算通常分别称为点加倍和点加法。这些运算的输出表示为 $(x_3, y_3)$。

对于曲线（1）在 $\mathbb{F}_{2^n}$ 上，

点加倍
$$
\lambda = \frac{x_1^2 + y_1}{x_1} \
x_3 = \lambda^2 + \lambda + a \
y_3 = (x_1 + x_3)\lambda + x_3 + y_1
\tag{3}
$$

点加法
$$
\lambda = \frac{y_1 + y_2}{x_1 + x_2} \
x_3 = \lambda^2 + \lambda + a + x_1 + x_2 \
y_3 = (x_1 + x_3)\lambda + x_3 + y_1.
\tag{4}
$$

对于曲线（2）在 $\mathbb{F}_p$ 上，

点加倍
$$
\lambda = \frac{3x_1^2 + a}{2y_1} \
x_3 = \lambda^2 - 2x_1 \
y_3 = \lambda(x_1 - x_3) - y_1
\tag{5}
$$

点加法
$$
\lambda = \frac{y_2 - y_1}{x_2 - x_1} \
x_3 = \lambda^2 - x_1 - x_2 \
y_3 = \lambda(x_1 - x_3) - y_1.
\tag{6}
$$

显然，计算 $2P+Q$ 有两种方法。一种是使用(3)或(5)计算 $2P$，然后使用(4)或(6)将 $Q$ 加到结果上。另一种是通过两次连续的点加法将 $2P+Q$ 计算为 $(P+Q)+P$。当然，结果是相等的，但验证计算并不简单。因此，在加速椭圆曲线密码学的研究中，通常认为这两种方式是不同的。在[19]中首次指出，$(P+Q)+P$ 方式比 $2P+Q$ 方式更高效。

他们的发现如下：对于曲线（2）上的 $\mathbb{F}_p$ 例如，$(x_5, y_5) = ((x_1, y_1)+(x_2, y_2)) + (x_1, y_1)$ 通过 $(P+Q)+P$ 方式计算，即
$$
\lambda = \frac{y_2 - y_1}{x_2 - x_1},\quad x_3 = \lambda^2 - x_1 - x_2,\quad y_3 = \lambda(x_1 - x_3) - y_1 \
\lambda_1 = \frac{y_3 - y_1}{x_3 - x_1},\quad x_5 = \lambda_1^2 - x_1 - x_3,\quad y_5 = \lambda_1(x_1 - x_5) - y_1.
\tag{7}
$$

然后，由于 $\lambda_1 = \frac{y_3 - y_1}{x_3 - x_1} = -\lambda + \frac{2y_1}{x_1 - x_3}$, $y_3$ 无需获得 $(x_5, y_5)$，因此在 $y_3$ 计算中可以节省一次域乘法。

在本工作中，我们将采用通用方法，将点运算的代价表示为域求逆（I）、乘法（M）和平方（S）的次数，为了简化起见，忽略域加法/减法以及与小常数的乘法的代价。

在 [16] 中表明，不仅 $y_3$，而且 $x_3$ 在计算 $(x_5, y_5)$ 时也可以忽略，并且当一次域求逆的代价不低于六次域乘法时，存在一种更高效的 $2P+Q$ 算法。以下是文献[16]中提出的关于 $2P+Q$ 在仿射坐标下的算法（$\mathbb{F}_p$ 上，此后简称为 A）：

[算法 1] ($2P+Q$ 在 $\mathbb{F}_p$ [16] 上)
输入 ：
$P=(x_1, y_1)$, $Q=(x_2, y_2) \in E_p$
输出 ：$(x_3, y_3) = 2P+Q$

$$
d = (x_2 - x_1)^2 \cdot (2x_1 + x_2) - (y_2 - y_1)^2, \
D = d \cdot (x_2 - x_1),\quad I = D^{-1}, \
\lambda = (y_2 - y_1) \cdot dI, \
\lambda_1 = -\lambda + 2y_1 \cdot (x_2 - x_1)^3 I, \
x_3 = (\lambda_1 - \lambda) \cdot (\lambda_1 + \lambda) + x_2, \
y_3 = (x_1 - x_3) \cdot \lambda_1 - y_1.
$$

如[16,20,21]、 $1I \geq 8M$ 所示，因此算法1优于[19]中提出的算法。

算法1中的相同思想也适用于二元域上的曲线(1)，将 $2P+Q$ 成本降低至 $1I + 9M + 2S$（见现[16]），我们将进一步降低这些成本。

3 加快仿射坐标下某些椭圆曲线运算

基于 [16] 中的 $2P+Q$, $3P$ 公式，我们提出了新的公式。

引理1 。 存在一种算法，在给定大特征域 $\mathbb{F}_p$ 上椭圆曲线 $y^2 = x^3 + ax + b$ 的点 $P=(x_1, y_1)$ 和 $Q=(x_2, y_2)$ 的情况下，以 $1I + 7M + 2S$ 的代价计算 $2P+Q = (x_3, y_3)$。

证明 。 设 $d = (x_2 - x_1)^2 \cdot (2x_1 + x_2) - (y_2 - y_1)^2$，$\lambda = (y_2 - y_1) \cdot (x_2 - x_1)^{-1}$。在算法1中，我们注意到
$$
x_3 = 2(-\lambda + y_1 \cdot (x_2 - x_1)^2 d^{-1}) \cdot (2y_1 \cdot (x_2 - x_1)^2 \cdot d^{-1}) + x_2 \
= 2(y_1 - y_2 + y_1 \cdot (x_2 - x_1)^3 d^{-1}) \cdot (2y_1 \cdot (x_2 - x_1) \cdot d^{-1}) + x_2.
$$
并且 $y_3$ 可以计算为
$$
y_3 = (x_1 - x_3) \cdot (-\lambda + 2y_1 \cdot (x_2 - x_1)^2 d^{-1}) - y_1 \
= [-2(y_1 - y_2 + y_1 \cdot (x_2 - x_1)^3 d^{-1}) \cdot (2y_1 \cdot (x_2 - x_1) \cdot d^{-1}) + (x_1 - x_2)] \
\cdot (-\lambda + 2y_1 \cdot (x_2 - x_1)^2 d^{-1}) - y_1 \
= [-2(y_1 - y_2 + y_1 \cdot (x_2 - x_1)^3 d^{-1}) \cdot (2y_1 \cdot d^{-1}) - 1] \
\cdot (x_2 - x_1)(-\lambda + 2y_1 \cdot (x_2 - x_1)^2 d^{-1}) - y_1 \
= [-2(y_1 - y_2 + y_1 \cdot (x_2 - x_1)^3 d^{-1}) \cdot (2y_1 \cdot d^{-1}) - 1] \
\cdot (-(y_2 - y_1) + 2y_1 \cdot (x_2 - x_1)^3 d^{-1}) - y_1.
$$

求逆仅需 $\frac{1}{d}$。然后 $2P+Q = (x_3, y_3)$ 在 $1I + 7M + 2S$ 的代价下显示于算法2中。

[算法 2] $2P+Q$ 在 $\mathbb{F}_p$ 上
输入 ： $P=(x_1, y_1)$, $Q=(x_2, y_2) \in E_p$
输出 ：$(x_3, y_3) = 2P+Q$

$$
d = (x_2 - x_1)^2 \cdot (2x_1 + x_2) - (y_2 - y_1)^2, \
I = d^{-1},\quad A = y_1 \cdot I, \
B = A \cdot (x_2 - x_1),\quad C = B \cdot (x_2 - x_1)^2, \
\lambda_2 = C - (y_2 - y_1), \
x_3 = x_2 + 4B \cdot \lambda_2, \
y_3 = -y_1 + (-1 - 4A \cdot \lambda_2) \cdot (\lambda_2 + C).
$$

令 $d = (x_2 - x_1)^2 \cdot (2x_1 + x_2) - (y_2 - y_1)^2$。在素域上的 $2P+Q$ 计算中，我们通过 $\frac{1}{d}$ 而非 $\frac{1}{d(x_2 - x_1)}$ 来获得 $(x_3, y_3)$。该技巧也可用于二进制域上的 $2P+Q$ 计算。

引理2 。 存在一种算法，可在给定二进制域 $\mathbb{F}_{2^n}$ 上非超奇异椭圆曲线 $y^2 + xy = x^3 + ax^2 + b$ 的点 $P=(x_1, y_1)$ 和 $Q=(x_2, y_2)$ 的情况下，以 $1I + 8M + 2S$ 的代价计算 $2P+Q = (x_3, y_3)$。

证明 。 可以直接验证，以下算法以 $1I + 8M + 2S$ 的代价计算 $(x_3, y_3)$。

[算法 3] ($2P+Q$ 在 $\mathbb{F}_{2^n}$ 上)
输入 ： $P=(x_1, y_1)$, $Q=(x_2, y_2) \in E_p$
输出 ： $(x_3, y_3) = 2P+Q$

$$
d = x_1 \cdot (y_2 + x_1^2) + x_2 \cdot (y_1 + x_1^2), \
I = d^{-1},\quad N = I \cdot (x_1^2 + x_1 \cdot x_2), \
M = (x_1 + x_2) \cdot N, \
x_3 = M^2 + M + x_2, \
y_3 = [1 + N \cdot (1 + M)] \cdot [y_1 + y_2 + (x_1 + x_2) \cdot (1 + M)] + x_3 + y_1.
$$

该技巧也可应用于 $P = Q$ 的情况，即 $3P$。据我们所知，迄今为止，在大特征和二进制域椭圆曲线上，$3P$ 计算的仿射坐标速度记录是 $1I + 7M + 4S$ 在 [16] 上。

引理3 。 存在一种算法，在给定大特征域 $\mathbb{F}_p$ 上椭圆曲线 $y^2 = x^3 + ax + b$ 的一个点 $P=(x, y)$ 的情况下，以 $1I + 7M + 3S$ 的代价计算 $3P = (x_3, y_3)$。

证明 。 与上述类似，我们将给出一个构造性证明。以下计算 $3P$ 的算法耗费 $1I + 7M + 3S$。

[算法 4] ($3P$ 在 $\mathbb{F}_p$ 上)
输入 ： $P=(x, y) \in E_p$
输出 ：$(x_3, y_3) = 3P$

$$
d = 12x \cdot y^2 - (3x^2 + a)^2, \
I = d^{-1},\quad A = 4y \cdot I, \
B = 8y^2 \cdot I,\quad C = y^2 \cdot B, \
x_3 = B \cdot (C - 3x^2 - a) + x, \
y_3 = A \cdot (3x^2 + a - C) \cdot (2C - 3x^2 - a) - y.
$$

应用于二进制域 $3P$ 计算的技巧可带来更多的节省。

引理4 。 存在一种算法，在给定二进制域椭圆曲线 $y^2 + xy = x^3 + ax^2 + b$ 上的点 $P=(x_1, y_1)$ 的情况下，以 $1I + 5M + 2S$ 的代价计算 $3P = (x_3, y_3)$。

证明 。 以下算法验证了该引理。

[算法 5] ($3P$ over $\mathbb{F}_{2^n}$)
输入 ： $P=(x_1, y_1) \in E_2$
输出 ：$(x_3, y_3) = 3P$

$$
d = x_1^2 \cdot (x_1^2 + x_1) + b, \
I = d^{-1},\quad M = I \cdot x_1^2, \
N = x_1 \cdot M, \
x_3 = N^2 + N + x_1, \
y_3 = (1 + N) \cdot [N^2 + M \cdot (x_1^2 + y_1)] + x_1 + y_1.
$$

4 比较

为了与其他研究进行比较，我们重点关注了点运算的成本及其在标量乘法中的应用。

4.1 点运算

本工作提出的算法2–5在素域和二元域上的仿射坐标下，改进了 $2P+Q$ 和 $3P$ 的运算速度记录。表1比较了本工作所提出算法的开销与之前速度记录之间的差异。

$2P+Q$ 的公式相比素域上的[16]节省了 $2M$，相比二元域上的[16]节省了 $1M$。而 $3P$ 的公式相比素域上的[16]节省了 $1S$，相比二元域上的[16]节省了 $2M+2S$。

接下来，我们将分析点运算的改进如何影响标量乘法。

表1. 不同有限域 $\mathbb{F}_p$ 上椭圆曲线各种算术算法的成本比较

	奇耶‐乔伊‐劳特‐蒙哥马利 [16]	本工作	节省成本
$2P+Q$	$1I + 9M + 2S$	$1I + 7M + 2S$	$2M$
$3P$	$1I + 7M + 4S$	$1I + 7M + 3S$	$1S$

表2. 在 $\mathbb{F}_{2^n}$ 上椭圆曲线各种算术算法的成本比较

	西埃‐乔伊‐劳特‐蒙哥马利 [16]	本工作	节省成本
$2P+Q$	$1I + 9M + 2S$	$1I + 8M + 2S$	$1M$
$3P$	$1I + 7M + 4S$	$1I + 5M + 2S$	$2M + 2S$

表3. 不同方法的点运算比率

	NAF [7]	二进制/三进制 [16] $0.2285 / 0.7996 = 0.286$	树形方法 [10] $0.3415 / 0.8364 = 0.408$
$2P$	$2/3$	$0.2284 / 0.7996 = 0.285$	$0.2154 / 0.8364 = 0.258$
$2P+Q$	$1/3$	$0.3427 / 0.7996 = 0.429$	$0.2795 / 0.8364 = 0.334$
$3P$	—	—	—

表4. 在素域上的椭圆曲线仿射坐标中每比特标量乘法的代价

	NAF [7]	二进制/三进制 [16]	树形方法 [10]
[6,16]	$I + \frac{13}{3}M + 2S$	$0.7996I + 4.9115M + 2.2846S$	$0.8364I + 4.5981M + 2.2318S$
本工作	$I + 11M + 2S$	$0.7996I + 4.4547M + 1.9419S$	$0.8364I + 4.1473M + 1.9523S$

表5. 在二元域上椭圆曲线仿射坐标中每比特标量乘法的代价

	NAF [7]	二进制/三进制 [16]	树形方法 [10]
[6,16]	$I + \frac{13}{3}M + 2S$	$0.7996I + 4.9115M + 2.2846S$	$0.8364I + 4.5981M + 2.2318S$
本工作	$I + 4M + 2S$	$0.7996I + 3.9977M + 1.5992S$	$0.8364I + 3.8037M + 1.6728S$

4.2 标量乘法

在[6]的成本上，$P+Q$ 和 $2P$ 分别为 $I + 2M + S$ 和 $I + 2M + 2S$，在素域和二元域上均如此。

不同标量乘法算法（如NAF、二进制/三进制方法和树形方法）中所使用的 $P+Q$、$2P$、$2P+Q$ 和 $3P$ 点运算的比例见表3。$2P+Q$ 在采用非相邻形式、二进制/三进制方法和树形方法计算标量乘法时，分别占所有点运算的 $1/3$、28.5% 和 25.8%。$3P$ 在采用返回双基表示来计算标量乘法的二进制/三进制方法和树形方法中，分别占所有点运算的 42.9% 和 33.4%。

使用NAF、二进制/三进制方法和树形方法在素域和二元域上每比特标量乘法的代价总结于表4和5中，这些数据是根据表1、2和3计算得出的。表4和表5中的结果表明，使用本工作提出的公式进行标量乘法明显快于使用[16]。

5 结论

我们提出了在素域和二元域上仿射坐标中进行点运算 $2P+Q$ 和 $3P$ 的新公式，其速度优于[16]中的相应公式。尽管当一次域求逆的代价 > 八次域乘法时，仿射坐标中的点运算通常比投影坐标中的慢，但它们通常占用更小的空间。本工作提出的点运算 $2P+Q$ 和 $3P$ 在一次域求逆的代价 < 四次域乘法的情况下，在[3,4]中比投影坐标中的点运算更快。