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

miqps_mosek

PURPOSE ^

MIQPS_MOSEK Mixed Integer Quadratic Program Solver based on MOSEK.

SYNOPSIS ^

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

DESCRIPTION ^

MIQPS_MOSEK  Mixed Integer Quadratic Program Solver based on MOSEK.
   [X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
       MIQPS_MOSEK(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT)
   [X, F, EXITFLAG, OUTPUT, LAMBDA] = MIQPS_MOSEK(PROBLEM)
   A wrapper function providing a standardized interface for using
   MOSEKOPT to solve 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).
       OPT : optional options structure with the following fields,
           all of which are also optional (default values shown in
           parentheses)
           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
           mosek_opt - options struct for MOSEK, value in verbose
                   overrides these options
       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 = success
             0 = terminated at maximum number of iterations
            -1 = primal or dual infeasible
           < 0 = the negative of the MOSEK return code
       OUTPUT : output struct with the following fields:
           r - MOSEK return code
           res - MOSEK result struct
       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_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt)

       x = miqps_mosek(H, c, A, l, u)
       x = miqps_mosek(H, c, A, l, u, xmin, xmax)
       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0)
       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype)
       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
       x = miqps_mosek(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_mosek(...)
       [x, f] = miqps_mosek(...)
       [x, f, exitflag] = miqps_mosek(...)
       [x, f, exitflag, output] = miqps_mosek(...)
       [x, f, exitflag, output, lambda] = miqps_mosek(...)

   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] = miqps_mosek(H, c, A, l, u, xmin, [], x0, vtype, opt);

   See also MIQPS_MASTER, MOSEKOPT.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [x, f, eflag, output, lambda] = miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
