path: root/src/algebra/gbintern.spad.pamphlet

\documentclass{article}
\usepackage{open-axiom}
\begin{document}
\title{\$SPAD/src/algebra gbintern.spad}
\author{The Axiom Team}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject
\section{package GBINTERN GroebnerInternalPackage}
<<package GBINTERN GroebnerInternalPackage>>=
)abbrev package GBINTERN GroebnerInternalPackage
++ Author: 
++ Date Created: 
++ Date Last Updated: 
++ Keywords: 
++ Description
++ This package provides low level tools for Groebner basis computations
GroebnerInternalPackage(Dom, Expon, VarSet, Dpol): T == C where
 Dom:   GcdDomain
 Expon: OrderedAbelianMonoidSup
 VarSet: OrderedSet
 Dpol:  PolynomialCategory(Dom, Expon, VarSet)
 NNI    ==> NonNegativeInteger
   ------  Definition of Record critPair and Prinp

 critPair ==> Record( lcmfij: Expon, totdeg: NonNegativeInteger,
                      poli: Dpol, polj: Dpol )
 sugarPol ==> Record( totdeg: NonNegativeInteger, pol : Dpol)
 Prinp    ==> Record( ci:Dpol,tci:Integer,cj:Dpol,tcj:Integer,c:Dpol,
                tc:Integer,rc:Dpol,trc:Integer,tF:Integer,tD:Integer)
 Prinpp   ==> Record( ci:Dpol,tci:Integer,cj:Dpol,tcj:Integer,c:Dpol,
                tc:Integer,rc:Dpol,trc:Integer,tF:Integer,tDD:Integer,
                 tDF:Integer)
 T== with

     credPol:  (Dpol, List(Dpol))  -> Dpol
	++ credPol \undocumented
     redPol:   (Dpol, List(Dpol))  -> Dpol
	++ redPol \undocumented
     gbasis:  (List(Dpol), Integer, Integer) -> List(Dpol)
	++ gbasis \undocumented
     critT:  critPair   -> Boolean
	++ critT \undocumented
     critM:  (Expon, Expon) -> Boolean
	++ critM \undocumented
     critB:  (Expon, Expon, Expon, Expon) -> Boolean
	++ critB \undocumented
     critBonD:  (Dpol, List(critPair)) -> List(critPair)
	++ critBonD \undocumented
     critMTonD1: (List(critPair)) -> List(critPair)
	++ critMTonD1 \undocumented
     critMonD1:  (Expon, List(critPair)) -> List(critPair)
	++ critMonD1 \undocumented
     redPo: (Dpol, List(Dpol) )  ->  Record(poly:Dpol, mult:Dom)
	++ redPo \undocumented
     hMonic:  Dpol  -> Dpol
	++ hMonic \undocumented
     updatF: (Dpol, NNI, List(sugarPol) ) -> List(sugarPol)
	++ updatF \undocumented
     sPol:  critPair  -> Dpol
	++ sPol \undocumented
     updatD: (List(critPair), List(critPair)) -> List(critPair)
	++ updatD \undocumented
     minGbasis: List(Dpol) -> List(Dpol)
	++ minGbasis \undocumented
     lepol: Dpol -> Integer
	++ lepol \undocumented
     prinshINFO : Dpol -> Void
	++ prinshINFO \undocumented
     prindINFO: (critPair, Dpol, Dpol,Integer,Integer,Integer) -> Integer
	++ prindINFO \undocumented
     fprindINFO: (critPair, Dpol, Dpol, Integer,Integer,Integer
                 ,Integer) ->  Integer
	++ fprindINFO \undocumented
     prinpolINFO: List(Dpol) -> Void
	++ prinpolINFO \undocumented
     prinb: Integer-> Void
	++ prinb \undocumented
     critpOrder: (critPair, critPair) -> Boolean
	++ critpOrder \undocumented
     makeCrit: (sugarPol, Dpol, NonNegativeInteger) -> critPair
	++ makeCrit \undocumented
     virtualDegree : Dpol -> NonNegativeInteger
	++ virtualDegree \undocumented

 C== add
   Ex ==> OutputForm
   import OutputForm

   ------  Definition of intermediate functions
   if Dpol has totalDegree: Dpol -> NonNegativeInteger then
     virtualDegree p == totalDegree p
   else
     virtualDegree p == 0

   ------  ordering of critpairs

   critpOrder(cp1,cp2) ==
     cp1.totdeg < cp2.totdeg => true
     cp2.totdeg < cp1.totdeg => false
     cp1.lcmfij < cp2.lcmfij

   ------    creating a critical pair

   makeCrit(sp1, p2, totdeg2) ==
     p1 := sp1.pol
     deg := sup(degree(p1), degree(p2))
     e1 := subtractIfCan(deg, degree(p1))::Expon
     e2 := subtractIfCan(deg, degree(p2))::Expon
     tdeg := max(sp1.totdeg + virtualDegree(monomial(1,e1)),
                 totdeg2 + virtualDegree(monomial(1,e2)))
     [deg, tdeg, p1, p2]$critPair

   ------    calculate basis

   gbasis(Pol: List(Dpol), xx1: Integer, xx2: Integer ) ==
     D, D1: List(critPair)
     ---------   create D and Pol

     Pol1:= sort(degree #1 > degree #2, Pol)
     basPols:= updatF(hMonic(first Pol1),virtualDegree(first Pol1),[])
     Pol1:= rest(Pol1)
     D:= nil
     while not null Pol1 repeat
        h:= hMonic(first(Pol1))
        Pol1:= rest(Pol1)
        toth := virtualDegree h
        D1:= [makeCrit(x,h,toth) for x in basPols]
        D:= updatD(critMTonD1(sort(critpOrder, D1)),
                   critBonD(h,D))
        basPols:= updatF(h,toth,basPols)
     D:= sort(critpOrder, D)
     xx:= xx2
     --------  loop

     redPols := [x.pol for x in basPols]
     while not null D repeat
         D0:= first D
         s:= hMonic(sPol(D0))
         D:= rest(D)
         h:= hMonic(redPol(s,redPols))
         if xx1 = 1  then
              prinshINFO(h)
         h = 0  =>
          if xx2 = 1 then
           prindINFO(D0,s,h,# basPols, # D,xx)
           xx:= 2
          " go to top of while "
         degree(h) = 0 =>
           D:= nil
           if xx2 = 1 then
            prindINFO(D0,s,h,# basPols, # D,xx)
            xx:= 2
           basPols:= updatF(h,0,[])
           leave "out of while"
         D1:= [makeCrit(x,h,D0.totdeg) for x in basPols]
         D:= updatD(critMTonD1(sort(critpOrder, D1)),
                   critBonD(h,D))
         basPols:= updatF(h,D0.totdeg,basPols)
         redPols := concat(redPols,h)
         if xx2 = 1 then
            prindINFO(D0,s,h,# basPols, # D,xx)
            xx:= 2
     Pol := [x.pol for x in basPols]
     if xx2 = 1 then
       prinpolINFO(Pol)
       messagePrint("    THE GROEBNER BASIS POLYNOMIALS")
     if xx1 = 1 and not one? xx2 then
       messagePrint("    THE GROEBNER BASIS POLYNOMIALS")
     Pol

             --------------------------------------

             --- erase multiple of e in D2 using crit M

   critMonD1(e: Expon, D2: List(critPair))==
      null D2 => nil
      x:= first(D2)
      critM(e, x.lcmfij) => critMonD1(e, rest(D2))
      cons(x, critMonD1(e, rest(D2)))

             ----------------------------

             --- reduce D1 using crit T and crit M

   critMTonD1(D1: List(critPair))==
           null D1 => nil
           f1:= first(D1)
           s1:= #(D1)
           cT1:= critT(f1)
           s1= 1 and cT1 => nil
           s1= 1 => D1
           e1:= f1.lcmfij
           r1:= rest(D1)
           e1 = (first r1).lcmfij  =>
              cT1 =>   critMTonD1(cons(f1, rest(r1)))
              critMTonD1(r1)
           D1 := critMonD1(e1, r1)
           cT1 => critMTonD1(D1)
           cons(f1, critMTonD1(D1))

             -----------------------------

             --- erase elements in D fullfilling crit B

   critBonD(h:Dpol, D: List(critPair))==
         null D => nil
         x:= first(D)
         critB(degree(h), x.lcmfij, degree(x.poli), degree(x.polj)) =>
           critBonD(h, rest(D))
         cons(x, critBonD(h, rest(D)))

             -----------------------------

             --- concat F and h and erase multiples of h in F

   updatF(h: Dpol, deg:NNI, F: List(sugarPol)) ==
       null F => [[deg,h]]
       f1:= first(F)
       critM(degree(h), degree(f1.pol))  => updatF(h, deg, rest(F))
       cons(f1, updatF(h, deg, rest(F)))

             -----------------------------

             --- concat ordered critical pair lists D1 and D2

   updatD(D1: List(critPair), D2: List(critPair)) ==
      null D1 => D2
      null D2 => D1
      dl1:= first(D1)
      dl2:= first(D2)
      critpOrder(dl1,dl2) => cons(dl1, updatD(D1.rest, D2))
      cons(dl2, updatD(D1, D2.rest))

            -----------------------------

            --- remove gcd from pair of coefficients

   gcdCo(c1:Dom, c2:Dom):Record(co1:Dom,co2:Dom) ==
      d:=gcd(c1,c2)
      [(c1 exquo d)::Dom, (c2 exquo d)::Dom]

            --- calculate S-polynomial of a critical pair

   sPol(p:critPair)==
      Tij := p.lcmfij
      fi := p.poli
      fj := p.polj
      cc := gcdCo(leadingCoefficient fi, leadingCoefficient fj)
      reductum(fi)*monomial(cc.co2,subtractIfCan(Tij, degree fi)::Expon) -
        reductum(fj)*monomial(cc.co1,subtractIfCan(Tij, degree fj)::Expon)

            ----------------------------

            --- reduce critpair polynomial mod F
            --- iterative version

   redPo(s: Dpol, F: List(Dpol)) ==
      m:Dom := 1
      Fh := F
      while not ( s = 0 or null F ) repeat
        f1:= first(F)
        s1:= degree(s)
        e: Union(Expon, "failed")
        (e:= subtractIfCan(s1, degree(f1))) case Expon  =>
           cc:=gcdCo(leadingCoefficient f1, leadingCoefficient s)
           s:=cc.co1*reductum(s) - monomial(cc.co2,e)*reductum(f1)
           m := m*cc.co1
           F:= Fh
        F:= rest F
      [s,m]

   redPol(s: Dpol, F: List(Dpol)) ==  credPol(redPo(s,F).poly,F)

            ----------------------------

            --- crit T  true, if e1 and e2 are disjoint

   critT(p: critPair) == p.lcmfij =  (degree(p.poli) + degree(p.polj))

            ----------------------------

            --- crit M - true, if lcm#2 multiple of lcm#1

   critM(e1: Expon, e2: Expon) ==
         en: Union(Expon, "failed")
         (en:=subtractIfCan(e2, e1)) case Expon

            ----------------------------

            --- crit B - true, if eik is a multiple of eh and eik not equal
            ---          lcm(eh,ei) and eik not equal lcm(eh,ek)

   critB(eh:Expon, eik:Expon, ei:Expon, ek:Expon) ==
       critM(eh, eik) and (eik ~= sup(eh, ei)) and (eik ~= sup(eh, ek))

            ----------------------------

            ---  make polynomial monic case Domain a Field

   hMonic(p: Dpol) ==
        p= 0 => p
        -- inv(leadingCoefficient(p))*p
        primitivePart p

            -----------------------------

            ---  reduce all terms of h mod F  (iterative version )

   credPol(h: Dpol, F: List(Dpol) ) ==
        null F => h
        h0:Dpol:= monomial(leadingCoefficient h, degree h)
        while (h:=reductum h) ~= 0 repeat
           hred:= redPo(h, F)
           h := hred.poly
           h0:=(hred.mult)*h0 + monomial(leadingCoefficient(h),degree h)
        h0

            -------------------------------

            ----  calculate minimal basis for ordered F

   minGbasis(F: List(Dpol)) ==
        null F => nil
        newbas := minGbasis rest F
        cons(hMonic credPol( first(F), newbas),newbas)

            -------------------------------

            ----  calculate number of terms of polynomial

   lepol(p1:Dpol)==
      n: Integer
      n:= 0
      while p1 ~= 0 repeat
         n:= n + 1
         p1:= reductum(p1)
      n

            ----  print blanc lines

   prinb(n: Integer)==
      for x in 1..n repeat
         messagePrint("    ")

            ----  print reduced critpair polynom

   prinshINFO(h: Dpol)==
           prinb(2)
           messagePrint(" reduced Critpair - Polynom :")
           prinb(2)
           print(h::Ex)
           prinb(2)

            -------------------------------

            ----  print info string

   prindINFO(cp: critPair, ps: Dpol, ph: Dpol, i1:Integer,
             i2:Integer, n:Integer) ==
       ll: List Prinp
       a: Dom
       cpi:= cp.poli
       cpj:= cp.polj
       if n = 1 then
        prinb(1)
        messagePrint("you choose option  -info-  ")
        messagePrint("abbrev. for the following information strings are")
        messagePrint("  ci  =>  Leading monomial  for critpair calculation")
        messagePrint("  tci =>  Number of terms of polynomial i")
        messagePrint("  cj  =>  Leading monomial  for critpair calculation")
        messagePrint("  tcj =>  Number of terms of polynomial j")
        messagePrint("  c   =>  Leading monomial of critpair polynomial")
        messagePrint("  tc  =>  Number of terms of critpair polynomial")
        messagePrint("  rc  =>  Leading monomial of redcritpair polynomial")
        messagePrint("  trc =>  Number of terms of redcritpair polynomial")
        messagePrint("  tF  =>  Number of polynomials in reduction list F")
        messagePrint("  tD  =>  Number of critpairs still to do")
        prinb(4)
        n:= 2
       prinb(1)
       a:= 1
       ph = 0  =>
          ps = 0 =>
            ll:= [[monomial(a,degree(cpi)),lepol(cpi),
                  monomial(a,degree(cpj)),
                   lepol(cpj),ps,0,ph,0,i1,i2]$Prinp]
            print(ll::Ex)
            prinb(1)
            n
          ll:= [[monomial(a,degree(cpi)),lepol(cpi),
             monomial(a,degree(cpj)),lepol(cpj),monomial(a,degree(ps)),
              lepol(ps), ph,0,i1,i2]$Prinp]
          print(ll::Ex)
          prinb(1)
          n
       ll:= [[monomial(a,degree(cpi)),lepol(cpi),
            monomial(a,degree(cpj)),lepol(cpj),monomial(a,degree(ps)),
             lepol(ps),monomial(a,degree(ph)),lepol(ph),i1,i2]$Prinp]
       print(ll::Ex)
       prinb(1)
       n

            -------------------------------

            ----  print the groebner basis polynomials

   prinpolINFO(pl: List(Dpol))==
       n:Integer
       n:= # pl
       prinb(1)
       n = 1 =>
         messagePrint("  There is 1  Groebner Basis Polynomial ")
         prinb(2)
       messagePrint("  There are ")
       prinb(1)
       print(n::Ex)
       prinb(1)
       messagePrint("  Groebner Basis Polynomials. ")
       prinb(2)

   fprindINFO(cp: critPair, ps: Dpol, ph: Dpol, i1:Integer,
             i2:Integer, i3:Integer, n: Integer) ==
       ll: List Prinpp
       a: Dom
       cpi:= cp.poli
       cpj:= cp.polj
       if n = 1 then
        prinb(1)
        messagePrint("you choose option  -info-  ")
        messagePrint("abbrev. for the following information strings are")
        messagePrint("  ci  =>  Leading monomial  for critpair calculation")
        messagePrint("  tci =>  Number of terms of polynomial i")
        messagePrint("  cj  =>  Leading monomial  for critpair calculation")
        messagePrint("  tcj =>  Number of terms of polynomial j")
        messagePrint("  c   =>  Leading monomial of critpair polynomial")
        messagePrint("  tc  =>  Number of terms of critpair polynomial")
        messagePrint("  rc  =>  Leading monomial of redcritpair polynomial")
        messagePrint("  trc =>  Number of terms of redcritpair polynomial")
        messagePrint("  tF  =>  Number of polynomials in reduction list F")
        messagePrint("  tD  =>  Number of critpairs still to do")
        messagePrint("  tDF =>  Number of subproblems still to do")
        prinb(4)
        n:= 2
       prinb(1)
       a:= 1
       ph = 0  =>
          ps = 0 =>
            ll:= [[monomial(a,degree(cpi)),lepol(cpi),
              monomial(a,degree(cpj)),
               lepol(cpj),ps,0,ph,0,i1,i2,i3]$Prinpp]
            print(ll::Ex)
            prinb(1)
            n
          ll:= [[monomial(a,degree(cpi)),lepol(cpi),
            monomial(a,degree(cpj)),lepol(cpj),monomial(a,degree(ps)),
             lepol(ps), ph,0,i1,i2,i3]$Prinpp]
          print(ll::Ex)
          prinb(1)
          n
       ll:= [[monomial(a,degree(cpi)),lepol(cpi),
            monomial(a,degree(cpj)),lepol(cpj),monomial(a,degree(ps)),
             lepol(ps),monomial(a,degree(ph)),lepol(ph),i1,i2,i3]$Prinpp]
       print(ll::Ex)
       prinb(1)
       n

@
\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 GBINTERN GroebnerInternalPackage>>
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}