Some Recent Advances in Mixed Integer Nonlinear Optimization IBM Analytics

IBM Analytics
Pierre Bonami (and many co-authors)
IBM ILOG CPLEX
CORS INFORMS 2015 - Montreal - June 17 2015
Some Recent Advances in Mixed Integer Nonlinear
Optimization
1
©2015 IBM Corporation
IBM Analytics
The mother of all deterministic optimization problems
[Lee, 2008]
2
min f (x)
s.t. gi (x) ≤ 0
x ∈X
xj ∈ Z
lj ≤ xj ≤ uj
1
i = 1, . . . , m
−2
(MINLP)
2
−1
j = 1, . . . , p
j = 1, . . . , p
−2
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
n
X ⊆ R polyhedral.
f and gi : X → R, i = 1, . . . , m,
continuous, differentiable.
2
1
−1
©2015 IBM Corporation
IBM Analytics
"Well solved" subproblems
Nonlinear Programming (NLP)
p = 0 : local optima.
+ f and gi convex ⇒ global optima.
Mixed-Integer linear programming (MILP)
f linear, m = 0, p > 0
3
©2015 IBM Corporation
IBM Analytics
MINLP
2
1
min
s.t.
f (x)
gi (x) ≤ 0
x ∈X
xj ∈ Z
lj ≤ xj ≤ uj
−2
i = 1, . . . , m
2
−1
(PNLM)
−2
j = 1, . . . , p
j = 1, . . . , p
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
To be solvable in general, lj , uj finite.
4
1
−1
©2015 IBM Corporation
IBM Analytics
Two main classes of MINLP
Mixed Integer Convex Program
Assume that the continuous relaxation is a convex optimization problem.
f is a convex function.
gi are convex functions.
Mixed Integer Nonlinear Program
Don’t assume any convexity on f or gi .
Continuous relaxation is NP-hard to solve in general.
Remark: if lj and uj are finite integers, variable xj can be seen as a
continuous variable satisfying:
(xj − lj )(xj − lj − 1)....(xj − uj ) = 0
5
©2015 IBM Corporation
IBM Analytics
A special class of convex MINLP: MISOCP
Second order cone constraints
Ld defines the (d + 1)-dimensional second order cone:
Ld := {(x, x0 ) ∈ Rd+1 :
d
X
xi2 ≤ x02 , x0 ≥ 0}
i=1
Through basic transformation can be used to model any convex quadratic
constraint. (note Ld describes a convex region)
MISOCP
MINLP’s where all nonlinear constraints are SOC.
Continuous relaxation solved efficiently by interior points.
convex MINLP algorithms work with some added technicality due to
non-differentiability [Drewes, 2009, Drewes and Ulbrich, 2012].
Supported by most commercial MIP solvers (CPLEX in particular).
6
©2015 IBM Corporation
IBM Analytics
Agenda
The convex case
Main algorithmic approaches
Glimpse of computations
Exploiting separability, disaggregations
Lift-and-project cuts
Application in CPLEX 12.6.2
A step into nonconvexity.
MIQP in CPLEX 12.6.
Basic setup of a spatial branch-and-bound.
BQP cuts.
Computational results CPLEX 12.6.2.
7
©2015 IBM Corporation
IBM Analytics
The mixed integer convex program
min c T x
s.t. gi (x) ≤ 0 i = 1, . . . , m
x ∈X
xj ∈ Z
j = 1, . . . , p
(MICP)
gi : X → R, i = 1, . . . , m, convex, differentiable.
Assume linear objective. If necessary, add var α ∈ R and min α with
f (x) ≤ α a constraint.
8
©2015 IBM Corporation
IBM Analytics
Main Algorithms for solving (MICP)
z
x
y
Fundamental property is convexity of the continuous relaxation, which can be
efficiently solved.
1 NLP Branch-and-bound [Gupta and Ravindran, 1985].
9
2
Outer Approximation Algorithm [Duran and Grossmann, 1986]. Builds an
MILP equivalent of (MICP)
3
LP/NLP branch-and-cut [Quesada and Grossmann, 1992].
©2015 IBM Corporation
IBM Analytics
NLP based branch-and-bound
Straightforward generalization of main MILP
algorithm:
solve an NLP at each node of the tree.
10
©2015 IBM Corporation
IBM Analytics
NLP based branch-and-bound
Straightforward generalization of main MILP
algorithm:
solve an NLP at each node of the tree.
Branch on variables with fractional value.
integer
feasible
infeasible
10
fathomed
by
bound
©2015
IBM Corporation
IBM Analytics
NLP based branch-and-bound
Straightforward generalization of main MILP
algorithm:
solve an NLP at each node of the tree.
Branch on variables with fractional value.
Prune by infeasibility, bounds and integer
feasibility.
integer
feasible
infeasible
10
fathomed
by
bound
©2015
IBM Corporation
IBM Analytics
NLP based branch-and-bound
Straightforward generalization of main MILP
algorithm:
solve an NLP at each node of the tree.
Branch on variables with fractional value.
Prune by infeasibility, bounds and integer
feasibility.
Main issues
integer
feasible
infeasible
Warm-starting of NLP solves.
Difficulty of reusing MILP technologies.
10
fathomed
by
bound
©2015
IBM Corporation
IBM Analytics
Outer Approximation [Duran and Grossmann, 1986]
min c T x
s.t.
gi (x) ≤ 0
i = 1, . . . , m,
xj ∈ Z,
j = 1, . . . , p.
Idea: Take first-order approximations of constraints at different points and
build an equivalent MILP.
11
©2015 IBM Corporation
IBM Analytics
Outer Approximation [Duran and Grossmann, 1986]
min c T x
s.t.
gi (x) ≤ 0
i = 1, . . . , m,
xj ∈ Z,
j = 1, . . . , p.
Idea: Take first-order approximations of constraints at different points and
build an equivalent MILP.
min
s.t.
cT x
gi (x k ) + ∇gi (x k )T (x − x k ) ≤ 0
xj ∈ Z,
11
i = 1, . . . , m, k = 1, . . . , K
j = 1, . . . , p.
©2015 IBM Corporation
IBM Analytics
Subproblems
Given x̂ ∈ Rp :
fixed NLP (NLP(x̂))
min c T x
s.t.
fixed feasibility subproblem
min
m
X
wi max{0, gi (x)}
i=1
gi (x) ≤ 0,
x ∈ X;
i = 1, . . . , m
(NLP(x̂))
s.t.
x ∈ X,
xj = x̂j ,
j = 1, . . . , p.
xj = x̂j , j = 1, . . . , p
(NLPF(x̂))
If x̂ ∈ Zp , and feasible: gives an upper
bound.
12
©2015 IBM Corporation
IBM Analytics
Subproblems
Given x̂ ∈ Rp :
fixed NLP (NLP(x̂))
min c T x
s.t.
fixed feasibility subproblem
min
m
X
wi max{0, gi (x)}
i=1
gi (x) ≤ 0,
x ∈ X;
i = 1, . . . , m
(NLP(x̂))
s.t.
x ∈ X,
xj = x̂j ,
j = 1, . . . , p.
xj = x̂j , j = 1, . . . , p
(NLPF(x̂))
If x̂ ∈ Zp , and feasible: gives an upper
bound.
Remark If (NLP(x̂)) is infeasible, NLP software will typically return a solution
to (NLPF(x̂)). By abuse, always say solution to (NLP(x̂))
12
©2015 IBM Corporation
IBM Analytics
Equivalent MILP formulation of convex MINLP
For each x̂ k ∈ K = Proj1,...,p (X ) ∩ Zp , let x k be an optimal solution to
(NLP(x̂)).
Theorem ([Duran and Grossmann, 1986])
If X 6= ∅, f and g are convex, continuously differentiable, and a constraint
qualification holds for each x k then
min
cT x
gi (x k ) + ∇gi (x k )T (x − x k ) ≤ 0 i = 1, . . . , m, x̂ k ∈ K ,
x ∈ X , xj ∈ Z, j = 1, . . . , p.
has the same optimal value as (MICP).
13
©2015 IBM Corporation
IBM Analytics
OA decomposition
Generate MILP equivalent by constraint generation.
Initialize with one set of linearizations.
min
cT x
s.t.
gi (x 0 ) + ∇gi (x 0 )T (x − x 0 ) ≤ 0,
i = 1, . . . , m,
,
(OA(K))
x ∈ X , xj ∈ Z, j = 1, . . . , p.
Where x 0 is the solution to the continuous relaxation:
min{c T x : x ∈ X , gi (x) ≤ 0, i = 1, . . . , m}
14
©2015 IBM Corporation
IBM Analytics
OA decomposition
Generate MILP equivalent by constraint generation.
Initialize with one set of linearizations.
Enrich iteratively the set of linearizations K.
min
cT x
s.t.
gi (x k ) + ∇gi (x k )T (x − x k ) ≤ 0,
i = 1, . . . , m,
,
x̂ k ∈ K
(OA(K))
x ∈ X , xj ∈ Z, j = 1, . . . , p.
Where x̂ k is a solution to (OA(K)) and, for k = 1, . . . , |K|, x k is the solution
to (NLP(x̂)).
14
©2015 IBM Corporation
IBM Analytics
OA decomposition
Generate MILP equivalent by constraint generation.
Initialize with one set of linearizations.
Enrich iteratively the set of linearizations K.
Convergence
At each iteration:
(OA(K)) gives a lower bound,
If feasible, (NLP(x̂)) gives an upper bound.
The OA Theorem guarantees that the two bounds converge in finite # of
iterations.
14
©2015 IBM Corporation
IBM Analytics
LP/NLP Branch-and-bound
OA can be embedded in a single tree search.
Start solving the same initial MILP by
branch-and-bound.
At each integer feasible node:
integer
feasible
15
©2015 IBM Corporation
IBM Analytics
LP/NLP Branch-and-bound
OA can be embedded in a single tree search.
Start solving the same initial MILP by
branch-and-bound.
At each integer feasible node:
1
2
3
15
solve (NLP(x̂)), and enrich the set of
linearizations.
Resolve the LP relaxation of the node with the
new cuts.
Repeat as long as node is integer feasible.
integer
feasible
©2015 IBM Corporation
IBM Analytics
LP/NLP Branch-and-bound
OA can be embedded in a single tree search.
Start solving the same initial MILP by
branch-and-bound.
At each integer feasible node:
1
2
3
solve (NLP(x̂)), and enrich the set of
linearizations.
Resolve the LP relaxation of the node with the
new cuts.
Repeat as long as node is integer feasible.
Never prune by integer feasibility.
15
integer
feasible
©2015 IBM Corporation
IBM Analytics
Solvers for Mixed Integer Convex Programs
Solver
Dicopt
MINLP_BB
SBB
α-ECP
Bonmin
FilMINT
KNITRO
SCIP
16
Reference
[Leyffer, 1998]
[Bussieck and Drud, 2001]
[Westerlund and Lundqvist, 2005]
[Bonami et al., 2008]
[Abhishek et al., 2010]
[Byrd et al., 2006]
[Vigerske, 2012]
Algorithm(s)
OA
NLP BB
NLP BB
ECP (variant of OA)
NLP BB, OA, LP/NLP
LP/NLP
NLP BB, LP/NLP
LP/NLP
©2015 IBM Corporation
IBM Analytics
Comparison of solvers in GAMS [Vigerske, 2013]
100
90
% of model solved
80
70
60
50
40
alpha-ECP
Bonmin BB
Bonmin Hyb
Bonmin OA/CPLEX
Bonmin OA/CBC
SBB
SCIP
30
20
10
0
1
17
10
100
time factor
1000
10000
©2015 IBM Corporation
IBM Analytics
Comparison of solvers in GAMS [Vigerske, 2013]
100
90
% of model solved
80
70
60
5 Best Algorithms are OA based
50
40
alpha-ECP
Bonmin BB
Bonmin Hyb
Bonmin OA/CPLEX
Bonmin OA/CBC
SBB
SCIP
30
20
10
0
1
17
10
100
time factor
1000
10000
©2015 IBM Corporation
IBM Analytics
Notes on results with Bonmin
Bonmin’s OA using CPLEX seems the best algorithm overall.
It is also the simplest: a loop calling CPLEX (MILP) and Ipopt (NLP)
alternatively as black boxes.
Improves with CPLEX.
Bonmin’s Hyb is in the pack of relatively good solvers
own variant of LP/NLP BB.
Reuse CBC infrastructure, LP solver, Cuts, MIP presolve.
Improves at a slower pace.
Bonmin’s BB clearly behind.
pure NLP based branch-and-bound. Doesn’t reuse much from Cbc.
Everything specifically tailored.
Better implementation exists that should be on par with Hyb.
Bonmin’s OA using CBC seems the worst algorithm overall.
18
©2015 IBM Corporation
IBM Analytics
Notes on results with Bonmin
Bonmin’s OA using CPLEX seems the best algorithm overall.
It is also the simplest: a loop calling CPLEX (MILP) and Ipopt (NLP)
alternatively as black boxes.
Improves with CPLEX.
Bonmin’s Hyb is in the pack of relatively good solvers
own variant of LP/NLP BB.
Reuse CBC infrastructure, LP solver, Cuts, MIP presolve.
Improves at a slower pace.
Bonmin’s BB clearly behind.
pure NLP based branch-and-bound. Doesn’t reuse much from Cbc.
Everything specifically tailored.
Better implementation exists that should be on par with Hyb.
Bonmin’s OA using CBC seems the worst algorithm overall.
18
©2015 IBM Corporation
IBM Analytics
The MIQCP/MISOCP solver in CPLEX
Implements the two main algorithms:
A branch-and-bound based on the continuous SOCP solver (barrier).
An outer approximation branch-and-cut algorithm.
Choice is controled by the parameter CPXPARAM_MIP_Strategy_MIQCPStrat.
Default is trying to decide which of the two algorithms to run in a “clever” way.
History of MIQCP with CPLEX
class
Convex QCP
convex MIQCP
–
19
p
0
>0
–
algorithm
barrier
barrier based B&B
Outer approximation B&C
V. (Year)
9.0 (2003)
9.0 (2003)
11.0 (2007)
©2015 IBM Corporation
IBM Analytics
A comparison of OA and SOCP-BB in CPLEX 12.6.1
1
Default strategy picked
100
OA 113 times
SOCP-BB 46 times
% model solved
80
56 models identical with
both
To be perfect should have
picked
14 more models with OA
60
40
20
1261_miqcpstrat_0
1261_miqcpstrat_1
1261_miqcpstrat_2
0
1
1 225
20
10
100
1000
time factor
10000
100000
36 more models with
SOCP-BB
models solved by at least one method and failed by none.
©2015 IBM Corporation
IBM Analytics
How bad can outer approximation be?
Consider the following convex MINLP:
min
s.t.
cT x
Pn
i=1
x ∈Z
xi −
n
1 2
2
≤
n−1
4
(1)
z
(1) is infeasible:
The ball is too small to contain
integer points.
It is large enough to touch every
edge of the hypercube.
x
21
y
©2015 IBM Corporation
IBM Analytics
Solving (1) with OA cuts
z
No OA constraint can cut 2
vertices of the hypercube.
If an inequality cuts two
vertices, it cuts the segment
joining them. This can not be:
the ball has non-empty
intersection with any such
segment.
22
x
y
©2015 IBM Corporation
IBM Analytics
Solving (1) with OA cuts
z
No OA constraint can cut 2
vertices of the hypercube.
If an inequality cuts two
vertices, it cuts the segment
joining them. This can not be:
the ball has non-empty
intersection with any such
segment.
OA decomposition takes at least
2n iterations (each iteration
requires solving a MILP).
22
x
y
©2015 IBM Corporation
IBM Analytics
Solving (1) with OA cuts
z
No OA constraint can cut 2
vertices of the hypercube.
If an inequality cuts two
vertices, it cuts the segment
joining them. This can not be:
the ball has non-empty
intersection with any such
segment.
OA decomposition takes at least
2n iterations (each iteration
requires solving a MILP).
An OA Based branch-and-cut
would take at least 2n nodes.
22
x
y
©2015 IBM Corporation
IBM Analytics
Solving (1) with OA cuts
z
No OA constraint can cut 2
vertices of the hypercube.
If an inequality cuts two
vertices, it cuts the segment
joining them. This can not be:
the ball has non-empty
intersection with any such
segment.
OA decomposition takes at least
2n iterations (each iteration
requires solving a MILP).
An OA Based branch-and-cut
would take at least 2n nodes.
22
x
y
Note: A basic NLP
branch-and-bound also
enumerates at least 2n integer
sols.
©2015 IBM Corporation
IBM Analytics
Experimental illustration
n
2n
CPLEX 12.4
nodes
SCIP 2.1
nodes
B-OA
OA it.
B-Hyb
nodes
IBM Analytics
Experimental illustration
n
10
2n
1,024
CPLEX 12.4
nodes
2,047
SCIP 2.1
nodes
720
B-OA
OA it.
1,105
B-Hyb
nodes
11,156
IBM Analytics
Experimental illustration
n
10
15
2n
1,024
32,768
CPLEX 12.4
nodes
2,047
65,535
SCIP 2.1
nodes
720
31,993
B-OA
OA it.
1,105
...
B-Hyb
nodes
11,156
947,014
IBM Analytics
Experimental illustration
n
10
15
20
23
2n
1,024
32,768
1,048,576
CPLEX 12.4
nodes
2,047
65,535
2,097,151
SCIP 2.1
nodes
720
31,993
1,216,354
B-OA
OA it.
1,105
...
...
B-Hyb
nodes
11,156
947,014
...
©2015 IBM Corporation
IBM Analytics
Experimental illustration
n
10
15
20
2n
1,024
32,768
1,048,576
CPLEX 12.4
nodes
2,047
65,535
2,097,151
SCIP 2.1
nodes
720
31,993
1,216,354
B-OA
OA it.
1,105
...
...
B-Hyb
nodes
11,156
947,014
...
Remark
Problem is simple for CPLEX/SCIP if variables are 0 − 1: replace xi2 by
xi , the contradiction n4 ≤ n−1
4 follows.
23
©2015 IBM Corporation
IBM Analytics
Experimental illustration
n
10
15
20
2n
1,024
32,768
1,048,576
CPLEX 12.4
nodes
2,047
65,535
2,097,151
SCIP 2.1
nodes
720
31,993
1,216,354
B-OA
OA it.
1,105
...
...
B-Hyb
nodes
11,156
947,014
...
Remark
Problem is simple for CPLEX/SCIP if variables are 0 − 1: replace xi2 by
xi , the contradiction n4 ≤ n−1
4 follows.
SCIP ≥ 3.0 and CPLEX ≥ 12.6.1 solve it in a blink.
23
©2015 IBM Corporation
IBM Analytics
Separable mixed integer convex programs
min c T x
s.t. gi (x) ≤ 0
x ∈X
xj ∈ Z
l ≤x ≤u
i = 1, . . . , m
(sMINLP)
j = 1, . . . , p
For i = 1, . . . , m, gi : X → R are convex separable:
gi (x) =
n
X
gij (xj )
j=1
with gij : [lj , uj ] → R convex.
24
©2015 IBM Corporation
IBM Analytics
Disaggregated formulation
Introduce one variable yij for each elementary function:
min c T x
n
P
s.t.
yij ≤ 0
i = 1, . . . , m,
x ∈ X,
xi ∈ Z
l ≤ x ≤ u.
i = 1, . . . , p,
j=1
gij (xj ) ≤ yij
25
i = 1, . . . , m,
j = 1, . . . , n,
(sMINLP∗ )
©2015 IBM Corporation
IBM Analytics
Application to (1) [Hijazi et al., 14]
z
Extended formulation of (1)
min c T x
n
X
yi ≤ (n − 1)/4
s.t.
i=1
2
(xi − 0.5) ≤ yi
x ∈ Zn .
i = 1, . . . , n
Its outer approximation
x
y
(2)
min c T x
n
X
yi ≤ (n − 1)/4
s.t.
i=1
2
2 x ki − 0.5 (xi − x ki ) + x ki − 0.5 ≤ yi
26
x ∈ Zn
i = 1, . . . , n
k = 1, . . . , K
©2015 IBM Corporation
IBM Analytics
Application to (1) [Hijazi et al., 14]
Extended formulation of (1)
min c T x
n
X
yi ≤ (n − 1)/4
s.t.
i=1
2
(xi − 0.5) ≤ yi
x ∈ Zn .
i = 1, . . . , n
(2)
Its outer approximation
min c T x
n
X
yi ≤ (n − 1)/4
s.t.
i=1
26
2
2 x ki − 0.5 (xi − x ki ) + x ki − 0.5 ≤ yi
x ∈ Zn
i = 1, . . . , n
k = 1, . . . , K
©2015 IBM Corporation
IBM Analytics
Application to (1) [Hijazi et al., 14]
Extended formulation of (1)
2 points suffice to make it infeasible x 1 = 0 and x 2 = 1:
T
min c x
n
− xi + 0.25 ≤ yi
i = 1, . . . , n
X
yi ≤ (n − 1)/4
s.t.
xi − 0.75+ ≤ yi
i = 1, . . . , n
i=1
2
(xi − 0.5) ≤ yi
x ∈ Zn .
i = 1, . . . , n
(2)
Its outer approximation
min c T x
n
X
yi ≤ (n − 1)/4
s.t.
i=1
26
2
2 x ki − 0.5 (xi − x ki ) + x ki − 0.5 ≤ yi
x ∈ Zn
i = 1, . . . , n
k = 1, . . . , K
©2015 IBM Corporation
IBM Analytics
Experimental Illustration
In the standard benchmark for MICP, out of > 100 instances, 8 are not
directly separable.
Constructing separated formulations on a subset of 47 instances gives a
3x speed up: [Hijazi et al., 14].
90
% of model solved
80
70
60
50
40
30
20
27
10
Sepa
Original
1
10
©2015
100 IBM Corporation
IBM Analytics
Cutting planes for MICP
Cuts are an essential component of MILP solvers.
Of course one can always apply MILP cuts to a linear OA of MICP.
How can we generate cut that also exploit nonlinear constraints?
Can we generate better cuts by looking directly at nonlinear functions?
A partial answer: as long as the cut generated is linear it could also have
been obtained from a linear outer approximation.
In the past three years, tremendous activity towards conic cuts for conic
programming but no general method yet, and no striking computational
results.
28
©2015 IBM Corporation
IBM Analytics
Cutting planes for MICP
Cuts are an essential component of MILP solvers.
Of course one can always apply MILP cuts to a linear OA of MICP.
How can we generate cut that also exploit nonlinear constraints?
Can we generate better cuts by looking directly at nonlinear functions?
A partial answer: as long as the cut generated is linear it could also have
been obtained from a linear outer approximation.
In the past three years, tremendous activity towards conic cuts for conic
programming but no general method yet, and no striking computational
results.
28
©2015 IBM Corporation
IBM Analytics
Split Relaxation
x2
Consider C and M := C ∩ (Zp × Rn−p ).
Let π ∈ Zp × {0}n−p , π0 ∈ Z and
(π,π0 )
C
:= conv C∩ x : π T x ≤ π0 ∪
.
x : π x ≥ π0 + 1
M
C
T
(clearly M ⊆ C(π,π0 ) ⊆ C).
x1
x1 = 0
29
x1 = 1
©2015 IBM Corporation
IBM Analytics
Split Relaxation
Consider C and M := C ∩ (Zp × Rn−p ).
Let π ∈ Zp × {0}n−p , π0 ∈ Z and
(π,π0 )
C
:= conv C∩ x : π T x ≤ π0 ∪
.
x : π x ≥ π0 + 1
T
(clearly M ⊆ C(π,π0 ) ⊆ C).
π T x ≤ π0
29
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Split Relaxation
Consider C and M := C ∩ (Zp × Rn−p ).
Let π ∈ Zp × {0}n−p , π0 ∈ Z and
(π,π0 )
C
:= conv C∩ x : π T x ≤ π0 ∪
.
x : π x ≥ π0 + 1
C (π,π0 )
T
(clearly M ⊆ C(π,π0 ) ⊆ C).
π T x ≤ π0
29
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Split Relaxation
Consider C and M := C ∩ (Zp × Rn−p ).
Let π ∈ Zp × {0}n−p , π0 ∈ Z and
(π,π0 )
C
:= conv C∩ x : π T x ≤ π0 ∪
.
x : π x ≥ π0 + 1
T
αT x = β
x̂
C (π,π0 )
(clearly M ⊆ C(π,π0 ) ⊆ C).
In the remainder, x̂ is the point to separate,
π = ek , x̂k ∈]0, 1[ (k ≤ p), and π0 = 0
29
©2015 IBM Corporation
IBM Analytics
MILP case
Consider a polyhedron P := {x : Ax = b, x ≥ 0}
Cut Generation LP [Balas, 1979, Balas et al., 1993]
x̂ ∈ P is separated from P (ek ,0) using the LP:
min αT x̂ − β
s.t. :
α = u T A − u 0 ek ,
α = v T A + v0 ek ,
(CGLP)
β = u T b,
β = v T b + v0 ,
n
m
α ∈ R , β ∈ R, u, v ∈ R , u0 , v0 ∈ R+
If x̂ 6∈ P (ek ,0) , αT x ≥ β cuts x̂; otherwise produces certificate that x̂ ∈ P (ek ,0)
with x 0 ∈ P ∩ {xk = 0}, x 1 ∈ P ∩ {xk = 1} such that x̂ = x̂k x 1 + (1 − x̂k )x 0 .
30
©2015 IBM Corporation
IBM Analytics
MILP case (primal view)
Membership LP [Bonami, 2012]
x̂ ∈ P wih 0 < x̂k < 1 also in P (ek ,0) if ∃ x 0 ∈ P ∩ {xk = 0} and
x 1 ∈ P ∩ {xk = 1} with x̂ = x̂k x 1 + (1 − x̂k )x 0 , or if
max yk
s.t.
Ay = bx̂k
(MLP)
0 ≤ y ≤ x̂,
y ∈ Rn .
has a solution with yk = x̂k otherwise can deduce a cut from dual optimal
solution.
x̂−y
is x 0 ).
(Hint x̂yk is x 1 , 1−x̂
k
31
©2015 IBM Corporation
IBM Analytics
Generalization to MICP: Two competing approaches
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
Using only LP [Kılınc et al., 2011].
Using NLP [Bonami, 2011]
1 Start with any linear OA of C
1 Solve a single NMLP that tells if
2 Solve CGLP. If a cut is found
x̂ is in the split relaxation.
STOP otherwise
2 If not, deduce from solution two
3 Deduce from dual of CGLP two
points such that
points such that
x̂ = λx 1 + (1 − λ)x 0 and closest
x̂ = λx 1 + (1 − λ)x 0 and
to the disjunction.
satisfying the disjunction.
3 Build OA around these two points.
4 If point(s) not in C generate new
4 Solve MLP and get the cut.
OA and goto 2, otherwise there is
no cut, STOP.
32
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x̂
π T x ≤ π0
33
x̂
π T x ≥ π0 + 1
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x̂
π T x ≤ π0
33
x̂
π T x ≥ π0 + 1
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x0
π T x ≤ π0
33
x̂
x1
π T x ≥ π0 + 1
x0
π T x ≤ π0
x̂
x1
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x0
π T x ≤ π0
33
x̂
x1
π T x ≥ π0 + 1
x0
π T x ≤ π0
x̂
x1
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x̂
π T x ≤ π0
33
x̂
π T x ≥ π0 + 1
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Generalization: Two competing approaches in pictures
Goal: build a linear OA from which a "best" cut can be deduced using CGLP.
x̂
π T x ≤ π0
33
x̂
π T x ≥ π0 + 1
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Snapshot of results
[Kılınc et al., 2011], report a speedup of 3 on a set of "hard" instances
with the NLP/LP Filimint.
[Bonami, 2011], report a speedup of 24 % on nontrivial instances with
NLP BB.
In both case, some instances not solved without these cuts are solved
within seconds.
Combination with disaggregation[Kılınç, 2011]
Even better results are obtained by combining the extended formulation trick
for separability and these cuts.
Batch
Markowitz
SLay
uflquad
34
n
10
10
14
15
Original
root gap
sol time
58.40
376.2
0.00 > 10 800
68.50
36
10.85
784
Extended
root gap sol time
68.77
58.7
98.07
1 262
86.08
5.0
96.25
145
©2015 IBM Corporation
IBM Analytics
New in CPLEX 12.6.2
Cone disaggregation for MISOCP
Lift-and-project cuts for MISOCP
Redesigned heuristic choice of most promising algorithm
Propagation of conic constraints (12.6.1).
Updated history of MIQCP with CPLEX
class
Convex QCP
convex MIQCP
–
—
–
35
p
0
>0
–
—
–
algorithm
barrier
barrier based B&B
Outer approximation B&C
Cone propagations
Major improvements
V. (Year)
9.0 (2003)
9.0 (2003)
11.0 (2007)
12.6.1 (2014)
12.6.2 (June 2015)
©2015 IBM Corporation
IBM Analytics
Cone disaggregation for MISOCP
In standard form the nonlinear constraint describing the second order cone is
not convex separable:
n
X
xi2 ≤ x02
i=1
Trick [Vielma et al., 2015], divide the constraint by x0 ≥ 0 to get a convex
separable constraint:
n
X
xi2
≤ x0 .
x0
i=1
Now introduce y1 , . . . , yn and rewrite as:
n
X
yi ≤ x0
i=1
xi2 ≤
36
x0 yi
©2015 IBM Corporation
IBM Analytics
Lift-and-project for MISOCP
Hybridizing the approaches of [Kılınc et al., 2011] and
[Bonami, 2011]
From [Kılınc et al., 2011]
x̂
Pure LP approach,
Dynamic generation of additional OA
constraints.
From [Bonami, 2011]
compact formulation using MLP,
penalization of the non-linear constraints.
37
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Lift-and-project for MISOCP
Hybridizing the approaches of [Kılınc et al., 2011] and
[Bonami, 2011]
From [Kılınc et al., 2011]
x̂
Pure LP approach,
Dynamic generation of additional OA
constraints.
From [Bonami, 2011]
compact formulation using MLP,
penalization of the non-linear constraints.
37
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Lift-and-project for MISOCP
Hybridizing the approaches of [Kılınc et al., 2011] and
[Bonami, 2011]
From [Kılınc et al., 2011]
Pure LP approach,
Dynamic generation of additional OA
constraints.
x1
x̂
x0
From [Bonami, 2011]
compact formulation using MLP,
penalization of the non-linear constraints.
37
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Lift-and-project for MISOCP
Hybridizing the approaches of [Kılınc et al., 2011] and
[Bonami, 2011]
From [Kılınc et al., 2011]
Pure LP approach,
Dynamic generation of additional OA
constraints.
x1
x̂
x0
From [Bonami, 2011]
compact formulation using MLP,
penalization of the non-linear constraints.
37
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Lift-and-project for MISOCP
Hybridizing the approaches of [Kılınc et al., 2011] and
[Bonami, 2011]
From [Kılınc et al., 2011]
Pure LP approach,
Dynamic generation of additional OA
constraints.
x1
x̂
x0
From [Bonami, 2011]
compact formulation using MLP,
penalization of the non-linear constraints.
37
π T x ≤ π0
π T x ≥ π0 + 1
©2015 IBM Corporation
IBM Analytics
Remarks
Cone disaggregation
Automatically applied by default during presolve,
Only interesting (and done):
1
2
If using the OA based algorithm,
If cones are long enough.
Lift-and-project
Also only done in OA based algorithm.
Still an expensive algorithm (may not speed up every “easy” model).
Also...
Redesigned heuristic to choose algorithm to apply in view of these changes.
38
©2015 IBM Corporation
IBM Analytics
Computational experiments MISOCP 12.6.22
Test bed
Cplex internal MIQCP test bed: 296 models
Benchmarks test set of CBLIB ( http://cblib.zib.de): 80 models
Root node experiments
Impact of combining cone disaggregation and L&P cuts by measuring the integrality
gap closed at the root node:
%gap = 100
Cut − Rel
,
Opt − Rel
Branch-and-bound results
Compare the new release Cplex 12.6.2 against Cplex 12.6.1: geo. mean of
branch-and-bound computing times.
2 All tests are carried on Linux machines: Intel X5650 @ 2.67 GHz, 24 GB RAM, 12
threads, deterministic.
39
©2015 IBM Corporation
IBM Analytics
Root node experiments in CPLEX 12.6.23
L&P cuts with all other CPLEX cuts at default
L&P cut level
Off
Default
3 All
40
Cplex test bed
No disag. Cone disag.
28.47
31.59
36.59
40.07
CBLIB
No disag. Cone disag.
12.64
23.19
34.78
45.30
entries in the table are arithmetic means of %gap closed at the root node.
©2015 IBM Corporation
IBM Analytics
Root node experiments in CPLEX 12.6.23
L&P cuts with all other CPLEX cuts at default
L&P cut level
Off
Default
Cplex test bed
No disag. Cone disag.
28.47
31.59
36.59
40.07
CBLIB
No disag. Cone disag.
12.64
23.19
34.78
45.30
L&P cuts with all other CPLEX cuts disabled
L&P cut level
Default
Aggressive
Very agg.
3 All
40
Cplex test bed
No disag. Cone disag.
23.64
26.49
32.59
34.97
36.74
39.02
CBLIB
No disag. Cone
33.47
39.92
40.70
disag.
43.13
48.96
49.47
entries in the table are arithmetic means of %gap closed at the root node.
©2015 IBM Corporation
IBM Analytics
CPLEX 12.6.1 vs 12.6.2
1.0
0.8
0.6
0.4
0.2
0
1.0
0.8
0.6
0.4
0.2
0
41
×2.86
×5.0
CPLEX test bed
CPLEX 12.6.1: 62 time limits
CPLEX 12.6.2: 38 time limits
> 0 sec
> 1 sec
231 Models 149 Models
×3.57
×6.25
CBLIB
CPLEX 12.6.1: 17 time limits
CPLEX 12.6.2: 8 time limits
> 0 sec
66 Models
> 1 sec
46 Models
©2015 IBM Corporation
IBM Analytics
A comparison of OA and SOCP-BB in CPLEX 12.6.2
4
Default strategy picked
100
OA 186 times
SOCP-BB 4 times
% model solved
80
55 models identical with
both
To be perfect should have
picked
2 more models with OA
60
40
20
1262_miqcpstrat_0
1262_miqcpstrat_1
1262_miqcpstrat_2
0
1
4 245
42
10
100
1000
time factor
10000
100000
9 more models with
SOCP-BB
models solved by at least one method and failed by none.
©2015 IBM Corporation
IBM Analytics
Reminder of CPLEX 12.6.1
5
Default strategy picked
100
OA 113 times
SOCP-BB 46 times
% model solved
80
56 models identical with
both
To be perfect should have
picked
14 more models with OA
60
40
20
1261_miqcpstrat_0
1261_miqcpstrat_1
1261_miqcpstrat_2
0
1
5 225
43
10
100
1000
time factor
10000
100000
36 more models with
SOCP-BB
models solved by at least one method and failed by none.
©2015 IBM Corporation
IBM Analytics
(MI)QP
min
1 T
x Qx + c T x
2
s.t.
Ax = b
xj ∈ Z
(MIQP)
j = 1, . . . , p
l ≤x ≤u
(with Q symmetric),
44
©2015 IBM Corporation
IBM Analytics
(MI)QP
min
1 T
x Qx + c T x
2
s.t.
(MIQP)
Ax = b
xj ∈ Z
j = 1, . . . , p
l ≤x ≤u
(with Q symmetric),
History of MIQP with CPLEX
class
Convex QP
–
convex MIQP
nonconvex QP
–
nonconvex MIQP
44
p
0
–
>0
0
–
>0
Q
0
–
0
6 0
–
6 0
algorithm
barrier
QP simplex
B&B
barrier (local)
spatial B&B (global)
spatial B&B (global)
V. (Year)
4.0 (1995)
8.0 (2002)
8.0 (2002)
12.3 (2011)
12.6 (2013)
12.6 (2013)
©2015 IBM Corporation
IBM Analytics
Example
Let G = (N, E ) be a graph and Q be the incidence matrix of G . The optimal
value of:
1
max x T Qx
2
s.t.
X
xj = 1
x ≥ 0.
1
where χ(G ) is the clique number of G [Motzkin and Straus,
is 21 1 − χ(G
)
1965],
⇒ QP is NP-hard
More generally QPs on the simplex (general Q) can be solved by a
nonlinear maximum clique algorithm [Scozzari and Tardella, 2008].
45
©2015 IBM Corporation
IBM Analytics
Local solver of nonconvex QP
Primal Dual Interior Point Algorithm.
Available since IBM CPLEX 12.3.
Not enabled by default, if Q is indefinite CPLEX will return
CPXERR_Q_NOT_POS_DEF.
Activated by setting the option optimality target to 2 (or
CPX_OPTIMALITYTARGET_FIRSTORDER).
Approach used by Ipopt but no need for
Feasibility restoration
Second order correction
Filter
Own implementation of indefinite factorization.
46
©2015 IBM Corporation
IBM Analytics
Global (MI)QP
Activated by setting optimality target to 3 (or
CPX_OPTIMALITYTARGET_OPTIMALGLOBAL ).
Note: previous versions could already solve some nonconvex MIQPs (pure
0-1 QPs, convex after presolve...)
Notes on complexity
Checking if a feasible solution is not a local minimum is coNP-Complete.
Checking if a nonconvex QP is unbounded is NP-complete.
Spatial B&B
Establish a convex (easily solvable) relaxation.
Establish branching rules on solutions of this relaxation.
47
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: Secant Approximation
The convex hull relaxations of a a square x12
y
x1 = l1
x1 = u1
{y ≤ x12 }
x1
48
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: Secant Approximation
The convex hull relaxations of a a square x12
x1 = u1
x1 = l1
Secant approximation
48
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: Secant Approximation
The convex hull relaxations of a a square x12
x12 ≤ yii+ := (l1 + u1 )x1 − l1 u1
48
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: McCormick formulas
The convex hull relaxations of a single product x1 x2 [McCormick, 1976]
x1 x2
x2
x1
49
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: McCormick formulas
The convex hull relaxations of a single product x1 x2 [McCormick, 1976]
x1 x2
(
)
u2 x1 + u1 x2 − u1 u2
x2
−
x1 x2 ≥ y12
:= max
l2 x1 + l1 x2 − l1 l2
x1
+
x1 x2 ≤ y12
:= min
49
(
u2 x1 + l1 x2 − l1 u2
l2 x1 + u1 x2 − u1 l2
)
©2015 IBM Corporation
IBM Analytics
Elementary relaxations: McCormick formulas
The convex hull relaxations of a single product x1 x2 [McCormick, 1976]
x1 x2
(
)
u2 x1 + u1 x2 − u1 u2
x2
−
x1 x2 ≥ y12
:= max
l2 x1 + l1 x2 − l1 l2
x1
+
x1 x2 ≤ y12
:= min
(
u2 x1 + l1 x2 − l1 u2
l2 x1 + u1 x2 − u1 l2
)
Depending on the sign of qij we only need y + or y − .
For simplicity, we assume we put all in the remainder.
49
©2015 IBM Corporation
IBM Analytics
Q-space reformulation and relaxation
Let Q = P + Q̃ with P the diagonal psd matrix containing qii > 0.
min
1 T
1
x Px + x T Q̃x + c T x
2
2
s.t.
(MIQP)
Ax = b
xj ∈ Z
j = 1, . . . , p
l ≤x ≤u
50
©2015 IBM Corporation
IBM Analytics
Q-space reformulation and relaxation
Let Q = P + Q̃ with P the diagonal psd matrix containing qii > 0.
Add one yij = xi xj variable for each non-zero entry qij of Q̃.
min
1 T
1
x Px + hQ̃, Y i + c T x
2
2
s.t.
Ax = b
xj ∈ Z
(MIQP)
j = 1, . . . , p
Y = xx T
l ≤x ≤u
( hQ, Y i =
50
P
i,j
qij yij )
©2015 IBM Corporation
IBM Analytics
Q-space reformulation and relaxation
Let Q = P + Q̃ with P the diagonal psd matrix containing qii > 0.
Add one yij = xi xj variable for each non-zero entry qij of Q̃.
Relax yij = xi xj using McCormick and Secant approximations.
min
1
1 T
x Px + hQ̃, Y i + c T x
2
2
s.t.
Ax = b
xj ∈ Z
j = 1, . . . , p
(q-MIQP)
yij− ≤ yij ≤ yij+
yii ≤ yii+
l ≤x ≤u
50
©2015 IBM Corporation
IBM Analytics
Factorizations of Q
Our block indefinite decomposition: M and B such that M 2-block
triangular and B 2-blocks diagonal with Q = M T BM
Reformulate x T Qx using additional variables z so that z T Dz = x T Bx
and D diagonal. Let L, D give the spectral decomposition of B, z = Lζ,
ζ = Mx.
(For simplicity assume z = Lx gives the system we want)
51
©2015 IBM Corporation
IBM Analytics
Factorized Eigenvector space reformulation and
relaxation
Use a decomposition to get z = Lx and z T Dz = x T Qx and do the same steps
as before (but more simple)....
min
1 T
z Dz + c T x
2
s.t.
Ax = b, Lx = z
xj ∈ Z
(MIQP)
j = 1, . . . , p
l ≤x ≤u
52
©2015 IBM Corporation
IBM Analytics
Factorized Eigenvector space reformulation and
relaxation
Use a decomposition to get z = Lx and z T Dz = x T Qx and do the same steps
as before (but more simple)....
Let D = D + − D − with D ± diagonal psd matrices.
min
1 T +
(z D z − z T D − z) + c T x
2
s.t.
Ax = b, Lx = z
xj ∈ Z
(MIQP)
j = 1, . . . , p
l ≤x ≤u
52
©2015 IBM Corporation
IBM Analytics
Factorized Eigenvector space reformulation and
relaxation
Use a decomposition to get z = Lx and z T Dz = x T Qx and do the same steps
as before (but more simple)....
Let D = D + − D − with D ± diagonal psd matrices.
Add yii ≤ z 2 variable for each non-zero of D − .
n
min
X dii
1 T +
z D z−
yii + c T x
2
2
i=1
s.t.
(MIQP)
Ax = b, Lx = z
xj ∈ Z
yii ≤
j = 1, . . . , p
zi2
l ≤x ≤u
52
©2015 IBM Corporation
IBM Analytics
Factorized Eigenvector space reformulation and
relaxation
Use a decomposition to get z = Lx and z T Dz = x T Qx and do the same steps
as before (but more simple)....
Let D = D + − D − with D ± diagonal psd matrices.
Add yii ≤ z 2 variable for each non-zero of D − .
Infer finite bounds, l z , u z for z and relax yii ≤ zi2 using Secant
approximations.
n
min
X dii
1 T +
z D z−
yii + c T x
2
2
i=1
s.t.
(ev-MIQP)
Ax = b, Lx = z
xj ∈ Z
j = 1, . . . , p
yii+
yii ≤
l ≤ x ≤ u, l z ≤ z ≤ u z
52
©2015 IBM Corporation
IBM Analytics
Notes on the two relaxations
The steps are almost the same.
If Q is diagonal the two relaxations are identical.
In general they are not comparable.
If Q 0, EV-space is better it preserves convexity.
Q-space gives a surpisingly good approximation [Luedtke et al., 2012]
show that, if Q has a 0 diagonal, for the box QP:
min{x T Qx : 0 ≤ x ≤ 1}:
if Q ≥ 0 the approximation is within a factor 2:
if Q 6≥ 0 the approximation is within a factor of # nnz in Q (conjecture it
is better)
Many ways to do different splittings of Q for eg. with SDP [Billionnet
et al., 2012].
CPLEX current strategy
By default, uses EV-space if problem looks almost convex.
Can be controled with parameter.
53
©2015 IBM Corporation
IBM Analytics
Branching
Let (x, y ) be the solution of the chosen QP relaxation after
presolve/cutting. And assume xj ∈ Z, j = 1, . . . , p.
If ∃y ij 6= x i x j , (x, y ) is not a solution of the problem and we need to
branch.
i
Pick such an index i, choose a value θ between li +u
and x i .
2
Branch by changing the bound to θ and updating all Secant and
McCormick approximations involving this bound.
x1 x2
x2
x1
54
©2015 IBM Corporation
IBM Analytics
Branching
Let (x, y ) be the solution of the chosen QP relaxation after
presolve/cutting. And assume xj ∈ Z, j = 1, . . . , p.
If ∃y ij 6= x i x j , (x, y ) is not a solution of the problem and we need to
branch.
i
Pick such an index i, choose a value θ between li +u
and x i .
2
Branch by changing the bound to θ and updating all Secant and
McCormick approximations involving this bound.
54
©2015 IBM Corporation
IBM Analytics
Branching
Let (x, y ) be the solution of the chosen QP relaxation after
presolve/cutting. And assume xj ∈ Z, j = 1, . . . , p.
If ∃y ij 6= x i x j , (x, y ) is not a solution of the problem and we need to
branch.
i
Pick such an index i, choose a value θ between li +u
and x i .
2
Branch by changing the bound to θ and updating all Secant and
McCormick approximations involving this bound.
x1 = θ
xi = θ
54
©2015 IBM Corporation
IBM Analytics
Branching
Let (x, y ) be the solution of the chosen QP relaxation after
presolve/cutting. And assume xj ∈ Z, j = 1, . . . , p.
If ∃y ij 6= x i x j , (x, y ) is not a solution of the problem and we need to
branch.
i
Pick such an index i, choose a value θ between li +u
and x i .
2
Branch by changing the bound to θ and updating all Secant and
McCormick approximations involving this bound.
x1 = θ
xi = θ
54
©2015 IBM Corporation
IBM Analytics
Notes on unbounded problems
Try to bound all auxiliary variables with a basic presolve.
If not possible, do it by solving LPs.
If there is an unbounded direction r look at its cost r T Qr :
If r T Qr < 0: problem is unbounded,
If r T Qr ≥ 0: relaxation is unbounded but can’t conclude on problem
status, return RELAXATION_UNBOUNDED.
(Very easy to construct examples where can’t conclude).
[Hu et al., 2012]
Propose a KKT system that detects unbounded problems correctly.
Use a combinatorial benders approach to solve it.
55
©2015 IBM Corporation
IBM Analytics
Other ingredients
Convex QP relaxation solved by a QP simplex.
Interior point solver for improving incumbents.
Bound strengthening based on the KKT system.
Linearize completely parts of the problem involving binary variables.
Heuristic detection of unbounded problems.
Multi-threaded.
BQP cuts (new in CPLEX 12.6.2)
56
©2015 IBM Corporation
IBM Analytics
Box QP
We consider the box constrained QP:
1
min x T Qx + c T x
2
s.t.
0≤x ≤1
(box-QP)
Bounds 0 and 1 are without loss of generality (every box QP can be
scaled to those bounds).
Still NP-Hard.
Has some academic interest [Vandenbussche and Nemhauser, 2005, Burer
and Vandenbussche, 2009, Chen and Burer, 2012]
Also some applications [Moré and Toraldo, 1989] (usually huge size).
57
©2015 IBM Corporation
IBM Analytics
Box QP and Boolean Quadratic Optimization
Proposition ([Burer and Letchford, 2009])
Assume that Q is without diagonal term (Qii = 0, i = 1, . . . , n).
Let Y Q be the set where variables y represent the products in Q:
Y Q = {(x, Y ) : yij = xi xj , ∀i, j such that i 6= j and qij 6= 0}.
We have
conv (x, Y ) ∈ Y Q : x ∈ [0, 1]
n
= conv (x, Y ) ∈ Y Q : x ∈ {0, 1}
n
.
Corollary
This set is the Boolean Quadratic Polytope (BQP) [Padberg, 1989]. BQP
gives a relaxation of box-QP. Every valid cut for BQP is valid for the box-QP.
58
©2015 IBM Corporation
IBM Analytics
Using BQP for nonconvex QPs [O. Günlük et al., 2015]
Boolean Quadratic Optimization in CPLEX
CPLEX has a lot of technology available to optimize over BQP.
In particular the 0 − 1/2 Chvátal Gomory cut separator finds strong
cutting planes for this polytope.
When used in the context of a non-convex QP we call these BQP-cuts.
Remarks
By scaling and shifting separation can always been reduced to 0 − 1 case.
In particular, after branching, rescaling using tighter bounds leads to
tighter local cuts.
BQP cuts arise from Box QP but can be used to strengthen any
non-convex QP.
Then the efficiency depends on how much the side constraints already
put (x, Y ) in BQP.
59
©2015 IBM Corporation
IBM Analytics
Related approaches
Separation of cuts for the box-QP (without BQP) have been developed in
Global Optimization
The McCormick formula gives the exact hull for 2 variable sets.
[Meyer and Floudas, 2005] give closed form formula for 3 variables sets.
Many approximation results on the McCormick Q-space formulation
[Coppersmith et al., 1999, Meyer and Floudas, 2005, Luedtke et al., 2012]
Exploit closed form formula for set with up to 6 variables [Misener and
Floudas, 2013].
60
©2015 IBM Corporation
IBM Analytics
Computational experiments global QP-MIQP 12.6.26
Test set
Internal nonconvex MIQP (with three variant original, 50% integer
relaxed, 100 % relaxed).
GAMS Globallib
minlp.org, Box-QP, Tardella instances...
CUTEr problems with flipped objective
Experiments
Compare CPLEX 12.6.1 to CPLEX 12.6.2.
Time limit of 3 hours, shifted geo means of computing time.
6 All tests are carried on Linux machines: Intel X5650 @ 2.67 GHz, 24 GB RAM, 12
threads, deterministic.
61
©2015 IBM Corporation
IBM Analytics
Nonconvex (MI)QP CPLEX 12.6.1 vs 12.6.2
1.0
0.8
0.6
0.4
0.2
0
1.0
0.8
0.6
0.4
0.2
0
62
×1.47
×2.33
CPLEX test bed
CPLEX 12.6.1: 270 time limits
CPLEX 12.6.2: 262 time limits
> 0 sec
> 1 sec
395 Models 134 Models
×68
×152
Box QP
CPLEX 12.6.1: 55 time limits
CPLEX 12.6.2: 19 time limits
> 0 sec
79 Models
> 1 sec
152 Models
©2015 IBM Corporation
IBM Analytics
Conclusion
CPLEX 12.6.2 is out with many good improvements for quadratic problems.
Nonconvex (MI)QP
MIQCP
1.0
0.8
0.6
0.4
0.2
0
×2.86
×5.0
> 0 sec
> 1 sec
231 Models 149 Models
1.0
0.8
0.6
0.4
0.2
0
×1.47
×2.33
> 0 sec
> 1 sec
395 Models 134 Models
We are looking for more new challenging MISOCPs and MIQPs.
63
©2015 IBM Corporation
IBM Analytics
Legal Disclaimer
© IBM Corporation 2015. All Rights Reserved.
The information contained in this publication is provided for informational purposes only. While efforts were
made to verify the completeness and accuracy of the information contained in this publication, it is provided AS
IS without warranty of any kind, express or implied. In addition, this information is based on IBM’s current
product plans and strategy, which are subject to change by IBM without notice. IBM shall not be responsible for
any damages arising out of the use of, or otherwise related to, this publication or any other materials. Nothing
contained in this publication is intended to, nor shall have the effect of, creating any warranties or
representations from IBM or its suppliers or licensors, or altering the terms and conditions of the applicable
license agreement governing the use of IBM software.
References in this presentation to IBM products, programs, or services do not imply that they will be available in
all countries in which IBM operates. Product release dates and/or capabilities referenced in this presentation
may change at any time at IBM’s sole discretion based on market opportunities or other factors, and are not
intended to be a commitment to future product or feature availability in any way. Nothing contained in these
materials is intended to, nor shall have the effect of, stating or implying that any activities undertaken by you
will result in any specific sales, revenue growth or other results.
Performance is based on measurements and projections using standard IBM benchmarks in a controlled
environment. The actual throughput or performance that any user will experience will vary depending upon
many factors, including considerations such as the amount of multiprogramming in the user’s job stream, the
I/O configuration, the storage configuration, and the workload processed. Therefore, no assurance can be given
that an individual user will achieve results similar to those stated here.
IBM, the IBM logo, CPLEX and SPSS are trademarks of International Business Machines Corporation in the
United States, other countries, or both.
Intel, Intel Xeon are trademarks or registered trademarks of Intel Corporation or its subsidiaries in the United
States and other countries.
Linux is a registered trademark of Linus Torvalds in the United States, other countries, or both. Other company,
product, or service names may be trademarks or service marks of others.
65
©2015 IBM Corporation
IBM Analytics
References I
K. Abhishek, S. Leyffer, and J. T. Linderoth. FilMINT: An outer-approximation-based solver for convex mixed-integer
nonlinear programs. INFORMS Journal on Computing, 2010. To appear, DOI: 10.1287/ijoc.1090.0373.
E. Balas. Disjunctive programming. Annals of Discrete Mathematics, 5:3–51, 1979.
E. Balas, S. Ceria, and G. Cornuéjols. A lift-and-project cutting plane algorithm for mixed 0-1 programs. Math.
Programming, 58:295–324, 1993.
A. Billionnet, S. Elloumi, and A. Lambert. Extending the qcr method to general mixed-integer programs. Mathematical
Programming, 131(1-2):381–401, 2012. doi: 10.1007/s10107-010-0381-7.
Bonami. On optimizing over lift-and-project closures. Mathematical Programming Computation, 4(2):151–179, 2012.
ISSN 1867-2949. doi: 10.1007/s12532-012-0037-0. URL http://dx.doi.org/10.1007/s12532-012-0037-0.
P. Bonami. Lift-and-project cuts for mixed integer convex programs. In IPCO, pages 52–64, 2011.
P. Bonami, L. T. Biegler, A. R. Conn, G. Cornuéjols, I. E. Grossmann, C. D. Laird, J. Lee, A. Lodi, F. Margot,
N. Sawaya, and A. Wächter. An algorithmic framework for convex mixed integer nonlinear programs. Discrete
Optimization, 5(2):186–204, 2008.
S. Burer and A. N. Letchford. On nonconvex quadratic programming with box constraints. SIAM Journal on
Optimization, 20(2):1073–1089, 2009. doi: 10.1137/080729529.
S. Burer and D. Vandenbussche. Globally solving box-constrained nonconvex quadratic programs with semidefinite-based
finite branch-and-bound. Comput Optim Appl, 43:181–195, 2009.
M. R. Bussieck and A. Drud. Sbb: A new solver for mixed integer nonlinear programming. Talk, OR 2001, Section
"Continuous Optimization", 2001.
R. H. Byrd, J. Nocedal, and R. A. Waltz. KNITRO: An integrated package for nonlinear optimization. In Large Scale
Nonlinear Optimization, 35–59, 2006, pages 35–59. Springer Verlag, 2006.
J. Chen and S. Burer. Globally solving nonconvex quadratic programming problems via completely positive programming.
Mathematical Programming Computation, 4(1):33–52, 2012.
66
©2015 IBM Corporation
IBM Analytics
References II
D. Coppersmith, O. Günlük, J. Lee, and J. Leung. A polyhedron for products of linear functions in 0/1 variables.
Technical Report RC21568, IBM Research Division, September 1999. Revised and published in: "Mixed Integer
Nonlinear Programming", S. Leyffer and J. Lee, Eds., The IMA Volumes in Mathematics and its Applications, Vol.
154. 1st Edition, 2012, pp. 513-529.
S. Drewes. Mixed Integer Second Order Cone Programming. PhD thesis, Technische Universität Darmstadt, 2009.
S. Drewes and S. Ulbrich. Subgradient based outer approximation for mixed integer second order cone programming. In
J. Lee and S. Leyffer, editors, Mixed Integer Nonlinear Programming, volume 154 of The IMA Volumes in
Mathematics and its Applications, pages 41–59. Springer New York, 2012. ISBN 978-1-4614-1926-6. doi:
10.1007/978-1-4614-1927-3_2. URL http://dx.doi.org/10.1007/978-1-4614-1927-3_2.
M. A. Duran and I. Grossmann. An outer-approximation algorithm for a class of mixed-integer nonlinear programs.
Mathematical Programming, 36:307–339, 1986.
O. K. Gupta and A. Ravindran. Branch and bound experiments in convex nonlinear integer programming. Management
Science, 31:1533–1546, 1985.
H. Hijazi, P. Bonami, and A. Ouorou. An outer-inner approximation for separable mixed-integer nonlinear programs.
INFORMS Journal on Computing, 26(1):null, 14. doi: 10.1287/ijoc.1120.0545.
J. Hu, J. E. Mitchell, and J.-S. Pang. An lpcc approach to nonconvex quadratic programs. Math. Program., 133(1-2):
243–277, 2012.
M. Kılınc, J. Linderoth, and J. Luedtke. Effective separation of disjunctive cuts for convex mixed integer nonlinear
programs. Technical Report 1681, 2011.
M. R. Kılınç. Disjunctive Cutting Planes ann Algorithms for Convex Mixed Integer Nonlinear Programming. PhD thesis,
University of Wisconsin-Madison, 2011.
J. Lee. How we participate in open source agreements, 2008. URL http://ibm.co/1nPMkAM.
S. Leyffer. Integrating SQP and branch-and-bound for mixed integer nonlinear programming. Technical report, University
of Dundee, 1998.
J. Luedtke, M. Namazifar, and J. Linderoth. Some results on the strength of relaxations of multilinear functions. Math.
Program., 136(2):325–351, 2012.
67
©2015 IBM Corporation
IBM Analytics
References III
G. P. McCormick. Computability of global solutions to factorable nonconvex programs: Part I—Convex underestimating
problems. Mathematical Programming, 10:147–175, 1976.
C. A. Meyer and C. A. Floudas. Convex envelopes for edge-concave functions. Math. Program., 103(2):207–224, 2005.
R. Misener and C. A. Floudas. GloMIQO: Global Mixed-Integer Quadratic Optimizer. J. Glob. Optim., 57:3–30, 2013.
J. Moré and G. Toraldo. Algorithms for bound constrained quadratic programming problems. Numerische Mathematik, 55
(4):377–400, 1989. ISSN 0029-599X. doi: 10.1007/BF01396045. URL http://dx.doi.org/10.1007/BF01396045.
T. S. Motzkin and E. G. Straus. Maxima for graphs and a new proof of a theorem of turán. Canadian Journal of
Mathematics, 17:533–540, 1965.
O. Günlük, J. Linderoth, and P. Bonami. In preparation, 2015.
M. Padberg. The boolean quadric polytope: Some characteristics, facets and relatives. Mathematical Programming, 45
(1-3):139–172, 1989. ISSN 0025-5610. doi: 10.1007/BF01589101. URL http://dx.doi.org/10.1007/BF01589101.
I. Quesada and I. E. Grossmann. An LP/NLP based branch–and–bound algorithm for convex MINLP optimization
problems. Computers and Chemical Engineering, 16:937–947, 1992.
A. Scozzari and F. Tardella. A clique algorithm for standard quadratic programming. Discrete Applied Mathematics, 156
(13):2439–2448, 2008.
D. Vandenbussche and G. L. Nemhauser. A polyhedral study of nonconvex quadratic programs with box constraints.
Mathematical Programming, 2005. to appear.
J. P. Vielma, I. Dunning, J. Huchette, and M. Lubin. Extended Formulations in Mixed Integer Conic Quadratic
Programming. Research Report, 2015.
S. Vigerske. Decomposition in multistage stochastic programming and a constraint integer programming approach to
mixed-integer nonlinear programming. PhD thesis, 2012.
S. Vigerske. Gams minlp benchmarks, 2013. URL http://bit.ly/1ftcys8.
T. Westerlund and K. Lundqvist. Alpha-ECP, version 5.101. an interactive minlp-solver based on the extended cutting
plane method. In Updated version of Report 01-178-A, Process Design Laboratory, Abo Akademi Univeristy, 2005.
68
©2015 IBM Corporation