I found how to convert from jacobian to standard coordinates, and I made some
minor modification to Cliff's code, there seems to be a minor bug when
repeatedly adding points when Z not equal to 1. So, it defeats the purpose of
using alternate coords when they need to be transformed on each call. (See
below for why bug is minor)
stdPoint =: 4 : 'x | (xx, yy)* x | 2 3 ^~ x invmod z [ ''xx yy z'' =. y '"0 1
NB. Ellitpic curve add (mod proj coord)
NB. P pab eca Q
eca=: 1 : 0("1)
:
n=. {. m
'x1 y1 z1'=. 3{.x,1
'x2 y2 z2'=. 3{.y,1
if. z1=0 do. 'x3 y3 z3'=.x2,y2,z2 else.
if. z2=0 do. 'x3 y3 z3'=.x1,y1,z1 else.
NB. genl case
u1=.x2*z12=.*:z1
u2=.x1*z22=.*:z2
s1=.y2*z1*z12
s2=.y1*z2*z22
w=.u1-u2
r=.s1-s2
if. w=0 do. 'x3 y3 z3' =. m ecd x1,y1,z1 else.
NB. continue genl
t=.u1+u2
m=.s1+s2
x3=.(*:r)-t*w2=.*:w
y3=.-:(r*(_2*x3)+t*w2)-m*w*w2
z3=.z1*z2*w
end. end. end.
n|x3,y3,z3
)
NB. Elliptic curve double point (mod proj coord)
NB. pab ecd Q
ecd=: 1 : 0
'n a '=. 2{.m
'x1 y1 z1'=.3{.y,1
if. 0 e. y1,z1 do. 'x2 y2 z2'=.0 1 0x else.
m1=.(3**:x1)+a**:*:z1
s=.4*x1*yy=.*:y1
x2=.(*:m1)+_2*s
y2=.(m1*s-x2)+_8**:yy
z2=.2*y1*z1
end.
n|x2,y2,z2
)
ecm=: 1 : 0
:
if. x=0 do. 0 1 0x return. end.
q=.y=.3{.y,1
for_mn. #.}.}:|:#:3 1*x do.
q=.m ecd q
select. mn
case. 2 do. q=.q m eca y
case. 1 do. q=.q m eca 1 _1 1* y
end.
end.
q
)
This only provides the correct addition results when the coordinates are
transformed (to z=1) between each call.
23&stdPoint@:( 23 1 1 eca ) each /\ 15 $ < 3 10
┌────┬────┬────┬────┬────┬────┬────┬─────┬───┬───┬─────┬───┬───┬───┬────┐
│3 10│7 12│19 5│17 3│9 16│12 4│11 3│13 16│0 1│6 4│18 20│5 4│1 7│4 0│1 16│
└────┴────┴────┴────┴────┴────┴────┴─────┴───┴───┴─────┴───┴───┴───┴────┘
but: ecd (double) does work:
23&stdPoint 23 1 1 ecd^:(i.5) 3 10 1
3 10
7 12
17 3
13 16
5 19
ecm also works which makes the addition weirdness irrelevant.
23 stdPoint ,. (2+ i.15) (23 1 1 ecm )"0 1] 3 10 1
7 12
19 0
17 3
9 16
12 0
18 0
13 16
0 1
6 4
16 0
5 0
1 7
1 15
5 9
5 19
ts '1943424321543223423423422342342234234234234343x (23 1 1 ecm2 )"0 1] 3 10
1x'
0.00761888 303232
23 1 1 ([ isCurvePoint stdPoint) 1943424321543223423423422342342343x (23 1 1
ecm )"0 1] 3 10 1x
1
where stdpoint and iscurvepoint
stdPoint =: 4 : 'x | (xx, yy)* x | 2 3 ^~ x invmod z [ x=. {.x[ ''xx yy z'' =.
y '"_ 1
isCurvePoint=: 4 : '0=p| (*:@{: - b + (a * {.) + ^&3@{.) y [ ''p a b''=. x'"1
----- Original Message -----
From: Pascal Jasmin <[email protected]>
To: "[email protected]" <[email protected]>
Cc:
Sent: Thursday, January 30, 2014 10:41:59 PM
Subject: Re: [Jprogramming] math requests
Hi Cliff,
I don't understand how to go from xyz back to xy coordinates.
At any rate, here is the affine (python) implementation. (I posted invmod
earlier):
Pointadd =: 1 : 0 NB. n is curve p a b
:
p =. {. m NB.'p a b' =. n
if. y -: p,0 do. x return. end.
if. x -: p,0 do. y return. end.
if. x -: y do. m Pointdouble y return. end.
'xx xy' =. x
'yx yy' =. y
if. (xx = yy) *. 0=p|xy+yy do. p,0 return. end.
l =. p | (p invmod yx - xx ) * yy - xy
(p| (l*xx -x3 ) - xy ) ,~ x3=. p | yx -~ xx -~ l * l
)
Pointdouble =: 4 : 0
'p a'=. 2{. x
if. y -: p,0 do. y return. end.
'xx xy' =. y
l =. p | (p invmod xy * 2 ) * a + 3 * *: xx
(p| (l*xx -x3 ) - xy ) ,~ x3=. p | (+: xx) -~ l * l
)
Pointmul =: 1 : 0 NB. sum of binary mask of repeated squares
:
m Pointadd/^:(1<#) |. bin # |. m Pointdouble^:(i. # bin =. 2 #. inv x) y
)
It passes the python tests, but it worries me that addition is not commutative.
I also don't know how to code the point at infinity (I put 0,p but that is
never reached).
3 10 (23 Pointadd) 9 7
17 20
(23 1 Pointadd) each /\ 18 # <3 10
┌────┬────┬────┬────┬────┬────┬────┬─────┬───┬─────┬───┬─────┬────┬─────┬────┬─────┬────┬────┐
│3 10│7 12│19 5│17 3│9 16│12 4│11 3│13 16│0 1│20 13│6 3│22 19│16 2│12 15│12
8│16 21│22 4│6 20│
└────┴────┴────┴────┴────┴────┴────┴─────┴───┴─────┴───┴─────┴────┴─────┴────┴─────┴────┴────┘
,. (2+ i.16) (23 1 Pointmul)"0 1 ] 3 10
7 12
19 5
17 3
9 16
12 4
11 3
13 16
0 1
6 4
18 20
16 20
5 15
13 21
2 21
5 19
18 3
These lists diverge after the item 0 1 is reached, which is the origin and a
good candidate for infinity? I don't seem to understand what order is.
----- Original Message -----
From: Cliff Reiter <[email protected]>
To: [email protected]
Cc:
Sent: Wednesday, January 29, 2014 3:32:21 PM
Subject: Re: [Jprogramming] math requests
Some elliptic curve stuff; I think there is a +1 error that Roger Hui
noticed in the factorization method.
http://archive.vector.org.uk/art10007270
http://archive.vector.org.uk/art10007280
Best, Cliff
On 1/29/2014 11:35 AM, Pascal Jasmin wrote:
>
> With all of the mathematicians on this list, these functions have likely been
> implemented before in J.
>
> elyptic curve point add, multiplication and double
> a python reference implementation:
> https://github.com/warner/python-ecdsa/blob/master/ecdsa/ellipticcurve.py
>
> the functions are: __add__ __mul__ and double
>
> if I may suggest J explicit signatures for the first 2 functions as:
>
> F =: 4 : 0
> 'yx yy yo' =. y
> 'xx xy xo' =. x
> )
>
> Some other methods than the python reference could be considered here:
>
> http://en.wikipedia.org/wiki/Elliptic_curve_point_multiplication
>
>
> also appreciated if you have in implementation of inverse_mod
> for reference function of same nate at:
> https://github.com/warner/python-ecdsa/blob/master/ecdsa/numbertheory.py
> ----------------------------------------------------------------------
> For information about J forums see http://www.jsoftware.com/forums.htm
>
--
Clifford A. Reiter
Lafayette College, Easton, PA 18042
http://webbox.lafayette.edu/~reiterc/
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
