The source files for all examples can be found in /examples.
Lovász Theta Function
The Lovász theta function is an important concept in graph theory. Consider an undirected graph $G = (V,E)$ with vertex set $V = \{1,\dots,n\}$ and edge set $E$ with $E_{ij} = 1$ if there exists an edge between the vertices $i$ and $j$. A stable set (or independent set) is a subset of $V$ such that the induced subgraph does not contain edges. The stability number $\alpha(G)$ of the graph is equal to the cardinality of the largest stable set. It is closely related to the Shannon capacity $\Theta(G)$ that models the amount of information that a noisy communication channel can carry if certain signal values can be confused with each other. Unfortunately, the determination of $\alpha(G)$ is an NP-hard problem and the computational complexity to determine $\Theta(G)$ remains unknown. The Lovász theta function $\vartheta(G)$ was introduced in [1], can be computed in polynomial time, and represents an upper bound on both the stability number and the Shannon capacity:
\[\alpha(G) \leq \Theta(G) \leq \vartheta(G).\]
The value of $\vartheta(G)$ can be determined by computing the optimal value $p^*$ of the following SDP:
\[\begin{array}{ll} \text{maximize} & \text{Tr}(JX)\\ \text{subject to} & \text{Tr}(X) = 1, \\ & X_{ij} = 0, \quad (i, j) \in E, \\ & X \succeq 0, \end{array}\]
with matrix variable $X$, matrix $J \in \mathbf{R}^{n \times n}$ of all ones and $\text{Tr}()$ denoting the matrix trace. In this simple example we will compute the value of the Lovász theta function for the Petersen Graph (which is known to be $4$).
using LinearAlgebra, COSMO, JuMP, Plots, GraphRecipes
# let's define the Petersen graph
n = 10
E = zeros(n, n)
E[1, 2] = E[1, 5] = E[1, 6] = 1.
E[2, 3] = E[2, 7] = 1.
E[3, 4] = E[3, 8] = 1.
E[4, 5] = E[4, 9] = 1.
E[5, 10] = 1.
E[6, 8] = E[6, 9] = 1.
E[7, 9] = E[7, 10] = 1.
E[8, 10] = 1.
# plot the graph
ri = 1.
ro = 2.
coordinates = []
for θ = 90:-72:-198
push!(coordinates, [ro * cosd(θ), ro * sind(θ)])
end
for θ = 90:-72:-198
push!(coordinates, [ri * cosd(θ), ri * sind(θ)])
end
graphplot(E, names = 1:n, x = getindex.(coordinates, 1), y = getindex.(coordinates, 2), fontsize = 10, nodesize = 1, nodeshape =:circle, curvature = 0.)
Let's solve the SDP with COSMO and JuMP:
model = JuMP.Model(COSMO.Optimizer);
@variable(model, X[1:n, 1:n], PSD)
x = vec(X)
@objective(model, Max, sum(x))
@constraint(model, tr(X)== 1.)
for j = 1:n
for i = 1:j-1
if E[i, j] == 1.
@constraint(model, X[i, j] == 0)
end
end
end
status = JuMP.optimize!(model)Maximize ScalarAffineFunction{Float64}:
0.0 + 1.0 X[1,1] + 2.0 X[1,2] + 2.0 X[1,3] + 2.0 X[1,4] + 2.0 X[1,5] + 2.0 X[1,6] + 2.0 X[1,7] + 2.0 X[1,8] + 2.0 X[1,9] + 2.0 X[1,10] + 1.0 X[2,2] + 2.0 X[2,3] + 2.0 X[2,4] + 2.0 X[2,5] + 2.0 X[2,6] + 2.0 X[2,7] + 2.0 X[2,8] + 2.0 X[2,9] + 2.0 X[2,10] + 1.0 X[3,3] + 2.0 X[3,4] + 2.0 X[3,5] + 2.0 X[3,6] + 2.0 X[3,7] + 2.0 X[3,8] + 2.0 X[3,9] + 2.0 X[3,10] + 1.0 X[4,4] + 2.0 X[4,5] + 2.0 X[4,6] + 2.0 X[4,7] + 2.0 X[4,8] + 2.0 X[4,9] + 2.0 X[4,10] + 1.0 X[5,5] + 2.0 X[5,6] + 2.0 X[5,7] + 2.0 X[5,8] + 2.0 X[5,9] + 2.0 X[5,10] + 1.0 X[6,6] + 2.0 X[6,7] + 2.0 X[6,8] + 2.0 X[6,9] + 2.0 X[6,10] + 1.0 X[7,7] + 2.0 X[7,8] + 2.0 X[7,9] + 2.0 X[7,10] + 1.0 X[8,8] + 2.0 X[8,9] + 2.0 X[8,10] + 1.0 X[9,9] + 2.0 X[9,10] + 1.0 X[10,10]
Subject to:
VectorAffineFunction{Float64}-in-Zeros
┌ ┐
│-1.0 + 1.0 X[1,1] + 1.0 X[2,2] + 1.0 X[3,3] + 1.0 X[4,4] + 1.0 X[5,5] + 1.0 X[6,6] + 1.0 X[7,7] + 1.0 X[8,8] + 1.0 X[9,9] + 1.0 X[10,10]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[1,2]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[2,3]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[3,4]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[1,5]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[4,5]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[1,6]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[2,7]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[3,8]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[6,8]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[4,9]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[6,9]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[7,9]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[5,10]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[7,10]│
└ ┘ ∈ Zeros(1)
┌ ┐
│0.0 + 1.0 X[8,10]│
└ ┘ ∈ Zeros(1)
VectorAffineFunction{Float64}-in-Scaled{PositiveSemidefiniteConeTriangle}
┌ ┐
│0.0 + 1.0 X[1,1] │
│0.0 + 1.4142135623730951 X[1,2] │
│0.0 + 1.0 X[2,2] │
│0.0 + 1.4142135623730951 X[1,3] │
│0.0 + 1.4142135623730951 X[2,3] │
│0.0 + 1.0 X[3,3] │
│0.0 + 1.4142135623730951 X[1,4] │
│0.0 + 1.4142135623730951 X[2,4] │
│0.0 + 1.4142135623730951 X[3,4] │
│0.0 + 1.0 X[4,4] │
│0.0 + 1.4142135623730951 X[1,5] │
│0.0 + 1.4142135623730951 X[2,5] │
│0.0 + 1.4142135623730951 X[3,5] │
│0.0 + 1.4142135623730951 X[4,5] │
│0.0 + 1.0 X[5,5] │
│0.0 + 1.4142135623730951 X[1,6] │
│0.0 + 1.4142135623730951 X[2,6] │
│0.0 + 1.4142135623730951 X[3,6] │
│0.0 + 1.4142135623730951 X[4,6] │
│0.0 + 1.4142135623730951 X[5,6] │
│0.0 + 1.0 X[6,6] │
│0.0 + 1.4142135623730951 X[1,7] │
│0.0 + 1.4142135623730951 X[2,7] │
│0.0 + 1.4142135623730951 X[3,7] │
│0.0 + 1.4142135623730951 X[4,7] │
│0.0 + 1.4142135623730951 X[5,7] │
│0.0 + 1.4142135623730951 X[6,7] │
│0.0 + 1.0 X[7,7] │
│0.0 + 1.4142135623730951 X[1,8] │
│0.0 + 1.4142135623730951 X[2,8] │
│0.0 + 1.4142135623730951 X[3,8] │
│0.0 + 1.4142135623730951 X[4,8] │
│0.0 + 1.4142135623730951 X[5,8] │
│0.0 + 1.4142135623730951 X[6,8] │
│0.0 + 1.4142135623730951 X[7,8] │
│0.0 + 1.0 X[8,8] │
│0.0 + 1.4142135623730951 X[1,9] │
│0.0 + 1.4142135623730951 X[2,9] │
│0.0 + 1.4142135623730951 X[3,9] │
│0.0 + 1.4142135623730951 X[4,9] │
│0.0 + 1.4142135623730951 X[5,9] │
│0.0 + 1.4142135623730951 X[6,9] │
│0.0 + 1.4142135623730951 X[7,9] │
│0.0 + 1.4142135623730951 X[8,9] │
│0.0 + 1.0 X[9,9] │
│0.0 + 1.4142135623730951 X[1,10]│
│0.0 + 1.4142135623730951 X[2,10]│
│0.0 + 1.4142135623730951 X[3,10]│
│0.0 + 1.4142135623730951 X[4,10]│
│0.0 + 1.4142135623730951 X[5,10]│
│0.0 + 1.4142135623730951 X[6,10]│
│0.0 + 1.4142135623730951 X[7,10]│
│0.0 + 1.4142135623730951 X[8,10]│
│0.0 + 1.4142135623730951 X[9,10]│
│0.0 + 1.0 X[10,10] │
└ ┘ ∈ Scaled{PositiveSemidefiniteConeTriangle}(PositiveSemidefiniteConeTriangle(10))
------------------------------------------------------------------
COSMO v0.8.10 - A Quadratic Objective Conic Solver
Michael Garstka
University of Oxford, 2017 - 2022
------------------------------------------------------------------
Problem: x ∈ R^{55},
constraints: A ∈ R^{71x55} (80 nnz),
matrix size to factor: 126x126,
Floating-point precision: Float64
Sets: DensePsdConeTriangle of dim: 55 (10x10)
ZeroSet of dim: 16
Settings: ϵ_abs = 1.0e-05, ϵ_rel = 1.0e-05,
ϵ_prim_inf = 1.0e-04, ϵ_dual_inf = 1.0e-04,
ρ = 0.1, σ = 1e-06, α = 1.6,
max_iter = 5000,
scaling iter = 10 (on),
check termination every 25 iter,
check infeasibility every 40 iter,
KKT system solver: QDLDL
Acc: Anderson Type2{QRDecomp},
Memory size = 15, RestartedMemory,
Safeguarded: true, tol: 2.0
Setup Time: 171.08ms
Iter: Objective: Primal Res: Dual Res: Rho:
1 -1.2711e+03 2.1055e+01 1.0052e+00 1.0000e-01
25 -4.0000e+00 1.7503e-10 3.0642e-14 1.0000e-01
------------------------------------------------------------------
>>> Results
Status: Solved
Iterations: 25
Optimal objective: -4
Runtime: 0.574s (573.72ms)The optimal objective is given by:
JuMP.objective_value(model)4.000000007423991Which is the correct known value for the Petersen Graph.
References
[1] Lovász - On the Shannon Capacity of a Graph, IEEE Transactions on Information Theory (1979)
This page was generated using Literate.jl.