椭圆曲线加法运算的数据仿真

刘增芳，韦性佳，芦殿军

(青海师范大学 数学与统计学院，青海 西宁，810008)

随着信息技术的高速发展，信息安全已经成为全社会所面临的重大挑战，椭圆曲线密码体制(Elliptic curve cryptosystem，简称ECC)已经广泛的应用于各类密码学方案的构造之中，作为保障信息安全性的核心工具之一，具有十分重要的研究价值。1985年Neal Koblitz和Victor Miller首次在密码学的研究之中启用了椭圆曲线公钥密码体制[1-2]。经过研究发现，ECC具有如下两方面优势：1)密钥具有较小的存储空间与运算量，这就极大的简化了密码体制的运算成本，具有更高的实用性。2)可定义群之间基于Weil对或是Tate对的双线性映射[3-4]，基于此ECC已成为理想的公钥密码体制之一。3)椭圆曲线的加法运算，在很大程度上加快了密码体制的运算速度。目前对于ECC的硬件实现已成为许多学者研究的重要课题之一。2013年，宋春玉[5]利用matlab软件，编译了椭圆曲线密码系统的基本运算与原理，并且第一次将其应用于汉字的加密和解密之中。2016年，谢天艺等人[6]利用ECC硬件加速器实现了硬件可配的点乘和素数检测。同年，王千喜等人[7]对于ECC算法的硬件的实现进行了优化，使得在主时钟为300M、密钥长度是160bit的条件下，ECC加密解密的效率提升至大约5k/s，这极大的提高了目前ECC加密和解密的速度。

文中讨论的是基于有限域上椭圆曲线的群中加法运算的实现方案。椭圆曲线加法运算的数据模拟并不是一件简单的事情，需要做如下几项工作：第一，判断两数是否互素。这个工作是为下一步模逆运算做准备的，因为只有两个互素的数才可求逆，否则不成立。第二，实现有限域上的模逆运算，这也是椭圆曲线加法群运算中的关键性算法。第三，对有限域上的椭圆曲线利用Euler准则测试给定一个x对应的y值是否是一个二次剩余，如果是二次剩余，文中将给出相应的点，并输出点的个数。第四，实现椭圆曲线加法运算。

1 预备知识

1.1 有限域Fp上的椭圆曲线

定义1 令p>3是一个奇素数。有限域Fp上的椭圆曲线E可以定义为等式[8]

y2=x3+ax+b，a,b∈Fp

并且4a3+27b2≠0(modp)，集合E(Fp)包括所有的点(x,y)，x∈Fp，y∈Fp，且点坐标(x,y)满足等式y2=x3+ax+b,集合E(Fp)还有一个特殊的点称作无穷远点。

1.2 Euclidean算法

给出两数a,b∈Z+，利用Euclidean算法求两数的最大公因子。令r0=a,r1=b，ri(i=2,…,m)：ri/ri+1所得余数，qj(j=1,2,…,m)：ri/ri+1所得商。具体算法执行如下(即辗转相除法)：

容易看到

gcd(r0,r1)=gcd(r1,r2)=…=gcd(rm-1,rm)=rm

因此，可以得到gcd(a,b)=rm。

扩展的Euclidean算法，它以两个整数a和b作为输入，计算出整数r,s和t使得r=gcd(a,b)且sa+tb=r。若r=1即gcd(a,b)=1，则有b-1moda=tmoda。

1.3 二次剩余

设E是Zp上的椭圆曲线y2=x3+ax+b。首先确定E的点。这可以通过对每个x∈Zp，计算x3+ax+b(modp)解方程

y2≡x3+ax+b(modp)

求y。对于给定的x，可用Euler准则测试是否z=x3+ax+b是一个二次剩余。对于素数p≡3(mod4)，利用公式z(p-1)/2modp判断z是否为一个二次剩余。若此值等于1，则说明z是一个二次剩余，否则z不是一个二次剩余。二次剩余z的平方根是±z(p+1)/4modp。

1.4 椭圆曲线的群的加法运算

椭圆曲线上的加法运算定义如下(这里的所有运算都在Zp中)。假设P=(x1,y1)并且Q=(x2,y2)是E上的点。

