Integer programming techniques for Polynomial
Optimization
Gonzalo Muñoz

Submitted in partial fulfillment of therequirements for the degreeof Doctor of Philosophy
in the Graduate School of Arts and Sciences

COLUMBIA UNIVERSITY2017c 2017
Gonzalo Muñoz
All Rights Reserved

ABSTRACT
Integer programming techniques for Polynomial
Optimization
Gonzalo Muñoz
Modern problems arising in many domains are driving a need for more capable, state-ofthe-art optimization tools. A sharp focus on performance and accuracy has appeared, forexample, in science and engineering applications. In particular, we have seen a growth instudies related to Polynomial Optimization: a field with beautiful and deep theory, offering
flexibility for modeling and high impact in diverse areas.
The understanding of structural aspects of the feasible sets in Polynomial Optimization,
mainly studied in Real Algebraic Geometry, has a long tradition in Mathematics and it has
recently acquired increased computational maturity, opening the gate for new Optimization
methodologies to be developed. The celebrated hierarchies due to Lasserre, for example,
emerged as good algorithmic templates. They allow the representation of semi-algebraic
sets, possibly non-convex, through convex sets in lifted spaces, thus enabling the use of longstudied Convex Optimization methods. Nonetheless, there are some computational drawbacks for these approaches: they often rely on possibly large semidefinite programs, and due
to scalability and numerical issues associated with SDPs, alternatives and complements are
arising.
In this dissertation, we will explore theoretical and practical Integer-Programming-based
techniques for Polynomial Optimization problems. We first present a Linear Programming
relaxation for the AC-OPF problem in Power Systems, a non-convex quadratic problem, andshow how such relaxation can be used to develop a tractable MIP-based algorithm for the AC
Transmission Switching problem. From a more theoretical perspective, and motivated by the
AC-OPF problem, we study how sparsity can be exploited as a tool for analysis of the fundamental complexity of a Polynomial Optimization problem, by showing LP formulations thatcan efficiently approximate sparse polynomial problems. Finally, we show a computationallypractical approach for constructing strong LP approximations on-the-fly, using cutting plane
approaches. We will show two different frameworks that can generate cutting planes, which
are based on classical methods used in Mixed-Integer Programming.
Our methods mainly rely on the maturity of current MIP technology; we believe thesecontributions are important for the development of manageable approaches to general Polynomial Optimization problems.

Table of Contents
List of Figuresv

List of Tablesvi

Acknowledgments

vii

1 Introduction

1

1.1

The main motivation: AC-OPF problem . . . . . . . . . . . . . . . . . . . . .

1

1.2

Network Polynomial Optimization Problems . . . . . . . . . . . . . . . . . . .

3

1.3

Structured Sparsity in General Polynomial Optimization . . . . . . . . . . . .

4

1.4

Cutting Plane Approaches to Polynomial Optimization . . . . . . . . . . . . .

7

1.5

Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8

1.6

Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

102 AC Optimal Power Flow Problem

12

2.1

Modeling Power Flows on the Grid . . . . . . . . . . . . . . . . . . . . . . . .

12

2.2

The AC-OPF Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

14

2.3

Description of our approach . . . . . . . . . . . . . . . . . . . . . . . . . . . .

16

2.4

Valid inequalities for AC-OPF . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

2.4.1

y sh = 0 and N = 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

2.4.2

General y sh but N = 1 . . . . . . . . . . . . . . . . . . . . . . . . . .

23

2.4.3

General y sh and N. . . . . . . . . . . . . . . . . . . . . . . . . . . .

25

2.5

Tightening inequalities through reference angle fixings . . . . . . . . . . . . .

27

2.6

Computational experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

i

2.7

AC transmission switching problem . . . . . . . . . . . . . . . . . . . . . . . .

30

2.8

Formulation and algorithm for ACTS . . . . . . . . . . . . . . . . . . . . . . .

31

2.8.1

Non-standard branching for ACTS . . . . . . . . . . . . . . . . . . . .

33

Computational experiments for ACTS . . . . . . . . . . . . . . . . . . . . . .

34

2.9

3 Binary Optimization with small tree-width

36

3.1

A brief tutorial on tree-width . . . . . . . . . . . . . . . . . . . . . . . . . . .

36

3.2

Problem description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

39

3.3

Examples of GB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

41

3.4

Reduction to the linear case . . . . . . . . . . . . . . . . . . . . . . . . . . . .

43

3.5

Preliminary definitions for LP reformulations . . . . . . . . . . . . . . . . . .

45

3.6

Lovász-Schrijver-based reformulation . . . . . . . . . . . . . . . . . . . . . . .

46

3.6.1

Correctness of formulation LPz . . . . . . . . . . . . . . . . . . . . . .

47

Alternative reformulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

48

3.7.1

Correctness of formulation LP-GB . . . . . . . . . . . . . . . . . . . .

49

3.7.2

Proof of Theorem 3.7.6 . . . . . . . . . . . . . . . . . . . . . . . . . .

52

3.7

4 Mixed-Integer Polynomial Optimization with small tree-width

55

4.1

Problem description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

55

4.2

Binary approximations of polynomial optimization problems . . . . . . . . . .

57

4.3

Linear reformulation of the binary approximation . . . . . . . . . . . . . . . .

61

4.3.1

Sparsity of the approximation . . . . . . . . . . . . . . . . . . . . . . .

61

4.3.2

From sparse PO to small LP approximations . . . . . . . . . . . . . .

62

Final comments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

63

4.4.1

Can the dependence on  be improved upon? . . . . . . . . . . . . . .

63

4.4.2

Example of LP approximation to PO . . . . . . . . . . . . . . . . . . .

63

4.4

5 Network Polynomial Optimization

67

5.1

Problem description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

67

5.2

Impact of node-degrees in an NPO . . . . . . . . . . . . . . . . . . . . . . . .

70

5.3

Example of NPO reformulation . . . . . . . . . . . . . . . . . . . . . . . . .

72

ii

5.4

NPO reformulation: general case . . . . . . . . . . . . . . . . . . . . . . . . .

75

5.4.1

Transforming the network and reformulating . . . . . . . . . . . . . .

76

5.4.2

Validity of the reformulation . . . . . . . . . . . . . . . . . . . . . . .

79

5.4.3

Implications to Theorem 5.1.3 . . . . . . . . . . . . . . . . . . . . . . .

82

5.4.4

Constructing reformulations on small tree-width networks . . . . . . .

83

6 Cutting planes for Polynomial Optimization
6.1

6.2

87

Digitization Cuts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

88

6.1.1

Approximations of Polynomial Optimization revisited . . . . . . . . .

88

6.1.2

Why we need to approximate . . . . . . . . . . . . . . . . . . . . . . .

94

6.1.3

Digitization-based Cuts . . . . . . . . . . . . . . . . . . . . . . . . . .

95

6.1.4

Computational Experiments . . . . . . . . . . . . . . . . . . . . . . . .

99

Intersection Cuts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
6.2.1

Review of S-free sets and Intersection cuts

. . . . . . . . . . . . . . . 101

6.2.2

S-free sets for Polynomial Optimization . . . . . . . . . . . . . . . . . 104

6.2.3

Cutting plane procedure . . . . . . . . . . . . . . . . . . . . . . . . . . 108

6.2.4

Computational Experiments . . . . . . . . . . . . . . . . . . . . . . . . 112

7 Conclusions

117

Bibliography

119

Appendices130

A Extra details for Chapter 2131

A.1 Derivation of General Case . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131
A.2 Proof of Theorem 2.4.4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135
B Extra details for Chapter 3138

B.1 Proof of part (b) of Theorem 3.2.3 . . . . . . . . . . . . . . . . . . . . . . . . 138
C Extra details for Chapter 4

140

C.1 Scaled feasiblity in Theorem 4.1.2 . . . . . . . . . . . . . . . . . . . . . . . . . 140

iii

C.2 LP size dependency on  of Theorem 4.1.2 . . . . . . . . . . . . . . . . . . . . 141
D Supplementary Experiments on Intersection Cuts

iv

143

List of Figures2.1

π-model of a line {k, m}, including transformer on the k side, and shunt admittance. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

14

3.1

Square grid formed by k colums and k rows of nodes. . . . . . . . . . . . . . .

38

3.2

Intersection Graph for system 3.2 . . . . . . . . . . . . . . . . . . . . . . . . .

40

5.1

Star network witn n + 1 nodes. . . . . . . . . . . . . . . . . . . . . . . . . . .

68

5.2

Binary tree replacement for star with 4 leaves. . . . . . . . . . . . . . . . . .

73

5.3

(a) Intersection graph for reformulation of Ex. 5.1.1. (b) A tree-decomposition. 74

5.4

Node splitting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

76

5.5

Complete node splitting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

77

5.6

Incorrect vertex splitting. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

84

5.7

Correct vertex splitting of graph in Figure 5.6(a). . . . . . . . . . . . . . . . .

84

6.1

Graphical representation of Digitization and the corresponding convex hull .

95

6.2

Graphical representation of Digitization Cut . . . . . . . . . . . . . . . . . . .

96v

List of Tables
2.1

Comparison of different relaxations of AC-OPF. . . . . . . . . . . . . . . . . .

30

2.2

ACTS algorithm on case2383wp . . . . . . . . . . . . . . . . . . . . . . . . . .

34

2.3

ACTS algorithm on case2746wp with ≥ 5 switched lines . . . . . . . . . . . .

35

6.1

Digitization Cuts Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . 100

6.2

Proposed Intersection Cuts . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109

6.3

Averages for GLOBALLib instances . . . . . . . . . . . . . . . . . . . . . . . 115

6.4

Averages for BoxQP instances . . . . . . . . . . . . . . . . . . . . . . . . . . . 115

6.5

Distribution of Closed Gaps for GLOBALLib and BoxQP . . . . . . . . . . . 115

6.6

Distribution of End Gaps for GLOBALLib and BoxQP . . . . . . . . . . . . . 116

D.1 Detailed results for 2x2 + OA on GLOBALLib instances . . . . . . . . . . . . 146
D.2 Detailed results for 2x2 + OA cuts on BoxQP instances . . . . . . . . . . . . . . . 151
D.3 Comparison with V2 on GLOBALLib instances . . . . . . . . . . . . . . . . . 152
D.4 Comparison with V2 on BoxQP instances . . . . . . . . . . . . . . . . . . . . 155vi

Acknowledgments
First of all, I would like to thank my advisor, Professor Daniel Bienstock. I came to Columbia
hoping to work with him, and everything turned out the best possible way. I encounteredan extremely smart and dedicated person, an important figure in development of the Optimization field, and an encouraging and friendly advisor. I feel privileged to have been able
to share so much with him. Working with Dan these 5 years has been a truly enjoyable
experience; he always kept a wise balance between hard work and having fun while doing
research.
I also want to thank the IEOR department for providing such a warm environment. The
faculty members were always available to give us a hand and support our research. They
constantly spread their enthusiasm and helped us mature as researchers. As for the IEOR
staff (past and current), they continuously provided help with all our extra-academic issues,
they pushed the social side of our PhD student life, and they did everything with a contagious
smile.
To my fellow PhD students: I am extremely grateful for the good friends I made here
at IEOR. To all of those who were interested in creating a strong bond, not constrained tothe boundaries of the Mudd building, thank you very much. A special mention for my friend
Antoine: I was very lucky that in such a small cohort I was able to find such a good friend.
I was lucky to see him form a beautiful family, lucky to see him start a successful academic
career and lucky to have been part of this process.
I also want to express my deepest gratitude to my family. They are my most faithful
supporters, and they have always made me believe I could do whatever I set my mind to. I
want to thank them for being next to me this whole time, even if we were 10.000 kms apart.

vii

I would like to thank my old friends from Chile as well, for confirming that time and distancedoes not wear out good friendships. And to all the new friends I made in New York, I would
like to thank them for sharing with me during this part of my life. I feel immense gratitude
to everyone who crossed paths with me at some point in this incredible city.
Last, but most importantly, to Fernanda for being next to me the whole time. Having
a partner capable of everything was the biggest support I had during these last years. She
is one of the strongest and most resourceful woman I know, and being by her side filled thisadventure with pure joy.

viii

CHAPTER 1. INTRODUCTION

1

Chapter 1

Introduction
1.1

The main motivation: AC-OPF problem

The main motivation for this work, which involves the development of efficient and accurate
Polynomial Optimization techniques, is drawn from a Power System problem: the so-called
AC-OPF problem.
The AC-OPF problem was introduced in [33], and it is a fundamental software component
in the operation of electrical power transmission systems. According to the Federal Energy
Regulatory Commission1 :Finding a good solution technique for the AC-OPF could potentially save tens ofbillions of dollars annually.
The goal of this problem is to determine where and how much power should be generated
in order to satisfy given demands in the transmission grid, in a way that optimizes generation
costs. The AC-OPF problem is solved every day, sometimes in time intervals of not morethan 10 minutes, thus efficient solution techniques are needed. Moreover, given its importantrole, solutions should be as accurate as possible.
The problem can be formulated as a non-convex, continuous optimization problem. Nonconvexities arise naturally from the physics laws governing the power flows in the transmission
grid. Mainly Ohm’s Law and Kirchhoff’s current law. Any reliable model or solution method
1

https://www.ferc.gov/industries/electric/indus-act/market-planning/opf-papers.asp

CHAPTER 1. INTRODUCTION

2accounts for these somehow, as they are not a simple modeling choice, but a physical phenomenon.
Due to its intrinsic non-convex nature, the problem is NP-Hard (see [107, 74]). In routine problem instances, solutions of excellent quality can be quickly obtained using a variety
of methodologies, including sequential linearization and interior point methods (e.g. MATPOWER [111]). Instances involving grids under stress or extreme conditions, however, canprove significantly more difficult, as the quality of the usual approximations and simplifications becomes modest. This, along with the recent advances in computation technology andoptimization methodologies, has motivated considerable efforts in solving the problem moreaccurately and efficiently. For more background on this problem, see [16].

Recently, in [73], a semidefinite programming (SDP) approach is explored, where relaxations of this type are used to obtain valid bounds on the AC-OPF problem, and thus certify
when a given feasible solution is close to optimal. These relaxations can provide good bounds,
nevertheless, SDP solvers in general can be numerically unstable and slow. It has also beenargued that, in general, the semidefinite approach can be inexact [62], even in small networks.
Other approaches typically involve convex quadratic relaxations to provide good bounds, relying in more mature and stable optimization routines. See [73, 86, 59, 80, 87, 57, 28, 93, 61]
and the citations therein for these and other cutting edge approaches.

In this dissertation we will present fast methods to obtain good bounds without having torely on SDP tools. We will exploit linear programming tools as much as possible, providing
elaborate linear inequalities that can yield valid and good quality bounds in short time. Thesecan be useful on their own, or to work as a complement of other state-of-the-art relaxations.
In addition, we will show these tools can be used to tackle the AC transmission switchingproblem, an extension of the AC-OPF problem where a planner can decide to turn off somelines of the grid.

CHAPTER 1. INTRODUCTION

1.2

3

Network Polynomial Optimization Problems

Building on the aforementioned motivation, the following question arises naturally: to what
extent can one exploit linear programming (LP) techniques in non-convex problems?. Since,
in general, the AC-OPF problem is strongly NP-Hard [107], for the question to have a nontrivial answer further characteristics need to be considered. In this case, we will rely on
the fact that transmission grids are typically sparse networks in a very specific way; theypresent a tree-like structure (see [81]). More precisely, they are graphs of low tree-width, agraph-theoretical parameter that measures, roughly, how tree-like a given graph is. We will
study what is the effect of this structure in the tractability of the problem.
In a more abstract level, we will introduce the concept of Network Polynomial Optimization (NPO) problems. These are optimization problems defined over networks where thedecision variables are associated to the nodes of the network. Additionally, there are “flows”
associated to each edge, which depend on the variables of the end-nodes. The expressions
defining these flows will be allowed to be non-linear. In the AC-OPF problem, for example,
the network is given by the transmission grid itself, the decision variables represent voltages
in each node of the network, and the physics laws governing the power flows yield non-linear
expressions for the power flowing through any line.
Mathematically, NPOs are defined as follows:• There is a network G, with vertices V (G) and edges E(G).
• For each u ∈ V (G) there is a set of variables Xu and a set Ku of constraints associatedwith u.
And the optimization problem is of the following form:

subject to:min cT xX(k)p(u,v) (Xu ∪ Xv ) ≥ 0,k ∈ Ku ,u ∈ V (G){u,v}∈δG (u)x ∈ {0, 1}p × [0, 1]n−p ,(k)where δG (u) is the set of edges of G incident with u and each p(u,v) is a polynomial that only
involves variables in Xu ∪ Xv .

CHAPTER 1. INTRODUCTION

4This class of problems not only generalizes the AC-OPF problem, but also other important
optimization problems. Standard network flow problems fit in category [1], capacitated fixedcharge network flow models [56], the unit-commitment problem [91] and some optimization
problems on gas networks [60]. The well-studied pooling problem [84] can also be cast as an
NPO.
In this dissertation we will study the effect of low tree-width of the underlying network G
defining an NPO. We will provide a way of obtaining LP approximations to it that can achieve
any desired tolerance, moreover, if G possesses bounded tree-width, our LP approximation will
be of moderate size. As a consequence, we will derive polynomial-size linear programs thatapproximate accurately the AC-OPF problem on networks that have bounded tree-width.

1.3

Structured Sparsity in General Polynomial Optimization

We will also study sparse mixed-integer polynomial optimization problems, that are not necessarily defined using a network as before. Here, sparsity will be measured using the tree-widthparameter as well, however, the graph we assume has low tree-width will be constructed froma given formulation. More precisely, we will study problems of the form
minsubject to:cT x
fi (x) ≥ 0

1≤i≤mx ∈ {0, 1}p × [0, 1]n−p .
where each fi is a polynomial. These are general Polynomial Optimization (PO) problems.
We will develop a reformulation operator which relies on the combinatorial structure of theconstraints to produce linear programming approximations which attain provable bounds. In
order to include the graph-theoretical parameter tree-width in this optimization context, we
will use the concept of intersection graph: a graph with vertices corresponding to variables,
and that has an edge whenever two variables appear together in a constraint. We will showthat a polynomial optimization problem whose intersection graph possesses low tree-width
has a “small” extended linear formulation that approximates it. We will also argue on the
difference between assuming sparsity of the intersection graph of a PO problem and assuming

CHAPTER 1. INTRODUCTION

5sparsity in the underlying network that defines an NPO.

Polynomial optimization has recently seen a considerable growth in studies related to it.
This is a long-studied field with beautiful theory, offering an appealing flexibility for modeling, and that recently has acquired much more maturity on the computational aspect. The
celebrated hierarchies by Lasserre [66], for example, emerged as a compelling approach forsolving polynomial optimization problems. Nevertheless, since they rely on SDP, scalability
issues and numerical instability associated to current SDP solvers are frequent. For thesereasons, some alternatives based on Linear Programming, or Second-order Cone Programming are arising.

From the tree-width-based sparsity perspective, there is broad literature dating from the
1980s on polynomial-time algorithms for combinatorial problems on graphs with bounded
tree-width. An early reference is [6]. Also see [4, 5, 30, 17, 24, 19] and from a very generalperspective, [27]. These algorithms rely on “nonserial dynamic programming”, i.e., dynamicprogramming on trees. See [3], [88], [18].
A parallel research stream concerns “constraint satisfaction problems”, or CSPs. One
can obtain efficient algorithms for CSPs, whenever the constraints present a sparse pattern
given by low tree-width. These algorithms rely on similar dynamic programming ideas as
the algorithms above, from the perspective of belief propagation on an appropriately defined
graphical model. Another central technique is the tree-junction theorem of [72], which shows
how a set of marginal probability distributions on the edges of a hypertree can be extended
to a joint distribution over the entire vertex set. Early references are [92, 47, 41]. Also see[108, 34, 109] (and references therein).

In the integer programming context, extended formulations for binary linear programs
whose constraints present intersection graphs of small tree-width have been developed in[20, 108, 70]. A different use of tree-width in integer programming is given in [38]. Analternative perspective on structural sparsity in optimization is taken in [29].
In this context, [21] (also see the PhD thesis [112]) develop extended formulations forCHAPTER 1. INTRODUCTION

6binary linear programs by considering the subset algebra of feasible solutions for individual
constraints or small groups of constraints; this entails a refinement of the cone of set-functions
approach of [77]. The method in [21] is similar to the one used here, in that here we rely ona similar algebra and on extended reformulations for 0/1 integer programs. The classical examples in this vein are the reformulation-linearization technique of [99], the cones of matrices
method [77], the lift-and-project method of [10], and the moment relaxation methodology of[66]. See [69] for a unifying analysis; another comparison is provided in [8].
In [108], binary polynomial optimization problems are considered, i.e problems asmin{ cT x : x ∈ {0, 1}n , fi (x) ≥ 0, 1 ≤ i ≤ m}where each fi (x) is a polynomial. They show that if the tree-width of the intersection graphof the constraints is ≤ ω, then the level-ω Sherali-Adams or Lasserre reformulation of the
problem is exact. Hence there is an LP formulation with O(nω+2 ) variables and O(nω+2 m)
constraints.
A comprehensive survey of results on polynomial optimization and related topics is provided in [70]. Section 8 of [70] builds on the work in [69], which provides a common framework
for the Sherali-Adams, Lovász-Schrijver and Lasserre reformulation operators. In addition to
the aforementioned results, [70] explicitly shows that the special case of the vertex-packing
problem on a graph with n vertices and tree-width ≤ ω has a formulation of size O(2ω n);
this is stronger than the implication from [20] discussed above. Similarly, it is shown in [70]that the max-cut problem on a graph with n vertices and tree-width ≤ ω has a formulationof size O(2ω n).

In the continuous variable polynomial optimization setting, [63, 110] present methods for
exploiting low tree-width of the intersection graph e.g. to speed-up the sum-of-squares or
moment relaxations of a problem. Also see [53] and Section 8 of [70]. [65] shows that whereω is the tree-width of the intersection graph of a polynomial optimization problem, there isa hierarchy of semidefinite relaxations where the rth relaxation (r = 1, 2, . . . ...) has O(nω 2r )variables and O(n + m) LMI constraints; further, as r → +∞ the value of the relaxation
converges to the optimum. Also see [89, 68].

CHAPTER 1. INTRODUCTION

7

Finally, there are a number of results on using lifted formulations for polynomial optimization problems, along the lines of the RLT methodology of [99]. See [101, 100] and references
therein.

1.4

Cutting Plane Approaches to Polynomial Optimization

Besides showing the existence of moderate-sized LPs that can approximate sparse polynomial
problems, and with a more pragmatic focus, we also aim to the development of computationally tractable techniques for Polynomial Optimization problems that can have an empiricallygood performance. For this purpose, we will make use of cutting planes algorithms, widelyused in Mixed-Integer Programming, in a Polynomial Optimization context.
Cutting planes applied to non-convex problems are typically derived using problemspecific structures; either using particular assumptions on the data that drives a given optimization model, or tackling single non-linear terms or single constraints separately. In thiswork we pursue the development of general cutting plane techniques that rely on minimal
assumptions on the problem structure.

We will explore two different families of cutting planes applicable to Polynomial Optimization: Digitization cuts and Intersection cuts. Digitization cuts make use of a discretization
technique that can be traced back to [51]; also see [22, 39, 54]. Using this discretization,
polynomials can be approximated accurately with a linear expression over binary variables,
thus allowing MIP technology to be of use in generating cuts. We will present theoreticalarguments for the validity of the inequalities and the quality of the cuts, and show heuristics
that can speed-up the cut-finding algorithms.
As for Intersection cuts, we will make use of a reformulation of Polynomial Optimization
problems that represents the feasible solutions (semi-algebraic sets) as P ∩ S, where P is apolyhedron and S is a closed set. These cuts will be generated from convex forbidden zones,or S-free sets, according to the Intersection cuts introduced by Balas [9]. We review the work
developed in [23], showing different families of S-free sets where cuts can be derived from.
We will also provide details on the computational efficacy of these cuts.

CHAPTER 1. INTRODUCTION

8There is a considerable amount of literature concerning cuts and strong linear relaxations
that can be applicable to Polynomial Optimization problems. These cuts typically tackle
single non-linear terms or single constraints, making use of substructures like edge-concavity[85, 103], multilinearity [79, 12, 94], or other special characteristics to derive convex envelopes
and cuts. Also see [76, 105, 104]. All cuts we will present can account for several (or all)variables and functions in the problem simultaneously. To the best of our knowledge, twopapers are similar in this regard. The disjunctive cuts of Saxena, Bonami, and Lee [97, 98]
and the lift-and-project-based cuts, using moment relaxations, proposed Ghaddar, Vera, andAnjos [50].

1.5

Overview

This dissertation is organized as follows. In Chapter 2 we discuss the AC-OPF problem andour contributions towards fast bounding procedures. We begin by showing the commonlyused models for Power Flows in Section 2.1, which are used in AC-OPF problem formulation
presented in Section 2.2. In Section 2.3 we describe our approach in general terms as an
introduction to Section 2.4, where we present in full detail the inequalities we developed forthis problem. A tightening procedure of these inequalities is presented in Section 2.5, and in
Section 2.6 we present experiments testing the performance of them.
In Section 2.7 we describe the AC transmission switching problem (ACTS), to then specify our proposed algorithm for it in Section 2.8, which makes use of our linear relaxation of
AC-OPF. Finally, in Section 2.9 we provide computational experiments on the ACTS problem.

In the subsequent chapters, and motivated by the sparse structure present in transmission
grids, we will provide a thorough study on tree-width-based sparsity in optimization problems.
In Chapter 3 we discuss a general class of optimization problems over binary variables andanalyze the effect of a low tree-width intersection graph. This will serve as the basic building
block for all subsequent results. Here, we show how to obtain a reformulation of these binary

CHAPTER 1. INTRODUCTION

9problems, where the size of the reformulation will be parametrized by the tree-width of thecorresponding intersection graph. In Section 3.1 we provide a brief tutorial on the tree-widthconcept, along with the key results we will use. In Sections 3.2 through 3.5 we provide
a preliminary discussion on the results we will prove, the necessary definitions and theirconsequences. Finally, in Sections 3.6 and 3.7 we provide two different LP reformulations
that attain the desired effect.
Chapter 4 concerns Mixed-integer Polynomial Optimization problems that present anintersection graph with small tree-width. In Section 4.1 we state the family of problems wewill tackle, as well as the main results of the chapter. In Section 4.2 we present the digitization
technique that will allow us to approximate polynomial problems with pure binary problems.
In Section 4.3 we show how to reformulate this pure binary approximation with an LP,
analyzing what is the effect of low tree-width in this process. We conclude the chapter with
Section 4.4, where we further discuss some aspects of our results and give a full detailed
example of the LP approximation we obtain on a concrete problem.
In Chapter 5 we move to approximability results for Network Polynomial Optimization
problems, which serve as a generalization of the AC-OPF problem since they have an underlying graph in the problem description, as opposed to the problems studied in Chapter
4, where a graph is drawn from a formulation. In Section 5.1 we formally define Network
Polynomial Optimization problems and argue on the difference between assuming a sparse
structure in the underlying graph and assuming the same structure on the intersection graphof a formulation. In Section 5.2 we outline the technique we will rely on to obtain tractability
for these problems. We provide a complete example of this proof technique in Section 5.3,
and then in Section 5.4 we provide the proof of the general case.

On a more pragmatic spirit, in Chapter 6 we discuss the two families of cutting planes
mentioned above, meant for generating strong LP relaxations of Polynomial Optimization in
a computationally effective manner. The first family, Digitization cuts, is presented in Section
6.1. These are cuts based on a digitization technique used for theoretical purposes in Chapter
4, and used here for generating cuts efficiently. The second family, Intersection cuts, whichwere were introduced in [23] in the polynomial optimization context, is discussed in Section

CHAPTER 1. INTRODUCTION

10

6.2. We review the key concepts behind them and show details on the computational aspects.

Finally, in Chapter 7 we present our concluding remarks, providing a general view on
the proposed techniques of this dissertation and their impact, as well as future extensions of
them.

1.6

Notation

We will denote as R and Z the set of real and integer numbers, respectively. We denote as
Rn and Zn the sets of n-dimensional vectors with coordinates in R and Z, respectively. We
n×n
use Sn×n to denote the space of symmetric n × n matrices, and S+
for the symmetric n × nmatrices that are positive semidefinite. In some occasions we will use [n] to denote the set{1, . . . , n}.
Given a set S, we let int(S) be its interior, bd(S) its boundary, conv(S) its convex hull,clconv(S) the closure of its convex hull, cone(S) its conic hull and clcone(S) the closure ofits conic hull. We use proj(S) to denote the projection of S onto a lower dimensional space,which will be specified using a subscript. For example, if V is a subspace, then projV (S)will be the projection of S onto V. In the case we distinguish different coordinates in S, forexample, if we use (x, y) to refer to vectors in S, we use projx (S) to denote the projection of
S onto the x coordinates.
We use | · | to denote the absolute value, or magnitude, of the argument. When applied
to a complex number (in the AC-OPF context), it denotes the complex number magnitude.
k · k is used to denote the euclidean norm in Rn and k · k1 for the 1-norm. The inner productbetween two vectors u and v will be denoted as uT v, where (·)T is the transpose operator.
We reserve ≤ (<) for component-wise (strict) inequality between vectors.
We let h·, ·i be the Frobenius inner product of matrices
hA, Bi = trace(AT B) =X

Aij Bij ,i,jand k·kF be the corresponding norm. We use X  0 (X  0) whenever a matrix X is positive(semi) definite, and we let rank(X) be its rank, i.e, the number of linearly independent rows

CHAPTER 1. INTRODUCTION

11(or columns) of X. X[i,j] represents the principal submatrix of X induced by indices i and j,i.e,
X[i,j] = 

Xii

Xij

Xij

Xjj


.

For a polynomial f (x), we will use the following representation. Given α ∈ Zn+ we writenxα. Y αj=xj .
j=1

Then, we represent f (x) as a sum of weighted monomials:
f (x) =

Xfα xα ,α∈I(f )

P
where each fα is rational and I(f ) ⊆ Zn+ is a finite set. We write kf k1 =α∈I(f ) |fα |. The
P
degree of f is defined as maxα∈I(f ) j αj and the support of f , denoted supp(f ), is definedas the set of variables that appear explicitly in f (x).

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

12Chapter 2

AC Optimal Power Flow Problem
In this chapter we formally describe the AC-OPF problem, a key component in the operations of the power grid. This problem can be modeled using a non-convex optimization
problem, and needs to be solved, in some form, every day, even every 10 minutes. Thus, ityields challenging optimization instances that are typically approximated in order to achieve
tractability. We present details on the formulation of the problem, along with novel linear relaxations that can provide valid strong bounds quickly. We will also discuss how this efficientrelaxation can be used on the AC Transmission Switching problem.

2.1

Modeling Power Flows on the Grid

The power grid can be represented as a network. We call the nodes of the network buses (thatcan be generators, loads or nodes where the power is redistributed) and the edges lines. Somebuses will have a demand of power that must be satisfied, some buses will be able to generate power, and the power will flow through lines. Even though this setting is reminiscent ofclassical network-flow problems, extra difficulties arise from the way current flows in the lines.

Each bus k will have a voltage associated to it, given by a complex number Vk . Twovariables ek and fk are used to represent the real and imaginary part of Vk , i.e,
Vk = ek + jfk

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM
where j =√13

−1. Alternatively, one can use a polar representation of Vk given by
Vk = |Vk |ejθkwhere θk is known as the phase angle. We will use both representations, depending on which
is better suited for arguing a given statement.

In order to represent the lines of the network, we use the traditional π-model (see Figure2.1). Consider a line {k, m} between buses k and m; its series impedance z{k,m} is defined asz{k,m} = r{k,m} + jx{k,m}where r{k,m} and x{k,m} are constants representing the line’s resistance and reactance, respectively. A line’s series admittance is given by. −1
y{k,m} = z{k,m}= g{k,m} + jb{k,m} , wherer{k,m}
x{k,m}
g{k,m} = 2
and b{k,m} = − 2.
2r{k,m} + x{k,m}r{k,m} + x2{k,m}

(2.1)
(2.2)

sh
sh
In addition, there will be a shunt admittance y{k,m}= g{k,m}+ jbsh{k,m} , and a transformer(assuming it is located at the k side of the line) with ratio.
Nk = τk ejσk ,where τk and σk are constants representing the transformer’s magnitude and phase shift
angle, respectively.
We let Ikm be the complex current injection of k in line {k, m}. Pkm will represent theactive power injected by k in line {k, m} and Qkm will represent the reactive power injectedby k in line {k, m}. Similarly, we define Imk , Pmk and Qmk . Defining
Skm = Pkm + jQkm
and using Ohm’s Law, we have
Ikm
Imk
Skm



1
1
1 sh
jσk
=
y{k,m}
Vk − e Vm + 2 y{k,m}
Vk
τk
τk
2τk


1
1 sh
= y{k,m} − jσ Vk + Vm + y{k,m}
Vm
k
τk e
2
∗
= Vk Ikm∗
Smk = Vm Imk
.

(2.3)
(2.4)
(2.5)
(2.6)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

I
k

km

V

N* I km

V

k

k

N

y= g + jb

Vm

I

mk

m

k

1 sh
_
y

1 sh
_
y

1

N:1
0

14

2

2

0

0

0

Figure 2.1: π-model of a line {k, m}, including transformer on the k side, and shunt admittance.

These are the so-called power flow equations. If we denote I the vector with components
Ikm , and V the vector with components Vk , equations (2.3) and (2.4) are typically written
in summarized form as
I = YV,
where Y is referred to as the branch admittance matrix. This matrix is defined as follows: let
Y{k,m} be the submatrix of Y given by rows km and mk, and columns k and m, thenY{k,m} 
shy{k,m}1
τk2
 y{k,m} + 2

= 


−y τ e1jσkk1−y τ e−jσk 
k

.


sh
y{k,m}
y{k,m} + 2(2.7)

All other components of Y (i.e, entries that do not appear in Y{k,m} for some line {k, m})
are 0. Now that we have defined the formulas for the power flow equations, we can move
forward in properly defining the AC-OPF problem.

2.2

The AC-OPF Problem

The AC-OPF problem consists on determining the power to be generated in the generator
buses and the appropriate voltages in each bus, in order to meet the demands of power
throughout the grid and in a way that minimizes the generating costs and satisfies engineering and operational constraints.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

15

The parameters of the model are given by• The branch admittance matrix Y defined above.
• Pkd (resp. Qdk ) the active (reactive) load, or demand, at bus k.
min the active and reactive generator output limits in bus k. If• Pkmax , Pkmin , Qmax
k , Qkbus k is not connected to a generator, then when we set
Pkmax = Pkmin = Qmax
= Qmin= 0.
kk• Vkmin and Vkmax the voltage magnitude limits in bus k.
And the variables in the model are:
• For each line {k, m} two (complex) variables associated to current Ikm and Imk , two(complex) variables associated to active power Pkm and Pmk , and two (complex) variables associated to reactive power Qkm and Qmk .
• For each bus k a (complex) variable for the voltage Vk , and additionally Pkg (resp., Qgk )the active (or reactive) generation at k (which would be fixed at zero by the limits
above, if no generator is connected to bus k).
For a bus k, we denote as δk the set of buses connected to k through a line. The AC-OPF
problem can be represented using the following form (see [73], equations (2)):

