Home > matpower6.0 > 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)
   A wrapper function providing a MATPOWER 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 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_mosek(H, c, A, l, u, xmin, [], x0, vtype, opt);

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

Generated on Fri 16-Dec-2016 12:45:37 by m2html © 2005