1)如果x1=x2且y1=-y2，则P+Q=O(无穷远点)

2)否则P+Q=(x3,y3)，其中

并且对λ有如下定义：

3)最后，对于所有的P∈E，有P+O=O+P=P

需要注意的是，Zp上的椭圆曲线没有像实数域上的椭圆曲线那样直观的几何解释。然而，定义加法运算亦可用同样的公式，(E,+)仍是一个阿贝尔群。

2 模素数的椭圆曲线加法

2.1 判断两数互素

利用Euclidean算法计算两数的最大公约数，若最大公约数为1，则说明两数互素。C语言程序如下：

#include

void main()

{int a,b,m;printf("please input two number:");

{scanf("%d%d",&a,&b);if(a==0‖b==0)

{printf("wrong input ");}

else{m=a%b;while(m!=0)

{ a=b;b=m;m=a%b;}

printf("%d ",b);}}}

2.2 模逆运算

设对于任意的a,b∈Zn，且gcd(a,b)=1。因为模逆运算所得到的值不唯一，文中取其最小值，利用扩展Euclidean算法可实现模逆运算。基本思想如下(分情况论述)：

1)当b>0时，计算模数a的正整数倍：1×a,2×a,3×a,…,n×a，将所得到的模数a的正整数倍加1，即1×a+1,2×a+1,3×a+1,…,n×a+1 。采用试除法，将模数a的倍数加1的每一项都除以b，直到余数为0即正好整除，则所得到的商即为b-1moda的值。

2)当b<0时，计算模数a的负整数倍：(-1)×a,(-2)×a,(-3)×a,…,(-n)×a，将所得到的模数a的负整数倍加1，即(-1)×a+1,(-2)×a+1,(-3)×a+1,…,(-n)×a+1。利用试除法，将模数a的倍数加1的每一项都除以b，直到余数为0即正好整除，这时将所得到的商是负数，因此文中再将商模a得到正值，故而此时的正值即为b-1moda的值。

根据上述算法，相应的C语言程序如下所示：

#include

#include

int gcd(int a,int b);

void main()

