MIQPS_MATPOWER Mixed Integer Quadratic Program Solver for MATPOWER.
[X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
MIQPS_MATPOWER(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT)
A common wrapper function for various QP solvers.
Solves the following QP (quadratic programming) problem:
min 1/2 X'*H*X + C'*X
X
subject to
L <= A*X <= U (linear constraints)
XMIN <= X <= XMAX (variable bounds)
Inputs (all optional except H, C, A and L):
H : matrix (possibly sparse) of quadratic cost coefficients
C : vector of linear cost coefficients
A, L, U : define the optional linear constraints. Default
values for the elements of L and U are -Inf and Inf,
respectively.
XMIN, XMAX : optional lower and upper bounds on the
X variables, defaults are -Inf and Inf, respectively.
X0 : optional starting value of optimization vector X
VTYPE : character string of length NX (number of elements in X),
or 1 (value applies to all variables in x),
allowed values are 'C' (continuous), 'B' (binary),
'I' (integer), 'S' (semi-continuous), or 'N' (semi-integer).
(MOSEK, GLPK, OT allow only 'C', 'B', or 'I')
OPT : optional options structure with the following fields,
all of which are also optional (default values shown in
parentheses)
alg ('DEFAULT') : determines which solver to use, can be either
a (new-style) string or an (old-style) numerical alg code
'DEFAULT' : (or 0) automatic, first available of CPLEX,
Gurobi, MOSEK, Opt Tbx (MILPs only), GLPK (MILPs only)
'CPLEX' : (or 500) CPLEX
'GLPK' : GLPK, (MILP problems only, i.e. empty H matrix)
'GUROBI' : (or 700) Gurobi
'MOSEK' : (or 600) MOSEK
'OT' : (or 300) Optimization Toolbox, INTLINPROG
(MILP problems only, i.e. empty H matrix)
verbose (0) - controls level of progress output displayed
0 = no progress output
1 = some progress output
2 = verbose progress output
skip_prices (0) - flag that specifies whether or not to
skip the price computation stage, in which the problem
is re-solved for only the continuous variables, with all
others being constrained to their solved values
cplex_opt - options struct for CPLEX
glpk_opt - options struct for GLPK
grb_opt - options struct for GBUROBI_MEX
intlinprog_opt - options struct for INTLINPROG
linprog_opt - options struct for LINPROG
mosek_opt - options struct for MOSEK
PROBLEM : The inputs can alternatively be supplied in a single
PROBLEM struct with fields corresponding to the input arguments
described above: H, c, A, l, u, xmin, xmax, x0, vtype, opt
Outputs:
X : solution vector
F : final objective function value
EXITFLAG : exit flag
1 = converged
0 or negative values = algorithm specific failure codes
OUTPUT : output struct with the following fields:
alg - algorithm code of solver used
(others) - algorithm specific fields
LAMBDA : struct containing the Langrange and Kuhn-Tucker
multipliers on the constraints, with fields:
mu_l - lower (left-hand) limit on linear constraints
mu_u - upper (right-hand) limit on linear constraints
lower - lower bound on optimization variables
upper - upper bound on optimization variables
Note the calling syntax is almost identical to that of QUADPROG
from MathWorks' Optimization Toolbox. The main difference is that
the linear constraints are specified with A, L, U instead of
A, B, Aeq, Beq.
Calling syntax options:
[x, f, exitflag, output, lambda] = ...
miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
x = miqps_matpower(H, c, A, l, u)
x = miqps_matpower(H, c, A, l, u, xmin, xmax)
x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0)
x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype)
x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
x = miqps_matpower(problem), where problem is a struct with fields:
H, c, A, l, u, xmin, xmax, x0, vtype, opt
all fields except 'c', 'A' and 'l' or 'u' are optional
x = miqps_matpower(...)
[x, f] = miqps_matpower(...)
[x, f, exitflag] = miqps_matpower(...)
[x, f, exitflag, output] = miqps_matpower(...)
[x, f, exitflag, output, lambda] = miqps_matpower(...)
Example: (problem from from http://www.jmu.edu/docs/sasdoc/sashtml/iml/chap8/sect12.htm)
H = [ 1003.1 4.3 6.3 5.9;
4.3 2.2 2.1 3.9;
6.3 2.1 3.5 4.8;
5.9 3.9 4.8 10 ];
c = zeros(4,1);
A = [ 1 1 1 1;
0.17 0.11 0.10 0.18 ];
l = [1; 0.10];
u = [1; Inf];
xmin = zeros(4,1);
x0 = [1; 0; 0; 1];
opt = struct('verbose', 2);
[x, f, s, out, lambda] = miqps_matpower(H, c, A, l, u, xmin, [], x0, vtype, opt);0001 function [x, f, eflag, output, lambda] = miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0002 %MIQPS_MATPOWER Mixed Integer Quadratic Program Solver for MATPOWER. 0003 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... 0004 % MIQPS_MATPOWER(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT) 0005 % A common wrapper function for various QP solvers. 0006 % Solves the following QP (quadratic programming) problem: 0007 % 0008 % min 1/2 X'*H*X + C'*X 0009 % X 0010 % 0011 % subject to 0012 % 0013 % L <= A*X <= U (linear constraints) 0014 % XMIN <= X <= XMAX (variable bounds) 0015 % 0016 % Inputs (all optional except H, C, A and L): 0017 % H : matrix (possibly sparse) of quadratic cost coefficients 0018 % C : vector of linear cost coefficients 0019 % A, L, U : define the optional linear constraints. Default 0020 % values for the elements of L and U are -Inf and Inf, 0021 % respectively. 0022 % XMIN, XMAX : optional lower and upper bounds on the 0023 % X variables, defaults are -Inf and Inf, respectively. 0024 % X0 : optional starting value of optimization vector X 0025 % VTYPE : character string of length NX (number of elements in X), 0026 % or 1 (value applies to all variables in x), 0027 % allowed values are 'C' (continuous), 'B' (binary), 0028 % 'I' (integer), 'S' (semi-continuous), or 'N' (semi-integer). 0029 % (MOSEK, GLPK, OT allow only 'C', 'B', or 'I') 0030 % OPT : optional options structure with the following fields, 0031 % all of which are also optional (default values shown in 0032 % parentheses) 0033 % alg ('DEFAULT') : determines which solver to use, can be either 0034 % a (new-style) string or an (old-style) numerical alg code 0035 % 'DEFAULT' : (or 0) automatic, first available of CPLEX, 0036 % Gurobi, MOSEK, Opt Tbx (MILPs only), GLPK (MILPs only) 0037 % 'CPLEX' : (or 500) CPLEX 0038 % 'GLPK' : GLPK, (MILP problems only, i.e. empty H matrix) 0039 % 'GUROBI' : (or 700) Gurobi 0040 % 'MOSEK' : (or 600) MOSEK 0041 % 'OT' : (or 300) Optimization Toolbox, INTLINPROG 0042 % (MILP problems only, i.e. empty H matrix) 0043 % verbose (0) - controls level of progress output displayed 0044 % 0 = no progress output 0045 % 1 = some progress output 0046 % 2 = verbose progress output 0047 % skip_prices (0) - flag that specifies whether or not to 0048 % skip the price computation stage, in which the problem 0049 % is re-solved for only the continuous variables, with all 0050 % others being constrained to their solved values 0051 % cplex_opt - options struct for CPLEX 0052 % glpk_opt - options struct for GLPK 0053 % grb_opt - options struct for GBUROBI_MEX 0054 % intlinprog_opt - options struct for INTLINPROG 0055 % linprog_opt - options struct for LINPROG 0056 % mosek_opt - options struct for MOSEK 0057 % PROBLEM : The inputs can alternatively be supplied in a single 0058 % PROBLEM struct with fields corresponding to the input arguments 0059 % described above: H, c, A, l, u, xmin, xmax, x0, vtype, opt 0060 % 0061 % Outputs: 0062 % X : solution vector 0063 % F : final objective function value 0064 % EXITFLAG : exit flag 0065 % 1 = converged 0066 % 0 or negative values = algorithm specific failure codes 0067 % OUTPUT : output struct with the following fields: 0068 % alg - algorithm code of solver used 0069 % (others) - algorithm specific fields 0070 % LAMBDA : struct containing the Langrange and Kuhn-Tucker 0071 % multipliers on the constraints, with fields: 0072 % mu_l - lower (left-hand) limit on linear constraints 0073 % mu_u - upper (right-hand) limit on linear constraints 0074 % lower - lower bound on optimization variables 0075 % upper - upper bound on optimization variables 0076 % 0077 % Note the calling syntax is almost identical to that of QUADPROG 0078 % from MathWorks' Optimization Toolbox. The main difference is that 0079 % the linear constraints are specified with A, L, U instead of 0080 % A, B, Aeq, Beq. 0081 % 0082 % Calling syntax options: 0083 % [x, f, exitflag, output, lambda] = ... 0084 % miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0085 % 0086 % x = miqps_matpower(H, c, A, l, u) 0087 % x = miqps_matpower(H, c, A, l, u, xmin, xmax) 0088 % x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0) 0089 % x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype) 0090 % x = miqps_matpower(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0091 % x = miqps_matpower(problem), where problem is a struct with fields: 0092 % H, c, A, l, u, xmin, xmax, x0, vtype, opt 0093 % all fields except 'c', 'A' and 'l' or 'u' are optional 0094 % x = miqps_matpower(...) 0095 % [x, f] = miqps_matpower(...) 0096 % [x, f, exitflag] = miqps_matpower(...) 0097 % [x, f, exitflag, output] = miqps_matpower(...) 0098 % [x, f, exitflag, output, lambda] = miqps_matpower(...) 0099 % 0100 % Example: (problem from from http://www.jmu.edu/docs/sasdoc/sashtml/iml/chap8/sect12.htm) 0101 % H = [ 1003.1 4.3 6.3 5.9; 0102 % 4.3 2.2 2.1 3.9; 0103 % 6.3 2.1 3.5 4.8; 0104 % 5.9 3.9 4.8 10 ]; 0105 % c = zeros(4,1); 0106 % A = [ 1 1 1 1; 0107 % 0.17 0.11 0.10 0.18 ]; 0108 % l = [1; 0.10]; 0109 % u = [1; Inf]; 0110 % xmin = zeros(4,1); 0111 % x0 = [1; 0; 0; 1]; 0112 % opt = struct('verbose', 2); 0113 % [x, f, s, out, lambda] = miqps_matpower(H, c, A, l, u, xmin, [], x0, vtype, opt); 0114 0115 % MATPOWER 0116 % Copyright (c) 2010-2015 by Power System Engineering Research Center (PSERC) 0117 % by Ray Zimmerman, PSERC Cornell 0118 % 0119 % $Id: miqps_matpower.m 2644 2015-03-11 19:34:22Z ray $ 0120 % 0121 % This file is part of MATPOWER. 0122 % Covered by the 3-clause BSD License (see LICENSE file for details). 0123 % See http://www.pserc.cornell.edu/matpower/ for more info. 0124 0125 %%----- input argument handling ----- 0126 %% gather inputs 0127 if nargin == 1 && isstruct(H) %% problem struct 0128 p = H; 0129 if isfield(p, 'opt'), opt = p.opt; else, opt = []; end 0130 if isfield(p, 'vtype'), vtype = p.vtype;else, vtype = []; end 0131 if isfield(p, 'x0'), x0 = p.x0; else, x0 = []; end 0132 if isfield(p, 'xmax'), xmax = p.xmax; else, xmax = []; end 0133 if isfield(p, 'xmin'), xmin = p.xmin; else, xmin = []; end 0134 if isfield(p, 'u'), u = p.u; else, u = []; end 0135 if isfield(p, 'l'), l = p.l; else, l = []; end 0136 if isfield(p, 'A'), A = p.A; else, A = []; end 0137 if isfield(p, 'c'), c = p.c; else, c = []; end 0138 if isfield(p, 'H'), H = p.H; else, H = []; end 0139 else %% individual args 0140 if nargin < 10 0141 opt = []; 0142 if nargin < 9 0143 vtype = []; 0144 if nargin < 8 0145 x0 = []; 0146 if nargin < 7 0147 xmax = []; 0148 if nargin < 6 0149 xmin = []; 0150 end 0151 end 0152 end 0153 end 0154 end 0155 end 0156 0157 %% default options 0158 if ~isempty(opt) && isfield(opt, 'alg') && ~isempty(opt.alg) 0159 alg = opt.alg; 0160 %% convert integer codes to string values 0161 if ~ischar(alg) 0162 switch alg 0163 case 0 0164 alg = 'DEFAULT'; 0165 case 300 0166 alg = 'OT'; 0167 case 500 0168 alg = 'CPLEX'; 0169 case 600 0170 alg = 'MOSEK'; 0171 case 700 0172 alg = 'GUROBI'; 0173 otherwise 0174 error('miqps_matpower: %d is not a valid algorithm code', alg); 0175 end 0176 end 0177 else 0178 alg = 'DEFAULT'; 0179 end 0180 if ~isempty(opt) && isfield(opt, 'verbose') && ~isempty(opt.verbose) 0181 verbose = opt.verbose; 0182 else 0183 verbose = 0; 0184 end 0185 if strcmp(alg, 'DEFAULT') 0186 if have_fcn('cplex') %% use CPLEX by default, if available 0187 alg = 'CPLEX'; 0188 elseif have_fcn('gurobi') %% if not, then Gurobi, if available 0189 alg = 'GUROBI'; 0190 elseif have_fcn('mosek') %% if not, then MOSEK, if available 0191 alg = 'MOSEK'; 0192 elseif isempty(H) || ~any(any(H)) %% if not, and linear objective 0193 if have_fcn('intlinprog') %% then Optimization Tbx, if available 0194 alg = 'OT'; 0195 elseif have_fcn('glpk') %% if not, and then GLPK, if available 0196 alg = 'GLPK'; 0197 end 0198 else 0199 error('miqps_matpower: no solvers available - requires CPLEX, Gurobi, MOSEK, INTLINPROG or GLPK'); 0200 end 0201 end 0202 0203 %%----- call the appropriate solver ----- 0204 switch alg 0205 case 'CPLEX' 0206 [x, f, eflag, output, lambda] = ... 0207 miqps_cplex(H, c, A, l, u, xmin, xmax, x0, vtype, opt); 0208 case 'GLPK' 0209 [x, f, eflag, output, lambda] = ... 0210 miqps_glpk(H, c, A, l, u, xmin, xmax, x0, vtype, opt); 0211 case 'GUROBI' 0212 [x, f, eflag, output, lambda] = ... 0213 miqps_gurobi(H, c, A, l, u, xmin, xmax, x0, vtype, opt); 0214 case 'MOSEK' 0215 [x, f, eflag, output, lambda] = ... 0216 miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt); 0217 case 'OT' 0218 [x, f, eflag, output, lambda] = ... 0219 miqps_ot(H, c, A, l, u, xmin, xmax, x0, vtype, opt); 0220 otherwise 0221 error('miqps_matpower: %d is not a valid algorithm code', alg); 0222 end 0223 if ~isfield(output, 'alg') || isempty(output.alg) 0224 output.alg = alg; 0225 end