[mathjax]Weierstrass標準形を併用しながらスカラー倍，点の加算にMontgomery curvesを使った話．
以前はPARI GPの機能を使って，Weierstrass標準形でECDSAを実装した．今回は新しくMontgomery curvesを用いたスカラー倍，点の加算を実装した．そのため，以前行ったWeierstrass標準形での実装と，新しくつくるMontgomery curvesでの実装を見比べてみる．
1.スカラー倍
作成したMontgomery ladderのラッパー関数としてscalarという関数を作成した．入力は楕円曲線のパラメータ\(A,B\)とアフィン座標の点\(P\)，スカラー\(k\)，有限体の位数\(p\)である．出力はアフィン座標の点\(Q=kP\)である．
scalar(A,B,P,k,p)={
local(P_M,pair,kP_M,k1P_M,x_k,x_k1,y_k);
P_M=[P[1],1];
pair=ladder(A,P_M,k,p);
kP_M=pair[2];k1P_M=pair[1];
if(kP_M!=[0]&&kP_M[2]==0,kP_M=[0]);
if(k1P_M!=[0]&&k1P_M[2]==0,k1P_M=[0]);
if(kP_M==[0],
return([0]);
);
if(kP_M!=[0],
x_k=Mod(kP_M[1]/kP_M[2],p);
if(k1P_M==[0]||k1P_M[2]==0,
x_k1=-1;,
x_k1=Mod(k1P_M[1]/k1P_M[2],p);
);
y_k=recover(A,B,P_M[1],P[2],x_k,x_k1,p);
);
return([x_k,y_k]);
}
2.点の加算
署名検証で点の加算が必要となるため，点の加算を実装する必要がある．アルゴリズムをそのままコードにして実装する．ただし，\(P_3=P_1+P_2\)とする際に\(P_1=-P_2\)となっている場合は\(P-P=O\)となるため無限遠点を返す．また，\(P_1=P_2\)の場合は\(P+P=2P\)となるため2倍算を行う．
mdbl(A,B,P,p)={
local(x,y);
x=Mod(((P[1]^2-1)^2)/(4*B*(P[2]^2)),p);
y=Mod(((2*P[1]+P[1]+A)*(3*(P[1]^2)+2*A*P[1]+1))/(2*B*P[2])-(B*(3*(((P[1]^2)+2*A*P[1]+1)^3)))/((2*B*P[2])^3)-P[2],p);
return([x,y]);
}
madd(A,B,P1,P2,p)={
local(x,y);
if(P1==[0],return(P2));
if(P2==[0],return(P1));
if(P1[1]==P2[1],
if(P1[2]==-P2[2],
return([0]);,
return(mdbl(A,B,P1,p));
);
);
x=Mod((B*((P2[2]-P1[2])^2))/((P2[1]-P1[1])^2)-A-P1[1]-P2[1],p);
y=Mod(((2*P1[1]+P2[1]+A)*(P2[2]-P1[2]))/(P2[1]-P1[1])-(B*((P2[2]-P1[2])^3))/((P2[1]-P1[1])^3)-P1[2],p);
return([x,y]);
}
3.動作確認
動作確認には\(F_{30011}\)上のMontgomery curves \(y^2=x^3+84x^2+x\)とWeierstrass標準形 \(y^2=x^3-2351x+43876\)を用いる．ベースポイントについてはWeierstrass標準形では\((7796,9454)\)，Mongtomery curvesでは\((7768,9454)\)を用いる．なお，点位数\(l=1249\)はとなる．基本的には以前作成したWeierstrass標準形のコードと同時にMontgomery curvesのコードを実行する．以下にコードと実行例を示す．
実行例からWeierstrass標準形の点\(P_W=(x_W,y_W)\)とMontgomery curvesの点\(P_M=(x_M,y_M)\)について\((x_M,y_M)=(x_W-28,y_W)\)が成り立っており，スカラー倍や点の加算が正しく行えていることが分かる．
\\m-n==1
add(P1,Pm,Pn,p)={
local(a,b,c,d,Xmn,Zmn);
if(Pm==0||Pm[2]==0,return(Pn));
if(Pn==0||Pn[2]==0,return(Pm));
a=Pm[1]+Pm[2];
b=Pm[1]-Pm[2];
c=Pn[1]+Pn[2];
d=Pn[1]-Pn[2];
Xmn=Mod(P1[2]*((b*c+a*d)^2),p);
Zmn=Mod(P1[1]*((b*c-a*d)^2),p);
return([Xmn,Zmn]);
}
dbl(A,P,p)={
local(tmp1,tmp2,tmp3,X2,Z2);
if(P==0||P[2]==0,return([0]));
tmp1=(P[1]+P[2])^2;
tmp2=(P[1]-P[2])^2;
tmp3=tmp1-tmp2;
X2=Mod(tmp1*tmp2,p);
Z2=Mod(tmp3*(tmp2+((tmp3*(A+2))/4)),p);
return([X2,Z2]);
}
recover(A,B,x1,y1,xn,xn1,p)={
local(up,down,res);
if(xn1==-1,
return(-y1);
);
up=(x1*xn+1)*(x1+xn+2*A)-2*A-((x1-xn)^2)*xn1;
down=2*B*y1;
if(down==0,
res=0;,
res=Mod(up/down,p);
);
return(res);
}
ladder(A,P,k,p)={
local(R0,R1,i,bin);
bin=binary(k);
R0=[0];
R1=P;
for(i=1,length(bin),
if(bin[i]==0,
R1=add(P,R1,R0,p);
R0=dbl(A,R0,p);,
R0=add(P,R1,R0,p);
R1=dbl(A,R1,p);
);
);
return([R1,R0]);
}
\\input and output is affine point P=(x,y)
scalar(A,B,P,k,p)={
local(P_M,pair,kP_M,k1P_M,x_k,x_k1,y_k);
P_M=[P[1],1];
pair=ladder(A,P_M,k,p);
kP_M=pair[2];k1P_M=pair[1];
if(kP_M!=[0]&&kP_M[2]==0,kP_M=[0]);
if(k1P_M!=[0]&&k1P_M[2]==0,k1P_M=[0]);
if(kP_M==[0],
return([0]);
);
if(kP_M!=[0],
x_k=Mod(kP_M[1]/kP_M[2],p);
if(k1P_M==[0]||k1P_M[2]==0,
x_k1=-1;,
x_k1=Mod(k1P_M[1]/k1P_M[2],p);
);
y_k=recover(A,B,P_M[1],P[2],x_k,x_k1,p);
);
return([x_k,y_k]);
}
mdbl(A,B,P,p)={
local(x,y);
x=Mod(((P[1]^2-1)^2)/(4*B*(P[2]^2)),p);
y=Mod(((2*P[1]+P[1]+A)*(3*(P[1]^2)+2*A*P[1]+1))/(2*B*P[2])-(B*(3*(((P[1]^2)+2*A*P[1]+1)^3)))/((2*B*P[2])^3)-P[2],p);
return([x,y]);
}
madd(A,B,P1,P2,p)={
local(x,y);
if(P1==[0],return(P2));
if(P2==[0],return(P1));
if(P1[1]==P2[1],
if(P1[2]==-P2[2],
return([0]);,
return(mdbl(A,B,P1,p));
);
);
x=Mod((B*((P2[2]-P1[2])^2))/((P2[1]-P1[1])^2)-A-P1[1]-P2[1],p);
y=Mod(((2*P1[1]+P2[1]+A)*(P2[2]-P1[2]))/(P2[1]-P1[1])-(B*((P2[2]-P1[2])^3))/((P2[1]-P1[1])^3)-P1[2],p);
return([x,y]);
}
{
\\準備
a=-2351;b=43876;A=84;B=1;p=30011;
e=ellinit([0,0,0,a,b],p);
P_W=[Mod(7796,p),Mod(9454,p)];
P=[P_W[1]-28,P_W[2]];
l=ellorder(e,P_W);
print("pripare");
print("\ty^2=x^3+",a,"x+",b," on F_",p);
print("\tbase point P_W=",P_W);
print("\tbase point P =",P);
print("\torder l=",l,"\n");
\\鍵生成
d_B=random(l-2)+2;
P_WB=ellmul(e,P_W,d_B);
P_B=scalar(A,B,P,d_B,p);
print("key generation");
print("\tsecret key d_B =",d_B);
print("\tpublic key P_B =",P_B);
print("\tpublic key P_WB=",P_WB,"\n");
\\署名生成
until(u!=0&&v!=0,
r=random(l-1)+1;
U_W=ellmul(e,P_W,r);
U=scalar(A,B,P,r,p);
m=random(p);
u=lift(lift(U[1])%l);
v=((m+u*d_B)/r)%l;
);
print("signature generation");
print("\trandom value r=",r);
print("\tU =",U);
print("\tU_W=",U_W);
print("\tm=",m);
print("\tsignature (u,v)=(",u,",",v,")","\n");
\\署名検証
d=(1/v)%l;
V_W=elladd(e,ellmul(e,P_W,d*m),ellmul(e,P_WB,d*u));
V=madd(A,B,scalar(A,B,P,d*m,p),scalar(A,B,P_B,d*u,p),p);
print("signature verification");
print("\td=",d);
print("\tV =",V);
print("\tV_W=",V_W,"\n");
\\結果表示
print("result");
print("\t(u,Vx)=(",u,",",lift(V[1])%l,")");
if(u==(lift(V[1])%l),
print("\tOK");,
print("\tNG");
);
}
(17:34) gp > \r m_ecdsa.gp
pripare
y^2=x^3+-2351x+43876 on F_30011
base point P_W=[Mod(7796, 30011), Mod(9454, 30011)]
base point P =[Mod(7768, 30011), Mod(9454, 30011)]
order l=1249
key generation
secret key d_B =1094
public key P_B =[Mod(12663, 30011), Mod(3408, 30011)]
public key P_WB=[Mod(12691, 30011), Mod(3408, 30011)]
signature generation
random value r=1112
U =[Mod(12043, 30011), Mod(16108, 30011)]
U_W=[Mod(12071, 30011), Mod(16108, 30011)]
m=18826
signature (u,v)=(802,305)
signature verification
d=86
V =[Mod(12043, 30011), Mod(16108, 30011)]
V_W=[Mod(12071, 30011), Mod(16108, 30011)]
result
(u,Vx)=(802,802)
OK