A Bundle Approach for SDPs with Exact Subgraph Constraints

TL;DR

This paper introduces a bundle method-based computational framework for solving large-scale semidefinite relaxations with exact subgraph constraints, improving efficiency in tackling NP-hard graph problems.

Contribution

It proposes a novel partial Lagrangian dual approach that decomposes the problem, enabling the use of bundle methods to efficiently solve SDPs with many subgraph constraints.

Findings

Efficiently solves SDPs with many subgraph constraints.

Demonstrates improved performance on Max-Cut, stable set, and coloring problems.

Shows the approach's scalability and effectiveness.

Abstract

The 'exact subgraph' approach was recently introduced as a hierarchical scheme to get increasingly tight semidefinite programming relaxations of several NP-hard graph optimization problems. Solving these relaxations is a computational challenge because of the potentially large number of violated subgraph constraints. We introduce a computational framework for these relaxations designed to cope with these difficulties. We suggest a partial Lagrangian dual, and exploit the fact that its evaluation decomposes into two independent subproblems. This opens the way to use the bundle method from non-smooth optimization to minimize the dual function. Computational experiments on the Max-Cut, stable set and coloring problem show the efficiency of this approach.

Tables

Table 1. Table 1: The running times for one Max-Cut and one stable set instance with different fixed sets of ESCs. The graphs of order n = 100 𝑛 100 n=100 are from the Erdős-Rényi model.

					interior point		our bundle
	#ESC of size			$b$	time (sec)		time (sec)		$\%$ of MOSEK
	$3$	$4$	$5$		MOSEK	SDPT3	oracle	overall	time	value
MC	2000		0	6000	18.37	49.22	1.01	6.05	32.93	97.20
	2000		300	9000	55.24	134.78	1.18	9.33	16.90	95.02
	4000		0	12000	104.56	289.78	1.71	11.13	10.64	93.66
	3000		600	15000	184.43	525.85	1.56	14.83	8.04	94.54
SS	1050	0	0	5914	23.54	79.25	7.86	10.65	45.22	98.25
	1050	212	63	8719	50.11	174.33	10.61	16.52	32.96	97.89
	2100	0	0	11780	126.40	388.07	7.43	12.27	9.71	93.65
	1575	318	212	14653	241.29	648.83	10.79	20.21	8.38	94.44

Table 2. Table 2: The deviation of the ESB to z m ​ c subscript 𝑧 𝑚 𝑐 z_{mc} for several Max-Cut instances.

name	3	5	7	$z_{mc}$
pw01-100.1	0.40	0.00	0.00	2060
pw05-100.1	0.90	0.51	0.39	8045
pw09-100.1	0.58	0.38	0.31	13417
w01-100.1	0.13	0.00	0.00	719
w05-100.1	3.91	1.41	0.85	1606
w09-100.1	8.06	5.66	5.09	2096

Table 3. Table 3: The gap of the ESB to z m ​ c subscript 𝑧 𝑚 𝑐 z_{mc} for two Max-Cut instances.

name	BBnodes	root gap	7-gap	$z_{mc}$
beas-250-6	223	1.02	0.00	41014
beas-250-8	4553	2.19	0.49	35726

Table 4. Table 4: Tighten ϑ ​ ( G ) italic-ϑ 𝐺 \vartheta(G) towards α ​ ( G ) 𝛼 𝐺 \alpha(G) for several instances for 10 cycles.

name	$n$	$m$	$\vartheta \left(G\right)$	$s_{ss}^2$	$s_{ss}^3$	$s_{ss}^4$	$s_{ss}^5$	$s_{ss}^6$	$\alpha \left(G\right)$
CubicVT26_5	26	39	11.82	11.82	11.00	10.98	10.54	10.46	10
Circulant47_030	47	282	14.30	14.30	13.61	13.21	13.24	13.14	13
G_50_0_5	50	308	13.56	13.46	13.13	12.96	12.82	12.67	12
hamming6_4	64	1312	5.33	4.00	4.00	4.00	4.00	4.00	4
spin5	125	375	55.90	55.90	50.42	50.17	50.00	50.00	50
keller4	171	5100	14.01	13.70	13.54	13.50	13.49	13.49	11
sanr200_0_9	200	2037	49.27	49.04	48.94	48.86	48.78	48.75	42
c_fat200_5	200	11427	60.35	60.34	58.00	58.00	58.00	58.00	58

Table 5. Table 5: Maximum found projection distance of X I subscript 𝑋 𝐼 X_{I} to STAB 2 ⁡ ( G I ) superscript STAB 2 subscript 𝐺 𝐼 \operatorname{STAB}^{2}(G_{I}) for the computations of Table 4 .

		beginning			end
name	$n$	$s_c^2$	$s_c^4$	$s_c^6$	$s_c^2$	$s_c^4$	$s_c^6$
CubicVT26_5	26	0.000	0.102	0.193	0.000	0.029	0.013
G_50_0_5	50	0.087	0.093	0.118	0.000	0.013	0.024
spin5	125	0.000	0.084	0.269	0.000	0.046	0.006
sanr200_0_9	200	0.044	0.062	0.107	0.072	0.028	0.020