{int a,b,i,s;printf("input two number:");scanf("%d%d",&a,&b);a=a%b;

if(gcd(a,b)==1){for(i=1;i<100000;i++)

{if(a>0){if((i*b+1)%a==0) s=(i*b+1)/a;s=s%b;}

else{if(((-1)*i*b+1)%a==0)

{s=((-1)*i*b+1)/a;s=s+b;s=(b+s%b)%b;}}}

printf("%d ",s);}

else {printf("waring ! waring ! wrong input: ");}

int gcd(int a,int b){int m;if(a>0) {m=a%b;while(m!=0) {a=b;b=m;m=a%b;} return b;}

else{a=a+b;m=a%b;while(m!=0) {a=b; b=m;m=a%b;} return b;}}

2.3 判断二次剩余

判断二次剩余的目的是为了选取生成元。对于椭圆曲线y2=x3+ax+b，令a=-1，b=8即y2=x3-x+8，如图1所示：

图1 y2=x3-x+8的函数图像Fig.1 The functional image of y2=x3-x+8

依据模数要求选取模数p为11。相应的C语言程序如下：

#include

#include

void main()

{int z,x,p,y1,y2,k=0;scanf("%d",&p);

if(p%4!=3) {printf("wrong input，please change the number: ");scanf("%d",&p);}

else{for(x=0;x<=p-1;x++) {z=(x*x*x-x+8)%p;if((int)(pow(z,(p-1)/2))%p==1)

{y1=(int)(pow(z,(p+1)/4))%p; y2=(p-(int)(pow(z,(p+1)/4))%p)%p;

printf("x=%d,y1=%d,y2=%d ",x,y1,y2);k=k+2;}}printf("k=%d ",k);}}

运行此程序，文中得到了相应的数据，如表1所示：

表1Z11上的椭圆曲线y2=x3-x+8的点

Table 1 The point of the elliptic curvey2=x3-x+8inZ11

xx3-x+8(mod11)是二次剩余吗y08否18否23是5,6310否42否57否69是3,873是5,686否92否108否

由上表可知产生6个点，还有一个点是无穷远点，所以E共有7个点。

2.4 椭圆曲线加法运算

对于模素数椭圆曲线y2≡x3-x+8(mod 11)加法运算的C语言程序如下：

#include

#include

int ECCadd(int x1,int y1,int x2,int y2,int k，int p)

{int x3,y3,r,m,i;for(i=2;i<=k;i++) {if(x1==x2&&y1==-y2) {printf("无穷远点 ");}

else if(x1==x2&&y1==y2)

{m=invert(2*y1,p);r=((3*x1*x1-1)*m)%p;x3=(p*p+(r*r-x1-x2))%p;

y3=(p*p+(r*(x1-x3)-y1))%p;printf("%da=(%d,%d) ",i,x3,y3);}

else{m=invert(x2-x1,p);r=(p*p+((y2-y1)*m)%p)%p;x3=(p*p+(r*r-x1-x2)%p)%p;

y3=(p*p+(r*(x1-x3)-y1)%p)%p;printf("%da=(%d,%d) ",i,x3,y3);} x2=x3;y2=y3;}}

int sy(int p) {int z,x,y1,y2,k=0;if(p%4!=3) {printf("please change the number: ");}

else{for(x=0;x<=p-1;x++) {z=(x*x*x-x+8)%p;if((int)(pow(z,(p-1)/2))%p==1)

{y1=(int)(pow(z,(p+1)/4))%p;y2=(p-(int)(pow(z,(p+1)/4))%p)%p;

printf("x=%d,y1=%d,y2=%d ",x,y1,y2);k=k+2;}} printf("k=%d ",k);}}

void main()

{int x1,y1,t,u,p,z,x,k=0;printf("p=");scanf("%d",&p);sy(p);

for(x=0;x<=p-1;x++) {z=(x*x*x-x+8)%p;if((int)(pow(z,(p-1)/2))%p==1){ k=k+2;}}

printf("x1=");scanf("%d",&x1);printf("y1=");scanf("%d",&y1);

printf("a=(%d,%d) ",x1,y1);t=x1;u=y1;ECCadd(x1,y1,t,u,k,p);}

运行程序后得到的结果如下：由判断二次剩余的表可知，任意选取一对生成元α=(2,5)。则相应的结果为：

α=(2,5) 2α=(7,6) 3α=(6,3)

4α=(6,8) 5α=(7,5) 6α=(2,6)

由此可以看出，α=(2,5)确是本原元。因此，完成了对椭圆曲线y2=x3-x+8，且模素数11的数据仿真。

3 总结

文中运用C语言完成了判断两数互素，利用扩展的Euclidean算法实现有限域上的模逆运算，对模素数椭圆曲线y2≡x3-x+8(mod 11)利用Euler准则判断二次剩余且给出列表。以及任意选取一个生成元，实现有限域上椭圆曲线的群的加法运算，并给出相应的结果。

参考文献:

[1]KOBLITZ N.Elliptic curve cryptosystems[J].Mathematics of Computation，1987，48(177):203-209.

[2]MILLER VS.Use of elliptic curve in cryptography[J].Lecture Notes in Computer Science,1986,218(01):417-426.

[3]MENEZES A，OKAMOTO T，VANSTONE S.Reducing elliptic curve logarithms to logarithms in a finite field[J].IEEE Transaction on Information Theory,1993，39(05):1639-1646.

[4]CHEN L，HARRISION K，MOSS A，et al.Certification of public keys within an identity based system[C]//Lecture Notes in Computer Science.5th International Conference on Information Security.Berlin Germany:Springer-Verlag，2002,2433:322-333.

[5]宋春玉.椭圆曲线密码系统的matlab实现[J].西北民族大学学报(自然科学版),2013,34(02):14-17.

[6]谢天艺,黄凯,修思文,等.素数域椭圆曲线密码加速器的VLSI实现[J].计算机工程与应用,2016,52(01):89-94.

[7]王千喜,刘海法,王占厚,等.椭圆曲线密码(ECC)算法的一种硬件实现方案[J].信息通信,2016,(08):87-88.

[8]户占良.椭圆曲线密码体制的研究与应用[J].山西师范大学学报(自然科学版)，2010，24(03)：13-17.