> restart: with(ListTools): with(LinearAlgebra): with(Groebner):

>

> #Input 11 matrices

> A:=Matrix([[-4, -4, 4, 0, 1, 4, -3, -4, 3, 4, 0, 2, 4, 3, -3, 1], [0, 0, 0, 1, 3, -2, 3, -1, -1, -1, -2, -3, -1, 2, 3, 2], [-1, 4, -2, 0, -4, 0, 0, -1, -1, -1, 0, 2, -1, 1, 2, 2], [-2, -2, 1, -3, -2, 0, -2, 3, 4, 2, -4, 4, 1, 4, 3, -2], [4, -4, 1, 3, 2, -3, -4, 0, -4, 0, 4, 2, -4, 0, -2, 0], [2, 2, 0, -1, 2, -4, -2, 0, 3, 3, 1, -1, 4, 1, -2, -4], [2, -1, 4, 0, 1, -3, -1, 3, 4, 0, -4, 0, 0, -1, 3, 4], [0, 4, -2, 3, 3, 2, -1, 0, -1, 1, 3, 2, 2, 1, 4, 3], [2, -2, -1, -3, -1, 2, 1, -4, 4, 3, 4, 4, -4, -1, -1, 3], [-4, 4, -1, -3, 2, 1, -3, 2, 0, 0, 4, 4, -1, 4, 1, -4], [1, 3, 2, 4, 1, 0, -4, 3, -2, -2, -2, 2, 0, -4, 4, -2]]):

>

> #Input 4 by 4 matrix Q

> Q:=Matrix([[x[1,1],x[1,2],x[1,3],x[1,4]],[x[2,1],x[2,2],x[2,3],x[2,4]],[x[3,1],x[3,2],x[3,3],x[3,4]],[x[4,1],x[4,2],x[4,3],x[4,4]]]):

> # Caculate the inner product <A_j,Q> (11 linear polynomials )

> VQ:=convert(Q,Vector):

> L:=[seq(simplify(Transpose(A[k,1..16]).VQ),k=1..11)]:

> # Caculate the minors of Q (16 homogeneous polynomials)

> S:=combinat[choose](4,3):

> mjk:=[seq(seq(Determinant(SubMatrix(Q,r,c)),r=S),c=S)]:

>

> Eq:=[op(L),op(mjk)]:

> # Caculate Grobner Basis to obtain f0

> Gbb:=Basis(Eq,plex(x[1,1],x[1,2],x[1,3],x[1,4],x[2,1],x[2,2],x[2,3],x[2,4],x[3,1],x[3,2],x[3,3],x[3,4],x[4,1],x[4,2],x[4,3],x[4,4])):

> f0:=subs({x[4,4]=1,x[4,3]=y},Gbb[1]):

> # Count the number of nonzero real roots of f_0

> s := sturmseq(f0,y):

> NumberZero:=sturm(s, y, -infinity, +infinity):

> printf(" The number of nonzero real roots of f_0 is %d \n", NumberZero):

>

> for j from 1 to 4 do

> for k from 1 to 4 do

> NewEq:=[op(Eq),x[j,k]-1,x[4,3],x[4,4]]:

> NewGb:=Basis(NewEq,plex(x[1,1],x[1,2],x[1,3],x[1,4],x[2,1],x[2,2],x[2,3],x[2,4],x[3,1],x[3,2],x[3,3],x[3,4],x[4,1],x[4,2],x[4,3],x[4,4])):

> if NewGb[1]=1 then

> printf(" If x[%d,%d]=1 then the polynomial system has no roots \n", j,k):

> end if

> end do:

> end do:

>

>

>

>

The number of nonzero real roots of f_0 is 0 If x[1,1]=1 then the polynomial system has no roots

If x[1,2]=1 then the polynomial system has no roots

If x[1,3]=1 then the polynomial system has no roots

If x[1,4]=1 then the polynomial system has no roots

If x[2,1]=1 then the polynomial system has no roots

If x[2,2]=1 then the polynomial system has no roots

If x[2,3]=1 then the polynomial system has no roots

If x[2,4]=1 then the polynomial system has no roots

If x[3,1]=1 then the polynomial system has no roots

If x[3,2]=1 then the polynomial system has no roots

If x[3,3]=1 then the polynomial system has no roots

If x[3,4]=1 then the polynomial system has no roots

If x[4,1]=1 then the polynomial system has no roots

If x[4,2]=1 then the polynomial system has no roots

If x[4,3]=1 then the polynomial system has no roots

If x[4,4]=1 then the polynomial system has no roots