17.1 Command Reference

The MOSEK toolbox provides a set of functions to interface to the MOSEK solver.

Main interface

mosekopt is the main interface to MOSEK.

Helper functions

These functions provide an easy-to-use but less flexible interface than the mosekopt function. In fact these procedures are just wrappers around the mosekopt interface and they are defined in MATLAB m-files.

  • msklpopt : Solves linear optimization problems.
  • mskqpopt : Solves quadratic optimization problems.
  • mskenopt : Solves entropy optimization problems.
  • mskgpopt : Solves geometric optimization problems.
  • mskscopt : Solves separable convex optimization problems.

I/O

  • mskgpwri : Write a geometric optimization problem to file.
  • mskgpread: Read a geometric optimization problem from file.

Options

MATLAB optimization toolbox compatible functions.

  • linprog : Solves linear optimization problems.
  • quadprog : Solves quadratic optimization problems.
  • intlinprog : Solves linear optimization problems with integer variables.
  • lsqlin : Solves least-squares with linear constraints.
  • lsqnonneg : Solves least-squares with non-negativity constraints.

17.1.1 Main Interface

rcode, res = mosekopt(cmd, prob, param, callback)

Solves an optimization problem. Data specifying the optimization problem can either be read from a file or be inputted directly from MATLAB. It also makes possible to write a file.

The following command strings are recognized for the cmd parameter:

  • anapro: Runs the problem analyzer.

  • echo(n): Controls how much information is echoed to the screen. n must be a nonnegative integer, where 0 means that no information is displayed and 3 means that all information is displayed.

  • info: Return the complete task information database in the field info of a res struct.

  • param: Return the complete parameter database in res.param.

  • primalrepair: Performs a primal feasibility repair. See Sec. 14 (Analyzing Infeasible Problems) for details.

  • maximize: Maximize the objective.

  • max : Sets the objective sense (similar to .maximize), without performing an optimization.

  • minimize: Minimize the objective.

  • min: Sets the objective sense (similar to .minimize), without performing an optimization.

  • nokeepenv: Delete the MOSEK environment after each run. This can increase the license checkout overhead significantly and is therefore only intended as a debug feature.

  • read(name): Request that data is read from a file name.

  • statuskeys(n): Controls the format of status keys (problem status, solution status etc.) in the returned problem:

    • statuskeys(0) – all the status keys are returned as strings,
    • statuskeys(1) – all the status keys are returned as numeric codes.
  • symbcon: Return the symbcon data structure in res.symbcon. See structure symbcon for details.

  • write(name): Write problem to the file name.

  • log(name): Write solver log-output to the file name.

  • version: Return the MOSEK version numbers in res.version.

Parameters:
 
  • cmd [in] (string) – The commands to be executed. By default it takes the value minimize.
  • prob [in] (prob) – [optional] a structure containing the problem data.
  • param [in] (struct) – [optional] a structure which is used to specify algorithmic parameters to MOSEK. The fields of param must be valid MOSEK parameter names. Moreover, the values corresponding to the fields must be of a valid type, i.e. the value of a string parameter must be a string, the value of an integer paramter must be an integer etc.
  • callback [in] (callback) – [optional] A MATLAB structure defining call-back data and functions.
Return:
 
  • rcode (rescode) – Return code. The interpretation of the value of the return code is listed in Sec. 17.6 (Enumerations).
  • res (res) – [optional] Solution obtained by the interior-point algorithm.

17.1.2 Helper Functions

res = msklpopt(c, a, blc, buc, blx, bux, param, cmd)

Solves a linear optimization problem of the form

\[\begin{split}\begin{array} {ll} min & c^T x \\ st. & \\ & blc \leq Ax \leq buc\\ & bux \leq x \leq bux. \end{array}\end{split}\]

Note

lc=[] and buc=[] means that the lower and upper bounds are plus and minus infinite respectively. The same interpretation is used for blx and bux. Note -inf is allowed in blc and blx. Similarly, inf is allowed in buc and bux.