Table 6. Table 6: Tighten ϑ ​ ( G ) italic-ϑ 𝐺 \vartheta(G) towards χ ​ ( G ) 𝜒 𝐺 \chi(G) for several instances for 10 cycles.

name	$n$	$m$	$\vartheta \left(G\right)$	$s_c^2$	$s_c^3$	$s_c^4$	$s_c^5$	$s_c^6$	$\chi \left(G\right)⩽$
myciel4	23	71	2.53	2.53	2.90	2.91	3.28	3.29	5
myciel5	47	236	2.64	2.64	3.05	3.09	3.45	3.45	6
mug88_1	88	146	3.00	3.00	3.00	3.00	3.00	3.00	4
1_FullIns_4	93	593	3.12	3.12	3.25	3.37	3.80	3.80	5
myciel6	95	755	2.73	2.73	3.02	3.09	3.57	3.51	7
myciel7	191	2360	2.82	2.82	3.02	3.08	3.63	3.50	8
2_FullIns_4	212	1621	4.06	4.06	4.32	4.38	4.66	4.68	6
flat300_26_0	300	21633	16.99	17.04	17.12	17.10	17.12	17.12	26

Table 7. Table 7: Maximum found projection distance of X I subscript 𝑋 𝐼 X_{I} to COL ⁡ ( G I ) COL subscript 𝐺 𝐼 \operatorname{COL}(G_{I}) for the computations of Table 6 .

		beginning			end
name	$n$	$s_c^2$	$s_c^4$	$s_c^6$	$s_c^2$	$s_c^4$	$s_c^6$
myciel4	23	0.000	0.365	0.760	0.000	0.000	0.000
1_FullIns_4	93	0.009	0.349	0.629	0.000	0.158	0.203
myciel7	191	0.000	0.356	0.621	0.000	0.207	0.272
flat300_26_0	300	0.127	0.279	0.360	0.143	0.142	0.091

Equations

$$z_{mc}=\max \{⟨L,X⟩:X\in \mathrm{CUT}\}.$$

$$\max \left\{⟨L,X⟩:X\in {\mathcal{X}}^E\right\}$$

$${\mathcal{X}}^{S}=\left\{X\in{\mathcal{S}}_{n}:~{}X_{ij}=0\quad\forall\{i,j\}\in E,~{}x=\operatorname{diag}(X),~{}\left(\begin{array}[]{cc}1&x^{{T}}\\ x&X\end{array}\right)\succcurlyeq 0\right\},$$

$$\vartheta (G)=\max \left\{\mathrm{trace}(X):X\in {\mathcal{X}}^S\right\}$$

$$\operatorname{COL}^{\varepsilon}(G)=\operatorname{conv}\left\{\left(\begin{array}[]{cc}k&\mathbbm{1}^{{T}}\\ \mathbbm{1}&X\end{array}\right)=\sum_{i=1}^{k}\binom{1}{s_{i}}\binom{1}{s_{i}}^{{T}}:\begin{array}[]{c}S=(s_{1},\ldots,s_{k})\text{ is a}\\ \text{SSPM of }G,~{}X=SS^{{T}}\end{array}\right\}.$$

$${\mathcal{X}}^{C}=\left\{\left(\begin{array}[]{cc}t&\mathbbm{1}^{{T}}\\ \mathbbm{1}&X\end{array}\right)\succcurlyeq 0:~{}\operatorname{diag}(X)=\mathbbm{1}_{n},~{}X_{ij}=0~{}\forall\{i,j\}\in E\right\}$$

$$t^{*}(G)=\min\left\{t:~{}\left(\begin{array}[]{cc}t&\mathbbm{1}^{{T}}\\ \mathbbm{1}&X\end{array}\right)\in{\mathcal{X}}^{C}\right\}.$$

$$X_I\in {\mathrm{CUT}}_{|I|}.$$

$$\max \{⟨L,X⟩:X\in {\mathcal{X}}^E,X_I\in {\mathrm{CUT}}_{|I|}\forall I\in J\}$$

$$\max \{\mathrm{trace}(X):X\in {\mathcal{X}}^S,X_I\in {\mathrm{STAB}}^2(G_I)\forall I\in J\}.$$

$$\min\left\{t:~{}\left(\begin{array}[]{cc}t&\mathbbm{1}^{{T}}\\ \mathbbm{1}&X\end{array}\right)\in{\mathcal{X}}^{C},~{}X_{I}\in\operatorname{COL}(G_{I})~{}\forall I\in J\right\}.$$

$$X_I=\sum _{i=1}^{t_I}{\lambda }_i{C}_i^I,\lambda \in {\Delta }_{t_I},$$

$${\mathcal{A}}_I^{\top }({\lambda }_I)-{\mathcal{P}}_I(X)=0,{\lambda }_I\in {\Delta }_{t_I}.$$

$${\mathcal{M}}_I^{\top }({\mathcal{A}}_I^{\top }({\lambda }_I)-{\mathcal{P}}_I(X))=0,{\lambda }_I\in {\Delta }_{t_I},$$

$$z=\max \{⟨C,\hat{X}⟩:\hat{X}\in \mathcal{X},{\lambda }_I\in {\Delta }_{t_I},{\mathcal{M}}_I^{\top }({\mathcal{A}}_I^{\top }({\lambda }_I)-{\mathcal{P}}_I(X))=0\forall I\in J\},$$

