A practical pilot-based channel estimation method for space-time coded ofdm system

1.   引言

椭圆曲线密码(Elliptic Curve Cryptography, ECC)已经被美国国家标准学会(American National Standards Institute, ANSI)^[1]、美国国家标准与技术研究所(National Institute of Standards and Technology, NIST)^[2]、电气电子工程师协会(Institute of Electrical and Electronics Engineers, IEEE)^[3]等国际标准组织广泛接受和部署^[4]，应用在各种各样的场景。在比特币的基础技术区块链中，广泛应用ECC中的secp256k1作为区块链的认证机制^[5]；我国自主研发的密码算法SM2密码标准，是国密体系中重要环节之一。ECC的应用日益增加，不同应用性能需求各异。例如在服务器和大型计算平台等，需要密码芯片具有较高的运算性能。在微控制器等设备中，一般要求密码芯片具有更小的面积。此外不同的标准通常使用不同的椭圆曲线，所以对多种曲线的兼容也一直是椭圆曲线的重点研究目标。

标量乘架构的研究成为热点，Hu等人^[6]提出SCA-256曲线的二阶段快速模约简方法，利用半字乘法器5个周期计算大数乘法实现高速标量乘的设计；文献[7]和文献[8]都对wNAF算法进行了优化，降低了标量乘的计算复杂度，提高了标量乘计算速度；Liu等人^[9]利用软硬结合的方式以较低的面积实现了支持任意曲线的双域ECC处理器；Hossain等人^[10]设计了同时适用于专用集成电路(Application Specific Integrated Circuit, ASIC)和现场可编程逻辑门阵列(Field Programmable Gate Array, FPGA)的标量乘操作步骤，适用于两条曲线；Liu等人^[11]优化了非相邻形式算法，以高资源消耗实现了支持P256的高速ECC标量乘设计；Hu等人^[12]和Choi等人^[13]都提出高度资源复用的设计，以极低的硬件资源实现标量乘法器；Zhang等人^[14]改进了蒙哥马利阶梯算法，提出了支持SCA-256的抗功耗攻击的高性能处理器。文献[15]改进了基2交叉模乘算法，提出了一种高面积时间效率的新型标量乘法器架构。这些设计的时间面积(Area·Time, AT)普遍较差，部分仅适用于特定曲线，ECC处理器需要兼顾速度和面积两者之间的关系，为解决此问题本文提出了全新的设计。

本文提出了一种基于比特重组快速模约简的高面积效率ECC标量乘法器，支持secp256k1，secp256r1和SCA-256 3条曲线。本文的主要工作为：

(1)设计了一种兼顾乘法和模逆两种功能的硬件复用运算单元，提高硬件资源使用率。同时采用Karatsuba-Ofman算法提高运算单元计算性能；

(2)设计了基于比特重组快速模约简算法，实现了支持secp256k1, secp256r1和SCA-256快速模约简计算的硬件架构；

(3)优化了点加和倍点的模运算操作调度，使乘法和快速模约简排布尽量紧密，减少标量乘计算需要的周期数量。

2.   研究背景

椭圆曲线有很多种，使用较多的是Weierstrass曲线^[16]。这一类曲线满足式(1)：

NIST提出的曲线secp256r1和secp256k1以及国密sm2密码算法提出的曲线SCA-256都属于Weierstrass曲线。

标量乘是ECC的核心运算，定义为$kP = P + \cdots + P$。其中$k$是整数，$P$是椭圆曲线上一点。标量乘算法较多，主要有二进制标量乘算法、非相邻形式(Non-Adjacent Form, NAF)标量乘算法等^[17]。NAF标量乘算法^[18]如算法1所示，通过对$k$进行重新编码，减少$k$中非0元素出现的次数以达到提高标量乘性能的目的。

算法1 NAF标量乘算法
输入：$k,P$
输出：$Q = kP$
步骤1：计算${\text{NAF} }\left( k \right) = \left\{ { {k_{m - 1} }, \cdots ,{k_0} } \right\}$
步骤2：$Q \leftarrow \infty$
步骤3：${\text{for }}i \leftarrow m - 1{\text{ down to }}0{\text{ do}}$
$Q \leftarrow 2Q$
${\text{if }}{k_i} = 1{\text{ then}}$ $Q \leftarrow Q + P$
${\text{if }}{k_i} = - 1{\text{ then}}$ $Q \leftarrow Q - P$
步骤4：返回$Q$

