Initial population.

1 files changed, 268 insertions, 0 deletions

diff --git a/src/algebra/listgcd.spad.pamphlet b/src/algebra/listgcd.spad.pamphlet
new file mode 100644
index 00000000..f626adbd
--- /dev/null
+++ b/src/algebra/listgcd.spad.pamphlet
@@ -0,0 +1,268 @@
+\documentclass{article}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/algebra listgcd.spad}
+\author{Patrizia Gianni}
+\maketitle
+\begin{abstract}
+\end{abstract}
+\eject
+\tableofcontents
+\eject
+\section{package HEUGCD HeuGcd}
+<<package HEUGCD HeuGcd>>=
+)abbrev package HEUGCD HeuGcd
+++ Author: P.Gianni
+++ Date Created:
+++ Date Last Updated: 13 September 94
+++ Basic Functions:
+++ Related Constructors:
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description:
+++ This package provides the functions for the heuristic integer gcd.
+++ Geddes's algorithm,for univariate polynomials with integer coefficients
+HeuGcd (BP):C == T
+ where
+  BP       :   UnivariatePolynomialCategory Integer
+  Z        ==> Integer
+  ContPrim ==> Record(cont:Z,prim:BP)
+
+
+  C == with
+     gcd          : List BP  -> BP
+       ++ gcd([f1,..,fk]) = gcd of the polynomials fi.
+     gcdprim      : List BP  -> BP
+       ++ gcdprim([f1,..,fk]) = gcd of k PRIMITIVE univariate polynomials
+     gcdcofact    : List BP  -> List BP
+       ++ gcdcofact([f1,..fk]) = gcd and cofactors of k univariate polynomials.
+     gcdcofactprim: List BP  -> List BP
+       ++ gcdcofactprim([f1,..fk]) = gcd and cofactors of k
+       ++ primitive polynomials.
+     content      : List BP  -> List Z
+       ++ content([f1,..,fk]) = content of a list of univariate polynonials
+     lintgcd      : List  Z  -> Z
+       ++ lintgcd([a1,..,ak]) = gcd of a list of integers
+
+  T == add
+
+    PI    ==> PositiveInteger
+    NNI   ==> NonNegativeInteger
+    Cases ==> Union("gcdprim","gcd","gcdcofactprim","gcdcofact")
+    import ModularDistinctDegreeFactorizer BP
+
+    --local functions
+    localgcd     :        List BP       -> List BP
+    constNotZero :           BP         -> Boolean
+    height       :           BP         -> PI
+    genpoly      :         (Z,PI)       -> BP
+    negShiftz    :         (Z,PI)       -> Z
+    internal     :     (Cases,List BP ) -> List BP
+    constcase    : (List NNI ,List BP ) -> List BP
+    lincase      : (List NNI ,List BP ) -> List BP
+    myNextPrime  :        ( Z , NNI )   -> Z
+
+    bigPrime:= prevPrime(2**26)$IntegerPrimesPackage(Integer)
+
+    myNextPrime(val:Z,bound:NNI) : Z == nextPrime(val)$IntegerPrimesPackage(Z)
+
+    constNotZero(f : BP ) : Boolean == (degree f = 0) and ^(zero? f)
+
+    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..  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 ==
+      empty? lval => 0$Z
+      member?(1,lval) => 1$Z
+      lval:=sort(#1<#2,lval)
+      val:=lval.first
+      for val1 in lval.rest while ^(val=1) repeat val:=gcd(val,val1)
+      val
+
+    --content for a list of univariate polynomials
+    content(listf:List BP ):List(Z) ==
+      [lintgcd coefficients f for f in listf]
+
+    --content of a list of polynomials with the relative primitive parts
+    contprim(listf:List BP ):List(ContPrim) ==
+       [[c:=lintgcd coefficients f,(f exquo c)::BP]$ContPrim  for f in listf]
+
+    -- one polynomial is constant, remark that they are primitive
+    -- but listf can contain the zero polynomial
+    constcase(listdeg:List NNI ,listf:List BP ): List BP  ==
+      lind:=select(constNotZero,listf)
+      empty? lind =>
+        member?(1,listdeg) => lincase(listdeg,listf)
+        localgcd listf
+      or/[n>0 for n in listdeg] => cons(1$BP,listf)
+      lclistf:List(Z):= [leadingCoefficient f for f in listf]
+      d:=lintgcd(lclistf)
+      d=1 =>  cons(1$BP,listf)
+      cons(d::BP,[(lcf quo d)::BP for lcf in lclistf])
+
+    testDivide(listf: List BP, g:BP):Union(List BP, "failed") ==
+      result:List BP := []
+      for f in listf repeat
+        if (f1:=f exquo g) case "failed" then return "failed"
+        result := cons(f1::BP,result)
+      reverse!(result)
+
+    --one polynomial is linear, remark that they are primitive
+    lincase(listdeg:List NNI ,listf:List BP ):List BP  ==
+      n:= position(1,listdeg)
+      g:=listf.n
+      result:=[g]
+      for f in listf repeat
+        if (f1:=f exquo g) case "failed" then return cons(1$BP,listf)
+        result := cons(f1::BP,result)
+      reverse(result)
+
+    IMG := InnerModularGcd(Z,BP,67108859,myNextPrime)
+
+    mindegpol(f:BP, g:BP):BP ==
+      degree(g) < degree (f) => g
+      f
+
+    --local function for the gcd among n PRIMITIVE univariate polynomials
+    localgcd(listf:List BP ):List BP  ==
+      hgt:="min"/[height(f) for f in listf|^zero? f]
+      answr:=2+2*hgt
+      minf := "mindegpol"/[f for f in listf|^zero? f]
+      (result := testDivide(listf, minf)) case List(BP) =>
+           cons(minf, result::List BP)
+      if degree minf < 100 then for k in 1..10 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 : Boolean := 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
+      gg:=modularGcdPrimitive(listf)$IMG
+      cons(gg,[(f exquo gg) :: BP for f in listf])
+
+    --internal function:it evaluates the gcd and avoids duplication of
+    --code.
+    internal(flag:Cases,listf:List BP ):List BP  ==
+      --special cases
+      listf=[] => [1$BP]
+      (nlf:=#listf)=1 => [first listf,1$BP]
+      minpol:=1$BP
+      -- extract a monomial gcd
+      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
+      -- make the polynomials primitive
+      Cgcd:List(Z):=[]
+      contgcd : Z := 1
+      if (flag case "gcd") or (flag case "gcdcofact") then
+        contlistf:List(ContPrim):=contprim(listf)
+        Cgcd:= [term.cont for term in contlistf]
+        contgcd:=lintgcd(Cgcd)
+        listf:List BP :=[term.prim for term in contlistf]
+        minpol:=contgcd*minpol
+      listdeg:=[degree f for f in listf ]
+      f:= first listf
+      for g in rest listf  repeat
+        f:=gcd(f,g,bigPrime)
+        if degree f = 0 then return cons(minpol,listf)
+      ans:List BP :=
+         --one polynomial is constant
+         member?(0,listdeg) => constcase(listdeg,listf)
+         --one polynomial is linear
+         member?(1,listdeg) => lincase(listdeg,listf)
+         localgcd(listf)
+      (result,ans):=(first ans*minpol,rest ans)
+      if (flag case "gcdcofact") then
+        ans:= [(p quo contgcd)*q for p in Cgcd for q in ans]
+      cons(result,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)
+
+@
+\section{License}
+<<license>>=
+--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.
+@
+<<*>>=
+<<license>>
+
+<<package HEUGCD HeuGcd>>
+@
+\eject
+\begin{thebibliography}{99}
+\bibitem{1} nothing
+\end{thebibliography}
+\end{document}