$$\mathcal{L}(X,\lambda ,y)=⟨C,X⟩+\sum _{I\in J}⟨y_I,{\mathcal{M}}_I^{\top }({\mathcal{A}}_I^{\top }({\lambda }_I)-{\mathcal{P}}_I(X))⟩$$

$$z=\underset{y}{\min }\underset{\begin{array}{c}X\in \mathcal{X}\\ {\lambda }_I\in {\Delta }_{t_I}\end{array}}{\max }\mathcal{L}(X,\lambda ,y).$$

$$\underset{\begin{array}{c}X\in \mathcal{X}\\ {\lambda }_I\in {\Delta }_{t_I}\end{array}}{\max }⟨C-\sum _{I\in J}{\mathcal{P}}_I^{\top }{\mathcal{M}}_I(y_I),X⟩+\sum _{I\in J}⟨{\mathcal{A}}_I{\mathcal{M}}_I(y_I),{\lambda }_I⟩.$$

$$\underset{{\lambda }_I\in {\Delta }_{t_I}}{\max }⟨{\mathcal{D}}_I(y_I),{\lambda }_I⟩=\underset{1⩽i⩽t_I}{\max }{\left[{\mathcal{D}}_I(y_I)\right]}_i.$$

$$h(y)=\underset{X\in \mathcal{X}}{\max }⟨C-\sum _{I\in J}{\mathcal{P}}_I^{\top }{\mathcal{M}}_I(y_I),X⟩=⟨C-\sum _{I\in J}{\mathcal{P}}_I^{\top }{\mathcal{M}}_I(y_I),X^{\ast }⟩,$$

$$z=\underset{y}{\min }\left\{h(y)+\sum _{I\in J}\underset{1⩽i⩽t_I}{\max }{\left[{\mathcal{D}}_I(y_I)\right]}_i\right\}.$$

$$h(y)⩾\underset{1⩽j⩽r}{\max }\left\{\overline{h}-e_j+⟨g_j,y-\overline{y}⟩\right\}.$$

$$\underset{y,w,v}{\min }w$$

$$stw$$

$$v_I$$

Full text

1 Introduction

The study of NP-hard problems has led to the introduction of various hierarchies of relaxations, which typically involve several levels. Moving from one level to the next the relaxations get increasingly tighter and ultimately the exact optimum may be reached, but the computational effort grows accordingly.

Among the most prominent hierarchies are the polyhedral ones from Boros, Crama and Hammer [3] as well as the ones from Sherali and Adams [20], Lovász and Schrijver [15] and Lasserre [13] which are based on semidefinite programming (SDP). Even though on the starting level they have a simple SDP relaxation, already the first nontrivial level in the hierarchy requires the solution of SDPs in matrices of order $\binom{n}{2}$ and on level $k$ the matrix order is $n^{O(k)}$. Hence they are considered mainly as theoretical tools and from a practical point of view these hierarchies are of limited use.

Not all hierarchies are of this type. In [3], a polyhedral hierarchy for the Max-Cut problem is introduced which maintains $\binom{n}{2}$ variables in all levels, with a growing number of constraints. More recently, Adams, Anjos, Rendl and Wiegele [1] introduced a hierarchy of SDP relaxations which act in the space of symmetric $n\times n$ matrices and at level $k$ of the hierarchy all submatrices of order $k$ have to be ’exact’ in a well-defined sense, i.e. they have to fulfill an exact subgraph constraint (ESC).

It is the main purpose of this paper to describe an efficient way to optimize over level $k$ of this hierarchy for small values of $k$, e.g. $k\leqslant 6$, and demonstrate the efficiency of our approach for the Max-Cut, stable set and coloring problem.

Maintaining $\binom{n}{k}$ possible ESCs in an SDP in matrices of order $n$ is computationally infeasible even for $k=2$ or $k=3$, because each ESC creates roughly $\binom{k}{2}$ additional equality constraints and at most $2^{k}$ additional linear variables.

We suggest the following ideas to overcome this difficulty. First we proceed iteratively, and in each iteration we include only (a few hundred of) the most violated ESCs. More importantly, we propose to solve the dual of the resulting SDP. The structure of this SDP with ESCs admits a reformulation of the dual in the form of a non-smooth convex minimization problem with attractive features. First, any dual solution yields a valid bound for our relaxations, so it is not necessary to carry out the minimization to optimality. Secondly, the dual function evaluation decomposes into two independent problems. The first one is simply a sum of max-terms (one for each subgraph constraint), and the second one consists in solving a ’basic’ SDP, independent of the ESCs. The optimizer for this second problem also yields a subgradient of the objective function. With this information at hand we suggest to use the bundle method from non-smooth convex optimization. It provides an effective machinery to get close to a minimizer in few iterations.

As a result we are able to get near optimal solutions where all ESCs for small values of $k$ ($k\leqslant 6$) are satisfied up to a small error tolerance. Our computational results demonstrate the practical potential of this approach.

