# Quadrics vanishing on orbits in Pluecker coordinates

As in Example 6.9, we consider the following nets:

$\begin{split}\mathcal{L}_1=\begin{pmatrix}x&0&0&0\\0&x&0&0\\0&0&y&0\\0&0&0&z\end{pmatrix},\mathcal{L}_2=\begin{pmatrix}x&y&0&0\\y&z&0&0\\0&0&x&y\\0&0&y&z\end{pmatrix},\mathcal{L}_3=\begin{pmatrix}0&0&x&y\\0&-2x&-y&z\\x&-y&-2z&0\\y&z&0&0\end{pmatrix}\end{split}$

In Example 6.10, we claim that 85 quadrics (in dual Plücker coordinates) vanish on the orbit of $$\mathcal{L}_1$$, 189 quadrics vanish on the orbit of $$\mathcal{L}_2$$, which also vanish on the orbit of $$\mathcal{L}_3$$, and that 918 quadrics vanish on the orbit of $$\mathcal{L}_3$$. We also give examples of such quadrics. The goal of this code is to verify these claims.

Our first goal is to parametrize the congruence orbits of $$\mathcal{L}_1,\mathcal{L}_2,\mathcal{L}_3$$ in the Grasmannian.

[1]:

# We generate variables for the action of GL(4)
R = PolynomialRing(QQ,['x'+str(i)+str(j) for i in range(1,5) for j in range(1,5)])
R.inject_variables()

# In stead of nets of symmetric matrices, we consider nets of quadrics.
# We generate their variables.
S.<a,b,c,d> = PolynomialRing(R,4)

# The congruence action replaces a,b,c,d by general linear combinations of themselves.
# We generate these.
[A,B,C,D] = [sum([R('x'+str(i)+str(j))*[a,b,c,d][i-1] for i in range(1,5)]) for j in range(1,5)]

# Next, we generate a basis for the net of quadrics
forms = [A^2+B^2,C^2,D^2]          # L1 has basis a^2+b^2, c^2, d^2
#forms = [A^2+C^2,A*B+C*D,B^2+D^2]  # L2 has basis a^2+b^2, a*b+c*d, b^2+d^2
#forms = [A*C-B^2,A*D-B*C,B*D-C^2]  # L3 has basis a*c-b^2, a*d-b*c, b*d-c^2

# We generate the matrix seen in the example.
M = Matrix([form.coefficients() for form in forms])
# Note that SAGE orders the indices differently:
# 11, 12, 22, 13, 23, 33, 14, 24, 34, 44
# Also, off-diagonal matrix entries are half of the corresponding quadric entry.
tmp = Matrix(QQ,[[0 for t1 in range(10)] for t2 in range(10)])
tmp[0,0] = 1
tmp[1,1] = 1/2
tmp[3,2] = 1/2
tmp[6,3] = 1/2
tmp[2,4] = 1
tmp[4,5] = 1/2
tmp[7,6] = 1/2
tmp[5,7] = 1
tmp[8,8] = 1/2
tmp[9,9] = 1
M = M*tmp

# We generate all dual Pluecker indices and compute the parametrisation in dual Pluecker coordinates.
plueckerindex = [[i,j,k] for i in range(10) for j in range(i+1,10) for k in range(j+1,10)]
plueckercoords = [Matrix([[M[s,t] for t in index] for s in range(3)]).det() for index in plueckerindex]

Defining x11, x12, x13, x14, x21, x22, x23, x24, x31, x32, x33, x34, x41, x42, x43, x44


Next, we need to generate all degree-$$2$$ monomials in the Plücker coordinates that are linearly independent modulo the Plücker relations.

[2]:

# We first compute the pair of indices these monomials correspond to. There are 4950 such pairs.
standardpairs = [[s,t] for s in range(120) for t in range(s,120) if plueckerindex[s][0]<=plueckerindex[t][0] if plueckerindex[s][1]<=plueckerindex[t][1] if plueckerindex[s][2]<=plueckerindex[t][2]]
standardproducts = [plueckercoords[s]*plueckercoords[t] for [s,t] in standardpairs]


We now generated a list of $$4950$$ homogeneous degree-$$12$$ polynomials in $$16$$ variables. We want to know the linear dependencies between them.

As this problem is too big to tackle directly, we consider in stead, for each polynomial, whether the polynomial is a linear combination of the previous ones, and if so how.

[3]:

def Dothings(polys,exprs,i):
# input:  a list of polynomials,
#         a list of expressions of those polynomials as linear combinations of standardproducts
#         an index
# output: IF the i-th polynomial of standardproducts if a linear combination of the previous ones,
#           we return [polys,exprs,expr] where expr is the list of coefficients of the dependence.
#         ELSE we return [polys+[f],exprs+[expr],[]] where f is the reduced version of
#           the i-th polynomial of standardproducts and expr expresses f as as a linear combination
#           of standardproducts.
#           NOTE: We sort polys such that the leading monomials are a descending chain.
#                 This ensures that we can reduce f fast!

# We take f the i-th polynomial of standardproducts.
f = standardproducts[i]
# It equals the i-th polynomial of standardproducts.
expr = [0 for k in range(4950)]
expr[i] = 1

