--------------------------------------------------------------------------- -- PURPOSE: Computations of implicit equation of bigraded rational surfaces -- PROGRAMMER : Nicolás Botbol -- UPDATE HISTORY : May 2011 --------------------------------------------------------------------------- newPackage("BigradedImplicit", Version => "1.4", Date => "29 June 2010 => 12 May 2011", Authors => { { Name => "Nicolás Botbol", Email => "nbotbol@dm.uba.ar", HomePage => "http://mate.dm.uba.ar/~nbotbol/"} }, Headline => "Implicit equations of bigraded rational surfaces by means of approximation complexes", DebuggingMode => true ) export { maxMinor, isGoodDegree, degreeImplicitEq, representationMatrix, implicitEq, symbol Numeric, symbol Exact } X := symbol X; lll := symbol lll; m := symbol m; M := symbol M; --------------------------------------------------------------------------- --------------------------------------------------------------------------- -- Functions that end with 'num' (rankNum, maxColNum, maxMinorNum, listDetComplexNum, minorsComplexNum and detComplexNum) are numerical functions, meaning that it does the same computation than the non-num function but after a random evaluation. -- It is perhaps recommendable to run them at least twice. --------------------------------------------------------------------------- --------------------------------------------------------------------------- getComputationStrategy := (options) -> ( exactStr := 0; numericStr := 1; strat := if options.?Strategy then options.Strategy else null; if strat === global Exact then exactStr else if strat === global Numeric then numericStr else if strat =!= null then ( error "'Strategy' keyword must be 'Numeric' or 'Exact'"; ) else exactStr ); --------------------------------------------------------------------------- --------------------------------------------------------------------------- -- PURPOSE : it returns the rank of a matrix obtained after specialization -- of the variables. With high probability, this is the rank of the input matrix. rankNE = method(Options => { Strategy => null }) rankNE(Matrix) := ZZ => options -> (M) -> ( if getComputationStrategy(options) == 0 then return(rank M) else return(rankNum M) ); rankNum = method(); rankNum(Matrix) := (M) -> ( Aring:=ring M; l:={}; i:=0; while i { Strategy => null }) maxCol(Matrix) := (Matrix,List) => options -> (M) -> ( rm:=rankNE(M,options); i:=0; c:={}; while #c { Strategy => null }) maxMinor(Matrix) := (Matrix,List) => options -> (M) -> ( --Apply two times maxCol to obtain a maximal minor --This is not very efficient... return (transpose (maxCol(transpose (maxCol(M,options))_0 ))_0 ) ); --------------------------------------------------------------------------- --------------------------------------------------------------------------- --------------------------------------------------------------------------- --------------------------------------------------------------------------- isGoodDegree = method(); isGoodDegree(List, List) := (polynomialList, nu) -> { F := matrix {polynomialList}; Z0 := (ring F)^1; Z1 := kernel koszul(1,F); Z2 := kernel koszul(2,F); Z3 := kernel koszul(3,F); d := degree (polynomialList_0); eulerChar := hilbertFunction(nu,Z0)-hilbertFunction(nu+d,Z1)+hilbertFunction(nu+2*d,Z2)-hilbertFunction(nu+3*d,Z3); return ((eulerChar == 0) and (nu_0 == nu_1+nu_2)) } --------------------------------------------------------------------------- --------------------------------------------------------------------------- -- We compute here the degree of the MacRae's invariant = degree of implicit equation degreeImplicitEq = method(); degreeImplicitEq(List, List) := (polynomialList, nu) -> { F := matrix {polynomialList}; Z0 := (ring F)^1; Z1 := kernel koszul(1,F); Z2 := kernel koszul(2,F); Z3 := kernel koszul(3,F); d := degree (polynomialList_0); return (hilbertFunction(nu+d,Z1)-2*hilbertFunction(nu+2*d,Z2)+3*hilbertFunction(nu+3*d,Z3)) } --------------------------------------------------------------------------- --------------------------------------------------------------------------- representationMatrix = method(); representationMatrix(List, List) := (polynomialList, nu) -> { l := length polynomialList -1; S := ring polynomialList_0; use S; F := matrix {polynomialList}; R := S[X_0..X_l]; G := sub(F,R); -- matrix presentation of Z1_nu in K1 Z1nu := super basis(nu+(degree (polynomialList_0)),kernel koszul(1,F)); Tnu := matrix{{X_0..X_l}}*substitute(Z1nu,R); --first map Z1_nu -> S lll=sub(matrix {gens S},R); gensSinR := apply(gens S, i -> sub(i,R)); (m,M)=coefficients(Tnu,Variables=>gensSinR,Monomials=>substitute(basis(nu,S),R)); T := QQ[X_0..X_l]; -- p has to be changed depending on l=3 (there should be l+1 zeros in the declaration of p) listofXand0 := toList(X_0..X_l); for i from 0 to l do { listofXand0 = append(listofXand0, 0) }; p := map(T,R,listofXand0); return (p(M)) } --------------------------------------------------------------------------- --------------------------------------------------------------------------- implicitEq = method(Options => { Strategy => null }); implicitEq(List, List) := (Matrix) => options -> (polynomialList, nu) -> { N := maxMinor(representationMatrix(polynomialList, nu), options); -- A multiple of the implicit equation (taking several minors we erase extra factors) return det(N) }; --------------------------------------------------------- --------------------------------------------------------- --------------------------------------------------------- beginDocumentation() --------------------------------------------------------- --------------------------------------------------------- -- Simple Doc information --------------------------------------------------------- --------------------------------------------------------- --******************************************************* -- DOCUMENTATION FOR PACKAGE --******************************************************* doc /// Key BigradedImplicit Headline A package for computing implicit equations of bigraded rational surfaces by means of approximation complexes Description Text @EM "BigradedImplicit"@ is a package for computing implicit equations of bigraded rational surfaces by means of approximation complexes. We provide a methd to compute a matrix representation (see @TO representationMatrix @) and the implicit equation (see @TO implicitEq @ ) by means of the method developed in [Bot11]. As it is probably the most interesting case from a practical point of view, we restrict our computations to parametrizations of bigraded surfaces. This implementation allows to compute small examples for the better understanding of the theory developed in [Bot11]. [Bot11] "The implicit equation of a multigraded hypersurface. arXiv:1007.3437v3 . Journal of Algebra. Vol 348, Issue 1 (2011), 381-401" {\bf Note:} We require all the polynomial ring to be multigraded, i.e. the ring should be of the form QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}], where {a+b,a,b} stands for the multidegree, where the first coordinate is the sum of the degrees, and the seecond and third coordinate are the degreen in each group of variables (s,u) and (t,v) respectively. /// document { Key => {isGoodDegree, (isGoodDegree, List, List)}, Headline => "verifies if the Z-complex is acyclic in the given degree", Usage => " isGoodDeg = isGoodDegree(polynomialList, nu)", Inputs => { "polynomialList" => List => {"polynomialList 'f={f0,...,fn}' defining the rational map"}, "nu" => List => {"the multidegree where ti take the homogeneous strand of the map f (i.e.: R_d)"} }, Outputs => { "bool" => Boolean => {"A boolean 'True' or 'False'"} }, PARA {" The list 'nu' needs to be a 'good degree' for the first parameter '{f0,...,fn}' that can be verified by doing isGoodDegree(polinomialList,nu) "}, PARA { "isGoodDegree({f0,...,fn}, nu) verifies if the approximation complex Z associated to the polynomials given is acyclic in degree 'nu' "}, PARA { "Given a list of polynomials {f_0,...,f_n}, the approximation complex of cycles is multigraded. This funtion verifies if the strand of degree 'nu' is acyclic."}, PARA { "Precisely, it computes the Euler characteristic of the nu-strand of the Z-complex, by computing the alternate sum of (-1)^i * hilbertFunction(nu+i*d,Z_i) "}, EXAMPLE { "A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}];", "f0=s*t*u*v+t*u^2*v; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v;", "isGoodDegree ({f0,f1,f2,f3},{2,1,1})" } } document { Key => {degreeImplicitEq, (degreeImplicitEq, List, List)}, Headline => "computes the degree of det((Z)_nu)", Usage => " theDeg = degreeImplicitEq(polynomialList, nu)", Inputs => { "polynomialList" => List => {"polynomialList 'f={f0,...,fn}' defining the rational map"}, "nu" => List => {"the multidegree where ti take the homogeneous strand of the map f (i.e.: R_d)"} }, Outputs => { "theDeg" => ZZ => {"the degree of the implicit equation"} }, PARA {" The list 'nu' needs to be a 'good degree' for the first parameter '{f0,...,fn}' that can be verified by doing isGoodDegree(polinomialList,nu) "}, PARA { "degreeImplicitEq({f0,...,fn}, nu) computes the degree of the gcd of the maximal minors of the matrix representation of {f0,...,fn} in degree 'nu', equivalently, the degree of determinant of the Z-complex associated to {f0,...,fn} in degree 'nu'. This is:"}, PARA { "degreeImplicitEq({f0,...,fn}, nu)=deg(representationMatrix({f0,...,fn},nu)=deg(det((Z(f0,...,fn))_nu))"}, PARA { "Given a list of polynomials {f_0,...,f_n}, it computes the degree of the gcd of the maximal minors of the right-most map of the strand of degree 'nu'."}, EXAMPLE { "A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}];", "f0=s*t*u*v+t*u^2*v; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v;", "degreeImplicitEq ({f0,f1,f2,f3},{2,1,1})" }, SeeAlso => { representationMatrix } } document { Key => {representationMatrix, (representationMatrix, List, List)}, Headline => "computes the right-most map of the Z-complex in degree nu)", Usage => " repMatrix = representationMatrix(polynomialList, nu)", Inputs => { "polynomialList" => List => {"polynomialList 'f={f0,...,fn}' defining the rational map"}, "nu" => List => {"the multidegree where ti take the homogeneous strand of the map f (i.e.: R_d)"} }, Outputs => { "repMatrix" => Matrix => {"a matrix of linear functions in QQ[X_0,...,X_n] defining a representation Matrix of f"} }, PARA {" The list 'nu' needs to be a 'good degree' for the first parameter '{f0,...,fn}' that can be verified by doing isGoodDegree(polinomialList,nu) "}, PARA { "representationMatrix({f0,...,fn},nu) computes the right-most map of the Z-complex in degree nu. Its determinant vanishes on the implicit equation of the image of the map given by (f0:...:fn) in P^n"}, EXAMPLE {"A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}];", "f0=s*t*u*v+t*u^2*v; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v;", "degreeImplicitEq ({f0,f1,f2,f3},{2,1,1})" }, SeeAlso => { implicitEq } } document { Key => {implicitEq, (implicitEq, List, List), [implicitEq, Strategy]}, Headline => "computes the gcd of the right-most map of the Z-complex in degree nu)", Usage => " iEq = implicitEq(polynomialList, nu)", Inputs => { "polynomialList" => List => {"polynomialList 'f={f0,...,fn}' defining the rational map"}, "nu" => List => {"the multidegree where ti take the homogeneous strand of the map f (i.e.: R_d)"} }, Outputs => { "iEq" => Matrix => {"An homogeneous polynomial in QQ[X_0,...,X_n] defining the implicit equation of f"} }, PARA {" The list 'nu' needs to be a 'good degree' for the first parameter '{f0,...,fn}' that can be verified by doing isGoodDegree(polinomialList,nu) "}, PARA { "implicitEq({f0,...,fn},nu) computes the determinant of the maximal minors of representationMatrix({f0,...,fn},nu). Equivalently, it computes the determinant of the Z-complex in degree 'nu'. This is:"}, PARA { "implicitEq({f0,...,fn}, nu)=det((Z(f0,...,fn))_nu)"}, PARA { "We also have that: degree(implicitEq({f0,...,fn}, nu))=degreeImplicitEq({f0,...,fn}, nu)"}, EXAMPLE {"A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}];", "f0=s*t*u*v+t*u^2*v; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v;", "implicitEq({f0,f1,f2,f3},{2,1,1})" }, SeeAlso => { representationMatrix, degreeImplicitEq, Exact, Numeric} } ------------------------ \ documentation maxMinor / ------------------------ document { Key => {maxMinor, (maxMinor, Matrix), [maxMinor, Strategy]}, Headline => "Returns a maximal minor of the matrix of full rank.", Usage => " MM = maxMinor(Mat)", Inputs => { "Mat" => Matrix }, Outputs => { "MM" => Matrix }, PARA {}, " From a given m x n - Matrix of rank r, ", TO maxMinor, " returns an r x r full rank Matrix. This method uses twice the method ", TO maxCol," by transposing twice.", EXAMPLE { " M = matrix {{1,2,3},{1,2,3},{4,5,6},{4,5,6}}", " maxMinor M"}, PARA {}, "NOTE: because of the necessity of ", TO rank," the base field need to be QQ for doing generic evaluation. If not, one gets the message: expected an affine ring (consider Generic=>true to work over QQ).", EXAMPLE { " R=QQ[a..g]", " M = matrix {{a,a,b},{c,c,d},{e,e,f},{g,g,g}}", " maxMinor M"}, SeeAlso => {rank, Exact, Numeric} } ------------------------ \ Strategy for functions on Complex / ------------------------ document { Key => [implicitEq, Strategy], Headline => "computes the gcd of the right-most map of the Z-complex in degree nu)", Usage => " iEq = implicitEq(polynomialList, nu)", Inputs => { "polynomialList" => List => {"polynomialList 'f={f0,...,fn}' defining the rational map"}, "nu" => List => {"the multidegree where ti take the homogeneous strand of the map f (i.e.: R_d)"} }, Outputs => { "iEq" => Matrix => {"An homogeneous polynomial in QQ[X_0,...,X_n] defining the implicit equation of f"} }, Consequences => {{ "If ", TT "s", " is ", TO "Exact", ", then the", TO "rank", " algortihms is used computing minors; if ", TT "s", " is ", TO "Numeric", ", then numerical rank computation is used, this is, all coefficients are evaluated in the ground field before computing ranks."}}, PARA{}, TO "Exact", "is the default Strategy.", SeeAlso => {Exact, Numeric} } document { Key => [maxMinor, Strategy], Headline => "Returns a maximal minor of the matrix of full rank.", Usage => " MM = maxMinor(Mat)", Inputs => { "Mat" => Matrix }, Outputs => { "MM" => Matrix }, Consequences => {{ "If ", TT "s", " is ", TO "Exact", ", then the", TO "rank", " algortihms is used computing minors; if ", TT "s", " is ", TO "Numeric", ", then numerical rank computation is used, this is, all coefficients are evaluated in the ground field before computing ranks."}}, PARA{}, TO "Exact", "is the default Strategy.", } document { Key => {Numeric}, Headline => "Strategy for functions that uses rank computation.", Usage => " Strategy => Numeric", Consequences => { "If 'Numeric' Strategy is used, then numerical rank computation is used instead of ", TO "Exact", " rank computation. Precisely, all coefficients are evaluated in the ground field before computing ranks."}, } document { Key => {Exact}, Headline => "Strategy for functions that uses rank computation.", Usage => " Strategy => Exact", Consequences => { "If 'Exact' Strategy is used, then exact rank computation is used. This could be notably slower that ", TO "Numeric", " rank computation when many generic elements are manipulated."}, } end --******************************************************* -- EXAMPLES FOR PACKAGE --******************************************************* ----------------------------------------------- -- example 1 ----------------------------------------------- A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}]; f0=s*t*u*v+t*u^2*v; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v; lista = {f0,f1,f2,f3}; nuu = {2,1,1} load "/home/Nicolas/Mis Notas y mis charlas/JSAG-Multigraded/codigo/BigradedImplicit.m2" representationMatrix (lista,nuu) isGoodDegree (lista,nuu) degreeImplicitEq (lista,nuu) implicitEq (lista,nuu) ----------------------------------------------- -- example 2 ----------------------------------------------- R=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}]; f1=1*s^2*t^3+2*s*u*t^3+3*u^2*t^3+4*s^2*t^2*v+5*s*u*t^2*v+6*u^2*t^2*v+7*s^2*t*v^2+8*s*u*t*v^2+9*u^2*t*v^2+10*s^2*v^3+1*s*u*v^3+2*u^2*v^3; f2=2*s^2*t^3-3*s^2*t^2*v-s^2*t*v^2+s*u*t^2*v+3*s*u*t*v^2-3*u^2*t^2*v+2*u^2*t*v^2-u^2*v^3; f3=2*s^2*t^3-3*s^2*t^2*v-2*s*u*t^3+s^2*t*v^2+5*s*u*t^2*v-3*s*u*t*v^2-3*u^2*t^2*v+4*u^2*t*v^2-u^2*v^3; f4=3*s^2*t^2*v-2*s*u*t^3-s^2*t*v^2+s*u*t^2*v-3*s*u*t*v^2-u^2*t^2*v+4*u^2*t*v^2-u^2*v^3; lista = {f1,f2,f3,f4}; nuu = {6,1,5} nuu = {5,3,2} load "/home/Nicolas/Mis Notas y mis charlas/JSAG-Multigraded/codigo/BigradedImplicit.m2" representationMatrix (lista,nuu) isGoodDegree (lista,nuu) degreeImplicitEq (lista,nuu) implicitEq (lista,nuu) nuu={2,1,1} Eq= (factor implicitEq (lista,nuu))#1#0 S=R[X_0..X_3] use S; Eq=sub(Eq,S) sub(Eq,{X_0=>sub(f1,S),X_1=>sub(f2,S),X_2=>sub(f3,S),X_3=>sub(f4,S)}) ----------------------------------------------- -- example 3 -- The theoretical bound can be improved !!!¿?¿?¿? ----------------------------------------------- A=QQ[s,u,t,v,Degrees=>{{1,1,0},{1,1,0},{1,0,1},{1,0,1}}]; f1=s*u*v^2; f2=u*t^2*s+u^2*t^2; f3=s^2*t*v+s*t*u*v; f4=s*t*u*v+t*u^2*v; lista = {f1,f2,f3,f4}; nuu = {2,1,1} load "/home/Nicolas/Mis Notas y mis charlas/JSAG-Multigraded/codigo/BigradedImplicit.m2" representationMatrix (lista,nuu) isGoodDegree (lista,nuu) degreeImplicitEq (lista,nuu) implicitEq (lista,nuu)