min F (P g , Qg )(2.8a)

subject to:∀ bus k : Pkmin ≤ Pkg ≤ Pkmax(2.8b)Qmin≤ Qgk ≤ Qmaxkk

(2.8c)

Vkmin ≤ kVk k ≤ Vkmax
X(Pkm + jQkm ) = (Pkg − Pkd ) + j(Qgk − Qdk )(2.8d)(2.8e)m∈δk∗∀ line {k, m} : Pkm + jQkm = Vk Ikm(2.8f)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

16∗
Pmk + jQmk = Vm Imk(2.8g)I = YV.

(2.8h)In (2.8a) F is a cost function, usually a sum of quadratics depending on the active power
generation. Constraints (2.8b) and (2.8c) set the generation limits in each bus. Constraints(2.8d) indicate line voltage magnitude limits in each bus. Constraints (2.8e) are active andreactive balance constraints, i.e, they account for the power demanded on a bus and/or
power generation limits on a bus. Constraints (2.8f)-(2.8h) capture the Power Flow equations
described before. Note that (2.8) can be formulated using only the voltage variables Vk ,nonetheless, we keep the other variables to simplify expressions.
Additionally, one can include limits on the flow over a line (given by the material properties of the line). These are typically convex constraints of the type kPkm + jQkm k ≤ U{k,m} .
We omit these to simplify the discussion, but they can easily be added to all results.

The non-convexity of this problem comes mainly from constraints (2.8f), (2.8g) and (2.8h)(also from (2.8d)), which we will see can be very complicated quadratic expressions. In thefollowing section we will describe different ways of obtaining convex relaxations to these
non-convex constraints.

2.3

Description of our approach

Here we focus on developing linear relaxations to AC-OPF problems, in lifted spaces, withthe primary goal of quickly proving lower bounds and enabling fast, standard optimization
methodologies, such as branching, to be used in this context. To motivate our approach,let (P, Q, V (2) ) be a vector that includes, for each line {k, m}, the real and reactive power
injections Pkm , Pmk , Qkm and Qmk , and for each bus k the squared bus voltage magnitude(2)|Vk |2 , denoted by Vk . Using these variables, we first write the AC-OPF (2.8) problem inthe following summarized form
min F (P, Q)
subject to:(2.9a)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

17

dL ≤ AP + BQ + CV (2) ≤ dU

(2.9b)(P, Q, V (2) ) ∈ Ω.

(2.9c)

Here,• In constraints (2.9b), A, B and C are matrices and dL and dU are vectors, all of
appropriate dimension. These constraints describe basic relationships such as generator
output limits, voltage limits, and active and reactive flow balance constraints. Clearly
constraints (2.8d) fit in this category. Constraints (2.8b), (2.8c) and (2.8e) can be
included here as well. Note that (2.8b) and (2.8e) can be combined to obtain a constraintof the form
PkL ≤

X

Pkm ≤ PkU ,m∈δkwhere PkL and PkU are given values. We can argue the same for constraints (2.8c).
• Constraints (2.9c) describe the underlying physics, e.g. Ohm’s law. In formulation (2.8)
these correspond to constraints (2.8f), (2.8g) and (2.8h). Note that in (2.9c) we do not
explicitly state the dependency on I, since we can either fully replace these variables
using (2.8h), or alternatively include them in the description of the set Ω.
• In standard OPF problem formulations, the objective F (P, Q) is typically the sumof active power generation costs (summed over the generators); a separable convex
quadratic function of the generator outputs.
Our basic approach will approximate (2.9c) with linear inequalities obtained by lifting
formulation (2.9) to a higher-dimensional space. By ‘lifting’ we mean a procedure that adds
new variables (with specific interpretations) and then writes inequalities that such variables,
together with (P, Q, V (2) ), must satisfy in any feasible solution to the AC-OPF problem.
To fix our language, we view the quantities Pkm , Pmk , Qkm , Qmk (for each line {k, m})and |Vk |2 (for each bus k) as foundational. All other variables, including those that arise
naturally from constraint (2.9c) as well as those that we introduce, will be called lifted 1 .

1Occasionally we may view the rectangular voltage coordinates as foundational.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

18

In Section 2.4 we introduce our lifted variables, as well as the inequalities that we derive
so as to obtain a convex relaxation of (2.9c). The inequalities described here will be of thefollowing types:
1. ∆ inequalities, (basic form given by (2.18))
2. (active power) loss inequalities, (basic form given by (2.22))
3. Circle inequalities, (basic form given by (2.28))
All these inequalities are convex; some linear and some conic. In the case of conic inequalities
we rely on outer approximation through tangent cutting planes so as to ultimately obtainlinear formulations as desired.

2.4

Valid inequalities for AC-OPF

In this section we derive valid inequalities, first for the simplest case (no shunt, no transformer) then for the case with shunts but no transformers, and finally for the most general
case2 .
For ease of notation, we will omit the {k, m} subscript in the line-related constants when
the dependency is evident.

2.4.1

y sh = 0 and N = 1

In this case the equations in (2.8h) take the form
Ikm = y(Vk − Vm ).

(2.10)

Using rectangular coordinates this can be expressed as
Ikm = g(ek − em ) − b(fk − fm ) + j[ b(ek − em ) + g(fk − fm ) ](2.11)

with a symmetric expression for Imk . This implies
Pkm = ek g(ek − em ) − ek b(fk − fm ) + fk g(fk − fm ) + fk b(ek − em )= (ek − em )(g , b)(efkk ) + (fk − fm )(−b , g)(efkk )2

Bus shunt admittances are omitted, but can be easily incorporated into our inequalities.

(2.12)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

19

with a symmetric expression for Pmk . Similarly,
Qkm = fk g(ek − em ) − fk b(fk − fm ) − ek g(fk − fm ) − ek b(ek − em )= (ek − em )(−b , g)(efkk ) + (fk − fm )(−g , −b)(efkk ).

(2.13)

To obtain expressions in polar coordinates we write the impedance and admittance inpolar coordinates:
z = |z|ej∠z ,y=1 −j∠ze.
|z|Then (see e.g. Bergen and Vittal [16], p. 104)∗
Skm = Vk Ikm= Vk (Vk∗ − Vm∗ )y ∗ =|Vk |2 j∠z |Vk ||Vm | j∠z jθkm
e
−
e e
,
|z|
|z|

(2.14)

where.
θkm = θk − θm .
We also can rewrite (2.14) as
Skm = |Vk |2 (g − jb) − |Vk ||Vm |(g − jb)(cos θkm + j sin θkm )= |Vk |2 g − |Vk ||Vm |g cos θkm − |Vk ||Vm |b sin θkm+j −|Vk |2 b + |Vk ||Vm |b cos θkm − |Vk ||Vm |b sin θkm .

(2.15)

Likewise, the power received at m (rather than injected), −Smk , satisfies−Smk = −|Vm |2 j∠z |Vk ||Vm | j∠z −jθkm
.
e
+
e e
|z|
|z|

(2.16)

We can also obtain an expression for Smk similar to (2.15) by switching the k and m symbols.
For these expressions in polar coordinates in the most general case, see the Appendix A.1.
2.4.1.1∆ and loss inequalities, 1

Let µkm and νkm denote known upper bounds on|(g , b)(efkk )| and |(−b , g)(efkk )|,respectively. For example, from the Cauchy-Schwarz inequality, both quantities are upperbounded by k(g, b)kVkmax , where Vkmax is an upper bound on |Vk |. Then, using (2.12) we
obtain:|Pkm | ≤ µkm |ek − em | + νkm |fk − fm |.

(2.17)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

20

Variables |ek − em | and |fk − fm | will be represented using (nonnegative) lifted variablesde,km , and df,km , thus obtaining the inequality|Pkm | ≤ µkm de,km + νkm df,km .

(2.18)

This is the basic ∆ inequality. Note that the vectors (gb ) and (−bg ) are of equal norm
and orthogonal, so further elaborations of the ∆ inequalities are possible. The upper boundk(g, b)kVkmax on µkm and νkm appears loose because, when (g, b) 6= 0 and kVk k > 0, it could
not be the case that both bounds are simultaneously tight. However,Lemma 2.4.1 The ∆ inequality|Pkm | ≤ k(g, b)kVkmax de,km + k(g, b)kVkmax df,km(2.19)

is the best possible.
Here, by “best possible” what is meant is that one can produce examples such that if thecoefficient of either de,km or df,km is tightened (i.e. decreased from k(g, b)kVkmax to a smaller
value) the resulting inequality becomes invalid – it cuts-off a feasible solution. However, theresult should not be interpreted as saying that the inequality can never be tightened. In
particular, if we consider a set of lines L and apply inequality (2.19) to each line {k, m} ∈ L,
obtaining a system of linear inequalities, it may well be the case that the system itself can betightened, that is to say, not all inequalities (2.19) for {k, m} ∈ L can simultaneously hold as
equations.

By adding the expression for Pkm in (2.15) and the corresponding expression for Pmk weobtainPkm + Pmk = g(|Vk |2 + |Vm |2 ) − 2g|Vk ||Vm | cos θkm = g|Vk − Vm |2 ,(2.20)

which can be relaxed asg(ek − em )2 + g(fk − fm )2 ≤ Pkm + Pmk ,(2.21)

or equivalently, using lifted variables,
gd2e,km + gd2f,km ≤ Pkm + Pmk ,(2.22)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

21

We term (2.22) the loss inequality. Note that by definition g ≥ 0 (unless by a modelingartifact we have r < 0), thus (2.22) is convex and a linear outer-approximation can be used
for it.
The quantity Pkm + Pmk represents the active power loss on line {k, m}. When g ≥ 0 (i.e.
r ≥ 0, the usual setting) it implies that losses are nonnegative. It is important to understandthe connection between the ∆ and the loss inequalities, which is highlighted by Theorem 2.4.4
given below, which may be stated in simplified form as “total active power generation equals
total active power loads plus total losses.” However, in order to obtain a precise statement(which is also valid in the cases where negative resistances occur) we proceed as follows.
Definition 2.4.2 Let G be an undirected graph. A pseudo-flow is a vector P that assigns toeach edge {k, m} of G two reals, Pkm and Pmk . For any node k of G we write.
δk = set of nodes of G adjacent to kand. X
ok (P ) =
Pkm .
m∈δk

We call ok (P ) the net output of k. We say that k is a source if ok (P ) > 0 and is a sink if
ok (P ) < 0. Likewise, an edge {k, m} is termed a sink-edge (or source-edge) if
Pkm + Pmk > 0

or, respectively

Pkm + Pmk < 0.

Remark 2.4.3
1. When P denotes the (standard) vector of active power flows in a transmission system,
ok (P ) = Gk −Dk where Gk is generation at k and Dk is load at k. In general, we expectthat generators will be sources and that sinks will be pure loads, though in principle one
could have a generator bus k with ok (P ) < 0.
2. In the normal case of a transmission line we will have r ≥ 0 and as discussed aboveline losses are non-negative. This means that no edge will ever be a source edge, but
could be if negative resistances are used.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

22

Theorem 2.4.4 Let G be a graph and P be a pseudo-flow on G. Then
X

X

ok (P ) =(−ok (P )) +k : ok (P )<0

k : ok (P )>0

X(Pkm + Pmk ).

{k,m}Furthermore, P can be decomposed into directed flow paths, each originating at a source orsource-edge and terminating at a sink or sink-edge.
A proof following standard flow conservation concepts is given in the Appendix A.2.

We can use this result to understand how the ∆ inequalities work. If a vector P satisfies
flow balance constraints, i.e. constraints at each bus k of the form
XPkm = Gk − Dkkmsome of the individual values Pkm must be large enough (if the active power loads are nonzero).
In the language of Definition 2.4.2, if some values ok (P ) are nonzero we will necessarily havethat some values Pkm are also nonzero. If inequality (2.18) is enforced, it will cause the liftedvariables δe,km and δf,km to be large enough. And in that case (2.22) implies that losses areappropriately large; hence by Theorem 2.4.4 total generation will (typically) have to be larger
than the sum of loads.
2.4.1.2

Circle inequalities, 1

We can rewrite equation (2.14) as
Skm = Ckm − Bkm ejθkm. |Vk |2 j∠z
Ckm =
e
and
|z|
. |Vk ||Vm | j∠z
Bkm =
e .
|z|

where

(2.23)(2.24)(2.25)

Note that Ckm and Bkm are obtained in the complex plane by rotating the real numbers|Vm |2|z|||Vm |
and |Vk|z|(respectively) by the same angle ∠z. As θkm varies, (2.14) indicates that

Skm describes a circle (the “sending circle”) with center Ckm and radius. |Vk ||Vm |
ρ=
.
|z|

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

23
2

m|
and radius
Likewise, −Smk describes a circle (the “receiving circle”) with center − |V|z|

ρ. Refer to Bergen and Vittal [16] for more details. Using either circle we can obtain valid
convex inequalities. For example, clearly we have
[Re(Skm − Ckm )]2 + [Im(Skm − Ckm )]2 ≤ ρ2 ,(2.26)

2
2

r|Vk |2
x|Vk |2|Vk |2 |Vm |2
Pkm − 2
.
+
Q
−
≤
km
2
2
2
r +x
r +x
r2 + x2(2.27)or in other words,

As discussed in Section 2.3, our formulation has variables used to represent Pkm , Qkm , |Vk |2
and |Vm |2 . Using these variables, from (2.27) we obtain a convex system by adding two liftedvariables αkm , βkm and the constraints(2)rV
Pkm − 2 k 2
r +x

= αkm

(2.28a)= βkm(2.28b)(2)xV
Qkm − 2 k 2
r +x

(2)

2
2
+ βkm
αkm

≤

(2)

Vk Vm
.
r2 + x2(2.28c)We term (2.28) the circle inequalities. Constraints (2.28a) and (2.28b) are linear, andconstraint (2.28c) is a rotated cone constraint.
Remark 2.4.5 Many other inequalities can be obtained, in particular so as to bound theratio Pkm /Qkm , using the geometry of sending- and receiving-circles.

2.4.2

General y sh but N = 1

In this case we have
Ikm = y(Vk − Vm ) +1 shy Vk ,
2(2.29)and so in rectangular coordinates
1
Ikm = g(ek − em ) − b(fk − fm ) + (g sh ek − bsh fk ) +21
j[ b(ek − em ) + g(fk − fm ) + (bsh ek + g sh fk ) ].
2

(2.30)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

24

This implies:g sh 2(e + fk2 )
2 k
bsh 2
= (ek − em )(−b , g)(efkk ) + (fk − fm )(−g , −b)(efkk ) −(e + fk2 ).
2 k

Pkm = (ek − em )(g , b)(efkk ) + (fk − fm )(−b , g)(efkk ) +(2.31)

Qkm(2.32)

Note that expressions in (2.31) and (2.32) are obtained from (2.12) and (2.13) by adding theshshterms g2 (e2k + fk2 ) and − b2 (e2k + fk2 ), respectively. In polar coordinates we will get
1∗
Skm = Vk Ikm= Vk (Vk∗ − Vm∗ )y ∗ + (g sh − jbsh )Vk Vk∗ .
2
2.4.2.1(2.33)∆ and loss inequalities, 2

Using (2.31), we obtain|Pkm | −g sh (2)V≤ µkm |ek − em | + νkm |fk − fm |.
2 k(2.34)

This is the second version of the ∆ inequality. Since the right-hand side of (2.31) isshobtained by adding g2 (e2k + fk2 ) to the right-hand side of (2.12), we have the following
analogue of (2.21):
g(ek − em )2 + g(fk − fm )2 ≤ Pkm + Pmk −g sh (2)(Vk + Vm(2) ),2(2.35)

the second version of our loss inequality, which again is a conic constraint in the space oflifted variables.
2.4.2.2

Circle inequalities, 2

From (2.33) we get
Skm

 j∠z

1 sh
|Vk ||Vm | j∠z jθkmesh
= |Vk |+ (g − jb ) −
e e
,
|z|
2
|z|
2

which again describes a circle, with center and radius, respectively,
 j∠z

1 sh
|Vk ||Vm |
2 e
sh
|Vk |
+ (g − jb )
and
.
|z|
2
|z|
Using (2.36), (2.37), and since

 j∠z



