椭圆曲线加密算法（ECC)示例（北太天元实现）
本文主要是使用北太天元演示椭圆曲线加密算法（ECC）的加密/解密过程，内容包括密钥、公钥生成，以及通过加密并解密一个简单数字的过程来描述其使用方法。我实际上就是把参考文献的文章中的代码移植到北太天元上来。这个移植工作非常容易，在最新版本的北太天元下几乎一行程序都不用改。 因此，我把原来的文章稍微精简一下，方便和我类似的人更容易理解吧。
按道理说，ECC加密还是挺难的，涉及到有限域(finite field, 又称为伽罗瓦域).但是实际上，我们只要这样想, 我们会取定一个素数p, 然后，我们所要考虑的数最后都会对p 取余数.
例如: 7 是一个素数，那么 3 和 10 就看成相等的，因为 3 除以7 所得余数 (北太天元用 mod(3,7) 来计算) 和 10 除以 7 所得余数是相等。 我们计算 3 * 5 之后得到了15， 然后回到7再取余数, 得到 1. 我们还可以定义 3 的逆， 因为 mod(3*5, 7) = 1, 因此 我们就说5是3的逆。
从而，我们还可以对 2/3 对 7 取余数进行定义了，因为 1/3 = 5 (mod 7)因此 2/3 = 2*5 (mod 7) 从而 2/3 对7 取余数就可以定义为mod(2*5, 7) , 计算结果得到3.
从而我们可以定义 2/3 对7的余数是 3.
椭圆曲线加密涉及到一个椭圆曲线
y^2 = x^3 + a*x + b
其中 a,b 是满足一定条件的常数.
工具
当指定x的值时，计算区间[0,p]内所有满足
y^2 = x^3 + a*x + b
的y值, 当然上面的等号成立也是在对p取余数的条件下。当然，这里我有一个疑问，为啥不是求{0,,1,...,p-1} 之内的数呢？对p取余数的话 0 == p (mod p) 和 0^2 == p^2 (mod p) 总是成立啊。不过，也许另有玄机啊，我这里保持了参考文献的做法。
这个工具功能很容易理解，通俗点说就是给定x求y,y的取值范围为[0,p]
%file:ECCCal.m
%a,b为椭圆参数，p为质数，x为给定的值
function [ y ] = ECCCal( a,b,p,x )
y=[];
mm = mod(x^3+a*x+b,p);
index = 1;
for yy = 0:1:p
if mod(yy^2,p) == mm
y(index)=yy;
index=index+1;
end
end
end
我们考虑 y^2 = x^3 +x +1 , 也就是椭圆曲线 y^2 = x^3 +a*x+b 中的 a = 1, b= 1
的情况。我们选定素数p =23, 计算x=3时的 y 值,
在北太天元的命令行输入
>> ECCCal(1,1,23,3)
我们会得到
ans =
10 13
我们可以得到可以得到两个结果 y =10和y = 13
（有多个结果是正常的不用在意，因为一元二次方程还有两个根，何况这是一个
更加奇怪的方程呢, 实际上如果不限定y的范围，理论上是有无数个解).
我们可以求出所有满足椭圆曲线方程
y^2 = x^3 + a*x + b
的解，
其中 a= 1, b= 1, 满足的意思是说方程左右两边在对p取余数之后是相等的，
自然，我们可以仅仅考虑x 和 y 都是属于 {0,1,....p} 这个集合的。 然后我把求的解作成一对(x,y), 作为平面上的点画出来。
计算(0,0)->(p,p)正方形范围内所有满足"椭圆曲线"的点
我们给出在北太天元上代码来完成这个工作, 使用这个工具可以方便我们将该范围内所有满足“椭圆曲线”画出来，
并且把数值也数出来。
%file:ECCPlot.m
%a,b分别为椭圆的参数，p为一个质数
function [x,y] = ECCPlot( a,b,p)
if (nargin == 0)
a = 1;
b = 1;
p = 23;
end
x=[];
y=[];
index = 1;
for xr = 0:1:p
mm = mod(xr^3+a*xr+b ,p);
for yr=0:1:p
if mod(yr^2,p) == mm
x(index)=xr;
y(index)=yr;
index = index+1;
end
end
end
plot(x,y,'.','MarkerSize',50)
hold on
grid on
end
假设椭圆曲线y^2 = x^3 + a*x+b 的 a= 1, b=1, 然后选取素素p =23,
在北太天元命令行输入
>>[x,y] = ECCPlot(1,1,23)
我们可以看到如下图片
因为椭圆曲线加密算法的特殊性我们需要自己实现几个操作方法分数的模运算
我们很容易知道 2对23取模的结果是2，但1/2对23的模是多少呢？北太天元暂时是没有这种计算的方法的，我们可以使用下面的modfrac函数来完成这个操作
%file:modfrac.m
% n 分子 d 分母 m 模数
function y = modfrac( n,d,m )
n=mod(n,m);
d=mod(d,m);
i=1;
while mod(d*i,m) ~=1
i=i+1;
end
y=mod(n*i,m);
end
我们还需要椭圆曲线加密而定义的一个特殊的加法, 具体规则此处不再赘述了，下面的北太天元代码得到点(resx, resy) 就是点(x1,y1)和(x2,y2)的和，当然这种和是用特殊的加法得到。
%file:Add.m
%a,b 椭圆参数 p 质数 x1,y1 第一个点的坐标 x2,y2 第二个点的坐标
function [ resx,resy ] = Add( a,b,p,x1,y1,x2,y2 )
if x1==x2 && y1==y2
k=modfrac(3*x1^2+a,2*y1,p);
resx = mod(k^2-x1-x2,p);
resy = mod(k*(x1-resx)-y1,p);
end
if x1==x2 && y1~=y2
resx = inf;
resy = inf;
end
if x1 ~= x2
k=modfrac(y2-y1,x2-x1,p);
resx = mod(k^2-x1-x2,p);
resy = mod(k*(x1-resx)-y1,p);
end
end
如两点P(3,10)，Q(9,7) 计算P+Q：我们在北太天元的命令行输入
[x,y]=Add(1,1,23,3,10,9,7)
得到结果是
x =
17
y =
20
也就是 P(3,10) + Q(9,7) 得到了 (17,20).
常量乘法
什么是常量乘法呢？实际上就是N（N为正整数）乘以点P。同时也需要说明的是下边北太天元的实现是来源于参考文献的MATLAB代码
%file:NP.m
%a,b 椭圆参数，p 质数，n表示 n个点P相加也就是n*P ,x,y 表示P点的横纵坐标
function [resx,resy] = NP( a,b,p,n,x,y )
if n ==1
resx = x;
resy = y;
return;
end
if n>=2
[xsub,ysub]=NP(a,b,p,n-1,x,y);
if xsub==Inf && ysub == Inf
resx=Inf;
resy=Inf;
else
[resx,resy]=Add(a,b,p,x,y,xsub,ysub);
end
end
end
如我们计算N=2, 点P(3,10), 计算N*P （特殊定义的加法把N个p相加的结果)
在北太天元上输入
>>[x,y]=NP(1,1,23,2,3,10)
x =
7
y =
12
N=2, P(3,10) 时N*P的结果是(7,12).
如果N=3, P(3,10) 相乘，
[x,y]=NP(1,1,23,3,3,10)
x =
19
y =
5
得到的结果(19,5).
加密解密过程
至此我们加密解密需要的所有算法模块就全部足够了。
流程以摘要中提到的第一篇文章中给的示例为例
Alice选定一条椭圆曲线E，并取椭圆曲线上一点作为基点G 假设选定的椭圆为a=4,b=20,p=29所示示的椭圆，基点G(13,23) , 基点G的阶数n=37（阶数的概念本文没提,你可以理解为在进行常量乘法运算的时候常量的最大值要小于n ）
Alice选择一个私有密钥k（k<n）并生成公开密钥K=kG 比如k=25, K= kG = 25G = (14,6）(使用我们的乘法运算函数去测试一下）
Alice将E和点K、G传给Bob （ 这没什么说的，E,K,G 共同组成了公钥，需要把公钥发给对方）
Bob收到信息后，将待传输的明文编码到上的一点M（编码方法略），并产生一个随机整数r（r<n,n为G的阶数） 假设r=6 要加密的信息为3,因为该信息也需要满足曲线方程E，所以我们很容易能够选取到一点(3,28)（其他的点如(3,1)也是可以的）
Bob计算点C1=M+rK和C2=rG
Bob将C1、C2传给Alice（容易理解，将加密信息传回给Bob，由Bob解密以知道Alice给他传的是啥）
Alice收到信息后，计算C1-kC2得到的结果应该是(3,28) （解密过程）
至此椭圆曲线加密、解密的整个流程就结束了。将其写成一个MATLAB脚本 如下：
%file:ECC.m
%演示曲线加密算法加/解密过程
a=4;
b=20;
p=29;
GX=13;
GY=23;
k=25;
[KX,KY]=NP(a,b,p,k,GX,GY);
r=6;
MX = 3;
MY = 28;
[rKX,rKY] = NP(a,b,p,r,KX,KY);
[C1X,C1Y]=Add(a,b,p,MX,MY,rKX,rKY);
[C2X,C2Y]=NP(a,b,p,r,GX,GY);
[kC2X,kC2Y]=NP(a,b,p,k,C2X,C2Y);
kC2Y=mod(-1*kC2Y,p);
[resx,resy]=Add(a,b,p,C1X,C1Y,kC2X,kC2Y)
在北太天元上执行m脚本ECC
>>ECC
resx =
3
resy =
28
我们执行一下可以看到输出结果 resx=3,resy=28 。
这里的3正是Bob给Alice传的数字内容，其他人因为没有私钥k所以无法像Alice
那样简单地还原出数字3
总结
椭圆加密算法看似简单, 仅仅通过一个加法运算和一个常量乘法运算便完成了加密解密算法, 当然背后的原来涉及到数学知识这里没有讨论。这个文章仅仅是起演示功能。
另外, 这里的代码是从参考文献的MATLAB代码移植到北太天元上的，一行不改也是可以
运行的。我修改了ECCPlot.m这么一行
plot(x,y,'*')
成了下面这么一行，这样画出的点更大
plot(x,y,'.','MarkerSize',50)
原来的代码在北太天元上画的图是
参考: https://blog.csdn.net/alphags/article/details/79660819
