Strong check matrices

1024 days ago by hthall

This worksheet provides functions for production verification matrices or checking various strong properties.

Included properties:

  • SAP (Strong Arnold Property)
  • SSP (Strong Spectral Property)
  • SMP (Strong Multiplicity Property)
  • SIP (Strong Interlacing Property) under deletion of the first row and column.
# Function definitions def make_zerolist(sym, zerolist=None): """ Given a symmetric matrix sym, create or validate a list of places, one per symmetric pair of off-diagonal entries, where the entries of the matrix are zero. If the list is non-empty, assume that missing items in the list correspond to +/-0 entries, double edges that are allowed to vary. If no list is passed in, create one and populate it. In any of the functions below, if you want to know the meaning of each column of the verification matrix, pass in an empty list and upon return the same list will have column labels. """ nn = sym.ncols() # size of matrix if zerolist is None: zerolist = [] full_list = [] for ii in range(nn): for jj in range(ii + 1, nn): # strict upper triangular part if not sym[ii, jj]: # entry is zero if sym[jj, ii]: raise ValueError('Matrix is not combinatorially symmetric.') full_list.append((ii, jj)) if not sym[jj, ii] and sym[ii, jj]: raise ValueError('Matrix is not combinatorially symmetric.') if zerolist is None: return full_list if not len(zerolist): # empty list was passed in for xx in full_list: zerolist.append(xx) return zerolist # Non-empty list just requires validation for xx in zerolist: if (xx[0], xx[1]) not in full_list and (xx[1], xx[0]) not in full_list: raise ValueError('zerolist value ' + str(xx) + ' is not a zero location.') return zerolist print("make_zerolist(sym, zerolist=None): Create column labels for the X matrix.") print("") def checkMatrixSAP(sym, zerolist=None): """ Given a symmetric matrix sym, produce a verification matrix for the Strong Arnold Property, whose columns are independent iff sym has SAP. """ zerolist = make_zerolist(sym, zerolist) nn = sym.ncols() # size of matrix # blank starting matrix outMatrix = Matrix(sym.base_ring(), nn*nn, len(zerolist)) col = 0 # one column for each for pair in zerolist: # zero pair in the matrix row = 0 # one row for each entry of for kk in range(nn): # the matrix AX, which should be zero. for ll in range(nn): # record the value of AX in this position. if ll == pair[0]: outMatrix[row, col] += sym[kk, pair[1]] if ll == pair[1]: outMatrix[row, col] += sym[pair[0], kk] row += 1 # advance a row in the check matrix col += 1 # advance a column in the check matrix return outMatrix # A has SAP if the rank is the number of columns print("checkMatrixSAP(sym, zerolist=None): Produce a check matrix for SAP.") print("") def hasSAP(sym): """ Return True if sym has SAP. """ check = checkMatrixSAP(sym) return rank(check) == check.ncols() print("hasSAP(sym): Return True if sym has SAP.") print("") def checkMatrixSSP(sym, zerolist=None): """ Given a symmetric matrix sym, produce a verification matrix for the Strong Spectral Property, whose columns are independent iff sym has SSP. """ zerolist = make_zerolist(sym, zerolist) nn = sym.ncols() # size of matrix # blank starting matrix outMatrix = Matrix(sym.base_ring(), nn*(nn-1)/2, len(zerolist)) col = 0 # one column for each for pair in zerolist: # zero pair in the matrix row = 0 # one row for each skew-symmetric for kk in range(nn): # basis element K with 1 in kk, ll for ll in range(kk+1, nn): # and -1 in ll, kk. # record the value of AK - KA in the zero pair place if ll == pair[1]: outMatrix[row, col] += sym[pair[0], kk] if kk == pair[1]: outMatrix[row, col] -= sym[pair[0], ll] if ll == pair[0]: outMatrix[row, col] += sym[kk, pair[1]] if kk == pair[0]: outMatrix[row, col] -= sym[ll, pair[1]] row += 1 # advance a row in the check matrix col += 1 # advance a column in the check matrix return outMatrix # A has SSP if the rank is the number of columns print("checkMatrixSSP(sym, zerolist=None): Produce a check matrix for SSP.") print("") def hasSSP(sym): """ Return True if sym has SSP. """ check = checkMatrixSSP(sym) return rank(check) == check.ncols() print("hasSSP(sym): Return True if sym has SSP.") print("") def checkMatrixSIP(sym, zerolist=None): """ Given a symmetric matrix sym, produce a verification matrix for the Strong Interlacing Property, whose columns are independent iff sym has SIP. """ zerolist = make_zerolist(sym, zerolist) nn = sym.ncols() # size of matrix # blank starting matrix outMatrix = Matrix(sym.base_ring(), (nn-1)*(nn-2)/2, len(zerolist)) col = 0 # one column for each for pair in zerolist: # zero pair in the matrix row = 0 # one row for each skew-symmetric for kk in range(1, nn): # basis element K with 1 in kk, ll for ll in range(kk+1, nn): # and -1 in ll, kk. # record the value of AK - KA in the zero pair place if ll == pair[1]: outMatrix[row, col] += sym[pair[0], kk] if kk == pair[1]: outMatrix[row, col] -= sym[pair[0], ll] if ll == pair[0]: outMatrix[row, col] += sym[kk, pair[1]] if kk == pair[0]: outMatrix[row, col] -= sym[ll, pair[1]] row += 1 # advance a row in the check matrix col += 1 # advance a column in the check matrix return outMatrix # A has SSP if the rank is the number of columns print("checkMatrixSIP(sym, zerolist=None): Produce a check matrix for SIP.") print("") def hasSIP(sym): """ Return True if sym has SIP. """ check = checkMatrixSIP(sym) return rank(check) == check.ncols() print("hasSIP(sym): Return True if sym has SIP.") print("") def checkMatrixSMP(sym, zerolist=None, q=None, find_q=True): """ Given a symmetric matrix sym, produce a verification matrix for the Strong Multiplicity Property, whose columns are independent iff sym has SSP. The number of rows depends on an upper bound q on the number of distinct eigenvalues of sym. If q is passed in, this is accepted as a known upper bound. If no value for q is passed in (or None is passed in), then calculate the degree of the minimal polynomial by default; otherwise use the universal upper bound of n. """ zerolist = make_zerolist(sym, zerolist) nn = sym.ncols() # size of matrix if q is None: if find_q: q = sym.minpoly().degree() else: q = nn # universal upper bound powerlist = [] # list of powers of sym power = sym for ii in range(2, q): # powers 2 through q - 1 power *= sym powerlist.append(power + 0) # append a copy, just to be safe # blank starting matrix outMatrix = Matrix(sym.base_ring(), nn*(nn-1)/2 + len(powerlist), len(zerolist)) col = 0 # one column for each for pair in zerolist: # zero pair in the matrix row = 0 # one row for each skew-symmetric for kk in range(nn): # basis element K with 1 in kk, ll for ll in range(kk+1, nn): # and -1 in ll, kk. # record the value of AK - KA in the zero pair place if ll == pair[0]: outMatrix[row, col] += sym[kk, pair[1]] if kk == pair[0]: outMatrix[row, col] -= sym[ll, pair[1]] if ll == pair[1]: outMatrix[row, col] += sym[kk, pair[0]] if kk == pair[1]: outMatrix[row, col] -= sym[ll, pair[0]] row += 1 # advance a row in the check matrix for ii in range(2, q): outMatrix[row, col] += powerlist[ii - 2][pair[0], pair[1]] row += 1 # progress through the extra rows for powers col += 1 # advance a column in the check matrix return outMatrix # A has SMP if the rank is the number of columns print("checkMatrixSMP(sym, zerolist=None): Produce a check matrix for SMP.") print("") def hasSMP(sym): """ Return True if sym has SMP """ check = checkMatrixSMP(sym) return rank(check) == check.ncols() print("hasSMP(sym): Return True if sym has SMP.") print("") def strongXList(sym, prop='SAP'): """ Return a basis for the space of X matrices which, if nonzero, show that sym does not have property prop, which is one of 'SAP', 'SSP', 'SMP', or 'SIP'. """ compedge = [] prop = prop.upper() if prop == 'SAP': checkmat = checkMatrixSAP(sym, compedge) elif prop == 'SMP': checkmat = checkMatrixSMP(sym, compedge) elif prop == 'SSP': checkmat = checkMatrixSSP(sym, compedge) elif prop == 'SIP': checkmat = checkMatrixSIP(sym, compedge) else: raise ValueError("Choices for prop are 'SAP', 'SMP', 'SSP', and 'SIP'.") matlist = [] ker = checkmat.transpose().kernel().basis_matrix() for ii in range(ker.nrows()): thismat = Matrix(ker.base_ring(), sym.ncols()) for jj in range(ker.ncols()): thismat[compedge[jj][0], compedge[jj][1]] = ker[ii,jj] thismat[compedge[jj][1], compedge[jj][0]] = ker[ii,jj] matlist.append(thismat) return matlist print("strongXList(sym, prop='SAP'): Produce a basis of X matrices making prop='SAP', 'SSP', 'SIP', or 'SMP' fail.") print("") 
       
