--* From postmaster%watson.vnet.ibm.com@yktvmv.watson.ibm.com Fri Aug 20 16:05:16 1993 --* Received: from yktvmv2.watson.ibm.com by radical.watson.ibm.com (AIX 3.2/UCB 5.64/900524) --* id AA13544; Fri, 20 Aug 1993 16:05:16 -0400 --* X-External-Networks: yes --* Received: from watson.vnet.ibm.com by yktvmv.watson.ibm.com (IBM VM SMTP V2R3) --* with BSMTP id 4401; Fri, 20 Aug 93 16:06:36 EDT --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id --* for asbugs@watson; Fri, 20 Aug 93 16:06:36 -0400 --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id 5353; Fri, 20 Aug 1993 16:06:34 EDT --* Received: from matthew.watson.ibm.com by yktvmv.watson.ibm.com --* (IBM VM SMTP V2R3) with TCP; Fri, 20 Aug 93 16:06:32 EDT --* Received: by matthew.watson.ibm.com (AIX 3.2/UCB 5.64/4.03) --* id AA21471; Fri, 20 Aug 1993 16:07:42 -0400 --* Date: Fri, 20 Aug 1993 16:07:42 -0400 --* From: gianni@matthew.watson.ibm.com --* Message-Id: <9308202007.AA21471@matthew.watson.ibm.com> --* To: asbugs@watson.ibm.com --* Subject: Bug: sstNext: too many simultaneous traversals [listgcd.as][27.9 (current)] --@ Fixed by: SMW Sat Oct 09 14:43:37 1993 --@ Tested by: reclist0.as --@ Summary: New sefo marking scheme. +++ This package provides the functions for the heuristic integer gcd, +++ Geddes's algorithm,for univariate polynomials with integer coefficients #include "aslib.as" #library DemoLib "asdem" import DemoLib --HeuGcd (BP):C == T (HeuGcd ():C == T) where -- BP : UnivariatePolynomialCategory(Integer) BP ==> Polynomial(Integer,Integer) ContPrim ==> Record(cont:Integer,prim:BP) Z ==> Integer C ==> with gcd : List(BP) -> BP gcdprim : List(BP) -> BP gcdcofact : List(BP) -> List(BP) gcdcofactprim: List(BP) -> List(BP) content : List(BP) -> List(Z) contprim : List(BP) -> List(ContPrim) lintgcd : List(Z) -> Z T ==> add PI ==> PositiveInteger NNI ==> NonNegativeInteger Cases ==> Union("gcdprim","gcd","gcdcofactprim","gcdcofact") --local functions negShiftz : (Z,PI) -> Z height : BP -> PI constcase : (BP,List(NNI),List(BP)) -> List(BP) lincase : (BP,List(NNI),List(BP)) -> List(BP) genpoly : (Z,PI) -> BP zerocase : (BP,List(BP)) -> List(BP) internal : (Cases,List(BP)) -> List(BP) negShiftz(n:Z,Modulus:PI):Z == n < 0 => n:= n+Modulus n > (Modulus quo 2) => n-Modulus n --compute the height of a polynomial height(f:BP):PI == k:PI:=1 while f^=0 repeat k:=max(k,abs(leadingCoefficient(f)@Z)::PI) f:=reductum f k --reconstruct the polynomial from the value-adic representation of --dval. genpoly(dval:Z,value:PI):BP == d:=0$BP val:=dval for i in 0..1000 while (val^=0) repeat val1:=negShiftz(val rem value,value) d:= d+monomial(val1,i) val:=(val-val1) quo value d --gcd of a list of integers lintgcd(lval:List(Z)):Z == member?(1,lval) => 1$Z val:=lval.first for val1 in lval.rest repeat val:=gcd(val,val1) val --content for a list of univariate polynomials content(listf:List(BP)):List(Z) == [content f for f in listf] --content of a list of polynomials with the relative primitive parts contprim(listf:List(BP)):List(ContPrim) == [[c:=content f,(f exquo c)::BP]$ContPrim for f in listf] --the list contains the zero polynomial zerocase(pol:BP,listf:List(BP)):List(BP)== n:=position(0$BP,listf) #listf=2 => g:=listf.(1-n) g^=0$BP => [pol*g,(listf.1 exquo g)::BP,(listf.2 exquo g)::BP] [0$BP,0$BP,0$BP] -- listgcd:=internal("gcdcofactprim",[:take(listf,n),:drop(listf,n+1)]) -- I hope I have this right .... JHD 2.1.90 listgcd:=internal("gcdcofactprim",delete(listf,n)) if pol^=1$BP then listgcd.1:=listgcd.0*pol -- [:take(listgcd,n+1),0$BP,:drop(listgcd,n+1)] insert(0$BP,listgcd,n) --one polynomial is constant, remark that they are primitive constcase(pol:BP,listdeg:List(NNI),listf:List(BP)): List(BP) == reduce(\/,[n>0 for n in listdeg]) => cons(pol,listf) lclistf:List(Z):= [leadingCoefficient f for f in listf] d:=lintgcd(lclistf) d=1 => cons(pol,listf) cons(d*pol,[(lcf quo d)::BP for lcf in lclistf]) --one polynomial is linear, remark that they are primitive lincase(pol:BP,listdeg:List(NNI),listf:List(BP)):List(BP) == n:= position(1,listdeg) g:=listf.n flag:=true result:=[g] for f in listf while flag repeat (f1:=f exquo g) case "failed" => flag:=false result := cons(f1::BP,result) ~flag => cons(pol,listf) result:=reverse(result) if pol^=1$BP then result.1:=pol*result.1 result --local function for the gcd among n PRIMITIVE univariate polynomials localgcd(listf:List(BP)):List(BP) == hgt:="min"/[height(f) for f in listf] answr:=2+2*hgt for k in 0..16 repeat listval:=[f answr for f in listf] dval:=lintgcd(listval) dd:=genpoly(dval,answr) contd:=content(dd) d:=(dd exquo contd)::BP result:List(BP):=[d] flag:=true for f in listf while flag repeat (f1:=f exquo d) case "failed" => flag:=false result := cons (f1::BP,result) if flag then return reverse(result) nvalue:= answr*832040 quo 317811 if ((nvalue + answr) rem 2) = 0 then nvalue:=nvalue+1 answr:=nvalue::PI error "try with other gcd" --internal function:it evaluates the gcd and avoids duplication of --code. internal(flag:Cases,listf:List(BP)):List(BP) == --special cases listf=[] => [1$BP] #listf=1 => [first listf,1$BP] minpol:=1$BP mdeg:= "min"/[minimumDegree f for f in listf] if mdeg>0 then minpol1:= monomial(1,mdeg) listf:= [(f exquo minpol1)::BP for f in listf] minpol:=minpol*minpol1 listdeg:=[degree f for f in listf ] --make f and g primitive Cgcd:List(Z):=[] if (flag case "gcd") or (flag case "gcdcofact") then contlistf:List(ContPrim):=contprim(listf) Cgcd:= [term.cont for term in contlistf] contgcd:=lintgcd(Cgcd) Cgcd:=[(n quo contgcd) for n in Cgcd] listf:List(BP):=[term.prim for term in contlistf] minpol:=contgcd*minpol ans:List(BP):= --zerocase, constcase, lincase rewrite without pol member?(0$BP,listf) => zerocase(1$BP,listf) --one polynomial is constant member?(0,listdeg) => constcase(1$BP,listdeg,listf) --one polynomial is linear member?(1,listdeg) => lincase(1$BP,listdeg,listf) localgcd(listf) (result,ans):=(first ans*minpol,rest ans) (flag case "gcdprim") or (flag case "gcd") => [result] flag case "gcdcofactprim" => cons(result,ans) cons(result,[p*q for p in Cgcd for q in ans]) --gcd among n PRIMITIVE univariate polynomials gcdprim (listf:List(BP)):BP == first internal("gcdprim",listf) --gcd and cofactors for n PRIMITIVE univariate polynomials gcdcofactprim(listf:List(BP)):List(BP) == internal("gcdcofactprim",listf) --gcd for n generic univariate polynomials. gcd(listf:List(BP)): BP == first internal("gcd",listf) --gcd and cofactors for n generic univariate polynomials. gcdcofact (listf:List(BP)):List(BP) == internal("gcdcofact",listf)