MIQPS_MASTER Mixed Integer Quadratic Program Solver wrapper function. [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... MIQPS_MASTER(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT) [X, F, EXITFLAG, OUTPUT, LAMBDA] = MIQPS_MASTER(PROBLEM) A common wrapper function for various QP solvers. Solves the following MIQP (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) X(i) is integer, for i in I (integer variable constraints) X(b) is binary, for b in B (binary variable constraints) 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 string (new-style) or a numerical alg code (old-style) 'DEFAULT' : (or 0) automatic, first available of Gurobi, CPLEX, 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 price_stage_warn_tol (1e-7) - tolerance on the objective fcn value and primal variable relative match required to avoid mis-match warning message cplex_opt - options struct for CPLEX glpk_opt - options struct for GLPK grb_opt - options struct for GUROBI 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 = solver 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_master(H, c, A, l, u, xmin, xmax, x0, vtype, opt) x = miqps_master(H, c, A, l, u) x = miqps_master(H, c, A, l, u, xmin, xmax) x = miqps_master(H, c, A, l, u, xmin, xmax, x0) x = miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype) x = miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype, opt) x = miqps_master(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_master(...) [x, f] = miqps_master(...) [x, f, exitflag] = miqps_master(...) [x, f, exitflag, output] = miqps_master(...) [x, f, exitflag, output, lambda] = miqps_master(...) Example: (problem from from %% from MOSEK 6.0 Guided Tour, section 7.13.1 https://docs.mosek.com/6.0/toolbox/node009.html) c = [-2; -3]; A = sparse([195 273; 4 40]); u = [1365; 140]; xmax = [4; Inf]; vtype = 'I'; opt = struct('verbose', 2); p = struct('c', c, 'A', A, 'u', u, 'xmax', xmax, 'vtype', vtype, 'opt', opt); [x, f, s, out, lam] = miqps_master(p);
0001 function [x, f, eflag, output, lambda] = miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0002 %MIQPS_MASTER Mixed Integer Quadratic Program Solver wrapper function. 0003 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... 0004 % MIQPS_MASTER(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT) 0005 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = MIQPS_MASTER(PROBLEM) 0006 % A common wrapper function for various QP solvers. 0007 % Solves the following MIQP (quadratic programming) problem: 0008 % 0009 % min 1/2 X'*H*X + C'*X 0010 % X 0011 % 0012 % subject to 0013 % 0014 % L <= A*X <= U (linear constraints) 0015 % XMIN <= X <= XMAX (variable bounds) 0016 % X(i) is integer, for i in I (integer variable constraints) 0017 % X(b) is binary, for b in B (binary variable constraints) 0018 % 0019 % Inputs (all optional except H, C, A and L): 0020 % H : matrix (possibly sparse) of quadratic cost coefficients 0021 % C : vector of linear cost coefficients 0022 % A, L, U : define the optional linear constraints. Default 0023 % values for the elements of L and U are -Inf and Inf, 0024 % respectively. 0025 % XMIN, XMAX : optional lower and upper bounds on the 0026 % X variables, defaults are -Inf and Inf, respectively. 0027 % X0 : optional starting value of optimization vector X 0028 % VTYPE : character string of length NX (number of elements in X), 0029 % or 1 (value applies to all variables in x), 0030 % allowed values are 'C' (continuous), 'B' (binary), 0031 % 'I' (integer), 'S' (semi-continuous), or 'N' (semi-integer). 0032 % (MOSEK, GLPK, OT allow only 'C', 'B', or 'I') 0033 % OPT : optional options structure with the following fields, 0034 % all of which are also optional (default values shown in 0035 % parentheses) 0036 % alg ('DEFAULT') : determines which solver to use, can be either 0037 % a string (new-style) or a numerical alg code (old-style) 0038 % 'DEFAULT' : (or 0) automatic, first available of Gurobi, 0039 % CPLEX, MOSEK, Opt Tbx (MILPs only), GLPK (MILPs only) 0040 % 'CPLEX' : (or 500) CPLEX 0041 % 'GLPK' : GLPK, (MILP problems only, i.e. empty H matrix) 0042 % 'GUROBI' : (or 700) Gurobi 0043 % 'MOSEK' : (or 600) MOSEK 0044 % 'OT' : (or 300) Optimization Toolbox, INTLINPROG 0045 % (MILP problems only, i.e. empty H matrix) 0046 % verbose (0) - controls level of progress output displayed 0047 % 0 = no progress output 0048 % 1 = some progress output 0049 % 2 = verbose progress output 0050 % skip_prices (0) - flag that specifies whether or not to 0051 % skip the price computation stage, in which the problem 0052 % is re-solved for only the continuous variables, with all 0053 % others being constrained to their solved values 0054 % price_stage_warn_tol (1e-7) - tolerance on the objective fcn 0055 % value and primal variable relative match required to avoid 0056 % mis-match warning message 0057 % cplex_opt - options struct for CPLEX 0058 % glpk_opt - options struct for GLPK 0059 % grb_opt - options struct for GUROBI 0060 % intlinprog_opt - options struct for INTLINPROG 0061 % linprog_opt - options struct for LINPROG 0062 % mosek_opt - options struct for MOSEK 0063 % PROBLEM : The inputs can alternatively be supplied in a single 0064 % PROBLEM struct with fields corresponding to the input arguments 0065 % described above: H, c, A, l, u, xmin, xmax, x0, vtype, opt 0066 % 0067 % Outputs: 0068 % X : solution vector 0069 % F : final objective function value 0070 % EXITFLAG : exit flag 0071 % 1 = converged 0072 % 0 or negative values = solver specific failure codes 0073 % OUTPUT : output struct with the following fields: 0074 % alg - algorithm code of solver used 0075 % (others) - algorithm specific fields 0076 % LAMBDA : struct containing the Langrange and Kuhn-Tucker 0077 % multipliers on the constraints, with fields: 0078 % mu_l - lower (left-hand) limit on linear constraints 0079 % mu_u - upper (right-hand) limit on linear constraints 0080 % lower - lower bound on optimization variables 0081 % upper - upper bound on optimization variables 0082 % 0083 % Note the calling syntax is almost identical to that of QUADPROG 0084 % from MathWorks' Optimization Toolbox. The main difference is that 0085 % the linear constraints are specified with A, L, U instead of 0086 % A, B, Aeq, Beq. 0087 % 0088 % Calling syntax options: 0089 % [x, f, exitflag, output, lambda] = ... 0090 % miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0091 % 0092 % x = miqps_master(H, c, A, l, u) 0093 % x = miqps_master(H, c, A, l, u, xmin, xmax) 0094 % x = miqps_master(H, c, A, l, u, xmin, xmax, x0) 0095 % x = miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype) 0096 % x = miqps_master(H, c, A, l, u, xmin, xmax, x0, vtype, opt) 0097 % x = miqps_master(problem), where problem is a struct with fields: 0098 % H, c, A, l, u, xmin, xmax, x0, vtype, opt 0099 % all fields except 'c', 'A' and 'l' or 'u' are optional 0100 % x = miqps_master(...) 0101 % [x, f] = miqps_master(...) 0102 % [x, f, exitflag] = miqps_master(...) 0103 % [x, f, exitflag, output] = miqps_master(...) 0104 % [x, f, exitflag, output, lambda] = miqps_master(...) 0105 % 0106 % Example: (problem from from %% from MOSEK 6.0 Guided Tour, section 7.13.1 0107 % https://docs.mosek.com/6.0/toolbox/node009.html) 0108 % c = [-2; -3]; 0109 % A = sparse([195 273; 4 40]); 0110 % u = [1365; 140]; 0111 % xmax = [4; Inf]; 0112 % vtype = 'I'; 0113 % opt = struct('verbose', 2); 0114 % p = struct('c', c, 'A', A, 'u', u, 'xmax', xmax, 'vtype', vtype, 'opt', opt); 0115 % [x, f, s, out, lam] = miqps_master(p); 0116 0117 % MP-Opt-Model 0118 % Copyright (c) 2010-2020, Power Systems Engineering Research Center (PSERC) 0119 % by Ray Zimmerman, PSERC Cornell 0120 % 0121 % This file is part of MP-Opt-Model. 0122 % Covered by the 3-clause BSD License (see LICENSE file for details). 0123 % See https://github.com/MATPOWER/mp-opt-model 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_master: %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_feature('gurobi') %% use Gurobi by default, if available 0187 alg = 'GUROBI'; 0188 elseif have_feature('cplex') %% if not, then CPLEX, if available 0189 alg = 'CPLEX'; 0190 elseif have_feature('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_feature('intlinprog') %% then Optimization Tbx, if available 0194 alg = 'OT'; 0195 elseif have_feature('glpk') %% if not, and then GLPK, if available 0196 alg = 'GLPK'; 0197 end 0198 else 0199 error('miqps_master: 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_master: ''%s'' is not a valid algorithm code', alg); 0222 end 0223 if ~isfield(output, 'alg') || isempty(output.alg) 0224 output.alg = alg; 0225 end