[1]:
using QuinticSpectrahedra
[2]:
using HomotopyContinuation, LinearAlgebra
Computing the combinatorial types of a quintic spectrahedra.
Given a linear space of \(5 \times 5\) real symmetric matrices of the form
The corresponding \(\textit{quintic spectrahedron}\) is the collection
The algebraic boundary of this spectrahedron is the variety
called a \(\textit{quintic symmetroid}\).
Quintic symmetroids generically have \(20\) singularities, some of which are real and some of which are nonreal. The real singularities may be partitioned into those which are on the spectrahedron, and those which are not. The \(\textit{combinatorial type}\) of a spectrahedron is \((\rho,\sigma)\) where \(\rho\) counts the number of real singularities on the symmetroid, and \(\sigma\) counts the real singularities on the spectrahedron.
Example: computing the singularities on 100 random quintic symmetroids
[3]:
@var x[1:4], a[1:3,1:5,1:5]
A=vcat([LinearAlgebra.I],[Symmetric(a[i,:,:]) for i in 1:3])
Ax=sum([x[i]*A[i] for i in 1:4])
eqs=[differentiate(det(Ax),x[i]) for i in 1:4]
affine_chart=sum(x.*randn(Float64,4))+randn()
eqs=vcat(eqs,affine_chart)
symm_parameters = unique(vcat(vec(A[2]),vec(A[3]),vec(A[4])))
quintic_system=System(eqs;parameters=symm_parameters)
init_params=randn(ComplexF64,45)
init_sols=solve(quintic_system,target_parameters=init_params)
samples=solve(quintic_system,
solutions(init_sols),
start_parameters=init_params,
target_parameters=[randn(Float64,45) for i in 1:100]);
Solving for 100 parameters... 100%|█████████████████████| Time: 0:00:06
# parameters solved: 100
# paths tracked: 2000
\(\texttt{samples}\) is now a list of \(100\) pairs:
\(\bullet\) the second element in the pair is a list of the \(15 \times 3 =45\) entries of the matrices \(A_2,A_3,A_4\)
\(\bullet\) the first element in the pair is a list of the \(20\) nodes on \(V(\text{det}(A(x)))\)
Computing the combinatorial type of a quintic spectrahedron
We can verify the combinatorial type of one such pair using the following function.
\(\textbf{Disclaimer}\): the computation of the combinatorial type below is numerical. Thus, it is subject to error.
(Next, we will introduce a function which certifies the combinatorial type)
[4]:
S=samples[1]
S_sols=S[1]
S_params=S[2]
C=combinatorial_type(S)
[4]:
(6, 2)
Certifying the combinatorial type of a quintic spectrahedron
The function Certify_Combinatorial_Type takes the entries of \(A_2,A_3,A_4\) as input.
It computes a square subsystem of the polynomial system above, solves for the \(64=20+44\) solutions, and then provably determines which \(20\) are solutions to the full system. For each of these \(20\) nodes, it also computes the values of its principal minors and then certifies their signs, thus determining the combinatorial type.
[5]:
Certify_Combinatorial_Type(S_params);
System solved
Certification performed
Every sign was certifiable
Certified 20 nodal solutions
Certified 6 real solutions
Certified 2 semidefinite solutions (0 PSD & 2 NSD)
Certification was successful
Certifying each type of quintic spectrahedra
The dictionary, Witness, encodes witnesses of each type of quintic spectrahedron, indexed by their combinatorial type.
Below, we certify that there exists a quintic spectrahedron with \(20\) nodes on its Euclidean boundary.
[6]:
W=Witness[20,20];
[7]:
Certify_Combinatorial_Type(W);
System solved
Certification performed
Every sign was certifiable
Certified 20 nodal solutions
Certified 20 real solutions
Certified 20 semidefinite solutions (0 PSD & 20 NSD)
Certification was successful
The following command will produce the defining polynomial for a spectrahedron of the given combinatorial type. You can copy the string into surf to produce a visualization of the symmetroid. Note the affine chart chosen to produce this polynomial is \(\text{tr}(A(x))=5\).
[8]:
str=jsurf((10,10))
[8]:
"1 - 1.11022302462516e-16*x + 1.11022302462516e-16*y + 13.5340331245054*x*y - 7.50082190675219*x*y^2 - 30.5947158533286*x*y^3 + 16.3556137746578*x*y^4 + 10.1250625635179*x*z - 13.5238469359924*x*z^2 - 34.7107588504292*x*z^3 + 46.2797918544673*x*z^4 - 11.8011128061791*x^2*y + 106.520994774996*x^2*y^2 - 36.1346631378685*x^2*y^3 - 0.912883555677102*x^2*z + 93.1119661438995*x^2*z^2 - 40.7493198322911*x^2*z^3 - 112.40944986166*x^3*y - 2.57851657765584*x^3*y^2 - 50.2183265153595*x^3*z + 34.7187817032043*x^3*z^2 + 87.3087978248684*x^4*y - 18.3629244608769*x^4*z - 4.88506621973035*y*z - 0.0839712657765846*y*z^2 + 24.6476526523819*y*z^3 - 17.589041046132*y*z^4 + 1.87109089280124*y^2*z + 51.980046097328*y^2*z^2 - 45.2540998481445*y^2*z^3 + 19.0583489603145*y^3*z - 34.1059841873924*y^3*z^2 - 13.6372255329615*y^4*z - 4.81380504464746*x*y*z - 103.884695693223*x*y*z^2 + 101.846614227372*x*y*z^3 - 86.5420161963045*x*y^2*z + 143.840185198655*x*y^2*z^2 + 88.8512079931257*x*y^3*z + 111.442737741989*x^2*y*z - 117.912907068058*x^2*y*z^2 - 153.635803466823*x^2*y^2*z + 109.985405146789*x^3*y*z - 18.7032651432321*x^2 + 13.3748754238749*x^3 + 61.4555517635348*x^4 - 46.351832498193*x^5 - 7.04515792259221*y^2 + 5.31535175043697*y^3 + 4.0935232861909*y^4 - 2.46438263066026*y^5 - 10.5285778845096*z^2 + 4.40736695602141*z^3 + 19.4756254479112*z^4 - 7.58970631225188*z^5"