下载: 导出CSV 

| 显示表格

标量乘运算是由一系列的点加(Point Addition, PA)和倍点(Point Doubling, PD)操作实现的。假设${P_1} = \left( {{X_1},{Y_1},{Z_1}} \right)$, ${P_2} = \left( {{X_2},{Y_2},{Z_2}} \right)$, ${P_3} = \left( {{X_3},{Y_3},{Z_3}} \right)$是椭圆曲线上的3个点，那么点加操作可以定义为${P_3} = {P_1} + {P_2}$，倍点操作可以定义为${P_3} = 2{P_1}$。为了避免模逆这一复杂的操作，通常使用转换坐标系的方式减少模逆的计算次数。仿射-雅可比混合坐标系可以产生更快的点加操作，雅可比坐标系可以产生更快的倍点操作^[19]。

仿射-雅可比混合坐标系下的点加公式如式(2)所示：

雅可比坐标系下的倍点公式如式(3)所示：

3.   标量乘法器设计

ECC应至少包含3个算术层：有限域运算层、点运算层和标量乘运算层，本节将按照标量乘法器的结构层次，自底向上介绍本文提出的标量乘法器设计。

3.1   资源共享的运算单元设计

本文使用计算速度较快的快速模约简模乘，该模乘分为整数乘法运算和快速模约简运算两部分，整数乘法运算性能对模乘影响较大。使用大位宽乘法器计算乘法只需要一个周期，但是需要消耗大量硬件资源，使用小位宽乘法器计算乘法虽然面积较小，但是需要较多周期计算乘法，速度较慢，因此选择大小合适的乘法器是设计该部分的关键。Karatsuba–Ofman算法^[20]可以将两个256 bit的大数乘法转化为3个128 bit的乘法和一些加减及移位操作，如式(4)所示。

该算法虽然可以减少乘法大小，但是同样带来了额外一些额外的计算。为了减少这些额外的计算带来的影响，优化了Karatsuba算法的操作调度方式，如算法2所示。

算法2 改进的Karatsuba算法操作调度
输入：$A={a}_{1}{2}^{128}+{a}_{0},\;B={b}_{1}{2}^{128}+{b}_{0}$
输出：$C = A \cdot B$
1：$s = {a_1}\cdot{b_1}$；$v\left[ {128:0} \right] = {\text{ }}{a_0} + {a_1}$； 　　 $v\left[ {257:129} \right] = {b_0} + {b_1}$；
2：$r = {a_0}\cdot{b_0}$；${c_1} = r + s$；$u = {\text{～}}{c_1} + 1$；
3：$v = v\left[ {128:0} \right] \cdot v\left[ {257:129} \right]$；${c_1} = v + u$； 　　 ${c_2} = \left\{ {s,r\left[ {255:128} \right]} \right\} + {c_1}$；
4：返回$C = \{ {c_2},r\left[ {127:0} \right]\}$

下载: 导出CSV 

| 显示表格

本文提出的调度方式将式(4)分为3个周期进行计算，前两个周期使用乘法器和加法器可以得到${a_1}{b_1}$, ${a_0}{b_0}$, ${a_0} + {\text{ }}{a_1}$, ${b_0} + {b_1}$以及$- ({a_1}{b_1} + {a_0}{b_0})$的结果，分别存入寄存器$s$, $v$, $r$和$u$中。第3个周期利用乘法器和加法器${c_1}$得到$({a_0} + {a_1})({b_0} + {b_1}) - ({a_1}{b_1} + {a_0}{b_0})$的结果。式(3)中可以看出$({a_1}{b_1}{2^{256}} + {a_0}{b_0})$可以通过位拼接直接实现，$[({a_0} + {a_1}) ({b_0} + {b_1}) - ({a_1}{b_1} + {a_0}{b_0})]{2^{128}}$可以通过${c_1}$左移128 bit实现。${c_1}$左移128 bit后低位都为0，$({a_1}{b_1}{2^{256}} + {a_0}{b_0})$的低128 bit即${a_0}{b_0}$的低128 bit，所以两者相加可以先舍去低128 bit，只相加剩余部分，通过加法器${c_2}$得到这一结果，然后将${c_2}$与${a_0}{b_0}$的低128 bit进行拼接，可以得到最终的乘法结果C。