0002 %MIQPS_MOSEK  Mixed Integer Quadratic Program Solver based on MOSEK.
0003 %   [X, F, EXITFLAG, OUTPUT, LAMBDA] = ...
0004 %       MIQPS_MOSEK(H, C, A, L, U, XMIN, XMAX, X0, VTYPE, OPT)
0005 %   [X, F, EXITFLAG, OUTPUT, LAMBDA] = MIQPS_MOSEK(PROBLEM)
0006 %   A wrapper function providing a standardized interface for using
0007 %   MOSEKOPT to solve 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 %       VTYPE : character string of length NX (number of elements in X),
0027 %               or 1 (value applies to all variables in x),
0028 %               allowed values are 'C' (continuous), 'B' (binary),
0029 %               'I' (integer).
0030 %       OPT : optional options structure with the following fields,
0031 %           all of which are also optional (default values shown in
0032 %           parentheses)
0033 %           verbose (0) - controls level of progress output displayed
0034 %               0 = no progress output
0035 %               1 = some progress output
0036 %               2 = verbose progress output
0037 %           skip_prices (0) - flag that specifies whether or not to
0038 %               skip the price computation stage, in which the problem
0039 %               is re-solved for only the continuous variables, with all
0040 %               others being constrained to their solved values
0041 %           price_stage_warn_tol (1e-7) - tolerance on the objective fcn
0042 %               value and primal variable relative match required to avoid
0043 %               mis-match warning message
0044 %           mosek_opt - options struct for MOSEK, value in verbose
0045 %                   overrides these options
0046 %       PROBLEM : The inputs can alternatively be supplied in a single
0047 %           PROBLEM struct with fields corresponding to the input arguments
0048 %           described above: H, c, A, l, u, xmin, xmax, x0, vtype, opt
0049 %
0050 %   Outputs:
0051 %       X : solution vector
0052 %       F : final objective function value
0053 %       EXITFLAG : exit flag
0054 %             1 = success
0055 %             0 = terminated at maximum number of iterations
0056 %            -1 = primal or dual infeasible
0057 %           < 0 = the negative of the MOSEK return code
0058 %       OUTPUT : output struct with the following fields:
0059 %           r - MOSEK return code
0060 %           res - MOSEK result struct
0061 %       LAMBDA : struct containing the Langrange and Kuhn-Tucker
0062 %           multipliers on the constraints, with fields:
0063 %           mu_l - lower (left-hand) limit on linear constraints
0064 %           mu_u - upper (right-hand) limit on linear constraints
0065 %           lower - lower bound on optimization variables
0066 %           upper - upper bound on optimization variables
0067 %
0068 %   Note the calling syntax is almost identical to that of QUADPROG
0069 %   from MathWorks' Optimization Toolbox. The main difference is that
0070 %   the linear constraints are specified with A, L, U instead of
0071 %   A, B, Aeq, Beq.
0072 %
0073 %   Calling syntax options:
0074 %       [x, f, exitflag, output, lambda] = ...
0075 %           miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
0076 %
0077 %       x = miqps_mosek(H, c, A, l, u)
0078 %       x = miqps_mosek(H, c, A, l, u, xmin, xmax)
0079 %       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0)
0080 %       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype)
0081 %       x = miqps_mosek(H, c, A, l, u, xmin, xmax, x0, vtype, opt)
0082 %       x = miqps_mosek(problem), where problem is a struct with fields:
0083 %                       H, c, A, l, u, xmin, xmax, x0, vtype, opt
0084 %                       all fields except 'c', 'A' and 'l' or 'u' are optional
0085 %       x = miqps_mosek(...)
0086 %       [x, f] = miqps_mosek(...)
0087 %       [x, f, exitflag] = miqps_mosek(...)
0088 %       [x, f, exitflag, output] = miqps_mosek(...)
0089 %       [x, f, exitflag, output, lambda] = miqps_mosek(...)
0090 %
0091 %   Example: (problem from from https://v8doc.sas.com/sashtml/iml/chap8/sect12.htm)
0092 %       H = [   1003.1  4.3     6.3     5.9;
0093 %               4.3     2.2     2.1     3.9;
0094 %               6.3     2.1     3.5     4.8;
0095 %               5.9     3.9     4.8     10  ];
0096 %       c = zeros(4,1);
0097 %       A = [   1       1       1       1;
0098 %               0.17    0.11    0.10    0.18    ];
0099 %       l = [1; 0.10];
0100 %       u = [1; Inf];
0101 %       xmin = zeros(4,1);
0102 %       x0 = [1; 0; 0; 1];
0103 %       opt = struct('verbose', 2);
0104 %       [x, f, s, out, lambda] = miqps_mosek(H, c, A, l, u, xmin, [], x0, vtype, opt);
0105 %
0106 %   See also MIQPS_MASTER, MOSEKOPT.
0107 
0108 %   MP-Opt-Model
0109 %   Copyright (c) 2010-2020, Power Systems Engineering Research Center (PSERC)
0110 %   by Ray Zimmerman, PSERC Cornell
0111 %
0112 %   This file is part of MP-Opt-Model.
0113 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0114 %   See https://github.com/MATPOWER/mp-opt-model for more info.
0115 
0116 %% check for Optimization Toolbox
0117 % if ~have_feature('mosek')
0118 %     error('miqps_mosek: requires MOSEK');
0119 % end
0120 
0121 %%----- input argument handling  -----
0122 %% gather inputs
0123 if nargin == 1 && isstruct(H)       %% problem struct
0124     p = H;
0125 else                                %% individual args
0126     p = struct('H', H, 'c', c, 'A', A, 'l', l, 'u', u);
0127     if nargin > 5
0128         p.xmin = xmin;
0129         if nargin > 6
0130             p.xmax = xmax;
0131             if nargin > 7
0132                 p.x0 = x0;
0133                 if nargin > 8
0134                     p.vtype = vtype;
0135                     if nargin > 9
0136                         p.opt = opt;
0137                     end
0138                 end
0139             end
0140         end
0141     end
0142 end
0143 
0144 %% define nx, set default values for H and c
0145 if ~isfield(p, 'H') || isempty(p.H) || ~any(any(p.H))
0146     if (~isfield(p, 'A') || isempty(p.A)) && ...
0147             (~isfield(p, 'xmin') || isempty(p.xmin)) && ...
0148             (~isfield(p, 'xmax') || isempty(p.xmax))
0149         error('miqps_mosek: LP problem must include constraints or variable bounds');
0150     else
0151         if isfield(p, 'A') && ~isempty(p.A)
0152             nx = size(p.A, 2);
0153         elseif isfield(p, 'xmin') && ~isempty(p.xmin)
0154             nx = length(p.xmin);
0155         else    % if isfield(p, 'xmax') && ~isempty(p.xmax)
0156             nx = length(p.xmax);
0157         end
0158     end
0159     p.H = sparse(nx, nx);
0160     qp = 0;
0161 else
0162     nx = size(p.H, 1);
0163     qp = 1;
0164 end
0165 if ~isfield(p, 'c') || isempty(p.c)
0166     p.c = zeros(nx, 1);
0167 end
0168 if ~isfield(p, 'x0') || isempty(p.x0)
0169     p.x0 = zeros(nx, 1);
0170 end
0171 if ~isfield(p, 'vtype') || isempty(p.vtype)
0172     p.vtype = '';
0173 end
0174 
0175 %% default options
0176 if ~isfield(p, 'opt')
0177     p.opt = [];
0178 end
0179 if ~isempty(p.opt) && isfield(p.opt, 'verbose') && ~isempty(p.opt.verbose)
0180     verbose = p.opt.verbose;
0181 else
0182     verbose = 0;
0183 end
0184 if ~isempty(p.opt) && isfield(p.opt, 'mosek_opt') && ~isempty(p.opt.mosek_opt)
0185     mosek_opt = mosek_options(p.opt.mosek_opt);
0186 else
0187     mosek_opt = mosek_options;
0188 end
0189 
0190 %% set up problem struct for MOSEK
0191 prob.c = p.c;
0192 if qp
0193    [prob.qosubi, prob.qosubj, prob.qoval] = find(tril(sparse(p.H)));
0194 end
0195 if isfield(p, 'A') && ~isempty(p.A)
0196     prob.a = sparse(p.A);
0197     nA = size(p.A, 1);
0198 else
0199     nA = 0;
0200 end
0201 if isfield(p, 'l') && ~isempty(p.A)
0202     prob.blc = p.l;
0203 end
0204 if isfield(p, 'u') && ~isempty(p.A)
0205     prob.buc = p.u;
0206 end
0207 if ~isempty(p.vtype)
0208     if length(p.vtype) == 1
0209         if p.vtype == 'I'
0210             prob.ints.sub = (1:nx);
0211         elseif p.vtype == 'B'
0212             prob.ints.sub = (1:nx);
0213             p.xmin = zeros(nx, 1);
0214             p.xmax = ones(nx, 1);
0215         end
0216     else
0217         k = find(p.vtype == 'B' | p.vtype == 'I');
0218         prob.ints.sub = k;
0219         k = find(p.vtype == 'B');
0220         if ~isempty(k)
0221             if isempty(p.xmin)
0222                 p.xmin = -Inf(nx, 1);
0223             end
0224             if isempty(p.xmax)
0225                 p.xmax = Inf(nx, 1);
0226             end
0227             p.xmin(k) = 0;
0228             p.xmax(k) = 1;
0229         end
0230     end
0231 end
0232 if isfield(p, 'xmin') && ~isempty(p.xmin)
0233     prob.blx = p.xmin;
0234 end
0235 if isfield(p, 'xmax') && ~isempty(p.xmax)
0236     prob.bux = p.xmax;
0237 end
0238 
0239 %% A is not allowed to be empty
0240 if ~isfield(prob, 'a') || isempty(prob.a)
0241     unconstrained = 1;
0242     prob.a = sparse(1, 1, 1, 1, nx);
0243     prob.blc = -Inf;
0244     prob.buc =  Inf;
0245 else
0246     unconstrained = 0;
0247 end
0248 sc = mosek_symbcon;
0249 s = have_feature('mosek', 'all');
0250 if isfield(prob, 'ints') && isfield(prob.ints, 'sub') && ~isempty(prob.ints.sub)
0251     mi = 1;
0252     if s.vnum >= 8
0253         mosek_opt.MSK_IPAR_OPTIMIZER = sc.MSK_OPTIMIZER_MIXED_INT;
0254 %     else
0255 %         mosek_opt.MSK_IPAR_OPTIMIZER = sc.MSK_OPTIMIZER_MIXED_INT_CONIC;
0256     end
0257 else
0258     mi = 0;
0259 end
0260 
0261 %%-----  run optimization  -----
0262 if verbose
0263     if s.vnum < 7
0264         alg_names = {           %% version 6.x
0265             'default',              %%  0 : MSK_OPTIMIZER_FREE
0266             'interior point',       %%  1 : MSK_OPTIMIZER_INTPNT
0267             '<conic>',              %%  2 : MSK_OPTIMIZER_CONIC
0268             '<qcone>',              %%  3 : MSK_OPTIMIZER_QCONE
0269             'primal simplex',       %%  4 : MSK_OPTIMIZER_PRIMAL_SIMPLEX
0270             'dual simplex',         %%  5 : MSK_OPTIMIZER_DUAL_SIMPLEX
0271             'primal dual simplex',  %%  6 : MSK_OPTIMIZER_PRIMAL_DUAL_SIMPLEX
0272             'automatic simplex',    %%  7 : MSK_OPTIMIZER_FREE_SIMPLEX
0273             '<mixed int>',          %%  8 : MSK_OPTIMIZER_MIXED_INT
0274             '<nonconvex>',          %%  9 : MSK_OPTIMIZER_NONCONVEX
0275             'concurrent'            %% 10 : MSK_OPTIMIZER_CONCURRENT
0276         };
0277     elseif s.vnum < 8
0278         alg_names = {           %% version 7.x
0279             'default',              %%  0 : MSK_OPTIMIZER_FREE
0280             'interior point',       %%  1 : MSK_OPTIMIZER_INTPNT
0281             '<conic>',              %%  2 : MSK_OPTIMIZER_CONIC
0282             'primal simplex',       %%  3 : MSK_OPTIMIZER_PRIMAL_SIMPLEX
0283             'dual simplex',         %%  4 : MSK_OPTIMIZER_DUAL_SIMPLEX
0284             'primal dual simplex',  %%  5 : MSK_OPTIMIZER_PRIMAL_DUAL_SIMPLEX
0285             'automatic simplex',    %%  6 : MSK_OPTIMIZER_FREE_SIMPLEX
0286             'network simplex',      %%  7 : MSK_OPTIMIZER_NETWORK_PRIMAL_SIMPLEX
0287             '<mixed int conic>',    %%  8 : MSK_OPTIMIZER_MIXED_INT_CONIC
0288             '<mixed int>',          %%  9 : MSK_OPTIMIZER_MIXED_INT
0289             'concurrent',           %% 10 : MSK_OPTIMIZER_CONCURRENT
0290             '<nonconvex>'           %% 11 : MSK_OPTIMIZER_NONCONVEX
0291         };
0292     else
0293         alg_names = {           %% version 8.x
0294             '<conic>',              %%  0 : MSK_OPTIMIZER_CONIC
0295             'dual simplex',         %%  1 : MSK_OPTIMIZER_DUAL_SIMPLEX
0296             'default',              %%  2 : MSK_OPTIMIZER_FREE
0297             'automatic simplex',    %%  3 : MSK_OPTIMIZER_FREE_SIMPLEX
0298             'interior point',       %%  4 : MSK_OPTIMIZER_INTPNT
0299             '<mixed int>',          %%  5 : MSK_OPTIMIZER_MIXED_INT
0300             'primal simplex'        %%  6 : MSK_OPTIMIZER_PRIMAL_SIMPLEX
0301         };
0302     end
0303     if qp
0304         lpqp = 'QP';
0305     else
0306         lpqp = 'LP';
0307     end
0308     if mi
0309         lpqp = ['MI' lpqp];
0310     end
0311     vn = have_feature('mosek', 'vstr');
0312     if isempty(vn)
0313         vn = '<unknown>';
0314     end
0315     fprintf('MOSEK Version %s -- %s %s solver\n', ...
0316             vn, alg_names{mosek_opt.MSK_IPAR_OPTIMIZER+1}, lpqp);
0317 end
0318 cmd = sprintf('minimize echo(%d)', verbose);
0319 [r, res] = mosekopt(cmd, prob, mosek_opt);
0320 
0321 %%-----  repackage results  -----
0322 if isfield(res, 'sol')
0323     if isfield(res.sol, 'int')
0324         sol = res.sol.int;
0325     elseif isfield(res.sol, 'bas')
0326         sol = res.sol.bas;
0327     else
0328         sol = res.sol.itr;
0329     end
0330     x = sol.xx;
0331 else
0332     sol = [];
0333     x = NaN(nx, 1);
0334 end
0335 
0336 %%-----  process return codes  -----
0337 eflag = -r;
0338 msg = '';
0339 switch (r)
0340     case sc.MSK_RES_OK
0341         if ~isempty(sol)
0342 %            if sol.solsta == sc.MSK_SOL_STA_OPTIMAL
0343             if strcmp(sol.solsta, 'OPTIMAL') || strcmp(sol.solsta, 'INTEGER_OPTIMAL')
0344                 msg = 'The solution is optimal.';
0345                 eflag = 1;
0346             else
0347                 eflag = -1;
0348 %                 if sol.prosta == sc.MSK_PRO_STA_PRIM_INFEAS
0349                 if strcmp(sol.prosta, 'PRIMAL_INFEASIBLE')
0350                     msg = 'The problem is primal infeasible.';
0351 %                 elseif sol.prosta == sc.MSK_PRO_STA_DUAL_INFEAS
0352                 elseif strcmp(sol.prosta, 'DUAL_INFEASIBLE')
0353                     msg = 'The problem is dual infeasible.';
0354                 else
0355                     msg = sol.solsta;
0356                 end
0357             end
0358         end
0359     case sc.MSK_RES_TRM_STALL
0360         if strcmp(sol.solsta, 'OPTIMAL') || strcmp(sol.solsta, 'INTEGER_OPTIMAL')
0361             msg = 'Stalled at or near optimal solution.';
0362             eflag = 1;
0363         else
0364             msg = 'Stalled.';
0365         end
0366     case sc.MSK_RES_TRM_MAX_ITERATIONS
0367         eflag = 0;
0368         msg = 'The optimizer terminated at the maximum number of iterations.';
0369     otherwise
0370         if isfield(res, 'rmsg') && isfield(res, 'rcodestr')
0371             msg = sprintf('%s : %s', res.rcodestr, res.rmsg);
0372         else
0373             msg = sprintf('MOSEK return code = %d', r);
0374         end
0375 end
0376 
0377 if (verbose || r == sc.MSK_RES_ERR_LICENSE || ...
0378         r == sc.MSK_RES_ERR_LICENSE_EXPIRED || ...
0379         r == sc.MSK_RES_ERR_LICENSE_VERSION || ...
0380         r == sc.MSK_RES_ERR_LICENSE_NO_SERVER_SUPPORT || ...
0381         r == sc.MSK_RES_ERR_LICENSE_FEATURE || ...
0382         r == sc.MSK_RES_ERR_LICENSE_INVALID_HOSTID || ...
0383         r == sc.MSK_RES_ERR_LICENSE_SERVER_VERSION || ...
0384         r == sc.MSK_RES_ERR_MISSING_LICENSE_FILE) ...
0385         && ~isempty(msg)  %% always alert user of license problems
0386     fprintf('%s\n', msg);
0387 end
0388 
0389 %%-----  repackage results  -----
0390 if nargout > 1
0391     if eflag == 1
0392         f = p.c' * x;
0393         if ~isempty(p.H)
0394             f = 0.5 * x' * p.H * x + f;
0395         end
0396     else
0397         f = [];
0398     end
0399     if nargout > 3
0400         output.r = r;
0401         output.res = res;
0402         if nargout > 4
0403             if ~isempty(sol)
0404                 if isfield(sol, 'slx')
0405                     lambda.lower = sol.slx;
0406                 else
0407                     lambda.lower = [];
0408                 end
0409                 if isfield(sol, 'sux')
0410                     lambda.upper = sol.sux;
0411                 else
0412                     lambda.upper = [];
0413                 end
0414                 if isfield(sol, 'slc')
0415                     lambda.mu_l  = sol.slc;
0416                 else
0417                     lambda.mu_l  = [];
0418                 end
0419                 if isfield(sol, 'suc')
0420                     lambda.mu_u  = sol.suc;
0421                 else
0422                     lambda.mu_u  = [];
0423                 end
0424             else
0425                 if isfield(p, 'xmin') && ~isempty(p.xmin)
0426                     lambda.lower = NaN(nx, 1);
0427                 else
0428                     lambda.lower = [];
0429                 end
0430                 if isfield(p, 'xmax') && ~isempty(p.xmax)
0431                     lambda.upper = NaN(nx, 1);
0432                 else
0433                     lambda.upper = [];
0434                 end
0435                 lambda.mu_l = NaN(nA, 1);
0436                 lambda.mu_u = NaN(nA, 1);
0437             end
0438             if unconstrained
0439                 lambda.mu_l  = [];
0440                 lambda.mu_u  = [];
0441             end
0442         end
0443     end
0444 end
0445 
0446 if mi && eflag == 1 && (~isfield(p.opt, 'skip_prices') || ~p.opt.skip_prices)
0447     if verbose
0448         fprintf('--- Integer stage complete, starting price computation stage ---\n');
0449     end
0450     if isfield(p.opt, 'price_stage_warn_tol') && ~isempty(p.opt.price_stage_warn_tol)
0451         tol = p.opt.price_stage_warn_tol;
0452     else
0453         tol = 1e-7;
0454     end
0455     pp = p;
0456     x(prob.ints.sub) = round(x(prob.ints.sub));
0457     pp.xmin(prob.ints.sub) = x(prob.ints.sub);
0458     pp.xmax(prob.ints.sub) = x(prob.ints.sub);
0459     pp.x0 = x;
0460     if qp
0461         pp.opt.mosek_opt.MSK_IPAR_OPTIMIZER = sc.MSK_OPTIMIZER_FREE;
0462     else
0463         pp.opt.mosek_opt.MSK_IPAR_OPTIMIZER = sc.MSK_OPTIMIZER_PRIMAL_SIMPLEX;
0464     end
0465     [x_, f_, eflag_, output_, lambda] = qps_mosek(pp);
0466     if eflag ~= eflag_
0467         error('miqps_mosek: EXITFLAG from price computation stage = %d', eflag_);
0468     end
0469     if abs(f - f_)/max(abs(f), 1) > tol
0470         warning('miqps_mosek: relative mismatch in objective function value from price computation stage = %g', abs(f - f_)/max(abs(f), 1));
0471     end
0472     xn = x;
0473     xn(abs(xn)<1) = 1;
0474     [mx, k] = max(abs(x - x_) ./ xn);
0475     if mx > tol
0476         warning('miqps_mosek: max relative mismatch in x from price computation stage = %g (%g)', mx, x(k));
0477     end
0478     output.price_stage = output_;
0479 end

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