------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --INVARIANT DEFORMATION THEORY OF AFFINE SCHEMES WITH REDUCTIVE GROUP ACTION --Christian LEHN and Ronan TERPEREAU --22/02/2014 (preliminary version) --Situation 1 of our paper ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Initialisation ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Now we initialise the ideal I \subset RX of X_0, and some rings that we will use later on --Let us mention that the rings kk and RR will be useful in the routine "Reynolds" only --the weights of the t_i have to be computed by hand once we know a basis of the tangent space T_{X_0}\HH g=9; --the dimension of the group G(=GL3 here) RX=QQ[x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33]; --RX=k[W] x={x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33}; --the variables of RX I=ideal{y11*x11+y12*x21+y13*x31,y11*x12+y12*x22+y13*x32,y11*x13+y12*x23+y13*x33,y21*x11+y22*x21+y23*x31,y21*x12+y22*x22+y23*x32,y21*x13+y22*x23+y23*x33,y31*x12+y32*x22+y33*x32,y31*x13+y32*x23+y33*x33,x11*y33,x21*y33,x11*y32,x31*y33+x11*y31,x31*y33+x21*y32,x21*y31,x21*y32+x11*y31,x31*y32,x31*y31,- x12*x21*y23 + x11*x22*y23,x12*x21*y21 - x11*x22*y21 + x22*x31*y23 - x21*x32*y23,x12*x21*y22 - x11*x22*y22 - x12*x31*y23 + x11*x32*y23,- x22*x31*y21 + x21*x32*y21,x12*x31*y21 - x11*x32*y21 - x22*x31*y22 + x21*x32*y22,x12*x31*y22 - x11*x32*y22, x12*y22*y33 - x12*y23*y32, -x12*y21*y33 + x12*y23*y31 + x22*y22*y33 - x22*y23*y32,- x22*y21*y33 + x22*y23*y31,x22*y21*y32 - x22*y22*y31 - x32*y21*y33 + x32*y23*y31,x12*y21*y32 - x12*y22*y31 + x32*y22*y33 - x32*y23*y32,x32*y21*y32 - x32*y22*y31}; --the ideal I with a minimal number of generators and such that the generators form bases of the underlying G-modules RT=QQ[t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,t11,t12,Degrees=>{2,2,2,1,2,1,1,1,2,2,1,1}]; kk=frac RT; RGm=kk[x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33,xx11,xx12,xx13,xx21,xx22,xx23,xx31,xx32,xx33,yy11,yy12,yy13,yy21,yy22,yy23,yy31,yy32,yy33,Degrees=>{1,1,1,1,1,1,1,1,1,-1,-1,-1,-1,-1,-1,-1,-1,-1,1,1,1,1,1,1,1,1,1,-1,-1,-1,-1,-1,-1,-1,-1,-1}]; RRrey=kk[x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33, MonomialOrder => Lex]; R=QQ[t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,t11,t12,Degrees=>{2,2,2,1,2,1,1,1,2,2,1,1}][x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33,MonomialOrder => Lex]; Rbis=QQ[t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,t11,t12,x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33,MonomialOrder => Lex]; use(R); t={t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,t11,t12}; --t1,...,t_12 is a basis of the cotangent space (T_{X_0}\HH)* and thus RT is the ring S of the article ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Some general routines ------------------------------------------------------------------------------------------------------------------------------------------------------------------ -------------------------------------------------------------------------------------------- --the next routine takes a polynomial P and a list Lfi of polynomials f1...fn --and it determines a n-uplet a1...an of rational numbers such that P=a1f1+...+anfn --However f1...fn are not linearly independant in general and thus the ai are non-unique --the output is the list Lai which contains the coefficients a1...an -------------------------------------------------------------------------------------------- --This routine is an implementation of Algorithm 4.5.19 of Derksen and Klemper and here we follow their notation combilinear=method(); combilinear(RingElement,List):= (List) => (P,Lfi) -> ( Lfi=Lfi|{P}; --Lfi is the list which contains the fi, and we add the polynomial P at the end of the list nf:=length Lfi; Lgi:={}; Lai:={}; --Lgi is the list which will contain the gi, Lai is the list which will contain the ai MatBli:=mutableMatrix(QQ,nf,nf); --MatBli is the matrix which contains the bli for l from 0 to nf-1 do ( --we are going to fill in the list Lgi, if P is linear combination of the fi then we will get Lgi_(nf-1)=0 gl=Lfi_l; MatBli_(l,l)=1; p:=0; while p (M) -> ( zeilen=numRows M; spalten=numColumns M; for i from 0 to zeilen-1 do ( for j from 0 to spalten-1 do ( if M_(i,j)!=0 then return (leadTerm M_(i,j),i,j); ););); leitMonom=method(); --return the leading monomial of a matrix M leitMonom(Matrix):= (RingElement) => (M) -> ( lt=leitTerm M; return (leadMonomial(lt_0),lt_1,lt_2); ); leitKoeffizient=method(); --return the leading coefficient of a matrix M leitKoeffizient(Matrix):= (RingElement) => (M) -> ( lt=leitTerm M; return (leadCoefficient(lt_0),lt_1,lt_2); ); -- given a list li of polynomials and a matrix M the next routine checks if there is an index such that LeadingMonomial(li_i)==LeadingMonomial(M) -- If this equality holds for some i, then the routine returns (true, the smallest such i), else it returns (false, 1+Card(li)) existequalLM=method(); existequalLM(List,Matrix):= (Boolean,ZZ) => (li,M) -> ( existequal=false; indexi=#li+1; if M==0 then return (existequal,indexi); for i from 0 to #li-1 do ( if (leitMonom(li_i)==leitMonom(M)) then (existequal=true;indexi = i;return(existequal,indexi););); return (existequal,indexi);); leitTerm2=method(); --return the leading term of a vector v and the corresponding coordinate i leitTerm2(Vector):=(RingElement,ZZ)=>(v)->( w=leadTerm v; n=#entries v -1; for i from 0 to n do (if w_i!=0 then return(w_i,i)); return (0,n+1);); leitKoeffizient2=method(); --return the leading coefficient of a vector v and the corresponding coordinate i leitKoeffizient2(Vector):=(RingElement,ZZ)=>(v)->( w=leadTerm v; n=#entries v -1; for i from 0 to n do (if w_i!=0 then return(leadCoefficient(w_i),i)); return (0,n+1);); -- the input is an integer indi and a list si of vectors with leading coefficient equals 1 -- if there exists a vector si_i in the list si such that LM(si_i)=LM(si_indi) with i (si,indi) ->( if indi<=0 then return(false,#si+1); if si=={} then return(false,#si+1); if #si<=indi then return(false,#si+1); if matrix(si_indi)==0 then return(false,#si+1); for i from 0 to indi-1 do (if leitTerm2(si_i)==leitTerm2(si_indi) then return(true,i);); return(false,indi+1);); -- Given a list li of non-zero vectors, the next routine compute a list "simple" which contains -- vectors with leading coefficient=1 et whose the leading terms are pairwise differents -- In particular the vectors of the list simple form a basis of the vector space generated by the vectors of li -- The output is the list "base" which is a sublist of li whose vectors form a basis of the vector space generated by the vectors of li lindeplist=method(); lindeplist(List):=(List)=>(li)->( base:={}; simple:=new MutableList from li; nn=#simple-1; for j from 0 to nn do (simple#j=substitute(1/((leitKoeffizient2(simple#j))_0),QQ) * simple#j); for i from 0 to nn do ( (existe,indexi)=existequalLM2(toList(simple),i); while existe do (simple#i=simple#i-simple#indexi; if (matrix(simple#i)!=0) then simple#i=substitute(1/((leitKoeffizient2(simple#i))_0),QQ) * simple#i; (existe,indexi)=existequalLM2(toList(simple),i);); if (matrix(simple#i)!=0) then base=base|{li_i}; ); return base;); ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Routines to make the Lie algebra Lie(G) act on RX=k[W], and on the matrices A_i and B_i --From now on we will denote by D1,D2...Dg the canonical basis of the Lie algebra g=Lie(G) ------------------------------------------------------------------------------------------------------------------------------------------------------------------ -- here we have the following correspondance between the Dp and the Eij: p=1,2,3,4... corresponds to (i,j)=(1,1),(1,2),(1,3),(2,1)... that is p=3(i-1)+j (actionLieAlgebra=method(); --Given a polynomial f in k[W] and an integer p between 1 and 3, we return the polynomial Dp.f actionLieAlgebra(RingElement,ZZ):= (RingElement) => (f,p) -> ( if f==0 then return f; if p==1 then return (x11*diff(x11,f)+x12*diff(x12,f)+x13*diff(x13,f) - (y11*diff(y11,f)+y21*diff(y21,f)+y31*diff(y31,f))); if p==2 then return (x21*diff(x11,f)+x22*diff(x12,f)+x23*diff(x13,f) - (y11*diff(y12,f)+y21*diff(y22,f)+y31*diff(y32,f))); if p==3 then return (x31*diff(x11,f)+x32*diff(x12,f)+x33*diff(x13,f) - (y11*diff(y13,f)+y21*diff(y23,f)+y31*diff(y33,f))); if p==4 then return (x11*diff(x21,f)+x12*diff(x22,f)+x13*diff(x23,f) - (y12*diff(y11,f)+y22*diff(y21,f)+y32*diff(y31,f))); if p==5 then return (x21*diff(x21,f)+x22*diff(x22,f)+x23*diff(x23,f) - (y12*diff(y12,f)+y22*diff(y22,f)+y32*diff(y32,f))); if p==6 then return (x31*diff(x21,f)+x32*diff(x22,f)+x33*diff(x23,f) - (y12*diff(y13,f)+y22*diff(y23,f)+y32*diff(y33,f))); if p==7 then return (x11*diff(x31,f)+x12*diff(x32,f)+x13*diff(x33,f) - (y13*diff(y11,f)+y23*diff(y21,f)+y33*diff(y31,f))); if p==8 then return (x21*diff(x31,f)+x22*diff(x32,f)+x23*diff(x33,f) - (y13*diff(y12,f)+y23*diff(y22,f)+y33*diff(y32,f))); if p==9 then return (x31*diff(x31,f)+x32*diff(x32,f)+x33*diff(x33,f) - (y13*diff(y13,f)+y23*diff(y23,f)+y33*diff(y33,f))); ); ) --The next routine is useful only for the forthcoming routine "actionN2" in which we have to consider relations x1 \otimes f1+x2 \otimes f2+... --that we represent as an element xx1f1+xx2f2+... of the ring Rb --and thus to make Lie(G) act on these relations, we will use this routine (actionLieAlgebraBis=method(); actionLieAlgebraBis(RingElement,ZZ):= (RingElement) => (f,p) -> ( if f==0 then return f; if p==1 then return (xx11*diff(xx11,f)+xx12*diff(xx12,f)+xx13*diff(xx13,f) - (yy11*diff(yy11,f)+yy21*diff(yy21,f)+yy31*diff(yy31,f))); if p==2 then return (xx21*diff(xx11,f)+xx22*diff(xx12,f)+xx23*diff(xx13,f) - (yy11*diff(yy12,f)+yy21*diff(yy22,f)+yy31*diff(yy32,f))); if p==3 then return (xx31*diff(xx11,f)+xx32*diff(xx12,f)+xx33*diff(xx13,f) - (yy11*diff(yy13,f)+yy21*diff(yy23,f)+yy31*diff(yy33,f))); if p==4 then return (xx11*diff(xx21,f)+xx12*diff(xx22,f)+xx13*diff(xx23,f) - (yy12*diff(yy11,f)+yy22*diff(yy21,f)+yy32*diff(yy31,f))); if p==5 then return (xx21*diff(xx21,f)+xx22*diff(xx22,f)+xx23*diff(xx23,f) - (yy12*diff(yy12,f)+yy22*diff(yy22,f)+yy32*diff(yy32,f))); if p==6 then return (xx31*diff(xx21,f)+xx32*diff(xx22,f)+xx33*diff(xx23,f) - (yy12*diff(yy13,f)+yy22*diff(yy23,f)+yy32*diff(yy33,f))); if p==7 then return (xx11*diff(xx31,f)+xx12*diff(xx32,f)+xx13*diff(xx33,f) - (yy13*diff(yy11,f)+yy23*diff(yy21,f)+yy33*diff(yy31,f))); if p==8 then return (xx21*diff(xx31,f)+xx22*diff(xx32,f)+xx23*diff(xx33,f) - (yy13*diff(yy12,f)+yy23*diff(yy22,f)+yy33*diff(yy32,f))); if p==9 then return (xx31*diff(xx31,f)+xx32*diff(xx32,f)+xx33*diff(xx33,f) - (yy13*diff(yy13,f)+yy23*diff(yy23,f)+yy33*diff(yy33,f))); ); ) ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of A_0, initialisation of n1 ------------------------------------------------------------------------------------------------------------------------------------------------------------------ A={}; --A is a list which will contain the matrices Ai if (length A)==0 then A=A|{substitute(gens I,R)}; --we construct A0 (if it does not already exist) n1:=numColumns A_0; --n1 is the number of generators of the ideal I --Given A0, the next routine determines a liste LLN1 of g matrices of size n1*n1 with rational coefficients such that Dp.gen_i=sum_j (LLN1_p)_(j,i) gen_j --where gen_i is the i-th generator of the ideal I, that is the i-th coefficient of A0 --In other words, we have Dp.A0=A0*LLN1_p for p=1...g --This routine encodes the action of Lie(G) on N1 and thus it will be useful to make Dp act on the matrices Ai (actionN1=method(); actionN1(Matrix):= (List) => (A0) -> ( use(RX); LLN1={}; LA0={}; LM={}; for i from 0 to n1-1 do LA0=LA0|{substitute(A0_(0,i),RX)}; --LA0 is a list which contains the coefficients of A0 Z=mutableMatrix(QQ,n1,n1); for p from 1 to g do (for i from 0 to n1-1 do ( LM=combilinear(actionLieAlgebra(LA0_i,p),LA0); --here we use the routine "combilinear" to express each Dp.A0_i as a linear combination of the A0_j for j from 0 to n1-1 do Z_(j,i)=LM_j; ); LLN1=LLN1|{matrix(Z)};); return LLN1; );) LN1:=actionN1(A_0); -- around 7 seconds ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of B_0, initialisation of n2 ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --The next routine takes a matrix M of size n1*b2 which contains relations between the generators of the ideal I --and for each columns Cl of M it constructs a matrix Z whose columns are D1.Cl, D2.Cl,...,Dg.Cl --then it makes M=M|Z for each column Cl and at the end it returns "mingens image M" --Hence if we repeat this procedure enough times, we can replace the matrix M by a matrix Mbis whose every column CMi verifies Dp.CMi is a linear combination of the CMj --Rk: The image of the matrix M generates the vector space N'_2 in the article, while Mbis generates the G-module N2 Gstabilization=method(); Gstabilization(Matrix):= (Matrix) => (M) -> ( b2=numColumns M; for c from 0 to b2-1 do ( Z=mutableMatrix(RX,n1,g); for p from 1 to g do ( for l from 0 to n1-1 do Z_(l,p-1)=actionLieAlgebra(M_(l,c),p)+sum(n1,k->M_(k,c)*substitute((LN1_(p-1))_(l,k),RX)); ); M=M|matrix(Z); ); return mingens image M ); use(RX); G2 = gb(I, Syzygies=>true); C=mingens image syz G2; --relations generated by the kernel of I3, I4 and I5 C1=submatrix(C,,{0..29}); --kernel of I3, already G-stable C2=submatrix(C,,{30..175}); --kernel of I4, C3=submatrix(C,,{176..187}); --kernel of I5 --now we have relations between the genrators of I but a priori the vector subspace of R \otimes N1 generated by this set of relations is not G-stable --so we have to make it G-stable by using the method "Gstabilization" above --In our situation, it is enough to apply it only once C2=Gstabilization(C2); --around 6 minutes C3=Gstabilization(C3); --around 30 seconds c1=numColumns C1; c2=numColumns C2; c3=numColumns C3; use(R); B={}; --the list B will contain matrices Bi if (length B)==0 then B=B|{substitute((C1|C2|C3),R)}; --we define the matrix B0 (if it has not been already defined) n2:=numColumns B_0; --Given B0, the next routine determines a liste LLN2 of g matrices of size n2*n2 with rational coefficients such that Dp.r_j=sum_i (LLN2_p)_(i,j) r_i --where r_i is the i-th relation between the generators of the ideal I, that is r_i corresponds to the i-th column of B0 --In other words, we have Dp.B0=B0*LLN2_p for p=1...g --This routine encodes the action of Lie(G) on N2 and thus it will be useful to make Dp act on the matrices Bi actionN2=method(); actionN2(Matrix):= (List) => (B0) -> ( b2=numColumns B0; LLN2={}; Lrel={}; LM={}; Rb=QQ[x11,x12,x13,x21,x22,x23,x31,x32,x33,y11,y12,y13,y21,y22,y23,y31,y32,y33,xx11,xx12,xx13,xx21,xx22,xx23,xx31,xx32,xx33,yy11,yy12,yy13,yy21,yy22,yy23,yy31,yy32,yy33]; use(Rb); B0b=substitute(B0,Rb); A0b=matrix{{yy11*xx11+yy12*xx21+yy13*xx31,yy11*xx12+yy12*xx22+yy13*xx32,yy11*xx13+yy12*xx23+yy13*xx33,yy21*xx11+yy22*xx21+yy23*xx31,yy21*xx12+yy22*xx22+yy23*xx32,yy21*xx13+yy22*xx23+yy23*xx33,yy31*xx12+yy32*xx22+yy33*xx32,yy31*xx13+yy32*xx23+yy33*xx33,xx11*yy33,xx21*yy33,xx11*yy32,xx31*yy33+xx11*yy31,xx31*yy33+xx21*yy32,xx21*yy31,xx21*yy32+xx11*yy31,xx31*yy32,xx31*yy31,- xx12*xx21*yy23 + xx11*xx22*yy23,xx12*xx21*yy21 - xx11*xx22*yy21 + xx22*xx31*yy23 - xx21*xx32*yy23,xx12*xx21*yy22 - xx11*xx22*yy22 - xx12*xx31*yy23 + xx11*xx32*yy23,- xx22*xx31*yy21 + xx21*xx32*yy21,xx12*xx31*yy21 - xx11*xx32*yy21 - xx22*xx31*yy22 + xx21*xx32*yy22,xx12*xx31*yy22 - xx11*xx32*yy22, xx12*yy22*yy33 - xx12*yy23*yy32, -xx12*yy21*yy33 + xx12*yy23*yy31 + xx22*yy22*yy33 - xx22*yy23*yy32,- xx22*yy21*yy33 + xx22*yy23*yy31,xx22*yy21*yy32 - xx22*yy22*yy31 - xx32*yy21*yy33 + xx32*yy23*yy31,xx12*yy21*yy32 - xx12*yy22*yy31 + xx32*yy22*yy33 - xx32*yy23*yy32,xx32*yy21*yy32 - xx32*yy22*yy31}}; Rel=mutableMatrix(Rb,1,b2); for i from 0 to b2-1 do ( for j from 0 to n1-1 do Rel_(0,i)=Rel_(0,i)+A0b_(0,j)*B0b_(j,i);); Rell=matrix(Rel); --the matrix Rell contains the relations r_1,...,r_n2 for i from 0 to b2-1 do Lrel=Lrel|{Rell_(0,i)}; --the lements of the list Lrel contains the relations r_1,...,r_n2 Z=mutableMatrix(QQ,b2,b2); for p from 1 to g do (for i from 0 to b2-1 do ( use(Rb); LM=combilinear(actionLieAlgebra(Lrel_i,p)+actionLieAlgebraBis(Lrel_i,p),Lrel); for j from 0 to b2-1 do Z_(j,i)=LM_j; ); LLN2=LLN2|{matrix(Z)}; ); return LLN2; ); --We have three types of relations (depending on the degree of the coefficients) --and the Lie algebra action preserves each of these type, so we can compute LN2 by blocks to compute faster LN21:=actionN2(C1); --around 7 seconds LN22:=actionN2(C2); --around 70 minutes LN23:=actionN2(C3); LN2={}; --now it remains to stick the blocks together to construct the matrices of the list LN2 for p from 1 to g do ( Z=mutableMatrix(QQ,n2,n2); for i from 0 to c1-1 do (for j from 0 to c1-1 do Z_(i,j)=(LN21_(p-1))_(i,j); ); for i from 0 to c2-1 do (for j from 0 to c2-1 do Z_(c1+i,c1+j)=(LN22_(p-1))_(i,j); ); for i from 0 to c3-1 do (for j from 0 to c3-1 do Z_(c1+c2+i,c1+c2+j)=(LN23_(p-1))_(i,j); ); LN2=LN2|{matrix Z} ); ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of the action of the Casimir (of Lie(G)) on the matrices Ai and Bi ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --If nb==1, then M is a matrix Ai, that is a matrix Rn1->R --If nb==2, then M is a matrix Bi, that is a matrix Rn2->Rn1 --The next routine take a matrix M and computes the matrix Dp.M (where Lie(G) acts on the coefficients but also on the basis vectors) actionLieAlgebraSurMat=method(); actionLieAlgebraSurMat(Matrix,ZZ,ZZ):= (Matrix) => (M,p,nb) -> ( if nb==1 then ( ZZ1=mutableMatrix(R,1,n1); for i from 0 to n1-1 do ZZ1_(0,i)=actionLieAlgebra(M_(0,i),p); ZZ2=matrix(ZZ1)-M*LN1_(p-1); ); if nb==2 then ( ZZ1=mutableMatrix(R,n1,n2); for i from 0 to n1-1 do ( for j from 0 to n2-1 do ZZ1_(i,j)=actionLieAlgebra(M_(i,j),p); ); ZZ2=matrix(ZZ1)+LN1_(p-1)*M-M*LN2_(p-1); ); return ZZ2; ); actionLieAlgebraSurMat = memoize actionLieAlgebraSurMat; --optional, cf doc MACAULAY2 --here we consider the Casimir operator for the SL3 action --recall that the Casimir is given (up to a scalar) by c= sum_ij Dij*D_ji - 1/3 * sum_ij Dii*Djj (where we denote by Dij the Dp such that p=3(i-1)+j) casimir=method(); casimir(Matrix,ZZ):= (Matrix) => (M,nb) -> ( if nb==1 then (ln=1; col=n1;); if nb==2 then (ln=n1; col=n2;); Z2=matrix(mutableMatrix(R,ln,col)); Z4=matrix(mutableMatrix(R,ln,col)); Z5=matrix(mutableMatrix(R,ln,col)); for i from 1 to 3 do ( for j from 1 to 3 do ( Z1=actionLieAlgebraSurMat(M,3*(j-1)+i,nb); Z1=actionLieAlgebraSurMat(Z1,3*(i-1)+j,nb); Z2=Z2+Z1; ); ); for j from 1 to 3 do ( Z3=actionLieAlgebraSurMat(M,4*j-3,nb); Z4=Z4+Z3; ); for i from 1 to 3 do ( Z3=actionLieAlgebraSurMat(Z4,4*i-3,nb); Z5=Z5+Z3; ); Z=Z2-(1/3)*Z5; return Z; ); casimir = memoize casimir; --optional, cf doc MACAULAY2 ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of the Reynolds operator applied on the matrices Ai and Bi --This routine is an implementation of Algorithm 4.5.19 of Derksen and Klemper and here we follow their notation ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --the list deglist1 and deglist2 will be useful in the next routines when we will compute the Reynolds operator use(RGm); deglist1={}; for k from 0 to n1-1 do deglist1=deglist1|degree (substitute(A_0,RGm))_(0,k); deglist2={}; B0b=substitute(B_0,RGm); A0b=matrix{{yy11*xx11+yy12*xx21+yy13*xx31,yy11*xx12+yy12*xx22+yy13*xx32,yy11*xx13+yy12*xx23+yy13*xx33,yy21*xx11+yy22*xx21+yy23*xx31,yy21*xx12+yy22*xx22+yy23*xx32,yy21*xx13+yy22*xx23+yy23*xx33,yy31*xx12+yy32*xx22+yy33*xx32,yy31*xx13+yy32*xx23+yy33*xx33,xx11*yy33,xx21*yy33,xx11*yy32,xx31*yy33+xx11*yy31,xx31*yy33+xx21*yy32,xx21*yy31,xx21*yy32+xx11*yy31,xx31*yy32,xx31*yy31,- xx12*xx21*yy23 + xx11*xx22*yy23,xx12*xx21*yy21 - xx11*xx22*yy21 + xx22*xx31*yy23 - xx21*xx32*yy23,xx12*xx21*yy22 - xx11*xx22*yy22 - xx12*xx31*yy23 + xx11*xx32*yy23,- xx22*xx31*yy21 + xx21*xx32*yy21,xx12*xx31*yy21 - xx11*xx32*yy21 - xx22*xx31*yy22 + xx21*xx32*yy22,xx12*xx31*yy22 - xx11*xx32*yy22, xx12*yy22*yy33 - xx12*yy23*yy32, -xx12*yy21*yy33 + xx12*yy23*yy31 + xx22*yy22*yy33 - xx22*yy23*yy32,- xx22*yy21*yy33 + xx22*yy23*yy31,xx22*yy21*yy32 - xx22*yy22*yy31 - xx32*yy21*yy33 + xx32*yy23*yy31,xx12*yy21*yy32 - xx12*yy22*yy31 + xx32*yy22*yy33 - xx32*yy23*yy32,xx32*yy21*yy32 - xx32*yy22*yy31}}; Rel=mutableMatrix(RGm,1,n2); for i from 0 to n2-1 do (for j from 0 to n1-1 do Rel_(0,i)=Rel_(0,i)+A0b_(0,j)*B0b_(j,i);); for k from 0 to n2-1 do deglist2=deglist2|degree Rel_(0,k); reynoldsGm=method(); --to compute Reynolds(M) for the action of the multiplicative group Gm reynoldsGm(Matrix,ZZ):= (Matrix) => (M,nb) -> ( use(RGm); M=substitute(M,RGm); if nb==1 then ( ReyM=mutableMatrix(RGm,1,n1); for i from 0 to n1-1 do ( if M_(0,i)!=0 then ( (M1,C1)=coefficients M_(0,i); for k from 0 to numColumns M1-1 do ( if (degree M1_(0,k))_0==deglist1_i then (ReyM_(0,i)=ReyM_(0,i)+M1_(0,k)*C1_(k,0)); ); ); ); ); if nb==2 then ( ReyM=mutableMatrix(RGm,n1,n2); for i from 0 to n1-1 do ( for j from 0 to n2-1 do ( if M_(i,j)!=0 then ( (M2,C2)=coefficients M_(i,j); for k from 0 to numColumns M2-1 do ( if (degree M2_(0,k))_0==(deglist2_j-deglist1_i) then (ReyM_(i,j)=ReyM_(i,j)+M2_(0,k)*C2_(k,0)); ); ); ); ); ); ReyM=matrix(ReyM); use(R); return substitute(ReyM,R); ); reynolds=method(); --to compute Reynolds(MR) for the action of the group G reynolds(Matrix,ZZ):= (Matrix) => (MR,nb) -> ( ff={MR}; use(RRrey); M=substitute(MR,RRrey); reyf=matrix mutableMatrix(RRrey,numRows(MR),numColumns(MR)); gg:={}; ell=0; bb:={}; a:=0; gell:=0; gellnonzero=(M!=0); while gellnonzero do ( use RRrey; bell={}; for i from 0 to ell-1 do bell=bell|{0}; bell=bell|{1}; bell=new MutableList from bell; gell=substitute(ff_ell,RRrey); if gell!=0 then a=leitKoeffizient(gell); (existe,indexi) = existequalLM(gg,gell); while existe do ( gell=gell-a_0*(gg_indexi); for j from 0 to indexi do bell#j=bell#j-a_0*bb_indexi_j; if gell!=0 then a=leitKoeffizient(gell); (existe,indexi) = existequalLM(gg,gell);); if gell==0 then ( gellnonzero=false; if bell#0==0 then (for k from 1 to ell do reyf=reyf+bell#k*substitute(ff_(k-1),RRrey);); return(1/bell#1*reyf);); gell=1/a_0 * gell; bell=1/a_0 * toList bell; bb=bb|{bell}; gg=gg|{gell}; ell= ell+1; use R; ff=ff|{casimir(ff_(ell-1),nb)};); use(R); return(substitute(reyf,R));); reynolds=memoize reynolds; --optional, cf doc MACAULAY2 ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of the tangent space of the invariant Hilbert scheme at the point I ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Here the computation of the tangent space is partially by hand to save computation time use(RX); H=Hom(I,RX/I); --H is a RX/I module G=gens H; MM=matrix(G_0)|matrix(G_19)|matrix(G_20)|matrix(G_21)|matrix(G_22)|matrix(G_23)|(y23*matrix(G_13))|(y13*matrix(G_11))|(y13*matrix(G_13))|(x33*matrix(G_9))|(x33*matrix(G_15))|(x32*matrix(G_9)); --the columns of the matrix MM are the vectors on which we have to apply the Reynolds operator in order to obtain a basis of the tangent space use(R); MM=substitute(lift(MM,RX),R); A1=sum(numColumns MM,l->substitute(t_l,R)*substitute(transpose matrix(MM_l),R)); A1=matrix(mutableMatrix A1); --trick to kill the multidegree A1=substitute(reynolds(A1,1),R); if length A==1 then (A=A|{A1}); ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of the obstruction map ------------------------------------------------------------------------------------------------------------------------------------------------------------------ kfinal=method(); --this routines takes and ideal Kn and an integer n and returns the subideal obtained by keeping the generators of degree at most n kfinal(Ideal,ZZ):= (Ideal) => (Kn,n) -> ( k=numgens(Kn)-1; liste=for i from 0 to k list (if (degree(Kn_i) > {0,n}) then continue; Kn_i); Kfinaln=ideal liste ; return substitute(Kfinaln,R); ); --here we construct some matrices wich will be useful in the routine "obs" Z=mutableMatrix(R,n1*n2,n2); for i from 0 to n1-1 do for j from 0 to n2-1 do Z_(plus(times(i,n2),j),j)=(A_0)_(0,i); T1= transpose matrix(Z); T2=transpose B_0; T3= T2|T1; Q2=image T3; Q1=super Q2; Q=Q1/Q2; h=inducedMap(Q,Q1); --here we have to be careful with the direction of the arrow: Q<--Q1, around 45 seconds --the next routine works with the lists of matrices A, B (defined above) and P (defined below) obs=method(); obs(ZZ):= (Ideal,Boolean) => (n) -> ( -- the n means that we are going to compute An, Bn and the ideal Kn bool=false; (q,r)=quotientRemainder(h*transpose P_(n-2),h*T3); rc=cover r; (mm,cc) = coefficients rc; Knp=ideal cc; K=kfinal(Knp,n); K=ideal mingens K; --to simplify the presentation of the ideal K Rnp=R/Knp; use(Rnp); PPn=transpose(P_(n-2))**Rnp; Sn=T3**Rnp; (X1,X2)=quotientRemainder(PPn,Sn); --X1 give us the matrices An and Bn, X2=0 is the rest X3=X1^{0..n1-1}; x4=mutableMatrix(Rnp,n1,n2); for i from 0 to n1-1 do for j from 0 to n2-1 do x4_(i,j)=X1_(plus(plus(n1,times(n2,i)),j),0); X5=matrix(x4); --X5 will be the matrix Bn if (length A)==n then A=A|{substitute(reynolds(lift(transpose X3,R),1),R)}; if (length B)==n then B=B|{substitute(reynolds(lift(X5,R),2),R)}; if ((sum(n+1,l->A_l)*sum(n+1,l->B_l))**(R/K))==0 then bool=true; return (K,bool); ); ------------------------------------------------------------------------------------------------------------------------------------------------------------------ --Computation of the universal invariant deformation ------------------------------------------------------------------------------------------------------------------------------------------------------------------ B1=(-(A_1*B_0)//A_0); if (length B)==1 then B=B|{substitute(reynolds(B1,2),R)}; P={}; --will contain the sum of products Ai*Bj required to compute the obstruction at each step if (length P)==0 then P=P|{-A_1*B_1}; n=2; --n is the integer such that (A0+...+An,B0+...+Bn) gives the universal deformation bb=false; if P_0==0 then bb=true; while not bb do ( (K,bb)=obs(n); if length P==n-1 then (P=P|{P_(n-2)-A_n*B_0-A_0*B_n-sum(n,l->A_(l+1)*B_(n-l))}); n=n+1; ); AA=sum(n,l->A_l); --the ideal K is such that RT/K is the base space of the universal invariant deformation --the image of the matrix AA in k[W] \otimes RT/K is the ideal of the universal subscheme ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------------------------------------------------------------------------------