优化后的调度方式需要使用1个乘法器和2个加法器，以及4个存储中间结果的寄存器，经过3个时钟周期可以得到512 bit的乘法结果C。

为进一步提高设计的面积效率，本文提出了一种乘法与模逆可以实现资源共享的运算单元设计，此设计尽可能地复用加法器和寄存器实现乘法和模逆两种运算。模逆使用RS-A模逆算法^[21]，结合算法2及RS-A模逆算法设计出的新的运算单元如图1所示。此运算单元通过控制信号切换功能，当使用乘法功能时，进位保存加法器(Carry Save Adder, CSA)相当于算法2中的加法器${c_1}$，只需将第3个输入置0即可，加法器1相当于算法2中的加法器${c_2}$。然后根据算法2中的调度方式控制乘法器和加法器选择不同的输入实现乘法计算的功能，最后512 bit乘法结果分为高256 bit和低256 bit分别存入4个寄存器中，方便后续快速模约简模块进行比特重组。

图  1  资源共享的运算单元

下载: 全尺寸图片 幻灯片

使用模逆功能时，CSA加法器负责RS-A模逆算法中加减模数$p$的相关计算，除2操作通过移位实现。加法器1负责RS-A模逆算法中的其他加减计算，除2也通过移位实现。然后根据RS-A模逆算法中的不同条件控制寄存器的输入输出选择和加法器的输入选择，最终实现模逆功能。

此运算单元在不同功能的各个工作阶段，通过数据选择器选择不同的数据来源，即可实现乘法功能或模逆功能。因为乘法和模逆运算过程中都需要使用加法器，所以两种运算不能同时进行。通过转换坐标系的方式，标量乘运算过程中模乘和模逆不会同时进行，所以不会影响整体性能。与文献[6]相比，我们提出的结构乘法速度更快，也不需要额外的资源实现单独的模逆；与文献[22]相比乘法运算消耗周期数相同，但本文资源共享的结构使得面积更有优势。

3.2   基于比特重组的快速模约简设计

模约简需要根据曲线特点进行单独设计，不能直接套用其他曲线模约简公式。对于secp256k1，$p = {2^{256}} - {2^{32}} - {2^9} - {2^8} - {2^7} - {2^6} - {2^4} - 1$ = ${2^{256}} - {2^{32}} - {2^{10}} + {2^6} - {2^4} - 1$，$p$的位宽为256 bit。乘法结果最高为512 bit，将乘法结果表示为$C = {c_{127}}{2^{508}} + {c_{126}}{2^{504}} + \cdots + {c_1}{2^4} + {c_0}$，其中${c_i}$位宽为4 bit。根据$p$的特点，secp256k1有如式(5)同余关系式：

模约简过程中，首先需要处理$C$的高256 bit，即${c_{127}}{2^{508}} + \cdots + {c_{64}}{2^{256}}$。结合式(5)，可以得到式(6)：

乘法结果的高256 bit可以通过式(6)右半部分的式子相加得到约简后的结果，然后再与乘法结果的低256 bit${c_{63}}{2^{252}} + \cdots + {c_0}$相加。将其按幂次由高到低排列，整理为和项和差项后，得到算法3。本文使用两个步骤完成secp256k1快速模约简的计算。首先第1步将乘法结果$C$分解为128个4位数据段，然后重组为${s_1}{\text{～}}{s_{16}}$，此过程为比特重组。然后将重组后的16个数据段进行计算得到${z}_{1}\in (-3p, 8p)$。第1步的结果可能产生进位或借位，所以第二阶段需要对其再进行处理。当${z_1}$大于0时，我们将${z_1}$高位补0然后重复一次第1步操作，约简结果${z}_{2}\in [0，2p)$，最多还需要一次减法就能得到最终约简结果。当${z_1}$小于0时，通过最多3次累加可以得到最终约简结果。