We finish this introductory section with some notation. We denote the vector of all-ones of size $n$ with $\mathbbm{1}_{n}$ and ${\Delta}_{n}=\{x\in{\mathbb{R}}^{n}_{+}:\sum_{i=1}^{n}x_{i}=1\}$. If the dimension is clear from the context we may omit the index and write $\mathbbm{1}$ and ${\Delta}$. Furthermore let $N=\{1,2,\dots,n\}$. A graph $G$ on $n$ vertices has vertex set $N$ and edge set $E$ and $\overline{G}$ is its complement graph. ${\mathcal{S}}_{n}$ is the set of $n$-dimensional symmetric matrices.

2 The Problems and their Semidefinite Relaxations

In the Max-Cut problem a symmetric matrix $L\in{\mathcal{S}}_{n}$ is given and $c\in\{-1,1\}^{n}$ which maximizes $c^{T}Lc$ should be determined. If the matrix $L$ corresponds to the Laplacian matrix of a (edge-weighted undirected) graph $G$, this is equivalent to finding a bisection of the vertices of $G$ such that the total weight of the edges joining the two bisection blocks is maximized. Such an edge set is also called a cut in $G$.

Bisections of $N$ can be expressed as $c\in\{-1,1\}^{n}$ where the two bisection blocks correspond to the entries in $c$ of the same sign. Given $c\in\{-1,1\}^{n}$ we call $C=cc^{{T}}$ a cut matrix. The convex hull of all cut matrices (of order $n$) is denoted by $\operatorname{CUT}_{n}$ or simply $\operatorname{CUT}$ if the dimension is clear. Since $c^{{T}}Lc=\langle L,cc^{{T}}\rangle$ Max-Cut can also be written as the following (intractable) linear program

[TABLE]

$\operatorname{CUT}$ is contained in the spectrahedron ${\mathcal{X}}^{E}=\left\{X\in{\mathcal{S}}_{n}:\operatorname{diag}(X)=\mathbbm{1}_{n},\,X\succcurlyeq 0\right\},$ hence

[TABLE]

is a basic semidefinite relaxation for Max-Cut. This model is well-known, attributed to Schrijver and was introduced in a dual form by Delorme and Poljak [4]. It can be solved in polynomial time to a fixed prescribed precision and solving this relaxation for $n=1000$ takes only a few seconds.

It is well-known that the Max-Cut problem is NP-hard. On the positive side, Goemans and Williamson [8] show that one can find a cut in a graph with nonnegative edge weights of value at least 0.878$z_{mc}$ in polynomial time.

In the stable set problem the input is an unweighted graph $G$. We call a set of vertices stable, if no two vertices are adjacent. Moreover we call a vector $s\in\{0,1\}^{n}$ a stable set vector if it is the incidence vector of a stable set. The convex hull of all stable set vectors of $G$ is denoted with $\operatorname{STAB}(G)$. In the stable set problem we want to determine the stability number $\alpha(G)$, which denotes the cardinality of a largest stable set in $G$, hence $\alpha(G)=\max\left\{\mathbbm{1}^{{T}}s:~{}s\in\operatorname{STAB}(G)\right\}.$ Furthermore we denote with $\operatorname{STAB}^{2}(G)=\operatorname{conv}\left\{ss^{{T}}:~{}s\in\operatorname{STAB}(G)\right\}$ the convex hull of all stable set matrices $ss^{{T}}$. Then with the arguments of Gaar [7] it is easy to check that $\alpha(G)=\max\{\operatorname{trace}(X):~{}X\in\operatorname{STAB}^{2}(G)\}.$ Furthermore $\operatorname{STAB}^{2}(G)$ is contained in the following spectrahedron

[TABLE]

which is known as the theta body in the literature. Therefore

[TABLE]

is a relaxation of the stable set problem. The Lovász theta function $\vartheta(G)$ was introduced in a seminal paper by Lovász [14]. We refer to Grötschel, Lovász and Schrijver [9] for a comprehensive analysis of $\vartheta(G)$.

Determining $\alpha(G)$ is again NP-hard. Contrary to Max-Cut, which has a polynomial time .878-approximation, for every $\varepsilon>0$ there can be no polynomial time algorithm that approximates $\alpha(G)$ within a factor better than $O(n^{1-\varepsilon})$ unless $P=NP$, see Håstad [11].

The coloring problem for a given graph $G$ consists in determining the chromatic number $\chi(G)$, which is the smallest $t$ such that $N$ can be partitioned into $t$ stable sets. Let $S=(s_{1},\ldots,s_{k})$ be a matrix where each column is a stable set vector and these stable sets partition $V$ into $k$ sets. Let us call such matrices $S$ stable-set partition matrices (SSPM). The $n\times n$ matrix $X=SS^{T}$ is called coloring matrix. The convex hull of the set of all coloring matrices of $G$ is denoted by $\operatorname{COL}(G)$. We also need the extended coloring polytope

[TABLE]

The difficult set $\operatorname{COL}^{\varepsilon}$ can be relaxed to the easier spectrahedron ${\mathcal{X}}^{C}$

[TABLE]

and we can consider the semidefinite program

[TABLE]

