|
|
last edited 11 years ago by test1 |
1 2 | ||
Editor: test1
Time: 2013/07/31 17:38:02 GMT+0 |
||
Note: |
changed: - A := u**2*X.3*Y.3-v**3-2*v**2*X.1*Y.3 A := u^2*X.3*Y.3-v^3-2*v^2*X.1*Y.3 changed: - y3 := (u*(v**2*X.1*Y.3-A)-v**3*X.2*Y.3) rem N - z3 := (v**3*X.3*Y.3) rem N y3 := (u*(v^2*X.1*Y.3-A)-v^3*X.2*Y.3) rem N z3 := (v^3*X.3*Y.3) rem N changed: - w := a*X.3**2+3*X.1**2 w := a*X.3^2+3*X.1^2 changed: - h := w**2 - 8*B h := w^2 - 8*B changed: - y3 := (w*(4*B-h)-8*X.2**2*s**2) rem N - z3 := (8*s**3) rem N y3 := (w*(4*B-h)-8*X.2^2*s^2) rem N z3 := (8*s^3) rem N
On Jul 10, 2007 9:08 AM Alasdair McAndrew? wrote:
I have written some highly unoptimized code for factoring integers using Lenstra's elliptic curve method, with the "birthday paradox" phase two developed by Richard Brent. Even at this stage, it can factor the seventh Fermat number 2^2^7+1 in 352 seconds, as opposed to 1877 seconds by the in-built factoring method. If anybody is interesting in developing this code further, do let me know.
All I really want to do is to extend the methods in intfact.spad to include ecm - which the file itself recommends.
-Alasdair
(1) -> mod(a:Integer,b:PositiveInteger):Integer == (a>0 => a rem b::Integer; b + (a rem b))
Function declaration mod : (Integer,PositiveInteger) -> Integer has been added to workspace.
a
is a parameter which describes the elliptic curve, and which is
the one used for addition and point doubling.
ecAdd(a:Integer,X:List Integer, Y:List Integer, N:Integer):List Integer == local u, v, A, x3, y3, z3 if (X=Y)::Boolean then return ecDouble(a, X, N) else u := Y.2*X.3-X.2*Y.3 v := Y.1*X.3-X.1*Y.3 A := u^2*X.3*Y.3-v^3-2*v^2*X.1*Y.3 x3 := (v*A) rem N y3 := (u*(v^2*X.1*Y.3-A)-v^3*X.2*Y.3) rem N z3 := (v^3*X.3*Y.3) rem N return [x3, y3, z3]
Function declaration ecAdd : (Integer,List(Integer), List(Integer) , Integer) -> List(Integer) has been added to workspace.
ecDouble(a:Integer,X:List Integer, N:Integer):List Integer == local w, s, B, h, x3, y3, z3 w := a*X.3^2+3*X.1^2 s := X.2*X.3 B := X.1*X.2*s h := w^2 - 8*B x3 := (2*h*s) rem N y3 := (w*(4*B-h)-8*X.2^2*s^2) rem N z3 := (8*s^3) rem N return [x3, y3, z3]
Function declaration ecDouble : (Integer,List(Integer), Integer) -> List(Integer) has been added to workspace.
ecMultiply(a:Integer,X:List Integer, n:Integer, N:Integer):List Integer == local lst, P lst := wholeRagits(n::RadixExpansion 2) P:=X for i in rest lst repeat P:=ecDouble(a, P, N) if i = 1 then P:=ecAdd(a, P, X, N) return P
Function declaration ecMultiply : (Integer,List(Integer), Integer, Integer) -> List(Integer) has been added to workspace.
Here is the main routine; you have to enter both the bounds B and r; for example:
ecFactor(2^101-1,10000,700)
ecFactor(N:Integer,B:Integer, r:Integer):List Integer == local a:Integer, x:Integer, y:Integer, b, g, curve, found, Q, f, i, tmp found:=false i:=0 while not(found) repeat curve:=false while not(curve) repeat a:=random(N) x:=random(N) y:=random(N) b:=y^2-x^3-a*x g:=gcd(4*a^3+27*b^2, N) if g = 1 then curve:=true i:=i+1 if (i rem 20)=0 then output i Q:=ecPhaseOne(a, [x, y, 1], B, N) f := gcd(Q.3, N) if f>1 and f<N then found:=true else tmp:=ecPhaseTwo(a, Q, N, r) if tmp>1 and tmp<N then found:=true return [f, N/f]
Function declaration ecFactor : (Integer,Integer, Integer) -> List( Integer) has been added to workspace.
ecPhaseOne(a:Integer,X:List Integer, B:PositiveInteger, N:PositiveInteger):List Integer == local P, i, e P := X for i in primes(2, B) repeat e := wholePart(log2(B::DoubleFloat)/log2(i::DoubleFloat))::PositiveInteger P := ecMultiply(a, P, i^e, N) return P
Function declaration ecPhaseOne : (Integer,List(Integer), PositiveInteger, PositiveInteger) -> List(Integer) has been added to workspace.
ecPhaseTwo(a:Integer,Q:List Integer, N:PositiveInteger, r:PositiveInteger):Integer == local randr, Qs, tmp, i, j, k, d, result result := 1 randr := rest wholeRagits(2^r+random(2^r)::RadixExpansion 2) Qs := [Q] for i in 2..r repeat tmp := ecDouble(a, Qs.(i-1), N) if randr.i = 1 then tmp := ecAdd(a, tmp, Q, N) Qs := append(Qs, [tmp]) for i in 2..r repeat if member?(Qs.i, [Qs.k for k in 1..i-1]) then --j := position(Qs.i, [Qs.k for k in 1..i-1]) d := reduce(*, [reduce(*, [Qs.i.3-Qs.j.3 for j in i+1..r]) for i in _ 1..r-1]) if gcd(N, d)>1 then result:=gcd(N, d) return result
Function declaration ecPhaseTwo : (Integer,List(Integer), PositiveInteger, PositiveInteger) -> Integer has been added to workspace.
See http://www.alpertron.com.ar/ECM.HTM for values of B.
For example:
)set message time on
factor 2^101-1
(1) |
Time: 1.32 (EV) + 0.03 (OT) = 1.35 sec
ecFactor(2^101-1,10000, 700)
Compiling function ecDouble with type (Integer,List(Integer), Integer) -> List(Integer)
Compiling function ecAdd with type (Integer,List(Integer), List( Integer), Integer) -> List(Integer)
Compiling function ecMultiply with type (Integer,List(Integer), Integer, Integer) -> List(Integer)
Compiling function ecPhaseOne with type (Integer,List(Integer), PositiveInteger, PositiveInteger) -> List(Integer)
Compiling function ecPhaseTwo with type (Integer,List(Integer), PositiveInteger, PositiveInteger) -> Integer
Compiling function ecFactor with type (Integer,Integer, Integer) -> List(Integer)
Compiling function G27 with type Integer -> Boolean
Compiling function G29 with type NonNegativeInteger -> Boolean 20
(2) |
Time: 0.02 (IN) + 3.79 (EV) + 0.05 (OT) + 0.01 (GC) = 3.88 sec
test( reduce(*,%)=2^101-1 )
(3) |
Time: 0 sec
)abbrev package EC EllipticCurve --#pile --#include "axiom" --import from Integer EllipticCurve(): with ecAdd:(Integer,List Integer, List Integer, Integer)->List Integer ecDouble:(Integer, List Integer, Integer) -> List Integer ecMultiply:(Integer, List Integer, Integer, Integer) -> List Integer == add ecAdd(a:Integer, X:List Integer, Y:List Integer, N:Integer):List Integer == if X=Y then ecDouble(a, X, N) else u := Y.2*X.3-X.2*Y.3 v := Y.1*X.3-X.1*Y.3 A := u^2*X.3*Y.3-v^3-2*v^2*X.1*Y.3 x3 := (v*A) rem N y3 := (u*(v^2*X.1*Y.3-A)-v^3*X.2*Y.3) rem N z3 := (v^3*X.3*Y.3) rem N [x3, y3, z3]
ecDouble(a:Integer,X:List Integer, N:Integer):List Integer == w := a*X.3^2+3*X.1^2 s := X.2*X.3 B := X.1*X.2*s h := w^2 - 8*B x3 := (2*h*s) rem N y3 := (w*(4*B-h)-8*X.2^2*s^2) rem N z3 := (8*s^3) rem N [x3, y3, z3]
ecMultiply(a:Integer,X:List Integer, n:Integer, N:Integer):List Integer == -- import from RadixExpansion 2, List Integer lst := wholeRagits(n::RadixExpansion 2) P:=X for i in rest lst repeat P:=ecDouble(a, P, N) if i = 1 then P:=ecAdd(a, P, X, N) P
Compiling FriCAS source code from file /var/lib/zope2.10/instance/axiom-wiki/var/LatexWiki/5435117267253907166-25px007.spad using old system compiler. EC abbreviates package EllipticCurve ------------------------------------------------------------------------ initializing NRLIB EC for EllipticCurve compiling into NRLIB EC compiling exported ecAdd : (Integer,List Integer, List Integer, Integer) -> List Integer Time: 0.34 SEC.
compiling exported ecDouble : (Integer,List Integer, Integer) -> List Integer Time: 0.07 SEC.
compiling exported ecMultiply : (Integer,List Integer, Integer, Integer) -> List Integer Time: 0 SEC.
(time taken in buildFunctor: 0)
;;; *** |EllipticCurve| REDEFINED
;;; *** |EllipticCurve| REDEFINED Time: 0 SEC.
Cumulative Statistics for Constructor EllipticCurve Time: 0.41 seconds
finalizing NRLIB EC Processing EllipticCurve for Browser database: --->-->EllipticCurve(constructor): Not documented!!!! --->-->EllipticCurve((ecAdd ((List (Integer)) (Integer) (List (Integer)) (List (Integer)) (Integer)))): Not documented!!!! --->-->EllipticCurve((ecDouble ((List (Integer)) (Integer) (List (Integer)) (Integer)))): Not documented!!!! --->-->EllipticCurve((ecMultiply ((List (Integer)) (Integer) (List (Integer)) (Integer) (Integer)))): Not documented!!!! --->-->EllipticCurve(): Missing Description ; compiling file "/var/aw/var/LatexWiki/EC.NRLIB/EC.lsp" (written 22 NOV 2024 01:07:36 AM):
; wrote /var/aw/var/LatexWiki/EC.NRLIB/EC.fasl ; compilation finished in 0:00:00.028 ------------------------------------------------------------------------ EllipticCurve is now explicitly exposed in frame initial EllipticCurve will be automatically loaded when needed from /var/aw/var/LatexWiki/EC.NRLIB/EC
factor(N:Integer,B:Integer, r:Integer):List Integer == local a:Integer, x:Integer, y:Integer, b, g, curve, found, Q, f, i, tmp found:=false i:=0 while not(found) repeat curve:=false while not(curve) repeat a:=random(N) x:=random(N) y:=random(N) b:=y^2-x^3-a*x g:=gcd(4*a^3+27*b^2, N) if g = 1 then curve:=true i:=i+1 if (i rem 20)=0 then output i Q:=ecPhase1(a, [x, y, 1], B, N) f := gcd(Q.3, N) if f>1 and f<N then found:=true else tmp:=ecPhase2(a, Q, N, r) if tmp>1 and tmp<N then found:=true return [f, N/f]
Function declaration factor : (Integer,Integer, Integer) -> List( Integer) has been added to workspace.
Time: 0 sec
ecPhase1(a:Integer,X:List Integer, B:PositiveInteger, N:PositiveInteger):List Integer == local P, i, e P := X for i in primes(2, B) repeat e := wholePart(log2(B::DoubleFloat)/log2(i::DoubleFloat))::PositiveInteger P := ecMultiply(a, P, i^e, N)$EC return P
Function declaration ecPhase1 : (Integer,List(Integer), PositiveInteger, PositiveInteger) -> List(Integer) has been added to workspace.
Time: 0 sec
ecPhase2(a:Integer,Q:List Integer, N:PositiveInteger, r:PositiveInteger):Integer == local randr, Qs, tmp, i, j, k, d, result result := 1 randr := rest wholeRagits(2^r+random(2^r)::RadixExpansion 2) Qs := [Q] for i in 2..r repeat tmp := ecDouble(a, Qs.(i-1), N)$EC if randr.i = 1 then tmp := ecAdd(a, tmp, Q, N)$EC Qs := append(Qs, [tmp]) for i in 2..r repeat if member?(Qs.i, [Qs.k for k in 1..i-1]) then --j := position(Qs.i, [Qs.k for k in 1..i-1]) d := reduce(*, [reduce(*, [Qs.i.3-Qs.j.3 for j in i+1..r]) for i in _ 1..r-1]) if gcd(N, d)>1 then result:=gcd(N, d) return result
Function declaration ecPhase2 : (Integer,List(Integer), PositiveInteger, PositiveInteger) -> Integer has been added to workspace.
Time: 0 sec
See http://www.alpertron.com.ar/ECM.HTM for values of B.
For example:
)set message time on
factor 2^101-1
(4) |
Time: 1.33 (EV) = 1.34 sec
factor(2^101-1,10000, 700)
Compiling function ecPhase1 with type (Integer,List(Integer), PositiveInteger, PositiveInteger) -> List(Integer)
Compiling function ecPhase2 with type (Integer,List(Integer), PositiveInteger, PositiveInteger) -> Integer
Compiling function factor with type (Integer,Integer, Integer) -> List(Integer) 20
(5) |
Time: 0.03 (IN) + 4.17 (EV) + 0.03 (OT) + 0.02 (GC) = 4.25 sec
test( reduce(*,%)=2^101-1 )
(6) |
Time: 0 sec