Parameters:
 
  • c [in] (double[]) – The objective function vector.
  • a [in] (double[][]) – A (preferably sparse) matrix.
  • blc [in] (double[]) – Constraints lower bounds.
  • buc [in] (double[]) – Constraints upper bounds.
  • blx [in] (double[]) – Variables lower bounds.
  • bux [in] (double[]) – Variables upper bounds.
  • param [in] (list) – New MOSEK parameters.
  • cmd [in] (list) – [optional] The command list. See mosekopt for a list of available commands.
Return:
 

res (res) – [optional] Solution information.

res = mskqpopt(q, c, a, blc, buc, blx, bux, param, cmd)

Solves the optimization problem

\[\begin{split}\begin{array}{lc} \min & \half x^T Q x + c^T x\\ st. & \\ & blc \leq Ax \leq buc\\ & blx\leq x \leq bux\\ \end{array}\end{split}\]

Note

blc=[] and buc=[] means that the lower and upper bounds are plus and minus infinite respectively. The same interpretation is used for blx and bux. Note -inf is allowed in blc and blx. Similarly, inf is allowed in buc and bux.

Parameters:
 
  • q (double[]) – It is assumed that q is a symmetric positive semi-definite matrix.
  • c [in] (double[]) – A vector.
  • a (float[][]) – A (preferably) sparse matrix.
  • blc [in] (double[]) – Constraints lower bounds.
  • buc [in] (double[]) – Constraints upper bounds
  • blx [in] (double[]) – Variables lower bounds
  • bux [in] (double[]) – Variables upper bounds
  • param [in] (list) – MOSEK parameters.
  • cmd [in] (string) – [optional] The command list. See mosekopt for a list of available commands.
Return:
 

res (res) – [optional] Solution information.

res = mskenopt(d, c, a, blc, buc, param, cmd)

Solves the entropy optimization problem

\[\begin{split}\begin{array}{lc} \min & d^T (\sum x_i ln(x_i))+c^T x \\ \text{s.t.} & blc \leq Ax \leq buc\\ & x\in \R^n_+ \end{array}\end{split}\]

It is required that \(d \geq 0.0\).

Parameters:
 
  • d [in] (double[]) – A vector of non negative values.
  • c [in] (double[]) – A vector.
  • a [in] (double[][]) – A (preferably) sparse matrix.
  • blc [in] (double[]) – [optional] Constraints lower bounds.
  • buc [in] (double[]) – [optional] Constraints Upper bounds.
  • param [in] (list) – [optional] MOSEK parameters.
  • cmd [in] (list) – [optional] MOSEK commands. See mosekopt for a list of available commands.
Return:
 

res (res) – [optional] Solution information.

res = mskgpopt(c, a, map, param, cmd)

Solves the posynomial version of the geometric optimization problem in exponential form:

(1)\[\begin{split}\begin{array}{ll} \min & \log(sum(k \in find(map==0),c(k)*exp(a(k,:)*x))\\ st. & \log(sum(k \in find(map==i),c(k)*exp(a(k,:)*x)) <= 0, for k=1,...,max(map) \end{array}\end{split}\]

It is required that \(c>0.0\). See Sec. 7.3 (Geometric Optimization).

Parameters:
 
  • c [in] (double[]) – A vector.
  • a [in] (double[][]) – A (preferably) sparse matrix.
  • map (int[][]) – Corresponds to the set \(J\) in Sec. 7.3 (Geometric Optimization).
  • param [in] (list) – [optional] MOSEK parameters.
  • cmd [in] (list) – [optional] The command list. See mosekopt for a list of available commands.
Return:
 

res (res) – [optional] Solution information.

See also:
 

mskgpwri, mskgpread

res = mskscopt(opr, opri, oprj, oprf, oprg, c, a, blc, buc, blx, bux, param, cmd)

Solves separable convex optimization problems on the form

\[\begin{split}\begin{array} {ll} \min & c^T x + \sum_j f_j(x_j)\\ \text{s.t.} & blc \leq a_ x + \sum_j g_{kj}(x_j) \leq buc(k), k=1,...,size(a)\\ & blx \leq x \leq bux\\ \end{array}\end{split}\]

The nonlinear functions \(f_j\) and \(g_{kj}\) are specified using opr, opri, oprj oprf, oprg as follows. For all k between 1 and length(opri) then following nonlinear expression

if opr(k,:)=='ent'
  oprf(k) * x(oprj(k)) * log(x(oprj(k)))
elseif if opr(k,:)=='exp'
  oprf(k) * exp(oprg(k)*x(oprj(k)))
elseif if opr(k,:)=='log'
  oprf(k) * log(x(oprj(k)))
elseif if opr(k,:)=='pow'
  oprf(k) * x(oprj(k))^oprg(k)
else
  An invalid operator has been specified.

Is added to the objective if opri(k)=0. Otherwise it is added to constraint opri(k).

Parameters:
 
  • c [in] (double[]) – Is a vector.
  • a [in] (double[][]) – Is a (preferably) sparse matrix.
  • blc [in] (double[]) – [optional] Lower bounds on constraints.
  • buc [in] (double[]) – [optional] Upper bounds on constraints.
  • blx [in] (double[]) – [optional] Lower bounds on variables.
  • bux [in] (double[]) – [optional] Upper bounds on variables.
  • param [in] (list) – [optional] MOSEK parameters.
  • cmd [in] (list) – [optional] The command list. See mosekopt for a list of available commands.
Return:
 

res (res) – [optional] Solution information.

17.1.3 I/O

c, a, map = mskgpread(filename)

This function reads a Geometric Programming (gp) problem from a file compatible with the mskenopt command tool.

Parameters:
 

filename [in] (string) – The name of the file to read.

Return:
 
  • c (double[]) – Objective function coefficients. See problem (1).
  • a (double[]) – Linear constraints coefficients. See problem (1).
  • map (struct) – Data in the same format accepted by mskgpopt.
mskgpwri(c, a, map, filename)

This function writes a Geometric Programming (gp) problem to a file in a format compatible with the mskenopt command tool.

Parameters:
 
  • c [in] (double[]) – Objective function coefficients. See problem (1).
  • a [in] (double[]) – Linear constraints coefficients. See problem (1).
  • map [in] (struct) – Data in the same format accepted by mskgpopt.
  • filename (string) – The output file name.

17.1.4 Options

val = mskoptimget(options, param, default)

Obtains a value of an optimization parameter. See the mskoptimset function for which parameters that can be set.

Parameters:
 
  • options [in] (struct) – The optimization options structure.
  • param [in] (string) – Name of the optimization parameter for which the value should be obtained.
  • default [in] (string) – [optional] If param is not defined, the value of default is returned instead.
Return:
 

val (list) – Value of the required option. If the option does not exist, then [] is returned unless the value default is defined in which case the default value is returned.

options = mskoptimset(arg1, arg2, param1, value1, param2, value2, ...)

Obtains and modifies the optimization options structure. Only a subset of the fields in the optimization structure recognized by the MATLAB optimization toolbox is recognized by MOSEK. In addition the optimization options structure can be used to modify all the MOSEK specific parameters defined in Sec. 17.4 (Parameters (alphabetical list sorted by type)).

  • .Diagnostics Used to control how much diagnostic information is printed. Following values are accepted:

    off No diagnostic information is printed.
    on Diagnostic information is printed.
  • .Display Defines what information is displayed. The following values are accepted:

    off No output is displayed.
    iter Some output is displayed for each iteration.
    final Only the final output is displayed.
  • .MaxIter Maximum number of iterations allowed.

  • .Write A filename to write the problem to. If equal to the empty string no file is written. E.g the option Write(myfile.opf) writes the file myfile.opf in the opf format.

Parameters:
 
  • arg1 [in] (None) –

    [optional] Is allowed to be any of the following two things [in]:

    • Any string The same as using no argument.
    • A structure The argument is assumed to be a structure containing options, which are copied to the return options.
  • param1 [in] (string) – [optional] A string containing the name of a parameter that should be modified.
  • value1 [in] (None) – [optional] The new value assigned to the parameter with the name param1.
  • param2 [in] (None) – [optional] See param1.
  • value2 [in] (None) – [optional] See value1.
Return:
 

options (struct) – The updated optimization options structure.

17.1.5 MATLAB Optimization Toolbox Compatible Functions.

x, fval, exitflag, output = intlinprog(f, A, b, B, c, x0, options)
x, fval, exitflag, output = intlinprog(problem)

Solves the binary linear optimization problem:

\[\begin{split}\begin{array} {lc} \mbox{minimize} & f^T x \\ \mbox{subject to} & A x \leq b,\\ & B x = c,\\ & x \in \{0,1\}^n \end{array}\end{split}\]
Parameters:
 
  • f [in] (double[]) – The objective function.
  • A [in] (double[][]) – Constraint matrix for the inequalities. Use A=[] if there are no inequalities.
  • b [in] (double[]) – Right-hand side for the inequalities. Use b=[] if there are no inequalities.
  • B [in] (double[][]) – [optional] Constraint matrix for the equalities. Use B=[] if there are no equalities.
  • c [in] (double[]) – [optional] Right-hand side for the equalities. Use c=[] if there are no equalities.
  • x0 [in] (double[]) – [optional] A feasible starting point.
  • options [in] (struct) –

    [optional] An optimization options structure. See the mskoptimset function for the definition of the optimization options structure. intlinprog uses the options

    • .Diagnostics
    • .Display
    • .MaxTime The maximum number of seconds in the solution-time
    • .MaxNodes The maximum number of branch-and-bounds allowed
    • .Write Filename of problem file to save.
  • problem [in] (struct) – A structure containing the fields f, A, b, B, c, x0 and options.
Return:
 
  • x (double[]) – The solution \(x\).
  • fval (double) – The objective \(f^T x\).
  • exitflag (int) –

    A number which has the interpretation:

    • \(1\) The function returned an integer feasible solution.
    • \(-2\) The problem is infeasible.
    • \(-4\) maxNodes reached without converging.
    • \(-5\) maxTime reached without converging.
x, fval, exitflag, output, lambda = linprog(f, A, b, B, c, l, u, x0, options)
x, fval, exitflag, output, lambda = linprog(problem)

Solves the linear optimization problem:

\[\begin{split}\begin{array} {lc} \mbox{minimize} & f^T x \\ \mbox{subject to} & A x \leq b, \\ & B x = c, \\ & l \leq x \leq u. \end{array}\end{split}\]
Parameters:
 
  • f [in] (double[]) – The objective function.
  • A [in] (double[][]) – Constraint matrix for the inequalities. Use \(A=[]\) if there are no inequalities.
  • b [in] (double[]) – Right-hand side for the inequalities. Use \(b=[]\) if there are no inequalities.
  • B [in] (double[][]) – [optional] Constraint matrix for the equalities. Use \(B=[]\) if there are no equalities.
  • c [in] (double[]) – [optional] Right-hand side for the equalities. Use \(c=[]\) if there are no equalities.
  • l [in] (double[]) – [optional] Lower bounds on the variables. Use \(-\infty\) to represent infinite lower bounds.
  • u [in] (double[]) – [optional] Upper bounds on the variables. Use \(\infty\) to represent infinite upper bounds.
  • x0 [in] (double[]) – [optional] An initial guess for the starting point. Only used for the primal simplex algorithm. For more advanced warm-starting (e.g., using dual simplex), use mosekopt directly.
  • options [in] (struct) –

    [optional] An optimization options structure. See the mskoptimset function for the definition of the optimization options structure. linprog uses the options

    • .Diagnostics
    • .Display
    • .MaxIter
    • .Simplex Valid values are 'on', 'primal' or 'dual'. If Simplex is 'on' then MOSEK will use either a primal or dual simplex solver (similar as specifying "MSK_OPTIMIZER_FREE_SIMPLEX" in mosekopt; otherwise either a primal or dual simplex algorithm is used. Note, that the 'primal' and 'dual' values are specific for the MOSEK interface, and not present in the standard MATLAB version.
    • .Write Filename of problem file (e.g., 'prob.opf') to be saved. This is useful for reporting bugs or problems.
  • problem [in] (struct) – structure containing the fields f, A, b, B, c, l, u, x0 and options.
  • output [in] (struct) –

    A struct with the following fields

    • .iterations Number of iterations spent to reach the optimum.
    • .algorithm Always defined as ’large-scale [in]: interior-point’.
  • lambda [in] (struct) –

    A struct with the following fields

    • .lower Lagrange multipliers for lower bounds \(l\).
    • .upper Lagrange multipliers for upper bounds \(u\).
    • .ineqlin Lagrange multipliers for the inequalities.
    • .eqlin Lagrange multipliers for the equalities.
Return:
 
  • x (double[]) – The optimal \(x\) solution.
  • fval (double) – The optimal objective value, i.e. \(f^T x\).
  • exitflag (int) –

    A number which has the interpretation [in]:

    • \(<0\) The problem is likely to be either primal or dual infeasible.
    • \(=0\) The maximum number of iterations was reached.
    • \(>0\) \(x\) is an optimal solution.
x, resnorm, residual, exitflag, output, lambda = lsqlin(C, d, A, b, B, c, l, u, x0, options)

Solves the linear least squares problem:

(2)\[\begin{split}\begin{array} {lc} \mbox{minimize} & \half \left\| C x - d\right\|_2^2 \\ \mbox{subject to} & A x \leq b, \\ & B x = c, \\ & l \leq x \leq u. \end{array}\end{split}\]
Parameters:
 
  • C [in] (double[][]) – A matrix. See problem (2) for the purpose of the argument.
  • d [in] (double[]) – A vector. See problem (2) for the purpose of the argument.
  • A [in] (double[][]) – Constraint matrix for the inequalities. Use \(A=[]\) if there are no inequalities.
  • b [in] (double[]) – Right-hand side for the inequalities. Use \(b=[]\) if there are no inequalities.
  • B [in] (double[][]) – [optional] Constraint matrix for the equalities. Use \(B=[]\) if there are no equalities.
  • c [in] (double[]) – [optional] Right-hand side for the equalities. Use \(c=[]\) if there are no equalities.
  • l [in] (double[]) – [optional] Lower bounds on the variables. Use \(-\infty\) to represent infinite lower bounds.
  • u [in] (double[]) – [optional] Upper bounds on the variables. Use \(\infty\) to represent infinite lower bounds.
  • x0 [in] (double[]) – [optional] An initial guess for the starting point. This information is ignored by MOSEK
  • options [in] (struct) –

    [optional] An optimization options structure. See the function mskoptimset function for the definition of the optimization options structure. lsqlin uses the options

    • .Diagnostics
    • .Display
    • .MaxIter
    • .Write
Return:
 
  • x (double[]) – The optimal \(x\) solution.
  • resnorm (double) – The squared norm of the optimal residuals, i.e. \(\left\| Cx-d \right\|^2\) evaluated at the optimal solution.
  • residual (double) – The residual \(C x - d\).
  • exitflag (int) –

    A scalar which has the interpretation:

    • \(<0\) The problem is likely to be either primal or dual infeasible.
    • \(=0\) The maximum number of iterations was reached.
    • \(>0\) \(x\) is the optimal solution.
  • output (struct) –
    • .iterations Number of iterations spent to reach the optimum.
    • .algorithm Always defined as ’large-scale: interior-point’.
  • lambda (struct) –
    • .lower Lagrange multipliers for lower bounds \(l\).
    • .upper Lagrange multipliers for upper bounds \(u\).
    • .ineqlin Lagrange multipliers for inequalities.
    • .eqlin Lagrange multipliers for equalities.
x, resnorm, residual, exitflag, output, lambda = lsqnonneg(C, d, x0, options)

Solves the linear least squares problem:

(3)\[\begin{split}\begin{array} {lc} \mbox{minimize} & \half \left\|C x - d \right\|_2^2 \\ \mbox{subject to} & x \geq 0. \end{array}\end{split}\]

This procedure just provides an easy interface to lsqlin. Indeed all the procedure does is to call lsqlin with the appropriate arguments.

Parameters:
 
  • C [in] (double[][]) – See problem (3).
  • d [in] (double[]) – See problem (3).
  • x0 [in] (double[]) – [optional] An initial guess for the starting point. This information is ignored by MOSEK
  • options [in] (struct) –

    [optional] An optimizations options structure. See the mskoptimset function for the definition of the optimization options structure. lsqlin uses the options

    • .Diagnostics
    • .Display
    • .MaxIter
    • .Write
Return:
 
  • x (double[]) – The \(x\) solution.
  • resnorm (double) –

    The squared norm of the optimal residuals, i.e. .. math:: left| Cx-d right|^2

    evaluated at the optimal solution.

  • residual (double) – The residual \(C x - d\).
  • exitflag (int) –

    A number which has the interpretation:

    • \(<0\) The problem is likely to be either primal or dual infeasible.
    • \(=0\) The maximum number of iterations was reached.
    • \(>0\) \(x\) is optimal solution.
  • output (struct) –
    • .iterations Number of iterations spend to reach the optimum.
    • .algorithm Always defined to be ’large-scale: interior-point’.
  • lambda (struct) –
    • .lower Lagrange multipliers for lower bounds \(l\).
    • .upper Lagrange multipliers for upper bounds \(u\).
    • .ineqlin Lagrange multipliers for inequalities.
    • .eqlin Lagrange multipliers for equalities.
x, fval, exitflag, output, lambda = quadprog(H, f, A, b, B, c, l, u, x0, options)

Solves the quadratic optimization problem:

(4)\[\begin{split}\begin{array} {lc} \mbox{minimize} & \half x^T H x + f^T x \\ \mbox{subject to} & A x \leq b, \\ & B x = c, \\ & l \leq x \leq u. \end{array}\end{split}\]
Parameters:
 
  • H [in] (double[][]) – Hessian of the objective function. \(H\) must be a symmetric matrix. Contrary to the MATLAB optimization toolbox, MOSEK handles only the cases where \(H\) is positive semidefinite. On the other hand MOSEK always computes a global optimum, i.e. the objective function has to be strictly convex.
  • f [in] (double[]) – See (4) for the definition.
  • A [in] (double[][]) – Constraint matrix for the inequalities. Use \(A=[]\) if there are no inequalities.
  • b [in] (double[]) – Right-hand side for the nequalities. Use \(b=[]\) if there are no inequalities.
  • B [in] (double[][]) – [optional] Constraint matrix for the equalities. Use \(B=[]\) if there are no equalities.
  • c [in] (double[]) – [optional] Right-hand side for the equalities. Use \(c=[]\) if there are no equalities.
  • l [in] (double[]) – [optional] Lower bounds on the variables. Use \(-\infty\) to represent infinite lower bounds.
  • u [in] (double[]) – [optional] Upper bounds on the variables. Use \(\infty\) to represent infinite upper bounds.
  • x0 [in] (double[]) – [optional] An initial guess for the starting point. This information is ignored by MOSEK
  • options [in] (struct) –

    [optional] An optimization options structure. See the mskoptimset function for the definition of the optimizations options structure. quadprog uses the options

    • .Diagnostics
    • .Display
    • .MaxIter
    • .Write
Return:
 
  • x (double[]) – The \(x\) solution.
  • fval (double) – The optimal objective value i.e. \(\half x^T H x + f^T x\).
  • exitflag (int) –

    A scalar which has the interpretation:

    • \(<0\) The problem is likely to be either primal or dual infeasible.
    • \(=0\) The maximum number of iterations was reached.
    • \(>0\) \(x\) is an optimal solution.
  • output (struct) –

    A structure with the following fields

    • .iterations Number of iterations spent to reach the optimum.
    • .algorithm Always defined as ’large-scale: interior-point’.
  • lambda (struct) –

    A structure with the following fields

    • .lower Lagrange multipliers for lower bounds \(l\).
    • .upper Lagrange multipliers for upper bounds \(u\).
    • .ineqlin Lagrange multipliers for inequalities.
    • .eqlin Lagrange multipliers for equalities.