--* From BMT%WATSON.vnet.ibm.com@yktvmh.watson.ibm.com Thu Aug 19 16:02:53 1993 --* Received: from yktvmh.watson.ibm.com by radical.watson.ibm.com (AIX 3.2/UCB 5.64/900524) --* id AA12665; Thu, 19 Aug 1993 16:02:53 -0400 --* Received: from watson.vnet.ibm.com by yktvmh.watson.ibm.com (IBM VM SMTP V2R3) --* with BSMTP id 0114; Thu, 19 Aug 93 16:04:11 EDT --* Received: from YKTVMH by watson.vnet.ibm.com with "VAGENT.V1.0" --* id --* for asbugs@watson; Thu, 19 Aug 93 16:04:10 -0400 --* Received: from YKTVMH by watson.vnet.ibm.com with "VAGENT.V1.0" --* id 2245; Thu, 19 Aug 1993 16:04:08 EDT --* Received: from cyst.watson.ibm.com by yktvmh.watson.ibm.com (IBM VM SMTP V2R3) --* with TCP; Thu, 19 Aug 93 16:04:07 EDT --* Received: from spadserv.watson.ibm.com by cyst.watson.ibm.com (AIX 3.2/UCB 5.64/900528) --* id AA55316; Thu, 19 Aug 1993 16:03:03 -0400 --* Received: by spadserv.watson.ibm.com (AIX 3.2/UCB 5.64/900524) --* id AA24277; Thu, 19 Aug 1993 16:08:38 -0400 --* Date: Thu, 19 Aug 1993 16:08:38 -0400 --* From: bmt@spadserv.watson.ibm.com --* X-External-Networks: yes --* Message-Id: <9308192008.AA24277@spadserv.watson.ibm.com> --* To: asbugs@watson.ibm.com --* Subject: no output is produced when running: asharp -DTestRun -lasdem -r intfact.as [intfact.as][29.7 (current)] --@ Fixed by: SMW Fri Aug 20 00:33:00 1993 --@ Tested by: triv4.as --@ Summary: Installed genfoam.c changes for functions returning 0 results. --> testcomp -l asdem --> testrun -l asdem --Copyright The Numerical Algorithms Group Limited 1991. #include "aslib.as" #library DemoLib "asdem" import DemoLib macro Boolean == Bit B == Bit N == Integer I == Integer min(aa,bb) == if aa < bb then aa else bb one? aa == aa = 1 facFlags ==> 'nil, sqfr, irred, prime' FFE(S:BasicType): BasicType with bracket: (facFlags, S, Integer) -> % apply: (%, 'xpnt') -> Integer apply: (%, 'fctr') -> S apply: (%, 'flg') -> facFlags set!: (%, 'xpnt', Integer) -> Integer == add Rep ==> Record(flg:facFlags, fctr:S, xpnt:Integer) import Rep import S import Integer (f:%) = (g:%):Bit == rep(f).fctr = rep(g).fctr and rep(f).xpnt = rep(g).xpnt (f:%) ~= (g:%):Bit == not(f=g) apply(f:%, tag:'xpnt'):Integer == rep(f).xpnt apply(f:%, tag:'fctr'):S == rep(f).fctr apply(f:%, tag:'flg'):facFlags == rep(f).flg set!(f:%, tag:'xpnt', e:Integer):Integer == set!(rep(f),xpnt,e) set!(f:%, tag:'fctr', s:S):S == set!(rep(f),fctr,s) set!(f:%, tag:'flg', flag:facFlags):facFlags == set!(rep(f),flg,flag) [flag:facFlags, s:S, e:Integer]:% == per [flag,s,e] apply(p: Outport, l: %): Outport == import String import S import Integer p("ffe(fctr: ")(rep(l).fctr)(" xpnt: ")(rep(l).xpnt)(")") Factored(S: BasicType with 1:% ) : with makeFR: (S, List FFE(S)) -> % factorList: % -> List FFE(S) unit: % -> S apply: (Outport, %) -> Outport == add Rep ==> Record(unit:S, faclist:List FFE(S)) import Rep makeFR(u:S, lffe:List FFE(S)):% == per [u, lffe] factorList(x:%):List FFE(S) == apply(rep(x),faclist) unit(x:%):S == apply(rep(x),unit) apply(p:Outport, f:%): Outport == import String import S import List FFE(S) import FFE(S) import 'xpnt' import 'fctr' import I if unit(f) ~= 1 then p:=p(unit f) for fac in factorList f repeat p:=p(" ")((fac.fctr)@S) if fac.xpnt ~= 1@I then p:=p("^")((fac.xpnt)@I) p --% IntegerFactorizationPackage -- recoded MBM Nov/87 +++ This Package contains basic methods for integer factorization. +++ The factor operation employs trial division up to 10,000. It +++ then tests to see if n is a perfect power before using Pollards +++ rho method. Because Pollards method may fail, the result +++ of factor may contain composite factors. We should also employ +++ Lenstra's elliptic curve method. B ==> Bit FF ==> Factored I LMI ==> ListMultiDictionary I IntegerFactorizationPackage: with factor : I -> FF ++ factor(n) returns the full factorization of integer n squareFree : I -> FF ++ squareFree(n) returns the square free factorization of integer n BasicMethod : I -> FF ++ BasicMethod(n) returns the factorization ++ of integer n by trial division PollardSmallFactor: I -> Partial(I) ++ PollardSmallFactor(n) returns a factor ++ of n or "failed" if no one is found == add import Segment I import IntegerRoots import IntegerPrimesPackage import Partial I import 'nil,sqfr,irred,prime' import FFE(I) import List(FFE(I)) import FF import 'xpnt' concat!(l:List FFE I, v:FFE I):List FFE I == concat!(l, cons(v,nil)) BasicSieve: (I, I) -> FF squareFree(n:I):FF == if n<0 then (m := -n; u:I := -1) else (m := n; u := 1) (m > 1) and (not(failed? (v := perfectSqrt m))) => l := factorList(sv := squareFree(retract v)) for rec in generator l repeat rec.xpnt := 2 * rec.xpnt makeFR(u * unit sv, l) lim := 1 + approxNthRoot(m,3) lim > 100000 => makeFR(u, factorList factor m) x := BasicSieve(m, lim) y := one? unit x => factorList x concat_!(factorList x, [sqfr,unit(x),1]$FFE(I)) makeFR(u, y) -- Pfun(y: I,n: I): I == (y**2 + 5) rem n PollardSmallFactor(n:I):Partial(I) == -- Use the Brent variation m := 100 y:I := random(n)$RandomNumberSource r:I := 1 q:I := 1 G:I := 1 while G <= 1 repeat x:I := y for i in 1..r repeat y := (y*y+5) rem n q := (q*abs(x-y)) rem n k:I := 0 while not( (k>=r) or (G>1)) repeat ys:I := y for i in 1..min(m,r-k) repeat y := (y*y+5) rem n q := q*abs(x-y) rem n G := gcd(q,n) k := k+m r := 2*r if G=n then while G<=1 repeat ys := (ys*ys+5) rem n G := gcd(abs(x-ys),n) G=n => failed G::Partial(I) BasicSieve(r:I, lim:I):FF == import List(I) l:List(I) := [1, 2, 2, 4, 2, 4, 2, 4, 6, 2, 6] concat_!(l, rest rest rest l) d := 2 n := r ls := empty()$List(FFE(I)) for s in l repeat d > lim => return makeFR(n, ls) if n1 then ls := concat_!(ls, [prime,n,1]$FFE(I)) return makeFR(1, ls) m:I:=0 for mm in 1.. while zero?(n rem d) repeat (n := n quo d; m:=mm) if m>0 then ls := concat_!(ls, [prime,d, m]$FFE(I)) d := d+s never BasicMethod(n:I):FF == local u:I if n<0 then (m := -n; u := -1) else (m := n; u := 1) x := BasicSieve(m, 1 + approxSqrt(m)$IntegerRoots) makeFR(u, factorList x) factor(m:I):FF == local u:I zero? m => makeFR(m,nil) if negative? m then (n := -m; u := -1) else (n := m; u := 1) bfac := BasicSieve(n, 10000) flb := factorList bfac import List(FFE I) import String one?(n := unit bfac) => makeFR(u, flb) a:LMI := dictionary() -- numbers yet to be factored b:LMI := dictionary() -- prime factors found f:LMI := dictionary() -- number which could not be factored insert_!(n, a) while not empty? a repeat n := inspect a; c:SingleInteger := count(n, a); remove_!(n, a) prime?(n)$IntegerPrimesPackage => insert_!(n, b, c) -- test for a perfect power (s :Record(base:I, exponent:I) := perfectNthRoot n).exponent > 1 => insert_!(s.base, a, c * retract s.exponent) not(failed? (dd := PollardSmallFactor n)) => d := retract dd for m:SingleInteger in 0$SingleInteger.. while zero?(n rem d) repeat n := n quo d insert_!(d, a, retract(m) * c) if n > 1 then insert_!(n, a, c) -- an elliptic curve factorization attempt should be made here insert_!(n, f, c) -- insert prime factors found while not empty? b repeat n := inspect b; c := count(n, b); remove_!(n, b) flb := concat_!(flb, [prime,n,c::I]$FFE(I)) -- insert non-prime factors found while not empty? f repeat n := inspect f; c := count(n, f); remove_!(n, f) flb := concat_!(flb, [nil,n,c::I]$FFE(I)) makeFR(u, flb) #if TestRun import I a:I:=3 print(a)() import IntegerPrimesPackage print(nextPrime a)() import List(I) print(primes(3,50))() import IntegerFactorizationPackage import Factored(I) a:=50 print(factor a)() print(factor 6)() print(factor 2349587234)() print(factor 23434534)() #endif