Obviously $t^{*}(G)\leqslant\chi(G)$ holds because the SSPM $S$ consisting of $\chi(G)$ stable sets yields a feasible coloring matrix $X=SS^{{T}}$ with objective function value $\chi(G)$. It is in fact a consequence of conic duality that $t^{*}(G)=\vartheta(\overline{G})$ holds.

It is NP-hard to find $\chi(G)$, to find a 4-coloring of a 3-colorable graph [10] and to color a $k$-colorable graph with $O(k^{\frac{\log k}{25}})$ colors for sufficiently large $k$, [12].

3 Exact Subgraph Hierarchy

In this section we will discuss how to systematically tighten the relaxations (1), (2) and (3) with ’exactness conditions’ imposed on small subgraphs. We obtained these relaxations by relaxing the feasible regions $\operatorname{CUT}$, $\operatorname{STAB}^{2}$ and $\operatorname{COL}$ of the integer problem to simple spectrahedral sets. Now we will use small subgraphs to get closer to original feasible regions again.

For $I\subseteq N$ we denote with $X_{I}$ the principal submatrix of $X$ corresponding to the rows and columns in $I$. Furthermore let $G_{I}$ be the induced subgraph of $G$ on the set of vertices $I$ and let $k_{I}=|I|$ be the cardinality of $I$.

We first look at the exact subgraph relaxations for Max-Cut. The exact subgraph constraint (ESC) on $I\subseteq N$, introduced in [1] by Adams, Anjos, Rendl and Wiegele, requires that the matrix $X_{I}$ corresponding to the subgraph $G_{I}$ lies in the convex hull of the cut matrices of $G_{I}$, that is

[TABLE]

In this case we say that $X$ is exact on $I$.

Now we want the ESCs to be fulfilled not only for one but for a certain selection of subgraphs. We denote with $J$ the set of subgraphs which we require to be exact and get the following SDP relaxation with ESCs for Max-Cut.

[TABLE]

We proceed analogously for the stable set problem in a graph $G$. The ESC of a subgraph $G_{I}$ for the stable set problem requires that $X_{I}\in\operatorname{STAB}^{2}(G_{I})$ holds and the SDP with ESCs for the stable set problem is

[TABLE]

Turning to the coloring problem, we analogously impose additional constraints of the form $X_{I}\in\operatorname{COL}(G_{I})$ to obtain the SDP with ESCs

[TABLE]

Note that in the case of the stable set and the coloring problem the polytopes $\operatorname{STAB}^{2}(G_{I})$ and $\operatorname{COL}(G_{I})$ depend on the subgraph $G_{I}$, whereas in Max-Cut the polytope $\operatorname{CUT}_{|I|}$ only depends on the number of vertices of the subgraph.

From a theoretical point of view, we obtain the $k$-th level of the exact subgraph hierarchy of [1] if we use $J=\{I\subseteq N:~{}|I|=k\}$ in the relaxations (4), (5) and (6) respectively. We denote the corresponding objective function values with $z_{mc}^{k}$, $z_{ss}^{k}$ and $z_{c}^{k}$. So the $k$-th level of the hierarchy is obtained by forcing all subgraphs on $k$ vertices to be exact in the basic SDP relaxation.

In the case of the stable set and the Max-Cut problem we have $z_{ss}^{n}=\alpha(G)$ (see [7]) and $z_{mc}^{n}=z_{mc}$. For coloring $z_{c}^{n}\leqslant\chi(G)$ holds. Let $z_{c\varepsilon}^{k}$ be the resulting value if we add the inequalities $t\geqslant\sum_{i=1}^{t_{I}}[\lambda_{I}]_{i}|S^{I}_{i}|$ where $|S^{I}_{i}|$ is the number of colors used for the SSPM $S^{I}_{i}$ and $\lambda_{I}\in{\Delta}_{t_{I}}$ is a variable for the convex combination for each subgraph $I$ to the SDP for $z_{c}^{k}$. Then $z_{c\varepsilon}^{n}=\chi(G)$ holds. Since the focus of this paper are computational results we are interested only in the computational results we omit the details and further theoretical investigations.

An important feature of this hierarchy is that the size of the matrix variable remains $n$ or $n+1$ on all levels of the hierarchy and only more linear variables and constraints (enforcing the ESCs, hence representing convex hull conditions) are added on higher levels. So it is possible to approximate $z_{mc}^{k}$, $z_{ss}^{k}$ and $z_{c}^{k}$ by forcing only some subgraphs of order $k$ to be exact. This is our key ingredient to computationally obtain tight bounds on $z_{mc}$, $\alpha(G)$ and $\chi(G)$.

From a practical point of view solving the relaxations (4), (5) and (6) with standard interior point (IP) solvers like SDPT3 [21] or MOSEK [16] is very time consuming. In Table 1 we list computation times (in seconds) for one specific Max-Cut and one specific stable set instance. We vary the number of ESCs for subgraphs of order $3$, $4$ and $5$, so we solve (4) and (5) for different $J$. We choose $J$ such that the total number of equality constraints induced by the convex hull formulation of the ESCs $b$ ranges between 6000 and 15000. Since the matrix order $n$ is fixed to $n=100$, the overall computation time depends essentially on the number of constraints, independent of the specific form of the objective function. Aside from the ESC constraints, we have $n$ additional equations for Max-Cut and $n+m+1$ additional equations for the stable set problem. Here $m$ denotes the number of edges of the graph. We have $m=722$ in the example graph. Clearly the running times get huge for a large number of ESC. Furthermore MATLAB requires 12 Gigabyte of memory for $b=15000$, showing also memory limitations.

