Home > matpower7.1 > mp-opt-model > lib > qps_master.m

qps_master

PURPOSE ^

QPS_MASTER Quadratic Program Solver wrapper function.

SYNOPSIS ^

function [x, f, eflag, output, lambda] = qps_master(H, c, A, l, u, xmin, xmax, x0, opt)

DESCRIPTION ^

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);

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

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

Generated on Fri 09-Oct-2020 11:21:31 by m2html © 2005