Dear Peter, *
I was waiting for your small Aldor patch and never see it. Did you send
it somewhere ? The main concern of this mail is if you are able to
compile ecfact.as (attached) and _execute_ the
LenstraEllipticMethod(Integer) function. On my system it fails, a bug
somewhere. Am I alone with this issue ? Martin ? I'm still using the
old build process and want to switch to the new one but this stops me. I
do not want to take your time but maybe you have an idea of what's going
on with this issue.
The output is :
Looking in OutputPackage() for ??349042727 with code 483270060
and an error message : "export not found"
Greg
PS : I precise that I applied one of your patchs, sorry I don't remember
which one, to be able to use recent versions of Aldor. The problem
mentioned here still remains on my machine.
Le jeudi 17 janvier 2008 à 09:32 +0000, Peter Broadbery a écrit :
> Finally got round to looking at this.
> The code attached can be used to build the axiom/aldor interopability
> library (libaxiom.al).
>
> New features:
> Less abuse of make features (though I've only used a recent version of
> make, and have used some newer features)
> Dependency on java removed.
> Cleaner build process - all the dependency analysis is done in the
> aldor code. It should be more data driven, but it isn't.
>
> There's a small aldor patch that I'll send along this evening.
>
> Peter
> _______________________________________________
> Axiom-developer mailing list
> Axiom-developer@nongnu.org
> http://lists.nongnu.org/mailman/listinfo/axiom-developer
-- Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
-- All rights reserved.
--
-- Redistribution and use in source and binary forms, with or without
-- modification, are permitted provided that the following conditions are
-- met:
--
-- - Redistributions of source code must retain the above copyright
-- notice, this list of conditions and the following disclaimer.
--
-- - Redistributions in binary form must reproduce the above copyright
-- notice, this list of conditions and the following disclaimer in
-- the documentation and/or other materials provided with the
-- distribution.
--
-- - Neither the name of The Numerical ALgorithms Group Ltd. nor the
-- names of its contributors may be used to endorse or promote products
-- derived from this software without specific prior written permission.
--
-- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
-- IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
-- TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
-- PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
-- OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
-- EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
-- PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
-- PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
-- LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
-- NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
-- SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
#include "axiom.as"
#pile
--% Elliptic curve method for integer factorization
-- This file implements Lenstra's algorithm for integer factorization.
-- A divisor of N is found by computing a large multiple of a rational
-- point on a randomly generated elliptic curve in P2 Z/NZ.
-- The Hessian model is used for the curve (1) to simplify the selection
-- of the initial point on the random curve and (2) to minimize the
-- cost of adding points.
-- Ref: IBM RC 11262, DV Chudnovsky & GV Chudnovsky
-- SMW Sept 86.
--% EllipticCurveRationalPoints
--)abbrev domain ECPTS EllipticCurveRationalPoints
EllipticCurveRationalPoints(x0:Integer, y0:Integer, z0:Integer, n:Integer): ECcat == ECdef where
Point ==> Record(x: Integer, y: Integer, z: Integer)
ECcat ==> AbelianGroup with
double: % -> %
p0: %
HessianCoordinates: % -> Point
ECdef ==> add
Rep == Point
import from Rep
import from List Integer
Ex == OutputForm
default u, v: %
apply(u:%,x:'x'):Integer == rep(u).x
apply(u:%,y:'y'):Integer == rep(u).y
apply(u:%,z:'z'):Integer == rep(u).z
import from 'x'
import from 'y'
import from 'z'
coerce(u:%): Ex == [u.x, u.y, u.z]$List(Integer) :: Ex
p0:% == per [x0 rem n, y0 rem n, z0 rem n]
HessianCoordinates(u:%):Point == rep u
0:% ==
per [1, (-1) rem n, 0]
-(u:%):% ==
per [u.y, u.x, u.z]
(u:%) = (v:%):Boolean ==
XuZv := u.x * v.z
XvZu := v.x * u.z
YuZv := u.y * v.z
YvZu := v.y * u.z
(XuZv-XvZu) rem n = 0 and (YuZv-YvZu) rem n = 0
(u:%) + (v:%): % ==
XuZv := u.x * v.z
XvZu := v.x * u.z
YuZv := u.y * v.z
YvZu := v.y * u.z
(XuZv-XvZu) rem n = 0 and (YuZv-YvZu) rem n = 0 => double u
XuYv := u.x * v.y
XvYu := v.x * u.y
Xw := XuZv*XuYv - XvZu*XvYu
Yw := YuZv*XvYu - YvZu*XuYv
Zw := XvZu*YvZu - XuZv*YuZv
per [Yw rem n, Xw rem n, Zw rem n]
double(u:%): % ==
import from PositiveInteger
X3 := u.x**([EMAIL PROTECTED])
Y3 := u.y**([EMAIL PROTECTED])
Z3 := u.z**([EMAIL PROTECTED])
Xw := u.x*(Y3 - Z3)
Yw := u.y*(Z3 - X3)
Zw := u.z*(X3 - Y3)
per [Yw rem n, Xw rem n, Zw rem n]
(n:Integer)*(u:%): % ==
n < 0 => (-n)*(-u)
v := 0
import from UniversalSegment Integer
for i in 0..length n - 1 repeat
if bit?(n,i) then v := u + v
u := double u
v
--% EllipticCurveFactorization
--)abbrev package ECFACT EllipticCurveFactorization
EllipticCurveFactorization: with
LenstraEllipticMethod: (Integer) -> Integer
LenstraEllipticMethod: (Integer, Float) -> Integer
LenstraEllipticMethod: (Integer, Integer, Integer) -> Integer
LenstraEllipticMethod: (Integer, Integer) -> Integer
lcmLimit: Integer -> Integer
lcmLimit: Float-> Integer
solveBound: Float -> Float
bfloor: Float -> Integer
primesTo: Integer -> List Integer
lcmTo: Integer -> Integer
== add
import from List Integer
Ex == OutputForm
import from Ex
import from String
import from Float
NNI==> NonNegativeInteger
import from OutputPackage
import from Integer, NonNegativeInteger
import from UniversalSegment Integer
blather:Boolean := true
--% Finding the multiplier
flabs (f: Float): Float == abs f
flsqrt(f: Float): Float == sqrt f
nthroot(f:Float,n:Integer):Float == exp(log f/n::Float)
bfloor(f: Float): Integer == wholePart floor f
lcmLimit(n: Integer):Integer ==
lcmLimit nthroot(n::Float, 3)
lcmLimit(divisorBound: Float):Integer ==
y := solveBound divisorBound
lcmLim := bfloor exp(log divisorBound/y)
if blather then
output("The divisor bound is", divisorBound::Ex)
output("The lcm Limit is", lcmLim::Ex)
lcmLim
-- Solve the bound equation using a Newton iteration.
--
-- f = y**2 - log(B)/log(y+1)
--
-- f/f' = fdf =
-- 2 2
-- y (y + 1)log(y + 1) - (y + 1)log(y + 1) logB
-- ---------------------------------------------
-- 2
-- 2y(y + 1)log(y + 1) + logB
--
fdf(y: Float, logB: Float): Float ==
logy := log(y + 1)
ylogy := (y + 1)*logy
ylogy2:= y*logy*ylogy
(y*ylogy2 - logB*ylogy)/(([EMAIL PROTECTED])*ylogy2 + logB)
solveBound(divisorBound:Float):Float ==
-- solve y**2 = log(B)/log(y + 1)
-- although it may be y**2 = log(B)/(log(y)+1)
relerr := (10::Float)**(-5)
logB := log divisorBound
y0 := flsqrt log10 divisorBound
y1 := y0 - fdf(y0, logB)
while flabs((y1 - y0)/y0) > relerr repeat
y0 := y1
y1 := y0 - fdf(y0, logB)
y1
-- maxpin(p, n, logn) is max d s.t. p**d <= n
maxpin(p:Integer,n:Integer,logn:Float): NonNegativeInteger ==
d: Integer := bfloor(logn/log(p::Float))
if d < 0 then d := 0
d::NonNegativeInteger
multiple?(i: Integer, plist: List Integer): Boolean ==
for p in plist repeat if i rem p = 0 then return true
false
primesTo(n:Integer):List Integer ==
n < 2 => []
n = 2 => [2]
plist := [3, 2]
i:Integer := 5
while i <= n repeat
if not multiple?(i, plist) then plist := cons(i, plist)
i := i + 2
if not multiple?(i, plist) then plist := cons(i, plist)
i := i + 4
plist
lcmTo(n:Integer):Integer ==
plist := primesTo n
m: Integer := 1
logn := log(n::Float)
for p in plist repeat m := m * p**maxpin(p,n,logn)
if blather then
output("The lcm of 1..", n::Ex)
output(" is", m::Ex)
m
LenstraEllipticMethod(n: Integer):Integer ==
LenstraEllipticMethod(n, flsqrt(n::Float))
LenstraEllipticMethod(n: Integer, divisorBound: Float):Integer ==
lcmLim0 := lcmLimit divisorBound
multer0 := lcmTo lcmLim0
LenstraEllipticMethod(n, lcmLim0, multer0)
InnerLenstraEllipticMethod(n:Integer, multer:Integer,
X0:Integer, Y0:Integer, Z0:Integer):Integer ==
import from EllipticCurveRationalPoints(X0,Y0,Z0,n)
import from Record(x: Integer, y: Integer, z: Integer)
p := p0
pn := multer * p
Zn := HessianCoordinates.pn.z
gcd(n, Zn)
LenstraEllipticMethod(n: Integer, multer: Integer):Integer ==
X0:Integer := random()
Y0:Integer := random()
Z0:Integer := random()
InnerLenstraEllipticMethod(n, multer, X0, Y0, Z0)
LenstraEllipticMethod(n:Integer, lcmLim0:Integer, multer0:Integer):Integer ==
nfact: Integer := 1
for i:Integer in 1.. while nfact = 1 repeat
output("Trying elliptic curve number", i::Ex)
X0:Integer := random()
Y0:Integer := random()
Z0:Integer := random()
nfact := InnerLenstraEllipticMethod(n, multer0, X0, Y0, Z0)
if nfact = n then
lcmLim := lcmLim0
while nfact = n repeat
output("Too many iterations... backing off")
lcmLim := bfloor(lcmLim * 0.6)
multer := lcmTo lcmLim
nfact := InnerLenstraEllipticMethod(n, multer0, X0, Y0, Z0)
nfact
_______________________________________________
Axiom-developer mailing list
Axiom-developer@nongnu.org
http://lists.nongnu.org/mailman/listinfo/axiom-developer