算法3 secp256k1快速模约简算法
输入：$C = {c_{127}}{2^{508}} + {c_{126}}{2^{504}} + \cdots + {c_1}{2^4} + {c_0}$, $0 \le C < {p^2}$
输出：$C\text{mod} p$
(1) ${s_1} = ({c_{63}}, \cdots ,{c_0})$；
${s_2} = ({c_{127}}, \cdots ,{c_{64}})$；
${s_3} = ({c_{126}}, \cdots ,{c_{64}},{c_{127}})$；
${s_4} = (2'h0,{c_{124}}, \cdots ,{c_{64}},2'h0,{c_{127}},4'h0)$；
${s_5} = ({c_{119}}, \cdots ,{c_{64}},{c_{127}}, \cdots ,{c_{120}})$；
${s_6} = (214'h0,{c_{127}},2'h0,{c_{127}}, \cdots ,{c_{120}},4'h0)$；
${s_7} = (214'h0,{c_{127}}, \cdots ,{c_{120}},{c_{126}},{c_{126}},2'h0)$；
${s_8} = (192'h0,{c_{127}}, \cdots ,{c_{120}},12'h0,{c_{127}},2'h0,$ 　　 　　 ${c_{127}},10'h0)$；
${s_9} = (218'h0,{c_{126}},18'h0,{c_{126}},{c_{127}},8'h0)$；
${s_{10}} = (220'h0,{c_{127}},18'h0,{c_{127}},10'h0)$；
${s_{11}} = ({c_{125}},254'h0)$；
${s_{12}} = (2'h0,{c_{125}}, \cdots ,{c_{64}},{c_{127}},2'h0)$；
${s_{13}} = (218'h0,{c_{127}}, \cdots ,{c_{120}},6'h0)$；
${s_{14}} = (218'h0,{c_{127}},18'h0,{c_{127}},2'h0,{c_{127}},6'h0)$；
${s_{15}} = (240'h0,{c_{127}},{c_{126}},8'h0)$；
${s_{16}} = ({c_{126}},254'h0)$；
${{\text{z}}_1} = {s_1} + {s_2} + {s_3} + {s_4} + {s_5} + {s_6} + {s_7} + {s_8} + {s_9} + {s_{10}} +$ 　　 ${s_{11}} - {s_{12}} - {s_{13}} - {s_{14}} - {s_{15}} - {s_{16}}$；
$= ({r_{65}}, \cdots ,{r_0})$
(2) ${t_1} = ({r_{63}}, \ldots ,{r_0})$；
${t_2} = \left( {216'h0,{\text{ }}{r_{65}},{r_{64}},20'h0,{\text{ }}{r_{65}},{r_{65}},{r_{64}}} \right)$；
${t_3} = \left( {235'h0,{\text{ }}{r_{65}},{r_{64}},2'h0,{\text{ }}{r_{64}},4'h0} \right)$；
${t_4} = \left( {242'h0,{\text{ }}{r_{65}},{r_{64}},6'h0} \right)$；
${z_2} = {\text{ }}{t_1} + {t_2} + {t_3} - {t_4}$
(3)若$({z_1} \ge 0)$，返回$\left( {{z_2}\text{mod} p} \right)$；
若$( - p < {z_1} < 0)$，返回$\left( {{z_1} + p} \right)$；
若$( - 2p < {z_1} \le - p)$，返回$\left( {{z_1} + 2p} \right)$；
若$( - 3p < {z_1} \le- 2p)$，返回$\left( {{z_1} + 3p} \right)$；

下载: 导出CSV 

| 显示表格

文献[14]和文献[23]提供了SCA-256和secp256r1的快速模约简公式，是将512 bit数据拆分为16个32 bit数据段，32 bit数据段可以进一步拆分为8个4 bit数据段，因此同样可以使用比特重组方法进行快速模约简算法的设计。

本文利用比特重组方法提出兼容3条曲线的快速模约简结构，如图2所示。比特重组模块根据不同曲线的快速模约简公式将乘法结果进行拆分与重组，华莱士树模块负责快速模约简算法中第1步和第2步中${z_1} \ge 0$时的计算，加法器阵列负责第2步${z_1} < 0$时的计算。为减少资源消耗和路径延迟，将两个步骤的快速模约简分为3个时钟周期计算。首先第1个时钟周期利用比特重组模块将输入进行重组，重组为${s_1}{\text{～}}{s_{16}}$，其中secp256r1只有${s_1}{\text{～}}{s_{11}}$为有效值，${s_{12}}{\text{～}}{s_{16}}$为全0。然后使用CSA加法器组计算快速模约简公式的第1步中的${s_1} + {s_2} + {s_3} + {s_4} + {s_5} + {s_6} + {s_7} - {s_{12}} - {s_{13}}$，减法通过补码加法实现，结果存入$w$中。第2个周期计算$w + {s_8} + {s_9} + {s_{10}} + {s_{11}} - {s_{14}} - {s_{15}} - {s_{16}}$，结果存入$w$中，此时已经计算完毕${z_1}$。第3个时钟周期先判断$w$的符号位，若$w \ge 0$则使用比特重组模块将$w$重组为${t_1}{\text{～}}{t_4}$，然后使用华莱士树计算得到${z_2}$。若${z_2}$超过$p - 1$则通过加法器完成减$p$操作。若$w < 0$，则使用加法器阵列进行累加，通过判断每一级加法器中间结果的符号位选择正确的约简结果。

图  2  兼容3条曲线快速模约简结构

下载: 全尺寸图片 幻灯片

虽然快速模约简相比常规方式需要周期数更多，但是本文通过对点加和倍点的模运算操作调度进行优化降低了其负面影响，将在下一节详细介绍。

3.3   高速点加和倍点设计

标量乘运算通过点加和倍点运算实现，因此对于点加和倍点的优化也十分重要。仿射-雅可比混合坐标系下的点加和雅可比坐标系下的倍点公式模乘计算次数较少^[22]，结合上文中提出的设计对点加和倍点的模运算操作调度进行了优化。为使乘法与模约简的排布尽量紧密，减少标量乘运算所需周期数量，本文提出了倍点和点加两种操作合一的模运算操作调度方案。

倍点和倍点点加合一的模运算操作调度如图3所示，图3(a)为倍点的模运算操作调度，图3(b)为倍点点加合一的模运算操作调度。其中雅可比坐标系输入点$({X_1},{Y_1},{Z_1})$用$({t_5},{t_6},{t_7})$表示，仿射坐标系输入点用$({x_2},{y_2})$表示，得到的雅可比坐标系下点$({X_2},{Y_2},{Z_2})$也用$({t_5},{t_6},{t_7})$，${t_1}{\text{～}}{t_7}$为寄存器。

图  3  高速点加和倍点的模运算操作调度

下载: 全尺寸图片 幻灯片

通过优化后的调度方法，1次模乘平均只需要3个周期，消除了因快速模约简计算需要3个周期带来的影响。标量乘运算过程中点加和倍点操作是连续进行的，我们提出的两种调度方式在最后3个周期可以开始下一次运算的第1次乘法，此时需要额外的周期计算第1次模加。故1次倍点平均需要27个周期，1次倍点点加需要60个周期。

3.4   总体设计

本文将secp256k1, secp256r1, SCA-2563条曲线的标量乘集成为1个单元进行设计，图4给出了标量乘法器的硬件架构。标量乘法器由标量乘控制模块、倍点和倍点点加控制模块、NAF编码器、模运算单元和寄存器组5部分组成。

图  4  标量乘法器

下载: 全尺寸图片 幻灯片

标量乘控制模块对标量乘法器整体进行控制，控制NAF编码器、点运算和有限域运算的开始和结束，使标量乘法器可以准确地完成标量乘运算。切换曲线也通过该模块实现，当切换不同曲线时，标量乘控制模块控制点运算和有限域运算选取不同的曲线参数。倍点控制和倍点点加控制模块负责点运算的顺利完成，控制底层模运算单元的开始和结束、模运算单元和寄存器组的输入与输出选择。模运算单元包括模加减模块、资源共享的运算单元和快速模约简模块，负责有限域运算。寄存器组存储模运算单元产生的所有中间结果。

4.   实验结果

表1是使用55 nm CMOS工艺综合得到标量乘法器中主要模块面积，可以看到资源共享的运算单元为151k等效二输入与非门，占整个标量乘法器面积的54%，是占比最大的模块。表2给出了每个椭圆曲线操作需要的周期数量。1次标量乘共需要10366个周期，在500 MHz的频率下计算1次标量乘需要0.02 ms，每秒可以计算48309次标量乘。

表  1  主要模块面积

模块	NAF编码器	资源共享 运算单元	模约简	标量乘法器
面积(K·Gates/LUTs)	16	151	45	275

下载: 导出CSV 

| 显示表格

表  2  椭圆曲线操作周期

	乘法	模约简	模加/减	模逆	NAF编码	倍点	倍点点加	标量乘
周期	3	3	1	360	256	27	60	10366

下载: 导出CSV 

| 显示表格

本文标量乘法器最终性能如表3所示，为了得到合理的对比，在ASIC平台使用55 nm和130 nm的工艺进行了前端逻辑综合，在FPGA平台使用Virtex-7进行了实现，并与类似工作进行了比较。

表  3  256位素数域下椭圆曲线标量乘法器的性能比较

方案	工艺(nm)	曲线	主频(MHz)	周期	面积(K·Gates/LUTs)	单次时间(μs)	吞吐量(次/s)	AT
本文	55	secp256k1/r1, SCA-256	500	10366	275	20.7	48309	5.7
文献[22]	55	secp256r1	300	7955	306	27	37712	8.1
文献[9]	55	256	136	458200	189	1450	690	274
文献[10]	65	256	546.5	397300	447	730	1370	326.3
文献[11]	65	secp256r1	188	2300	3500	12.5	80000	43.7
本文	130	secp256k1/r1, SCA-256	243	10366	314	42.5	23529	13.3
文献[6]	130	SCA-256	250	14240	280	57	17544	15.9
文献[12]	130	256	150	610000	57.1	4070	246	232
文献[14]	130	SCA-256	214	45500	208	211	4740	44.2
文献[13]	180	secp256r1	185	36800	144.8	71	14084	28.9
本文	Virtex-7	secp256k1/r1, SCA-256	37	10366	47.4	279	3584	13.2
文献[24]	Virtex-7	256	177.7	262.7	32.7	1471	679	48.5
文献[25]	Virtex-7	256	72.9	215.9	96.9	2957	338	286.8

下载: 导出CSV 

| 显示表格

文献[22]使用快速模约简方式模乘，支持单曲线secp256r1，使用滑动窗口NAF标量乘算法，标量乘消耗周期数少，但是主频较低而且面积更大，所以本文AT值更有优势。文献[9]使用基4交错模乘，标量乘法器面积较小，但一次模乘需要大量周期，导致吞吐量和AT较差。文献[10]提出新的倍点和点加并行技术减少操作步骤，主频是所有设计中最高的，但是面积大，而且标量乘计算速度慢，与之相比本文的设计在AT值和吞吐量都很有优势。文献[11]采用128 bit的分段模乘，一次标量乘仅需0.0125 ms，每秒可实现80000次标量乘，是所有设计中吞吐量最高的，但是标量乘法器面积是本文的12倍，AT值远大于本文提出的设计。文献[6]设计面积比本文小，但一次乘法需要5个周期计算，计算标量乘需要的周期比本文多，所以标量乘计算速度比本文慢，AT值比本文高。文献[12]和文献[13]都以硬件资源消耗最少为目标，所以速度都较慢，AT值都远远大于本文。文献[14]通过改进蒙哥马利阶梯算法，使得点加和倍点同时进行，虽然面积较小，但标量乘速度慢，吞吐量和AT值都差于本文。文献[24]采用蒙哥马利阶梯算法实现标量乘，虽然采用了转换坐标系方式提高标量乘计算速度，但是基2交错模乘速度慢，吞吐量和AT都不高。文献[25]提出了一种基于RNS的ECC处理器，可同时处理多个标量乘操作，但AT差于本文设计。

5.   结论

本文提出了一种基于比特重组快速模约简的高面积效率ECC标量乘法器设计。利用Karatsuba-Ofman算法设计出仅需3个时钟周期的乘法器，提升了乘法的计算速度。基于此乘法器，又提出乘法与模逆资源共享的运算单元设计，在原乘法器结构基础上只增加一个CSA加法器可实现乘法或模逆的计算，在不影响两者速度的情况下进一步降低硬件资源消耗。推导了secp256k1曲线的快速模约简公式，基于比特重组方法与secp256r1和SCA-256两条曲线融合，设计了兼容3条曲线的快速模约简结构。最后根据所设计的运算单元和快速模约简结构提出了新的点加和倍点的操作调度方法，提高了乘法与快速模约简的并行性，从而获得了更高的吞吐量。实验结果表明，该标量乘法器的AT值比现有的大部分设计更低，此外在吞吐量方面也具有优越性。