QPS_MASTER Quadratic Program Solver wrapper function. [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... QPS_MASTER(H, C, A, L, U, XMIN, XMAX, X0, OPT) [X, F, EXITFLAG, OUTPUT, LAMBDA] = QPS_MASTER(PROBLEM) 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 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 (if MATLAB), GLPK (LPs only), BPMPD, MIPS 'MIPS' : (or 200) MIPS, MATPOWER Interior Point Solver pure MATLAB implementation of a primal-dual interior point method, if mips_opt.step_control = 1 (or alg=250) it uses MIPS-sc, a step controlled variant of MIPS 'BPMPD' : (or 100) BPMPD_MEX 'CLP' : CLP 'CPLEX' : (or 500) CPLEX 'GLPK' : GLPK, (LP problems only, i.e. empty H matrix) 'GUROBI' : (or 700) Gurobi 'IPOPT' : (or 400) IPOPT, requires MEX interface to IPOPT solver, https://github.com/coin-or/Ipopt 'MOSEK' : (or 600) MOSEK 'OSQP' : OSQP, https://osqp.org 'OT' : (or 300) Optimization Toolbox, QUADPROG or LINPROG verbose (0) - controls level of progress output displayed 0 = no progress output 1 = some progress output 2 = verbose progress output bp_opt - options vector for BP clp_opt - options vector for CLP cplex_opt - options struct for CPLEX glpk_opt - options struct for GLPK grb_opt - options struct for GUROBI ipopt_opt - options struct for IPOPT linprog_opt - options struct for LINPROG mips_opt - options struct for QPS_MIPS mosek_opt - options struct for MOSEK osqp_opt - options struct for OSQP quadprog_opt - options struct for QUADPROG 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, 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] = ... qps_master(H, c, A, l, u, xmin, xmax, x0, opt) x = qps_master(H, c, A, l, u) x = qps_master(H, c, A, l, u, xmin, xmax) x = qps_master(H, c, A, l, u, xmin, xmax, x0) x = qps_master(H, c, A, l, u, xmin, xmax, x0, opt) x = qps_master(problem), where problem is a struct with fields: H, c, A, l, u, xmin, xmax, x0, opt all fields except 'c', 'A' and 'l' or 'u' are optional x = qps_master(...) [x, f] = qps_master(...) [x, f, exitflag] = qps_master(...) [x, f, exitflag, output] = qps_master(...) [x, f, exitflag, output, lambda] = qps_master(...) Example: (problem from from https://v8doc.sas.com/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] = qps_master(H, c, A, l, u, xmin, [], x0, opt);
0001 function [x, f, eflag, output, lambda] = qps_master(H, c, A, l, u, xmin, xmax, x0, opt) 0002 %QPS_MASTER Quadratic Program Solver wrapper function. 0003 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = ... 0004 % QPS_MASTER(H, C, A, L, U, XMIN, XMAX, X0, OPT) 0005 % [X, F, EXITFLAG, OUTPUT, LAMBDA] = QPS_MASTER(PROBLEM) 0006 % A common wrapper function for various QP solvers. 0007 % Solves the following QP (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 % 0017 % Inputs (all optional except H, C, A and L): 0018 % H : matrix (possibly sparse) of quadratic cost coefficients 0019 % C : vector of linear cost coefficients 0020 % A, L, U : define the optional linear constraints. Default 0021 % values for the elements of L and U are -Inf and Inf, 0022 % respectively. 0023 % XMIN, XMAX : optional lower and upper bounds on the 0024 % X variables, defaults are -Inf and Inf, respectively. 0025 % X0 : optional starting value of optimization vector X 0026 % OPT : optional options structure with the following fields, 0027 % all of which are also optional (default values shown in 0028 % parentheses) 0029 % alg ('DEFAULT') : determines which solver to use, can be either 0030 % a string (new-style) or a numerical alg code (old-style) 0031 % 'DEFAULT' : (or 0) automatic, first available of Gurobi, 0032 % CPLEX, MOSEK, Opt Tbx (if MATLAB), GLPK (LPs only), 0033 % BPMPD, MIPS 0034 % 'MIPS' : (or 200) MIPS, MATPOWER Interior Point Solver 0035 % pure MATLAB implementation of a primal-dual 0036 % interior point method, if mips_opt.step_control = 1 0037 % (or alg=250) it uses MIPS-sc, a step controlled 0038 % variant of MIPS 0039 % 'BPMPD' : (or 100) BPMPD_MEX 0040 % 'CLP' : CLP 0041 % 'CPLEX' : (or 500) CPLEX 0042 % 'GLPK' : GLPK, (LP problems only, i.e. empty H matrix) 0043 % 'GUROBI' : (or 700) Gurobi 0044 % 'IPOPT' : (or 400) IPOPT, requires MEX interface to IPOPT 0045 % solver, https://github.com/coin-or/Ipopt 0046 % 'MOSEK' : (or 600) MOSEK 0047 % 'OSQP' : OSQP, https://osqp.org 0048 % 'OT' : (or 300) Optimization Toolbox, QUADPROG or LINPROG 0049 % verbose (0) - controls level of progress output displayed 0050 % 0 = no progress output 0051 % 1 = some progress output 0052 % 2 = verbose progress output 0053 % bp_opt - options vector for BP 0054 % clp_opt - options vector for CLP 0055 % cplex_opt - options struct for CPLEX 0056 % glpk_opt - options struct for GLPK 0057 % grb_opt - options struct for GUROBI 0058 % ipopt_opt - options struct for IPOPT 0059 % linprog_opt - options struct for LINPROG 0060 % mips_opt - options struct for QPS_MIPS 0061 % mosek_opt - options struct for MOSEK 0062 % osqp_opt - options struct for OSQP 0063 % quadprog_opt - options struct for QUADPROG 0064 % PROBLEM : The inputs can alternatively be supplied in a single 0065 % PROBLEM struct with fields corresponding to the input arguments 0066 % described above: H, c, A, l, u, xmin, xmax, x0, opt 0067 % 0068 % Outputs: 0069 % X : solution vector 0070 % F : final objective function value 0071 % EXITFLAG : exit flag 0072 % 1 = converged 0073 % 0 or negative values = solver specific failure codes 0074 % OUTPUT : output struct with the following fields: 0075 % alg - algorithm code of solver used 0076 % (others) - algorithm specific fields 0077 % LAMBDA : struct containing the Langrange and Kuhn-Tucker 0078 % multipliers on the constraints, with fields: 0079 % mu_l - lower (left-hand) limit on linear constraints 0080 % mu_u - upper (right-hand) limit on linear constraints 0081 % lower - lower bound on optimization variables 0082 % upper - upper bound on optimization variables 0083 % 0084 % Note the calling syntax is almost identical to that of QUADPROG 0085 % from MathWorks' Optimization Toolbox. The main difference is that 0086 % the linear constraints are specified with A, L, U instead of 0087 % A, B, Aeq, Beq. 0088 % 0089 % Calling syntax options: 0090 % [x, f, exitflag, output, lambda] = ... 0091 % qps_master(H, c, A, l, u, xmin, xmax, x0, opt) 0092 % 0093 % x = qps_master(H, c, A, l, u) 0094 % x = qps_master(H, c, A, l, u, xmin, xmax) 0095 % x = qps_master(H, c, A, l, u, xmin, xmax, x0) 0096 % x = qps_master(H, c, A, l, u, xmin, xmax, x0, opt) 0097 % x = qps_master(problem), where problem is a struct with fields: 0098 % H, c, A, l, u, xmin, xmax, x0, opt 0099 % all fields except 'c', 'A' and 'l' or 'u' are optional 0100 % x = qps_master(...) 0101 % [x, f] = qps_master(...) 0102 % [x, f, exitflag] = qps_master(...) 0103 % [x, f, exitflag, output] = qps_master(...) 0104 % [x, f, exitflag, output, lambda] = qps_master(...) 0105 % 0106 % Example: (problem from from https://v8doc.sas.com/sashtml/iml/chap8/sect12.htm) 0107 % H = [ 1003.1 4.3 6.3 5.9; 0108 % 4.3 2.2 2.1 3.9; 0109 % 6.3 2.1 3.5 4.8; 0110 % 5.9 3.9 4.8 10 ]; 0111 % c = zeros(4,1); 0112 % A = [ 1 1 1 1; 0113 % 0.17 0.11 0.10 0.18 ]; 0114 % l = [1; 0.10]; 0115 % u = [1; Inf]; 0116 % xmin = zeros(4,1); 0117 % x0 = [1; 0; 0; 1]; 0118 % opt = struct('verbose', 2); 0119 % [x, f, s, out, lambda] = qps_master(H, c, A, l, u, xmin, [], x0, opt); 0120 0121 % MP-Opt-Model 0122 % Copyright (c) 2010-2020, Power Systems Engineering Research Center (PSERC) 0123 % by Ray Zimmerman, PSERC Cornell 0124 % 0125 % This file is part of MP-Opt-Model. 0126 % Covered by the 3-clause BSD License (see LICENSE file for details). 0127 % See https://github.com/MATPOWER/mp-opt-model for more info. 0128 0129 %%----- input argument handling ----- 0130 %% gather inputs 0131 if nargin == 1 && isstruct(H) %% problem struct 0132 p = H; 0133 if isfield(p, 'opt'), opt = p.opt; else, opt = []; end 0134 if isfield(p, 'x0'), x0 = p.x0; else, x0 = []; end 0135 if isfield(p, 'xmax'), xmax = p.xmax; else, xmax = []; end 0136 if isfield(p, 'xmin'), xmin = p.xmin; else, xmin = []; end 0137 if isfield(p, 'u'), u = p.u; else, u = []; end 0138 if isfield(p, 'l'), l = p.l; else, l = []; end 0139 if isfield(p, 'A'), A = p.A; else, A = []; end 0140 if isfield(p, 'c'), c = p.c; else, c = []; end 0141 if isfield(p, 'H'), H = p.H; else, H = []; end 0142 else %% individual args 0143 if nargin < 9 0144 opt = []; 0145 if nargin < 8 0146 x0 = []; 0147 if nargin < 7 0148 xmax = []; 0149 if nargin < 6 0150 xmin = []; 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 100 0166 alg = 'BPMPD'; 0167 case 200 0168 alg = 'MIPS'; 0169 opt.mips_opt.step_control = 0; 0170 case 250 0171 alg = 'MIPS'; 0172 opt.mips_opt.step_control = 1; 0173 case 300 0174 alg = 'OT'; 0175 case 400 0176 alg = 'IPOPT'; 0177 case 500 0178 alg = 'CPLEX'; 0179 case 600 0180 alg = 'MOSEK'; 0181 case 700 0182 alg = 'GUROBI'; 0183 otherwise 0184 error('qps_master: %d is not a valid algorithm code', alg); 0185 end 0186 end 0187 else 0188 alg = 'DEFAULT'; 0189 end 0190 if ~isempty(opt) && isfield(opt, 'verbose') && ~isempty(opt.verbose) 0191 verbose = opt.verbose; 0192 else 0193 verbose = 0; 0194 end 0195 if strcmp(alg, 'DEFAULT') 0196 if have_feature('gurobi') %% use Gurobi by default, if available 0197 alg = 'GUROBI'; 0198 elseif have_feature('cplex') %% if not, then CPLEX, if available 0199 alg = 'CPLEX'; 0200 elseif have_feature('mosek') %% if not, then MOSEK, if available 0201 alg = 'MOSEK'; 0202 elseif have_feature('quadprog') && have_feature('matlab') %% if not, then Opt Tbx, if available in MATLAB 0203 alg = 'OT'; 0204 elseif (isempty(H) || ~any(any(H))) && have_feature('glpk') %% if not, and 0205 alg = 'GLPK'; %% prob is LP (not QP), then GLPK, if available 0206 elseif have_feature('bpmpd') %% if not, then BPMPD_MEX, if available 0207 alg = 'BPMPD'; 0208 else %% otherwise MIPS 0209 alg = 'MIPS'; 0210 end 0211 end 0212 0213 %%----- call the appropriate solver ----- 0214 switch alg 0215 case 'BPMPD' %% use BPMPD_MEX 0216 [x, f, eflag, output, lambda] = ... 0217 qps_bpmpd(H, c, A, l, u, xmin, xmax, x0, opt); 0218 0219 if eflag == -99 0220 if verbose 0221 fprintf(' Retrying with QPS_MIPS solver ...\n\n'); 0222 end 0223 %% save (incorrect) solution from BPMPD 0224 bpmpd = struct('x', x, 'f', f, 'eflag', eflag, ... 0225 'output', output, 'lambda', lambda); 0226 opt.alg = 'MIPS'; 0227 [x, f, eflag, output, lambda] = ... 0228 qps_master(H, c, A, l, u, xmin, xmax, x0, opt); 0229 output.bpmpd = bpmpd; 0230 end 0231 case 'CLP' 0232 [x, f, eflag, output, lambda] = ... 0233 qps_clp(H, c, A, l, u, xmin, xmax, x0, opt); 0234 case 'CPLEX' 0235 [x, f, eflag, output, lambda] = ... 0236 qps_cplex(H, c, A, l, u, xmin, xmax, x0, opt); 0237 case 'GLPK' 0238 [x, f, eflag, output, lambda] = ... 0239 qps_glpk(H, c, A, l, u, xmin, xmax, x0, opt); 0240 case 'GUROBI' 0241 [x, f, eflag, output, lambda] = ... 0242 qps_gurobi(H, c, A, l, u, xmin, xmax, x0, opt); 0243 case 'IPOPT' 0244 [x, f, eflag, output, lambda] = ... 0245 qps_ipopt(H, c, A, l, u, xmin, xmax, x0, opt); 0246 case 'MIPS' 0247 %% set up options 0248 if ~isempty(opt) && isfield(opt, 'mips_opt') && ~isempty(opt.mips_opt) 0249 mips_opt = opt.mips_opt; 0250 else 0251 mips_opt = []; 0252 end 0253 mips_opt.verbose = verbose; 0254 0255 %% call solver 0256 [x, f, eflag, output, lambda] = ... 0257 qps_mips(H, c, A, l, u, xmin, xmax, x0, mips_opt); 0258 case 'MOSEK' 0259 [x, f, eflag, output, lambda] = ... 0260 qps_mosek(H, c, A, l, u, xmin, xmax, x0, opt); 0261 case 'OSQP' 0262 [x, f, eflag, output, lambda] = ... 0263 qps_osqp(H, c, A, l, u, xmin, xmax, x0, opt); 0264 case 'OT' %% use QUADPROG or LINPROG from Opt Tbx ver 2.x+ 0265 [x, f, eflag, output, lambda] = ... 0266 qps_ot(H, c, A, l, u, xmin, xmax, x0, opt); 0267 otherwise 0268 error('qps_master: ''%s'' is not a valid algorithm code', alg); 0269 end 0270 if ~isfield(output, 'alg') || isempty(output.alg) 0271 output.alg = alg; 0272 end