Note that it is argued in [1] that $z_{mc}^{4}=z_{mc}^{3}$, so we omit subgraphs of order $k_{I}=4$ for Max-Cut. This is because in the back of our minds our final algorithm to determine the best possible bounds first includes ESCs of size $k$, starting for example with $k=3$. As soon as we do not find violated ESCs of size $k$ anymore, we repeat this for size $k+1$.

4 Partial Lagrangian Dual

To summarize we are interested in solving relaxations (4), (5) and (6) with a potentially large number of ESCs, where using interior point solvers is too time consuming. In this section we will first establish a unified formulation of the relaxations (4), (5) and (6). Then we will build the partial Lagrangian dual of this formulation, where only the ESCs are dualized. This model will be particularly amenable for the bundle method, because it will be straightforward to obtain a subgradient of the model when evaluating it at a certain point.

In order to unify the notation for the three problems observe that the ESCs $X_{I}\in\operatorname{CUT}_{|I|}$, $X_{I}\in\operatorname{STAB}^{2}(G_{I})$ and $X_{I}\in\operatorname{COL}(G_{I})$ can be represented as

[TABLE]

where $C^{I}_{i}$ is the $i$-th cut, stable set or coloring matrix of the subgraph $G_{I}$ and $t_{I}$ is their total number.

A formal description of ESC in (7) requires some additional notation. First we introduce the projection ${\mathcal{P}}_{I}\colon{\mathcal{S}}_{n}\mapsto{\mathcal{S}}_{k_{I}}$, mapping $X$ to the submatrix $X_{I}$. Second we define a map ${\mathcal{A}}_{I}\colon{\mathcal{S}}_{k_{I}}\mapsto{\mathbb{R}}^{t_{I}}$, such that its adjoint map ${\mathcal{A}}_{I}^{\top}\colon{\mathbb{R}}^{t_{I}}\mapsto{\mathcal{S}}_{k_{I}}$ is given by ${\mathcal{A}}_{I}^{\top}(\lambda)=\sum_{i=1}^{t_{I}}{\lambda_{i}C_{i}^{I}}$ and produces a linear combination of the cut, stable set or coloring matrices. Thus we can rewrite (7) as

[TABLE]

The left-hand side of the matrix equality is a symmetric matrix, of which some entries (depending on which problem we consider) are zero for sure, so we do not have to include all $k_{I}\times k_{I}$ equality constraints into the SDP. Let $b_{I}$ be the number of equality constraints we have to include. Note that $b_{I}=\binom{k_{I}}{2}$, $b_{I}=\binom{k_{I}+1}{2}-m_{I}$ and $b_{I}=\binom{k_{I}}{2}-m_{I}$ for the Max-Cut, stable set and coloring problem respectively, if $m_{I}$ denotes the number of edges of $G_{I}$. This is because in the case of the stable set problem we also have to include equations for the entries of the main diagonal contrary to Max-Cut and the coloring problem. Then we define a linear map ${\mathcal{M}}_{I}\colon{\mathbb{R}}^{b_{I}}\mapsto{\mathcal{S}}_{k_{I}}$ such that the adjoint operator ${\mathcal{M}}_{I}^{\top}\colon{\mathcal{S}}_{k_{I}}\mapsto{\mathbb{R}}^{b_{I}}$ extracts the $b_{I}$ positions, for which we have to include the equality constraints, into a vector. So eventually we can rephrase (8) equivalently as

[TABLE]

which are $b_{I}+1$ equalities and $t_{I}$ inequalities. In consequence all three relaxations (4), (5) and (6) have the generic form

[TABLE]

where $C$, ${\mathcal{X}}$, ${\mathcal{A}}_{I}$, ${\mathcal{M}}_{I}$ and $b_{I}$ have to be defined problem specific. Furthermore $\widehat{X}=X$ in the case of Max-Cut and stable set and \widehat{X}=\left(\begin{array}[]{cc}t&\mathbbm{1}^{{T}}\\ \mathbbm{1}&X\end{array}\right) for coloring, but for the sake of understandability we will just use $X$ in the following.

The key idea to get a handle on problem (9) is to consider the partial Lagrangian dual where the ESCs (without the constrains $\lambda_{I}\in{\Delta}_{t_{I}}$) are dualized. We introduce a vector of multipliers $y_{I}$ of size $b_{I}$ for each $I$ and collect them in $y=(y_{I})_{I\in J}$ and also collect $\lambda=(\lambda_{I})_{I\in J}$. The Lagrangian function becomes

[TABLE]

and standard duality arguments (Rockafellar [19, Corollary 37.3.2]) yield

[TABLE]

For a fixed set of multipliers $y$ the inner maximization becomes

[TABLE]

