一种基于余数系统的椭圆曲线密码算法

一种基于余数系统的椭圆曲线密码算法

1椭圆曲线点乘的计算椭圆曲线的密码是一个重要的函数密码，其基本计算方法是椭圆曲线的点幂。定义在NIST素域上的椭圆曲线点乘较易计算,而在一般的素域上计算椭圆曲线点乘则较为困难。其原因在于:前者的阶是特殊模数,模约简较为方便。目前,有若干种方法实现素域椭圆曲线点乘。一般来说,首先需要构建一个多精度的模乘器,然后顺序计算椭圆曲线点乘中的模乘运算。在文献中,利用商流水线搭建了一个大的蒙哥马利模乘器,进而计算椭圆曲线点乘。在文献中,则基于心动阵列构建了一个蒙哥马利模乘器,并用它来计算全部模乘。另一种计算方法利用了余数系统[7~9]。在文献[7~9]中,椭圆曲线点乘分成一些列的倍点和点加。每次倍点和点加前,先将坐标转化到余数系统中,将它们所对应的模乘在余数系统中用乘法计算。这种方法仅用到了一个余数系统,但是其动态范围很大。然后,将上述计算结果转换到二进制系统中,利用Horner法则进行模约简。约简完成后,再将结果转换到余数系统中,以进行后续倍点和点加。这种策略充分利用了余数系统的并行性,却包含了余数系统和二进制间的频繁转换,因而影响了它的实施效率。此外,在文献[10,11]中,同样应用余数系统来计算椭圆曲线点乘,但是余数系统和二进制间的相互转换只在首尾各计算一次。该方法的特点是,利用余数系统下的蒙哥马利模乘[12~17]计算所有的模乘。而余数系统蒙哥马利模乘中,通常用到两个近乎对称的余数系统,它们的模数个数和动态范围都很接近。本文提出了一种新的余数系统下快速计算素域椭圆曲线点乘的方法,它充分汲取了上述两种余数系统算法的优点,而避免了相应的缺点。首先,本文采用了非对称的余数系统,其中由两个模数{2L,2L-1}组成余数系统Γ的模数组,另外八个模数组成余数系统Γ的模数组,它们都是广义梅森数(GeneralizedMersennenumbers)2L-2Ki+1(i=1,2,…,8)。这些特殊模数十分方便模约简运算。其次,利用该非对称的余数系统,椭圆曲线中每次点加、倍点的模乘绝大部分在余数系统中通过乘法实现,仅有四个模乘利用余数系统下的蒙哥马利模乘完成。后者在运算中进行了模约简,将之前乘法所增大的中间值压缩到了开始时较小的范围,避免了转换到二进制的步骤,使得点加和倍点可以在余数系统中连续计算。与此相应,在计算椭圆曲线点乘时,用到了蒙哥马利梯子,它每次仅计算x坐标,点加和倍点可以并行运算。上述方法结合在一起,大幅减少了素域椭圆曲线点乘的周期数,加快了运算速度。2动态范围与乘法余数系统是一种无权重的整数表示系统,它由一组互质的模数定义:{m1,m2,…,mn}。余数系统所能表示的整数范围[0,M-1]称为动态范围,其中,M=∑ni=1mi。在余数系统中,加法、减法、乘法可以并行运算。设A=(a1,a2,…,an),B=(b1,b2,…,bn),用ue0c9表示+、-或者×,则有:2.1余数ls-ls设余数系统Ω的模数组为{m1,m2,…,mn},动态范围是[0,M-1],M为模数的连乘积,且Mi=M/mi。又设整数R处在[0,M-1]内,它的余数表示为(r1,r2,…,rn),其中,ri=Rmodmi,i=1,2,…,n。此时,中国剩余定理表示如下:在文献中,式(1)中最后模M的操作被减去M的若干倍数代替。设这个倍数因子为α,则它必定处在[0,n-1]内,且有其中,σi=(ri·Mi-1)modmi。为了计算α,需要额外的信息,它可以通过冗余的模数mn+1及相应的余数rn+1实现,也可以利用近似法求出。2.2蒙哥马利模乘算法1的运算余数系统下的蒙哥马利模乘算法,可以通过在两个余数系统间进行基底扩展实现,如算法1所示。算法1余数系统下的蒙哥马利模乘算法1的主要运算发生在步骤2和步骤4,前者是从余数系统Ω向余数系统Γ扩展基底,后者刚好相反。以步骤2为例,其运算可以用向量表示如下:其中,σi=ai·bi·(-pi-1)modmi,Mi,j=(M/mi)modmn+j,μj=Mmodmn+j。第j个分量的运算需要模mn+j。3椭圆曲线点乘通常应用的素域椭圆曲线方程形式为y2=x3+ax+b,它的运算定义在FP中,且有4a3+27b2≠0(modp),p>3。可以利用蒙哥马利梯子计算椭圆曲线点乘,如算法2所示。算法2利用蒙哥马利梯子计算椭圆曲线点乘在算法2中,每一行的运算都并行执行,运算过程中始终保持V-U=G。其中的2倍表示倍点,加法表示点加。如果U、V、V-U的x坐标都已知,则2U、U+V的x坐标可以仅用上述三点的x坐标表示出来。其y坐标可以用U-V的x坐标以及U和V点的椭圆曲线方程消去。在文献[21,23]中,进一步将仿射坐标x用投影坐标表示,即xue039X/Z,以避免椭圆曲线点加和倍点中模逆。椭圆曲线点乘中某一步中U+V和2U的公式如下:上面的表达式中,x是G点的横坐标,下标u、v分别表示U点和V点。此外,只需要将2U的公式下标中u置换为v,即可得到2V的计算公式。参数a、b分别是椭圆曲线方程中一次项的系数和常数项。4新的方法是使用余数系统快速计算椭圆曲线的点幂4.1进位信号cf的产生本文用到了两个余数系统Ω和Γ。Ω的模数组为{2L-2K1+1,…,2L-2K8+1},Γ的模数组为{2L,2L-1}。其中,1≤Ki<L/2。与算法1比较,这两个余数系统的设定显著不同。这里Γ是基本的余数系统,而Ω是引入的冗余系统;前者仅有两个模数,以减少余数系统蒙哥马利模乘中的乘法的累积次数;后者有八个模数,以保证两个余数系统合在一起有足够的动态范围。形如2L-2Ki+1的数通常称为广义的梅森数,对它的模算术比较简单。设mi=2L-2Ki+1=2L-2k+1,a和b是两个L位的二进制数,且T=a·b。则根据文献中的讨论,有:T=(tL-1..0)2=2L·Th+Tl+cf·2L,Th1=(tL-1..L-k)2<2k.Th2=(tL-k-1..0)2<2L-k,即Th=2L-kTh1+Th2。进而:T≡C″(modmi),其中C″=Tl+(2k-1-2L-k)Th1+(2k-1)Th2+cf(2k-1)。进位信号cf的产生过程如图1所示。图1中将特殊模数模乘器的乘法和模约简综合考虑,并不计算完整的乘积,而是得到两个L位的二进制数以及一个中间的进位信号。它避免了一个2L位的长进位。上述乘法的关键路径长度为:其中,下标‘mult’表示乘法,‘FA’表示全加器(FullAdder)。图1中的乘积即为C″,需要进一步约简。它的表达式中合并为六项,可以记为C″=T1+T2+…+T6。如图2所示,其关键路径为:在余数系统Γ中,只需要对{2L,2L-1}进行模约简,更为简单。Tmod2L=Tl,Tmod(2L-1)=(Th+Tl+cf)mod(2L-1)。设c′f·2L+S=Th+Tl+cf,则Tmod(2L-1)=S+c′f。4.2算法3的并行步骤算法3非对称余数系统中的蒙哥马利模乘算法3采用了大量预计算,以节省串行的操作步骤。可以利用5n(即10)个处理器核并行执行,如图3所示。一般情况下,左边八个处理器核处理余数系统Ω中的运算,右边两个处理器核则处理余数系统Γ中的运算。但是,在算法3的步骤7~步骤9中,左边八个处理器核其实是用来计算Γ中的模乘,n步完成。算法的步骤1由左侧八个处理器计算,步骤2由右侧两个处理器计算;因此这两步可以并行执行。类似地,步骤3和步骤4、步骤17和步骤18都可以并行执行。步骤10可以和步骤7~步骤9并行,在右侧处理器核中完成。因此,算法3总的串行步骤为e=2n+3=7。图3中的buffer表示缓冲器,可以由寄存器、锁存器、存储器等构成。4.3蒙哥马利模乘在本文的非对称余数系统中,椭圆曲线点乘中大部分的模乘可以直接用乘法代替,少部分的余数系统蒙哥马利模乘用来压缩输出值,保证结果不溢出余数系统动态范围。据此,式(4)所示的椭圆曲线并行倍点、点加公式重写如下:在图4中,假设余数系统Γ扩展k个模数,使得其动态范围任意大,则经过算法3中的余数系统蒙哥马利模乘之后,其输出值为(c4n+1,c4n+2,…,c4n+k)。它的前两个余数与未扩展的值相同。由于上述两个值就能够确定输出值,因此其它模数是多余的,即不必扩充余数系统Γ。因此,可以按照式(7)同时计算椭圆曲线点乘中的点加和倍点,其具体步骤如式(8)所示。其中,符号·表示该行当前进行的乘法所在位置,其它连写在一起的变量则表示已经求出的乘积。而仅压缩了Y输出范围,并不改变它在Fp中的对应值。根据式(8),计算一次并行的点加和倍点,需要23+4×(7-1)=47次串行的乘法。4.4x,l-1,x转化在利用余数系统进行并行计算之前,需要将输入坐标值由二进制转换为余数表示。求得最终结果后,还需要将它转换回到二进制。假设x需要从二进制转换到余数表示,根据前述模数组的选择,0≤x<22L。可以直接在图1中令,Tl=xmod2L,cf=0。余数表示到二进制表示的转换可以在Γ中完成。不妨设最后所得的结果在Γ中的表示为(x1,x2),对应模数组{2L,2L-1}。根据中国剩余定理,它所对应的二进制数x为:当x2≥x1时,x=2L(x2-x1)+x1。当x2<x1时,x=2L[(2L-1-x1)+x2]+x1。其中,(2L-1-x1)可以用对x1的二进制表示求反得到。4.5模装置的组合利用式(8)中的投影坐标公式,可以在余数系统中计算椭圆曲线点乘。其结果需要转换回仿射坐标,此时需要对Z求逆。可以利用费马小定理,即Z-1modP=ZP-2modP。如此,对于N位的椭圆曲线点乘,则需要计算一个N位的模幂。可以将幂指数(P-2)的二进制表示平均分为N/2个部分,则每个部分中最多有两次模平方、两次模乘。在本文的余数系统中,上述运算可以用四次连续的乘法完成,然后再与|M|P作一次蒙哥马利模乘即可。这样一来,连续进行约(N-1)/2次上述操作,即可完成模幂。通常N为偶数,则总共需要的乘法次数约为:(N/2)·(7+4)=11N/2。5蒙哥马利土尼并无意义的椭圆点乘大体而言,假设按照传统余数系统蒙哥马利方法进行设置,在单个模数大小相当的情形下,它只需要四个模数。本文扩展模数组到10个模数,若不考虑每步约简的运算量区别,相当于计算量增长了1.5倍;但是,实际每次模乘约简的消耗远远小于原值,故实际运算大为提速。具体而言,由于本文采用了蒙哥马利梯子,N位的素域椭圆曲线点乘由N次并行的点加和倍点以及最后从投影坐标到仿射坐标的转换完成。文献[10,21]同样采用蒙哥马利梯子的方法计算素域椭圆曲线点乘,表1比较了本文和它们在并行点加、倍点所需要的串行乘法次数。从表1可见,在采用蒙哥马利梯子进行并行倍点和点加时,本文比参考文献[10,21]所需要的串行乘法次数大为减少。例如,在参考文献中,对于256位的椭圆曲线点乘,取n=15,此时,本文的步骤数不到它的1/9。总体来看,如果采用费马小定理将结果从X/Z→x,本文中要完成一次N位的椭圆曲线点乘,需在余数系统中计算的串行乘法次数约为:余数系统中的乘法其实是对模数组的模乘。如图1和图2所示,它可以分为两级实现,则一级的周期为Iperiod=max{Imult,Ireduc}。当L较大时,有Imult>Ireduc,则从式(8)可知:如果采用两级流水,注意到算法3的步骤7和步骤14会有数据相关性问题。式(11)中前三次出现ue3c1的位置处,可以将其前出现的两次乘法插入到该次ue3c1的步骤7、步骤14分别执行。第四次出现ue3c1的位置处会增加三个周期。则总的周期数变为:则完成一次N位的素域椭圆曲线点乘需要的时间约为:以256为素域椭圆曲线点乘为例,N=256,可以选择L=136,则:通常,周期长度和具体的实现方式有关。如果假设(I68-bitmult+IFA+I68-bitFA)=10ns,则总的计算时间约为140μs。5.1存储于保护作用下的存储函数上述算法若采用VLSI实现,需要的硬件包括:(1)10个L位的模乘器。相当于约40个L/2位的乘法器,以及(6×8+2)=50个L位的加法器。(2)55nL=110L位的寄存器,包括六组5nL位的寄存器用于保存式(8)中的变量,15nL位用于模乘器,10nL位用于模加-模减器。其中,式(8)中的变量可以按如下方式分别访问各自的寄存器:(3)约(8n2+12n)L=56L位的ROM。用于存储余数系统蒙哥马利模乘中的预计算值以及椭圆曲线阶P的余数表示,以及|M|P的余数表示。(4)nL=2L位的RAM,存储出发点G的横坐标。(5)控制逻辑,包括式(8)中的状态选择、寄存器复用控制、模乘器复用控制、模加-模减器复用控制等。6生成多次蒙哥马利模乘本文提出了一种新的在余数系统下快速计算素域椭圆曲线点乘的方法。我们所选取的余数系统由两个不对称的余数系统构成,其中一个包含两个模数{2L,2L-1},另