Gurobi Mex: A MATLAB interface for Gurobi, Presolver

From Wikimization

(Difference between revisions)
Jump to: navigation, search
m (→Install Gurobi Mex in MATLAB)
Current revision (21:23, 3 October 2016) (edit) (undo)
(→Assumptions)
 
Line 1: Line 1:
-
[[Image:Gurobi_matlab_logos.jpg|right]] Gurobi Mex is a free MATLAB interface for '''Gurobi 2''', '''3''', and '''4''' written by [http://www.caam.rice.edu Wotao Yin] with contributions by Jon Dattorro, Imre Polik, and Tomáš Strnad. It calls [http://www.gurobi.com Gurobi] to solve linear/quadratic/mixed-integer optimization problems. Gurobi offers [http://www.gurobi.com/html/freetrial.html free trial] and [http://www.gurobi.com/html/academic.html free academic license].
+
==Introduction==
 +
''Presolving'' conventionally means quick elimination of some variables and constraints prior to numerical solution of an optimization problem.
 +
Presented with constraints <math>a^{\rm T}x_{\!}=_{\!}0\,,\,\,x\succeq0</math> for example,
 +
a presolver is likely to check whether constant vector <math>a\,</math> is positive; for if so,
 +
variable <math>\,x</math> can have only the trivial solution.
 +
The effect of such tests is to reduce the problem dimensions.
-
This interface is open source and subject to [http://creativecommons.org/licenses/by-sa/3.0/us Creative Commons Attribution-Share Alike 3.0 United States License]. It is a tool for MATLAB users to quickly call Gurobi. Its source code may also serve as a starting point for those who want to develop a customized MATLAB interface for Gurobi.
+
Most commercial optimization problem solvers incorporate presolving.
 +
Particular reductions can be proprietary or invisible,
 +
while some control or selection may be given to a user.
 +
But all presolvers have the same motivation:
 +
to make an optimization problem smaller and (ideally) easier to solve.
 +
There is profit potential because
 +
a solver can then compete more effectively in the marketplace for large-scale problems.
-
Current version 1.55 was published July 23, 2011. Versions 1.45 and higher support only '''Gurobi 4'''.
+
We present a method for reducing variable dimension based upon geometry of constraints in the problem statement:
-
[[http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi#Download.2C_Installation.2C_and_Limitations Skip to Download Release]]
+
<center><math>
-
----
+
\begin{array}{rl}
 +
\mbox{minimize}_{x\in_{}\mathbb{R}^n} & f(x)
 +
\\ \mbox{subject to} & A_{}x=b
 +
\\ & g(x)_{\!}\in\mathcal{S}
 +
\\ & x\succeq0
 +
\\ & x_{j\!}\in\mathbb{Z}\,,\qquad j\in\mathcal{J}
 +
\end{array}
 +
</math></center>
 +
where <math>A\,</math> is a real <math>m\!\times_{\!}n</math> matrix,
 +
<math>\mathbb{Z}</math> represents the integers,
 +
<math>\mathbb{R}</math> the real numbers,
 +
<math>\mathcal{S}</math> is some predetermined set,
 +
and <math>\mathcal{J}</math> is some possibly empty index set.
 +
A ''caveat'' to use of our proposed method for presolving is that it is not fast.
 +
'''('''One would incorporate this method only when a problem is too big to be solved;
 +
that is, when solver software chronically exits with error or hangs.''')'''
 +
That it can be decomposed into <math>n\,</math> independent parallel subproblems is its saving grace;
 +
speedup becomes proportional to number of parallel processors.
-
== Model ==
+
===Geometry of Constraints===
-
<pre>
+
[[Image:pyramid.jpg|thumb|right|234px|Euclidean pyramid]]
-
min/max x'Cx + c'x
+
The idea central to our presolve method is best understood geometrically.
 +
Constraints
 +
<center><math>
 +
\begin{array}{l}
 +
\\ A_{}x=b
 +
\\ x\succeq0
 +
\end{array}
 +
</math></center>
 +
suggest that a ''polyhedral cone'' comes into play.
 +
Geometers in ''convex analysis'' regard cones as convex Euclidean bodies semiinfinite in extent.
 +
Finite circular cones hold ice cream and block road traffic in daily life.
 +
Each of the great Pyramids of Egypt is a finite polyhedral cone.
 +
The geometer defines a polyhedral (semiinfinite) cone <math>\,\mathcal{K}</math> in <math>\,\mathbb{R}^m</math> as a set
 +
<center><math>\mathcal{K}=\{A_{}x~|~x\succeq0\}</math></center>
 +
that is closed and convex but not necessarily pointed (might not have a vertex).
-
subject to
+
[[Image:Polycone2.jpg|thumb|right|200px|A pointed polyhedral cone (truncated)]]
-
Ax >= / <= / = b,
+
To visualize a pointed polyhedral cone in three dimensions,
-
lb <= x <= ub,
+
imagine one Egyptian Pyramid continuing into the ground
-
SOS constraints,
+
and then indefinitely out into Space from the opposite side of Earth.
-
x(i) is continuous / binary / integer / semi-continuous / semi-integer.
+
Its four edges correspond to four columns from matrix <math>A\,</math>.
-
</pre>
+
Those four columns completely describe the semiinfinite Pyramid
 +
per definition of <math>\mathcal{K}\,</math>.
 +
But <math>A\,</math> can have more than four columns and still describe the same Pyramid.
 +
For such a fat <math>A\,</math>, each additional column resides anywhere in <math>\mathcal{K}\,</math>:
 +
either interior to the cone <math>(\{A_{}x~|~x\succ0\})</math> or on one of its ''faces'' (the vertex, an edge, or a facet).
-
== Syntax ==
+
Polyhedral cones have infinite variety. Most are not so regularly shaped as a Pyramid,
 +
and they can have any number of edges and facets.
 +
We assume a ''pointed'' polyhedral cone throughout, so there can be only one vertex
 +
(which resides at the origin in Euclidean space by definition).
-
:x = gurobi_mex(c, objtype, A, b, contypes, lb, ub, vartypes);
+
====Assumptions====
-
:x = gurobi_mex(c, objtype, A, b, contypes, lb, ub, vartypes, opts);
+
'''1) Pointedness of cones.'''&nbsp;
-
:[x,val] = gurobi_mex(...);
+
Our presolve method applies only to pointed polyhedral cones,
-
:[x,val,flag] = gurobi_mex(...);
+
namely cones that contain no ''line''.
-
:[x,val,flag,output] = gurobi_mex(...);
+
(Recall that a line is a 1-dimensional Euclidean body that extends infinitely in two opposite directions.)
-
:[x,val,flag,output,lambda] = gurobi_mex(...);
+
A cone <math>\mathcal{K}</math> contains no line if and only if the cone
 +
and its ''polar cone'' <math>-\mathcal{K}</math> intersect only at the origin:
 +
<center><math>\mathcal{K}\,\cap-\mathcal{K}=\{\mathbf{0}\}.</math></center>
 +
This condition is necessary and sufficient to establish pointedness.
 +
'''Halfspace figure as example of cone that contains line...'''
-
== Input Description ==
+
It is now prudent to ask whether a given set of ''generators'' for <math>\mathcal{K}\,</math>
-
* '''c''': objective coefficient vector, double type.
+
(the columns of <math>A</math>)
-
:[] (empty array) means uniformly <font face="Times">0</font> coefficients, and ''scalar'' means all coefficients equal to ''scalar''.
+
describe a pointed cone.
-
* '''objtype''': 1 (minimization) or -1 (maximization).
+
A practical test for pointedness therefore attempts to establish existence of a line.
-
* '''A''': constraint coefficient matrix, double type, sparse.
+
The test we choose is proof by counterexample:
-
* '''b''': constraint right-hand side vector, double type.
+
given matrix <math>A\,</math> having no zero columns,
-
:Gurobi takes a dense vector for this input. If a ''sparse'' vector is specified, it will be converted to ''full'' by Gurobi Mex.
+
<center><math>\begin{array}{rl}{\text find}&x\\
-
* '''contypes''': constraint types. Char array of '>', '<', '='.
+
\mbox{subject to}&A_{}x=\mathbf{0}\\
-
:Warning: '>=' specifies two constraints, not one.
+
&\mathbf{1}^{\rm T}x=1\\
-
:Example: '>><=' specifies four constraints. The first two constraints have ''greater or equal to'' signs, the third has ''less than or equal to'' sign, and the last is an equality constraint.
+
&x\succeq0.\end{array}
-
:If a single character is specified, all constraints have the corresponding type uniformly.
+
</math></center>
-
* '''lb''': variable lower bound, double type, either scalar or vector.
+
This test says that if <math>\mathcal{K}</math> contains a line, then some conic combination of its generators can be made to point opposite in direction to another <math>(A_{}x\!=\mathbf{0})</math>. Constraint <math>\mathbf{1}^{\rm T}x\!=\!1</math> precludes the trivial solution.
-
:[] (empty array) means <font face="Times">0</font> lower bound. -inf means no lower bound. If a scalar is specified, it will be the uniform lower bound for all variables.
+
If there is an <math>x\,</math> solving this convex feasibility problem, then the cone described by <math>A\,</math> cannot be pointed.
-
* '''ub''': variable upper bound, double type, either scalar or vector.
+
-
:[] (empty array) means no (or infinity) upper bound. If a scalar is specified, it will be the uniform upper bound for all variables.
+
-
* '''vartypes''': variable types, char array of 'C', 'B', 'I', 'S', 'N'. 'C' for continuous; 'B' for binary; 'I' for integer; 'S' for semi-continuous; 'N' for semi-integer. [] (empty array) means all variables are continuous.
+
-
:Example: 'CCCCC' specifies five continuous variables.
+
-
:Note that semi-continuous variables are variables that must take a value between their minimum and maximum or zero.
+
-
Semi-integer variables are similarly defined.
+
-
:If a single character is specified, all variables will be signed to the corresponding type uniformly.
+
-
=== [http://www.gurobi.com/doc/40/refman/node572.html Gurobi Parameters] ===
+
-
* '''opts''': optional structure that may have any number of following parameters.
+
-
** '''opts.[any Gurobi parameter]''': See [http://www.gurobi.com/doc/40/refman/node572.html Gurobi Parameters] for their allowed values.
+
-
** '''opts.QP''': Quadratic objective terms. See [http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi#Quadratic_Programming here].
+
-
** '''opts.SOS''': Special ordered set constraints. See [http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi#SOS_Constraints here].
+
-
** '''opts.Start''': MIP start vector, or Gurobi's attribute 'Start'. See [http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi#MIP_Start_Vector here].
+
-
** '''opts.TrapCtrlC''': true (to trap Ctrl-C) or false (not to trap Ctrl-C).
+
-
** '''opts.Display''': screen output level. <font face="Times">0</font> (no output); 1 (error only); 2 (default output). For complete silence, set '''opts.DisplayInterval=<font face="Times">0</font>''' and '''opts.OutputFlag=<font face="Times">0</font>''' (to silence Gurobi) and '''opts.Display=<font face="Times">0</font>''' (to make Gurobi Mex silent).
+
-
** '''opts.WriteToFile''': char array of the name of the file to which optimization data is written. See Gurobi C-Reference entry [http://www.gurobi.com/doc/40/refman/node39.html GRBwrite] for supported formats. This option helps one verify whether the model is correctly passed to Gurobi.
+
-
== Output Description ==
+
'''2) <math>\,b</math> on Boundary.'''&nbsp;
-
* '''x''': primal solution vector; empty if Gurobi encounters errors or stops early (in this case, check output flag).
+
Success of our method relies on point <math>b\,</math> belonging
-
* '''val''': optimal objective value; empty if Gurobi encounters errors or stops early.
+
to the boundary of cone <math>\mathcal{K}\,</math>.
-
* '''flag''': value meanings:
+
One way to think of a pointed closed convex cone is as a kind
-
** [ ] (empty 1x1 array) general failure
+
of coordinate system whose basis is generally nonorthogonal and overcomplete;
-
** 1 for not started
+
a conic system, very much like the familiar Cartesian system
-
** 2 for optimal
+
whose analogous cone is the first quadrant (the ''nonnegative orthant'')
-
** 3 for infeasible
+
and whose generators are the Cartesian coordinate axes.
-
** 4 for infeasible or unbounded
+
-
** 5 for unbounded
+
-
** 6 for objective worse than user-specified cutoff
+
-
** 7 for reaching iteration limit
+
-
** 8 for reaching node limit
+
-
** 9 for reaching time limit
+
-
** 10 for reaching solution limit
+
-
** 11 for user interruption
+
-
** 12 for numerical difficulties
+
-
** 13 for suboptimal solution ('''Gurobi 3''' and later)
+
-
* '''output''': structure contains the following fields
+
-
** '''output.IterCount''': number of Simplex iterations
+
-
** '''output.Runtime''': running time in seconds
+
-
** '''output.ErrorMsg''': contains Gurobi error message, if any
+
-
* '''lambda''': constraint dual variables (up to v1.50) or structure contains the following fields (v1.55 and later)
+
-
** '''lambda.RC''': reduced cost
+
-
** '''lambda.Pi''': constraint dual variables
+
-
** '''lambda.Slack''': constraint slack variables
+
 +
In the context of pointed polyhedral cones,
 +
''overcomplete'' means an abundance of linearly dependent generators in <math>A\,</math>.
 +
Consequently, nonnegative coordinates <math>x\,</math> of <math>\,b</math> are not unique;
 +
there are, generally, an infinite number of coordinates <math>x\,</math>.
 +
We can discern a smaller set of coordinate axes among the columns of <math>A\,</math> by finding
 +
the smallest face of cone <math>\mathcal{K}\,</math> to which <math>b\,</math> belongs.
 +
Then the generators of that smallest face constitute a subset of conic coordinate axes
 +
that can describe <math>b\,</math>. Once found, remaining coordinates may be assumed zero.
-
== Notes ==
+
This assumption, that <math>b\,</math> resides on the boundary <math>\partial\,</math> of cone
-
=== Quadratic Programming ===
+
<math>\mathcal{K}\,</math>, is realistic. This phenomenon occurs whenever
-
'''Gurobi 4''' and later solve quadratic programs. The quadratic terms in the objective function should be specified by '''opts.QP.qrow''', '''opts.QP.qcol''', and '''opts.QP.qval''', which correspond to the input arguments ''qrow'', ''qcol'', and ''qval'' of [http://gurobi.com/doc/40/refman/node11.html function GRBaddqpterms]. They are all 1D arrays. The first two arguments, ''qrow'' and ''qcol'', specify the row and column indices (starting from 0) of 2nd-order terms such as <math>x_1^2</math> and <math>x_1 x_2\,</math>. The third argument, ''qval'', gives their coefficients. An example is given [http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi#Example_5._Quadratic_programming below].
+
<math>b\,</math> has sparse representation in its coordinates.
-
 
+
It is therefore prudent to check, before presolving, whether <math>b\,</math> resides on the boundary:
-
To solve a quadratic program with no constraints, set A = [], b = [], and contypes = [].
+
For full-dimensional pointed closed convex cone <math>\mathcal{K}\,</math>
-
 
+
described by given matrix <math>A\,</math>
-
----
+
<center><math>\begin{array}{rl}{\text find}_{x,y}&y\not=\mathbf{0}\\
-
 
+
\mbox{subject to}&b^{\rm T}y=0\\
-
=== SOS Constraints ===
+
&A^{\rm T}y\succeq0\\
-
SOS stands for [http://www.google.com/search?q=special+ordered+set+SOS Special Ordered Sets]. This mex program uses opts.SOS.weights and opts.SOS.types to pass SOS constraints to Gurobi. '''opts.SOS.weights''' is a ''sparse'' matrix describing the weights of SOS variables, and '''opts.SOS.types''' a 1D array of type int32 or int64 (if sizeof(int) is 4 for your system, then you should use int32; if 8, use int64), which specifies the constraint types. Here is an example with 4 variables and 3 SOS constraints:
+
&A^{}x=b\\
-
 
+
&x\succeq0\end{array}
-
SOS1: x1 = 0 or x2 = 0
+
</math></center>
-
SOS1: x1 = 0 or x3 = 0
+
this linear feasibility problem has a solution iff <math>b\!\in\partial\mathcal{K}\,</math>.
-
SOS2: (x1, x3, x2, x4)
+
The requirement that <math>\mathcal{K}</math> be full-dimensional is equivalent to the condition that fat matrix <math>A\,</math> be full rank.
-
 
+
This feasibility problem is a statement of boundary membership relation from convex analysis []...
-
The corresponding code for a 32-bit system is
+
-
opts.SOS.weights = sparse([
+
-
1 1 1;
+
-
2 0 3;
+
-
0 2 2;
+
-
0 0 4]);
+
-
opts.SOS.types = int32([1 1 2]);
+
-
The ''i''th column of '''opts.SOS.weights''' specifies the weights (i.e., orders) of the variables in the ''i''th SOS constraint.
+
-
 
+
-
----
+
-
 
+
-
=== MIP Start Vector ===
+
-
To specify an MIP start vector (supported since v1.45), say x = [1 0 3 2], one can use one of the following two ways:
+
-
 
+
-
opts.Start = [1 0 3 2];
+
-
 
+
-
or
+
-
 
+
-
opts.Start.pos = int32([0 1 2 3]); % use int64 if sizeof(int) is 8 for your system
+
-
opts.Start.val = [1 0 3 2];
+
-
 
+
-
To specify start values for a subset of variables, for example to set x = [? ? -1 2], where ? means undefined, one can choose either one of the following two ways
+
-
 
+
-
GRB_UNDEFINED = 1e101;
+
-
opts.Start = [GRB_UNDEFINED GRB_UNDEFINED -1 2];
+
-
 
+
-
or
+
-
 
+
-
opts.Start.pos = int32([2 3]); % use int64 if sizeof(int) is 8 for your system
+
-
opts.Start.val = [-1 2];
+
-
 
+
-
----
+
-
 
+
-
=== How to pass a parameter from MATLAB to Gurobi? ===
+
-
v1.35 and v1.45 support all parameters of '''Gurobi 3''' and '''4''', respectively. However, if you want to specify a new or undocumented Gurobi parameter of your interest, you can DIY very easily. Suppose that we want to set a double-type parameter called '''SecretPara''' in MATLAB and pass it through this mex interface to Gurobi. Because the parameter '''TimeLimit''' has the same (double) type and it is already supported by this mex program, we can copy the code for '''TimeLimit''', replace TimeLimit by SecretPara in the code, and paste it at Line 1250 of v1.35 (or Line 510 of v1.30), just before the mex program checks unrecognized input option fields. After compiling ''gurobi.c'', the modified mex will let you assign a double value to '''opts.SecretPara'''. We compare the code for '''TimeLimit''' and '''SecretPara''' below where the differences are italicized:
+
-
 
+
-
/* Option ''TimeLimit'' */
+
-
field_n = mxGetFieldNumber(IN_OPTS, "''TimeLimit''");
+
-
if (field_n != -1) {
+
-
field = mxGetFieldByNumber(IN_OPTS,0,field_n);
+
-
bOpts[field_n] = true;
+
-
if (!mxIsDouble(field) || mxIsComplex(field) || mxIsEmpty(field)) {
+
-
mexPrintf("''Option TimeLimit must be real positive double (0 to inf).''");
+
-
goto QUIT;
+
-
}
+
-
error = GRBsetdblparam(env, "''TimeLimit''", mxGetScalar(field));
+
-
if (error) goto QUIT;
+
-
}
+
-
 
+
-
/* Option ''SecretPara'' */
+
-
field_n = mxGetFieldNumber(IN_OPTS, "''SecretPara''");
+
-
if (field_n != -1) {
+
-
field = mxGetFieldByNumber(IN_OPTS,0,field_n);
+
-
bOpts[field_n] = true;
+
-
if (!mxIsDouble(field) || mxIsComplex(field) || mxIsEmpty(field)) {
+
-
mexPrintf("''Option SecretPara must real double (?? through ??).''");
+
-
goto QUIT;
+
-
}
+
-
error = GRBsetdblparam(env, "''SecretPara''", mxGetScalar(field));
+
-
if (error) goto QUIT;
+
-
}
+
-
 
+
-
Note that you must start from a parameter of the ''same'' type (int, double, or string). In case memory allocation is needed, use ''mxCalloc'' and make sure that ''mxFree'' has been called whenever the mex program exits, normally or not.
+
-
 
+
-
== Callbacks ==
+
-
Callbacks are useful to obtain the progress of Gurobi (e.g., by calling GRBcbget) and to modify its behavior during runtime (e.g., by calling GRBcbcut and GRBcbsolution).
+
-
 
+
-
Gurobi Mex implements a callback function ''mycallback'' to obtain Gurobi's progress messages and print them on the MATLAB screen. The print frequency is set in '''opts.DisplayInterval''' (in seconds). The same function is also used to detect user input Ctrl-C.
+
-
 
+
-
Information for Gurobi callbacks can be found [http://www.gurobi.com/doc/40/refman/node84.html here] in Gurobi's help. An example can be found [http://www.gurobi.com/doc/40/examples/node8.html here].
+
-
 
+
-
== Examples ==
+
-
 
+
-
=== Example 1. Linear programming ===
+
-
This example is borrowed from MATLAB's linprog help.
+
-
 
+
-
Problem:
+
-
<pre>
+
-
min –5 x1 – 4 x2 –6 x3,
+
-
 
+
-
subject to
+
-
x1 – x2 + x3 ≤ 20
+
-
3 x1 + 2 x2 + 4 x3 ≤ 42
+
-
3 x1 + 2 x2 ≤ 30
+
-
0 ≤ x1, 0 ≤ x2, 0 ≤ x3.
+
-
</pre>
+
-
MATLAB code:
+
-
<pre>
+
-
c = [-5; -4; -6];
+
-
objtype = 1;
+
-
A = sparse([1 -1 1; 3 2 4; 3 2 0]);
+
-
b = [20; 42; 30];
+
-
lb = zeros(3,1); % same as lb = [];
+
-
ub = [];
+
-
contypes = '<<<';
+
-
vtypes = []; % same as vtypes = 'CCC'; [] means 'C...C'
+
-
 
+
-
clear opts
+
-
opts.IterationLimit = 20;
+
-
opts.FeasibilityTol = 1e-6;
+
-
opts.IntFeasTol = 1e-5;
+
-
opts.OptimalityTol = 1e-6;
+
-
opts.Method = 1; % 0 - primal, 1 - dual
+
-
opts.Presolve = -1; % -1 - auto, 0 - no, 1 - conserv, 2 - aggressive
+
-
opts.Display = 1;
+
-
opts.LogFile = 'test_gurobi_mex_LP.log'; % optional
+
-
opts.WriteToFile = 'test_gurobi_mex_LP.mps'; % optional; it can cause a long delay if problem is large
+
-
 
+
-
[x,val,exitflag,output,lambda] = gurobi_mex(c,objtype,A,b,contypes,lb,ub,vtypes,opts);
+
-
</pre>
+
-
Results:
+
-
<pre>
+
-
x' =
+
-
0 15 3
+
-
 
+
-
val =
+
-
-78
+
-
 
+
-
exitflag =
+
-
2
+
-
 
+
-
output =
+
-
IterCount: 2
+
-
Runtime: 0
+
-
ErrorMsg: []
+
-
 
+
-
lambda' =
+
-
0 -1.5000 -0.5000
+
-
</pre>
+
-
Log file: ''test_gurobi_mex_LP.log''. MPS file: ''test_gurobi_mex_LP.mps''.
+
-
 
+
-
=== Example 2. Integer programming ===
+
-
This example is borrowed from ''mip1_c.c'' of '''Gurobi 2'''.
+
-
 
+
-
Problem:
+
-
<pre>
+
-
max x + y + 2z,
+
-
 
+
-
subject to
+
-
x + 2 y + 3 z <= 4
+
-
x + y >= 1
+
-
x, y, z binary.
+
-
</pre>
+
-
MATLAB code:
+
-
<pre>
+
-
c = [1; 1; 2];
+
-
objtype = -1; % 1 for minimize, -1 for maximize
+
-
A = sparse([1 2 3; 1 1 0]);
+
-
b = [4; 1];
+
-
lb = [];
+
-
ub = [];
+
-
contypes = '<>';
+
-
vtypes = 'BBB';
+
-
 
+
-
clear opts
+
-
opts.IterationLimit = 20;
+
-
opts.FeasibilityTol = 1e-6;
+
-
opts.IntFeasTol = 1e-5;
+
-
opts.OptimalityTol = 1e-6;
+
-
opts.Method = 1; % 0 - primal, 1 - dual
+
-
opts.Presolve = -1; % -1 - auto, 0 - no, 1 - conserv, 2 - aggressive
+
-
opts.Display = 1;
+
-
opts.LogFile = 'test_gurobi_mex_MIP.log';
+
-
opts.WriteToFile = 'test_gurobi_mex_MIP.mps'; % this option can cause a long delay if problem is large
+
-
 
+
-
[x,val,exitflag,output] = gurobi_mex(c,objtype,A,b,contypes,lb,ub,vtypes,opts);
+
-
</pre>
+
-
Gurobi does not give lambda (Pi, or Lagrange multipliers) for MIPs, unless model fix is called.
+
-
 
+
-
Results:
+
-
<pre>
+
-
disp('Solution:');disp(x')
+
-
disp('Optimal obj value:');disp(val)
+
-
disp('Exit flag:');disp(exitflag)
+
-
disp('Optimization info:');disp(output)
+
-
 
+
-
Solution:
+
-
1 0 1
+
-
 
+
-
Optimal obj value:
+
-
3
+
-
 
+
-
Exit flag:
+
-
2
+
-
 
+
-
Optimization info:
+
-
IterCount: 0
+
-
Runtime: 0
+
-
ErrorMsg: []
+
-
</pre>
+
-
Log file: ''test_gurobi_mex_MIP.log''. MPS file: ''test_gurobi_mex_MIP.mps''.
+
-
 
+
-
=== Example 3. Feasibility test ===
+
-
Problem:
+
-
<pre>
+
-
Find a solution or report infeasibility of
+
-
 
+
-
5 x1 + 4 x2 + 5 x4 >= -21
+
-
5 x1 + 3 x2 + 1 x3 + 4 x4 = -14
+
-
3 x1 + 5 x2 + 2 x3 - 5 x4 = 11
+
-
x1,x2,x3,x4 >= 0.
+
-
</pre>
+
-
MATLAB code:
+
-
<pre>
+
-
c = []; % use [] or 0 for null objective
+
-
objtype = -1; % 1 for minimize, -1 for maximize
+
-
A = sparse([5 4 0 5; 5 3 1 4; 3 5 2 -5]);
+
-
b = [-21; -14; 11];
+
-
lb = []; % stands for 0 lower bound uniformly
+
-
ub = []; % stands for inf upper bound uniformly
+
-
contypes = '>==';
+
-
vtypes = []; % same as vtypes = 'CCCC'
+
-
 
+
-
clear opts
+
-
opts.FeasibilityTol = 1e-6;
+
-
opts.Presolve = -1; % -1 - auto, 0 - no, 1 - conserv, 2 - aggressive
+
-
opts.Display = 1;
+
-
opts.LogFile = 'test_gurobi_mex_Feasibility.log';
+
-
opts.WriteToFile = 'test_gurobi_mex_Feasibility.mps'; % this option can cause a long delay if problem is large
+
-
 
+
-
[x,val,exitflag,output] = gurobi_mex(c,objtype,A,b,contypes,lb,ub,vtypes,opts);
+
-
</pre>
+
-
 
+
-
Results:
+
-
<pre>
+
-
disp('Solution:');disp(x')
+
-
disp('Optimal obj value:');disp(val)
+
-
disp('Exit flag:');disp(exitflag)
+
-
 
+
-
Model is infeasible. No solution will be returned.
+
-
 
+
-
Solution:
+
-
Optimal obj value:
+
-
Exit flag:
+
-
3
+
-
</pre>
+
-
Log file: ''test_gurobi_mex_Feasibility.log''. MPS file: ''test_gurobi_mex_Feasibility.mps''.
+
-
 
+
-
=== Example 4. SOS constraint test ===
+
-
Problem:
+
-
<pre>
+
-
min –3 x1 – 1 x2 – 1 x3,
+
-
 
+
-
subject to
+
-
x1 + x2 + x3 <= 2,
+
-
0 <= x1 <= 1, 0 <= x2 <= 1, 0 <= x3 <= 2,
+
-
SOS type 1: x1 = 0 or x2 = 0,
+
-
SOS type 1: x1 = 0 or x3 = 0.
+
-
</pre>
+
-
MATLAB code:
+
-
<pre>
+
-
c = [-3; -1; -1];
+
-
objtype = 1;
+
-
A = sparse([1 1 1]);
+
-
b = 2;
+
-
lb = []; % means 0 lower bound
+
-
ub = [1 1 2];
+
-
contypes = '<';
+
-
vtypes = []; % same as vtypes = 'CCC'
+
-
sos.weights = sparse([1 1; 2 0; 0 2]);
+
-
sos.types = int32([1 1]); % Both the SOS constraints are of type 1
+
-
 
+
-
clear opts
+
-
opts.FeasibilityTol = 1e-6;
+
-
opts.Presolve = -1; % -1 - auto, 0 - no, 1 - conserv, 2 - aggressive
+
-
opts.Display = 1;
+
-
opts.LogFile = 'test_gurobi_mex_Feasibility.log';
+
-
opts.WriteToFile = 'test_gurobi_mex_Feasibility.mps';
+
-
 
+
-
[x,val,exitflag,output] = gurobi_mex(c,objtype,A,b,contypes,lb,ub,vtypes,opts);
+
-
</pre>
+
-
 
+
-
Results:
+
-
<pre>
+
-
Optimal solution found (tolerance 1.00e-04)
+
-
Best objective -3.0000000000e+00, best bound -3.0000000000e+00, gap 0.0%
+
-
+
-
Solution:
+
-
1 0 0
+
-
</pre>
+
-
Log file: ''test_gurobi_mex_SOS.log''. MPS file: ''test_gurobi_mex_SOS.mps''.
+
-
 
+
-
=== Example 5. Quadratic programming ===
+
-
This example appears in MATLAB Help entry for ''quadprog''.
+
-
 
+
-
Problem:
+
-
<pre>
+
-
min 0.5 x0^2 - x0*x1 + x1^2 - 2*x0 - 6*x1,
+
-
 
+
-
subject to
+
-
x0 + x1 <= 2,
+
-
-x0 + 2x1 <= 2,
+
-
2x0 + x1 <= 3,
+
-
x0 >= 0, x1 >= 0.
+
-
</pre>
+
-
MATLAB code:
+
-
<pre>
+
-
clear opts
+
-
 
+
-
c = [-2 -6]; % objective linear term
+
-
objtype = 1; % minimization
+
-
A = sparse([1 1; -1 2; 2 1]); % constraint coefficients
+
-
b = [2; 2; 3]; % constraint right-hand side
+
-
lb = []; % [] means 0 lower bound
+
-
ub = []; % [] means inf upper bound
+
-
contypes = '<<<'; % three <= inequalities
+
-
vtypes = []; % [] means all variables are continuous
+
-
% note that 0.5 x0^2 - x0*x1 + x1^2 = 0.5 x0*x0 - x0*x1 + x1*x1
+
-
opts.QP.qrow = int32([0 0 1]); % indices of (x0), (x0), (x1) in 0.5 (x0)*x0 - (x0)*x1 + (x1)*x1; use int64 if sizeof(int) is 8 for your system
+
-
opts.QP.qcol = int32([0 1 1]); % indices of (x0), (x1), (x1) in 0.5 x0*(x0) - x0*(x1) + x1*(x1); use int64 if sizeof(int) is 8 for your system
+
-
opts.QP.qval = [0.5 -1 1]; % coefficients of 0.5 x0^2 - x0*x1 + x1^2
+
-
 
+
-
opts.IterationLimit = 200;
+
-
opts.FeasibilityTol = 1e-6;
+
-
opts.IntFeasTol = 1e-5;
+
-
opts.OptimalityTol = 1e-6;
+
-
 
+
-
[x,val,exitflag,output,lambda] = gurobi_mex(c,objtype,A,b,contypes,lb,ub,vtypes,opts);
+
-
 
+
-
disp('Solution:');disp(x')
+
-
disp('Optimal obj value:');disp(val)
+
-
disp('Exit flag:');disp(exitflag)
+
-
disp('Optimization info:');disp(output)
+
-
disp('Lagrange multiplers:');disp(lambda')
+
-
</pre>
+
-
 
+
-
Results:
+
-
<pre>
+
-
Solution:
+
-
0.6667 1.3333
+
-
Optimal obj value:
+
-
-8.2222
+
-
Exit flag:
+
-
2
+
-
Optimization info:
+
-
IterCount: 4
+
-
Runtime: 0.0630
+
-
ErrorMsg: []
+
-
Lagrange multiplers:
+
-
-3.1111 -0.4444 0
+
-
</pre>
+
-
 
+
-
=== Example 6. Compressive sensing ===
+
-
See example m-file ''test_gurobi_mex_CS.m''.
+
-
 
+
-
== Feedback ==
+
-
 
+
-
[http://www.convexoptimization.com/wikimization/index.php/Special:Emailuser/Wotao.yin Wotao Yin] would be delighted to hear from you if you find Gurobi Mex useful, or if you have any suggestions, contributions, or bug reports.
+
-
 
+
-
== How to cite ==
+
-
 
+
-
Wotao Yin. ''Gurobi Mex: A MATLAB interface for Gurobi'', URL: http://convexoptimization.com/wikimization/index.php/gurobi_mex, 2009-2011.
+
-
 
+
-
== Download, Installation, and Limitations ==
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/download_request.html v1.55] ('''Gurobi 4''') latest version C source code and MATLAB examples.
+
-
 
+
-
New features: support the new options introduced in Gurobi v4.5; support the output of reduced costs and constraint slacks.
+
-
 
+
-
=== History ===
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.50.zip v1.50] ('''Gurobi 4''') New features: quadratic programming with no linear constraints. Fixed bugs on handling SOS constraints on 64-bit systems.
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.45.zip v1.45] ('''Gurobi 4''') New features: quadratic programming, MIP start vector. Fixed a bug on reporting unsupported options. Dattorro added support for all Gurobi options. Gurobi 4 changed parameter name "LPMethod" to "Method".
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.35.zip v1.35] '''(Gurobi 2&3)''' New features: support of [http://www.google.com/search?q=special+ordered+set+SOS Special Ordered Sets (SOS)] constraints of types 1 and 2; support all Gurobi parameters and a new option TrapCtrlC; detection of unrecognized options. Fixed minor display issues. Dattorro added support for all Gurobi options.
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.30.zip v1.30] Short release for SOS support.
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.20.zip v1.20] New features: [http://www.caam.rice.edu/~wy1/links/mex_ctrl_c_trick/ Ctrl-C detection], '''Gurobi 3''' support. Improved error and exception handling. Empty array [] is returned if an output argument is not available. Fixed the display interval option.
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.10.zip v1.10] New features: callback, runtime progress output, flexible log file, flexible input types, more options. Part of the code was contributed by ''Tomáš Strnad''.
+
-
Known bug: print an empty line even if options '''DisplayInterval''' and '''Display''' are both set to
+
-
<font face="Times">0</font>. Fix: remove Line 736 of ''gurobi_mex.c'': mexPrintf("\n");
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.05.zip v1.05] Major bug fix: char array of constraint sense has been fixed
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.04.zip v1.04] support writing model to files in various formats such as MPS, REW, LP, ...
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.03.zip v1.03] support log file
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.02.zip v1.02] fixed a memory leak issue
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.01.zip v1.01] update: support output dual solution lambda; allow vartypes to be empty (for all continuous variables).
+
-
 
+
-
[http://www.caam.rice.edu/~wy1/gurobi_mex/gurobi_mex_v1.0.zip v1.00] initial version.
+
-
 
+
-
=== Future Release ===
+
-
Please send your suggested features to [http://www.convexoptimization.com/wikimization/index.php/Special:Emailuser/Wotao.yin Wotao Yin].
+
-
 
+
-
=== Install Gurobi Mex in MATLAB ===
+
-
 
+
-
==== Step 1. Install Gurobi ====
+
-
Download and install Gurobi. Refer to Gurobi's installation guide. Make sure that (i) proper environment variables are set, and (ii) your copy of Gurobi has a valid license.
+
-
 
+
-
==== Step 2. Install C++ Compiler ====
+
-
Download and install a supported C/C++ compiler for your copy of MATLAB. Do not use the built-in compiler ''lcc'', which cannot link with Gurobi's library.
+
-
 
+
-
This step is typically easy for Linux systems but not as easy for Windows 32/64-bit systems. For Windows, one often uses Microsoft's Visual C++ compiler. If the compiler is not installed, you can download and install the free Microsoft Visual C++ Express. On 64-bit Windows, Windows SDK is also needed. Google for "Windows SDK" and you will find a Microsoft webpage from where you can download and install a SDK. Note that only the compiler part of the SDK is needed.
+
-
 
+
-
Once a compiler is installed, run ''mex -setup'' in MATLAB and it shall automatically locates the compiler and generates a configuration file.
+
-
 
+
-
==== Step 3. Compiling Gurobi Mex ====
+
-
 
+
-
Preparation: If you are going to use '''Gurobi 2''', uncomment Line 20 of gurobi_mex.c regarding GRB_SUBOPTIMAL. '''opts.Method''' in '''Gurobi 4''' corresponds to parameter '''opts.LPMethod''' in '''Gurobi 2''' and '''3'''.
+
-
 
+
-
In MATLAB, go to the folder where ''gurobi_mex.c'' is saved and call ''mex'' as follows:
+
-
 
+
-
===== Step 3a. Windows =====
+
-
 
+
-
mex -O -I"<gurobi include path>" "<gurobi_mex.c>" "<gurobi C library file>" "<MATLAB libut.lib>"
+
-
 
+
-
For 64-bit MATLAB, add option "-largeArrayDims".
+
-
 
+
-
Example with '''Gurobi 4.51''', MSVC2010 Express, MATLAB 2011a, and 64-bit Windows
+
-
 
+
-
mex -O -largeArrayDims -I"C:\Gurobi451\win64\include" "C:\folder\gurobi_mex.c" "C:\Gurobi451\win64\lib\gurobi45.lib" "C:\Program Files\MATLAB\R2011a\extern\lib\win64\microsoft\libut.lib"
+
-
 
+
-
===== Missing libut.lib? =====
+
-
Ctrl-C detection requires libut.lib. If it is not found under your copy of MATLAB, you can download one for [http://www.caam.rice.edu/~wy1/gurobi_mex/libut_32bit/libut.lib 32-bit Windows] and [http://www.caam.rice.edu/~wy1/gurobi_mex/libut_64bit/libut.lib 64-bit Windows] (courtesy of Imre Polik).
+
-
 
+
-
Alternatively, libut.lib can be manually generated by creating a .def text file including the following five lines
+
-
LIBRARY libut.dll
+
-
EXPORTS
+
-
utIsInterruptPending
+
-
utSetInterruptPending
+
-
[empty line]
+
-
and then calling lib.exe (included in MSVC) like
+
-
"C:\Program Files\Microsoft Visual Studio 9.0\VC\bin\lib" /def:libut.def /out:libut.lib /machine:x86
+
-
where /machine:x86 should be replaced by /machine:x64 for 64-bit Windows.
+
-
 
+
-
===== Step 2b. Linux/Unix =====
+
-
 
+
-
mex -O -I"<gurobi include path>" "<gurobi_mex.c>" -L"<gurobi lib path>" -l<gurobi C library file> -lut
+
-
 
+
-
For 64-bit MATLAB, add option "-largeArrayDims".
+
-
 
+
-
Example with '''Gurobi 3''', GCC, MATLAB 2009B, and 64-bit Linux
+
-
 
+
-
mex -O -I"/opt/gurobi300/linux64/include" "/home/wotao/gurobi mex/gurobi_mex.c" -L"/opt/gurobi300/linux64/lib" -lgurobi30 -lut -largeArrayDims
+
-
 
+
-
==== Tested platforms ====
+
-
 
+
-
* 64-bit Ubuntu Linux 9.10, 64-bit MATLAB, and gcc.
+
-
* 32-bit Windows, 32-bit MATLAB, and gcc (part of free [http://gnumex.sourceforge.net Mingw/GnuMex], alternatively [http://tdm-gcc.tdragon.net TDM-GCC]).
+
-
* 32-bit Windows, 32-bit MATLAB, and MSVC 2008 SP1 (the express Edition is free).
+
-
* 64-bit Windows, 64-bit MATLAB, and MSVC 2008 SP1 (the express Edition is free). Courtesy of Imre Polik.
+
-
* 64-bit Windows, 64-bit MATLAB, and MSVC 2010 Express.
+
-
 
+
-
For 64-bit MATLAB, Jon Dattorro suggests Microsoft SDK and a [http://www.mathworks.com/support/solutions/en/data/1-8FJXQE/index.html?solution=1-8FJXQE bug fix] for the linker.
+
-
 
+
-
== FAQs ==
+
-
=== compiling is successful, but gurobi_mex cannot run ===
+
-
 
+
-
'''Solution''': check and correct Gurobi license and environment variables
+
-
 
+
-
Step 1. Check and validate [http://www.gurobi.com/doc/40/quickstart/node2.html Gurobi license]
+
-
 
+
-
Step 2. Check [http://gurobi.com/doc/40/quickstart/node1.html system environment variables for Gurobi]
+
-
 
+
-
Step 3. Verify MATLAB knows the correct system environment variables by running
+
-
>> getenv('GUROBI_HOME')
+
-
>> getenv('GRB_LICENSE_FILE')
+
-
>> getenv('PATH')
+
-
>> getenv('LD_LIBRARY_PATH') % on Unix/Linux/Mac
+
-
You may need to restart MATLAB '''from the terminal''' to get all environment variables loaded to MATLAB.
+
-
 
+
-
=== "int32" or "int64" errors ===
+
-
 
+
-
'''Solution''': use int32 if sizeof(int) is 4 for your system; use int64 if sizeof(int) is 8. To determine sizeof(int), take the following steps
+
-
 
+
-
Step 1. create "check_sizeof_int.c" with the following lines
+
-
 
+
-
#include "mex.h"
+
-
void mexFunction(
+
-
int nlhs, mxArray *plhs[],
+
-
int nrhs, const mxArray *prhs[]
+
-
)
+
-
{
+
-
mexPrintf("Size of int is %d\n", sizeof(int));
+
-
return;
+
-
}
+
-
 
+
-
Step 2. Launch Matlab, run '''mex check_sizeof_int.c''', and then run '''check_sizeof_int'''
+
-
 
+
-
=== MATLAB reports "out of memory" ===
+
-
 
+
-
'''Solution''': run ''clear mex'' after each call to gurobi_mex
+
-
 
+
-
Running out of memory is often the result of memory leaks. However, the interface has been checked numerous times for memory leaks. If there still appears to be a leak, we are not sure if it is with the interface, Gurobi, or MATLAB itself.
+
-
 
+
-
=== Gurobi Mex in a loop returns incorrect solutions ===
+
-
 
+
-
'''Solution''': run ''clear mex'' after each call to gurobi_mex
+
-
 
+
-
Even though gurobi_mex releases all memory allocated before it exits, something seems to stay in memory since the first call and gets carried over to all the subsequent calls. We are aware of this problem but we cannot find a better solution than executing ''clear mex'' after each call to gurobi_mex.
+
-
 
+
-
== License ==
+
-
 
+
-
Creative Commons Attribution-Share Alike 3.0 United States License.
+
-
Allow commercial use of this work. Permit others to copy, distribute, display, and perform the work, including for commercial purposes.
+
-
Allow modification, as long as others share alike. Permit others to distribute derivative works only under the same license or one compatible with the one that governs the licensor's work.
+

Current revision

Introduction

Presolving conventionally means quick elimination of some variables and constraints prior to numerical solution of an optimization problem. Presented with constraints LaTeX: a^{\rm T}x_{\!}=_{\!}0\,,\,\,x\succeq0 for example, a presolver is likely to check whether constant vector LaTeX: a\, is positive; for if so, variable LaTeX: \,x can have only the trivial solution. The effect of such tests is to reduce the problem dimensions.

Most commercial optimization problem solvers incorporate presolving. Particular reductions can be proprietary or invisible, while some control or selection may be given to a user. But all presolvers have the same motivation: to make an optimization problem smaller and (ideally) easier to solve. There is profit potential because a solver can then compete more effectively in the marketplace for large-scale problems.

We present a method for reducing variable dimension based upon geometry of constraints in the problem statement:

LaTeX: <pre> \begin{array}{rl} \mbox{minimize}_{x\in_{}\mathbb{R}^n} & f(x) \\ \mbox{subject to} & A_{}x=b \\ & g(x)_{\!}\in\mathcal{S} \\ & x\succeq0 \\ & x_{j\!}\in\mathbb{Z}\,,\qquad j\in\mathcal{J} \end{array} </pre>

where LaTeX: A\, is a real LaTeX: m\!\times_{\!}n matrix, LaTeX: \mathbb{Z} represents the integers, LaTeX: \mathbb{R} the real numbers, LaTeX: \mathcal{S} is some predetermined set, and LaTeX: \mathcal{J} is some possibly empty index set.

A caveat to use of our proposed method for presolving is that it is not fast. (One would incorporate this method only when a problem is too big to be solved; that is, when solver software chronically exits with error or hangs.) That it can be decomposed into LaTeX: n\, independent parallel subproblems is its saving grace; speedup becomes proportional to number of parallel processors.

Geometry of Constraints

Euclidean pyramid
Euclidean pyramid

The idea central to our presolve method is best understood geometrically. Constraints

LaTeX: <pre> \begin{array}{l} \\ A_{}x=b \\ x\succeq0 \end{array} </pre>

suggest that a polyhedral cone comes into play. Geometers in convex analysis regard cones as convex Euclidean bodies semiinfinite in extent. Finite circular cones hold ice cream and block road traffic in daily life. Each of the great Pyramids of Egypt is a finite polyhedral cone. The geometer defines a polyhedral (semiinfinite) cone LaTeX: \,\mathcal{K} in LaTeX: \,\mathbb{R}^m as a set

LaTeX: \mathcal{K}=\{A_{}x~|~x\succeq0\}

that is closed and convex but not necessarily pointed (might not have a vertex).

A pointed polyhedral cone (truncated)
A pointed polyhedral cone (truncated)

To visualize a pointed polyhedral cone in three dimensions, imagine one Egyptian Pyramid continuing into the ground and then indefinitely out into Space from the opposite side of Earth. Its four edges correspond to four columns from matrix LaTeX: A\,. Those four columns completely describe the semiinfinite Pyramid per definition of LaTeX: \mathcal{K}\,. But LaTeX: A\, can have more than four columns and still describe the same Pyramid. For such a fat LaTeX: A\,, each additional column resides anywhere in LaTeX: \mathcal{K}\,: either interior to the cone LaTeX: (\{A_{}x~|~x\succ0\}) or on one of its faces (the vertex, an edge, or a facet).

Polyhedral cones have infinite variety. Most are not so regularly shaped as a Pyramid, and they can have any number of edges and facets. We assume a pointed polyhedral cone throughout, so there can be only one vertex (which resides at the origin in Euclidean space by definition).

Assumptions

1) Pointedness of cones.  Our presolve method applies only to pointed polyhedral cones, namely cones that contain no line. (Recall that a line is a 1-dimensional Euclidean body that extends infinitely in two opposite directions.) A cone LaTeX: \mathcal{K} contains no line if and only if the cone and its polar cone LaTeX: -\mathcal{K} intersect only at the origin:

LaTeX: \mathcal{K}\,\cap-\mathcal{K}=\{\mathbf{0}\}.

This condition is necessary and sufficient to establish pointedness. Halfspace figure as example of cone that contains line...

It is now prudent to ask whether a given set of generators for LaTeX: \mathcal{K}\, (the columns of LaTeX: A) describe a pointed cone. A practical test for pointedness therefore attempts to establish existence of a line. The test we choose is proof by counterexample: given matrix LaTeX: A\, having no zero columns,

LaTeX: \begin{array}{rl}{\text find}&x\\ <p>\mbox{subject to}&A_{}x=\mathbf{0}\\ &\mathbf{1}^{\rm T}x=1\\ &x\succeq0.\end{array} </p>

This test says that if LaTeX: \mathcal{K} contains a line, then some conic combination of its generators can be made to point opposite in direction to another LaTeX: (A_{}x\!=\mathbf{0}). Constraint LaTeX: \mathbf{1}^{\rm T}x\!=\!1 precludes the trivial solution. If there is an LaTeX: x\, solving this convex feasibility problem, then the cone described by LaTeX: A\, cannot be pointed.

2) LaTeX: \,b on Boundary.  Success of our method relies on point LaTeX: b\, belonging to the boundary of cone LaTeX: \mathcal{K}\,. One way to think of a pointed closed convex cone is as a kind of coordinate system whose basis is generally nonorthogonal and overcomplete; a conic system, very much like the familiar Cartesian system whose analogous cone is the first quadrant (the nonnegative orthant) and whose generators are the Cartesian coordinate axes.

In the context of pointed polyhedral cones, overcomplete means an abundance of linearly dependent generators in LaTeX: A\,. Consequently, nonnegative coordinates LaTeX: x\, of LaTeX: \,b are not unique; there are, generally, an infinite number of coordinates LaTeX: x\,. We can discern a smaller set of coordinate axes among the columns of LaTeX: A\, by finding the smallest face of cone LaTeX: \mathcal{K}\, to which LaTeX: b\, belongs. Then the generators of that smallest face constitute a subset of conic coordinate axes that can describe LaTeX: b\,. Once found, remaining coordinates may be assumed zero.

This assumption, that LaTeX: b\, resides on the boundary LaTeX: \partial\, of cone LaTeX: \mathcal{K}\,, is realistic. This phenomenon occurs whenever LaTeX: b\, has sparse representation in its coordinates. It is therefore prudent to check, before presolving, whether LaTeX: b\, resides on the boundary: For full-dimensional pointed closed convex cone LaTeX: \mathcal{K}\, described by given matrix LaTeX: A\,

LaTeX: \begin{array}{rl}{\text find}_{x,y}&y\not=\mathbf{0}\\ <p>\mbox{subject to}&b^{\rm T}y=0\\ &A^{\rm T}y\succeq0\\ &A^{}x=b\\ &x\succeq0\end{array} </p>

this linear feasibility problem has a solution iff LaTeX: b\!\in\partial\mathcal{K}\,. The requirement that LaTeX: \mathcal{K} be full-dimensional is equivalent to the condition that fat matrix LaTeX: A\, be full rank. This feasibility problem is a statement of boundary membership relation from convex analysis []...

Retrieved from "http://www.convexoptimization.com/wikimization/index.php/Gurobi_Mex:_A_MATLAB_interface_for_Gurobi"
Personal tools