This maximization is interesting in at least two aspects. First, it is separable in the sense that the first term depends only on $X$ and the second one only on the separate $\lambda_{I}$. Moreover, if we denote the linear map ${\mathcal{A}}_{I}{\mathcal{M}}_{I}(y_{I})\colon{\mathbb{R}}^{b_{I}}\mapsto{\mathbb{R}}^{t_{I}}$ with ${\mathcal{D}}_{I}$, the second term has an explicit solution, namely

[TABLE]

In order to consider the first term in more detail, we define the following function. Let $b=\sum_{I\in J}b_{I}$ be the dimension of $y$. Then $h\colon{\mathbb{R}}^{b}\to{\mathbb{R}}$ is defined as

[TABLE]

where $X^{\ast}$ is a maximizer over the set ${\mathcal{X}}$ for $y$ fixed. Note that $h(y)$ is convex but non-smooth, but (12) shows that $g_{I}=-{\mathcal{M}}_{I}^{T}{\mathcal{P}}_{I}(X^{\ast})$ is a subgradient of $h$ with respect to $y_{I}$. By combining (11) and (12) we can reformulate the partial Lagrangian dual (10) to

[TABLE]

The formulation (13) of the original relaxations (4), (5) and (6) fits perfectly into the bundle method setting described by Frangioni and Gorgone in [6], hence we suggest to approach this problem using the bundle method.

5 Solving (13) with the Bundle Method

The bundle method is an iterative procedure for minimizing a convex non-smooth function and firstly maintains the current center $\overline{y}$, which represents the current estimate to the optimal solution, throughout the iterations. Secondly it maintains the bundle of the form $\mathcal{B}=\{(y_{1},h_{1},g_{1},X_{1}),\dots,(y_{r},h_{r},g_{r},X_{r})\}$. Here $y_{1},\ldots,y_{r}$ are the points which we use to set up our subgradient model. Moreover $h_{i}=h(y_{i})$, $g_{i}$ is a subgradient of $h$ at $y_{i}$ and $X_{i}$ is a maximizer of $h$ at $y_{i}$ as in (12).

At the start we select $y_{1}=\overline{y}=0$ and evaluate $h$ at $\overline{y}$, which yields the bundle $\mathcal{B}=\{(y_{1},g_{1},h_{1},X_{1})\}$. A general iteration consists of the two steps determining the new trial point and evaluating the oracle. For determining a new trial point $\widetilde{y}$ the subgradient information of the bundle $\mathcal{B}$ translates into the subgradient model $h(y)\geqslant h_{j}+\langle g_{j},y-y_{j}\rangle$ for all $j=1$, …, $r$. It is common to introduce $e_{j}=h(\overline{y})-h_{j}-\langle g_{j},\overline{y}-y_{j}\rangle$ for $j=1,$ …, $r$ and with $\overline{h}=h(\overline{y})$ the subgradient model becomes

[TABLE]

The right-hand side above is convex, piecewise linear and minorizes $h$. In each iteration of the bundle method we minimize the right-hand side of (14) instead of $h$, but ensure that we do not move too far from $\overline{y}$ by adding a penalty term of the form $\frac{1}{2}\mu\left\lVert y-\overline{y}\right\rVert^{2}$ for a parameter $\mu\in{\mathbb{R}}_{+}$ to the objective function. With the auxiliary variables $w\in{\mathbb{R}}$ and $v_{I}\in{\mathbb{R}}$ for all $I\in J$ to model the maximum terms and with $v=(v_{I})_{I\in J}\in{\mathbb{R}}^{q}$ and $q=|J|$ we end up with

[TABLE]

This is a convex quadratic problem in $1+q+b$ variables with $r+\sum_{I\in J}t_{I}$ linear inequality constraints. Its solution $(\widetilde{y},\widetilde{w},\widetilde{v})$ includes the new trial point $\widetilde{y}$. Problems of this type can be solved efficiently in various ways, see [7] for further details. In our implementation we view (15) as a rotated second order cone program with one second-order cone constraint and solve it with MOSEK.

The second step in each bundle iteration is to evaluate the dual function $h$ at $\widetilde{y}$. In our case determining $h(\widetilde{y})$ means solving the basic SDP relaxation as introduced in Section 2 with a modified objective function. Hence in the case of Max-Cut the oracle can be evaluated very quickly, whereas evaluating the oracle is computationally more expensive for the stable set and the coloring problem.

The bundle iteration finishes by deciding whether $\widetilde{y}$ becomes the new center (serious step, roughly speaking if the increase of the objective function is good) or not (null step). In either case the new point is included in the bundle, some other elements of the bundle are possibly removed, the bundle parameter $\mu$ is updated and a new iteration starts.

6 Computational Results and Conclusions

We close with a small sample of computational results and start with comparing our bundle method with interior point methods. In our context we are mostly interested to improve the upper bounds quickly, so we do not run the bundle method described in Section 5 until we reach a minimizer, but stop after a fixed number of iterations, say $30$. In Table 1 one sees that the running times decrease drastically if we use the bundle method. For $b\approx 15000$ it takes the bundle method only around $8\%$ of the MOSEK running time to get as close as $95\%$ to the optimal value, which is sufficient for our purposes. One sees that our bundle method scales much better for increasing $|J|$.