make_zerolist(sym, zerolist=None): Create column labels for the X
matrix.

checkMatrixSAP(sym, zerolist=None): Produce a check matrix for SAP.

hasSAP(sym): Return True if sym has SAP.

checkMatrixSSP(sym, zerolist=None): Produce a check matrix for SSP.

hasSSP(sym): Return True if sym has SSP.

checkMatrixSIP(sym, zerolist=None): Produce a check matrix for SIP.

hasSIP(sym): Return True if sym has SIP.

checkMatrixSMP(sym, zerolist=None): Produce a check matrix for SMP.

hasSMP(sym): Return True if sym has SMP.

strongXList(sym, prop='SAP'): Produce a basis of X matrices making
prop='SAP', 'SSP', 'SIP', or 'SMP' fail.
make_zerolist(sym, zerolist=None): Create column labels for the X matrix.

checkMatrixSAP(sym, zerolist=None): Produce a check matrix for SAP.

hasSAP(sym): Return True if sym has SAP.

checkMatrixSSP(sym, zerolist=None): Produce a check matrix for SSP.

hasSSP(sym): Return True if sym has SSP.

checkMatrixSIP(sym, zerolist=None): Produce a check matrix for SIP.

hasSIP(sym): Return True if sym has SIP.

checkMatrixSMP(sym, zerolist=None): Produce a check matrix for SMP.