e
1
r
g sh
Re (|Vk |2+ (g sh − jbsh )= |Vk |2
+
and
|z|
2
r 2 + x2
2



 j∠z

1 sh
x
bsh
2 e
sh
2
Im (|Vk |+ + (g − jb )= |Vk |−,|z|2
r 2 + x22(2.36)(2.37)(2.38)(2.39)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

25

we obtain the following generalization of the convex lifted system (2.28):


r
g sh
(2)
Vk
= αkm
Pkm −
+
r 2 + x2
2


x
bsh
(2)
Qkm −
−
Vk
= βkm
2
2
r +x
2
(2)

2
2
αkm
+ βkm
≤

2.4.3

(2.40a)(2.40b)

(2)

Vk Vm
.
r 2 + x2

(2.40c)General y sh and N

In this case two approaches are possible.
• Indirect approach. Here we break up line {k, m} into two separate lines, i.e. line {k, k1 }and line {k1 , m} (see Figure 2.1). We have that
Vk1

=

1
1
(ek cos σ + fk sin σ) + j (fk cos σ − ek sin σ).
ττ(2.41)

which is an explicit linear inequality. Line {k1 , m} can be separately handled using theapproach in Section 2.4.2 so as to obtain ∆, loss and circle inequalities.
• Direct approach. This is the approach followed next.
From the general formula (2.5) we have (see the Appendix A.1 for details)
Pkm =i
1 h ek− em cos σ + fm sin σ (g , b)(efkk )
τ τ

1 1
+
fk − fm cos σ − em sin σ (−b , g)(efkk )τ τg sh+ 2 (e2k + fk2 ).
2τ

(2.42)

Similarly,
Pmk



1
1
= em − ek cos σ − fk sin σ (g , b)(efm
)
m
τ
τ


1
1
+ fm − fk cos σ + ek sin σ (−b , g)(efm)
mττg sh 22(e + fm).
+
2 m(2.43)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM
2.4.3.1

26∆ and loss inequalities, 3

In the transformer case there will be two ∆ inequalities. The first is obtained by from (2.42)
by taking absolute values:|Pkm | −g sh|Vk |2 ≤
2τ 2

µkm 1
νkm 1
ek − em cos σ + fm sin σ +fk − fm cos σ − em sin σ .
τ ττ τ(2.44)

Here as before µkm and νkm are known upper bounds on |(g , b)(efkk )| and |(−b , g)(efkk )|,respectively. Similarly, we obtain a second ∆ inequality from (2.43):
g sh|Vm |2 ≤
2τ 2
1
1
1
1
µmk em − ek cos σ − fk sin σ + νmk fm − fk cos σ + ek sin σ .
ττττ|Pmk | −(2.45)

Thus in order to represent these inequalities we need to introduce additional lifted variables used to model the expressions
1ek − em cos σ + fm sin στ1fk − fm cos σ − em sin σ
τ

,
,

1
1
em − ek cos σ − fk sin σ
ττ1
1
fm − fk cos σ + ek sin σ .
ττIn the no-transformer case the first two variables are equal to |em − ek | and the last two(2)

are equal to |fm − fk |. Replacing, in (2.44) and (2.45), |Vk |2 and |Vm |2 with Vk(2)

and Vmrespectively, we obtain the most general form of the ∆ inequalities.
Since all losses are incurred in the section of the line between k1 and m, applying (2.35)
and (2.41) we obtain:

2

2
1
1
1
1
g em − ek cos σ − fk sin σ + g fm − fk cos σ + ek sin σ
τ
τ
τ
τ
!
(2)
g sh Vk≤ Pkm + Pmk −+ Vm(2) .
2
τ2(2.46)In this form we obtain a convex inequality that employs the lifted variables introduced in(2.45). A similar construction yields an inequality using the lifted variables in (2.44), since

2 
2
1
1
1
1
em − ek cos σ − fk sin σ + fm − fk cos σ + ek sin σ
τ
τ
τ
τ

2 
2
1
1
=
ek − em cos σ + fm sin σ +fk − fm cos σ − em sin σ .
ττCHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

27

These inequalities represent a generalized version of the loss inequalities. For a complete
and detailed derivation of (2.46) see Appendix A.1.
2.4.3.2

Circle inequalities, 3

In the transformer case the structure of the circle inequalities differs due to the asymmetrycaused by the transformer. First, the system (2.40) applied at m is unchanged (i.e. system(2.40) with k and m interchanged). To obtain a system at k we again consider point k1 in
Figure 2.1 and we now obtain:

Pkm −

Qkm −

r
g sh
+
r2 + x2
2



x
bsh
−
r 2 + x2
2



(2)

Vk
τ2

= αkm

(2.47a)= βkm(2.47b)

(2)

Vk
τ2

(2)

2
2
≤
+ βkm
αkm

2.5

(2)

Vk Vm
.
2
τ (r2 + x2 )(2.47c)Tightening inequalities through reference angle fixings

Above we introduced a family of inequalities for each line of the underlying network. Herewe will describe a tightening procedure that can render significant improvements.
Recall the discussion in Section 2.4 regarding foundational and lifted variables. The liftedvariables include e.g. a variable used to represent em − τ1 ek cos σ − τ1 fk sin σ introduced inequation (2.45).
We can express these facts in compact form as follows. As in Section 2.4, let (P, Q, V (2) )(2)indicate the vector of all foundational variables. Here, for each bus k variable Vk

is used torepresent the quantity |Vk |2 . If N and M indicate the number of buses and lines, respectively,
then (P, Q, V 2 ) ∈ R2M +N . Let W indicate the vector of all lifted variables, say with H
components, and let K ⊆ R2M +N +H indicate the convex set described by all inequalities
introduced in Section 2.4. Then we can represent (P, Q, V (2) , W ) ∈ K more compactly bystating that.
(P, Q, V (2) ) ∈ K̂ = projR2M +N K

(2.48)

where projR2M +N K is the projection of K to the subspace of the first 2M + N variables.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

28

We now describe a procedure for tightening (2.48). As is well known, fixing an arbitrary
bus at an arbitrary angle does not change the set of feasible solutions to a standard AC-OPF
problem. Thus, let k̂ be a particular bus, and let θ̂k̂ be a particular angle; we can therefore
without loss of generality fix θk̂ = θ̂k̂ . How can we take advantage of this fact so as to obtainstronger constraints? Trivially, we can of course enforce fk̂ = tan θ̂k̂ ek̂ .
Moreover, consider for example the ∆ inequality (2.18) for a line k̂m (for simplicity weassume the line has zero shunt admittance and no transformer). We repeat the constraint
here for convenience:|Pk̂m | ≤ µk̂m |ek̂ − em | + νk̂m |fk̂ − fm |,(2.49)

where µk̂m and νk̂m are valid upper bounds onee|(g , b)(fk̂ )| and |(−b , g)(fk̂ )|,k̂k̂

respectively. As stated above, both b and g depend on the line but we omit the dependency
for simplicity of notation. Given that we know θk̂ = θ̂k̂ we can tighten the estimates on µk̂m
and νk̂m , thereby obtaining a tighter inequality from (2.49). We can likewise tighten manyof the inequalities introduced above.
More generally, suppose that rather than fixing θk̂ to a fixed value, we insist that it iscontained in a known set I(k̂) (in particular an interval), i.e.
θk̂ ∈ I (k̂)
As just argued we can therefore without loss of generality, tighten the valid inequalities
we described in previous section. This tightening is easiest in the case where the set is in fact
an interval.

Let K k̂, I (k̂) ⊆ R2M +N +H denote the resulting convex body, and let




.
Π k̂, I (k̂) = projR2M +N K k̂, I (k̂) .
As a consequence of the above observations, we now formally have:
Lemma 2.5.1 Suppose (P̃ , Q̃, Ṽ (2) ) is feasible for the AC-OPF problem. Then for any bus
k̂, and any set I (k̂) ,(P̃ , Q̃, Ṽ (2) ) ∈ Π k̂, I (k̂) .

(2.50)

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

29

Of course one can simply enforce (2.50) by explicitly writing down all the lifted variablesand all the constraints used to describe the set K k̂, I (k̂) . Alternatively, one can separatefrom the convex set Π k̂, I (k̂) and use such cuts as cutting planes. From this perspective,
the following result is important:Corollary 2.5.2 Suppose (P̃ , Q̃, Ṽ (2) ) is feasible for the AC-OPF problem. Then for any
family of buses ki (i ∈ F ) and sets I (ki ) we have(P̃ , Q̃, Ṽ (2) ) ∈\Π ki , I (ki ) .

(2.51)i∈F

In other words, in particular, we can separate a given vector (P̃ , Q̃, Ṽ (2) ) from sets obtained from our original family of valid inequalities by e.g. fixing one arbitrary bus to anarbitrary angle, and tightening.

2.6

Computational experiments

In the experiments reported here, we implemented the ∆, loss and circle inequalities in theirmost general form. To solve conic and linear programs, we used Gurobi 5.6.3 [58]. To solve
semidefinite programs, we used the system due to Lavaei et al. [80], which also includes aprocedure for extracting a feasible rank-one solution from the SDP. All runs were performed
on a current workstation with ample physical memory. All running times are in seconds
unless indicated.
In Table 2.1, “SDP time” is the time taken to solve the SDP relaxation of the OPF
problem, “SDP gap” is the percentage gap between the value of the SDP relaxation andthe upper bound (value of feasible solution) obtained by the SDP system. “SOCP time” and“LP time”, are, respectively, the time required to solve our conic relaxation and its first-order(outer) relaxation through a cutting-plane algorithm. “SOCP gap” and “LP gap” are thepercentage gaps relative to the SDP upper bound.
From Table 2.1 we can see the big potential of using LP technology. We are able to obtaincompetitive gaps in reduced time. In the smaller instances we are not able to improve on
the SOCP time; this behavior was expected, as commercial optimization software sometimes
follow the same approach we used for dealing with SOCP constraints (using linear approx-

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

30

SDP time

SDP gap

SOCP time

SOCP gap

LP timeLP gap

case9

1.04

0.0002 %0.05

0.7899 %0.04

0.7899 %case303.40

0.0185 %0.23

1.3808 %0.35

1.3964 %case574.23

0.0000 %0.62

0.9954 %1.41

0.9954 %case1188.73

0.0045 %0.98

1.4645 %5.12

1.4642 %case300

20.29

0.0018 %4.62

1.0585 %49.61

1.0559 %case2383wp

13 min0.6836 %2 min

3.6134 %1.63

5.6489 %case2746wp

16 min0.0375 %79.10

1.8593 %1.88

3.1235 %Table 2.1: Comparison of different relaxations of AC-OPF.

imations). For large instances we are able to obtain a much smaller running time, to theexpense of a slightly weaker gap.
The results are nevertheless encouraging, as they provide indications that LPs can still
be used in challenging non-convex settings such as the AC-OPF problem.

2.7

AC transmission switching problem

We now turn to a different power-grid-based problem, where the aforementioned relaxations
can be used with promising results. This is the so-called AC transmission switching problem.
In the AC transmission switching problem (ACTS for short) a planner seeks to switch off
transmission lines with the goal of reducing transmission cost, improving congestion, carrying
out line maintenance, or a number of other reasons. In addition, the planner may seek toenforce additional constraints on the set of switched off lines, such as only allowing a specificsubset of lines to be switched off and placing upper or lower bounds on the quantity ofswitched off lines. See [45], [55] for results and background. Network modeling relies on an
AC power flow model; thus ACTS is a nonlinear, non-convex, mixed-integer optimization
problem.
In the following sections we describe a methodology for addressing ACTS problems thatborrows from ideas in traditional mixed-integer optimization methods. In designing an algorithm for ACTS along these lines, one would start with an effective convex relaxation

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

31

for ACTS that yields good (lower) bounds. However, care needs to be exercised in designing such a relaxation so that it can be leveraged by the standard set of tools for solving
non-continuous optimization problems: disjunctions, branching, formulation tightening (e.g.
dynamically developed cutting planes) and repeated solutions. In this context, it is important to stress that when addressing truly difficult non-continuous optimization problems
one often needs iterative algorithms that require repeated solutions of progressively modified relaxations, as opposed to feeding a single, static problem formulation to a generic solver.

In what follows, we introduce an iterative method for tackling ACTS, at the core of whichis the lightweight relaxation of AC-OPF presented above. In order to handle a simplified
version of our relaxation to the AC-OPF problem for ease of notation, we will assume notransformer is present and that y sh = 0. The general case follows directly.

2.8

Formulation and algorithm for ACTS

The ACTS problem has a similar structure as AC-OPF; however for each line will we have abinary variable used to model the decision to switch off that line. We model line switching
using a binary variable s{k,m} for each line {k, m}. We include the constraints
kPkm + jQkm k ≤ M{k,m} (1 − s{k,m} )where M{k,m} is either the line’s limit U{k,m} (if present), or a large number. We also modify
the last constraint of the circle system (2.28) to state−[Vkmax ]4 − [Vkmin Vmmin ]222
s{k,m} + αkm+ βkm≤
r 2 + x2

(2)

(2)

Vk Vm
.
r 2 + x2

(2.52)

When s{k,m} = 0 (line not switched) (2.52) coincides with (2.28c). When s{k,m} = 1,
equations (2.28a)-(2.28b) and Pkm = Qkm = 0 give

2
2
αkm
+ βkm
=


(2) 2

Vk

r 2 + x2

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

32

which shows the validity of (2.52). The following is a valid relaxation for ACTS:
min F (P, Q)(2.53a)subject to:∀ bus k : PkL ≤X

Pkm ≤ PkU

(2.53b)Qkm ≤ QUk

(2.53c)

m∈δk

X

QL
k ≤

m∈δk

Vkmin

2

(2)

≤ Vk

≤ (Vkmax )2∀ line {k, m} : s{k,m} ∈ {0, 1}(2.53d)(2.53e)kPkm + jQkm k ≤ Mkm (1 − s{k,m} )(2.53f)∆ inequality (2.19), loss inequality (2.22)(2.53g)(2.28a), (2.28b), (2.52),(2.53h)side-constraints on the s{k,m} variables.

(2.53i)Constraint (2.53i) will be defined from the extra requirements of the line switching instance, e.g, it can represent the maximum number of lines that can be switched off, or aparticular set of lines from where to choose from. Note that the ∆ and loss constraints are
enforced even when line {k, m} is switched off. For convenience of the reader, we restatethese constraints:|Pkm | ≤ k(g, b)kVkmax de,km + k(g, b)kVkmax df,kmg d2e,km + g d2f,km≤

Pkm + Pmk .

(2.54)(2.55)

The critical observation is that, even though in our lifted relaxation for AC-OPF thevariables de,km and df,km modeled the rectangular voltage deviations |ek − em | and |fk − fm |,variables de,km and df,km appear only in the two inequalities listed above. Hence, subject
to satisfying (2.54) and (2.55), de,km and df,km can take arbitrary values. In particular,when s{k,m} = 1 (line switched off), setting de,km = df,km = 0 will satisfy (2.54)-(2.55). In
summary:
Theorem 2.8.1 Formulation (2.53) is a valid relaxation for ACTS, and as a consequenceits value is a lower bound on the value of ACTS.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

2.8.1

33

Non-standard branching for ACTS

Formulation (2.53) only proves a lower bound on the value of ACTS. To obtain upper bounds,
we use a non-standard branch-and-bound procedure which relies on any fast algorithm forcomputing upper bounds for AC-OPF; see e.g. [32]. The root node for our branching method
is given by formulation (2.53). A typical node will be endowed with an extension of formulation (2.53), and is processed by the following template (L is the set of all lines):

PROCESSING NODE v OF BRANCH-AND-BOUND
1. Solve the formulation at v. Let K v the set of lines switched off in the solution.
2. Run the AC-OPF upper bounding procedure with set K v switched off. If
feasible, we obtain an upper bound for ACTS as well.
3. If not, add to the branch-and-bound tree two nodes:→ The first extends the formulation at v with:X

s{k,m} ≤ |K v | − 1.

{k,m}∈K v→ The second extends the formulation at v with:X

s{k,m} = |K v |.

{k,m}∈K v

X

s{k,m} ≥ |K v | + 1.

{k,m}∈L

The two nodes added to the branch-and-bound tree in Step 3 guarantee correct enumeration;
the second one is justified by the observation that switching off a set K of lines may proveinfeasible, but a super-set of K could prove feasible. As branch-and-bound iterates it willproduce both upper bounds and rigorous lower bounds for ACTS3 . Any feasible solution(Step 2) will furnish an upper bound for ACTS and a corresponding set of lines to switch off.
3

Subject to round-off errors by the solver.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

34

H

MIP time

MIP value

Nodes

Feas. value

5

192 s

1804930.39

11868524.14

6

265 s

1804981.21

1

1868585.69

10

190 s

1805024.07

1

1868760.04

Table 2.2: ACTS algorithm on case2383wp

2.9Computational experiments for ACTS

In our numerical experiments with the above algorithm we consider ACTS with a sideconstraint (2.53i) stating that for some given value H, a minimum of H lines must be switched
off, indicating a maintenance or testing schedule. This constraint can be expressed asXs{k,m} ≥ H.

{k,m}∈L

We ran branch-and-bound until the first feasible solution was found.
Table 2.2 shows results using case2383wp, using Gurobi to solve the mixed-integer programin Step 1. These experiments required a single node of branch-and-bound to terminate. In
the Table, “MIP time” and “MIP value” indicate the time needed to solve the problem in
Step 1 and its value (i.e., K 1 ). “Nodes” is the number of iterations taken by the cutting-plane
method. “Feas. value” is the value of the AC-OPF problem reported by MATPOWER [111],
which we used to handle Step 2 of the algorithm, typically taking just a few seconds of CPU
time.
Note that case2383wp has more than 2800 lines; and thus, even in the case H = 5,complete enumeration of all subsets with exactly |H| lines is impractical.
Table 2.3 concerns case2746wp with H = 5. The column headed “Status” displays INF
when the AC-OPF problem in Step 2 of the template above is infeasible; and otherwise itshows the solution value.
In this case two nodes of branch-and-bound are needed. Some round-off error from the
MIP solver is noted in the sixth digit of the solution obtained in the second node.

CHAPTER 2. AC OPTIMAL POWER FLOW PROBLEM

35

Node 1 (the root)
MIP time

MIP value

Lines

Status114 s1601139.20

270, 1246, 1262

INF

1517, 3016
Node 2
MIP time

MIP value

Lines

Status135 s

1601085.68

3163, 3374, 3439

1631760.80

3492, 3500
Table 2.3: ACTS algorithm on case2746wp with ≥ 5 switched lines

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

36

Chapter 3

Binary Optimization with small
tree-width
In the following chapters we will focus on obtaining tractable methods to structurally sparse
optimization problems. Motivated by the AC-OPF problem discussed in the previous chapter
and the tree-like topology typically observed in transmission networks, tractability will bemeasured by tree-width; a graph theoretical parameter used to roughly measure how “treelike” a given graph is.
In this chapter we begin by introducing the tree-width concept to the reader, along with
key notions surrounding it in order to cement the foundations for this chapter. Next, we
will study a general class of binary optimization problems and prove how tree-width-based
sparsity can be used to tackle them effectively. The pure binary case will be the building
block toward the other more general results in the subsequent chapters.

3.1

A brief tutorial on tree-width

In what follows, given an undirected graph H, we will use V (H) and E(H) to denote thevertex-set and edge-set of H, respectively, and δH (u) will be the set of edges incident withvertex u.
Definition 3.1.1 Let G be an undirected graph. A tree-decomposition [95, 96] of G is apair (T, Q) where T is a tree and Q = {Qt : t ∈ V (T )} is a family of subsets of V (G) (theCHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

37

vertices of G) such that(i) For all v ∈ V (G), the set {t ∈ V (T ) : v ∈ Qt } forms a subtree Tv of T , and(ii) For each {u, v} ∈ E(G) there is a t ∈ V (T ) such that {u, v} ⊆ Qt , i.e. t ∈ Tu ∩ Tv .
(iii)

S

t∈V (T ) Qt

= V (G).

The width of the decomposition is defined as max {|Qt | : t ∈ V (T )} − 1. The tree-width of
G is the minimum width over all tree-decompositions of G.
Example 3.1.2 Let G be defined as
1

2

3

4

5

6

7

8

A valid tree-decomposition, with the sets Qt indicated inside each node of the tree, is asfollows:
12

46

4

7
24
7
25
7

23

57

5

8

The width of this tree-decomposition is 2.
In some cases we will refer to the sets Qt as “bags”, as they consists on sets of nodes ofthe original graph.

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH
1

···

2

···
..
.

..
.

..
.

1

2

3

..
.

..
.

k−1

k

38

···

k

Figure 3.1: Square grid formed by k colums and k rows of nodes.
Tree-width, roughly speaking, indicates how “tree-like” a graph is. Trees are the graphswith tree-width 1, cycles have tree-width 2 and cliques on n nodes have tree-width n − 1, tomention a few. It can be shown that a graph with tree-width ω and n vertices has O(ωn)
edges, and so graphs of small tree-width are sparse. However, not all sparse graphs havesmall tree-width. One well known example of this is the k × k grid given in Figure 3.1. Thisgraph is sparse in the usual sense, i.e, it posseses a low number of edges (they scale linearlywith the number of nodes), however its tree-width is k.
An alternative definition of tree-width the reader might find useful is the following:
Definition 3.1.3 A graph G is said to have tree-width at most ω if and only if G has achordal super-graph with clique number ω + 1.
A chordal super-graph of G is sometimes referred to as chordal completion of G.
Determining if a given graph has tree-width at most ω, with ω variable, is NP-hard [7].
For ω fixed, however, tree-width ω can be recognized in linear time [25]. In terms of finding
approximations to the tree-width of a graph, please see the recent work [26] and references
therein. We also refer the reader to [19, 24, 5, 44] for additional background.

Besides its width, another important feature of a tree-decomposition (T, Q) we will use
is the size of the tree-decomposition, given by |V (T )|. It was recently proven that, given
a graph H of width at most ω, computing a tree-decomposition (T, Q) of width ω thatminimizes |V (T )| is NP-hard in the class of graphs with tree-width at most ω [75]. However,for our purposes the following well known result will suffice.

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

39

Proposition 3.1.4 Given a graph G with tree-width at most ω, then there exists a treedecomposition (T, Q) of G of width ω such that
|V (T )| = O(|V (G)|)
Another important tree-decomposition result we will use is given in the following Remark.

Remark 3.1.5 Suppose (T, Q) is a tree-decomposition of some graph G, and let H be aconnected subgraph of G. Then the set of vertices t of T such that Qt intersects V (H) forms
a subtree of T .
Remark 3.1.5 follows directly from property (i) in Definition 3.1.1. And finally, a key
property relating cliques to the bags in a given tree-decomposition is given by:
Proposition 3.1.6 Consider a graph G and a tree-decomposition (T, Q) of G. Then forevery clique K ∈ V (G), there exists t ∈ T such thatK ⊆ Qt .
Proposition 3.1.6 is a standard result of graph theory, and not only it provides a condition
bags must satisfy, it also shows a direct lower bound for the tree-width of a graph.

We now move to the optimization context, where the graph-theoretical tools we justintroduced will be used to measure how sparse a problem is.

3.2

Problem description

We will study “general” binary problems, or GB for short, defined as follows.
(GB):
subject to:min cT x
xKi ∈ S ix ∈ {0, 1}n .

(3.1a)
1≤i≤m(3.1b)(3.1c)

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

40

2
1

6

5

3

4

Figure 3.2: Intersection Graph for system 3.2
For 1 ≤ i ≤ m, constraint i is characterized by a subset Ki ⊆ {1, . . . , n} and a set
S i ⊆ {0, 1}Ki . Set S i is implicitly given by a membership oracle, that is to say a mechanismthat upon input y ∈ {0, 1}Ki , truthfully reports whether y ∈ S i or not.
Any linear-objective, binary optimization problem whose constraints are explicitly stated
can be recast in the form GB; e.g., each set S i could be described by a system of algebraic
equations in the variables xj for j ∈ Ki . GB problems are related to classical constraint
satisfaction problems, however the terminology above will prove useful later.
The link between problems of the form GB and the sparsity structure defined in Section
3.1, i.e tree-width, will be given by the concept of intersection graph, which provides a wayof using this graph-theoretical parameter in an optimization context.
Definition 3.2.1 The intersection graph [48] for a system of constraints is the undirected
graph which has a vertex for each variable and an edge for each pair of variables that appear
in any common constraint.
Example 3.2.2 Consider the following system of constraints on variables x1 , . . . , x6 :3x21 − x2 ≥ 0,

−2x22 + x33 ≥ 1,x2 + x6 = 1,

x4 − x35 ≤ 0,

(3.2a)

x33 − x24 ≤ 2,

x1 + x4 ≤ 0,x2 + x5 ≥ 0,x25 − x24 = 0.

(3.2b)The intersection graph is shown in Figure 3.2, where vertex j represents xj for 1 ≤ j ≤ 6.
The intersection graph depicts the complexity of relationships among variables. If theintersection graph is dense then, potentially, problem GB could prove difficult. However, as
we will see in what follows, when the intersection graph presents low tree-width (hence, itCHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

41

is sparse) there exists tractable ways to approach GB. The main theorem we prove in thischapter is the following:
Theorem 3.2.3 Consider a GB problem whose intersection graph has tree-width ≤ ω.
(a) There is an exact linear programming formulation with O(2ω n) variables and constraints, with {0, 1, −1}-valued constraint coefficients.
(b) The formulation can be constructed by performing O(2ω m) oracle queries and withadditional workload Õ(ωn2ω (m + ω)), where the “Õ” notation indicates logarithmic
factors in m or n.
Note that the size of the formulation is independent of the number constraints in thegiven instance of GB. And even though we use the general setting of membership oracles,
this theorem gives an exact reformulation.
A proof of part (a) in Theorem 3.2.3 can be obtained using techniques in [70] (Section 8)
although not explicitly stated there. We will outline this proof, which relies on the “cone of
set-functions” approach of [77] and also present a new proof.
Regarding part (b) of the theorem, it can be easily seen that 2ω m is a lower bound on
the number of oracle queries that any algorithm for solving GB must perform.

Of course, Theorem 3.2.3 also implies the existence of an algorithm for solving GB in timepolynomial in (2ω , n, m). However one can also derive a direct algorithm of similar complexity
using well-known, prior ideas on polynomial-time methods for combinatorial problems ongraphs of bounded tree-width.
In the rest of this chapter we work on the context of Theorem 3.2.3. Before proving theTheorem, we begin by showing some examples for problem GB in order to emphasize its
potential and further analyze different angles of it.

3.3

Examples of GB

Example 3.3.1 (Linear binary integer programming). Let A be an m × n matrix, and consider a problem min{cT x : Ax ≥ b, x ∈ {0, 1}n }. To view this problem as a special case

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

42

of GB, we set for 1 ≤ i ≤ m, Ki = {1 ≤ j ≤ n : aij 6= 0} and S i = {x ∈ {0, 1}Ki :P
j∈Ki aij xj ≥ bi }.
In this special case, problem GB can be addressed by a variety of methods. Of particularinterest in this work are the reformulation or lifting methods of [77] and [99]. Next we consider
a more complex example, chosen to highlight the general nature of the problem.
Example 3.3.2 Let d, r, p be positive integers. Consider a constrained semidefinite program
over binary variables of the form
minsubject to:

r X
d X
d
X

k
ckij Xi,j

k=1 i=1 j=1k
M • X k = bk ,(3.3a)

1 ≤ k ≤ r,

(3.3b)

X k ∈ Sd×d
1 ≤ k ≤ r,
+ ,
X
k
Xi,j
≡ 0 mod p,
1 ≤ k ≤ r,

(3.3c)(3.3d)

i,j
k−1
k
,
Xi,1
= Xi,d1 ≤ i ≤ d,kXi,j
∈ {0, 1}, ∀i, j, k.

2 ≤ k ≤ r,

(3.3e)(3.3f)

Here M1 , . . . , Mr are symmetric d×d matrices, and b and c are vectors. Constraint (3.3e)
states that the first column of matrix X k is identical to the last column of matrix X k−1 .
We obtain an instance of problem GB with m = 2r − 1, as follows. First, for each1 ≤ k ≤ r we let Kk be the set of triples (i, j, k) with 1 ≤ i, j ≤ r, and S k to be the set ofk that satisfy (3.3b)-(3.3d). Next, for each 2 ≤ k ≤ r we let K
binary values Xi,jr+k−1 be theset of all triples (i, 1, k − 1) and all triples (i, d, k) and S r+k−1 to be the set of binary values(indexed by Kr+k−1 ) such that (3.3e) holds.
In the case of this last example, a direct application of standard integer programming
methods appears difficult. Moreover, we stress that the sets S i in problem GB are completely
generic and that the membership oracle perspective can prove useful as we discuss below.
Theorem 3.2.3 concerns the tree-width of the intersection graph of a problem of type GB.
Recall that as per Definition 3.2.1, given a problem instance I of GB, the intersection graphfor I has a vertex for each 1 ≤ j ≤ n, and an edge {j, k} whenever there exists 1 ≤ i ≤ m

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

43

such that {j, k} ⊆ Ki , that is to say, j and k appear in a common constraint in problem
GB.
Example 3.3.3 (Example 3.3.2, continued). Here the set of variables is given by{(i, j, k) : 1 ≤ k ≤ r and 1 ≤ i, j ≤ d}.
The intersection graph of the problem will have(a) the edge {(i, j, k), (i0 , j 0 , k)} for all 1 ≤ k ≤ r and 1 ≤ i, j, i0 , j 0 ≤ d, arising fromconstraints (3.3b)-(3.3d)(b) the edge {(i, 1, k), (i, d, k−1)} for each 1 ≤ k < r and 1 ≤ i ≤ d, arising from constraints(3.3e).
A tree-decomposition (T, Q) of the intersection graph, of width O(d2 ), is obtained as follows.
Define T as a path with vertices v1 , u2 , v2 , u3 , . . . , vr−1 , ur , vr . For 1 ≤ k ≤ r we set Qvk ={(i, j, k) : 1 ≤ i, j ≤ d} and for 2 ≤ k ≤ r we set Quk = {(i, 1, k), (i, d, k − 1) : 1 ≤ i ≤ d}.
Sets Qvk account for all edges of type (a), whereas the sets Quk cover all edges of type (b).
2

Thus Theorem 3.2.3 states that there is an LP formulation for problem (3.3) with O(2d d2 r)
variables and constraints.

3.4

Reduction to the linear case

Consider a problem instance of GB. An apparently simpler alternative to the general approach we follow would be to construct, for 1 ≤ i ≤ m, the polyhedron.
Pi = conv x ∈ {0, 1}Ki : x ∈ S i⊆ RKi .

Thus we can write Pi as the projection onto RKi of a polyhedron {x ∈ [0, 1]n : Ai x ≥ bi }where each row of Ai has zero entries on any column not in Ki . Thus, the GB problem canbe restated as the equivalent linear integer program
min cT x
subject to: Ai x ≥ bi ,(3.4a)1≤i≤mx ∈ {0, 1}n .

(3.4b)(3.4c)

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

44

Switching to this formulation makes it possible to apply general integer programming
methods to problem GB. However, this analysis ignores the size of formulation (3.4). In
particular, for any integer d ≥ 1 large enough there exist examples of 0/1-polytopes in Rd
with at least


d
log d

d/4

facets (up to constants). See [13], [49], [64]. Using this observation, one can constructexamples of problem GB where the tree-width of the intersection graph is ω = d − 1 andeach of the matrices Ai has more than ω ω/4 rows (see Example 3.4.1, below). This dependenceon ω makes any classical integer programming method more computationally expensive thanusing the method we will present.
Example 3.4.1 Choose d ≥ 2 large enough so that there is a 0/1-polyhedron P ⊆ Rd withmore than (cd/ log d)d/4 facets for some c. Let P be given by the system Ax ≥ b, where A
is M × d (M ≥ (cd/ log d)d/4 ). Choose N ≥ 1, and consider the system of inequalities overbinary variables xij , for 1 ≤ i ≤ N and 1 ≤ j ≤ d:
Axi ≥ b,1 ≤ i ≤ N,

(3.5a)

x1j = xij

2 ≤ i ≤ N, 1 ≤ j ≤ bd/2c.

(3.5b)xij ∈ {0, 1}∀ i, j.

(3.5c)

Constraint (3.5a) indicates that this system includes N copies of polyhedron P , with eachcopy described using a different coordinate system. Constraint (3.5b) states that the first
bd/2c coordinates take equal value across all such systems.
Any linear program over (3.5) is can be viewed as an example of problem GB with m =2N − 1; for 1 ≤ i ≤ N , Ki is used to represent the d variables xij (1 ≤ j ≤ d) and S i is acopy of the set of binary points contained in P (i.e. the extreme points of P ).
The intersection graph of this instance of GB will be the union of N cliques (one for each
set of variables xi ) plus the set of edges {x11 , xi1 } for 2 ≤ i ≤ N . A tree-decomposition (T, Q)
of this graph, of width d − 1, is as follows:• T has vertices u(0), as well as u(i) and v(i), for 1 ≤ i ≤ N .

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH45no• Let Qu(0) = x1j : 1 ≤ j ≤ bd/2c , for 1 ≤ i ≤ NQu(i) = Qu(0) ∪ xij : 1 ≤ j ≤ bd/2c
and Qv(i) = {xij , 1 ≤ j ≤ d}.
Thus, ω = d−1 and Theorem 3.2.3 states that any linear objective problem over constraints(3.5) can be solved as a continuous LP with O(2d dN ) variables and constraints. In contrast,
system (3.5) has more than (cd/ log d)d/4 N constraints.
As the example shows, formulation (3.4) may be exponentially larger than the linear program
stated in Theorem 3.2.3.

We will now proceed with the proof of Theorem 3.2.3, which we recall for convenience ofthe reader:
Theorem 3.2.3 Consider a GB problem whose intersection graph has tree-width ≤ ω.
(a) There is an exact linear programming formulation with O(2ω n) variables and constraints, with {0, 1, −1}-valued constraint coefficients.
(b) The formulation can be constructed by performing O(2ω m) oracle queries and withadditional workload Õ(ωn2ω (m + ω)), where the “Õ” notation indicates logarithmic
factors in m or n.
We will provide 2 formulations that prove part (a) of the Theorem. The first formulation
is obtained using techniques in [70]. We outline this proof in Section 3.6, which relies on the“cone of set-functions” approach of [77]. The second and new formulation is presented in
Section 3.7, and it is based on techniques presented in [21]. For part (b) please see Appendix
B.1.

3.5

Preliminary definitions for LP reformulations

For the next sections we consider a problem of the form GB and we assume that we have a tree
decomposition (T, Q) of width ω of the intersection graph of constraints (3.1b). Furthermore,
we recall and define useful notation:

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

46• As per the definition of GB, for 1 ≤ i ≤ m we let Ki be the set of variables xj that
explicitly appear in i-th constraint.
• Whenever we have values vj for all j in some superset of Ki , we will write vKi to indicate
the sub-vector of v with indices in Ki .
Definition 3.5.1 Let t ∈ V (T ).
(a) We say that v ∈ {0, 1}Qt is Qt -feasible if vKi ∈ S i for every 1 ≤ i ≤ m such thatKi ⊆ Qt .
(b) Write Ft = {v ∈ {0, 1}Qt : v is Qt -feasible}.
Roughly speaking, the set Ft collects all vectors that are partially or locally feasible. Notethat such vectors do not posses full support, and, moreover, they might not even sub-vectorsof fully feasible vectors. We will see, however, that this local information can be merged
through an LP in order to form fully feasible solutions using the tree-decomposition of theintersection graph.

3.6

Lovász-Schrijver-based reformulation

In this first formulation, which proves Theorem 3.2.3 (a), we use:
• A variable λtv , for each t ∈ V (T ) and each vector v ∈ Ft .
• A variable ZS , for each S such that S ⊆ Qt for some t ∈ T .
The formulation is as follows:(LPz) : minnX

cj Z{j}(3.6a)

j=1

s.t.

∀t ∈ V (T ) :X
ZS ={λtv : v ∈ Ft , vj = 1 ∀j ∈ S}
X
λtv = 1, λt ≥ 0.
v∈Ft

∀S ⊆ Qt(3.6b)(3.6c)

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

47

Constraints (3.6b) enforce consistency across different t ∈ V (T ) in the following sense.
When a set S is a subset of at least two different sets Qt , Qt0 then we will have two constraints(3.6b) used to define the same value ZS but with different right-hand sides; thus imposing a0

relationship between the λt and the λt variables. Additionally, when S is a singleton thenZS appears in the objective of LPz. But in all remaining cases for variable ZS (i.e., |S| =6 1
and S contained in a single set Qt ) the constraint (3.6b) is redundant because ZS does notappear anywhere else in the formulation.
Formulation LPz was built using the given tree-decomposition (T, Q) as an input. We
know that this is a tree-decomposition for a graph on n vertices and as a result, it can beassumed that T has at most O(n) vertices (see Proposition 3.1.4). Using this observation,
since each set Qt has cardinality at most ω + 1, it is easily derived that formulation LPz hasO(2ω n) variables and constraints.

3.6.1

Correctness of formulation LPz

Constraint (3.6b) can be restated in a more familiar way. Given t ∈ T , (3.6b) states:
ZS =X

supp(v,Qt )λtv ζS,∀S ∈ 2Qt .

(3.7)

v∈F (t)Here, given a set Y and a vector w ∈ RY , supp(w, Y ) = {j ∈ Y : wj 6= 0}, and for any
Yset Y and p ⊆ Y the vector ζ p ∈ {0, 1}2 is defined by setting, for each q ⊆ Y , 1 if q ⊆ p
ζqp = 0 otherwise.
Constraints (3.6b)-(3.6c) then describe the Lovász-Schrijver approach to lifted formulations, restricted to a given set Qt . It is clear that LPz amounts to a relaxation for the given
problem GB, in the sense that given x̂ feasible for GB then there is a vector (Ẑ, λ̂) feasible
for LPz where Ẑ{j} = x̂j for 1 ≤ j ≤ n. To do so, let t ∈ V (T ) and denote by x̂t the
restriction of x̂ to Qt . Then by definition we have that x̂t ∈ Ft . Thus we can set λ̂tx̂t = 1 andsupp(x̂t ,Qt )λ̂tv = 0 for any other v ∈ Ft , and for any S ⊆ Qt ẐS = ζS

. The last equation simply

states that ẐS = 1 iff S ⊆ supp(x̂, [n]), a consistent definition across t ∈ V (T ). Hence indeed(Ẑ, λ̂) is feasible for LPz and attains Ẑ{j} = x̂j for each j = 1, . . . , n, as desired. Note that,

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

48

Q

effectively, we have argued that the restriction of ζ supp(x̂,[n]) to R2 t for t ∈ V (T ) yields a
feasible solution to LPz.
Next we argue that (3.6b)-(3.6c) defines an integral polyhedron. This is a consequence ofthe following result, which can be obtained from Lemma 8.18 of [70].
Lemma 3.6.1 Suppose that (Z, λ) is a feasible solution to (3.6b)-(3.6c). Then there existsna vector W ∈ R2 , nonnegative values θ1 , . . . , θk and vectors y 1 , . . . , y k in Rn+ such that:(1) ZS = WS , for all S ∈ ∪t∈V (T ) 2Qt .
(2)Pk

i=1 θi= 1.

(3) y i is feasible for GB, for 1 ≤ i ≤ k.
(4) W =

Pk

i=1 θi ζ

supp(y i ,[n]) .
i

As a consequence of (2)-(4), the vector W is a convex combination of the vectors ζ supp(y ,[n])which as argued above yield feasible solutions to LPz, thus proving the desired result. We
remark that the proof of Lemma 8.18 of [70] is related to that of the tree-junction theorem;
this technique, evocative of dynamic programming, was also used in [20] in a closely related
setting.

3.7

Alternative reformulation

We now construct a formulation that yields Theorem 3.2.3 (a), as an alternative to the proof
outlined in [70] and developed in Section 3.6.
As before, we assume a tree-decomposition (T, Q) with width ω of the intersection graphfor formulation (3.1). We use the notation introduced in Definition 3.5.1, and additionally
Definition 3.7.1 Let t ∈ V (T ). We let Ωt denote the set of pairs (Y, N ) with Y ∩ N =∅, Y ∪ N ⊆ Qt , and such that
1. |Y | ≤ 1 and |N | = 0, or
2. (Y, N ) partition Qt ∩ Qt0 , for some t0 ∈ V (T ) with {t, t0 } ∈ E(T ).

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH49

The variables of this formulation are:
• A variable λtv , for each t ∈ V (T ) and each vector v ∈ Ft .
• A variable X[Y, N ], for each pair (Y, N ) ∈ 2n × 2n with (Y, N ) ∈ Ωt for some t ∈ V (T ).
And the formulation is as follows.
(LP-GB) : minnX

cj X[{j}, ∅](3.8a)

j=1

s.t. ∀t ∈ V (T ) :X

X[Y, N ] =

λtv

v∈Ft

X

λtv = 1,

Y
j∈Y

vj

Y

(1 − vj ) ∀ (Y, N ) ∈ Ωt

(3.8b)j∈N

λt ≥ 0.

(3.8c)v∈Ft

A similar counting argument to the one above can be used to prove that LP-GB has sizeO(2ω n).

3.7.1

Correctness of formulation LP-GB

In this section we will show that LP-GB is a relaxation of GB, that the relaxation is exact
and that the polyhedron defined by (3.8b)-(3.8c) is integral in the X[Y, N ] variables.
Remark 3.7.2(f.1) When (Y, N ) partition Qt ∩ Qt0 for some edge {t, t0 } then variable X[Y, N ] will appear
in the constraint (3.8b) arising from t and also that corresponding to t0 . This implies
0

an equation involving λt and λt .
(f.2) The sum on the right-hand side of constraint (3.8b) could be empty. This will be thecase if for any v ∈ {0, 1}Qt with vj = 1 for all j ∈ Y and vj = 0 for all j ∈ N there
exists 1 ≤ i ≤ m with Ki ⊆ Qt and yet vKi ∈/ S i . Then (3.8b) states X[Y, N ] = 0.
(f.3) When Y = N = ∅ the right-hand side of (3.8b) is
X[∅, ∅] = 1.

tv∈Ft λv .

PHence we will have

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

50(f.4) The λ variables are the same as those in formulation LPz. For any edge {t, t0 } ∈ E(T )and Y ⊆ Qt ∩ Qt0 the terms on the right-hand side of the row (3.8b) are a subsetof the terms on the right-hand side of the row (3.6b) corresponding to Y (additional
statements are possible).
First we show that LP-GB is a relaxation for GB, in a strong sense.
Lemma 3.7.3 Let x̃ be a feasible solution to an instance for GB.
(i) There is a feasible, 0/1-valued solution (X̃, λ̃) to LP-GB such that for each variable
Q
Q
X[Y, N ] in LP-GB we have X̃[Y, N ] = j∈Y x̃j j∈N (1 − x̃j ).
(ii) As a corollaryPn

j=1 cj X̃[{j}, ∅]= cT x̃.

Proof:(i) For each variable X[Y, N ] in problem LP-GB we set X̃[Y, N ] =Qj∈Y x̃jQj∈N (1− x̃j ).

Further, for each t ∈ V (T ) let ṽ(t) ∈ {0, 1}Qt be the restriction of x̃ to Qt , i.e. ṽ(t)j = x̃j
for each j ∈ Qt . Since x̃ is feasible, ṽ(t) ∈ Ft . Then we set λ̃tṽ(t) = 1 and λ̃tv = 0 for
every vector v ∈ Ft with v 6= ṽ(t). By construction for every t ∈ V (T ) and (Y, N ) ∈ Ωt
we have X̃[Y, N ] = 1 iff ṽ(t)j = 1 for all j ∈ Y and ṽ(t)j = 0 for all j ∈ N ; in otherwords (3.8b) is satisfied.
(ii) This follows from (i).

As a consequence of Lemma 3.7.3, Theorem 3.2.3 will follow if we can prove that theconstraint matrix in LP-GB defines an integral polyhedron in the X[Y, N ] variables. This
will be done in Lemma 3.7.7 given below. In what follows, we will view T as rooted, i.e. alledges are directed so that T contains a directed path from an arbitrarily chosen leaf vertexr (the root of T ) to every other vertex. If (v, u) is an edge thus directed, then we say that v
is the parent of u and u is a child of v.
Definition 3.7.4 A rooted subtree T̃ is a subtree of T , such that there exists a vertex u ofT̃ so that T̃ contains a directed path from u to every other vertex of T̃ . We then say that T̃
is rooted at u.

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

51

Definition 3.7.5 Let T̃ be a rooted subtree of T .
S(a) We denote by Ω(T̃ ) the sett∈T̃ Ωt .

(b) We denote by V(T̃ ) the set {j : j ∈ Qt for some t ∈ T̃ }.
Below we will prove the following result:
Theorem 3.7.6 Let (X̂, λ̂) be a feasible solution to the LP-GB problem. Then for every
rooted subtree T̃ there is a family of vectors
pk,T̃ ∈ {0, 1}Ω(T̃ ) ,vectors
xk,T̃ ∈ {0, 1}V(T̃ )and reals
0 < µk,T̃ ≤ 1,
(k = 1, 2, . . . , n(T̃ )) satisfying the following properties:(a) For each 1 ≤ k ≤ n(T̃ ) and each constraint 1 ≤ i ≤ m of problem GB, if Ki ⊆ Qt forsome t ∈ T̃ , then xk,T̃ ∈ S i .
(b) For 1 ≤ k ≤ n(T̃ ) and each pair (Y, N ) ∈ Ω(T̃ ),pk,T̃ [Y, N ] =

Y

T̃
xk,
j

j∈Y

Y


T̃
1 − xk,
.
j

j∈N

T̃
As a result, for each 1 ≤ k ≤ n(T̃ ) and j ∈ V(T̃ ), xk,= pk,T̃ [{j}, ∅].
j(c)Pn(T̃ )

k=1 µ

k,T̃

= 1.

(d) For each (Y, N ) ∈ Ω(T̃ ),n(T̃ )X̂[Y, N ] =Xµk,T̃ pk,T̃ [Y, N ].

k=1

The family of vectors pk,T̃ and reals µk,T̃ will be called a decomposition of (X̂, λ̂) over T̃ .

Pending a proof of Theorem 3.7.6, we can show that the polyhedron defined by theconstraints in LP-GB is integral.

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

52

Lemma 3.7.7 The polyhedron defined by (3.8b)-(3.8c) is integral in the X[Y, N ] variables
and problems GB and LP-GB have the same value.
Proof: Let (X̂, λ̂) be a feasible solution to LP-GB. We apply Theorem 3.7.6 with T̃ = T obtaining a family of vectors pk ∈ {0, 1}Ω(r) , vectors xk ∈ {0, 1}n and reals µk , for 1 ≤ k ≤ n(r),satisfying conditions (a)-(d) of the theorem. By (a) each vector xk is feasible for GB, sinceeach Ki induces a clique in the corresponding intersection graph and by Proposition 3.1.6 wecan conclude. By (d), the vector X̂ is a convex combination of the vectors pk . From this weconclude.

This result completes the proof of Theorem 3.2.3, pending Theorem 3.7.6.

3.7.2

Proof of Theorem 3.7.6

Assume we have a feasible solution (X̂, λ̂) to LP-GB. The proof of Theorem 3.7.6 will bedone by induction on the size of T̃ . First we handle the base case.
Lemma 3.7.8 If T̃ consists of a single vertex u, there is a decomposition of (X̂, λ̂) over T̃ .
Proof: We have that Ω(T̃ ) = Ωu (see Definition 3.7.5). By (3.8c) we haveuv∈Fu λ̂v = 1.

P

Let n(T̃ ) > 0 be the number of elements v ∈ Fu with λ̂uv > 0 and denote these vectors by{w(1), . . . , w(n(T̃ ))}. Then, for 1 ≤ k ≤ n(T̃ ) let xk,T̃ = w(k) and µk,T̃ = λ̂uw(k) . Finally, for1 ≤ k ≤ n(T̃ ) we define the vector pk,T̃ ∈ Ωu by setting
Y k,T̃ Y T̃
pk,T̃ [Y, N ] =xj
1 − xk,jj∈Yj∈Nfor each pair (Y, N ) ∈ Ωu . Now we will verify that conditions (a)-(d) of Theorem 3.7.6 hold.
Clearly (a)-(c) hold by construction. To see that (d) holds, note that (X̂, λ̂) satisfies (3.8b):
X̂[Y, N ] =

X

λ̂uv

v∈Fu

Y

vj

j∈Y

Y

(1 − vj )j∈Nn(T̃ )

=

X
k=1

µk,T̃

Y
j∈Y

T̃
xk,
j

Y
j∈N

n(T̃ )=

X
k=1µk,T̃ pk,T̃ [Y, N ]T̃
1 − xk,jCHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

53

which is condition (d), as desired.

Next we prove the general inductive step needed to establish Theorem 3.7.6. The technique used here is related to the junction tree theorem, is similar to one used in [20] and isreminiscent of Lemma L of [70].
Consider a vertex u of T and a subtree T̃ rooted at u with more than one vertex. Letv be a child of u. We will apply induction by partitioning T̃ into two subtrees: the subtree
L consisting of v and all its descendants in T̃ , and the subtree H = T̃ \ L. Consider adecomposition of (X̂, λ̂) over L given by the vectors pk,L ∈ {0, 1}Ω(L) and the positive realsµk,L for k = 1, 2, . . . , n(L), and a decomposition of (X̂, λ̂) over H given by the vectors
pk,H ∈ {0, 1}Ω(H) and the positive reals µk,H for k = 1, 2, . . . , n(H).
Denote by P̃ the set of partitions of Qu ∩ Qv into two sets. Thus, by Definition 3.7.1, foreach (α, β) ∈ P̃ we have a variable X[α, β]. Note that Ω(T̃ ) = Ω(H) ∪ Ω(L). We construct afamily of vectors and reals satisfying (a)-(d) Theorem 3.7.6 for T̃ , as follows.
For each (α, β) ∈ P̃ such that X̂[α, β] > 0, and each pair i, h such that 1 ≤ i ≤ n(L),
α,β
α,β
using
and a real γih
1 ≤ h ≤ n(H), and ph,H [α, β] = pi,L [α, β] = 1 we create a vector qih

the following rule. For any vertex t in T̃ and (Y, N ) ∈ Ωt :α,β[Y, N ] = pi,L [Y, N ].
(r.1) If t ∈ V (L) we set qihα,β(r.2) If t ∈ V (H) we set qih[Y, N ] = ph,H [Y, N ].

Further, we set
α,β
γih
=

µi,L µh,HX̂[α, β].

To argue that this construction is valid we note that since X̂[α, β] > 0, pairs of indicesi, h as listed above must exist, by (d) of the inductive assumption applied to H and L.
α,βFurthermore, we have γih> 0.

Now we will prove that the qih and the γih provide a decomposition of (X̂, λ̂) over T̃ . Let
i and h be given. Since the restriction of pi,L (and ph,H ) to L (resp., H) satisfy (a) and (b)of the inductive assumption, so will qih . Thus, there remains to prove (c) and (d).

CHAPTER 3. BINARY OPTIMIZATION WITH SMALL TREE-WIDTH

54

First, consider (d). Let (Y, N ) ∈ Ω(T̃ ), say (Y, N ) ∈ Ω(H). We claim thatXα,β α,βγihqih [Y, N ]α,β,i,hn(L) n(H)X X µi,L µh,HX=(α,β)∈P̃ : X̂[α,β]>0 i=1 h=1X̂[α, β]pi,L [α, β]ph,H [α, β]ph,H [Y, N ].

(3.9)

This equation holds because in any nonzero term in either expression we must have pi,L [α, β] =α,βph,H [α, β] = 1 and since (Y, N ) ∈ Ω(H) we also have that qih[Y, N ] = ph,H [Y, N ].

Now the right-hand side of (3.9) equals
n(L) i,L i,Ln(H)X
X µ p [α, β]
X


i=1

(α,β)∈P̃ : X̂[α,β]>0X̂[α, β]µh,H ph,H [α, β]ph,H [Y, N ]

h=1



n(H)



X

µh,H ph,H [α, β]ph,H [Y, N ] ,

X

=



(3.10)

h=1(α,β)∈P̃ : X̂[α,β]>0

by the inductive assumption (d) applied to subtree L. The expression in (3.10) equals



n(H)
X
X

ph,H [α, β] µh,H ph,H [Y, N ] .
(3.11)
h=1(α,β)∈P̃ : X̂[α,β]>0

By inductive property (b) applied to subtree H, given 1 ≤ h ≤ n(H) we have that ph,H [α, β] =1 for exactly one partition (α, β) ∈ P̃, and so expression (3.11) equalsn(H)Xµh,H ph,H [Y, N ].

(3.12)

h=1

In summary,n(H)Xα,β,i,hα,β α,βγihqih [Y, N ]=Xµh,H ph,H [Y, N ].

h=1

and by induction applied to the subtree H this quantity equals X̂[Y, N ]. Thus property (d)does indeed hold.

Finally we turn to (c). Inductively, (c) and (d) hold for trees L and H. Thus, as noted
in Remark 3.7.2 (f.3), X[∅, ∅] = 1. But we have just shown that (d) holds for T̃ , and inparticular that it holds for Y = N = ∅. Using Remark 3.7.2 (f.3) we obtain that (c) holds
for T̃ , as desired.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

55

Chapter 4

Mixed-Integer Polynomial
Optimization with small tree-width
In Chapter 3 we discussed how to exploit tree-width-based sparsity in pure binary optimization problems. We showed how to obtain an LP reformulation of these problems whose size
is parametrized by the tree-width of the problem.
In this chapter we will show how to use these results in order to obtain tractable approaches to mixed-integer polynomial optimization problems. Due to the inclusion of continuous variables, in this case we will obtain linear programming approximations instead ofexact reformulations. We begin by stating the setting of this chapter, as well as the main
theorem to be proved.

4.1

Problem description

Our main goal will be the study of optimization problems of the form(PO) : min cT x
subject to :fi (x) ≥ 0(4.1a)1≤i≤mx ∈ {0, 1}p × [0, 1]n−p ,(4.1b)(4.1c)

where each fi is a polynomial. Any linear-objective polynomial optimization problem where
explicit upper and lower bounds are known for all variables can be brought into this form by

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

56

appropriately translating and scaling variables, and restating any equation as two inequalities. We focus on obtaining linear programming approximations to PO that attain any
desired tolerance (both feasibility and optimality). Our goal is to obtain linear programming
approximation of polynomial size when problem PO is appropriately sparse. As in the previous chapter, we will make use of tree-width-based sparsity of the intersection graph of PO,
defined in 3.2.1.

We represent the polynomial of the i-th constraint of PO, fi (x) ≥ 0, as
fi (x) =X

fi,α xα ,(4.2)

α∈I(i)

where I(i) ⊆ Zn+ . Recall that, given α ∈ Zn+n. Y αj
xα =
xj .
j=1

and kfi k1 =P

α∈I(i) |fi,α |.

Definition 4.1.1 Given a set of polynomial constraints fi (x) ≥ 0 i = 1, . . . , m, and  ≥ 0,
we say a vector x̂ is -scaled feasible if fi (x̂) ≥ −kfi k1 ∀ i = 1, . . . , m.
The main result we will show in this chapter is as follows:
Theorem 4.1.2 Consider an instance of problem PO and any 0 <  < 1. Let ρ themaximum degree of any of the polynomials fi . Given a tree-decomposition of the intersection graph of the constraints of width ω, there is a linear programming formulation withO (2ρ/)ω+1 n log(ρ/) variables and constraints, such that any optimal solution x̂ for the
LP satisfies:
1. x̂ is -scaled feasible for (4.1b),
2. cT x̂ ≤ P ∗ + kck1 , where P ∗ is the value of PO.
The statement in Theorem 4.1.2 is indicative of the fact that as  → 0 we converge to an
optimal solution.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

57

In what follows we will provide a proof for Theorem 4.1.2 together with some of our core
constructions. Toward this goal, we begin by formally developing an approximation techniquefor PO.

4.2

Binary approximations of polynomial optimization problems

The general approximation scheme we will use for PO is based on a discretization method
originally proposed in [51]; also see [22, 39, 54] and citations therein for more insight. Let r
be a real with 0 ≤ r ≤ 1. Then, given 0 < γ < 1 we can approximate r as a sum of inverse
powers of 2, within additive error γ, by making use of the following. Let.
Lγ = dlog2 γ −1 e,then there exist 0/1-values zh , 1 ≤ h ≤ Lγ , so that
Lγ
X
h=1

2−h zh ≤

r

≤

Lγ
X

2−h zh + 2−Lγ ≤

h=1

Lγ
X

2−h zh + γ ≤ 1.

(4.3)

h=1

To apply this idea to problem PO, let 0 <  < 1 and as before let ρ denote the maximum
degree of any polynomial in PO. Choose γ so that = 1 − (1 − γ)ρ ,(4.4)

then for each j = p + 1, . . . , n (i.e, the indices of continuous variables) we will approximately
PLγ −hrepresent xj as h=1
2 zj,h where each zj,h is a (new) binary variable.

For each 1 ≤ i ≤ m and α ∈ I(i) we write
Z(i, α) = {j : αj 6= 0, 1 ≤ j ≤ p},in other words, set Z(i, α) is the set given by the indices of the binary variables for PO thatappear explicitly in monomial xα . We consider the following replacement for problem PO:

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH(GB(γ)) :

p
X

min

cj xj +

j=1

s.t.

fi,α α∈I(i)
Lγ
X

xj =

cj 

j=p+1


X



n
X

Y

j∈Z(i,α)
2−h zj,h 

h=1



n
Y

xj

Lγ
X

58

Lγ
X


j=p+1

αj 
2−h zj,h   ≥ −kfi k1 , 1 ≤ i ≤ m(4.5a)

h=1

2−h zj,h , p + 1 ≤ j ≤ n(4.5b)

h=1

xj ∈ {0, 1}, 1 ≤ j ≤ p(4.5c)

zj,h ∈ {0, 1}, p + 1 ≤ j ≤ n, 1 ≤ h ≤ Lγ

(4.5d)

where  = 1 − (1 − γ)ρ . Note that the dependency of  and γ can be swapped, as they arelinked through a bijection.

Next we prove a series of results describing the quality of the approximation to problem
PO provided by GB(γ). Toward this goal first we establish a technical property.
Lemma 4.2.1 Suppose that for j = 1, . . . s we have values uj ≥ 0, vj ≥ 0, qj ∈ Z+ withuj + vj ≤ 1. Then

s
Y

qj

−

(uj + vj )

j=1

s
Y

q
uj j

≤ 1−

j=1

s
Y

(1 − vj )qj .

j=1

Proof: Take any fixed index 1 ≤ i ≤ s. The expressionsY(uj + vj )qj −

j=1

s
Y

(uj )qj

j=1

is a non-decreasing function of ui when all uj and vj are non-negative, and so in the range
0 ≤ ui ≤ 1 − vi it is maximized when ui = 1 − vi .

Next we establish our approximation result.
Lemma 4.2.2(a) Suppose x̄ is feasible for PO. Then there is a feasible solution (x̃, z̃) for

GB(γ) with objective value at most cT x̄ + kck1 .
(b) Suppose (x̂, ẑ) is feasible for GB(γ). Then x̂ is -scaled feasible for (4.1b) and


Lγ
p
n
X
X
X
cT x̂ =
cj x̂j +cj 
2−h ẑj,h  .
j=1

j=p+1

h=1

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

59

Proof:(a) For each j ∈ {1, . . . , p} define x̃j = x̄j . For j ∈ {p + 1, . . . , n} choose binary valuesz̃j,h = z̃j,h (x̄j ) h = 1, . . . , Lγ so as to attain the approximation for x̄j as in (4.3). We
define x̃j for j ∈ {p + 1, . . . , n} from z̃ according to (4.5b).
We claim that for each 1 ≤ i ≤ m and α ∈ I(i)

αj

αj
Lγ
Lγ
n
n
n
Y
X
Y
Y
X
α


2−h z̃j,h ≤x̄j j ≤
2−h z̃j,h  + .
j=p+1

j=p+1

h=1

j=p+1

(4.6)

h=1

The left-hand inequality follows from (4.3). To obtain the right-hand inequality wePLγ −happly Lemma 4.2.1 with s = n and for p + 1 ≤ j ≤ n, qj = αj , uj = h=1
2 z̃j,h , andvj = x̄j − uj . Using (4.3) and the definition (4.4) of , Lemma 4.2.1 yields

αj
Lγ
n
n
n
Y
Y
X
Y
α

x̄j j −
2−h z̃j,h ≤1−(1 − vj )αj
j=p+1

j=p+1

j=p+1

h=1

≤1−

n
Y

(1 − γ)αj

j=p+1≤ 1 − (1 − γ)ρ
=
as desired. Therefore, (x̃, z̃) is feasible for GB(γ) and the second assertion in (a) issimilarly proved.
(b) By definition of x̂, for each 1 ≤ i ≤ m we havefi (x̂) =

X

fi,α

n
Y

α

x̂j j

j=1

α∈I(i)

=Xfi,α α∈I(i)Yj∈Z(i,α)

x̂j

n
Y



Lγ
X


j=p+1

αj 
2−h ẑj,h  

h=1≥ − kfi k1
as (x̂, ẑ) is feasible for GB(γ). The second assertion holds by definition of x̂.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

60

Corollary 4.2.3 Let P ∗ be the optimal value of problem PO, and let (x̂(γ), ẑ(γ)) be optimalfor GB(γ). Then cT x̂(γ) ≤ P ∗ + kck1 and x̂(γ) is -scaled feasible for (4.1b).
Remark 4.2.4 As per the corollary, x̂ achieves feasibility and optimality tolerance proportional to  for problem PO. However, x̂ may actually be super-optimal for PO. Neverthelessfor any sequence γk → 0+ , the vectors x̂(γk ) will have an accumulation point x∗ , and thispoint necessarily must be feasible and thus optimal for PO.
The reader may also notice in the proof of Lemma 4.2.2 (b) that the error bound −kfi k1
is conservative, because an approximation error is only incurred on terms of fi (x) that are
non-constant and that involve at least one continuous variable. Accounting for this fact will
provide a tighter approximation estimate, which we skip to simplify notation.
Additionally, note that in the proof of Lemma 4.2.1, the only property that each qj is
actually required to satisfy is that either qj = 0 or qj ≥ 1. One can use this observation to
generalize Theorem 4.1.2 beyond polynomials. However this produces additional complications such as the approximation of rational powers. We will take up these and related issues
in future work.
Corollary 4.2.3 implies that given any  > 0, problem GB(γ) with
γ = 1 − (1 − )1/ρ
provides an “-scaled” approximation to PO. Replacing the equality constraints (4.5b), problem GB(γ) becomes a pure binary problem, with m constraints and at mostnLγ = O(n log2 (ρ/))variables.

Henceforth we refer to GB(γ) as a pure binary problem. However, for the sake of notation,in the solution vectors of GB(γ) we will include the continuous components when referringto the x variables.
In the next section we will reformulate this problem as a linear program, thus yielding an
LP approximation to PO and completing the proof of Theorem 4.1.2.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

4.3

Linear reformulation of the binary approximation

4.3.1

Sparsity of the approximation

61

We already established that the pure binary problems GB(γ) provides an approximation to
PO within guaranteed tolerance. This result holds for all problems PO regardless of sparsity.
In this section we will show how to use structural sparsity to reformulate the pure binary
problems as equivalent LPs of moderate size.
To this effect, suppose that (T, Q) is a tree-decomposition for the intersection graph of an
instance of problem PO and choose γ as in (4.4) in order to obtain the desired tolerance .
We will now construct a tree-decomposition for the intersection graph for the corresponding
instance of problem GB(γ). This tree-decomposition will be of the form (T, Q0 ) (note: same
tree T ) where for any vertex t of T we set
Q0t = { xj : xj ∈ Qt , j = 1, . . . , p} ∪ { zj,h : xj ∈ Qt , j = p + 1, . . . , n and 1 ≤ h ≤ Lγ }.
Lemma 4.3.1(a) (T, Q0 ) is a tree-decomposition of the intersection graph for problem GB(γ).
(b) Further, if (T, Q) has width w then (T, Q0 ) has width at most Lγ (w + 1) − 1.
Proof:(a) We need to establish properties (i)-(iii) of Definition (3.1.1); (iii) is clear. Note thatzj,h ∈ Q0t if and only if xj ∈ Qt . This proves (i). To see that (ii) holds, note that anyedge of the intersection graph for GB(γ) arises from some constraint (4.1b).
(b) The width statement is true by construction of the sets Q0t and the fact that every|Qt | ≤ w + 1.

Remark 4.3.2 As a consequence of this result, not only does problem GB(γ) provide a close
approximation to problem PO, but when PO is structurally sparse (small tree-width of theintersection graph) then so is GB(γ), modulo the O(Lγ ) multiplicative increase in tree-width.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

62

Recall that Theorem 3.2.3 yields an LP reformulation of pure binary problems whose size
is parametrized by the tree-width of the intersection graph. Therefore, an application ofTheorem 3.2.3 along with Lemma 4.3.1 to GB(γ) will provide the desired result. We put
together these parts next.

4.3.2

From sparse PO to small LP approximations

Here we use the above results to prove Theorem 4.1.2. Consider an instance of PO on
n variables. As before, let ρ represent the maximum degree of any of the polynomials fi .
Suppose we have a tree-decomposition of width ω of the intersection graph of the constraints.
Given 0 <  < 1 we proceed as follows:
1. We choose
γ = 1 − (1 − )1/ρ ,so that  = 1 − (1 − γ)ρ as per (4.4). Note that without loss of generality  is small
enough so that γ > 2ρ.

2. We apply Theorem 3.2.3 to construct the linear programming reformulation of the
all-binary problem GB(γ) (formulation (4.5)). Let us call this linear program LP(γ).
3. As per Corollary 4.2.3, GB(γ) and thus, LP(γ), yields a vector x̂(γ) that is -scaled
feasible for (4.1b) and cT x̂(γ) ≤ P ∗ + kck1 .
Next we analyze the size of LP(γ).
1. By Lemma 4.3.1, there is a tree-decomposition of the intersection graph for GB(γ) of
width at most Lγ (ω + 1) − 1, where
Lγ = dlog2 γ −1 e < log2 (2ρ/) + 1
for  small enough.
2. Further, GB(γ) has at most nLγ = O(n log2 (ρ/)) variables.
3. Thus, since the tree-width of the intersection graph of GB(γ) is at most Lγ (ω + 1) − 1,by Theorem 3.2.3, the number of variables and constraints in LP(γ) isO((2ρ/)ω+1 n log2 (ρ/))CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

63

The proof of Theorem 4.1.2 is now complete.

4.4

Final comments

To conclude this chapter, we discuss additional aspects of the approximations we just introduced. We will discuss if an improvement on the dependency of  in Theorem 4.1.2 is possible,
and provide a complete example on the resulting LP approximation we obtain applying thetechnique outlined in Section 4.3.2.

4.4.1

Can the dependence on  be improved upon?

The approximation scheme given by Theorem 4.1.2 has two characteristics: first, it allows aviolation of each constraint by  times the 1-norm of the constraint, and second, the running
time is pseudo-polynomial in −1 . One may wonder if either characteristic can be improved.
For example, one might ask for constraint violations that are at most , independent of the
1-norm of the constraints. However this is not possible even for a fixed value of , unless
P=NP. For completeness, we include a detailed analysis of this fact in Section C.1 of the
Appendix. Intuitively, if we were allowed to approximately satisfy every constraint with anerror that does not depend on the data, we could appropriately scale constraint coefficientsso as to obtain exact solutions to NP-hard problems.
Similarly, it is not possible to reduce the pseudo-polynomial dependency on −1 in general.
The precise statement is given in Section C.2 of the Appendix, and the intuitive reasoning issimilar: if there was a formulation of size polynomially dependent on log(−1 ) (and not on−1 ) we could again solve NP-hard problems in polynomial time.

4.4.2

Example of LP approximation to PO

We now give a full detailed example with the explicit formulation we will obtain when approximating a problem PO with the pure binary problem GB(γ), and then reformulating
the latter as an LP. Consider the following low-dimensional example for problem PO:min{ x1 + x2 + x3 : x21 + x22 ≥ 1.95, x23 + x22 ≥ 1.95, x ∈ [0, 1]3 }.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

64

Suppose we build our binary powers approximation using γ = 1/8 = 0.125 (and as per(4.4)  = 0.234375), thus Lγ = 3. We obtain the pure binary approximation with N = 9variables:GB( 18 ) :

min

3
X

2−h z1,h +

h=1

s.t.

3
X

2−h z2,h +h=1!2
2−h z1,h

+

h=1
3
X

3
X

3
X

3
X

2−h z3,h(4.7a)h=1!22−h z2,h≥ 1.95 − 3.95 ≈ 1.024(4.7b)≥ 1.95 − 3.95

(4.7c)h=1!22−h z3,hh=1

zj,h ∈ {0, 1},

+

3
X

!2
2−h z2,hh=1

1 ≤ j ≤ 3,

1 ≤ h ≤ 3.

(4.7d)

A tree-decomposition of the intersection graph is (T, Q) where• T is the tree with vertices α and β and the single edge {α, β}, and• Qα = {z1,1 , z1,2 , z1,3 , z2,1 , z2,2 , z2,3 }• Qβ = {z3,1 , z3,2 , z3,3 , z2,1 , z2,2 , z2,3 }Next, in order to construct an LP formulation of (4.7) using Theorem 3.2.3, we construct
the set Fα , i.e. the set of assignments of binary values to the members of Qα such that
constraint (4.7b) holds. We will argue that there are six such assignments, which we indicate
as binary vectors using the same ordering as in the definition of Qα :(a) (1, 0, 1, 1, 1, 1), (b) (1, 1, 0, 1, 1, 0), (c) (1, 1, 0, 1, 1, 1)(d) (1, 1, 1, 1, 0, 1), (e) (1, 1, 1, 1, 1, 0), (f) (1, 1, 1, 1, 1, 1).
To show that (a)-(d) are the only feasible vectors, we first note that if either z1,1 = 0 or
z2,1 = 0 then (4.7b) is violated. Also, if z1,u = 0 and z2,v = 0 for some u, v then u = v = 3 andall other z1,w and z1,w must equal 1. Likewise there are six similar assignments corresponding
to the members of Qβ ; we denote these (g), (h), (i), (j), (k), (l). We assume these areconstructed in the same way as (a)-(f).
Now we turn to the construction of LPz (formulation (3.6)). As discussed there, we needto add constraints (3.6b) in just two kinds of cases: singletons S, and S which are common

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH

65

subsets of Qα and Qβ . In this latter case we have all S ⊆ {z2,1 , z2,2 , z2,3 }. We can now write
the LP (3.6) for our example. In this LP we will abbreviate a variable zj,h as the pair (j, h).
3 X
3X

2−h Z{(j,h)}(4.8a)

s.t. Z{(1,1)} = λαa + λαb + λαc + λαd + λαe + λαf(4.8b)

Z{(1,2)} = λαb + λαc + λαd + λαe + λαf(4.8c)

Z{(1,3)} = λαa + λαd + λαe + λαf(4.8d)

Z{(2,1)} = λαa + λαb + λαc + λαd + λαe + λαf(4.8e)

Z{(2,2)} = λαa + λαb + λαc + λαe + λαf(4.8f)

Z{(2,3)} = λαa + λαc + λαd + λαf(4.8g)

Z{(2,1)} = λβg + λβh + λβi + λβj + λβk + λβl(4.8h)

Z{(2,2)} = λβg + λβh + λβi + λβk + λβl(4.8i)

Z{(2,3)} = λβg + λβi + λβk + λβl(4.8j)

Z{(3,1)} = λβg + λβh + λβi + λβj + λβk + λβl(4.8k)

Z{(3,2)} = λβh + λβi + λβj + λβk + λβl(4.8l)

Z{(3,3)} = λβg + λβj + λβk + λβl(4.8m)

Z{(2,1),(2,2)} = λαa + λαb + λαc + λαe + λαf(4.8n)

Z{(2,1),(2,2)} = λβg + λβh + λβi + λβk + λβl(4.8o)

Z{(2,1),(2,3)} = λαa + λαc + λαd + λαf(4.8p)

Z{(2,1),(2,3)} = λβg + λβi + λβj + λβl(4.8q)

Z{(2,2),(2,3)} = λαa + λαc + λαf(4.8r)

Z{(2,2),(2,3)} = λβg + λβi + λβl(4.8s)

Z{(2,1),(2,2),(2,3)} = λαa + λαc + λαf(4.8t)

Z{(2,1),(2,2),(2,3)} = λβg + λβi + λβl(4.8u)

λαa + λαb + λαc + λαd + λαe + λαf = λβg + λβh + λβi + λβj + λβk + λβl = 1(4.8v)min

j=1 h=1

λ ≥ 0.
The above LP (4.8) produces the solution x1 = x3 = .625 and x2 = .875 with value 2.15.

CHAPTER 4. MIXED-INTEGER POLYNOMIAL OPTIMIZATION WITH SMALL
TREE-WIDTH66

The optimum solution for the original problem is on the other hand attained by x1 = x3 ≈0.9747, x2 = 1.0, with objective value ≈ 2.9494.
Remark 4.4.1 The goal of this example is to illustrate our techniques with an instance thatyields a small LP formulation. Should we use γ = 1/16 (Lγ = 4) for example, we would
obtain  ≈ 0.12 which would reduce our approximation error by a factor of approximately twoat the cost of a larger LP.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

67

Chapter 5

Network Polynomial Optimization
In this chapter we focus on polynomial optimization problems PO where there is network
structure as an explicit element of the problem description. In these problems we are given
a network1 , as in the AC-OPF problem, and both variables and constraints are associatedwith the network structure directly. We will study how low tree-width in this underlying
network affects tractability, and show how different this is from assuming low tree-width in
the intersection graph, as in the previous chapters.

5.1

Problem description

We consider problems we term Network Polynomial Optimization problems (NPO); these are
PO problems defined using the structure of a network. As was mentioned in the Introduction,in an NPO over a network G, for each node u ∈ V (G) there is a set Xu of variables associatedwith u and a set Ku of constraints associated with u. The optimization problem is of thefollowing form:(NPO):
subject to:min cT x
X(5.1a)(k)p(u,v) (Xu ∪ Xv ) ≥ 0,k ∈ Ku ,u ∈ V (G)(5.1b){u,v}∈δG (u)x ∈ {0, 1}p × [0, 1]n−p ,1(5.1c)We use the term ‘network’ to contrast with ‘graph’ which we reserve for intersection graphs. Likewise wewill typically use ‘node’ to refer to a vertex of a network.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

2
1

..

68

.
n

0

Figure 5.1: Star network witn n + 1 nodes.
(k)where δG (u) is the set of edges of G incident with u and each p(u,v) (Xu ∪ Xv ) is a polynomial(k)associated with the ordered pair (u, v). The notation p(u,v) (Xu ∪ Xv ) stresses the ordered pair(k)(u, v) and indicates that the polynomial only involves variables in2 Xu ∪ Xv . The p(u,v) will betermed the arc polynomials of the problem. We do not assume that the sets Xu are pairwise
disjoint. Rather (for a technical reason) we allow intersections but we assume that for any
variable xj , the set of u ∈ V (G) with xj ∈ Xu induces a connected subgraph of G. For futurereference, we will refer to this assumption as the connectedness requirement.
The definition of NPO can be seen as a direct generalization of the AC-OPF problem,
where the sets Xu correspond to the voltage variables in each node u, and the arc polynomials
are defined using the power flow equations.
We are interested in achieving tractability of NPO when G is appropriately sparse, i.e,
exploting its tree-width. In devising an algorithm for problem NPO parametrized by thetree-width of G, a direct reliance on Theorem 4.1.2 runs into a difficulty as detailed next.
Example 5.1.1 Consider an NPO where(i) The network G is a star on n + 1 nodes 0, 1, . . . , n and edge set {0, j} for 1 ≤ j ≤ n.
See Figure 5.1.
(ii) X0 = ∅ and Xj = {xj } for 1 ≤ j ≤ n.
(iii) There is a single constraint (5.1b) and it is associated with node 0. In this constraint(1)

we have p(0,j) (X0 ∪ Xj ) = jx3j − 1/n for 1 ≤ j ≤ n.
2

Using this notation we allow cases where only variables in Xu or in Xv are actually involved. Without
(k)

(k)

loss of generality we assume that for a given k, no two p(u,v) and p(u,v0 ) use a common monomial.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION
Thus the (single) constraint (5.1b) for this example is:Pn3
j=1 jxj69≥ 1. This is a dense

constraint, even though the input network G is a tree. In fact, the intersection graph of thisconstraint is a clique on n vertices (with tree-width n − 1), and, as a consequence, if we
directly apply Theorem 4.1.2 we will obtain a formulation of size exponential in n.
Another example of this issue can be given using the well-known kanpsack problem:
Example 5.1.2 Consider a knapsack problem min{ cT x : aT x ≥ b, x ∈ {0, 1}n }. Thiscan be cast as NPO, using a star network on n + 1 nodes as in Example 5.1.1, which has
tree-width 1. Yet, if aj 6= 0 for all j, the intersection graph is a clique of size n.
In summary: even if an NPO problem arises from a network with small tree-width,
the problem, if viewed directly as an instance of PO, may yield a very dense formulation.
However, we will argue that one can always reformulate any instance of NPO over a small
tree-width network G so as to obtain an equivalent instance which, when viewed as a general
problem PO gives rise to an intersection graph with (still) small tree-width.
(k)

For simplicity, in some cases we abbreviate constraint (5.1b) as pu (x) ≥ 0. The main
result we will prove in this chapter is as follows.
Theorem 5.1.3 Consider an instance of NPO and any 0 <  < 1. Let D, ∆ and ρ be suchthat:
• The network G has at most D edges incident with any node,• The number of variables plus the number of constraints associated with any node of G
is at most ∆,(k)• Every polynomial pu has maximum degree ≤ ρThen there is a linear program of size O((Dρ/)O(∆ω) n log(ρ/)), such that any optimal
solution x̂ for the LP satisfies:
1. x̂ is -scaled feasible for (5.1b),
2. cT x̂ ≤ P ∗ + kck1 , where P ∗ is the value of NPO.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION70

In Example 5.1.1 ω = 1, D = n, ρ = 3 and ∆ = 1. Thus Theorem 5.1.3 yields an LP
formulation of size O((n/)O(1) log(1/)).

5.2

Impact of node-degrees in an NPO

Let us continue analyzing Example 5.1.1, which illustrates how it is not possible to directly
obtain Theorem 5.1.3 as a consequence of Theorem 4.1.2. In this example the network G
was a “star” with node set 0, 1, . . . , n and center node 0. The NPO had a single constraint(1)(5.1b), associated with node 0, with arc polynomials p(0,j) (X0 ∪ Xj ) = jx3j − n1 for 1 ≤ j ≤ n.
Pn
3
The corresponding constraint (5.1b) reads:
j=1 jxj ≥ 1. This is a dense constraint and a
direct application of Theorem 4.1.2 will produce a formulation of size exponential in n.
This example illustrates the point that if, in an NPO, a node has high degree, theintersection graph of the NPO will likely have high tree-width. This observation suggeststhat we should try to reformulate an NPO into an equivalent NPO on a network where
every node has small degree. The following construction is a counterpoint to Example 5.1.1
and is a central component in our approach to proving Theorem 5.1.3.
Consider an NPO, P̄, on a network Ḡ all of whose nodes have degree at most three, andwith associated sets of variables X̄u . Let (T̄ , Q̄) be a tree decomposition of Ḡ. Form the pair(T̄ , Q̄0 ), where for each t ∈ V (T̄ ) we defineQ̄0t =[{X̄v : v ∈ Q̄t or {u, v} ∈ E(Ḡ) for some u ∈ Q̄t }.

(5.2)

We have that:
Theorem 5.2.1 (T̄ , Q̄0 ) is a tree-decomposition for the intersection graph of NPO P̄. If
the width of (T̄ , Q̄) is w̄, then the width of (T̄ , Q̄0 ) is at most 3(w̄ + 1)k̄ − 1, where k̄ =maxu∈V (Ḡ) |X̄u |.
We will prove this result shortly. The result suggests a way to obtain Theorem 5.1.3.
Namely, given an NPO, P, on a general network of “small” tree-width we reformulate it as
an equivalent NPO, P̄, on a network Ḡ with nodes of degree at most three and such thatḠ also has “small” tree-width. If, in addition the parameter k̄ in the statement of Theorem
5.2.1 is also “small”, then our reformulation will be an NPO which, as a general polynomial

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

71

optimization problem, has an intersection graph with small tree-width and thus can be handled using Theorem 4.1.2. Note that, additionally, Theorem 4.1.2 only yields approximately
feasible solutions, and so part of our task will be to guarantee that approximate feasibility
for P̄ carries over to P.

Proof: (Theorem 5.2.1) The width statement follows directly from definition (5.2). Thus
our main task is to show that (a) any pair of variables of P̄ that occur in a common constraintof P̄ are also found in some common set Q̄0t , and (b) for any variable xj of P̄ the set of vertices
t of T̄ such that Q̄0t contains xj induces a subtree of T̄ .
(a) Consider a pair of variables {xi , xj } that occur a common constraint of P̄, say in aconstraint associated with node u. Then{xi , xj } ⊆ X̄u ∪ X̄v ∪ X̄w ,where {u, v} and {u, w} are edges of Ḡ (possibly v = w). But in that case for any t ∈ T̄such that u ∈ Q̄t we will have that {xi , xj } ⊆ Q̄0t by construction of the set Q̄0t in (5.2).
(b) Here we use Remark 3.1.5. Consider any variable xj . Let us define.
Cj = {u ∈ V (Ḡ) : xj ∈ X̄u }and.
Nj = Cj ∪ {v ∈ V (Ḡ) : {u, v} ∈ E(Ḡ) for some u ∈ Cj }.
These sets are relevant because we defined Q̄0t such that xj ∈ Q̄0t iff Qt intersects Nj .
But we note that Cj induces a connected subgraph of Ḡ – in Section 5.1 this attribute
was called the connectedness requirement for an NPO. As a result, Nj also induces aconnected subgraph of Ḡ. If we apply the Remark 3.1.5 to Nj we obtain that, indeed,
the set of t such that xj ∈ Q̄0t induces a subtree of T̄ , as desired.

In the following sections we will follow the approach to proving Theorem 5.1.3 suggested
above, beginning with the reformulation of an NPO into an equivalent NPO on a networkwith degree ≤ 3 node, aiming toward the use of Theorem 5.2.1.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

5.3

72

Example of NPO reformulation

As motivation for the general recipe we just outlined, we will argue that the problem in
Example 5.1.1 can be reformulated as an equivalent NPO on network of maximum degree
three (and with small tree-width). And then, using the strategy in Theorem 5.2.1 we willargue that the intersection graph of such reformulation has small tree-width.
To fix ideas, let us consider the case n = 4. We will first produce the equivalent NPO byP
first constructing an extended formulation equivalent to 4j=1 jx3j ≥ 1. Later we will showthat this extended formulation amounts to a new NPO.
The extended formulation has additional variables yj (1 ≤ j ≤ 7), and the followingPsystem of constraints whose sum yields 4j=1 jx3j ≥ 1.
1
yj ≤ jx3j − ,41 ≤ j ≤ 4,y5 ≤ y1 + y2 ,y6 ≤ y3 + y4 ,(5.3a)y7 ≤ y5 + y6 ,y7 ≥ 0.
Effectively, this system splits the sum(5.3b)(5.3c)

P4

3
j=1 jxj into partial sums with two terms each; y5

and y6 are stand-ins for these partial sums and y7 represents the complete sum3 .
P
System (5.3) is equivalent to 4j=1 jx3j ≥ 1 in the sense that the projection to x-spaceP4
3of the set of solutions to (5.3) equals {x ∈ R4 :j=1 jxj ≥ 1}. Let us put aside, for the
moment, the issue that this equivalence might require some of the yj to take values outsideof the range [0, 1], which is not allowed in our formal definition (5.1) for an NPO. Modulo
this point, we can argue that (5.3) is the system of constraints for an NPO. To constructthis NPO we use a binary tree with nodes 1, . . . , 7 shown in Figure 5.2.
The variables for this NPO will be all the xj and yj . We associate with node 5 a family ofarc polynomials that will yield (5.3a) for j = 1, 2 and the first inequality in (5.3b). Namely, we
associate with node 5 the set of variables {x1 , . . . , x4 }∪{y1 , y2 , y5 }, define the arc polynomials
1(1)
p(5,1) = −y1 + x31 − ,
4
1
(2)
3
p(5,2) = −y2 + 2x2 − ,
4
3

(1)

p(5,2) = 0,(2)

p(5,1) = 0,This “splitting” technique is reminiscent of sparsification techniques for linear systems [52].

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

73

1
5
2
7
3
6
4

Figure 5.2: Binary tree replacement for star with 4 leaves.
1(3)
p(5,1) = − y5 + y1 ,2

1
(3)
p(5,2) = − y5 + y2 ,2and impose the NPO constraints (associated with node 5):
(k)

(k)

p(5,1) + p(5,2) ≥ 0,

k = 1, 2, 3.

(5.4)

Likewise, we associate with node 6 the set {x1 , . . . , x4 } ∪ {y3 , y4 , y6 }, define1(1)p(6,3) = −y3 + 3x33 − ,4
1(2)
p(6,4) = −y4 + 4x34 − ,41(3)
p(6,3) = − y6 + y3 ,2(1)

p(6,4) = 0,(2)

p(6,3) = 0,
1
(3)
p(6,4) = − y6 + y4 ,2and impose
(k)

(k)

p(6,3) + p(6,4) ≥ 0,

k = 1, 2, 3,(5.5)

which yields (5.3a) for j = 3, 4 and the second inequality in (5.3b). Finally we associate with
node 7 the variables {x1 , . . . , x4 } ∪ {y7 }, define1(1)p(7,5) = y5 − y7 ,2

1
(1)
p(7,6) = y6 − y7 ,
2

(5.6)

and set the constraint
(1)

(1)

p(7,5) + p(7,6) ≥ 0.

(5.7)

This yields the third inequality in (5.3b). The bound y7 ≥ 0 is already implicit in the definition of an NPO. Thus, indeed, system (5.3) arises as the constraint set for an NPO.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

x1

x2

x3

x4

y1

y2

y3

y4

x1 , y174

x2 , y2x3 , y3

y1 , y2y5y6x4 , y4y3 , y4y5y6y5 , y6

y7

(a)

(b)

y7

Figure 5.3: (a) Intersection graph for reformulation of Ex. 5.1.1. (b) A tree-decomposition.
Next we claim that as a general polynomial optimization problem, the problem withconstraints (5.4)-(5.7) has an intersection graph of tree-width 2.
The intersection graph is shown in Figure 5.3 (a) and a tree decomposition, of width 2, isgiven in Figure 5.3 (b). Note that the tree in Figure 5.3 (b) is isomorphic to that in Figure 5.2– this is not an accidental event and, rather, reflects the structure of constraints (5.4)-(5.7).
P
It is clear that the above process can be applied to the general case of nj=1 jx3j ≥ 1 so as
to always yield a reformulation as an NPO on a binary tree with n leaves, a yj variable foreach internal node (and so less than n yj variables), and such that the intersection graph has
tree-width 2.
In summary, system (5.3) represents a natural way to rewrite our original NPO as aproblem with structured sparsity. To complete this argument, however, we return to thePaforementioned issue that the equivalence between system (5.3) and 4j=1 jx3j ≥ 1 might
require that some yj falls outside [0, 1]. We handle this issue through a further reformulation,using the familiar trick of representing a real variable as the difference between two nonnegative variables and scaling. System (5.8) implements these ideas; we discuss the choice ofscaling factors below:(j + 1/4)(πj − µj ) ≤ jx3j − 1/4,

1≤j≤4

(5.8a)(3 + 1/2)(π5 − µ5 ) ≤ (1 + 1/4)(π1 − µ1 ) + (2 + 1/4)(π2 − µ2 ),(5.8b)(7 + 1/2)(π6 − µ6 ) ≤ (3 + 1/4)(π3 − µ3 ) + (4 + 1/4)(π4 − µ4 ),(5.8c)

11 π7 ≤ (3 + 1/2)(π5 − µ5 ) + (7 + 1/2)(π6 − µ6 ),(5.8d)x ∈ [0, 1]4 , π, µ ∈ [0, 1]6 , p7 ∈ [0, 1].

(5.8e)

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION75

The scaling constant j + 1/4 in (5.8a) can intuitively be justified by noting that if 0 ≤ xj ≤ 1then (5.8a) allows us to assume that πj ≤ 1 and µj ≤ 1. Indeed, we could have used j − 1/4
for this purpose; but the choice of j + 1/4 will simplify language in our proofs later. The
j + 1/4 constants are then carried into the right-hand side of (5.8b) and (5.8c) with e.g. thecoefficients 3 + 1/2 chosen so that we can argue that without loss of generality 0 ≤ π5 , µ5 ≤ 1,
and similarly with π6 , µ6 . Finally in (5.8d) we use the scaling multiplier of 11 so that π7 ≤ 1.
We can verify that the projection of the feasible set for (5.8) to x-space is {x ∈ [0, 1]4 :P4

3j=1 jxj ≥ 1}:

first, summing the inequalities in (5.8a)-(5.8d), and using π7 ≥ 0 yields
P4
3
3
j=1 jxj − 1 ≥ 0 we can choose π, µ so thatj=1 jxj − 1 ≥ 0. Conversely, given x such thatP44(5.3) holds using the following procedure. For 1 ≤ j ≤ 4 we set πj = 4j+1
max{jx3j − 14 , 0} ∈4
[0, 1] and µj = 4j+1
max{−jx3j + 41 , 0} ∈ [0, 1], π5 = 27 max{(1 + 1/4)(π1 − µ1 ) + (2 + 1/4)(π2 −µ2 ), 0} ∈ [0, 1], and so on.
To complete the switch from the y to the π, µ variables, in the tree-decomposition in
Figure 5.3 (b) we replace each yj by the pair πj , µj for 1 ≤ j ≤ 6 (and y7 with π7 ), obtaining
a tree-decomposition for the intersection graph of constraints (5.8) of width 5.
In summary, we have verified that the NPO in Example 5.1.1 can be reformulated as anequivalent NPO whose intersection graph has small tree-width. As a result, Theorem 4.1.2
yields a polynomial-size linear program that approximates the latter NPO.

5.4

NPO reformulation: general caseIn the following sections we will construct a proof for Theorem 5.1.3. As Example 5.1.1
shows, given an NPO on a network G it may be difficult to directly apply Theorem 4.1.2
so as to obtain Theorem 5.1.3 if G has nodes of high (graph-theoretic) degree. Our method
for handling Example 5.1.1 effectively splits the high-degree node 0 into a binary tree andappropriately reformulate the NPO.
Our general methodology builds on these ideas using the following strategy:(s.1) We will present a general technique for reformulating any NPO as an equivalent NPO
on a network where every node has degree at most three, obtained through a sequenceof node splitting operations. In general, several such reformulations will exist.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

..

76

A1 .
..

.
−→

. .A2
.
u1

u

u2
r

Figure 5.4: Node splitting(s.2) If G has small tree-width, we will show that there is a specific reformulation as in(s.1) where the resulting NPO, as a general polynomial optimization problem, has anintersection graph of small tree-width.
In Sections 5.4.1 and 5.4.2 we will focus on (s.1); here we will essentially repeat thetechnique used above in the reformulation of Example 5.1.1. Section 5.4.3 outlines how the
proof of Theorem 5.1.3 is obtained if (s.2) is assumed, and Section 5.4.4 focuses on (s.2) andwill complete the proof.

5.4.1

Transforming the network and reformulating

In this section we discuss a procedure that, given an NPO on a network G, yields an equivalent
NPO on a new network where every node has degree at most three. This procedure will
operate sequentially, at each step modifying the current network so as to decrease by one thenumber of nodes of degree larger than three, while at the same time producing an NPO on
the modified network which is equivalent to the original NPO.
We use the idea of node splitting, which we first discuss informally in a more generalcontext than Example 5.1.1. Suppose u ∈ V (G) has degree larger than three and consider apartition of δ(u) into two sets A1 , A2 . We obtain a new network G 0 from G by replacing uwith three new nodes, u1 , u2 and r, introducing the edges {u1 , r} {r, u2 } and replacing eachedge {u, v} ∈ Ai (for i = 1, 2) with {ui , v}. See Figure 5.4.
Iterating this procedure, given u ∈ V (G) with degG (u) > 3, we can replace u and the setof edges δ(u) with a tree, where each internal node will have degree equal to three except forthe special node r. See Figure 5.5.
In this figure, a degree-5 node u and the set of edges {(u, vi ) : 1 ≤ i ≤ 5} is converted intoa tree with four internal nodes (r, a, b and c). There is a one-to-one correspondence between

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

v1

v2

v3

v4

v5

−→

77

c

r

b

v1

a

v4

v5

u
v2

v3

Figure 5.5: Complete node splitting
the edges {u, vi } in the original network, and edges of the tree incident with leaves. Also
note that this operation does not change the degree of the vi .
In general there are several possible complete node splittings, and alternate sequences ofnode splittings can produce networks with dramatically different tree-width, an issue that wewill tackle in Section 5.4.4.
Formally, the splitting operation goes as follows. Consider u ∈ V (G) with degG (u) > 3.
Let T̂u be an arbitrary tree where(t.1) The set of non-leaf or internal nodes of T̂u is disjoint from V (G). Further, the set ofleaves of T̂u is {v ∈ V (G) : {u, v} ∈ E(G)}.
(t.2) Each internal node of T̂u has degree equal to three, except for a special node r (theroot), of degree two.
Completely splitting u using T̂u yields a new network, G 0 where V (G 0 ) consists of V (G)\{u}together the internal nodes of T̂u , and E(G 0 ) = (E(G) \ δG (u)) ∪ E(T̂u ). For convenience,
below we will think of T̂u as rooted at r (i.e., with all edges oriented away from r). Note thatthe degree of any v ∈ V (G) with v 6= u is unchanged in this process, thus we can iteratively
apply this process to any node with degree four or more. Such algorithm will terminate witha network where each node has degree at most three. Each complete splitting operation
will be accompanied by a reformulation of the NPO that yields an equivalent NPO. This
reformulation is described next.
Assume a complete splitting at u ∈ V (G) using T̂u . First, at all nodes v 6= u of G theoriginal NPO remains unchanged. For each internal node i of T̂u , in our reformulation we
first have:CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

78(v.0) All variables in Xu are associated with i in the new NPO.
Let k ∈ Ku be the index representing a constraint (5.1b) that is associated with u. We
will replace this constraint with a family of constraints associated with the internal nodes ofT̂u . This will be done in steps (v.k), (c.1.k) and (c.2.k) given next. Define
Ni,k =oX n (k)kp(u,v) k1 : {u, v} ∈ δ(u), v a descendant of i in T̂u

The quantities 5/4, 5/2 and 5 that we saw in (5.8) are precisely the Ni,k for Example 5.1.1.
Note that if i is an internal node of T̂u with children j, l, then our definition implies thatNi,k = Nj,k + Nl,k .

We now associate the following additional variables with node i:(v.k) If i 6= r we additionally associate two variables, πi,k , µi,k ∈ [0, 1]. If i = r we associate
a single additional variable πr,k ∈ [0, 1] with r.

If i has a child v that is a leaf, and hence {u, v} ∈ δG (u), then we associate with i two
additional variables πv,k , µv,k ∈ [0, 1]. We emphasize that these variables are associatedwith i but not with v.
Next we describe the arc polynomials. We omit polynomials that are identically zero.
(c.1.k) Suppose v is a child of i which is a leaf, i.e. {u, v} ∈ δG (u). Then we define the arc(k)polynomial kp(u,v) k1 (πv,k − µv,k ) and write the constraint (also associated with i)(k)kp(u,v) k1 (πv,k − µv,k )≤(k)p(u,v) (x).

(5.9)[See constraint (5.8a).] We note that this constraint is well-defined as an NPO constraint, because by (v.0) and (v.k) all variables in (5.9) are associated with node i orwith node v.
(c.2.k) Let j, l be the children of i. Then we define the (linear) arc polynomials Nj,k (πj,k −µj,k )and Nl,k (πl,k − µl,k ), associated with edges {i, j} and {i, k}, respectively.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

79

If i 6= r, then we define the additional arc polynomial Ni,k (πi,k − µi,k ) (which can beassociated with any edge incident with i) and write the constraint
Ni,k (πi,k − µi,k ) ≤ Nj,k (πj,k − µj,k ) + Nl,k (πl,k − µl,k ).

(5.10)[See constraints (5.8b), (5.8c).]

If on the other hand i = r, then we write the constraint [see (5.8d)]Nr,k πr,k ≤X

Ns,k (πs,k − µs,k ).

(5.11)

s=j,l

We have just described how the variables and constraints associated with node u in theoriginal NPO are to be reformulated in the new NPO. To complete the description of the
new NPO we consider a generic node v ∈ V (G) with v 6= u. Formally,(v’) We keep the same set Xv of variables associated with v.
We also keep the same set of constraints Kv associated with v. This is done by, effectively,
keeping the same set of arc polynomials p(v,w) . Specifically, consider any edge {v, w} ∈ E(G 0 ):(c’.1) If w 6= u, the original arc polynomials p(v,w) are arc polynomials in the reformulated
problem because of (v’) above.
(c’.2) If w = u then as discussed in item (t.1) above, there is one edge {i, v} of T̂u . By
item (v.0), all variables in Xu are associated with node i. Hence, we any arc polynomialp(v,u) in the original NPO can be viewed as an arc polynomial p(v,i) in the reformulated
NPO.

5.4.2

Validity of the reformulation

Here we consider a single node u and and index k ∈ K(u) as in Section 5.4.1. Lemma 5.4.1
proves that the original and reformulated NPOs are equivalent, and Lemma 5.4.3 discusses
the quality of the approximation to the first NPO provided by an approximate solution to
the second.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

80

Lemma 5.4.1(a) Suppose (π̂, µ̂, x̂) satisfies constraints (5.10), (5.11), (5.9) and additionP(k)ally π̂r,k ≥ 0. Then x satisfies the constraint {u,v}∈δG (u) p(u,v) (x) ≥ 0.

(b) Conversely, if(k){u,v}∈δG (u) p(u,v) (x̂) ≥ 0 with x̂jP∈ [0, 1] for all j, then we can find(π̂, µ̂, π̂r,k ) that satisfy (5.10), (5.11), (5.9) with all variables in the range [0, 1].
Proof:(a) This fact follows because summing (5.10), (5.11), (5.9) yields
X(k)p(u,v) (x̂) ≥ Nr,k πr,k .
{u,v}∈δG (u)(b) This fact is obtained as in the discussion following equations (5.8). We proceed, inductively, up the tree, choosing values π̂ and ν̂ so that all constraints (5.10) and (5.9) holdas equality, and so that at each internal node i, π̂i,k ≥ 0, µ̂i,k ≥ 0 and π̂i,k µ̂i,k = 0.
Inductively, by definition of the constants Ni,k , no variable will ever exceed 1. At
termination we will be able to likewise set π̂r,k so that (5.11) also holds an equation.

Part (a) of Lemma 5.4.1 shows that a feasible solution for the updated NPO gives rise
to a feasible solution to the original NPO. If, however, we were to apply Theorem 4.1.2 tothe new NPO we would only obtain an approximately feasible solution (π̂, µ̂, x̂) and we needto address the infeasibility of x̂ with respect to the original NPO. This issue is handled by
Lemma 5.4.3. We first need a technical result.
Lemma 5.4.2 The sum of 1-norms of all constraints (5.10), (5.11), (5.9) is at mostX(k)
4 degG (u)kp(u,v) k1 .
(u,v)∈δG (u)Proof: Write d = degG (u). The 1-norm of a constraint (5.10) at an internal node i (or(5.11) if i = r) equals 2Ni,k . Since all internal nodes of T̂u , except for r, have degree three
in T̂u and T̂u has d leaves, the total number of internal nodes equals d − 1. Thus the sumP(k)of 1-norms of all constraints (5.10) or (5.11) is at most 2d (u,v)∈δG (u) kp(u,v) k. A similar
observation regarding constraints (5.9) yields the result.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

81

Lemma 5.4.3 Let γ > 0. Suppose (π̂, µ̂, x̂) is γ-scaled feasible for constraints (5.10), (5.11),(5.9), and satisfies π̂r,k ≥ 0. Then
(k)

X

X

p(u,v) (x̂) ≥ −4 degG (u) γ(k)kp(u,v) k1 ,{u,v}∈δG (u){u,v}∈δG (u)i.e. x̂ is (4 degG (u) γ)-scaled feasible for constraint(k){u,v}∈δG (u) p(u,v) (x) ≥ 0.

P

Proof: Writing a generic constraint (5.10), (5.11), (5.9) as f (π, µ, x) ≥ 0 we have byassumption that f (π̂, µ̂, ŷr , x̂) ≥ −kf k1 γ. Since the sum of all constraints (5.10), (5.11), (5.9)
yields(k)Xp(u,v) (x) ≥ Nr,k πr,k(u,v)∈δG (u)we therefore obtain(k)Xp(u,v) (x̂) ≥ Nr,k π̂r,k − Sγ,(u,v)∈δG (u)where S is the sum of L1 norms of all constraints (5.10), (5.11), (5.9). Applying Lemma 5.4.2
yields the desired result.

We can now comment on the quality of an approximate solution to the reformulation, in
terms of the initial problem. Let us suppose that given an NPO on a network G (we will term
this the initial NPO) we apply a sequence of complete node splittings and reformulations so
as to obtain a final NPO. In what follows, we will denote P f is the final NPO, and G f the
final network.
Lemma 5.4.4 P f is equivalent to the initial NPO and a solution that is γ-scaled feasible
for the constraints in P f yields a (4Dγ)-scaled feasible solution to the initial NPO, whereD = max degV (G) (v).
v∈V (G)

Proof: Equivalence follows by applying Lemma 5.4.1 inductively. Applying Lemma 5.4.3
inductively proves the approximation error bound.

For future reference we also state the following result concerning the size of the reformulated NPO:

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

82

Lemma 5.4.5 Let ∆ = maxv∈V (G) {|Kv |}. Suppose the original NPO has n variables intotal, and |E(G)| = m. Then G f has O(m) nodes and edges. Further, the number of variables
and constraints in P f

n + O

X

degG (u)|Ku | = n + O(D∆).

(5.12)u∈V (G)

Proof: SincePv∈V (G) degG (v) = 2m we obtain that G

f has the desired size. The variables

for P f are the original variables, plus all π and µ variables, which are accounted for in thesecond term of (5.12).

In the next section we group all arguments together, including the approximate LP formulations that arise from applying Theorem 4.1.2 to P f .

5.4.3

Implications to Theorem 5.1.3

Let us follow the above recipe to obtain a final reformulated NPO P f on a network G f . We
will now apply Theorem 4.1.2 by treating P f as a general polynomial optimization problemof the form PO. In this section we address the size of the resulting LP and its approximation
guarantees. In particular, suppose that we want to produce a solution for the initial NPO
within tolerance , i.e. -scaled feasible and with objective value exceeding the optimum byat most kck1 .
1. Using Lemma 5.4.4 we see that we need to find a solution to P f within tolerance γ,
where γ = /(4D).
2. Suppose we have a tree-decomposition of G f of width W ≥ 1, say. Then by Theorem 5.2.1, there is a tree-decomposition of the intersection network of P f of width≤ O(∆W ).
3. The number of variables in P f , by Lemma 5.4.5, is n + O(D∆).
4. Let ρ be the maximum degree of any of the arc polynomials in the initial NPO (alsothe maximum degree of arc polynomials in P f ). We apply Theorem 4.1.2 to find aCHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

83

solution to P f as per item 1 above. This step will produce an LP of sizeO (2ρ/γ)O(∆W )+1 (n + D∆) log(Dρ/γ) = O (Dρ/)O(∆W )+1 n log(ρ/) .
Thus the construction almost yields Theorem 5.1.3. The key issue is the parameter W inthis estimation. This is the width of a tree-decomposition of the final network G f rather
than the initial network G. We are thus faced with a purely graph-theoretic issue: given an
arbitrary network, is there a sequence of complete node splittings that render a network withmaximum degree at most three, and without a large increase in tree-width?
This is precisely what will be done in the next section. We will prove that if the initial
network has a tree-decomposition of width ω then there is a (polynomial-time constructible)
network G f with a tree-decomposition of width O(ω). Combining this fact with the above
observations yields a proof of Theorem 5.1.3.

5.4.4

Constructing reformulations on small tree-width networks

Consider the example given in Figure 5.6. Here the vertices of the graph shown in (a) arearranged into k > 1 columns. The odd-numbered columns have k vertices each, which induce
a path, while the even-numbered columns have a single vertex which is adjacent to all vertices
in the preceding and following columns. It can be shown that this graph has tree-width 3.
In (b) we show the outcome after splitting vertices so that the maximum degree is four. Thissecond graph has tree-width k, and further splitting the degree-four vertices will not change
this fact.
In contrast to this situation, suppose that we split the graph in Figure 5.6(a) in two steps
as shown in Figure 5.7. The tree-width of the final graph is also 3. The difference between
Figures 5.6 and 5.7 is explained by the fact that the splitting initiated by the “first step” inFigure 5.7 exploits the tree-decomposition of width 3 of the graph Figure 5.6(a).
Next we turn to a formal approach that produces the desired outcome in the general
setting. Given a graph G, a simplification of G will be graph Ḡ obtained by a sequenceof complete vertex splittings, such that the maximum degree of a vertex in Ḡ is ≤ 3. The
following Lemma will show how to obtain a simplification of a graph via a vertex splittingthat maintains tree-width up to a constant factor. In the proof, the trees T̂u that yield theCHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

84

1
1

···

2

···

2
..
.

2

..
.

..
.

···
4

−→

k−1

..
.

..
.

..
.

1

2

3

1

3

k

..
.

k−1

k

···

k

k

..
.

(a)

(b)

Figure 5.6: Incorrect vertex splitting.

1
2
..
.

..
.

···

.. −→
.

k−1

..
.

k
1

3
First step

..
.

..
.

..
.

..
.

k
Second step

Figure 5.7: Correct vertex splitting of graph in Figure 5.6(a).

..
.

..
.

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

85

general splitting procedure described in Section 5.4.4 will be explicitly defined.
Lemma 5.4.6 Let G be a connected undirected graph and (T, Q) a tree-decomposition of G
of width Z. Then there is a simplification Ḡ of G and a tree-decomposition (T̄ , Q̄) of Ḡ ofwidth at most 2Z + 1.
Proof: We first modify (T, Q) in a sequence of steps.

Step 1. For any edge e = {u, v} ∈ E(G), choose an arbitrary t ∈ V (T ) with e ⊆ Qt . Then
we modify T by adding to T a new vertex, te and the edge {te , t}. Further, we set Qte = {u, v}.

Step 2. Without loss of generality, every vertex of T has degree at most 3. To attain thiscondition, consider any t ∈ V (T ) with δT (t) = {s1 , . . . , sd } (say) where d > 3. Then we alter
T by replacing t with two adjacent vertices t1 and t2 , such that t1 is also adjacent to s1 ands2 and t2 is adjacent to s3 , . . . , sd . Finally, we set Qt1 = Qt2 = Qt . Continuing inductively
we will attain the desired condition.

Step 3. For any vertex u ∈ V (G) let Tu be the subtree of T consisting of vertices t withu ∈ Qt , and T̆u be the subtree of Tu that spans {te : e ∈ δG (u)} (which is a subset of the leaves
of Tu ). Then we modify (T, Q) by replacing Tu with T̆u , yielding a new tree-decompositionof same or smaller width. In other words, without loss of generality every leaf of Tu is of theform te for some e ∈ δG (u).

We can now describe our vertex splitting scheme. Consider u ∈ V (G) with degG (u) > 3.
We say that a vertex of Tu is blue if it is either a leaf or of degree three in Tu . Now we form
the tree T̂u whose vertex-set is the set of blue vertices of Tu , and whose edge-set is obtained
as follows. By construction, E(Tu ) can be partitioned into a set of paths whose endpoints
are blue and which contain no other blue vertices. For each such path, with endpoints a andb (say), the tree T̂u will contain the edge {a, b} (in other words, Tu can be obtained from T̂u
by subdividing some edges and so Tu and T̂u are topologically equivalent). Note that T̂u has
degG (u) leaves, each internal vertex with degree 3, and for each edge {u, v} ∈ E(G) there is

CHAPTER 5. NETWORK POLYNOMIAL OPTIMIZATION

86

one pendant edge, as needed.
Let Ğ be the graph obtained by the complete splitting of u using T̂u . For each internalvertex t ∈ V (T̂u ) we name ut the corresponding new vertex in Ğ, to emphasize that each
non-leaf vertex in T̂u will create a vertex in the complete splitting (recall that the leaves inT̂u will correspond to the neighbors of u). This operation does not change the degree of any
vertex v ∈ V (G) with v 6= u. The eventual graph Ḡ in the proof will be obtained by applying
a complete splitting of this type at every vertex of degree > 3 in G.
Returning to Ğ, we construct, for it, a tree decomposition (T, Q̆) as follows. First, let us
regard the tree Tu as rooted at some internal blue vertex r(u). For a vertex t ∈ V (Tu ) let
Au (t) be the closest blue ancestor of t in Tu ; we write Au (r(u)) = r(u). Then, for t ∈ V (T ),we set(Qt − u) ∪ {uAu (t) },if t ∈ V (Tu ) and degTu (t) = 2 ,


Q̆t =
(Qt − u) ∪ {ut , uAu (t) }, if t ∈ V (Tu ) and degTu (t) 6= 2,



 Q,
if t ∈/ V (Tu ).
t(5.13)

Now we argue that (T, Q̆) is a tree-decomposition of Ğ. To see this, note that if t ∈ V (T̂u )then ut ∈ Q̆s iff s = t or s is a child of t in T̂u , thus the endpoints of any edge {ut , us },where t is the parent of s, will be contained in Q̆s . Further, for any edge of e = {u, v} of G,
by Step 3 above there will be a leaf te of Tu such that the edge {te , Au (te )} ∈ E(T̂u ). Thisecorresponds to a pendant edge {uAu (t ) , v} ∈ E(Ğ) and by construction both v ∈ Q̆te andeuAu (t ) ∈ Q̆te . The fact that every vertex in Ğ induces a connected subgraph in T can beeasily verified. This completes the argument (T, Q̆) is a tree-decomposition of Ğ.
Notice that for v ∈ V (G) with v 6= u, the subtree Tv is the same in (T, Q) and (T, Q̆),and that Ğ has one less node of degree greater than three than G. Thus, iteratively applying
the complete splitting of every vertex of Ğ of degree greater than three, and modifying thetree-decomposition as in (5.13) will produce a tree-decomposition (T, Q̄) of the final graphḠ.
By construction, for each t ∈ V (T ) we obtain Q̄t from Qt by replacing each elementwith (at most) two new elements. Thus, since |Qt | ≤ Z + 1, the width of (T, Q̆) is at most2(Z + 1) − 1.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

87

Chapter 6

Cutting planes for Polynomial
Optimization
In the previous chapter we introduced a theoretical framework, along with tractability results, involving linear programming techniques used to approximate polynomial optimization(PO) problems. Tractability was parametrized by tree-width, a graph theoretical measureof structured sparsity. In this chapter we follow on the idea of using linear programming toapproximate PO problems using the well-studied technique of cutting planes, widely used in
the context of mixed-integer programming. As before, we will be concerned on PO problems
which we recall for convenience of the reader:(PO) : min cT x
subject to :fi (x) ≥ 0(4.1a)1≤i≤mx ∈ {0, 1}p × [0, 1]n−p ,(4.1b)(4.1c)

where each fi is a polynomial. In a different, or complementary, spirit from Chapter 4 we will
not assume structured sparsity on PO. Instead, we will aim at developing general techniquesthat can be computationally efficient and stable. To this purpose, we discuss 2 families of
cutting planes for PO problems: digitization and intersection cuts.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

88

Digitization cuts make use of the maturity of Mixed-Integer Programming (MIP) solvers inorder to find valid cuts to PO problems. Using the digitization technique already introducedin Section 4.2, and used there to obtain theoretical tractability results, we will devise aseparation procedure for PO that uses MIP solvers as a sub-routine. As for Intersection cuts,
we will provide a review and details on the computational aspects of the cuts introduced in[23]. There, a generalization of classical Intersection cuts for IP is developed, along with apractical procedure for finding such cuts in the PO context.

6.1

Digitization Cuts

6.1.1

Approximations of Polynomial Optimization revisitedGiven a PO problem, in Section 4.2 we discussed the following pure binary approximation:(GB(γ)) :

min

p
X

cj xj +

j=1

s.t.

fi,α α∈I(i)

xj =

Y

Lγ
X


2−h zj,h 

h=1

n
Y

xj



Lγ
X


j=p+1

j∈Z(i,α)

Lγ
X

cj 

j=p+1


X



n
X

αj 
−h

2

zj,h   ≥ −kfi k1 , 1 ≤ i ≤ m(4.5a)

h=1

2−h zj,h , p + 1 ≤ j ≤ n(4.5b)

h=1

xj ∈ {0, 1}, 1 ≤ j ≤ p(4.5c)

zj,h ∈ {0, 1}, p + 1 ≤ j ≤ n, 1 ≤ h ≤ Lγ

(4.5d)

where  = 1 − (1 − γ)ρ . This approximation can be used to obtain solutions to PO that are-scaled feasible (Definition 4.1.1) and it is based on the fact that, for any 0 ≤ r ≤ 1, given
0 < γ < 1 we can approximate r as a sum of inverse powers of 2:
Lγ
X
h=1

−h

2

zh ≤

r

≤

Lγ
X

2

−h

−Lγ

zh + 2h=1

for some zh ∈ {0, 1}, 1 ≤ h ≤ Lγ and where.
Lγ = dlog2 γ −1 e.

≤

Lγ
X
h=1

2−h zh + γ ≤ 1,

(6.1)

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

89

What GB(γ) effectively does is, starting from PO, it replaces each continuous variable
by a digitized version of it, given by the leftmost expression in (6.1). Using the binary natureof the variables, all non-linearities in GB(γ) can be expressed as linear expressions in a lifted
space. This follows from the fact that, given x1 , x2 ∈ {0, 1}:
y = x1 · x2 ⇐⇒ max{x1 + x2 − 1, 0} ≤ y ≤ min{x1 , x2 }.
Thus, with a slight abuse of notation, we will consider GB(γ) as a binary linear problem. The
following result was proved in Section 4.2, which justifies the usefulness of this approximation.
Lemma 4.2.2(a) Suppose x̄ is feasible for PO. Then there is a feasible solution (x̃, z̃) for

GB(γ) with objective value at most cT x̄ + kck1 .
(b) Suppose (x̂, ẑ) is feasible for GB(γ). Then x̂ is -scaled feasible for (4.1b) and


Lγ
p
n
X
X
X
cT x̂ =
cj x̂j +2−h ẑj,h  .
cj 
j=1

j=p+1

h=1

In this chapter we will aim at obtaining valid inequalities to PO from GB(γ). However,even though Lemma 4.2.2 indicates that GB(γ) acts as a form of relaxation, a valid inequality
for GB(γ) might not be valid to PO. We justify this in the next example.
Example 6.1.1 Consider the following simple feasible set for PO
5
S = (x1 , x2 ) : x21 + x22 ≤ , (x1 , x2 ) ∈ [0, 1]2 .
4Choose γ = 1/4, thus Lγ = 2 and  = 0.4375. Formulating GB(γ) renders the following
feasible set:11
z1,1 + z1,2
2
4

2


+

1
1
z2,1 + z2,2
2
4

2

5 13
+  = 2.671875
441
1x1 = z1,1 + z1,2241
1x2 = z2,1 + z2,22
4

(6.2b)

zi,j ∈ {0, 1}.

(6.2d)≤

Constraint (6.2a) is redundant, hence the projection onto x-space of (6.2) is simply
1 1 3 20
S = 0, , ,.
4 2 4(6.2a)(6.2c)

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

90

p
p
p
On the other side, note that ( 5/8, 5/8) ∈ S. However, as 5/8 ≈ 0.79,pp( 5/8, 5/8) 6∈ conv(S 0 ).
This shows that a valid inequality for conv(S 0 ) (e.g, x1 ≤ 3/4) might be invalid for S.
The issue with GB(γ) is that, even though it provides an approximation with provableguarantees, it only uses one side of approximation (6.1). Therefore, there might be feasible
points of PO outside the convex hull of the feasible region of GB(γ). We fix this problem
by introducing a more elaborate version of this approximation, which considers both of theleftmost inequalities in (6.1) simultaneously. We call this new approximation GB(γ)+ :(GB(γ)+ ) :

min

p
X

n
X

cj x j +

j=1

s.t.

Lγ +1



cj 

X

dh zj,h 

j=p+1


X



fi,α α∈I(i)

Y

Lγ +1

n
Y

xj

j∈Z(i,α)

h=1



X


j=p+1

αj 
dh zj,h   ≥ −kfi k1 , 1 ≤ i ≤ m (6.3a)

h=1

Lγ +1

xj =

X

dh zj,h , p + 1 ≤ j ≤ n(6.3b)

h=1

xj ∈ {0, 1}, 1 ≤ j ≤ p(6.3c)

zj,h ∈ {0, 1}, p + 1 ≤ j ≤ n, 1 ≤ h ≤ Lγ + 1(6.3d)

where dh = 2−h for h = 1, . . . , Lγ and dLγ +1 = 2−Lγ . GB(γ)+ includes one extra digit in
the digitization of each continuous variable, which yields a relaxation with better properties
than GB(γ). As a starting point, we can easily prove the same result as in Lemma 4.2.2 forGB(γ)+ :Lemma 6.1.2(a) Suppose x̄ is feasible for PO. Then there is a feasible solution (x̃, z̃) for GB(γ)+ with
objective value at most cT x̄ + kck1 .
(b) Suppose (x̂, ẑ) is feasible for GB(γ)+ . Then x̂ is -scaled feasible for (4.1b) and


Lγ +1
p
n
X
X
X
cT x̂ =
cj x̂j +cj 
dh ẑj,h  .
j=1

j=p+1

h=1

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

91

Proof: For part (a), we use the same proof as in Lemma 4.2.2 and set the zj,Lγ +1 variables
to 0. Part (b) is direct.

The importance of the extra digit included in GB(γ)+ will become evident in the next
Lemma, and it will justify why GB(γ)+ can be used to derive valid cuts for PO.
Lemma 6.1.3 Suppose x̃ is feasible for PO. Then there exists (xi , z i ) i = 0, . . . , 2n−p − 1
feasible solutions for GB(γ)+ such thatx̃ ∈ conv({xi | i = 0, . . . , 2n−p − 1}).
Proof: We follow a similar procedure as in the proof of Lemma 4.2.2. We start by constructing (x0 , z 0 ) to fix ideas. For each j ∈ {1, . . . , p} define x0j = x̃j , for j ∈ {p + 1, . . . , n} choose
0 h = 1, . . . , L so as to attain the approximation for x̃ as in (6.1), and set z 0zj,hγ
jj,Lγ +1 = 0.

We define x0j for j ∈ {p + 1, . . . , n} from z 0 according to (6.3b). From the analysis in Lemma
4.2.2 we know (x0 , z 0 ) is feasible for GB(γ)+ . Building (x0 , z 0 ) can be seen as a “rounding”
of x̃ to the nearest smaller inverse power of 2, but we can also “round up” using the last
digit.
i
For i = 1, . . . , 2n−p − 1 we define (xi , z i ) the same way as (x0 , z 0 ), except for zp+j,L,
γ +1

j = 1, . . . , n − p, which we set to 1 if and only if the j-th digit of the binary encoding of i is
1. The idea is to have every possible combination of 0’s and 1’s in the last digit. From (6.1),
for every i:
LγXi
2−h zp+j,h

≤

x̃p+j

≤

Lγ
X

i
2−h zp+j,h+ 2−Lγ .

h=1

h=1

We claim that for every i:

 αj

 αj
n−p
n−p
n−p
γ +1
γ +1
Y LX
Y αj
Y LX
i
i

 − ≤
 + .
dh zp+j,hx̃p+j ≤
dh zp+j,h
j=1

j=1

h=1

j=1

(6.4)

h=1

For each i we denote Yi ⊆ {1, . . . , n − p} the digits with 1 in the binary encoding of i.
i
Similarly we denote Ni the digits with 0 in the binary encoding of i. Thus zp+j,L= 1 iff
γ +1

j ∈ Yi . From the proof of Lemma 4.2.2, we have

Y
j∈Ni

α

j
x̄p+j
≤

Y

Lγ +1

X


j∈Ni

h=1

 αj
i

dh zp+j,h+ .

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

92

and ∀j ∈ Yi
Lγ +1

X

x̃p+j ≤

idh zp+j,h
≤ 1.

h=1

This implies

n−p
Y


α

j
x̃p+j
=

j=1


j∈Ni

j 
x̄p+j

j∈Yi

 αj



Lγ +1

Y

X


j∈Ni
n−p
Y

≤

α

Y

j 
x̄p+j


≤



α

Y

i

dh zp+j,h

+ 

X

αj

Lγ +1

X


i

dh zp+j,h

h=1

αj

Lγ +1




Y
j∈Yi

h=1



j=1



i

dh zp+j,h

+ ,

h=1

which proves the right-hand side inequality of (6.4). To obtain the left-hand side inequality, we apply Lemma 4.2.1 for the indices j ∈ Yi .
PLγ +1ih=1 dh zp+j,h − x̃p+j ≥ 0. Clearly, since j ∈ Yi
−Lγ

vj = 2+Lγ
X

We identify uj = x̃p+j and vj =i− x̃p+j ≤ 2−Lγ ≤ γ
dh zp+j,hh=1

therefore,

 αj



Lγ +1

Y

X


j∈Yi

i

dh zp+j,h

−

Y

α

j
x̄p+j
≤1−

j∈Yi

h=1

Y

(1 − vj )αj

j∈Yi

≤1−

Y

(1 − γ)αjj∈Yi≤ 1 − (1 − γ)ρ = .
This implies
n−p
Y
j=1


α

j
+ = 
x̃p+j


Y

α

j 
x̄p+j

j∈Ni

α

j 
x̄p+j
+

j∈Yi


≥ 


Y


Y

j∈Ni

α

j 
x̄p+j


Y
j∈Yi

α

j
x̄p+j
+ 

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

≥ 


Y


j∈Ni

=

n−p
Y

Lγ
X



i
 
dh zp+j,h

Lγ +1

X


i
 
dh zp+j,h

h=1

 αj

Lγ +1

X

αj 


Y

j∈Yi

h=1


j=1

αj  93i
dh zp+j,h,h=1

which completes the proof of (6.4). It immediately follows that ∀i ∈ {0, . . . , 2n−p − 1} thevector (xi , z i ) is feasible for GB(γ)+ .
i
By definition, for each fixed j ∈ {1, . . . , n − p} all zp+j,h, h = 1, . . . , Lγ are equal. The

difference only exists in the last digit. Thus we can define
Cj =

Lγ
X

i
2−h zp+j,h.

h=1

Now define the polyhedronPγ (x̃) = {x ∈ [0, 1]n | xj = x̃j j = 1, . . . , p, Ck ≤ xp+k ≤ Ck + 2−Lγ k = 1, . . . n − p}.
We clearly see that xi i = 0, . . . , 2n−p − 1 is the set of all extreme points in Pγ (x̃), and sincex̃ ∈ Pγ (x̃), the convex combination claim follows.

In this case, we can also provide a tighter result than the one in Corollary 4.2.3:
Corollary 6.1.4 Let P ∗ be the optimal value of problem PO, and let (x̂(γ), ẑ(γ)) be optimalfor GB(γ)+ . Then cT x̂(γ) ≤ P ∗ and x̂(γ) is -scaled feasible for (4.1b).
Proof: Follows directly from Lemma 6.1.3.

Example 6.1.5 Let us continue with Example 6.1.1 to illustrate how GB(γ)+ avoids theaforementioned issue of GB(γ). In this example,

5
S = (x1 , x2 ) : x21 + x22 ≤ , (x1 , x2 ) ∈ [0, 1]2 .
4
We choose γ = 1/4 as before, thus Lγ = 2 and  = 0.4375. Instead of the feasible set definedCHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATIONby (6.2), GB(γ)+ defines the following one:
2 
2

1
1
1
1
1
1
z1,1 + z1,2 + z1,3 +z2,1 + z2,2 + + z2,3 ≤ 2.671875
2
4
4
2
4
4
1
1
1
x1 = z1,1 + z1,2 + z1,3
2
4
4
1
1
1
x2 = z2,1 + z2,2 + z2,3
2
4
4
zi,j ∈ {0, 1}.

94(6.5a)(6.5b)(6.5c)(6.5d)As before, constraint (6.5a) is redundant. We can discard it and see that the feasible set
defined in the x-space is now
00

S =



2
1 1 3
0, , , , 1 .
4 2 4Clearly, conv(S 00 ) = [0, 1]2 and thus(pp
5/8, 5/8) ∈ conv(S 00 ).

Even though these examples used inequalities that become redundant after the approximation is introduced, this might not be the case in general. This choice was made only tosimplify the discussion and argue the need of Lemma 6.1.3 even in easy cases.

6.1.2

Why we need to approximate

A reader may wonder why or if “exact” feasibility (or optimality) for PO cannot be guaranteed. From a trivial perspective, we point out that there exist simple instances of PO (infact convex, quadratically constrained problems) where all feasible solutions have irrational
coordinates. Should that be the case, if any algorithm outputs an explicit numerical solution
in finite time, such a solution will be infeasible.
Example 6.1.6 A simple 1-dimensional example evidencing the need for approximation is
given by the following problem:
max
subject to:x
12x ∈ [0, 1]
x2 ≤where the optimal solution is an irrational number.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION.
.
.
.. ... ... ... ... ... ..
.. .. .. . . . .. .. ..
.. .. .
. .. ..
..

95

Figure 6.1: Graphical representation of Digitization and the corresponding convex hull
Example 6.1.7 A more elaborate example is given by the following region in R2 :


1
2
2
2
(x, y) : x + y ≤ 1, xy ≥ , (x, y) ∈ [0, 1] .
2√√
The set is expressed only using rational data, and the unique feasible point is (1/ 2, 1/ 2).
Another question one might ask, is whether it is possible to achieve constraint violations
that are at most , independent of the 1-norm of the constraints. Intuitively, this should
not be possible since in such case one could augment the coefficients on the constraints andobtain arbitrarily high-quality solutions. For a more rigorous analysis of this see Section C.1of the Appendix, where we argue why the scaled feasibility feature is needed, even in thebounded tree-width case.

6.1.3

Digitization-based Cuts

In Figure 6.1 we show a graphical representation on what GB(γ)+ is doing with respect to
PO: it replaces each continuous variable with a weighted sum of binary variables. Moreover,by Lemma 6.1.3, the digitization is such that the convex hull contains the original feasible
set. The cuts we will discuss in this section will aim to separate from such convex hull. In
what follows, we assume  (and thus γ) are predetermined values.
To make the analysis precise, we label the different feasible sets we are considering theCHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION.
.
.
.. ... ... ... ... ... ..
.. .. .. . . . .. .. ..
.. .. .
. .. ..
.
..

96

Figure 6.2: Graphical representation of Digitization Cut
following way:F = {x ∈ {0, 1}p × [0, 1]n−p | x is feasible for PO}
F = {x ∈ {0, 1}p × [0, 1]n−p | x is -scaled feasible for PO}
Fγ = {(x, z) ∈ {0, 1}p × [0, 1]n−p × {0, 1}(n−p)Lγ | (x, z) is feasible for GB(γ)+ }Lemma 6.1.8 conv(F ) ⊆ conv(projx (Fγ )) ⊆ conv(F )Proof: This is direct from Lemma 6.1.3 and from the definition of GB(γ)+ .

Our goal will be to use Fγ to produce -feasible solutions to PO by separating from its
convex hull. We use the following algorithmic template:
1. Find a polyhedral relaxation of F and an extreme point solution x̄.
2. If x̄ ∈ F , then STOP.
3. Otherwise, x̄ 6∈ conv(projx (Fγ )) (by Lemma 6.1.8). We look for a cutting plane separating x̄ from conv(projx (Fγ )) and add it to the polyhedral relaxation of F . See Figure
6.2, where the black point represents x̄.
Separating from Fγ can be computationally challenging, as it can be as hard as optimizing
over GB(γ)+ in general. Thus, in order to speed up the cut-finding procedure, we digitize
only a subset of the continuous variables.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

97

Let J ⊆ {p + 1, . . . , n}, we define
I(J) = {i | supp(fi ) ⊆ {1, . . . , p} ∪ J}
i.e, the indices of the constraints that involve only continuous variables contained in J andbinary variables. Given I ⊆ {1, . . . , m}, we denote as Poly(I) a polyhedral relaxation (potentially using lifted variables) of the set{x ∈ {0, 1}p × [0, 1]n−p | fi (x) ≥ 0, i ∈ I}.
Given η ∈ Rn , we definenX(GB(γ, J, η)+ ): max

ηj xj

j=1

subject to: x ∈ Poly(I(J)c )XYfi,α α∈I(i)(6.6a)xj

j∈Z(i,α)

Y

Lγ +1

X


j∈J

αj 
dh zj,h   ≥ −kfi k1 , i ∈ I(J)

(6.6b)

h=1

zj,h ∈ {0, 1}, j ∈ J, 1 ≤ h ≤ Lγ + 1(6.6c)

Lγ +1

xj =

X

dh zj,h , j ∈ J(6.6d)

h=1

xj ∈ [0, 1], j 6∈ J

(6.6e)

xj ∈ {0, 1}, 1 ≤ j ≤ p(6.6f)

In simple words, we are just digitizing the continuous variables in J instead of all ofthem. As for the constraints, we are only imposing constraints in I(J) explicitly, as thesecan be represented using only binary variables and thus can be reformulated using linear
expressions. The rest of the constraints, i.e, the ones with indices in I(J)c are relaxed using
some polyhedron. This renders problem GB(γ, J, η)+ a mixed-integer linear problem withp + Lγ |J| binary variables and n − p continuous variables.
Following the same notation as before, we denote as Fγ (J) the feasible region of problem
GB(γ, J, η)+ . ClearlyFγ = Fγ ({p + 1, . . . , n})and using the same argument as in Lemma 6.1.3 we seeconv(F ) ⊆ conv(projx (Fγ (J))) ∀J ⊆ {p + 1, . . . , n}.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

98

Therefore, every valid inequality derived from Fγ (J) will be valid for F , and if we keep Jof moderate size, GB(γ, J, η)+ will be a MIP with a small number of binary variables. Thiscan make finding such valid inequality a manageable task.

Let x̄ 6∈ F . In order compute a cut we solve the following problem
SEPγ (J):max δ

subject to:αT (x̄ − x) ≥ δ∀x ∈ projx (Fγ (J))

kαk ≤ 1δ≥0with k · k some norm. We solve problem SEPγ (J) using a row generation approach, i.e, at
any given step we count with a finite number of solutions {xk }k∈K ⊆ projx (Fγ (J)) for some
K finite. Then we solve the following relaxation of SEPγ (J):SEPγrel (J):max δ

subject to:αT (x̄ − xk ) ≥ δ∀xk , k ∈ K

kαk ≤ 1
δ ≥ 0.
If the norm k · k is linearly representable (e.g, 1-norm, infinity norm), problem SEPγrel (J) isan LP. Given (α̃, δ̃) optimal for SEPγrel (J), we solve GB(γ, J, α̃)+ . Let (x0 , z 0 ) be its optimal
solution.
• If α̃T x̄ − δ̃ < α̃T x0 , then (α̃, δ̃) is not feasible for SEPγ (J), and we add x0 to K in order
to separate (α̃, δ̃).
• If α̃T x̄− δ̃ ≥ α̃T x0 , then by optimality of (x0 , z 0 ) the vector (α̃, δ̃) is feasible for SEPγ (J),and thus optimal.
Using this procedure we can directly look for cutting planes. The performance of it will
depend heavily on the size of J, thus the choice of J should aim for both a strong relaxation,and for efficiency. Below we use a natural choice of J that can leverage these notions if theconstraints are sparse.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

6.1.4

99

Computational Experiments

In the previous section we discussed a cutting plane approach that separates from Fγ (J), asa way of generating valid inequalities for PO. We implemented an iterative cutting plane
method that, in each iteration, chooses a different set J to attempt to separate the current
solution x̄ 6∈ F from Fγ (J). The following considerations were taken in the implementation:• The algorithm is stopped if no cut is found, or if some numerical instability is detected.
• At each iteration, given an extreme point solution x̄, we choose J to be the support ofthe most violated constraint by x̄.
• We used  = 10−12 for the approximation quality.
• When solving SEPγrel (J), the 1-norm was used.
Computations were run on a 4-core laptop computer with an Intel Core i7-5600U CPU
and 8 GB of RAM. The code was written in Python 2.7.13. The linear programming solver
is Gurobi 7.0.1 with default settings.
Test instances are taken from Floudas et al. [46] (available via GLOBALLib [82]). Even
though all cuts can handle arbitrary polynomial terms, for implementation purposes werestricted ourselves to problems with quadratic objective and constraints. Whenever an
objective function is non-linear, we use the usual epigraph reformulation that brings the nonlinearities to the constraints. To evaluate the quality of the cuts obtained, we compare with
[97, 98].
We choose the initial LP relaxation to be the standard RLT relaxation of QCQP, i.e, we
include a new variable representing each bilinear term and add the McCormick estimators
for them (see e.g. [2]). At each subsequent iteration we construct the polyhedra Poly(I(J)c )using the same RLT technique.
The results are summarized in Table 6.1. Let OP T denote the optimal value of an
instance, RLT the optimal value of the standard RLT relaxation, and GLB the objective
value obtained after applying the cutting plane procedure. We use the following metric fordetermining the quality of the solutions:
Gap Closed =GLB − RLT.
OP T − RLT

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

SBL-Gap

SBL-Time

100

Digitization

Instance

Gap Closed

Time (s)Gap Closed

Time (s)Gap Closed

Time (s)Cuts

Ex2 1 599.98%

0.173

77.09%

0.049

98.6%

2.3

8

Ex2 1 699.97%

0.385

99.96%

0.265

77.8%

0.12

3

Ex2 1 788.09%

221.363

69.74%

11.179100%224.545

Ex2 1 899.90%

1.692

99.64%

0.645

9.94%

49.7

100

Ex2 1 993.96%

0.985

93.96%

0.985

37.9%

10.8

3

Ex3 1 1

15.94%

3600

0.35%

0.186

87.76%

63.63

7

Ex5 2 2

0.22%

0.150

0.22%

0.150

100%

0.53

1

Ex5 2 479.31%

68.92773.30%

1.448

94.72%

235.244

Ex5 3 27.27%

245.821

0.0%

2.665

0%

5.3

1

Ex5 4 227.57%

3614

0.52%

0.177

30%

92.02

11

Ex8 4 191.84%

3659

89.05%

3.050

0%

Num

-

Ex9 1 4

0.0%

0.057

0.0%

0.057

39.53%

5.43

36

Ex9 2 168.74%

0.563

68.74%

0.563

0%No Cut-Ex9 2 761.34%

0.899

59.30%

0.468

0%No Cut-Table 6.1: Digitization Cuts Experiments
SBL-Gap shows the values obtained in [97, 98] that achieve the best Gap Closed, regardlessof the time spent in the cutting plane algorithm. SBL-Time shows the values obtained in thesame articles, but restricted to results obtained with a running time in order of seconds.
Table 6.1 shows encouraging results. First of all, only using a few cuts we are able to close
considerable gap in many instances, with moderate running times. This shows potential, as
fewer cuts yield algorithms that are numerically more stable. An interesting case is given by
instance Ex8 4 1, which presents numerical instability while generating the very first cut. In
the last 2 instances, Ex9 2 1 and Ex9 2 7, the algorithm was not able to find a cutting planewith our choice of J.
Overall these results show how useful this approach is. Digitization cuts are able toprovide good solutions in short times, only relying in MIP solvers, and the flexibility on the
choice of J opens a gate for more elaborate heuristics that can improve on these results.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

6.2

101

Intersection Cuts

In this section we review the family of cuts introduced in [23], and show the computational
aspects in detail, in order to analyze the practical performance of them in challenging PO
instances. At the core of this family there is the concept of S-free sets, which we review next.

6.2.1

Review of S-free sets and Intersection cuts

Consider a generic optimization problem
minsubject to:cT x(6.7a)x ∈ S ∩ P,

(6.7b)

where P ⊆ Rn is a polyhedron and S ⊆ Rn is a closed set. In order to address this problem,
we will follow a cutting plane approach as before, i.e, we construct a polyhedral relaxationof S ∩ P , find and extreme point solution x̄, and if x̄ 6∈ S we aim to find (π, π0 ) ∈ Rn+1 suchthat
π T x ≤ π0 ∀x ∈ S ∩ P∧

π T x̄ > π0 .

A typical example is given by integer programming, where S is the integer lattice Zn . In
our case, we will use a set S drawn from an equivalent formulation to PO, and then we will
derive cuts from convex forbidden zones, or S-free sets, which are defined as follows:
Definition 6.2.1 A set C ⊂ Rn is S-free if int(C) ∩ S = ∅ and C is convex.
Remark 6.2.2 For any S-free set C we have S ∩ P ⊆ clconv(P \ int(C)), and so any valid
inequalities for clconv(P \ int(C)) are valid for S ∩ P .
To generate cuts, we will adopt the standard procedure of finding a simplicial cone P 0
containing P and apply Balas’ intersection cut [9] for P 0 \ int(C). For relevant literature
related to this well studied family of cuts see [35, 37, 36]. Balas and Margot [11] propose
strengthened intersection cuts by using a tighter relaxation of P .
Definition 6.2.3 We say P 0 ⊇ P is a simplicial conic relaxation of P when P 0 is a displaced polyhedral cone with apex x̄ and defined by the intersection of n linearly independent
halfspaces.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION
A simplicial cone may be written as follows:


n


X
0
(j)
P = x|x = x̄ +λj r , λ ≥ 0 .



102

(6.8)

j=1

Remark 6.2.4 Any n linearly independent constraints describing P can be used to define a
simplicial conic relaxation, using the corresponding extreme point of P as the apex of P 0 .
Each extreme ray is of the form x̄ + λj r(j) , where each r(j) ∈ Rn is a direction andλj ∈ R+ .
As mentioned before, we will assume x̄ 6∈ S, with the goal of separating x̄ from S. Let
C be an S-free set with x̄ ∈ int(C), and P 0 a simplicial cone relaxation with apex x̄. Eachextreme ray of P 0 must satisfy one of the following:• Either it is entirely contained in C:x̄ + λj r(j) ∈ int(C) ∀λj ≥ 0• Or there is an intersection point with the boundary:∃λ∗j : x̄ + λ∗j r(j) ∈ bd(C)Then the intersection cut is defined as follows:
Definition 6.2.5 Let x̄ 6∈ S be an extreme point of a polyhedral relaxation of S ∩ P , C bean S-free set with x̄ ∈ int(C), and P 0 a simplicial cone relaxation of P with apex x̄. Theintersection cut is the halfspace whose boundary contains each intersection point (given byλ∗j ) and that is parallel to all extreme rays of P 0 fully contained in C.
In [9] a closed-form expression for the cut π T x ≤ π0 is provided. Consider P 0 in inequality
form, P 0 = {x|Ax ≤ b}, where A is an n × n invertible matrix and x̄ = A−1 b. Then thecoefficients of the intersection cut defined in Definition 6.2.5 can be expressed as

π0 =

n
X
i=1

(1/λ∗i )bi − 1, πj =nX
(1/λ∗i )aij .

(6.9)

i=1

Here, 1/λ∗j is treated as zero if the step-length is infinite. Note that, without loss ofgenerality, the rays are assumed to be aligned with the inequality form, in the sense that forCHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

103

each j, −r(j) is the j-th column of the inverse of A. A quick calculation shows the hyperplane
π T x = π0 intersects all intersection points. Let j be such that λ∗j < ∞:
π T x̄ + λ∗j r(j) = π T x̄ + λ∗j π T r(j)=nXπk x̄k + λ∗j

k=1

=
=
=

n
X
k=1
n
X
i=1
n
X

n
X
k=1

πk (A−1 b)k + λ∗j

=

n X
n
X

(j)

(1/λ∗i )aik rk

k=1 i=1(1/λ∗i )nX

aik (A−1 b)k + λ∗j(1/λ∗i )bi − λ∗jnX(1/λ∗i )

i=1

k=1

i=1
n
X

(j)

πk rk

n
X

n
X

(j)

aik rkk=1(1/λ∗i )1{i=j}i=1(1/λ∗i )bi − 1 = π0i=1and clearly the halfspace is parallel to all r(j) such that λ∗j = ∞.
The final step on constructing the intersection cut is to determine the right direction of
the inequality. To this end, let β0 := sgn(π T x̄ − π0 ), where sgn(·) denotes the sign function.
The intersection cut formula is then given by:
β0 (π T x − π0 ) ≤ 0.

(6.10)

The following proposition concludes the argument for the correctness of the cut:
Proposition 6.2.6 Let V := {x | β0 (π T x − π0 ) ≤ 0} be the halfspace defined by the intersection cut. ThenV ⊇ P 0 \ int(C).
Furthermore, if all step lengths are finite, then V ∩ P 0 = conv(P 0 \ int(C)).
For a proof of this statement see [9, Theorem 1].

From this section we see that the bottleneck of constructing intersection cuts efficiently
is given by determining the step-lengths λ∗ . These will depend on the S-free set constructed
for each x̄, and in general they might not be easily computable.
In the next section we will review how we can construct S-free sets for PO problems, aswell as showing how to compute the step-lengths for each one of the proposed S-free sets.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

6.2.2

104

S-free sets for Polynomial Optimization

We now turn to the construction of S-free sets for PO problems. We begin by casting PO
as a problem of the type (6.7) so as to identify S precisely.
6.2.2.1

Moment-based reformulation for PO problems

We will make use of the moment-based reformulation of PO problems (see [67, 70]). Let ρ
be the maximum degree among polynomials fi defined in (4.1b). Let mr (x) = [1, x1 , ..., xn ,x1 x2 , ..., x2n , ...xrn ] be a vector of all monomials up to degree r. Then it is easy to see thatwhenever 2r ≥ ρ any polynomial fi can be written asfi (x) = mTr (x)Ai mr (x)where Ai is an appropriately defined symmetric matrix obtained from fi . Then,fi (x) ≥ 0 ⇐⇒ mTr (x)Ai mr (x) ≥ 0 ⇐⇒ hmr (x)mr (x)T , Ai i
and we can apply this transformation to PO to obtain the following reformulation:(LPO) : min hA0 , Xi
subject to: hAi , Xi ≥ 0, i = 1, ..., m,
X = mr (x)mTr (x).

(6.11a)(6.11b)Remark 6.2.7 A problem of the form LPO can be obtained even when not assuming theobjective function in PO is linear.
Condition (6.11b) can be equivalently stated as X  0, rank(X) = 1, and consistency
constraints among the expressions that represent the same monomial. Constraint X  0 is
a convex constraint, thus the non-convexity of the original PO problem is now absorbed in
the rank constraint. Dropping the latter yields the standard semidefinite relaxation of PO[102]. In some cases, (e.g. [78, 67, 71, 90]) it can be proven that there exists sufficiently large
r to ensure an exact relaxation. However, there is a rapid increase in the size of LPO and its
semidefinite relaxation with respect to r.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

105

Note that LPO was written in a convenient way: constraints (6.11a) are polyhedral
constraints, and constraints (6.11b) represent a closed set. Hence, we are left with a naturaldefinition of S:
S = {xxT : x ∈ RN (r) }(6.12)

with N (r) appropriately defined. Furthermore, we let P correspond to the linear constraints(6.11a) together with the consistency constraints and X11 = 1. This brings a PO problem
to a problem of the form (6.7), and justifies the study of outer-product-free sets, i.e, sets ofmatrices that are not symmetric outer-products of vectors. In what follows, we will drop the
dependency on r in the definition N for alleviating notation.
For X ∈ SN ×N consider the vectorized matrix vec(X) ∈ RN (N +1)/2 , with entries fromthe upper triangular part of the matrix. We will work in this vector space, thus all notionsof the interior, convexity, and so forth are with respect to it. However, for simplicity of the
discussion, we drop the explicit vectorization where obvious, as translating the results inmatrix form to vector form is straight-forward.
6.2.2.2

Oracle-based S-free sets for PO

When there is an oracle that can provide the distance of any point to the set S, a simple
S-free set can be constructed. Say x̄ 6∈ S, and let d > 0 be the distance from x̄ to S, thenclearly the ball of center x̄ and radius d, denoted Boracle (x̄, d) is S-free and an intersection
cut can be easily derived from it. In [23], an in-depth study of these oracle-based cuts is
performed, providing theoretical guarantees when using them in an iterative fashion.
In the context of PO, given a matrix X̄ that is not an outer-product, we would liketo compute the distance of X̄ to the nearest symmetric outer product, or equivalently thenearest positive semidefinite matrix with rank at most one. Fortunately, this distance can beobtained using a low-rank approximation scheme. Given an integer r > 0, we define(PMA) :minYkX̄ − Y kF

s.t. rank(Y ) ≤ r,
Y  0.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

106

In what follows, we denote as {λi }ni=1 the eigenvalues of X̄, in decreasing order, and{di }ni=1 the corresponding eigenvectors. Dax [40] proves the following:
Theorem 6.2.8 (Dax’s Theorem [40]) Let k be the number of nonnegative eigenvalues ofPmin{k,r}
X̄. For 1 ≤ r ≤ N − 1, an optimal solution to PMA is given by Y = i=1
λi di dTi .
When X̄ is not negative semidefinite, the solution from Dax’s theorem coincides with
Eckart-Young-Mirsky [83, 43] solution to PMA without the positive semidefinite constraint.
With this, we can construct an outer-product-free ball: B(X̄, kX̄kF ),if X̄ is NSD,
Boracle (X̄) :=P B(X̄, k n λ d dT k ), otherwise.
i=2 i i i F
And directly from Dax’s Theorem:
Corollary 6.2.9 Boracle (X̄) is outer-product-free.
This construction centers the oracle ball around around X̄, however, for LPO we canobtain a larger ball containing the original one by re-centering it. Consider a ball B(X, r).
Let s > 0 and suppose Q is in the boundary of the ball. We construct the “shifted” ballB(Q + (s/r)(X − Q), s). This ball has radius s and its center is located on the ray that goes
through X emanating from Q, and moreover, whenever s > r the shifted ball contains theoriginal ball.
This allows us to construct Bshift (X̄, s), a shifted oracle ball. We will use the nearest
symmetric outer product as the boundary point in our construction: B(sX̄/kX̄kF , s),if X̄ is NSD,
Bshift (X̄, s) :=s B λ1 d1 dT +(X̄ − λ1 d1 dT1 ) , s , otherwise.
1kX̄−λ d dT k
1 1 1

F

As the reader might notice, the definition depends on the value of s. This will determine
when the shifted ball is outer-product-free. In [23], the following is proved:
Proposition 6.2.10 Suppose X̄ is not an outer product. If λ2 ≤ 0 then Bshift (X̄, kX̄kF + )is outer-product-free and strictly contains Boracle (X̄) for  > 0. If 0 < λ2 < λ1 , then for
Pn
λ1 Pn
T
T
i=2 λi di di F < s ≤ λ2i=2 λi di di F , Bshift (X̄, s) is outer-product-free and strictly
contains Boracle (X̄).

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

107

Additionally, the following Theorem allows to further improve outer-product-free sets wehave discussed so far.
Theorem 6.2.11 Let C ⊂ SN ×N be a full-dimensional outer-product-free set. Then clcone(C)is outer-product-free.
This implies that, provided s is chosen appropriately so Bshift (X̄, s) is outer-product-free,
the closure of the conic hull, clcone(Bshift (X̄, s)), is also outer-product-free.
6.2.2.3

Maximal S-free sets for PO

While in principle any S-free set can be used to obtain intersection cuts, larger S-free sets
can be useful for generating deeper cuts [35]. Thus, we aim to find maximal S-free sets.
Definition 6.2.12 An S-free set C is maximal if V 6⊃ C for all S-free V .
The following results of [23] yield two families of maximal outer-product-free sets.
Theorem 6.2.13 A halfspace hA, Xi ≥ 0 is maximal outer-product-free iff A is negative
semidefinite.
Remark 6.2.14 It is not hard to see that whenever the S-free set is of the form hA, Xi ≥ 0,then the intersection cut associated to it is given by hA, Xi ≤ 0.
Theorem 6.2.13 gives an interesting interpretation to the well known outer-approximationof the PSD cone, which can be used to solve SDPs. Say X̄ 6 0, then there exists c such thatcT X̄c < 0
and the set{X ∈ SN ×N | h−ccT , Xi ≥ 0}
is maximal outer-product free. By Remark 6.2.14 the intersection cut will be given byh−ccT , Xi ≤ 0which will cut off X̄.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

108

Remark 6.2.15 We can assume the matrix A in Theorem 6.2.13 is of rank 1. If not, letr > 1 be its rank, λi < 0 i = 1, . . . , r its eigenvalues and di i = 1, . . . , r its eigenvectors.
Then, all sets
{X ∈ SN ×N | h−di dTi , Xi ≥ 0}

i = 1, . . . , rwill be outer-product free, and they will generate the intersection cutsh−di dTi , Xi ≤ 0.
These cuts imply* rX+
λi di dTi , X

≤0

i=1

which is the cut hA, Xi ≥ 0 induces.
Another family of maximal outer-product-free sets is defined next.
Definition 6.2.16 A 2 × 2 PSD cone is of the form C i,j := {X ∈ SN ×N |X[i,j]  0}, for somepair 1 ≤ i 6= j ≤ N .
Remark 6.2.17 The boundary of C i,j is the set of X ∈ SN ×N such that X[i,j] has rank one
and is PSD.
Theorem 6.2.18 Every 2 × 2 PSD cone is maximal outer-product-free.
In [23] it is proven that for N = 2 every maximal outer-product-free sets is either ahalfspace, or the PSD cone (which corresponds to the single 2 × 2 PSD cone in this case).
However, for N ≥ 3 there exist other maximal outer-product-free sets.

6.2.3

Cutting plane procedure

From the previous section we can extract 4 families of intersection cuts, which are summarizedin Table 6.2. Since the 2 × 2 cut will separate any rank two or greater PSD matrix, and theouter approximation cuts will separate NSD or indefinite matrices, then the two families ofcuts can be considered complementary. Bshift provides a stronger cut than Boracle and can beused alone or in combination with the other cuts.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

Cut NameS-free set

Separation Condition

2 × 2 Cut

X[i,j]  0X̄[i,j]  0Outer Approximation Cut

cT Xc ≤ 0

X̄ NSD or indefinite

Oracle Ball Cut

BoracleX̄ not an outer product

Strengthened Oracle Cut

clcone(Bshift )X̄ not an outer product109

Table 6.2: Proposed Intersection Cuts
In order to test the efficacy of these cuts, we construct a pure cutting plane procedurethat makes use of these intersection cuts. From equations (6.10) and (6.9) we observe most of
the complexity of finding the cuts of a given S-free set is to find the associated step-lengths
λ∗ . In the case of the 4 families described above, the step-lengths can be computed using
closed-form expressions. This is one of the strongest features of these families, as finding
the corresponding cuts should not present a heavy computational burden. We specify how
to compute these next. In all cases we assume X̄ is matrix which is not an outer-product,
and as before we denote as λi and di its eigenvalues and eigenvectors, respectively, withthe eigenvalues sorted in decreasing order. We also assume we computed the rays D(k)k = 1, . . . , N of the simplicial cone relaxation in matrix form.
6.2.3.1

Oracle Ball step-length

The outer-product-free set in this case is given by B(X̄, kX̄kF ),if X̄ is NSD,
Boracle (X̄) :=P B(X̄, k n λ d dT k ), otherwise.
i=2 i i i F
Clearly the step-length in all directions will be given by the radius of the ball. SuchP
computation involves computing either kX̄kF or k ni=2 λi di dTi kF .

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION
6.2.3.2

110

Strengthened Oracle Ball step-length

The analysis to obtain the step-lengths associated with clcone(Bshift (X̄, s)) is more complicated, and requires separating into several different cases according to Proposition 6.2.10.
We summarize these cases next. For details we refer the reader to the full-length paper [23].
• If X̄ is NSD, it can be proven that the conic extension of Bshift (X̄, kX̄kF + ) is ahalfspace tangent to the shifted ball, at the origin. The best possible cut is then thehalfspace:hX̄/kXkF , Xi ≤ 0.
Hence, we do not require a step-length computation.
• If X̄ is not NSD, and λ2 is non-positive, then we may use the halfspace that contains
λ1 d1 dT1 on its boundary and that is perpendicular to the vector from X̄ to λ1 d1 dT1 . The
best possible cut is a halfspace:hX̄ − λ1 d1 dT1 , X − λ1 d1 dT1 i ≤ 0.
Again, we do not require a step-length computation.
• If X̄ is not NSD and λ2 is positive, we use the maximum shift of Proposition 6.2.10:PnTs = λλ12
i=2 λi di di F . This gives us a shifted ball with center
XC := λ1 d1 dT1 +λ1(X̄ − λ1 d1 dT1 )λ2

and radiusq :=λ1X̄ − λ1 d1 dT1
λ2

F

.

Given a ray D(k) , a simple computation can distinguish whether it is fully containedon the cone or not. This is determined using the projection of the direction vector onto
the axis of the cone: hD(k) , XC i/kXC kF .
If the step length is finite, using geometric arguments we can prove that the step λ∗krequired to reach the boundary of cone(B(XC , q)) using ray D(k) satisfies the quadratic

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

111

equation aλ2k + bλk + c = 0, wherea := q 2 hD(k) , XC i2 /(kXC k2F − q 2 ) − kZ4 k2F ,b := 2q 2 hX̄, XC ihD(k) , XC i/(kXC k2F − q 2 ) − 2hZ3 , Z4 i,
c := q 2 hX̄, XC i2 /(kXC k2F − q 2 ) − kZ3 k2F ,Z3 := kXC kF X̄ − hX̄, XC iXC /kXC kF ,Z4 := kXC kF D(k) − hD(k) , XC iXC /kXC kF .
λ∗k is equal to the positive root of this quadratic equation.
6.2.3.3

Outer Approximation step-length

From the discussion below Theorem 6.2.13, given a matrix X̄ 6 0, a maximal outer-productfree set can be computed directly using any negative eigenvalue of X̄ (and the corresponding
eigenvector) of X̄. This outer-product-free set will be a halfspace, and the corresponding
intersection cut can be obtained directly. Computing the step-lengths in this case is again
unnecessary, and the intersection cut coefficients computation rely only on the eigenvalue
decomposition of X̄.
6.2.3.4

2 × 2 Cut step-length

Given X̄, in order to find a 2 × 2 PSD cone containing X̄ in its interior we first cycle throughall pairs 1 ≤ i 6= j ≤ N and check whether X̄[i,j]  0 or not. This can be done by checkingX̄ii , X̄jj > 0 and X̄ii X̄jj > X̄ij2 .
If some pair i, j is such that X̄[i,j]  0, then X̄ ∈ int(C i,j ) and the step-lengths can be(k)found using the 2 × 2 principal sub-matrices X̄[i,j] and D[i,j] corresponding to X̄ and someextreme ray direction D(k) , respectively.
(k)

(k)

First suppose D[i,j] 6 0. Then we need to compute λk ≥ 0 such that X̄[i,j] + λk D[i,j] lieson the boundary of the 2 × 2 cone, i.e. its minimum eigenvalue is zero. For a 2 × 2 symmetricmatrix




a b
b c




CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION112

the minimum eigenvalue is given by(a + c −pa2 − 2ac + 4b2 + c2 )/2.
(k)

The eigenvalue is zero when ac = b2 , a + c ≥ 0. Applying this to X̄[i,j] + λk D[i,j] we obtainthe following equation:
(k)

(k)

(k)

(k)

(k)

(k)

((Dij )2 − Dii Djj )λ2k + (2Dij X̄ij − Dii X̄jj − Djj X̄ii )λk + X̄ij2 − X̄ii X̄jj = 0.
The desired step length is given by the greater root of this quadratic equation with respect
to λk (the lesser root sets the maximum eigenvalue to zero).
(k)

(k)

In the case D[i,j]  0, X̄[i,j] + λk D[i,j] is positive semidefinite for all nonnegative λk , whichyields an intersection at infinity.
Remark 6.2.19 If some extreme rays are fully contained in an outer-product-free set, thenthe intersection cut can be strengthened by using negative step lengths λ∗j [42, 15, 14]. The
strengthening procedure involves changing one step-length at a time, which rotates the intersection cut in the axis defined by the other (n − 1) intersection points. See [23] for details on
how we can apply this procedure to our proposed families of outer-product-free sets.

6.2.4

Computational Experiments

We present numerical experiments using a pure cutting plane algorithm implementation of
the cuts described above. These experiments are meant to provide a sense on how useful thecuts are, and not to be used as a complete polynomial optimization solver.
The cutting plane procedure follows the same structure as usual: we find a polyhedralrelaxation of S ∩ P and compute an extreme point optimal solution x̄. If x̄ 6∈ S ∩ P , we addintersection cuts to the relaxation separating x̄, and repeat. In our implementation we used
the following stopping conditions:• Time limit of 600 seconds• Objective value not improving for 10 iterations• Maximum violation of a cut not more than 10−8 . Here, if π T x ≤ π0 is the cut, i.e, avalid inequality such that π T x∗ > π0 , we define the violation as(π T x∗ − π0 )/kπk1

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

113• Linear program becomes numerically unstable
The machine where the tests were made is a 20-core server with an Intel Xeon E5-2687W
v33.10GHz CPU and 264 GB of RAM. Even though these tests were performed in a computerwith high specifications, we confirmed that similar performance can be obtained with a laptopof the same specifications as the Digitization cuts. The code is written in Python 2.7.13 using
the Anaconda 4.3 distribution. The linear programming solver is Gurobi 7.0.1 with default
settings. The code and full experimental data is available at https://github.com/g-munoz/
poly_cuts.
Test instances are taken from two sources. First, as in Digitization cuts, we include problem instances from GLOBALLib [82] that have quadratic objective and constraints. Second,
we include all 99 instances of BoxQP developed by several authors [106, 31]. As before, wechoose the initial LP relaxation to be the standard RLT relaxation of QCQP. Note that, unlike the Digitization cuts, Intersection cuts do not need the epigraph reformulation to bring
non-linearities of the objective function to the constraints. To obtain variable bounds forsome of the GLOBALLib instances we apply a simple bound tightening procedure: minimize/maximize a given variable subject to the RLT relaxation.

Results averaged over all included instances are shown in Table 6.3 for GLOBALLib andTable 6.4 for BoxQP. Cut family indicates one of 5 cut configurations: Oracle Ball cuts (OB),Strengthened Oracle cuts (SO), Outer Approximation cuts (OA), 2x2 cuts together with OA(2x2 + OA), and finally 2x2, OA and SO cuts together (SO + 2x2 + OA). As before, let
OP T denote the optimal value of an instance (with objective minimized), RLT the optimalvalue of the standard RLT relaxation, and GLB the objective value obtained after applying
the cutting plane procedure. Then the optimality gaps (per instance) are

Initial Gap =End Gap =OP T − RLT,|OP T | + 

OP T − GLB,|OP T | + 

Gap Closed =GLB − RLT.
OP T − RLT

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

114

where  := 1 is a parameter to avoid division by zero. The Initial Gap values are the sameirrespective of the cut configuration, so only one entry in the corresponding column is neededon Table 6.3 and Table 6.4. Iters is the number of cutting plane algorithm iterations per
problem instance. Time is the total time in seconds spent by the algorithm. LPTime is thepercentage of total time spent solving linear programs.
Large initial gaps are observed due to our choice of a simple initial relaxation. It is also
important to mention that initial and end gaps averages are negatively affected by a few
instances with high gap, thus we also provide a distribution of outcomes in Table 6.6. Table
6.5 shows a distribution on the closed gaps, which we believe is a more appropriate measure
for comparing averages, as it is normalized to between 0% and 100%.
In the oracle-based cuts (OB and SO), we first observe what was expected: the strengthened oracle cuts produce better results than the basic oracle ball cuts. However, both families
exhibit tailing off behavior, which results in modest closed gaps. The 2×2+ OA configuration
provides substantially better performance compared to OB, SO, and OA, on GLOBALLibinstances, indicating the benefits of separating over the 2 × 2 family of cuts. OA produces
the best results in terms of average end gap and closed gap on BoxQP instances.
The small LPTime values indicate a relatively high percent of time spent in producing
cuts. This was counter-intuitive, as generating the cuts do not require complex computations. A closer inspection showed that the majority of the time was spent in linear algebra
operations (eigenvalue computations, norms, solving linear systems, etc). For these operations we rely on the NumPy library, thus our simple Python implementation is bound tobe as efficient as Python and NumPy, which are mainly built for scripting and not for highperformance computing. We strongly believe using a lower-level language should increase
efficiency drastically.
Per-instance data for 2 × 2+OA can be found in Appendix D; for other configurations werefer the reader to the online supplement mentioned at the beginning of this section. The
appendix also provides data showing that 2×2+OA offers comparable bounds at substantiallyfaster computation times compared to the V2 configuration of Saxena, Bonami, and Lee [97].
However, we emphasize that cuts are complementary and not competitive.

CHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATION

Cut Family

115

Initial Gap

End Gap

Closed Gap# CutsItersTime (s)LPTime (%)

1387.92%

1387.85%

1.00%

16.48

17.20

2.59

2.06%

SO

1387.83%

8.77%

18.56

19.52

4.14

2.29%

OA

1001.81%8.61%

353.40

83.76

33.25

7.51%

2x2 + OA

1003.33%

32.61%

284.98118.0830.40

15.03%

SO+2x2+OA

1069.59%

31.91%

174.79107.16

29.55

12.56%OB

Table 6.3: Averages for GLOBALLib instances

Cut Family

Initial Gap

End Gap

Closed Gap# CutsItersTime (s)LPTime (%)

103.59%103.56%

0.04%

12.84

13.62

127.150.40%

SO103.33%

0.34%

14.34

15.45

132.070.49%

OA

30.88%

75.55%

676.90137.52459.2831.80%

2x2 + OA

32.84%

74.52%

349.21140.40473.1828.76%

SO+2x2+OA

33.43%

74.03%

227.39136.93475.3826.59%OB

Table 6.4: Averages for BoxQP instances

GLOBALLib

BoxQP
SO+2x2SO+2x2

Closed Gap

OB

SO

OA

2x2 + OA

+OA

OB

SO

OA

2x2 + OA

+OA

>98%

0

1

0

3

3

0

0

37

37

37

90-98 %

0

1

1

2

1

0

0

16

15

15

75-90 %

0

0

0

0

0

0

0

14

11

10

50-75 %

0

0

1

3

4

0

0

13

15

15

25-50 %

0

0

1

4

4

0

0

5

7

8

<25 %

25

23

22

13

13

99

99

14

14

14

Table 6.5: Distribution of Closed Gaps for GLOBALLib and BoxQPCHAPTER 6. CUTTING PLANES FOR POLYNOMIAL OPTIMIZATIONGLOBALLib
Initial

116

BoxQP

2x2

SO+2x2

Initial2x2

SO+2x2End Gap

Gap

OB

SO

OA

+ OA

+OA

Gap

OB

SO

OA

+ OA

+OA

<1 %

2

2

2

2

4

4

0

0

0

35

34

34

1 - 25 %

4

4

4

4

4

4

4

4

4

32

32

31

25-50 %

6

6

6

6

4

4

14

14

17

11

9

9

50-75 %

3

3

3

3

3

3

19

19

16

11

8

8

75-100 %

3

3

3

4

5

5

18

18

18

2

8

9

100-500 %

3

3

3

2

2

2

44

44

44

8

8

8

>500%

4

4

4

4

3

3

0

0

0

0

0

0

Table 6.6: Distribution of End Gaps for GLOBALLib and BoxQPCHAPTER 7. CONCLUSIONS

117

Chapter 7

Conclusions
The development of new polynomial optimization techniques is a cornerstone in the expansion of the current reach computational optimization has. A vast number of articles haveshowed a significant shift towards handling non-linear and non-convex expressions not only
from a theoretical standpoint, but also from the need of solving important real-word application such as the AC-OPF problem. The gap between theory and computations related topolynomial optimization problems seems to be decreasing, and we expect this work can help
in this regard by boosting a synergy between mixed-integer programming and polynomial
optimization.

We showed how useful this approach can be in the AC-OPF and ACTS problems. These
are difficult nonlinear, non-convex optimization problems with an important impact in ourdaily lives. By introducing mixed-integer programming techniques that yield simple relaxations, we are able to address problems quickly and effectively. Even though there is a loss
in expressiveness by only limiting ourselves to LP formulations, there is a large gain in terms
on numerical reliability and speed, and furthermore, these methodologies can also be of use
alongside other frameworks such as SDP or SOCP.
We expect this approach can also be of use in Power Systems problems that involve
renewable energy sources. One can envision models that merge AC-OPF models with elaborate uncertainty models arising from the presence of alternative energy sources, however,the current computational boundaries in optimization can hinder the development of them.

CHAPTER 7. CONCLUSIONS

118

Counting with quick procedures in the classical problems can make the transition to complex
models involving renewable sources more manageable.
From a more abstract perspective, our thorough study on the effect on structured sparsity
on polynomial problems given by tree-width provided strong results regarding how far we cango using LP formulations, and even though some extra machinery is needed to bring this topractice, we provided an important contribution in terms of analyzing much can tree-width be
exploited in a very general setting. There are two main questions that arise from this. First,
is there another graph-theoretical structure we can grasp on in general? More specifically,
if there is a family of instances whose intersection graphs have unbounded tree-width, canwe expect to solve them efficiently? The work of [34] shows that in the pure binary case theanswer is negative, and we strongly believe this is also the case for our approximation resultsin the continuous setting. The second question is related to how tight the sizes of our LP
approximations are. It would be interesting to show some examples where the sizes given in
our results cannot be improved.
The introduction of Network Polynomial Optimization problems, and the study of theirtractability with respect to the structured sparsity of the underlying network, also presents
a significant contribution. This concept encompasses many classical optimization problem,and can provide a unified framework for the study of them; our tractability result being animportant example of this. This result provides a glimpse on how much can be achieved,
even in a fairly generic setting.
Finally, the proposed cutting plane procedures for Polynomial Optimization we studied
provide an interesting approach to quickly compute bounds of challenging problems without
assuming much about their structure. These approaches are flexible enough so we can easily
extend an build upon them. We expect they can be further exploited and blended into moreelaborate solution methods such as branching.
In general, we were able to give an in-depth study on the reach of mixed-integer programming techniques in Polynomial Optimization. And we believe this is only a starting point.
From both a theoretical and computational standpoint we see an important contribution,
providing novel insights into the Polynomial Optimization field, and hopefully giving a boost
into the computational aspects of it.

BIBLIOGRAPHY

119

Bibliography[1]

R.K. Ahuja, T.L. Magnanti, and J.B. Orlin. Network flows: theory, algorithms, andapplications. Prentice Hall, 1993.

[2]

K M Anstreicher. Semidefinite programming versus the reformulation linearization
technique for nonconvex quadratically constrained quadratic programming. Journal of
Global Optimization, 43:471–484, 2009.

[3]

R. Aris, G.L. Nemhauser, and D.J. Wilde. Optimization of multistage cycle and branching systems by serial procedures. Journal of American Institute of Chemical Engineers,
10:913–919, 1964.

[4]

S. Arnborg. Efficient algorithms for combinatorial problems on graphs with bounded
decomposability. A survey. BIT Numerical Mathematics, 25:2–23, 1985.

[5]

S. Arnborg, D. Corneil, and A. Proskurowski. Complexity of finding embeddings in ak-tree. SIAM Journal on Algebraic Discrete Methods, 8(2):277–284, 1987.

[6]

S. Arnborg and A. Proskurowski. Linear time algorithms for NP-hard problems ongraphs embedded in k-trees. Technical Report TRITA-NA-8404, The Royal Instituteof Technology, Stockholm, 1984.

[7]

Stefan Arnborg, Derek G Corneil, and Andrzej Proskurowski. Complexity of finding
embeddings in ak-tree. SIAM Journal on Algebraic Discrete Methods, 8(2):277–284,
1987.

[8]

Y. H. Au and L. Tunçel. A comprehensive analysis of polyhedral lift-and-project methods, December 2013. arXiv:1312.5972.

BIBLIOGRAPHY
[9]

120

Egon Balas. Intersection cuts—a new type of cutting planes for integer programming.
Operations Research, 19(1):19–39, 1971.

[10]

Egon Balas, Sebastián Ceria, and Gérard Cornuéjols. A lift-and-project cutting planealgorithm for mixed 0-1 programs. Math. Program., 58(3):295–324, 1993.

[11]

Egon Balas and François Margot. Generalized intersection cuts and a new cut generating paradigm. Mathematical Programming, 137(1-2):19–35, 2013.

[12]

Xiaowei Bao, Nikolaos V Sahinidis, and Mohit Tawarmalani. Multiterm polyhedralrelaxations for nonconvex, quadratically constrained quadratic programs. Optimization
Methods & Software, 24(4-5):485–504, 2009.

[13]

Imre Bárány and Attila Pór. 0-1 polytopes with many facets. Advances Math., 161:209–
228, 2001.

[14]

Amitabh Basu, Michele Conforti, Gérard Cornuéjols, and Giacomo Zambelli. Maximal lattice-free convex sets in linear subspaces. Mathematics of Operations Research,
35(3):704–720, 2010.

[15]

Amitabh Basu, Gérard Cornuéjols, and Giacomo Zambelli. Convex sets and minimal
sublinear functions. Journal of Convex Analysis, 18(2):427–432, 2011.

[16]

A.R. Bergen and V. Vittal. Power Systems Analysis. Prentice Hall, 1999.

[17]

M.W Bern, E.L Lawler, and A.L Wong. Linear-time computation of optimal subgraphsof decomposable graphs. Journal of Algorithms, 8(2):216 – 235, 1987.

[18]

U. Bertele and F. Brioschi. Nonserial Dynamic Programming. Academic Press, 1972.

[19]

D. Bienstock and M. A. Langston. Chapter 8 algorithmic implications of the graphminor theorem. In C.L. Monma M.O. Ball, T.L. Magnanti and G.L. Nemhauser, editors, Network Models, volume 7 of Handbooks in Operations Research and Management
Science, pages 481 – 502. Elsevier, 1995.

[20]

D. Bienstock and N. Özbay. Tree-width and the Sherali-Adams operator. Discrete
Optimization, 1(1):13–21, June 2004.

BIBLIOGRAPHY
[21]

121

D. Bienstock and M. Zuckerberg. Subset algebra lift operators for 0-1 integer programming. SIAM Journal on Optimization, 15(1):63–95, 2005.

[22]

Daniel Bienstock. Histogram models for robust portfolio optimization. Journal ofComputational Finance, 11:1–64, 2007.

[23]

Daniel Bienstock, Chen Chen, and Gonzalo Muñoz. Outer-product-free sets for polynomial optimization and oracle-based cuts. arXiv preprint arXiv:1610.04604, 2016.

[24]

H. L. Bodlaender. Dynamic programming on graphs with bounded treewidth. In Timo
Lepistö and Arto Salomaa, editors, Automata, Languages and Programming, volume
317 of Lecture Notes in Computer Science, pages 105–118. Springer Berlin Heidelberg,
1988.

[25]

Hans L Bodlaender. A linear-time algorithm for finding tree-decompositions of small
treewidth. SIAM Journal on computing, 25(6):1305–1317, 1996.

[26]

Hans L Bodlaender, Pål Grnås Drange, Markus S Dregi, Fedor V Fomin, Daniel Lokshtanov, and Michal Pilipczuk. A cˆkn 5-approximation algorithm for treewidth. SIAM
Journal on Computing, 45(2):317–378, 2016.

[27]

R.B. Borie, R.G. Parker, and C.A. Tovey. Automatic generation of linear algorithms
on recursive graphs from a predicate calculus. Algorithmica, 7:555–581, 1992.

[28]

S. Bose, S.H. Low, T. Teeraratkul, and B. Hassibi. Equivalent relaxations of optimal
power flow. IEEE Transactions on Automatic Control, 60:729–742, 2014.

[29]

G. Braun, S. Pokutta, and A. Roy. Strong Reductions for Extended Formulations, pages
350–361. Springer International Publishing, Cham, 2016.

[30]

D. J. Brown, M. R. Fellows, and M. A. Langston. Polynomial-time self-reducibility:
Theoretical motivations and practical results. International Journal of Computer Mathematics, 31:19, 1989.

[31]

S. Burer. Optimizing a polyhedral-semidefinite relaxation of completely positive programs. Mathematical Programming Computation, 2(1):1–19, 2010.

BIBLIOGRAPHY
[32]

122

M.B. Cain, R.P. ONeill, and A. Castillo. The IV formulation and linear approximationsof the ac optimal power flow problem. FERC staff technical paper, 2012.

[33]

J. Carpentier. Contribution á l’étude du dispatching économique. Bulletin de la SociétéFranćaise des Électriciens, 3:431–447, 1962.

[34]

V. Chandrasekaran, N. Srebro, and P. Harsha. Complexity of Inference in GraphicalModels. In Proc. 24th Conference on Uncertainty in Artificial Intelligence, 2008.

[35]

Michele Conforti, Gérard Cornuéjols, Aris Daniilidis, Claude Lemaréchal, and Jérôme
Malick. Cut-generating functions and s-free sets. Mathematics of Operations Research,
40(2):276–391, 2014.

[36]

Michele Conforti, Gérard Cornuéjols, and Giacomo Zambelli. Equivalence between
intersection cuts and the corner polyhedron. Operations Research Letters, 38(3):153–
155, 2010.

[37]

Gérard Cornuéjols, Laurence Wolsey, and Sercan Yıldız. Sufficiency of cut-generating
functions. Mathematical Programming, pages 1–9, 2013.

[38]

W. H. Cunningham and Jim Geelen. On integer programming and the branch-widthof the constraint matrix. In Matteo Fischetti and DavidP. Williamson, editors, Integer Programming and Combinatorial Optimization, volume 4513 of Lecture Notes in
Computer Science, pages 158–166. 2007.

[39]

S. Dash, O. Günlük, and A. Lodi. On the mir closure of polyhedra. In Integer Programming and Combinatorial Optimization, pages 337–351. Springer Berlin Heidelberg,
2007.

[40]

Achiya Dax. Low-rank positive approximants of symmetric matrices. Advances in
Linear Algebra & Matrix Theory, 4(3):172–185, 2014.

[41]

R. Dechter and J. Pearl. Tree clustering for constraint networks (research note). Artif.
Intell., 38:353–366, 1989.

[42]

Santanu S Dey and Laurence A Wolsey. Constrained infinite group relaxations of MIPs.
SIAM Journal on Optimization, 20(6):2890–2912, 2010.

BIBLIOGRAPHY
[43]

123

Carl Eckart and Gale Young. The approximation of one matrix by another of lower
rank. Psychometrika, 1(3):211–218, 1936.

[44]

M. R. Fellows and M. A. Langston. Nonconstructive advances in polynomial-timecomplexity. Information Processing Letters, 26(3):157162, 1987.

[45]

E.B. Fisher, R.P. O’Neill, and M.C. Ferris. Optimal transmission switching. IEEE
Transactions on Power Systems, 23:1346–1355, 2008.

[46]

Christodoulos A Floudas, Panos M Pardalos, Claire Adjiman, William R Esposito,
Zeynep H Gümüs, Stephen T Harding, John L Klepeis, Clifford A Meyer, and Carl A
Schweiger. Handbook of test problems in local and global optimization, volume 33.
Springer Science & Business Media, 2013.

[47]

E. Freuder. A sufficient condition for backtrack-bounded search. J. ACM, 32(4):755–
761, October 1985.

[48]

D. R. Fulkerson and O. Gross. Incidence matrices and interval graphs. Pacific Journalof Mathematics, 15:835–855, 1965.

[49]

Dimitris Gatzouras, Apostolos Giannopoulos, and Nilolaos Markoulakis. Lower boundfor the maximal number of facets of a 0/1 polytope. Discrete Comput. Geometry,
34:331–349, 2005.

[50]

Bissan Ghaddar, Juan C Vera, and Miguel F Anjos. A dynamic inequality generation
scheme for polynomial programming. Mathematical Programming, pages 1–37, 2011.

[51]

F. Glover. Improved linear integer programming formulations of nonlinear integer problems. Management Science, 22(4):455–460, 1975.

[52]

J. Gondzio. Splitting dense columns of constraint matrix in interior point methods forlarge scale linear programming. Optimization, 24:285–297, 1992.

[53]

D. Grimm, T. Netzer, and M. Schweighofer. A note on the representation of positive
polynomials with structured sparsity. Archiv der Mathematik, 89:399–403, 2007.

BIBLIOGRAPHY
[54]

124

A. Gupte, S. Ahmed, M.S. Cheon, and S.S. Dey. Solving mixed integer bilinear problems
using MIP formulations. SIAM Journal on Optimization, 23:721–744, 2013.

[55]

K.W. Hedman, S.S. Oren, and R.P. O’Neill. A review of transmission switching andnetwork topology optimization. 2011.

[56]

M. Hewitt, G. L. Nemhauser, and M. W. P. Savelsbergh. Combining exact and heuristic
approaches for the capacitated fixed-charge network flow problem. INFORMS Journal
on Computing, 22(2):314–325, 2010.

[57]

Hassan Hijazi, Carleton Coffrin, and Pascal Van Hentenryck. Convex quadratic relaxations for mixed-integer nonlinear programs in power systems. Mathematical Programming Computation, pages 1–47, 2016.

[58]

Gurobi Optimization Inc. Gurobi optimizer reference manual, 2016.

[59]

R. Jabr. Exploiting sparsity in sdp relaxations of the opf problem. IEEE Transactions
on Power Systems, 27:1138–1139, 2012.

[60]

T. Koch, B. Hiller, M. Pfetsch, and L. Schewe, editors. Evaluating Gas Network Capacities. SIAM, Philadelphia, 2015.

[61]

Burak Kocuk, Santanu S Dey, and X Andy Sun. Strong socp relaxations for the optimal
power flow problem. Operations Research, 64(6):1177–1196, 2016.

[62]

Burak Kocuk, Santanu S Dey, and Xu Andy Sun. Inexactness of sdp relaxation and valid
inequalities for optimal power flow. IEEE Transactions on Power Systems, 31(1):642–
651, 2016.

[63]

M. Kojima, S. Kim, and H. Waki. Sparsity in sums of squares of polynomials. Mathematical Programming, 103:4562, 2005.

[64]

Ulrich H. Kortenkamp, Jürgen Richter-Gebert, A. Sarangarajan, and Günter M.
Ziegler. Extremal properties of 0/1-polytopes. Discrete & Computational Geometry,
17:439–448, 1997.

BIBLIOGRAPHY
[65]

125

J. Lasserre. Convergent SDP relaxations in polynomial optimization with sparsity.
SIAM Journal on Optimization, 17(3):822–843, 2006.

[66]

J. B. Lasserre. An explicit exact SDP relaxation for nonlinear 0-1 programs. In Karen
Aardal and Bert Gerards, editors, Integer Programming and Combinatorial Optimization, volume 2081 of Lecture Notes in Computer Science, pages 293–303. 2001.

[67]

J B Lasserre. Global optimization with polynomials and the problem of moments.
SIAM Journal on Optimization, 11:796–817, 2001.

[68]

J.B. Lasserre, K.-C. Toh, and S. Yang. A bounded degree SOS hierarchy for polynomial
optimization, January 2015. arXiv:1501.006126.

[69]

M. Laurent. A comparison of the Sherali-Adams, Lovász-Schrijver and Lasserre relaxations for 0-1 programming. Mathematics of Operations Research, 28:470–496, 2001.

[70]

M. Laurent. Sum of squares, moment matrices and optimization over polynomials.
IMA, pages 1–147, 2010.

[71]

Monique Laurent. Sums of squares, moment matrices and optimization over polynomials. In Emerging applications of algebraic geometry, pages 157–270. Springer, 2009.

[72]

S. L. Lauritzen. Graphical Models. Oxford University Press, 1996.

[73]

J. Lavaei and S. H. Low. Zero duality gap in optimal power flow problem. IEEE Trans.
Power Systems, 27(1):92–107, 2012.

[74]

Karsten Lehmann, Alban Grastien, and Pascal Van Hentenryck. Ac-feasibility on tree
networks is np-hard. IEEE Transactions on Power Systems, 31(1):798–801, 2016.

[75]

Bi Li, Fatima Zahra Moataz, Nicolas Nisse, and Karol Suchan. Minimum size treedecompositions. Electronic Notes in Discrete Mathematics, 50:21–27, 2015.

[76]

M Locatelli and F Schoen. On convex envelopes for bivariate functions over polytopes.
Mathematical Programming, pages 1–27, 2013.

[77]

L. Lovász and A. Schrijver. Cones of matrices and set-functions and 0-1 optimization.
SIAM J. on Optimization, 1:166–190, 1991.

BIBLIOGRAPHY
[78]

126

L Lovász and A Schrijver. Cones of matrices and set-functions and 0-1 optimization.
SIAM Journal on Optimization, 1:166–190, 1991.

[79]

James Luedtke, Mahdi Namazifar, and Jeff Linderoth. Some results on the strength ofrelaxations of multilinear functions. Mathematical programming, 136(2):325–351, 2012.

[80]

R. Madani, M. Ashraphijuo, and J. Lavaei. OPF Solver. http://www.ee.columbia.
edu/lavaei/Software.html, 2014.

[81]

Ramtin Madani, Morteza Ashraphijuo, and Javad Lavaei. Promises of conic relaxation
for contingency-constrained optimal power flow problem. IEEE Transactions on Power
Systems, 31(2):1297–1307, 2016.

[82]

A. Meeraus. GLOBALLib. http:/http://www.gamsworld.org/global/globallib.
htm.

[83]

Leon Mirsky. Symmetric gauge functions and unitarily invariant norms. The quarterly
journal of mathematics, 11(1):50–59, 1960.

[84]

Ruth Misener and Christodoulos A Floudas. Advances for the pooling problem: Modeling, global optimization, and computational studies. Applied and Computational Mathematics, 8(1):3–22, 2009.

[85]

Ruth Misener and Christodoulos A Floudas. Global optimization of mixed-integer
quadratically-constrained quadratic programs (miqcqp) through piecewise-linear andedge-concave relaxations. Mathematical Programming, 136(1):155–182, 2012.

[86]

D.K. Molzahn and I.A. Hiskens. Sparsity-exploiting moment-based relaxations of the
optimal power flow problem. IEEE Transactions on Power Systems, 30:3168–3180,
2014.

[87]

D.K. Molzahn, J.T. Holzer, B.C. Lesieutre, and C.L. DeMarco. Implementation ofa large-scale optimal power flow solver based on semidefinite programming. IEEE
Transactions on Power Systems, 28(4):3987–3998, 2013.

[88]

G.L. Nemhauser. Introduction to Dynamic Programming. Wiley, 1966.

BIBLIOGRAPHY
[89]

127

J. Nie. Optimality conditions and finite convergence of Lasserre’s hierarchy. Mathematical Programming, 146:97–121, 2014.

[90]

Jiawang Nie. Optimality conditions and finite convergence of lasserres hierarchy. Mathematical programming, 146(1-2):97–121, 2014.

[91]

Narayana Prasad Padhy. Unit commitment-a bibliographical survey. IEEE Transactions on power systems, 19(2):1196–1205, 2004.

[92]

J. Pearl. Reverend bayes on inference engines: a distributed hierarchical approach.
In in Proceedings of the National Conference on Artificial Intelligence, pages 133–136,
1982.

[93]

Dzung T. Phan. Lagrangian duality and branch-and-bound algorithms for optimal
power flow. Operations Research, 60(2):275–285, 2012.

[94]

Anatoliy D Rikun. A convex envelope formula for multilinear functions. Journal of
Global Optimization, 10(4):425–437, 1997.

[95]

N. Robertson and P.D. Seymour. Graph minors. III. Planar tree-width. Journal of
Combinatorial Theory, Series B, 36(1):49 – 64, 1984.

[96]

N. Robertson and P.D. Seymour. Graph minors II: Algorithmic aspects of tree-width.
Journal of Algorithms, 7:309 – 322, 1986.

[97]

A Saxena, P Bonami, and J Lee. Convex relaxations of non-convex mixed integer
quadratically constrained programs: extended formulations. Mathematical programming, 124:383–411, 2010.

[98]

A Saxena, P Bonami, and J Lee. Convex relaxations of non-convex mixed integer
quadratically constrained programs: projected formulations. Mathematical programming, 130:359–413, 2011.

[99]

H. Sherali and W. Adams. A hierarchy of relaxations between the continuous andconvex hull representations for zero-one programming problems. SIAM Journal on
Discrete Mathematics, 3(3):411–430, 1990.

BIBLIOGRAPHY

128[100] H.D. Sherali, E. Dalkiran, and L. Liberti. Reduced RLT representations for nonconvex
polynomial programs. Journal of Global Optimization, 52, 2012.
[101] H.D. Sherali and C.H. Tuncbilek. New reformulation linearization/convexification relaxations for univariate and multivariate polynomial programming problems. Operations Research Letters, 21:1–9, 1997.
[102] N Z Shor. Quadratic optimization problems. Soviet Journal of Circuits and Systems
Sciences, 25:6, 1987.
[103] F Tardella. Existence and sum decomposition of vertex polyhedral convex envelopes.
Optimization Letters, 2:363–375, 2008.
[104] M Tawarmalani, J-P P Richard, and C Xiong. Explicit convex and concave envelopes
through polyhedral subdivisions. Mathematical Programming, pages 1–47, 2013.
[105] Mohit Tawarmalani and Nikolaos V Sahinidis. Convex extensions and envelopes oflower semi-continuous functions. Mathematical Programming, 93:247–263, 2002.
[106] D. Vandenbussche and G. Nemhauser. A branch-and-cut algorithm for nonconvex
quadratic programs with box constraints. Mathematical Programming, 102(3):559–575,
2005.
[107] A. Verma.

Power Grid Security Analysis: an Optimization Approach.

disserta-tion, Columbia University, 2007. http://www.columbia.edu/~dano/theses/theses.
html/.
[108] M.J. Wainwright and M.I. Jordan. Treewidth-Based conditions for exactness of theSherali-Adams and Lasserre relaxations. Technical Report 671, University of California,
September 2004.
[109] M.J. Wainwright and M.I. Jordan. Graphical models, exponential families, and variational inference. Found. Trends Mach. Learn., 1:1–305, 2008.
[110] H. Waki, S. Kim, M. Kojima, and M. Muramatsu. Sums of squares and semidefinite programming relaxations for polynomial optimization problems with structured
sparsity. SIAM Journal on Optimization, 17:218–242, 2006.

BIBLIOGRAPHY

129[111] R. D. Zimmerman, C. E. Murillo-Sánchez, and R. J. Thomas. MATPOWER: Steadystate operations, planning, and analysis tools for power systems research and education.
IEEE Transactions on Power Systems, 26:12–19, 2011.
[112] M. Zuckerberg. A Set Theoretic Approch to Lifting Procedures for 0,1 Integer Programming. Ph.D. dissertation, Columbia University, 2009. http://www.columbia.
edu/~dano/theses/theses.html/.

130

Appendices

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2

131

Appendix A

Extra details for Chapter 2
A.1

Derivation of General Case

Here we specify the derivation of the formulas (2.42) and (2.43) for general case of a linewith a transformer and/or nonzero shunts. We assume the transformer is on the k side ofline {k, m} and we drop the subscripts to simplify notation. In this case we have


1 sh
1 1
y Vk
Ikm = y Vk − ejσ Vm +
τ
τ
2τ 2


1 1
1
= y Vk − (cos σ + j sin σ)Vm +(g sh + jbsh )Vk .
ττ2τ 2
In rectangular coordinates this can be further expanded as11
Ikm = (g + jb) (ek + jfk ) − (cos σ + j sin σ)(em + jfm )ττ1+ 2 (g sh + jbsh )(ek + jfk )
2τ


1
1
= (g + jb) (ek + jfk ) − em cos σ + fm sin σ − j(em sin σ + fm cos σ)
τ
τ
h
i
1
+ 2 g sh ek − bsh fk + j(bsh ek + g sh fk )
2τ


1
1
1
= (g + jb) ek − em cos σ + fm sin σ + j( fk − em sin σ − fm cos σ)
τ
τ
τ
h
i
1
+ 2 g sh ek − bsh fk + j(bsh ek + g sh fk ) .
2τAPPENDIX A. EXTRA DETAILS FOR CHAPTER 2

132

From this expression we obtain:i b 1g h ek− em cos σ + fm sin σ −fk − em sin σ − fm cos σ
Re Ikm =
τ τ
τ τ
i
1 h
+ 2 g sh ek − bsh fk ,2τ
i g 1b h ek− em cos σ + fm sin σ +fk − em sin σ − fm cos σ
Im Ikm =
τ τ
τ τ
i
1 h
+ 2 bsh ek + g sh fk .
2τ

(A.1)

(A.2)

We can easily verify that in the “no-transformer” case, i.e. τ = 1 and σ = 0, (A.1) and(A.2) match the expansion (2.30) for Ikm , as desired. As a result,∗
Pkm = Re Vk Ikm
 h
i b 1g ek= ek− em cos σ + fm sin σ −fk − em sin σ − fm cos σ
τ τ
τ τ
 h

i g 1b ek+ fk− em cos σ + fm sin σ+ +fk − em sin σ − fm cos σ
τ τ
τ τ
g sh
+ 2 (e2k + fk2 ).
2τThis can be rewritten asi
1 h ek1 1
ek
Pkm =− em cos σ (g , b)(fk ) +fk − fm cos σ (−b , g)(efkk )τ ττ τ
g sh 2(e + fk2 )+2τ 2 kgek fm + bek em + bfk fm − gfk em+
sin σ.
τ(A.3)

Note that in (A.3) the third term vanishes when there is no transformer, and the second
term vanishes when there is no shunt conductance. Thus, in the no-transformer case thisexpression evaluates to(ek − em )(g , b)(efkk ) + (fk − fm )(−b , g)(efkk ) +g sh 2(e + fk2 )2 k

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2

133

which is the same as (2.31), as desired. We can further rewrite (A.3) as
Pkm =i
1 h ek− em cos σ + fm sin σ (g , b)(efkk )
τ τ

1 1
+
fk − fm cos σ − em sin σ (−b , g)(efkk )τ τ
g sh 2+(e + fk2 ).
2τ 2 k(A.4)

This is equation (2.42) given above.

Next we compute an expression for Pmk . If there is a transformer present, the expressions
are not symmetric and we first need to compute Imk .


1
y sh
Imk = − jσ yVk + y +
Vm
τe
2
1
= − (cos σ − j sin σ)(g + jb)(ek + jfk )τ+ g + g sh /2 + j(b + bsh /2) (em + jfm )

(A.5)

1
= − (g + jb) [ek cos σ + fk sin σ + j(−ek sin σ + fk cos σ)]τhi+ (g + g sh /2)em − (b + bsh /2)fm + j (b + bsh /2)em + (g + g sh /2)fm .

(A.6)

Therefore,
Re Imk = −bg[ek cos σ + fk sin σ] + [−ek sin σ + fk cos σ]ττ+ (g + g sh /2)em − (b + bsh /2)fm ,Im Imk = −

and

(A.7)

g
b
[−ek sin σ + fk cos σ] − [ek cos σ + fk sin σ]ττ+ (b + bsh /2)em + (g + g sh /2)fm .
Finally, following the power flow equations we obtain∗
Pmk = Re Vm Imkgb= em − [ek cos σ + fk sin σ] + [−ek sin σ + fk cos σ]
τ
τ


g
b
+ fm − [−ek sin σ + fk cos σ] − [ek cos σ + fk sin σ]ττ2+ (g + g sh /2)(e2m + fm
)

(A.8)

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2

134




1
1
em
)
= em − ek cos σ (g , b)(fm ) + fm − fk cos σ (−b , g)(efmmττg sh 22+(e + fm)
2 m−gem fk − bfm fk − bem ek + gfm ek
+
sin σ
τ




1
1
1
1
)
+
f
−
)
f
cos
σ
+
e
sin
σ
(−b , g)(efm= em − ek cos σ − fk sin σ (g , b)(efm
m
k
k
m
m
τ
τ
τ
τ
g sh 22+(e + fm
),
(A.9)
2 m


which is equation (2.43). We now turn to the representation of Pkm and Pmk using polar
coordinates.
1
1 ∗ 1 −jσ ∗∗∗−(A.10)VeV(g sh − jbsh )Vk Vk∗
Skm = Vk Ikm
= Vk
m y +
τ2 k
τ
2τ 2


1 ∗ 1 −jσ ∗1= Vk
Vk − e Vm (g − jb) + 2 (g sh − jbsh )Vk Vk∗
2
τ
τ
2τ
g
b
g sh
g
= |Vk |2 2 − |Vk ||Vm | cos(θkm − σ) − |Vk ||Vm | sin(θkm − σ) + 2 |Vk |2
τ
τ
τ
2τ


b
b
g
bsh
2
2
+ j −|Vk | 2 + |Vk ||Vm | cos(θkm − σ) − |Vk ||Vm | sin(θkm − σ) − 2 |Vk | .
τττ2τSimilarly,gbg sh
Smk = |Vm |2 g − |Vk ||Vm | cos(θmk + σ) − |Vk ||Vm | sin(θmk + σ) +|Vm |2
τ
τ
2


g
bsh
b
2
2
|Vm | .
+ j −|Vm | b + |Vk ||Vm | cos(θmk + σ) − |Vk ||Vm | sin(θmk + σ) −ττ2Hence the active power loss equals
Pkm + Pmk
|Vk |2
g sh
|Vk |
g sh
2
2
=
+
|V
|
g
−
2g
|V
|
cos(θ
−
σ)
+
|V
|
+
|Vm |2 .
m
m
km
k
τ2
τ
2τ 2
2

(A.11)(2)

Relaxing this to an inequality, and using the lifted variables Vk , we obtain|Vk |2|Vk |g sh
2
Pkm + Pmk ≥
+
|V
|
g
−
2g
|V
|
cos(θ
−
σ)
+
m
m
km
τ2
τ
2


|Vk |2|Vk |=+ |Vm |2 g − 2g
|Vm | cos(θkm ) cos(σ)
τ2
τ
!
(2)
|Vk |
g sh Vk+ 2g|Vm | sin(θkm ) sin(σ) ++ Vm(2)
τ
2
τ2


(2)

Vk
+ Vm(2)
τ2!

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2

135

Thus, in order to prove the validity of (2.46) it suffices to show
2 
2

1
1
1
1
em − ek cos σ − fk sin σ + fm − fk cos σ + ek sin σ(A.12)
ττττ|Vk |2
|Vk |
|Vk |
2
=
+ |Vm | − 2|Vm | cos(θkm ) cos(σ) + 2|Vm | sin(θkm ) sin(σ)
2
τ
τ
τ
2

As e2k + fk2 = |Vk |2 , it can be easily seen that the term |Vτk2| + |Vm |2 corresponds to thepure quadratics terms in the expansion of (A.12). Thus, we only need to argue that
|Vk |
|Vk |
|Vm | cos(θkm ) cos(σ) +|Vm | sin(θkm ) sin(σ)
τ
τ
1
1
1
= − em ek cos σ − em fk sin σ + 2 ek fk cos σ sin σ
τ
τ
τ
1
1
1
− fm fk cos σ + fm ek sin σ − ek fk cos σ sin στττ−Which can be easily verified to be true. This concludes the proof of (2.46).

The details of the proof mentioned in the main body are as follows: consider point k1
in Figure 2.1. The power injection into the line, at k1 , is equal to Pkm (i.e. it equals
Vk∗∗
N (N ∗ Ikm ) = Vk Ikm = Pkm ). Moreover by construction the voltage magnitude at k1 equals|Vk |/τ and the phase angle difference from k1 to m equals θkm − σ. We can now recover(A.11) from (2.20), with the last two terms account for shunts, as when deriving (2.35).

A.2

Proof of Theorem 2.4.4

Here we provide a proof of Theorem 2.4.4. For convenience of the reader we restate theTheorem and the relevant definitions:Definition 2.4.2 Let G be an undirected graph. A pseudo-flow is a vector P that assigns toeach edge {k, m} of G two reals, Pkm and Pmk . For any node k of G we write.
δk = set of nodes of G adjacent to kand. X
ok (P ) =
Pkm .
m∈δk

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2

136

We call ok (P ) the net output of k. We say that k is a source if ok (P ) > 0 and is a sink if
ok (P ) < 0. Likewise, an edge {k, m} is termed a sink-edge (or source-edge) if
Pkm + Pmk > 0

or, respectively

Pkm + Pmk < 0.

Theorem 2.4.4 Let G be a graph and P be a pseudo-flow on G. Then
X

X

ok (P ) =(Pkm + Pmk ).

{k,m}k : ok (P )<0

k : ok (P )>0

X(−ok (P )) +Furthermore, P can be decomposed into directed flow paths, each originating at a source orsource-edge and terminating at a sink or sink-edge.
Proof: We obtain a directed graph D and a flow f on D as follows. First, we subdivide
each edge e = {k, m} of G by introducing a new node ve . Moreover• If Pkm ≥ 0 the edge between k and e is directed from k to e and we write fke = Pkm .
• If Pkm < 0 we direct the edge from e to k and we set fek = −Pkm .
It follows that f ≥ 0. Moreover, for any node v of D write
δv+ = set of arcs of D of the form (v, v 0 ),
δv− = set of arcs of D of the form (v 0 , v).
Then for any node k of GX

ok (P ) =

fa −

X

fa

a∈δ − (k)a∈δ + (k)while for any edge e = {k, m} of G,X

Pkm + Pmk =fa −a∈δ + (e)Xfa .

a∈δ − (e)Thus k is a source (sink) of P in G if and only ifXa∈δ + (k)fa −Xa∈δ − (k)fa > 0(< 0 resp.)

APPENDIX A. EXTRA DETAILS FOR CHAPTER 2
and edge e = {k, m} of G is a sink-edge (source-edge) for P if and only ifXa∈δ + (e)fa −Xfa > 0(< 0 resp.)a∈δ − (e)The theorem now follows from standard network flow concepts [1].

137

APPENDIX B. EXTRA DETAILS FOR CHAPTER 3

138

Appendix B

Extra details for Chapter 3
B.1

Proof of part (b) of Theorem 3.2.3

Here we describe a procedure that constructs formulation LPz or LP-GB which requires
O(2ω m) oracle queries and with additional workload Õ(ωn2ω (m + ω)), as per Theorem 3.2.3(b). Here, as per the formulation, we have a tree-decomposition (T, Q) of the intersection
graph of a problem GB, of width ω. The critical element in the procedure is the constructionof the sets Ft of Definition 3.5.1, and we remind the reader that for 1 ≤ i ≤ m constraint i
has support Ki and the set of feasible solutions for constraint i is indicated by S i ⊆ {0, 1}Ki .
Note that |Ki | ≤ ω + 1 for all i. The procedure operates as follows:
1. For each constraint i, enumerate each partition of Ki . Given a partition (A1 , A0 ) if thevector y ∈ {0, 1}Ki defined by yj = k if j ∈ Ak (k = 0, 1) is such that y ∈/ S i (i.e., notfeasible) then we record the triple (i, A1 , A0 ) as a vector of length |Ki | + 1 (with someabuse of notation). This process requires 2|Ki | oracle queries. This sum of all thesequantities is O(2ω m) but the more precise estimate will be needed.
2. Let L be the list of all vectors recorded in Step 1, sorted lexicographically; first by theindex i, then by A1 and then by A0 . After post-processing if necessary, we can assume
that L contains no duplicates. These can be performed in timeO(ω|L| log |L|) = O(2ω m(ω + log m)).

APPENDIX B. EXTRA DETAILS FOR CHAPTER 3

139

3. For each t ∈ V (T ) construct a list of all constraints i such that Ki ⊆ Qt . This can bedone in time O(ωmn).
4. For each t ∈ V (T ) we form the sub-list of L consisting of all vectors (i, A1 , A0 ) (constructed in Step 1) such that Ki ⊆ Qt . Note that for any such i the total number ofsuch vectors is at most 2|Ki | . Given a vector (i, A1 , A0 ) thus enumerated, we form allvectors of the form (A01 , A00 ) such that A01 ∪ A00 = Qt and Ak ⊆ A0k for k = 0, 1. Let
Lt be the list of all vectors obtained this way. Clearly, |Lt | ≤ 2ω | {i : Ki ⊆ Qt } |. We
lexicographically sort Lt .
5. For each t ∈ V (T ) we enumerate all vectors y ∈ {0, 1}Qt . For any such vector y, we havethat y ∈ Ft if and only if y is not found in the list Lt ; and this test can be performed
in time O(ω log |Lt |) after lexicographically sorting the list.
The total amount of work entailed in Step 4, using |Qt | ≤ ω + 1 for each t ∈ V (T ), is
!
!
X
X
O
ω|Lt | log |Lt | = Õ 2ω ω
| {i : Ki ⊆ Qt } | .
tLikewise, Step 5 requires O(ω2ωtt2 ωt log |L |) = Õ(nω 2 ). This completes the proof.

PAPPENDIX C. EXTRA DETAILS FOR CHAPTER 4

140

Appendix C

Extra details for Chapter 4
C.1

Scaled feasiblity in Theorem 4.1.2

In this section we will prove that the -scaled characteristic of Theorem 4.1.2 cannot be
avoided in general.
Suppose that there is an algorithm A such that any PO, whose intersection graph has
tree-width ≤ 2, can be solved in polynomial time to some given feasibility tolerance  < 1.
That is to say, the algorithm guarantees fi (x) ≥ − for any constraint fi (x) ≥ 0 (as opposedto the -scaled feasibility notion). Note that since  is fixed in this case, the formulation in
Theorem 4.1.2 yields an algorithm that runs in polynomial time (see the result on Theorem
3.2.3 for the time it takes to build the LP formulation) but with a weaker approximation
guarantee than the hypothetical algorithm A.
We claim that the existence of algorithm A implies P = NP. Consider the subset-sum
Pproblem: given n ≥ 2 positive integers a1 , . . . , an find I ⊆ {1, . . . , n} such that j∈I aj =Pj ∈I/ aj . Denotingn.
. 1X
aj and M = 4nS,
S=
2
j=1

the subset-sum problem can be cast as the following (pure feasibility) PO:
M Sy1 = M a1 x1 ,M Syi = M ai xi + M Syi−1 ,M yn = M,(C.1a)
2 ≤ i ≤ n,

(C.1b)(C.1c)

APPENDIX C. EXTRA DETAILS FOR CHAPTER 4

141

M xj (1 − xj ) = 0,

1 ≤ j ≤ n,

(C.1d)

0 ≤ yj ≤ 1,

1 ≤ j ≤ n,

(C.1e)

0 ≤ xj ≤ 1,

1 ≤ j ≤ n.

(C.1f)

Given a solution (x, y) to (C.1) it is clear that x ∈ {0, 1}n and thatP

1 Pn
j aj xj = 2
j=1 aj .

Moreover, the intersection graph of (C.1) has tree-width 2.
By assumption, algorithm A will produce a solution (x̂, ŷ) that violates each of the constraints (C.1a)-(C.1d) by at most  and that satisfies (C.1e)-(C.1f). Then adding (C.1a)(C.1c) yieldsnX

aj x̂j − S

≤

j=1

n
.
M

(C.2)

Moreover, by (C.1d) and (C.1f) for each 1 ≤ j ≤ n,22or 1 − M≤ x̂j ≤ 1.
either 0 ≤ x̂j ≤ MThis follows from the fact g(x) = x(1 − x) is strictly increasing in [0, 1/2), strictly decreasingin (1/2, 1], and g(2/M ) > /M . Thus, suppose we round each x̂j to the nearest integer,
obtaining binary values x̃j for 1 ≤ j ≤ n. Using (C.2) we obtain


n
n
X
X
n 4S
n 
2
aj x̃j − S ≤
+
aj 
=
+
M
M
M
M
j=1

j=1

and thereforenX

aj x̃j − S

< 1.

j=1

Since the left hand side of the inequality must be an integer, we conclude thatnX

aj x̃j = S.

j=1

This proves that, unless P = N P , algorithm A does not exist.

C.2

LP size dependency on  of Theorem 4.1.2

Now we will prove that the dependency on −1 of the resulting LP size cannot be improvedin general. Suppose that there is an algorithm A that, for any  < 1, solves PO problems toAPPENDIX C. EXTRA DETAILS FOR CHAPTER 4

142

-scaled tolerance (i.e. the violation of any constraint fi (x) ≥ 0 is at most kf k1 ) but whose
running time is polynomial, i.e. in particular it depends polynomially on log −1 . This is incontrast with the formulation in Theorem 4.1.2 yields an algorithm that runs time polynomial
on n, m and −1 . Consider an unscaled version of the previous formulation of the subset-sum
problem, i.e:
Sy1 = a1 x1 ,Syi = ai xi + Syi−1 ,(C.3a)
2 ≤ i ≤ n,yn = 1,(C.3b)(C.3c)

xj (1 − xj ) = 0,

1 ≤ j ≤ n,

(C.3d)

0 ≤ yj ≤ 1,

1 ≤ j ≤ n,

(C.3e)

0 ≤ xj ≤ 1,

1 ≤ j ≤ n.

(C.3f)

Define  = 1/(3SM ) and use algorithm A to find a solution (x̂, ŷ) that is -scaled feasible.
Since the 1-norm of any polynomial in constraints (C.3) is at most 2S + 1, we get that foreach constraint fi (x, y) ≥ 0
fi (x̂, ŷ) ≥ −kfi k1 ≥ −(2S + 1) ≥1
MThis way we can reuse the same argument as before to obtain a solution to the subset-sum
problem. Since we assume the running time depends on log(−1 ) we get a running time that
depends polynomially on log(nS) yielding the same contradiction as before.

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

143

Appendix D

Supplementary Experiments on
Intersection Cuts
Table D.1 details results for GLOBALLib instances using the 2 × 2+OA configuration; Table
D.2 shows results for BoxQP instances. OPT is the best known primal solution. RLT is thestandard RLT relaxation optimal value. RLT-BT is the RLT optimal value after applying
bound tightening (only applied to GLOBALLib instances). Final LB is the final lower boundobtained by the cutting plane algorithm. OA and 2×2 are the number of outer approximationand 2 × 2 cuts added in total, respectively. All other columns are as described in Table 6.2.4.
Table D.3 and Table D.1 compare the 2 × 2+OA configuration with the reported values ofalgorithm V2 by Saxena, Bonami, and Lee [97, 98]. For certain instances of GLOBALLib we
did not obtain the same initial bound and thus excluded these from comparison. For BoxQP
we start with a weak RLT relaxation (with optimal values in the wRLT column of Table D.1)
to match the setup of Saxena et al. Furthermore, result for V2 were reported only for thesmaller instances of BoxQP. V2 uses a RLT relaxation for QCQP problems and applies two
types of cuts in standard cutting plane algorithm fashion. The first is an outer-approximationof the PSD constraint of Shor’s relaxation; these are convex quadratic cuts unlike the linear
cuts we are adding with OA and hence provide stronger but more expensive approximationsof SDP+RLT. The second is disjunctive cuts for which the separation procedure involves amixed-integer linear program.

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

144

On GLOBALLib instances our algorithm terminates with slightly smaller gap closed on
average, with much smaller running times. In Table D.1 we run the 2 × 2+OA configuration
with both 600 and 3600 second time limits. Due to the strong performance of RLT+SDP
on these BoxQP instances, much of the difference could be attributed to differences in theouter-approximation methods used for the PSD constraint.

RLT

RLT-BT

Final LB

Initial Gap

End Gap

Gap Closed

OA

2x2

Iters

Time

LPTimeEx2 1 1-17.00-18.90-18.90-17.89

10.6%

4.9%

53.2%

8

228

67

0.41

2.4%

Ex2 1 5-268.01-269.45-269.08-268.01

0.4%

0.0%

99.6%

2

12

4

0.14

0.0%

Ex2 1 6-39.00-44.40-40.94-39.04

5.1%

0.1%

97.8%

5

132

39

0.64

1.6%

Ex2 1 7-4150.41-13698.36-5820.01-5085.63

40.2%

22.5%

44.0%

90

127

52

7

1.0%

Ex2 1 815639.00-82460.00

14439.00

15626.96

7.7%

0.1%

99.0%

222

73

59

13.82

7.2%

Ex2 1 9

-0.38

-2.20

-2.20

-1.66

131.9%

92.5%

29.9%

1400

891

474

36.9

69.8%

Ex3 1 17049.252533.20

2766.732788.77

60.7%

60.4%

0.5%

149

141

58

1.11

37.8%

Ex3 1 2-30665.54-30802.76-30709.10-30665.55

0.1%

0.0%

100.0%

4

8

5

0.03

0.0%

Ex3 1 3-310.00-310.00-310.00

-

-

-

-

-

-

-

-

-

Ex3 1 4

-4.00

-6.00

-6.00

-5.41

40.0%

28.2%

29.5%

39

228

69

0.26

7.7%

Ex5 2 2 case1-400.00-599.90-599.90-598.97

49.9%

49.6%

0.5%

417

353

154

4.65

47.7%

Ex5 2 2 case2-600.00-1200.00-1200.00-1200.00

99.8%

99.8%

0.0%

27

28

11

0.3

3.3%

Ex5 2 2 case3-750.00-875.00-874.80-873.60

16.6%

16.5%

1.0%

188

182

74

1.67

37.7%

Ex5 2 4-450.00-2933.33-2933.33-2224.97

550.6%

393.6%

28.5%

402

744

256

3.66

15.6%

Ex5 2 5

-3500.00-9700.00-9700.00-9700.00

177.1%

177.1%

0.0%

38

17

11

5.64

2.8%

Ex5 3 2

1.86

1.00

1.00

1.00

30.1%

30.1%

0.0%

38

17

11

1.97

1.5%

Ex5 3 3

3.23

1.63

1.63

1.68

37.9%

36.7%

3.1%

52

3

11

57.97

2.0%

Ex5 4 27512.23

2598.253010.673022.14

59.9%

59.8%

0.3%

94

81

35

0.57

22.8%

Ex8 4 1

0.62

-4.75

-3.84

-0.84

275.1%

90.0%

67.3%

5357

28

1077

600.6963.8%

Ex9 1 4-37.00-63.00-62.96-62.85

68.3%

68.0%

0.4%

389

421

162

3.84

21.1%

Ex9 2 1

17.00

-16.00

1.00

1.00

88.9%

88.9%

0.0%

21

34

11

0.25

0.0%

Ex9 2 2100.00-50.00

66.67

89.79

33.0%

10.1%

69.4%

297

348

129

2.62

11.8%

Ex9 2 30.00

-30.00-30.00-30.003000.0%3000.0%

0.0%

35

19

11

0.88

1.1%

Ex9 2 40.50

-597.00-296.50-296.50

19800.0%19800.0%

0.0%

25

30

11

0.22

0.0%

Ex9 2 6-1.00

-406.00-201.50-18.3110025.2%

865.5%

91.4%

623

127

150

14.37

16.8%

145

OPT

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

Instance

17.00-16.00

1.00

1.00

88.9%

88.9%

0.0%

21

34

11

0.3

0.0%

Ex9 2 8

1.50

0.50

1.50

-

-

-

-

-

-

-

-

-

Table D.1: Detailed results for 2x2 + OA on GLOBALLib instances

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

Ex9 2 7146

RLT

Initial Gap

End Gap

Gap Closed

OA

2x2

Iters

Time

LPTime %spar020-100-1-706.50-1066.00

50.81%

0.01%

99.97%

121

63

37

5.64

6.91%

spar020-100-2-856.50-1289.00

50.44%

0.10%

99.80%

360

146

102

16.87

14.82%

spar020-100-3-772.00-1168.50

51.29%

0.00%

100.00%

30

5

7

1.23

5.69%

spar030-060-1-706.00-1454.75

105.91%

1.13%

98.94%

3450

15

693

600.3558.32%

spar030-060-2-1377.17-1699.50

23.39%

0.00%

99.99%

93

77

34

13.41

9.55%

spar030-060-3-1293.50-2047.00

58.21%

0.37%

99.36%

2074

470

521

557.8768.14%

spar030-070-1-654.00-1569.00

139.69%2.87%

97.94%

3402

18

684

601.3060.23%

spar030-070-2-1313.00-1940.25

47.74%

0.00%

99.99%

160

73

47

19.61

11.32%

spar030-070-3-1657.40-2302.75

38.91%

0.01%

99.97%

459

203

135

57.93

20.40%

spar030-080-1-952.73-2107.50

121.08%1.31%

98.92%

3274

66

668

601.3860.68%

spar030-080-2-1597.00-2178.25

36.37%

0.00%

100.00%

65

20

17

7.18

6.55%

spar030-080-3-1809.78-2403.50

32.79%

0.00%

99.99%

71

26

20

8.42

6.18%

spar030-090-1-1296.50-2423.50

86.86%

0.01%

99.99%

355

137

102

42.96

20.25%

spar030-090-2-1466.84-2667.00

81.76%

0.01%

99.99%

335

90

89

39.40

20.76%

spar030-090-3-1494.00-2538.25

69.85%

0.00%

99.99%

118

74

40

16.25

11.45%

spar030-100-1-1227.13-2602.00

111.95%0.01%

99.99%

756

60

165

79.32

27.84%

spar030-100-2-1260.50-2729.25

116.43%0.02%

99.99%

1808

222

410

257.5444.54%

spar030-100-3-1511.05-2751.75

82.05%

0.15%

99.82%

793

354

234

121.1233.78%

spar040-030-1-839.50-1088.00

29.57%

0.00%

99.99%

945

77

205

270.4845.74%

spar040-030-2-1429.00-1635.00

14.41%

0.00%

99.99%

693

110

165

188.5937.53%

147

OPT

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

Instance

-1086.00-1303.25

19.99%

0.00%

99.99%

986

114

220

305.4546.79%

spar040-040-1-837.00-1606.25

91.80%

6.16%

93.29%

1669

26

339

600.7360.30%

spar040-040-2-1428.00-1920.75

34.48%

0.00%

99.99%

418

112

111

115.3832.72%

spar040-040-3-1173.50-2039.75

73.75%

2.44%

96.69%

1591

24

323

602.7362.16%

spar040-050-1-1154.50-2146.25

85.83%

2.07%

97.59%

1659

16

335

602.5560.50%

spar040-050-2-1430.98-2357.25

64.68%

0.63%

99.03%

1775

10

357

602.2357.78%

spar040-050-3-1653.63-2616.00

58.16%

0.41%

99.30%

1658

7

333

600.21

61.12%

spar040-060-1-1322.67-2872.00

117.05%4.77%

95.92%

1667

18

337

602.2661.10%

spar040-060-2-2004.23-2917.50

45.54%

0.00%

100.00%

984

48

207

268.0545.09%

spar040-060-3-2454.50-3434.00

39.89%

0.00%

100.00%

165

64

47

42.73

13.95%

spar040-070-1-1605.00-3144.00

95.83%

0.00%

100.00%

1353

79

288

426.2651.88%

spar040-070-2-1867.50-3369.25

80.37%

0.00%

100.00%

853

60

187

232.3541.41%

spar040-070-3-2436.50-3760.25

54.31%

0.00%

99.99%

1368

159

313

448.9250.06%

spar040-080-1-1838.50-3846.50

109.16%0.06%

99.94%

1649

11

332

602.6261.57%

spar040-080-2-1952.50-3833.00

96.26%

0.00%

100.00%

985

85

215

277.9444.40%

spar040-080-3-2545.50-4361.50

71.31%

0.03%

99.96%

1776

179

391

601.6954.47%

spar040-090-1-2135.50-4376.75

104.90%0.07%

99.93%

1737

8

349

602.5659.09%

spar040-090-2-2113.00-4357.75

106.18%0.15%

99.86%

1789

1

358

601.2258.01%

spar040-090-3-2535.00-4516.75

78.14%

0.00%

100.00%

441

117

117

118.2427.11%

spar040-100-1-2476.38-5009.75

102.26%0.00%

100.00%

808

92

183

228.6642.57%

spar040-100-2-2102.50-4902.75

133.12%0.64%

99.52%

1687

3

338

600.0260.07%

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar040-030-3

148

-1866.07-5075.75

171.91%

4.35%

97.47%

1617

13

326

600.1161.83%

spar050-030-1-1324.50-1858.25

40.27%

2.31%

94.27%

925

15

188

603.0752.48%

spar050-030-2-1668.00-2334.00

39.90%

4.45%

88.84%

936

14

190

601.4752.07%

spar050-030-3-1453.61-2107.25

44.94%

7.22%

83.94%

911

14

185

600.6753.98%

spar050-040-1-1411.00-2632.00

86.47%

7.56%

91.26%

906

19

185

602.2051.35%

spar050-040-2-1745.76-2923.25

67.41%

5.23%

92.24%

919

11

186

602.8853.03%

spar050-040-3-2094.50-3273.50

56.26%

1.09%

98.07%

890

5

179

603.9350.85%

spar050-050-1-1198.41-3536.00

194.89%50.26%

74.21%

920

65

197

601.7049.24%

spar050-050-2-1776.00-3500.50

97.05%

8.48%

91.26%

968

22

198

602.6350.32%

spar050-050-3-2106.10-4119.75

95.56%

6.83%

92.86%

888

17

181

603.6054.05%

spar060-020-1-1212.00-1757.25

44.95%

44.95%

0.00%

54

1

11

43.39

2.72%

spar060-020-2-1925.50-2238.25

16.23%

16.23%

0.00%

52

3

11

42.36

3.05%

spar060-020-3-1483.00-2098.75

41.49%

16.44%

60.38%

551

24

115

600.1745.57%

spar070-025-1-2538.91-3832.75

50.94%

17.53%

65.58%

354

36

78

603.1431.70%

spar070-025-2-1888.00-3248.00

72.00%

34.35%

52.29%

355

35

78

600.6330.92%

spar070-025-3-2812.28-4167.25

48.16%

13.91%

71.12%

346

39

77

609.4532.52%

spar070-050-1-3252.50-7210.75

121.66%18.93%

84.44%

350

35

77

609.6230.01%

spar070-050-2-3296.00-6620.00

100.82%14.28%

85.84%

355

25

76

607.5530.33%

spar070-050-3-4306.50-7522.00

74.65%

3.60%

95.17%

364

16

76

605.9629.81%

spar070-075-1-4655.50-11647.75

150.16%12.75%

91.51%

357

13

74

603.9430.52%

spar070-075-2-3865.15-10884.75

181.57%19.27%

89.39%

350

15

73

600.1729.43%

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar040-100-3

149

-4329.40-11262.25

160.10%15.57%

90.27%

343

37

76

601.3127.61%

spar080-025-1-3157.00-4840.75

53.32%

25.75%

51.71%

224

21

49

603.3920.23%

spar080-025-2-2312.34-4378.50

89.32%

54.23%

39.28%

240

15

51

601.91

20.41%

spar080-025-3-3090.88-5130.25

65.96%

31.26%

52.60%

237

18

51

610.5822.29%

spar080-050-1-3448.10-9783.25

183.68%48.85%

73.41%

220

30

50

609.1017.83%

spar080-050-2-4449.20-9270.00

108.33%15.77%

85.45%

229

16

49

607.9520.13%

spar080-050-3-4886.00-10029.75

105.25%15.37%

85.39%

229

26

51

608.6419.83%

spar080-075-1-5896.00-15250.75

158.64%14.65%

90.77%

220

20

48

610.9718.21%

spar080-075-2-5341.00-14246.50

166.71%17.20%

89.68%

222

18

48

605.5117.49%

spar080-075-3-5980.50-14961.50

150.15%16.46%

89.04%

217

28

49

604.1317.64%

spar090-025-1-3372.50-6171.50

82.97%

55.94%

32.58%

147

28

35

611.7513.68%

spar090-025-2-3500.29-6015.00

71.82%

50.93%

29.09%

145

30

35

608.43

12.53%

spar090-025-3-4299.00-6698.25

55.80%

32.76%

41.29%

144

31

35

616.6915.60%

spar090-050-1-5152.00-12584.00

144.23%42.44%

70.57%

133

37

34

605.2313.35%

spar090-050-2-5386.50-11920.50

121.28%33.83%

72.11%

131

44

35

613.62

12.83%

spar090-050-3-6151.00-12514.00

103.43%23.61%

77.17%

137

33

34

619.2715.37%

spar090-075-1-6267.45-19054.25

203.99%52.30%

74.36%

113

47

32

610.769.80%

spar090-075-2-5647.50-18245.50

223.03%56.16%

74.82%

116

44

32

615.87

10.41%

spar090-075-3-6450.00-18929.50

193.45%40.61%

79.01%

130

35

33

608.4512.21%

spar100-025-1-4027.50-7660.75

90.19%

76.91%

14.72%

100

20

24

608.357.84%

spar100-025-2-3892.56-7338.50

88.50%

76.29%

13.79%

108

12

24

610.508.82%

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar070-075-3

150

-4453.50-7942.25

78.32%

64.89%

17.15%

109

11

24

604.639.01%

spar100-050-1-5490.00-15415.75

180.76%98.98%

45.25%

8134

23

610.73

6.89%

spar100-050-2-5866.00-14920.50

154.33%79.56%

48.45%

8926

23

606.938.30%

spar100-050-3-6485.00-15564.25

139.98%71.68%

48.79%

87

33

24

620.598.42%

spar100-075-1-7384.20-23387.50

216.69%78.70%

63.68%

7926

21

622.82

6.96%

spar100-075-2-6755.50-22440.00

232.14%93.17%

59.86%

7124

19

607.02

6.05%

spar100-075-3-7554.00-23243.50

207.67%88.07%

57.59%

7020

18

607.93

5.51%

spar125-025-1-5572.00-12251.00

119.85%

119.76%

0.07%

43

2

9

632.220.98%

spar125-025-2-6156.06-12732.00

106.80%106.80%

0.00%

44

1

9

633.901.08%

spar125-025-3-6815.50-12650.75

85.60%

85.61%

0.00%

45

0

9

640.701.07%

spar125-050-1-9308.38-24993.00

168.48%156.65%

7.03%

31

4

7

612.132.84%

spar125-050-2-8395.00-24810.50

195.52%190.86%

2.38%

22

8

6

632.911.65%

spar125-050-3-8343.91-24424.50

192.70%178.70%

7.27%

31

9

8

648.462.89%

spar125-075-1-12330.00-38202.00

209.81%190.70%

9.11%

14

6

4

609.091.01%

spar125-075-2-10382.47-37466.75

260.84%244.29%

6.34%

15

0

3

616.681.02%

spar125-075-3-9635.50-36202.25

275.69%272.41%

1.19%

7

3

2

607.65

0.46%

Table D.2: Detailed results for 2x2 + OA cuts on BoxQP instances

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar100-025-3

151

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

Instance

V2 Gap

V2 TimeGap ClosedTime

Ex2 1 172.62%

704.4053.21%

0.41

Ex2 1 599.98%

0.17

99.68%

0.13

Ex2 1 699.95%

3397.65

93.87%

0.95

Ex2 1 884.70%

3632.28

73.23%

19.13

Ex2 1 998.79%

1587.9429.87%

36.9

Ex3 1 1

15.94%

3600.27

0.34%

0.55

Ex3 1 299.99%

0.08

99.98%

0.04

Ex3 1 486.31%

21.26

29.49%

0.26

Ex5 2 2 case1

0.00%

0.02

2.05%

0.47

Ex5 2 2 case2

0.00%

0.05

0.00%

0.26

Ex5 2 2 case3

0.36%

0.36

0.00%

0.16

Ex5 2 479.31%

68.93

29.04%

5.69

Ex5 3 27.27%

245.82

0.00%

2.33

Ex5 4 227.57%

3614.38

0.24%

0.59

Ex9 1 4

0.00%

0.60

0.00%

0.34

Ex9 2 160.04%

2372.6454.17%

28.37

Ex9 2 288.29%

3606.36

77.90%

30.84

Ex9 2 687.93%

2619.0290.45%

0.12

Ex9 2 8--83.27%

0.12

Table D.3: Comparison with V2 on GLOBALLib instances152

V2 Gap

V2 TimeGap ClosedTime

Gap Closed

Time

spar020-100-1-706.5-1137

95.40%

3638.278.32%

600.4578.34%

3600.97

spar020-100-2-856.5-1328.5

93.08%

3636.665

78.61%

600.6878.61%

1083.09

spar020-100-3

-772

-1224

97.47%

3632.56

83.68%

600.52

83.76%

3600.58

spar030-060-1-706-1472.5

60.00%

3823.051

58.38%

600.9758.51%

3603.02

spar030-060-2-1377.17-1741

91.16%

3715.979

88.55%

601.2288.76%

3601.78

spar030-060-3-1293.5-2073.5

77.41%

3696.495

77.14%

600.9877.22%

3602.11

spar030-070-1

-654

-1647

57.39%

3786.025

53.10%

600.6453.29%

3610.13

spar030-070-2-1313-1989.5

86.60%

3708.212

81.52%

601.1481.62%

3600.68

spar030-070-3-1657.4-2367.5

88.66%

3744.044

86.29%

601.0186.38%

3600.67

spar030-080-1-952.729-2189

69.67%

3600.777

56.96%

600.4

57.04%

3602.02

spar030-080-2

-1597

-2316

86.25%

3627.132

73.92%

601.1173.97%

3600.21

spar030-080-3-1809.78-2504.5

91.42%

3666.392

85.14%

601.3485.28%

3601.93

spar030-090-1-1296.5-2521

81.15%

3676.815

70.02%

600.8570.09%

3600.91

spar030-090-2-1466.84-2755

82.66%

3646.756

71.56%

601.0371.62%

3600.09

spar030-090-3-1494-2619.5

86.37%

3701.849

75.33%

600.5775.56%

3602.01

spar030-100-1-1227.13-2683.5

81.10%

3692.504

69.68%

601.1769.78%

3602.26

spar030-100-2-1260.5-2870.5

72.87%

3697.329

63.60%

600.0563.66%

3601.47

spar030-100-3-1511.05-2831.5

84.10%

3606.496

75.06%

600.4275.13%

3602.89

spar040-030-1-839.5-1162

31.05%

3719.223

50.80%

600.6757.21%

3605.46

spar040-030-2

-1429

-1695

27.74%

3937.898

17.58%

15.57

17.58%

15.49

153

wRLT

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

OPT-1086-1322

28.00%

3798.683

7.73%

15.88

7.73%

16.97

spar040-040-1

-837

-1641

33.31%

3817.844

39.77%

600.5342.76%

3602.06

spar040-040-2-1428-1967.5

35.19%

3968.111

56.73%

600.0861.35%

3600.24

spar040-040-3-1173.5-2089

26.71%

3972.902

43.30%

600.4247.36%

3603.8spar040-050-1-1154.5-2204

36.72%

3819.72

46.05%

600.5449.94%

3602.77

spar040-050-2-1430.98-2403.5

40.87%

3610.64

50.39%

601.6553.90%

3602.41

spar040-050-3-1653.63-2715

33.95%

3639.97746.20%

600.3548.93%

3604.85

spar040-060-1-1322.67-2934

47.75%

3760.964

51.73%

601.6454.33%

3604.92

spar040-060-2-2004.23-3011

55.79%

3707.992

68.04%

601.1

71.04%

3600.39

spar040-060-3-2454.5-3532

72.63%

3764.079

74.36%

601.5776.45%

3600.44

spar040-070-1-1605-3194.5

64.03%

3642.681

67.37%

601.4769.39%

3602

spar040-070-2-1867.5-3446.5

57.91%

3756.377

62.79%

601.5164.77%

3605.11

spar040-070-3-2436.5-3833.5

62.94%

3693.666

68.40%

600.4570.59%

3604.16

spar040-080-1-1838.5-3969

58.37%

3808.258

59.04%

601.5961.09%

3600.31

spar040-080-2-1952.5-3902.5

66.96%

4062.433

63.85%

600.6965.25%

3602.2

spar040-080-3-2545.5-4440

72.31%

4057.149

73.35%

600.4274.94%

3601.66

spar040-090-1-2135.5-4490

66.64%

3781.044

66.03%

600.9267.48%

3603.5spar040-090-2

-2113

-4474

66.46%

3931.349

65.74%

600.3367.12%

3603.76

spar040-090-3

-2535

-4641

73.49%

4003.706

73.80%

601.3275.18%

3600.81

spar040-100-1-2476.38-5118

76.24%

3853.573

74.81%

601.4576.31%

3605.11

spar040-100-2-2102.5-5043

63.89%

3658.261

64.80%

601.7

66.14%

3603.14APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar040-030-3

154

-1866.07-5196.5

59.92%

3842.685

61.46%

602.01

Table D.4: Comparison with V2 on BoxQP instances

63.58%

3600.51

APPENDIX D. SUPPLEMENTARY EXPERIMENTS ON INTERSECTION CUTS

spar040-100-3

155