If we are given a graph and want to get an approximation on $z_{mc}^{k}$, $z_{ss}^{k}$ and $z_{c}^{k}$, then we iteratively perform a fixed number, say 30, iterations of the bundle method and then update the set $J$. We denote the exact subgraph bounds (ESB) obtained in this way with $s_{mc}^{k}$, $s_{ss}^{k}$ and $s_{c}^{k}$.

For the sake of brevity we will only outline how to determine $J$ heuristically, see [7] for details. Let $X^{\ast}$ be the current solution of (4), (5) or (6). We use the fact that the inner product of $X_{I}^{\ast}$ and particular matrices of size $k_{I}$ is potentially small whenever $X_{I}^{\ast}$ is not in $\operatorname{STAB}^{2}(G_{I})$. Minimizing this inner product over all subgraphs of order $k_{I}$ would yield a quadratic assignment problem, so we repeatedly use a local search heuristic for fixed particular matrices in order to obtain potential subgraphs. Then we calculate the projection distances from $X_{I}^{\ast}$ to $\operatorname{STAB}^{2}(G_{I})$ for all these subgraphs and include those in $J$ which have the largest distances and hence are violated most.

Finally we present several computational results for obtained ESBs. Note that we refrain from comparing the running times of our bundle method with the running time of inter point methods, because interior point methods would reach their limit very soon. Hence the bounds presented can only be obtained with our methods in reasonable time.

When considering Max-Cut the graphs in Table 2 are from the Biq Mac library [2] with $n=100$ vertices. The edge density is 10%, 50% and 90%. The first 3 instances have positive weights and the remaining 3 have also some negative weights. The column labeled $3$ provides the deviation (in %) of the ESB with $k=3$ from $z_{mc}$. Thus if $p$ is the value in the column labeled $3$, then $s_{mc}^{3}=(1+p/100)z_{mc}$. The columns labeled 5 and 7 are to be understood in a similar way for $k=5$ and $k=7$. We note that the improvement of the bound from column 3 to column 7 is quite substantial in all cases. We also point out that the relative gap is much larger if also negative edge weights are present.

In Table 3 we look at graphs from the Beasley collection [2] with $n=250$. These instances were used by Rendl, Rinaldi and Wiegele [18] in a Branch-and-Bound setting. We only consider the ’hardest’ instances from [18] where the Branch-and-Bound tree has more than 200 nodes. The table provides the gap at the root node and also the number of nodes in the Branch-and-Bound tree as reported in [18]. The column 7-gap contains the gap after solving our new relaxation with ESCs up to size $k=7$. We find it remarkable that the first instance is solved to optimality and the gap in the second instance is reduced by 75 % compared to the original gap. This implies that using our ESBs would expectedly reduce the very high number of required Branch-And-Bound nodes tremendously.

We conclude that for Max-Cut our ESB constitute a substantial improvement compared to the previously used strongest bounds based on SDP with triangle inequalities. These correspond to the column 3-gap.

For the calculations for the stable set and the coloring problem all instances are chosen in such a way that $\vartheta(G)$ does not coincide and is not very close to $\alpha(G)$ and $\chi(G)$ respectively.

The instances for the stable set problem are taken partly from the DIMACS challenge [5] with some additional instances from [7] with $n$ ranging from 26 to 200. Table 4 contains the new bounds. Here the starting point is the relaxation $\vartheta(G)$. We carry out 10 cycles of adding ESCs. In each cycle we add at most 200 ESCs, so in the final round we have no more than 2000 ESCs. The column heading indicates the order of the subgraphs. Here the improvement of the bounds is smaller than in the Max-Cut case, but we see that including larger subgraphs leads to much tighter bounds. In Table 5 we show that our approach also reduces the largest found projection distance over all subgraphs $G_{I}$ of $X_{I}$ to the corresponding $\operatorname{STAB}^{2}(G_{I})$ in the course of the cycles. This indicates that the violation of the subgraphs decreases over the cycles and less and less subgraphs do not fulfill the ESCs. For example the value $0.000$ for the graph spin5 for $s_{ss}^{2}$ at the end of the cycles means that we did not find a violated subgraph of order $2$ anymore.

Results for a selection of coloring instances from [17] are provided in Tables 6 and 7. As in the stable set case there is only little improvement using small subgraphs ($k=2$ or 3). The inclusion of larger subgraphs ($k=6$) shows the potential of the exact subgraph approach.

Summarizing, we offer the following conclusions from these preliminary computational results.

$\bullet$ Our computational approach based on the partial Lagrangian dual is very efficient in handling also a large number of ESCs. The dual function evaluation separates the SDP part from the ESCs and therefore opens the way for large-scale computations. The minimization of the dual function is carried out as a convex quadratic optimization problem without any SDP constraints, and therefore is also suitable for a large number of ESCs.

$\bullet$ On the practical side we consider the small ESCs for Max-Cut a promising new way to tighten bounds for this problem. It will be a promising new project to explore these bounds also in a Branch-and-Bound setting.

$\bullet$ Our computational results for stable set and coloring confirm the theoretical hardness results for these problems. Here the improvement of the relaxations is small for $k\leqslant 3$ but including larger subgraphs yields a noticeable improvement of the bounds. It will be a challenge to extend our approach to larger subgraphs.

Appendix 0.A Tables