hasSMP(sym): Return True if sym has SMP.

strongXList(sym, prop='SAP'): Produce a basis of X matrices making prop='SAP', 'SSP', 'SIP', or 'SMP' fail.

As an example, we use a verification matrix to show that $P_4$ with a marked vertex of degree 2 cannot achieve $\mu_1=\mu_2$ with the Strong Interlacing Property.

First we define a polynomial ring with enough variables to parametrize $A \in \mathcal{S}(P_4)$.

P = PolynomialRing(QQ, ['d0', 'd1', 'd2', 'd3', 'a01', 'a02', 'a23']) P.inject_variables() print(P) 
       
Defining d0, d1, d2, d3, a01, a02, a23
Multivariate Polynomial Ring in d0, d1, d2, d3, a01, a02, a23 over
Rational Field
Defining d0, d1, d2, d3, a01, a02, a23
Multivariate Polynomial Ring in d0, d1, d2, d3, a01, a02, a23 over Rational Field
A = Matrix(P, [[d0, a01, a02, 0], [a01, d1, 0, 0], [a02, 0, d2, a23], [0, 0, a23, d3]]) show(A) 
       

                                
                            

                                

The graph looks correct.

G = Graph(A) G 
       

Deleting vertex 0 produces a disconnected graph: an isolated vertex plus a path on two vertices.

show(A[[1,2, 3], [1, 2, 3]]) 
       

                                
                            

                                
G.delete_vertex(0) G 
       

The path of length $1$ must have distinct eigenvalues so $A(1)$ has a repeated eigenvalue if and only if one of the eigenvalues is shared, which happens if and only if $A[\{2, 3\}] - d_1$ is a singular matrix. We solve for this determinant $= 0$.

(A[[2,3],[2,3]] - d1).determinant() 
       
d1^2 - d1*d2 - d1*d3 + d2*d3 - a23^2
d1^2 - d1*d2 - d1*d3 + d2*d3 - a23^2

