The Elliptic Curve Factoring Method
Roman Pearce, Simon Fraser University, Canada
roman_pearce@hotmail.com
This worksheet demonstrates an integer factorization method based on elliptic curves modulo n. I will present only a brief introduction to the theory of elliptic curves before writing the procedures for adding and multiplying points, finsihing with the general factorization algorithm. For a more thorough development of the ideas behind this method, I refer the reader to the excellent book; A Course in Number Theory and Cryptography , by Neal Koblitz. Please use this worksheet for good only. Do not use for evil.
> ### WARNING: persistent store makes one-argument readlib obsolete restart: with(numtheory): with(plots): readlib(ilog):
Warning, the protected name order has been redefined and unprotected
Warning, the name changecoords has been redefined
Introduction :
An elliptic curve is a set of points ( ) satisfying . We also require that .
Here's an example :
> implicitplot(y^2=x^3-10*x-7,x=-3..6,y=-10..10,grid=[60,60]);
>
In addition to the points on the curve, we will also require a special point, known as the point at infinity. You can think of the point at infinity as having an infinite y-value. The important idea is that a line drawn through this point and any other point on the curve will be exactly vertical. We will denote the point at infinity by 'O'.
The set of points on an elliptic curve, together with the point at infinity, form a group under the following operation :
Given two points on the curve P and Q, we draw a line through P and Q. This line will pass through a third point on the curve, which will be denoted R. We now draw a second line through R and O, the point at infinity. This second line will also intersect the curve at a third point, which we will denote P+Q.
All of that sounds a little crazy (or maybe a lot), so let's look at an example :
> e := x-> evalf(sqrt(x^3-10*x-7)):
> curve:=implicitplot(y^2=x^3-10*x-7,x=-10..10,y=-10..10,grid=[72,72]):
> points:=pointplot([[-2.5,e(-2.5)],[-1.5,e(-1.5)],[4.36,e(4.36)],[4.36,-e(4.36)]],color=black):
> lines:=implicitplot({y-e(-2.5)=(e(-1.5)-e(-2.5))/(1)*(x+2.5),x=4.35},x=-3..6,y=-10..10,color=blue,grid=[72,72]):
> text:=textplot({[-2.7,2.5,`P`],[-1.4,3.3,`Q`],[4.0,6.3,`R`],[3.8,-6,`P+Q`],[4.1,10,`O`]},font=[TIMES,ITALIC,12]):
> display(curve,points,lines,text);
Here you can see as a line is drawn through and , it passes through a third point R. We then draw a vertical line through R, and the point where it intersects the curve again is defined to be . If by chance , we instead take the tangent to the curve at P, (the only line through P with multiplicity 2).
You may want to verify the following properties of elliptic addition :
- The point is the group identity, ie: for any point on the curve. We therefore define .
- If and have the same coordinant, but , then .
Let , . On the elliptic curve we have
where denotes the coordinant of .
and for
Elliptic Curves modulo n :
In order to factor numbers, we will be working with elliptic curves modulo a number n. This is basically the same beast with the same equations, only this time the equations are solved mod n and instead of dividing we must calculate inverses mod n.
The following program generates a random elliptic curve mod n, and a point [0,y] on the curve. If it finds a divisor of n in the process, it outputs that instead.
You may notice that we actually use the coordinant of the point to construct the constant term of the curve. The reason for this is that to calculate a random point any other way would involve a modular square root. For large n, this is simply too costly.
> generate_curve := proc(n, p::name, bound) local a,b,y,g,random;
> if nargs = 3 then random := rand(bound):
> else random := rand(50); fi;
> a := 0; y := 0; g := n;
> while ( a = 0 or y = 0 or g >= n ) do
> a := random() mod n; y := random() mod n; b := y^2 mod n;
> g := igcd(4*a^3 + 27*b^2,n);
> od;
> if (g = 1) then
> if (nargs >= 2) then p := [0,y]; fi;
> [a,b,n];
> else g; fi;
> end:
> generate_curve(41);
> generate_curve(43,'p');
> p;
Here is a procedure to add the two points and on an elliptic curve mod n. The program accepts as parameters , , and . It returns .
If for whatever reason the program can not compute an inverse mod n it instead returns the gcd, which will be a divisor of n.
> elliptic_add := proc(p1,p2,E) local g,inv,x1,x2,x3,y1,y2,y3,a,b,n;
> x1 := p1[1]; x2 := p2[1]; y1 := p1[2]; y2 := p2[2];
> a := E[1]; b := E[2]; n := E[3];
> if (x1=infinity and y1=infinity) then [x2,y2];
> elif (x2=infinity and y2=infinity) then [x1,y1];
> elif (x1=x2 and y1 <> y2) then [infinity,infinity];
> elif (x1=x2 and y1=y2) then
> g := igcdex(2*y1,n,inv);
> if g = 1 then
> x3 := ((3*x1^2+a)*inv)^2-2*x1 mod n;
> y3 := ((3*x1^2+a)*inv)*(x1-x3)-y1 mod n;
> [x3,y3];
> else
> g := igcdex(x2-x1,n,inv);
> x3 := ((y2-y1)*inv)^2-x1-x2 mod n;
> y3 := (y2-y1)*inv*(x1-x3)-y1 mod n;
> fi;
Example :
> E := generate_curve(41,'p');
> p1 := elliptic_add(p,p,E);
> p2 := elliptic_add(p1,p,E);
Another helpful operation will be multiplication of a point by a natural number constant.
I present the following algorithm for this, which is analagous to binary powering.
> elliptic_mul := proc(p1,k,E) local pn,ps,d,r;
> pn := p1; r := irem(k,2);
> if r = 1 then ps := p1; else ps := [infinity,infinity]; fi;
> d := (k-r)/2;
> while (d > 0 and (not (type(ps,integer) or type(pn,integer)) and (pn <> [infinity,infinity]))) do
> r := irem(d,2);
> pn := elliptic_add(pn,pn,E);
> if r = 1 then ps := elliptic_add(ps,pn,E); fi;
> d := (d-r)/2;
> if type(pn,integer) then pn;
> else ps; fi;
> E := generate_curve(31,'p'); p;
> p2 := elliptic_add(p,p,E);
> elliptic_add(p,p2,E);
> elliptic_mul(p,3,E);
The Elliptic Curve Factoring Algorithm :
We are now ready to present the elliptic curve factorization method. The basic idea is to construct a random elliptic curve modulo and a point on . We then compute for some large integer . If at any point the elliptic addition formula fails, it will output a divisor of . Otherwise, we try again with a new elliptic curve.
> ecf := proc(n,boundB,boundC) local B,C,E,p,k,d,i,q;
> if isprime(n) then RETURN(n); fi;
> if nargs = 2 then B := boundB; else B := 2*length(n); fi;
> if nargs = 3 then C := boundC; else C := ilog[2](n); fi;
> q := 2; k := 1;
> while q < B do
> k := k*q^ilog[q](C);
> q := nextprime(q); od;
> d := 1;
> while (d=1 or d=n) do
> E := generate_curve(n,'p');
> if type(E,integer) then d := E;
> d := elliptic_mul(p,k,E);
> if not type(d,integer) then d := 1; fi;
> d;
And there you have it. Now let's try it out...
> ecf(1763);
> ecf(54354653);
> ecf(65873986798549);
> ifactor(65873986798549);
> ecf(7598437895743975894379);
> ecf(5766865804051782016322656112299);
> ecf(3607948903983773);
> ifactor(3607948903983773);
That's all for today folks, I hope you enjoyed the show. Email me if you have any suggestions or find any mistakes. Future topics are likely to include : elliptic curve cryptosystems, other integer factoring methods, and anything else that's interesting and fun to program. Feel free to send me suggestions.
undergraduate at Simon Fraser Universisty, Canada
Disclaimer: While every effort has been made to validate the solutions in this worksheet, Waterloo Maple Inc. and the contributors are not responsible for any errors contained and are not liable for any damages resulting from the use of this material (especially evil use).
This worksheet was originally written in Maple V R5 and completed August 11, 2000