# We reduce f with elements of polys.
# We change expr accordingly.
for j in range(len(polys)):
g = polys[j]
if g.lm()==f.lm():
[f,expr] = [g.lc()*f-f.lc()*g,[g.lc()*expr[k]-f.lc()*exprs[j][k] for k in range(4950)]]
tmp = gcd(expr)
f = f/tmp
expr = [expr[k]/tmp for k in range(4950)]
if f==0:
return [polys,exprs,[expr]]
elif (f+g).lm()==f.lm():
return [polys[:j]+[f]+polys[j:],exprs[:j]+[expr]+exprs[j:],[]]
return [polys+[f],exprs+[expr],[]]

[4]:

# We add the first polynomial to polys (since it is not zero).
polys = [standardproducts[0]]
# The first polynomial equals the first polynomial.
exprs = [[1]+[0 for k in range(1,4950)]]
# We record the linear dependencies.
res = []

# For each polynomial in the list, we check if it is a linear combination of the others.
# If so, add how to res. If not, reduce the polynomial, add it to polys (sorted) and
# keep track of what the reduced polynomial is in terms of the orginal ones.
for i in range(1,4950):
[polys,exprs,tmp] = Dothings(polys,exprs,i)
res = res + tmp
# This takes long... print how far we are and how many dependecies we found.
#print (i,len(res))

# Save results
f = open('L1res.txt', 'w') # change name to match L2, L3.
f.write(str(res))
f.close()

[5]:

def toString(eqn):
# input:  a depencency
# output: the dependency as a string
tmp = gcd(eqn)
eqn = [eqn[k]/tmp for k in range(4950)]
string = ''
for l in range(4950):
if eqn[l]!=0:
sign = (eqn[l]>0)
value = str(abs(eqn[l]))
if value=='1':
value=''
[a,b] = standardpairs[l]
[i,j,k] = plueckerindex[a]
[i2,j2,k2] = plueckerindex[b]
if string=='':
if sign==False:
string='-'
string = string + value +'p_{'+str(i)+str(j)+str(k)+'}p_{'+str(i2)+str(j2)+str(k2)+'}'
else:
if sign==False:
string = string + '-'
else:
string = string + '+'
string = string + value +'p_{'+str(i)+str(j)+str(k)+'}p_{'+str(i2)+str(j2)+str(k2)+'}'
return string


We verify the quadric vansihing on the orbit of $$\mathcal{L}_1$$.

[6]:

print toString(res[60])

p_{124}p_{789}-p_{125}p_{589}+p_{125}p_{679}-p_{126}p_{678}+p_{127}p_{569}-p_{128}p_{568}-p_{234}p_{579}+p_{234}p_{678}+p_{235}p_{479}+p_{235}p_{568}-p_{236}p_{478}-p_{236}p_{567}-p_{237}p_{459}+p_{238}p_{458}

[7]:

# Import the relations for L1
with open('./L1res.txt', 'r') as f:
res1 = eval(L[0])
X = Matrix(res1).transpose()
# NOTE: the matrix is in echelon form if we transpose.

# Import the relations for L2
with open('./L2res.txt', 'r') as f:
res2 = eval(L[0])
Y = Matrix(res2).transpose()
# NOTE: the matrix is in echelon form if we transpose.

# Import the relations for L3
with open('./L3res.txt', 'r') as f:
res3 = eval(L[0])
Z = Matrix(res3).transpose()
# NOTE: the matrix is in echelon form if we transpose.


We find the numbers in the paper.

[8]:

print len(res1)
print len(res2)
print len(res3)

85
189
918

[9]:

# Turn a degree-2 monomial into our format
def stdpair(i,j,k,i2,j2,k2):
res = [0 for tmp in range(4950)]
a = plueckerindex.index([i,j,k])
b = plueckerindex.index([i2,j2,k2])
res[standardpairs.index([a,b])] = 1
return Matrix(res).transpose()


We verify the quadric vansihing on the orbit of $$\mathcal{L}_2$$.

[10]:

eqn = stdpair(0,1,2,0,1,8)-stdpair(0,1,2,0,2,6)-stdpair(0,1,2,0,3,5)-2*stdpair(0,1,2,1,2,3)-stdpair(0,1,3,0,1,7)+2*stdpair(0,1,3,0,2,5)-stdpair(0,2,3,0,2,4)
# Check that the quadric is a linear combination of the ones we found.
block_matrix([[Y,eqn]]).rank() == Y.rank()

[10]:

True


We verify the quadric vansihing on the orbit of $$\mathcal{L}_3$$.

[11]:

# Generate the quadric on L3 from the paper.
eqn = 8*stdpair(0,1,2,4,5,7)-4*stdpair(0,4,5,0,5,7)+stdpair(0,4,7,0,4,7)
# Check that the quadric is a linear combination of the ones we found.
block_matrix([[Z,eqn]]).rank() == Z.rank()

[11]:

True


Finally, we check that the quadrics that vanish on the orbit of $$\mathcal{L}_2$$ also vanish on the orbit of $$\mathcal{L}_3$$.

[12]:

# Check that the equations on L2 are among the equation on L3.
block_matrix([[Y,Z]]).rank() == Z.rank()

[12]:

True