\documentclass{article} \usepackage{open-axiom} \begin{document} \title{\$SPAD/src/algebra idecomp.spad} \author{Patrizia Gianni} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{package IDECOMP IdealDecompositionPackage} <>= )abbrev package IDECOMP IdealDecompositionPackage ++ Author: P. Gianni ++ Date Created: summer 1986 ++ Date Last Updated: ++ Basic Functions: ++ Related Constructors: PolynomialIdeals ++ Also See: ++ AMS Classifications: ++ Keywords: ++ References: ++ Description: ++ This package provides functions for the primary decomposition of ++ polynomial ideals over the rational numbers. The ideals are members ++ of the \spadtype{PolynomialIdeals} domain, and the polynomial generators are ++ required to be from the \spadtype{DistributedMultivariatePolynomial} domain. IdealDecompositionPackage(vl,nv) : C == T -- take away nv, now doesn't -- compile if it isn't there where vl : List Symbol nv : NonNegativeInteger Z ==> Integer -- substitute with PFE cat Q ==> Fraction Z F ==> Fraction P P ==> Polynomial Z UP ==> SparseUnivariatePolynomial P Expon ==> DirectProduct(nv,NNI) OV ==> OrderedVariableList(vl) SE ==> Symbol SUP ==> SparseUnivariatePolynomial(DPoly) DPoly1 ==> DistributedMultivariatePolynomial(vl,Q) DPoly ==> DistributedMultivariatePolynomial(vl,F) NNI ==> NonNegativeInteger Ideal == PolynomialIdeals(Q,Expon,OV,DPoly1) FIdeal == PolynomialIdeals(F,Expon,OV,DPoly) Fun0 == Union("zeroPrimDecomp","zeroRadComp") GenPos == Record(changeval:List Z,genideal:FIdeal) C == with zeroDimPrime? : Ideal -> Boolean ++ zeroDimPrime?(I) tests if the ideal I is a 0-dimensional prime. zeroDimPrimary? : Ideal -> Boolean ++ zeroDimPrimary?(I) tests if the ideal I is 0-dimensional primary. prime? : Ideal -> Boolean ++ prime?(I) tests if the ideal I is prime. radical : Ideal -> Ideal ++ radical(I) returns the radical of the ideal I. primaryDecomp : Ideal -> List(Ideal) ++ primaryDecomp(I) returns a list of primary ideals such that their ++ intersection is the ideal I. contract : (Ideal,List OV ) -> Ideal ++ contract(I,lvar) contracts the ideal I to the polynomial ring ++ \spad{F[lvar]}. T == add import MPolyCatRationalFunctionFactorizer(Expon,OV,Z,DPoly) import GroebnerPackage(F,Expon,OV,DPoly) import GroebnerPackage(Q,Expon,OV,DPoly1) ---- Local Functions ----- genPosLastVar : (FIdeal,List OV) -> GenPos zeroPrimDecomp : (FIdeal,List OV) -> List(FIdeal) zeroRadComp : (FIdeal,List OV) -> FIdeal zerodimcase : (FIdeal,List OV) -> Boolean is0dimprimary : (FIdeal,List OV) -> Boolean backGenPos : (FIdeal,List Z,List OV) -> FIdeal reduceDim : (Fun0,FIdeal,List OV) -> List FIdeal findvar : (FIdeal,List OV) -> OV testPower : (SUP,OV,FIdeal) -> Boolean goodPower : (DPoly,FIdeal) -> Record(spol:DPoly,id:FIdeal) pushdown : (DPoly,OV) -> DPoly pushdterm : (DPoly,OV,Z) -> DPoly pushup : (DPoly,OV) -> DPoly pushuterm : (DPoly,SE,OV) -> DPoly pushucoef : (UP,OV) -> DPoly trueden : (P,SE) -> P rearrange : (List OV) -> List OV deleteunit : List FIdeal -> List FIdeal ismonic : (DPoly,OV) -> Boolean MPCFQF ==> MPolyCatFunctions2(OV,Expon,Expon,Q,F,DPoly1,DPoly) MPCFFQ ==> MPolyCatFunctions2(OV,Expon,Expon,F,Q,DPoly,DPoly1) convertQF(a:Q) : F == ((numer a):: F)/((denom a)::F) convertFQ(a:F) : Q == (ground numer a)/(ground denom a) internalForm(I:Ideal) : FIdeal == Id:=generators I nId:=[map(convertQF,poly)$MPCFQF for poly in Id] groebner? I => groebnerIdeal nId ideal nId externalForm(I:FIdeal) : Ideal == Id:=generators I nId:=[map(convertFQ,poly)$MPCFFQ for poly in Id] groebner? I => groebnerIdeal nId ideal nId lvint:=[variable(xx)::OV for xx in vl] nvint1:=(#lvint-1)::NNI deleteunit(lI: List FIdeal) : List FIdeal == [I for I in lI | not element?(1$DPoly,I)] rearrange(vlist:List OV) :List OV == vlist=[] => vlist sort(#1>#2,setDifference(lvint,setDifference(lvint,vlist))) ---- radical of a 0-dimensional ideal ---- zeroRadComp(I:FIdeal,truelist:List OV) : FIdeal == truelist=[] => I Id:=generators I x:OV:=truelist.last #Id=1 => f:=Id.first g:= (f exquo (gcd (f,differentiate(f,x))))::DPoly groebnerIdeal([g]) y:=truelist.first px:DPoly:=x::DPoly py:DPoly:=y::DPoly f:=Id.last g:= (f exquo (gcd (f,differentiate(f,x))))::DPoly Id:=groebner(cons(g,remove(f,Id))) lf:=Id.first pv:DPoly:=0 pw:DPoly:=0 while not one? degree(lf,y) repeat val:=random(23)$Z pv:=px+val*py pw:=px-val*py Id:=groebner([(univariate(h,x)).pv for h in Id]) lf:=Id.first ris:= generators(zeroRadComp(groebnerIdeal(Id.rest),truelist.rest)) ris:=cons(lf,ris) if not zero? pv then ris:=[(univariate(h,x)).pw for h in ris] groebnerIdeal(groebner ris) ---- find the power that stabilizes (I:s) ---- goodPower(s:DPoly,I:FIdeal) : Record(spol:DPoly,id:FIdeal) == f:DPoly:=s I:=groebner I J:=generators(JJ:= (saturate(I,s))) while not in?(ideal([f*g for g in J]),I) repeat f:=s*f [f,JJ] ---- is the ideal zerodimensional? ---- ---- the "true variables" are in truelist ---- zerodimcase(J:FIdeal,truelist:List OV) : Boolean == element?(1,J) => true truelist=[] => true n:=#truelist Jd:=groebner generators J for x in truelist while Jd~=[] repeat f := Jd.first Jd:=Jd.rest if ((y:=mainVariable f) case "failed") or (y::OV ~=x ) or not (ismonic (f,x)) then return false while Jd~=[] and (mainVariable Jd.first)::OV=x repeat Jd:=Jd.rest if Jd=[] and position(x,truelist) [J] zerodimcase(J,truelist) => (flag case "zeroPrimDecomp") => zeroPrimDecomp(J,truelist) (flag case "zeroRadComp") => [zeroRadComp(J,truelist)] x:OV:=findvar(J,truelist) Jnew:=[pushdown(f,x) for f in generators J] Jc: List FIdeal :=[] Jc:=reduceDim(flag,groebnerIdeal Jnew,remove(x,truelist)) res1:=[ideal([pushup(f,x) for f in generators idp]) for idp in Jc] s:=pushup((_*/[leadingCoefficient f for f in Jnew])::DPoly,x) degree(s,x)=0 => res1 res1:=[saturate(II,s) for II in res1] good:=goodPower(s,J) sideal := groebnerIdeal(groebner(cons(good.spol,generators J))) in?(good.id, sideal) => res1 sresult:=reduceDim(flag,sideal,truelist) for JJ in sresult repeat if not(in?(good.id,JJ)) then res1:=cons(JJ,res1) res1 ---- Primary Decomposition for 0-dimensional ideals ---- zeroPrimDecomp(I:FIdeal,truelist:List OV): List(FIdeal) == truelist=[] => list I newJ:=genPosLastVar(I,truelist);lval:=newJ.changeval; J:=groebner newJ.genideal x:=truelist.last Jd:=generators J g:=Jd.last lfact:= factors factor(g) ris:List FIdeal:=[] for ef in lfact repeat g:DPoly:=(ef.factor)**(ef.exponent::NNI) J1:= groebnerIdeal(groebner cons(g,Jd)) if not (is0dimprimary (J1,truelist)) then return zeroPrimDecomp(I,truelist) ris:=cons(groebner backGenPos(J1,lval,truelist),ris) ris ---- radical of an Ideal ---- radical(I:Ideal) : Ideal == J:=groebner(internalForm I) truelist:=rearrange("setUnion"/[variables f for f in generators J]) truelist=[] => externalForm J externalForm("intersect"/reduceDim("zeroRadComp",J,truelist)) -- the following functions are used to "push" x in the coefficient ring - ---- push x in the coefficient domain for a polynomial ---- pushdown(g:DPoly,x:OV) : DPoly == rf:DPoly:=0$DPoly i:=position(x,lvint) while not zero? g repeat g1:=reductum g rf:=rf+pushdterm(g-g1,x,i) g := g1 rf ---- push x in the coefficient domain for a term ---- pushdterm(t:DPoly,x:OV,i:Z):DPoly == n:=degree(t,x) xp:=convert(x)@SE cf:=monomial(1,xp,n)$P :: F newt := t exquo monomial(1,x,n)$DPoly cf * newt::DPoly ---- push back the variable ---- pushup(f:DPoly,x:OV) :DPoly == h:=1$P rf:DPoly:=0$DPoly g := f xp := convert(x)@SE while not zero? g repeat h:=lcm(trueden(denom leadingCoefficient g,xp),h) g:=reductum g f:=(h::F)*f while not zero? f repeat g:=reductum f rf:=rf+pushuterm(f-g,xp,x) f:=g rf trueden(c:P,x:SE) : P == degree(c,x) = 0 => 1 c ---- push x back from the coefficient domain for a term ---- pushuterm(t:DPoly,xp:SE,x:OV):DPoly == pushucoef((univariate(numer leadingCoefficient t,xp)$P), x)* monomial(inv((denom leadingCoefficient t)::F),degree t)$DPoly pushucoef(c:UP,x:OV):DPoly == c = 0 => 0 monomial((leadingCoefficient c)::F::DPoly,x,degree c) + pushucoef(reductum c,x) -- is the 0-dimensional ideal I primary ? -- ---- internal function ---- is0dimprimary(J:FIdeal,truelist:List OV) : Boolean == element?(1,J) => true Jd:=generators(groebner J) not one?(#(factors factor Jd.last)) => return false i:=subtractIfCan(#truelist,1) (i case nothing) => return true k: NNI := i JR:=(reverse Jd);JM:=groebnerIdeal([JR.first]);JP:List(DPoly):=[] for f in JR.rest repeat if not ismonic(f,truelist.k) then if not inRadical?(f,JM) then return false JP:=cons(f,JP) else x:=truelist.k k := (k-1)::NNI if not testPower(univariate(f,x),x,JM) then return false JM :=groebnerIdeal(append(cons(f,JP),generators JM)) true ---- Functions for the General Position step ---- ---- put the ideal in general position ---- genPosLastVar(J:FIdeal,truelist:List OV):GenPos == x := last truelist ;lv1:List OV :=remove(x,truelist) ranvals:List(Z):=[random(23)$Z for vv in lv1] val:=+/[rv*(vv::DPoly) for vv in lv1 for rv in ranvals] val:=val+(x::DPoly) [ranvals,groebnerIdeal(groebner([(univariate(p,x)).val for p in generators J]))]$GenPos ---- convert back the ideal ---- backGenPos(I:FIdeal,lval:List Z,truelist:List OV) : FIdeal == lval=[] => I x := last truelist ;lv1:List OV:=remove(x,truelist) val:=-(+/[rv*(vv::DPoly) for vv in lv1 for rv in lval]) val:=val+(x::DPoly) groebnerIdeal (groebner([(univariate(p,x)).val for p in generators I ])) ismonic(f:DPoly,x:OV) : Boolean == ground? leadingCoefficient(univariate(f,x)) ---- test if f is power of a linear mod (rad J) ---- ---- f is monic ---- testPower(uf:SUP,x:OV,J:FIdeal) : Boolean == df:=degree(uf) trailp:DPoly := inv(df:Z ::F) *coefficient(uf,(df-1)::NNI) linp:SUP:=(monomial(1$DPoly,1$NNI)$SUP + monomial(trailp,0$NNI)$SUP)**df g:DPoly:=multivariate(uf-linp,x) inRadical?(g,J) ---- Exported Functions ---- -- is the 0-dimensional ideal I prime ? -- zeroDimPrime?(I:Ideal) : Boolean == J:=groebner((genPosLastVar(internalForm I,lvint)).genideal) element?(1,J) => true n:NNI:=#vl;i:NNI:=1 Jd:=generators J #Jd~=n => false for f in Jd repeat if not ismonic(f,lvint.i) then return false if i1 => false lfact.1.exponent =1 -- is the 0-dimensional ideal I primary ? -- zeroDimPrimary?(J:Ideal):Boolean == is0dimprimary(internalForm J,lvint) ---- Primary Decomposition of I ----- primaryDecomp(I:Ideal) : List(Ideal) == J:=groebner(internalForm I) truelist:=rearrange("setUnion"/[variables f for f in generators J]) truelist=[] => [externalForm J] [externalForm II for II in reduceDim("zeroPrimDecomp",J,truelist)] ---- contract I to the ring with lvar variables ---- contract(I:Ideal,lvar: List OV) : Ideal == Id:= generators(groebner I) empty?(Id) => I fullVars:= "setUnion"/[variables g for g in Id] fullVars = lvar => I n:= # lvar #fullVars < n => error "wrong vars" n=0 => I newVars:= append([vv for vv in fullVars | not member?(vv,lvar)]$List(OV),lvar) subsVars := [monomial(1,vv,1)$DPoly1 for vv in newVars] lJ:= [eval(g,fullVars,subsVars) for g in Id] J := groebner(lJ) J=[1] => groebnerIdeal J J=[0] => groebnerIdeal empty() J:=[f for f in J| member?(mainVariable(f)::OV,newVars)] fullPol :=[monomial(1,vv,1)$DPoly1 for vv in fullVars] groebnerIdeal([eval(gg,newVars,fullPol) for gg in J]) @ \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}