It looks like solving for $d_3$ is a reasonable choice. To use the Sage solve routine requires casting values from our polynomial ring to the Symbolic Ring SR.

solve(SR((A[[2,3],[2,3]] - d1).determinant()), SR(d3)) 
       
[d3 == -(a23^2 - d1^2 + d1*d2)/(d1 - d2)]
[d3 == -(a23^2 - d1^2 + d1*d2)/(d1 - d2)]
A = Matrix(P, [[d0, a01, a02, 0], [a01, d1, 0, 0], [a02, 0, d1 + d2, a23], [0, 0, a23, d3]]) show(A) 
       

                                
                            

                                
solve(SR((A[[2,3],[2,3]] - d1).determinant()), SR(d3)) 
       
[d3 == (a23^2 + d1*d2)/d2]
[d3 == (a23^2 + d1*d2)/d2]

This is almost a polynomial. Reparametrizing again to replace $a_{23}$ by $a_{23}d_2$ would make it one.

A = Matrix(P, [[d0, a01, a02, 0], [a01, d1, 0, 0], [a02, 0, d1 + d2, a23*d2], [0, 0, a23*d2, d3]]) show(A) 
       

                                
                            

                                
solve(SR((A[[2,3],[2,3]] - d1).determinant()), SR(d3)) 
       
[d3 == a23^2*d2 + d1]
[d3 == a23^2*d2 + d1]

Replacing $d_3$ by this quantity gives us a full parametrization for the matrices of interest, with the requirement that $a_{01}, a_{02}, a_{23},$ and $d_2$ are all nonzero.

A = Matrix(P, [[d0, a01, a02, 0], [a01, d1, 0, 0], [a02, 0, d1 + d2, a23*d2], [0, 0, a23*d2, a23^2*d2 + d1]]) show(A) 
       

                                
                            

                                

We can check that the submatrix on indices $\{1,2, 3\}$ has a repeated eigenvalue.

show(A[[1, 2, 3], [1, 2, 3]].characteristic_polynomial().factor()) 
       

                                
                            

                                

Now we are ready to produce an SIP verification matrix for $A$. In order to know which zero entry each column of the verification matrix represents, we create a named empty list and feed it to the routine.

zlist = [] VM = checkMatrixSIP(A, zlist) show(VM) print(zlist) 
       
[(0, 3), (1, 2), (1, 3)]
[(0, 3), (1, 2), (1, 3)]

The three columns are not independent, and so this matrix does not have SIP.

VM.determinant() 
       
0
0

The routine strongXList produces a basis for the space of matrices $X$ that make a strong property fail.

This tends to complain for matrices over polynomial rings, so first we cast $A$ to a matrix over the fraction field.

PQ = (d0/d0).parent() PQ 
       
Fraction Field of Multivariate Polynomial Ring in d0, d1, d2, d3, a01,
a02, a23 over Rational Field
Fraction Field of Multivariate Polynomial Ring in d0, d1, d2, d3, a01, a02, a23 over Rational Field
A = Matrix(PQ, A) show(A) 
       

                                
                            

                                
Xlist = strongXList(A, 'SIP') show(Xlist) 
       

                                
                            

                                

The matrix $X$ that shows that $A$ fails to have SIP is unique up to scaling. Clearing the denominator produces an example over the original polynomial ring.

X = Matrix(P, Xlist[0] * a23) show(X) 
       

                                
                            

                                

The commutant of $A$ and $X$ is not zero, but taking the commutant and deleting the first row and column does produce a zero matrix.

show(A*X - X*A) 
       

                                
                            

                                
show((A*X - X*A)[[1, 2, 3], [1, 2, 3]]) 
       

                                
                            

                                

The matrix $A$ does, however, have the SSP, as can be seen from the SSP verification matrix.

zlist = [] VSSP = checkMatrixSSP(A, zlist) show(VSSP) show(zlist) 
       

                                
                            

                                

These columns are independent because deleting rows 0, 3, and 4 produces a square matrix that is invertible as long as $a_{01}$ and $a_{02}$ are both nonzero.

show(VSSP[[1, 2, 5], :]) show(VSSP[[1, 2, 5], :].determinant())