\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} <>= )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 := subtractIfCan(s1, degree(f1)) e 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 := subtractIfCan(e2, e1) en 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} <>= --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. @ <<*>>= <> <> @ \eject \begin{thebibliography}{99} \bibitem{1} nothing \end{thebibliography} \end{document}