17 List of examples¶
List of examples shipped in the distribution of Optimizer API for Julia:
File |
Description |
|---|---|
A simple problem with one affine conic constraint (ACC) |
|
A simple problem with two affine conic constraints (ACC) |
|
An example of data/progress callback |
|
A simple conic exponential problem |
|
Implementation of a concurrent optimizer for linear and mixed-integer problems |
|
A simple conic quadratic problem |
|
A simple problem with disjunctive constraints (DJC) |
|
A simple example of how to repair an infeasible problem |
|
A simple geometric program (GP) in conic form |
|
A Hello World example |
|
A simple linear problem |
|
A simple linear problem |
|
Implements logistic regression and simple log-sum-exp (CEO) |
|
A simple mixed-integer conic problem |
|
A simple mixed-integer linear problem |
|
A simple mixed-integer linear problem with an initial guess |
|
Uses MOSEK OptServer to solve an optimization problem asynchronously |
|
Uses MOSEK OptServer to solve an optimization problem synchronously |
|
Demonstrates parallel optimization using a batch method in MOSEK |
|
Shows how to set optimizer parameters and read information items |
|
Shows how to obtain and analyze a primal infeasibility certificate |
|
Portfolio optimization - basic Markowitz model |
|
Portfolio optimization - efficient frontier |
|
Portfolio optimization - market impact costs |
|
Portfolio optimization - transaction costs |
|
Portfolio optimization - cardinality constraints |
|
Portfolio optimization - factor model |
|
A simple power cone problem |
|
A simple quadratically constrained quadratic problem |
|
A simple quadratic problem |
|
Demonstrate how to modify and re-optimize a linear problem |
|
Demonstrates proper response handling |
|
A simple semidefinite problem with one matrix variable and a quadratic cone |
|
A simple semidefinite problem with two matrix variables |
|
A simple semidefinite problem with an LMI using the SVEC domain. |
|
Sensitivity analysis performed on a small linear problem |
|
A simple I/O example: read problem from a file, solve and write solutions |
|
Demonstrates how to examine the quality of a solution |
|
Demonstrates solving a linear system with the basis matrix |
|
Demonstrates solving a general linear system |
|
Shows how to find a Cholesky factorization of a sparse matrix |
Additional examples can be found on the MOSEK website and in other MOSEK publications.
acc1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : acc1.jl
#
# Purpose : Tutorial example for affine conic constraints.
# Models the problem:
#
# maximize c^T x
# subject to sum(x) = 1
# gamma >= |Gx+h|_2
##
using Mosek
# Define problem data
n, k = 3, 2
let Gsubi = Int64[1, 1, 2, 2],
Gsubj = Int32[1, 2, 1, 3],
Gval = Float64[1.5, 0.1, 0.3, 2.1]
# Make a MOSEK environment
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Attach a printer to the task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Create n free variables
appendvars(task,n)
putvarboundsliceconst(task,1, n+1, MSK_BK_FR, -Inf, Inf)
# Set up the objective
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
putclist(task,[1,2,3],[2.0,3.0,-1.0])
# One linear constraint - sum(x) = 1
appendcons(task,1)
putarow(task,1, [1,2,3], [1.0,1.0,1.0])
putconbound(task,1, MSK_BK_FX, 1.0, 1.0)
# Append empty AFE rows for affine expression storage
appendafes(task,k + 1)
# G matrix in sparse form
# Other data
h = Float64[0, 0.1]
gamma = 0.03
# Construct F matrix in sparse form
Fsubi = Int32[i + 1 for i in Gsubi] # G will be placed from row number 1 in F
Fsubj = Gsubj
Fval = Gval
# Fill in F storage
putafefentrylist(task,Fsubi, Fsubj, Fval)
# Fill in g storage
putafeg(task,1, gamma)
putafegslice(task,2, k+2, h)
# Define a conic quadratic domain
quadDom = appendquadraticconedomain(task,k + 1)
# Create the ACC
appendaccseq(task,
quadDom, # Domain index
1, # Indices of AFE rows [0,...,k]
nothing) # Ignored
# Solve and retrieve solution
optimize(task)
writedata(task,"acc1.ptf")
xx = getxx(task,MSK_SOL_ITR)
println("Solution: $xx")
# Demonstrate retrieving activity of ACC
activity = evaluateacc(task,MSK_SOL_ITR,
1) # ACC index
println("Activity of ACC:: $activity")
# Demonstrate retrieving the dual of ACC
doty = getaccdoty(task,MSK_SOL_ITR,
1) # ACC index
println("Dual of ACC: $doty")
end
end
acc2.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : acc2.jl
#
# Purpose : Tutorial example for affine conic constraints.
# Models the problem:
#
# maximize c^T x
# subject to sum(x) = 1
# gamma >= |Gx+h|_2
#
# This version inputs the linear constraint as an affine conic constraint.
##
using Mosek
# Define problem data
n = 3
k = 2
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Attach a printer to the task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Create n free variables
appendvars(task,n)
putvarboundsliceconst(task,1, n+1, MSK_BK_FR, -Inf, Inf)
# Set up the objective
c = Float64[2, 3, -1]
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
putclist(task,[1:n...], c)
# Set AFE rows representing the linear constraint
appendafes(task,1)
putafefrow(task,1, [1:n...], ones(n))
putafeg(task,1, -1.0)
# Set AFE rows representing the quadratic constraint
appendafes(task,k + 1)
putafefrow(task,3, # afeidx, row number
[1, 2], # varidx, column numbers
[1.5, 0.1]) # values
putafefrow(task,4, # afeidx, row number
[1, 3], # varidx, column numbers
[0.3, 2.1]) # values
h = [0, 0.1]
gamma = 0.03
putafeg(task,2, gamma)
putafegslice(task,3, k+2+1, h)
# Define domains
zeroDom = appendrzerodomain(task,1)
quadDom = appendquadraticconedomain(task,k + 1)
# Append affine conic constraints
appendacc(task,zeroDom, # Domain index
[1], # Indices of AFE rows
nothing) # Ignored
appendacc(task,quadDom, # Domain index
[2,3,4], # Indices of AFE rows
nothing) # Ignored
# Solve and retrieve solution
optimize(task)
writedata(task,"acc2.ptf")
@assert getsolsta(task,MSK_SOL_ITR) == MSK_SOL_STA_OPTIMAL
xx = getxx(task,MSK_SOL_ITR)
if getsolsta(task,MSK_SOL_ITR) == MSK_SOL_STA_OPTIMAL
println("Solution: $xx")
end
# Demonstrate retrieving activity of ACC
activity = evaluateacc(task,MSK_SOL_ITR,2)
println("Activity of quadratic ACC:: $activity")
# Demonstrate retrieving the dual of ACC
doty = getaccdoty(task,MSK_SOL_ITR,2)
println("Dual of quadratic ACC:: $doty")
end
callback.jl
using Mosek
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : callback.jl
#
# Purpose : To demonstrate how to use the progress
# callback.
#
# Use this script as follows:
#
# callback.py psim 25fv47.mps
# callback.py dsim 25fv47.mps
# callback.py intpnt 25fv47.mps
#
# The first argument tells which optimizer to use
# i.e. psim is primal simplex, dsim is dual simplex
# and intpnt is interior-point.
##
using Mosek
using Printf
function callback(task :: Mosek.Task,
maxtime :: Float64,
caller :: Callbackcode,
douinf :: Vector{Float64},
intinf :: Vector{Int32},
lintinf :: Vector{Int64})
opttime = 0.0
if caller == MSK_CALLBACK_BEGIN_INTPNT
println("Starting interior-point optimizer")
elseif caller == MSK_CALLBACK_INTPNT
itrn = intinf[MSK_IINF_INTPNT_ITER]
pobj = douinf[MSK_DINF_INTPNT_PRIMAL_OBJ]
dobj = douinf[MSK_DINF_INTPNT_DUAL_OBJ]
stime = douinf[MSK_DINF_INTPNT_TIME]
opttime = douinf[MSK_DINF_OPTIMIZER_TIME]
println("Iterations: $itrn")
@printf(" Elapsed time: %6.2f(%.2f) ",opttime, stime)
@printf(" Primal obj.: %-18.6e Dual obj.: %-18.6e",pobj, dobj)
elseif caller == MSK_CALLBACK_END_INTPNT
println("Interior-point optimizer finished.")
elseif caller == MSK_CALLBACK_BEGIN_PRIMAL_SIMPLEX
println("Primal simplex optimizer started.")
elseif caller == MSK_CALLBACK_UPDATE_PRIMAL_SIMPLEX
itrn = intinf[MSK_IINF_SIM_PRIMAL_ITER]
pobj = douinf[MSK_DINF_SIM_OBJ]
stime = douinf[MSK_DINF_SIM_TIME]
opttime = douinf[MSK_DINF_OPTIMIZER_TIME]
println("Iterations: %-3d", itrn)
println(" Elapsed time: %6.2f(%.2f)",opttime, stime)
println(" Obj.: %-18.6e", pobj)
elseif caller == MSK_CALLBACK_END_PRIMAL_SIMPLEX
println("Primal simplex optimizer finished.")
elseif caller == MSK_CALLBACK_BEGIN_DUAL_SIMPLEX
println("Dual simplex optimizer started.")
elseif caller == MSK_CALLBACK_UPDATE_DUAL_SIMPLEX
itrn = intinf[MSK_IINF_SIM_DUAL_iter]
pobj = douinf[MSK_DINF_SIM_OBJ]
stime = douinf[MSK_DINF_SIM_TIME]
opttime = douinf[MSK_DINF_OPTIMIZER_TIME]
println("Iterations: %-3d",itrn)
println(" Elapsed time: %6.2f(%.2f)",opttime, stime)
println(" Obj.: %-18.6e",pobj)
elseif caller == MSK_CALLBACK_END_DUAL_SIMPLEX
println("Dual simplex optimizer finished.")
elseif caller == MSK_CALLBACK_NEW_INT_MIO
println("New integer solution has been located.")
xx = getxx(task,MSK_SOL_ITG)
println(" x = $xx")
println("Obj.: %f", douinf[MSK_DINF_MIO_OBJ_INT])
end
if opttime >= maxtime
# mosek is spending too much time. Terminate it.
println("Terminating.")
1
else
0
end
end
# To run a continuous problem example
if length(ARGS) < 2
println("Usage: callback ( psim | dsim | intpnt ) filename")
else
slvr = ARGS[1]
filename = ARGS[2]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
readdata(task,filename)
if slvr == "psim"
putintparam(task,MSK_IPAR_OPTIMIZER, MSK_OPTIMIZER_PRIMAL_SIMPLEX)
elseif slvr == "dsim"
putintparam(task,MSK_IPAR_OPTIMIZER, MSK_OPTIMIZER_DUAL_SIMPLEX)
elseif slvr == "intpnt"
putintparam(task,MSK_IPAR_OPTIMIZER, MSK_OPTIMIZER_INTPNT)
end
putcallbackfunc(task,(caller,dinf,iinf,linf) -> callback(task,0.05,caller,dinf,iinf,linf))
optimize(task)
end
end
ceo1.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : ceo1.jl
#
# Purpose : Demonstrates how to solve small conic exponential
# optimization problem using the MOSEK Python API.
##
using Mosek
using Printf, SparseArrays
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Attach a printer to the task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
c = [1.0, 1.0, 0.0]
a = [1.0, 1.0, 1.0]
numvar = 3
numcon = 1
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
# Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Set up the linear part of the problem
putcslice(task,1, numvar+1, c)
putarow(task,1, [1, 2, 3], a)
putvarboundsliceconst(task,1, numvar+1, MSK_BK_FR, -Inf, Inf)
putconbound(task,1, MSK_BK_FX, 1.0, 1.0)
# Add a conic constraint
# Create a 3x3 identity matrix F
appendafes(task,3)
putafefentrylist(task,
[1, 2, 3], # Rows
[1, 2, 3], # Columns
ones(3))
# Exponential cone (x(0),x(1),x(2)) \in EXP
expdomain = appendprimalexpconedomain(task)
appendacc(task,
expdomain, # Domain
[1, 2, 3], # Rows from F
nothing) # Unused
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
# Optimize the task
optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
# Output a solution
xx = getxx(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
println("Optimal solution: $xx")
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility.")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
concurrent1.jl
#
# File: concurrent1.jl
#
# Purpose: Demonstrates a simple implementation of a concurrent optimizer.
#
# The concurrent optimizer starts a few parallel optimizations
# of the same problem using different algorithms, and reports
# a solution when the first optimizer is ready.
#
# This example also demonstrates how to define a simple callback handler
# that stops the optimizer when requested.
using Mosek
# Defines a Mosek callback function whose only function
# is to indicate if the optimizer should be stopped.
stop = false
firstStop = 0
function callback(caller :: Callbackcode,
douinf :: Vector{Float64},
intinf :: Vector{Int32},
lintinf :: Vector{Int64})
if stop
1
else
0
end
end
# firstOK, res, trm = optimize(tasks)
#
# Takes a list of tasks and optimizes then in parallel. The
# response code and termination code from each optimization is
# returned in ``res`` and ``trm``.
#
# When one task completes with rescode == ok, others are terminated.
#
# firstOK is the index of the first task that returned
# with rescode == ok. Whether or not this task contains the
# most valuable answer, is for the caller to decide. If none
# completed without error returns -1.
function runTask(num, task)
global stop
global firstStop
ttrm =
try
optimize(task)
catch e
if isa(e,MosekError)
return r.rcode,MSK_RES_ERR_UNKNOWN
else
rethrow()
end
end
# If this finished with success, inform other tasks to interrupt
# Note that data races around stop/firstStop are irrelevant
if ! stop
stop = true
firstStop = num
end
return MSK_RES_OK,ttrm
end
function optimizeconcurrent(tasks::Vector{Mosek.Task})
res = [ MSK_RES_ERR_UNKNOWN for t in tasks ]
trm = [ MSK_RES_ERR_UNKNOWN for t in tasks ]
# Set a callback function
for t in tasks
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putcallbackfunc(t, callback)
end
# Start parallel optimizations, one per task
Threads.@threads for i in 1:length(tasks)
(tres,ttrm) = runTask(i,tasks[i])
res[i] = tres
trm[i] = ttrm
end
# For debugging, print res and trm codes for all optimizers
for (i,(tres,ttrm)) in enumerate(zip(res,trm))
println("Optimizer $i res $tres trm $ttrm")
end
return firstStop, res, trm
end
#
# idx, winTask, winTrm, winRes = optimizeconcurrent(task, optimizers)
#
# Given a continuous task, set up jobs to optimize it
# with a list of different solvers.
#
# Returns an index, corresponding to the optimization
# task that is returned as winTask. This is the task
# with the best possible status of those that finished.
# If none task is considered successful returns -1.
function optimizeconcurrent(task, optimizers)
# Choose various optimizers for cloned tasks
tasks = Mosek.Task[ let t = maketask(task)
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putintparam(t,MSK_IPAR_OPTIMIZER, opt)
t
end for opt in optimizers ]
# Solve tasks in parallel
firstOK, res, trm = optimizeconcurrent(tasks)
if firstOK > 0
return firstOK, tasks[firstOK], trm[firstOK], res[firstOK]
else
return 0, Nothing, Nothing, Nothing
end
end
#
# idx, winTask, winTrm, winRes = optimizeconcurrent(task, optimizers)
#
# Given a mixed-integer task, set up jobs to optimize it
# with different values of seed. That will lead to
# different execution paths of the optimizer.
#
# Returns an index, corresponding to the optimization
# task that is returned as winTask. This is the task
# with the best value of the objective function.
# If none task is considered successful returns -1.
#
# Typically, the input task would contain a time limit. The two
# major scenarios are:
# 1. Some clone ends before time limit - then it has optimum.
# 2. All clones reach time limit - pick the one with best objective.
function optimizeconcurrentMIO(task, seeds)
# Choose various seeds for cloned tasks
tasks = Mosek.Task[ let t = maketask(task)
putintparam(MSK_IPAR_MIO_SEED, seed)
t
end for seed in seeds ]
# Solve tasks in parallel
(firstOK, res, trm) = optimizeconcurrent(tasks)
sense = getobjsense(task)
bestObj = if sense == MSK_OBJECTIVE_SENSE_MINIMIZE 1.0e+10 else -1.0e+10 end
bestPos = -1
if firstOK >= 0
# Pick the task that ended with res = ok
# and contains an integer solution with best objective value
for (i,t) in enumerate(tasks)
pobj = getprimalobj(t,MSK_SOL_ITG)
print("$i $pobj")
end
for (i,(tres,ttrm,t)) in enumerate(zip(res,trm,tasks))
solsta = getsolsta(t,MSK_SOL_ITG)
if tres == MSK_RES_OK &&
( solsta == MSK_SOL_STA_PRIM_FEAS ||
solsta == MSK_SOL_STA_INTEGER_OPTIMAL)
pobj = getprimalobj(t,MSK_SOL_ITG)
if ( ( sense == MSK_OBJECTIVE_SENSE_MINIMIZE &&
getprimalobj(t,MSK_SOL_ITG) < bestObj ) ||
( sense == MSK_OBJECTIVE_SENSE_MAXIMIZE &&
getprimalobj(t,MSK_SOL_ITG) > bestObj ) )
bestObj = pobj
bestPos = i
end
end
end
end
if bestPos > 0
return bestPos, tasks[bestPos], trm[bestPos], res[bestPos]
else
return 0, Nothing, Nothing, Nothing
end
end
# This is an example of how one can use the methods
# optimizeconcurrent
# optimizeconcurrentMIO
#
# argv[0] : name of file with input problem
# argv[1]: (optional) time limit
function main(fname::String,tlimit)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
if fname != "-"
readdata(task,fname)
else
readptfstring(task,lo1_ptf)
end
putobjname(task,"concurrent1")
# Optional time limit
if tlimit !== Nothing
putdouparam(task,MSK_DPAR_OPTIMIZER_MAX_TIME, tlimit)
end
let (idx,t,trm,res) = (
if getnumintvar(task) == 0
# If the problem is continuous
# optimize it with three continuous optimizers.
# (Simplex will fail for non-linear problems)
optimizers = [ MSK_OPTIMIZER_CONIC,
MSK_OPTIMIZER_DUAL_SIMPLEX,
MSK_OPTIMIZER_PRIMAL_SIMPLEX ]
optimizeconcurrent(task, optimizers)
else
# Mixed-integer problem.
# Try various seeds.
seeds = [ 42, 13, 71749373 ]
optimizeconcurrentMIO(task, seeds)
end )
# Check results and print the best answer
if idx > 0
println("Result from optimizer with index $idx: res $res trm $trm")
putstreamfunc(t,MSK_STREAM_LOG,msg -> print(msg))
optimizersummary(t,MSK_STREAM_LOG)
solutionsummary(t,MSK_STREAM_LOG)
else
println("All optimizers failed.")
end
end
end
end
const lo1_ptf = """Task
Objective
Maximize + 3 @x0 + @x1 + 5 @x2 + @x3
Constraints
@c0 [30] + 3 @x0 + @x1 + 2 @x2
@c1 [15;+inf] + 2 @x0 + @x1 + 3 @x2 + @x3
@c2 [-inf;25] + 2 @x1 + 3 @x3
Variables
@x0 [0;+inf]
@x1 [0;10]
@x2 [0;+inf]
@x3 [0;+inf]
"""
let fname = if length(ARGS) < 1 "../data/25fv47.mps" else ARGS[1] end,
tlimit = if length(ARGS) < 2 Nothing else parse(Float64,ARGS[2]) end
optimize(maketask()) # Just to initialize the thread pool from the main thread
main(fname,tlimit)
end
cqo1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : cqo1.jl
#
# Purpose: Demonstrates how to solve small conic
# optimization problem using the MOSEK Python API.
##
using Mosek
printstream(msg::AbstractString) = print(msg)
callback(where,dinf,iinf,liinf) = 0
# Since the actual value of Infinity is ignores, we define it solely
# for symbolic purposes:
bkc = [ MSK_BK_FX ]
blc = [ 1.0 ]
buc = [ 1.0 ]
c = [ 0.0, 0.0, 0.0,
1.0, 1.0, 1.0 ]
bkx = [ MSK_BK_LO,MSK_BK_LO,MSK_BK_LO,
MSK_BK_FR,MSK_BK_FR,MSK_BK_FR ]
blx = [ 0.0, 0.0, 0.0,
-Inf, -Inf, -Inf ]
bux = [ Inf, Inf, Inf,
Inf, Inf, Inf ]
asub = [ 1 ,2, 3 ]
aval = [ 1.0, 1.0, 2.0 ]
numvar = length(bkx)
numcon = length(bkc)
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,printstream)
putcallbackfunc(task,callback)
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
#Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Set the linear term c_j in the objective.
putclist(task,[1:6;],c)
# Set the bounds on variable j
# blx[j] <= x_j <= bux[j]
putvarboundslice(task,1,numvar+1,bkx,blx,bux)
putarow(task,1,asub,aval)
putconbound(task,1,bkc[1],blc[1],buc[1])
# Input the cones
appendafes(task,6)
putafefentrylist(task,
[1, 2, 3, 4, 5, 6], # Rows
[4, 1, 2, 5, 6, 3], # Columns */
[1.0, 1.0, 1.0, 1.0, 1.0, 1.0])
# Quadratic cone (x(3),x(0),x(1)) \in QUAD_3
quadcone = appendquadraticconedomain(task,3)
appendacc(task,
quadcone, # Domain
[1, 2, 3], # Rows from F
nothing)
# Rotated quadratic cone (x(4),x(5),x(2)) \in RQUAD_3
rquadcone = appendrquadraticconedomain(task,3)
appendacc(task,
rquadcone, # Domain
[4, 5, 6], # Rows from F
nothing);
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
# Optimize the task
#optimize(task,"mosek://solve.mosek.com:30080")
optimize(task)
writedata(task,"cqo1.ptf")
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Output a solution
xx = getxx(task,MSK_SOL_ITR)
println("Optimal solution: $xx")
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
djc1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : djc1.jl
#
# Purpose: Demonstrates how to solve the problem with two disjunctions:
#
# minimize 2x0 + x1 + 3x2 + x3
# subject to x0 + x1 + x2 + x3 >= -10
# (x0-2x1<=-1 and x2=x3=0) or (x2-3x3<=-2 and x1=x2=0)
# x0=2.5 or x1=2.5 or x2=2.5 or x3=2.5
#
using Mosek
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Append free variables
numvar = 4
appendvars(task,numvar)
putvarboundsliceconst(task,1, numvar+1, MSK_BK_FR, -Inf, Inf)
# The linear part: the linear constraint
appendcons(task,1)
putarow(task,1,[1:numvar...], ones(numvar))
putconbound(task,1,MSK_BK_LO, -10.0, -10.0)
# The linear part: objective
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
putclist(task,[1:numvar...], Float64[2, 1, 3, 1])
# Fill in the affine expression storage F, g
numafe = 10
appendafes(task,numafe)
fafeidx = Int64[1, 1, 2, 2, 3, 4, 5, 6, 7, 8, 9, 10]
fvaridx = Int32[1, 2, 3, 4, 1, 2, 3, 4, 1, 2, 3, 4]
fval = Float64[1.0, -2.0, 1.0, -3.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]
g = Float64[1.0, 2.0, 0.0, 0.0, 0.0, 0.0, -2.5, -2.5, -2.5, -2.5]
putafefentrylist(task,fafeidx, fvaridx, fval)
putafegslice(task,1, numafe+1, g)
# Create domains
zero1 = appendrzerodomain(task,1)
zero2 = appendrzerodomain(task,2)
rminus1 = appendrminusdomain(task,1)
# Append disjunctive constraints
numdjc = 2
appenddjcs(task,numdjc)
# First disjunctive constraint
putdjc(task,1, # DJC index
[rminus1, zero2, rminus1, zero2], # Domains (domidxlist)
[1, 5, 6, 2, 3, 4], # AFE indices (afeidxlist)
nothing, # Unused
[2, 2] ) # Term sizes (termsizelist)
# Second disjunctive constraint
putdjc(task,2, # DJC index
[zero1, zero1, zero1, zero1], # Domains (domidxlist)
[7, 8, 9, 10], # AFE indices (afeidxlist)
nothing, # Unused
[1, 1, 1, 1] ) # Term sizes (termidxlist)
# Useful for debugging
writedata(task,"djc1.ptf") # Write file in human-readable format
# Solve the problem
optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
# Get status information about the solution
sta = getsolsta(task,MSK_SOL_ITG)
println(sta)
@assert sta == MSK_SOL_STA_INTEGER_OPTIMAL
xx = getxx(task,MSK_SOL_ITG)
println("Optimal solution:")
for i in 1:numvar
println("x[$i]=$(xx[i])")
end
@assert maximum(abs.(xx - [0.0, 0.0, -12.5, 2.5])) < 1e-7
end
feasrepairex1.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : feasrepairex1.jl
#
# Purpose : To demonstrate how to use the MSK_relaxprimal function to
# locate the cause of an infeasibility.
#
# Syntax : On command line
#
# feasrepairex1.jl [ feasrepair.lp | - ]
#
# feasrepair.lp is located in mosek\<version>\tools\examples.
feasrepair_lp = """
minimize
obj: - 10 x1 - 9 x2
st
c1: + 7e-01 x1 + x2 <= 630
c2: + 5e-01 x1 + 8.333333333e-01 x2 <= 600
c3: + x1 + 6.6666667e-01 x2 <= 708
c4: + 1e-01 x1 + 2.5e-01 x2 <= 135
bounds
x2 >= 650
end
"""
using Mosek
if length(ARGS) < 1
println("Syntax: feasrepairex1 [ FILENAME | - ]")
else
filename = ARGS[1]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
if filename != "-"
readdata(task,filename)
else
readlpstring(task,feasrepair_lp)
end
putintparam(task,MSK_IPAR_LOG_FEAS_REPAIR, 3)
numvar = getnumvar(task)
numcon = getnumcon(task)
println(numvar,",",numcon)
primalrepair(task,zeros(numcon),zeros(numcon),zeros(numvar),zeros(numvar))
sum_viol = getdouinf(task,MSK_DINF_PRIMAL_REPAIR_PENALTY_OBJ);
println("Minimized sum of violations = $(sum_viol)");
optimize(task)
solutionsummary(task,MSK_STREAM_LOG)
end
end
gp1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : gp1.jl
#
# Purpose: Demonstrates how to solve a simple Geometric Program (GP)
# cast into conic form with exponential cones and log-sum-exp.
#
# Example from
# https:#gpkit.readthedocs.io/en/latest/examples.html#maximizing-the-volume-of-a-box
#
using Mosek
function max_volume_box(Aw :: Float64,
Af :: Float64,
alpha :: Float64,
beta :: Float64,
gamma :: Float64,
delta :: Float64)
numvar = 3 # Variables in original problem
# Create the optimization task.
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Add variables and constraints
appendvars(task,numvar)
x = 1
y = 2
z = 3
# Objective is the sum of three first variables
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
putcslice(task,1, numvar+1, [1.0,1.0,1.0])
putvarboundsliceconst(task,1, numvar+1, MSK_BK_FR, -Inf, Inf)
appendcons(task,3)
# s0+s1 < 1 <=> log(s0+s1) < 0
putaijlist(task,
[1,1,2,2,3,3],
[y, z, x, y, z, y],
[1.0, 1.0, 1.0, -1.0, 1.0, -1.0])
putconbound(task,1,MSK_BK_UP,-Inf,log(Af))
putconbound(task,2,MSK_BK_RA,log(alpha),log(beta))
putconbound(task,3,MSK_BK_RA,log(gamma),log(delta))
let afei = getnumafe(task)+1,
u1 = getnumvar(task)+1,
u2 = u1+1,
afeidx = [1, 2, 3, 3, 4, 4, 6, 6],
varidx = [u1, u2, x, y, x, z, u1, u2],
fval = [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],
gfull = [0.0, 0.0, log(2.0/Aw), log(2.0/Aw), 1.0, -1.0]
appendvars(task,2)
appendafes(task,6)
putvarboundsliceconst(task,u1, u1+2, MSK_BK_FR, -Inf, Inf)
# Affine expressions appearing in affine conic constraints
# in this order:
# u1, u2, x+y+log(2/Awall), x+z+log(2/Awall), 1.0, u1+u2-1.0
putafefentrylist(task,afeidx, varidx, fval)
putafegslice(task,afei, afei+6, gfull)
let dom = appendprimalexpconedomain(task)
# (u1, 1, x+y+log(2/Awall)) \in EXP
appendacc(task,dom, [1, 5, 3], nothing)
# (u2, 1, x+z+log(2/Awall)) \in EXP
appendacc(task,dom, [2, 5, 4], nothing)
end
let dom = appendrzerodomain(task,1)
# The constraint u1+u2-1 \in \ZERO is added also as an ACC
appendacc(task,dom, [6], nothing)
end
end
optimize(task)
writedata(task,"gp1.ptf")
exp.(getxxslice(task,MSK_SOL_ITR, 1, numvar+1))
end # maketask
end # max_volume_box
# maximize h*w*d
# subjecto to 2*(h*w + h*d) <= Awall
# w*d <= Afloor
# alpha <= h/w <= beta
# gamma <= d/w <= delta
#
# Variable substitutions: h = exp(x), w = exp(y), d = exp(z).
#
# maximize x+y+z
# subject log( exp(x+y+log(2/Awall)) + exp(x+z+log(2/Awall)) ) <= 0
# y+z <= log(Afloor)
# log( alpha ) <= x-y <= log( beta )
# log( gamma ) <= z-y <= log( delta )
#
# Finally, the model we will implement:
#
# maximize x+y+z
# subject to s0 > exp(x+y+log(2/Awall); (s0,1,x+y+log(2/Awall)) in PEXP
# s1 > exp(x+z+log(2/Awall); (s1,1,x+z+log(2/Awall)) in PEXP
# s0+s1 < 1
#
# y+z < log Afloor
#
# x-y in [log alpha; log beta]
# z-y in [log gamma; log delta]
#
# (x,y,z) in pexp : x0 > x1 * exp(x2/x1)
hwd = let Aw = 200.0,
Af = 50.0,
alpha = 2.0,
beta = 10.0,
gamma = 2.0,
delta = 10.0
max_volume_box(Aw, Af, alpha, beta, gamma, delta)
end
println("h=$(hwd[1]) w=$(hwd[2]) d=$(hwd[3])\n");
helloworld.jl
##
# Copyright: Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File: helloworld.jl
#
# The most basic example of how to get started with MOSEK.
##
using Mosek
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
appendvars(task, 1) # 1 variable x
putcj(task, 1, 1.0) # c_0 = 1.0
putvarbound(task, 1, MSK_BK_RA, 2.0, 3.0) # 2.0 <= x <= 3.0
putobjsense(task, MSK_OBJECTIVE_SENSE_MINIMIZE) # minimize
optimize(task) # Optimize
x = getxx(task, MSK_SOL_ITR) # Get solution
println("Solution x = $(x[1])") # Print solution
end
lo1.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : lo1.jl
#
# Purpose : Demonstrates how to solve small linear
# optimization problem using the MOSEK Python API.
##
using Mosek
using Printf, SparseArrays
############################
## Define problem data
bkc = [MSK_BK_FX
MSK_BK_LO
MSK_BK_UP]
# Bound values for constraints
blc = [30.0, 15.0, -Inf]
buc = [30.0, +Inf, 25.0]
# Bound keys for variables
bkx = [ MSK_BK_LO
MSK_BK_RA
MSK_BK_LO
MSK_BK_LO ]
# Bound values for variables
blx = [ 0.0, 0.0, 0.0, 0.0]
bux = [+Inf, 10.0, +Inf, +Inf]
numvar = length(bkx)
numcon = length(bkc)
# Objective coefficients
c = [ 3.0, 1.0, 5.0, 1.0 ]
# Below is the sparse representation of the A
# matrix stored by column.
A = sparse([1, 2, 1, 2, 3, 1, 2, 2, 3],
[1, 1, 2, 2, 2, 3, 3, 4, 4],
[3.0, 2.0, 1.0, 1.0, 2.0, 2.0, 3.0, 1.0, 3.0 ],
numcon,numvar)
############################
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
putobjname(task,"lo1")
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
for i=1:numcon
putconname(task,i,@sprintf("c%02d",i))
end
# Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
for j=1:numvar
putvarname(task,j,@sprintf("x%02d",j))
end
putclist(task,[1,2,3,4], c)
putacolslice(task,1,numvar+1,A)
putvarboundslice(task, 1, numvar+1, bkx,blx,bux)
# Set the bounds on constraints.
# blc[i] <= constraint_i <= buc[i]
putconboundslice(task,1,numcon+1,bkc,blc,buc)
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Solve the problem
optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
# Get status information about the solution
solsta = getsolsta(task,MSK_SOL_BAS)
if solsta == MSK_SOL_STA_OPTIMAL
xx = getxx(task,MSK_SOL_BAS)
print("Optimal solution:")
println(xx)
elseif solsta in [ MSK_SOL_STA_DUAL_INFEAS_CER,
MSK_SOL_STA_PRIM_INFEAS_CER ]
println("Primal or dual infeasibility certificate found.\n")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
@printf("Other solution status (%d)\n",solsta)
end
end
lo2.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : lo2.jl
#
# Purpose : Demonstrates how to solve a small linear
# optimization problem using the MOSEK Java API.
using Mosek
let numcon = 3,
numvar = 4,
NUMANZ = 9,
c = [3.0, 1.0, 5.0, 1.0],
aptrb = [1,4,8],
aptre = [4,8,10],
asub = [ 1, 2, 3,
1, 2, 3, 4,
2, 4 ],
aval = [ 3.0, 1.0, 2.0 ,
2.0, 1.0, 3.0, 1.0,
2.0, 3.0 ],
bkc = [ MSK_BK_FX,
MSK_BK_LO,
MSK_BK_UP ],
blc = [30.0,
15.0,
-Inf
],
buc = [30.0,
+Inf,
25.0
],
bkx = [ MSK_BK_LO,
MSK_BK_RA,
MSK_BK_LO,
MSK_BK_LO ],
blx = [ 0.0,
0.0,
0.0,
0.0 ],
bux = [ Inf,
10.0,
Inf,
Inf ]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
# Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
for j in 1:numvar
# Set the linear term c_j in the objective.
putcj(task,j,c[j])
# Set the bounds on variable j.
# blx[j] <= x_j <= bux[j]
putvarbound(task,j, bkx[j], blx[j], bux[j])
end
# Set the bounds on constraints.
# for i=1, ...,numcon : blc[i] <= constraint i <= buc[i]
for i in 1:numcon
putconbound(task,i, bkc[i], blc[i], buc[i])
end
#Input rows of A #
putarowslice(task,
1,numcon+1,
aptrb,aptre,
asub,
aval)
# A maximization problem
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Solve the problem
trm = optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
# Get status information about the solution
solsta = getsolsta(task,MSK_SOL_BAS)
x = getxx(task,MSK_SOL_BAS) # Basic solution.
if solsta == MSK_SOL_STA_OPTIMAL
println("Optimal primal solution")
println(" x = $x")
elseif ( solsta == MSK_SOL_STA_DUAL_INFEAS_CER ||
solsta == MSK_SOL_STA_PRIM_INFEAS_CER)
println("Primal or dual infeasibility.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status.")
else
println("Other solution status")
end
end
end
logistic.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : logistic.jl
#
# Purpose: Implements logistic regression with regulatization.
#
# Demonstrates using the exponential cone and log-sum-exp in Optimizer API.
using Mosek
"""
Adds ACCs for t_i >= log ( 1 + exp((1-2*y[i]) * theta' * X[i]) )
Adds auxiliary variables, AFE rows and constraints
"""
function softplus(task :: Mosek.Task, d::Int, n::Int, theta::Int, t::Int, X::Matrix{Float64}, y::Vector{Bool})
nvar = getnumvar(task)
ncon = getnumcon(task)
nafe = getnumafe(task)
z1 = nvar
z2 = z1+n
zcon = ncon
thetaafe = nafe
tafe = thetaafe+n
z1afe = tafe+n
z2afe = z1afe+n
appendvars(task,2*n) # z1, z2
appendcons(task,n) # z1 + z2 = 1
appendafes(task,4*n) #theta * X[i] - t[i], -t[i], z1[i], z2[i]
# Linear constraints
let subi = Vector{Int32}(undef,2*n),
subj = Vector{Int32}(undef,2*n),
aval = Vector{Float64}(undef,2*n),
k = 1
for i in 1:n
# z1 + z2 = 1
subi[k] = zcon+i; subj[k] = z1+i; aval[k] = 1; k += 1
subi[k] = zcon+i; subj[k] = z2+i; aval[k] = 1; k += 1
end
putaijlist(task,subi, subj, aval)
end
putconboundsliceconst(task,zcon+1, zcon+n+1, MSK_BK_FX, 1, 1)
putvarboundsliceconst(task,nvar+1, nvar+2*n+1, MSK_BK_FR, -Inf, Inf)
# Affine conic expressions
let afeidx = Vector{Int64}(undef,d*n+4*n),
varidx = Vector{Int32}(undef,d*n+4*n),
fval = Vector{Float64}(undef,d*n+4*n),
k = 1
# Thetas
for i in 1:n
for j in 1:d
afeidx[k] = thetaafe + i; varidx[k] = theta + j
fval[k] = if y[i] X[i,j]*(-1.0) else X[i,j] end
k += 1
end
end
# -t[i]
for i in 1:n
afeidx[k] = thetaafe + i; varidx[k] = t + i; fval[k] = -1; k += 1
afeidx[k] = tafe + i; varidx[k] = t + i; fval[k] = -1; k += 1
end
# z1, z2
for i in 1:n
afeidx[k] = z1afe + i; varidx[k] = z1 + i; fval[k] = 1; k += 1
afeidx[k] = z2afe + i; varidx[k] = z2 + i; fval[k] = 1; k += 1
end
# Add the expressions
putafefentrylist(task,afeidx, varidx, fval)
end
# Add a single row with the constant expression "1.0"
let oneafe = getnumafe(task)+1,
# Add an exponential cone domain
expdomain = appendprimalexpconedomain(task)
appendafes(task,1)
putafeg(task,oneafe,1.0)
# Conic constraints
for i in 1:n
appendacc(task,expdomain, [z1afe+i, oneafe, thetaafe+i], nothing)
appendacc(task,expdomain, [z2afe+i, oneafe, tafe+i], nothing)
end
end
end # softplus
"""
Model logistic regression (regularized with full 2-norm of theta)
X - n x d matrix of data points
y - length n vector classifying training points
lamb - regularization parameter
"""
function logisticRegression(X :: Matrix{Float64},
y :: Vector{Bool},
lamb :: Float64)
(n,d) = size(X)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Variables [r; theta; t]
nvar = 1+d+n
appendvars(task,nvar)
putvarboundsliceconst(task,1, nvar+1, MSK_BK_FR, -Inf, Inf)
r = 0
theta = r+1
t = theta+d
# Objective lambda*r + sum(t)
putcj(task,r+1,lamb)
for i in 1:n
putcj(task,t+i, 1.0)
end
# Softplus function constraints
softplus(task, d, n, theta, t, X, y)
# Regularization
# Append a sequence of linear expressions (r, theta) to F
numafe = getnumafe(task)
appendafes(task,1+d)
putafefentry(task,numafe+1, r+1, 1.0)
for i in 1:d
putafefentry(task,numafe+i+1, theta+i, 1.0)
end
# Add the constraint
appendaccseq(task,appendquadraticconedomain(task,1+d), numafe+1, zeros(d+1))
# Solution
optimize(task)
getxxslice(task,MSK_SOL_ITR, theta+1, theta+d+1)
end
end # logisticRegression
# Test: detect and approximate a circle using degree 2 polynomials
let n = 30
X = Matrix{Float64}(undef,(n*n,6))
Y = Vector{Bool}(undef,n*n)
k = 1
for i in 1:n
for j in 1:n
x = -1 + 2.0*i/(n-1)
y = -1 + 2.0*j/(n-1)
X[k,:] = [1.0, x, y, x*y, x*x, y*y]
# X[k,0] = 1.0; X[k,1] = x; X[k,2] = y; X[k,3] = x*y
# X[k,4] = x*x; X[k,5] = y*y
Y[k] = x*x+y*y>=0.69
k += 1
end
end
theta = logisticRegression(X, Y, 0.1)
println("theta = $theta")
end
mico1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : mico1.jl
#
# Purpose : Demonstrates how to solve a small mixed
# integer conic optimization problem.
# minimize x^2 + y^2
# subject to x >= e^y + 3.8
# x, y - integer
using Mosek
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream to the user specified
# method task_msg_obj.stream
appendvars(task,3); # x, y, t
x=1
y=2
t=3
putvarboundsliceconst(task,1, 4, MSK_BK_FR, -Inf, Inf)
# Integrality constraints for x, y
putvartypelist(task,
[x,y],
[MSK_VAR_TYPE_INT,MSK_VAR_TYPE_INT])
# Set up the affine expressions
# x, x-3.8, y, t, 1.0
appendafes(task,5)
putafefentrylist(task,
[1,2,3,4],
[x,x,y,t],
[1,1,1,1])
putafegslice(task,1, 6, Float64[0, -3.8, 0, 0, 1.0])
# Add constraint (x-3.8, 1, y) \in \EXP
appendacc(task,appendprimalexpconedomain(task), Int64[2, 5, 3], nothing)
# Add constraint (t, x, y) \in \QUAD
appendacc(task,appendquadraticconedomain(task,3), Int64[4, 1, 3], nothing)
# Objective
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
putcj(task,t, 1.0)
# Optimize the task
optimize(task)
solutionsummary(task,MSK_STREAM_MSG)
xx = getxxslice(task,MSK_SOL_ITG, 1, 3)
println("x = $(xx[1]) y = $(xx[2])")
end
milo1.jl
##
# File: milo1.jl
#
# Purpose: Demonstrates how to solve a small mixed
# integer linear optimization problem using the MOSEK Python API.
##
using Mosek
using Printf, SparseArrays
# Define a stream printer to grab output from MOSEK
printstream(msg::String) = print(msg)
bkc = [ MSK_BK_UP, MSK_BK_LO ]
blc = [ -Inf, -4.0 ]
buc = [ 250.0, Inf ]
bkx = [ MSK_BK_LO, MSK_BK_LO ]
blx = [ 0.0, 0.0 ]
bux = [ Inf, Inf ]
c = [ 1.0, 0.64 ]
A = sparse( [ 1, 1, 2, 2],
[ 1, 2, 1, 2],
[ 50.0, 31.0,
3.0, -2.0] )
numvar = length(bkx)
numcon = length(bkc)
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Attach a printer to the task
putstreamfunc(task,MSK_STREAM_LOG,printstream)
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
#Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Set the linear term c_j in the objective.
putclist(task,[1:numvar;],c)
# Set the bounds on variables
# blx[j] <= x_j <= bux[j]
putvarboundslice(task,1,numvar+1,bkx,blx,bux)
# Input columns of A
putacolslice(task,1,numvar+1, A.colptr[1:numvar],A.colptr[2:numvar+1],A.rowval,A.nzval)
putconboundslice(task,1,numcon+1,bkc,blc,buc)
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Define variables to be integers
putvartypelist(task,[ 1, 2 ],
[ MSK_VAR_TYPE_INT, MSK_VAR_TYPE_INT ])
# Set max solution time
putdouparam(task,MSK_DPAR_MIO_MAX_TIME, 60.0)
# Optimize the task
optimize(task)
writedata(task,"milo1.ptf")
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITG)
solsta = getsolsta(task,MSK_SOL_ITG)
if solsta == MSK_SOL_STA_INTEGER_OPTIMAL
# Output a solution
xx = getxx(task,MSK_SOL_ITG)
@printf("Optimal solution: %s\n", xx')
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
mioinitsol.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : mioinitsol.jl
#
# Purpose : Demonstrates how to solve a MIP with a start guess.
using Mosek
let numvar = 4,
numcon = 1,
c = [ 7.0, 10.0, 1.0, 5.0 ],
bkc = [MSK_BK_UP],
blc = [ -Inf],
buc = [2.5],
bkx = [ MSK_BK_LO,
MSK_BK_LO,
MSK_BK_LO,
MSK_BK_LO],
blx = [ 0.0,
0.0,
0.0,
0.0 ],
bux = [ Inf,
Inf,
Inf,
Inf ],
ptrb = Int64[1, 2, 3, 4],
ptre = Int64[2, 3, 4, 5],
aval = Float64[1.0, 1.0, 1.0, 1.0],
asub = Int32[1, 1, 1, 1 ],
intsub = Int32[1, 2, 3],
inttype = [ MSK_VAR_TYPE_INT,
MSK_VAR_TYPE_INT,
MSK_VAR_TYPE_INT ]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
inputdata(task,numcon, numvar,
c, 0.0,
ptrb, ptre,
asub, aval,
bkc, blc, buc,
bkx, blx, bux)
putvartypelist(task,intsub, inttype)
# A maximization problem
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Assign values to integer variables
# We only set that slice of xx
putxxslice(task,MSK_SOL_ITG, 1, 4, [1.0, 1.0, 0.0])
# Request constructing the solution from integer variable values
putintparam(task,MSK_IPAR_MIO_CONSTRUCT_SOL, MSK_ON)
# solve
optimize(task)
writedata(task,"mioinitsol.ptf")
solutionsummary(task,MSK_STREAM_LOG)
# Read and print solution
if solutiondef(task,MSK_SOL_ITG)
# Output a solution
xx = getxx(task,MSK_SOL_ITG)
println("Integer optimal solution:")
for i in 1:numvar
println(" x[$i] = $(xx[i])")
end
# Was the initial guess used?
constr = getintinf(task,MSK_IINF_MIO_CONSTRUCT_SOLUTION)
constrVal = getdouinf(task,MSK_DINF_MIO_CONSTRUCT_SOLUTION_OBJ)
println("Construct solution utilization: $constr")
println("Construct solution objective: $constrVal")
else
println("No integer solution is available.")
end
# Was the initial solution used?
constr = getintinf(task,MSK_IINF_MIO_CONSTRUCT_SOLUTION)
constrVal = getdouinf(task,MSK_DINF_MIO_CONSTRUCT_SOLUTION_OBJ)
println("Construct solution utilization: $constr")
println("Construct solution objective: $constrVal")
end
end
opt_server_async.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : opt_server_async.jl
#
# Purpose : Demonstrates how to use MOSEK OptServer
# to solve optimization problem asynchronously
##
using Mosek
if length(ARGS) < 2
println("Missing argument, syntax is:")
println(" opt_server_async inputfile host:port numpolls [cert]")
else
filename = ARGS[1]
serveraddr = ARGS[2]
numpolls = parse(Int,ARGS[3])
cert = (if length(ARGS) > 3
ARGS[4]
else
Nothing
end)
token = maketask() do task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
token = Nothing
print("reading task from file")
readdata(task,filename)
if cert !== Nothing
putstrparam(task,MSK_SPAR_REMOTE_TLS_CERT_PATH,cert)
end
println("Solve the problem remotely (async)")
asyncoptimize(task,serveraddr,"")
end
println("Task token: '$token'")
maketask() do task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
readdata(task,filename)
if cert !== Nothing
putstrparam(task,MSK_SPAR_REMOTE_TLS_CERT_PATH,cert)
end
i = 0
while i < numpolls
sleep(0.1)
println("poll $i...")
respavailable, res, trm = asyncpoll(task, serveraddr, "", token)
println("done!")
if respavailable
println("solution available!")
respavailable, res, trm = asyncgetresult(task, serveraddr, "", token)
solutionsummary(task,MSK_STREAM_LOG)
break
end
i = i + 1
if i == numpolls
println("max number of polls reached, stopping host.")
asyncstop(task,serveraddr,"", token)
end
end
end
end
opt_server_sync.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : opt_server_sync.jl
#
# Purpose : Demonstrates how to use MOSEK OptServer
# to solve optimization problem synchronously
##
using Mosek
if length(ARGS) < 2
println("Missing argument, syntax is:")
println(" opt_server_sync inputfile addr [certpath]")
else
maketask() do task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
inputfile = ARGS[1]
serveraddr = ARGS[2]
tlscert = if length(ARGS) < 3 Nothing else ARGS[3] end
# We assume that a problem file was given as the first command
# line argument (received in `argv')
readdata(task,inputfile)
# Set OptServer URL
putoptserverhost(task,serveraddr)
# Path to certificate, if any
if tlscert !== Nothing
putstrparam(task,MSK_SPAR_REMOTE_TLS_CERT_PATH, tlscert)
# Solve the problem remotely, no access token
trm = optimize(task)
# Print a summary of the solution
solutionsummary(task,MSK_STREAM_LOG)
end
end
end
parallel.jl
# Copyright: Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File: parallel.jl
#
# Purpose: Demonstrates parallel optimization using optimizebatch()
using Mosek
# Example of how to use env.optimizebatch().
# Optimizes tasks whose names were read from command line.
if length(ARGS) < 2
println("Usage: parallel FILENAME FILENAME [ FILENAME ... ]")
else
n = length(ARGS)
makeenv() do env
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
tasks = [ maketask(filename=f) for f in ARGS ]
# Size of thread pool available for all tasks
threadpoolsize = 6
for t in tasks
putintparam(t,MSK_IPAR_NUM_THREADS, 2)
end
# Optimize all the given tasks in parallel
(trm,res) = optimizebatch(env,
false, # No race
-1.0, # No time limit
threadpoolsize,
tasks) # Array of tasks to optimize
for (i,t) in enumerate(tasks)
println("Task $i res $(res[i]) trm $(trm[i]) obj_val $(getdouinf(t,MSK_DINF_INTPNT_PRIMAL_OBJ)) time $(getdouinf(t,MSK_DINF_OPTIMIZER_TIME))")
end
end
end
parameters.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : parameters.jl
#
# Purpose : Demonstrates a very simple example about how to get/set
# parameters with MOSEK Julia API
#
using Mosek
maketask() do task
println("Test MOSEK parameter get/set functions");
# Set log level (integer parameter)
putintparam(task,MSK_IPAR_LOG, 1)
# Select interior-point optimizer... (integer parameter)
putintparam(task,MSK_IPAR_OPTIMIZER, MSK_OPTIMIZER_INTPNT)
# ... without basis identification (integer parameter)
putintparam(task,MSK_IPAR_INTPNT_BASIS,MSK_BI_NEVER)
# Set relative gap tolerance (double parameter)
putdouparam(task,MSK_DPAR_INTPNT_CO_TOL_REL_GAP, 1.0e-7)
# The same using explicit string names
putparam(task,"MSK_DPAR_INTPNT_CO_TOL_REL_GAP", "1.0e-7");
putnadouparam(task,"MSK_DPAR_INTPNT_CO_TOL_REL_GAP", 1.0e-7 )
# Incorrect value
try
putdouparam(task,MSK_DPAR_INTPNT_CO_TOL_REL_GAP, -1.0)
catch
println("Wrong parameter value")
end
param = getdouparam(task,MSK_DPAR_INTPNT_CO_TOL_REL_GAP)
println("Current value for parameter intpnt_co_tol_rel_gap = $param")
# Define and solve an optimization problem here
# optimize(task,)
# After optimization:
println("Get MOSEK information items")
tm = getdouinf(task,MSK_DINF_OPTIMIZER_TIME)
iter = getintinf(task,MSK_IINF_INTPNT_ITER)
println("Time: $tm");
println("Iterations: $iter");
end
pinfeas.jl
# File : pinfeas.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# Purpose: Demonstrates how to fetch a primal infeasibility certificate
# for a linear problem
#
using Mosek
# Set up a simple linear problem from the manual for test purposes
function testProblem(func :: Function)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
appendvars(task,7)
appendcons(task,7);
putclist(task,
Int32[1,2,3,4,5,6,7],
Float64[1,2,5,2,1,2,1])
putaijlist(task,
Int32[1,1,2,2,3,3,3,4,4,5,6,6,7,7],
Int32[1,2,3,4,5,6,7,1,5,2,3,6,4,7],
Float64[1,1,1,1,1,1,1,1,1,1,1,1,1,1])
putconboundslice(task,
1, 8,
[MSK_BK_UP,MSK_BK_UP,MSK_BK_UP,MSK_BK_FX,MSK_BK_FX,MSK_BK_FX,MSK_BK_FX],
Float64[-Inf, -Inf, -Inf, 1100, 200, 500, 500],
Float64[200, 1000, 1000, 1100, 200, 500, 500])
putvarboundsliceconst(task,1, 8, MSK_BK_LO, 0.0, +Inf)
func(task)
end
end
"""
Analyzes and prints infeasibility contributing elements
sl - dual values for lower bounds
su - dual values for upper bounds
eps - tolerance for when a nunzero dual value is significant
"""
function analyzeCertificate(sl :: Vector{Float64}, su :: Vector{Float64}, eps :: Float64)
for i in 1:length(sl)
if abs(sl[i]) > eps
println("#$i, lower, dual = $(sl[i])")
end
if abs(su[i]) > eps
println("#$i, upper, dual = $(su[i])")
end
end
end
# In this example we set up a simple problem
# One could use any task or a task read from a file
testProblem() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Useful for debugging
writedata(task,"pinfeas.ptf"); # Write file in human-readable format
# Attach a log stream printer to the task
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Perform the optimization.
optimize(task)
solutionsummary(task,MSK_STREAM_LOG)
# Check problem status, we use the interior point solution
if getprosta(task,MSK_SOL_ITR) == MSK_PRO_STA_PRIM_INFEAS
# Set the tolerance at which we consider a dual value as essential
eps = 1e-7
println("Variable bounds important for infeasibility: ");
analyzeCertificate(getslx(task,MSK_SOL_ITR), getsux(task,MSK_SOL_ITR), eps)
println("Constraint bounds important for infeasibility: ")
analyzeCertificate(getslc(task,MSK_SOL_ITR), getsuc(task,MSK_SOL_ITR), eps)
else
println("The problem is not primal infeasible, no certificate to show")
end
end
portfolio_1_basic.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_1_basic.jl
#
# Description : Implements a basic portfolio optimization model.
using Mosek
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
gamma :: Float64,
GT :: Array{Float64,2})
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
# Constraints offsets
budget_ofs = 0
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,n)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]");
# No short-selling - x^l = 0, x^u = inf
putvarbound(task,x_ofs+j, MSK_BK_LO, 0.0, Inf);
end
# One linear constraint: total budget
appendcons(task,1);
putconname(task,1,"budget");
for j in 1:n
# Coefficients in the first row of A
putaij(task,budget_ofs+1, x_ofs+j, 1.0)
end
putconbound(task, budget_ofs+1, MSK_BK_FX, totalBudget, totalBudget)
# Input (gamma, GTx) in the AFE (affine expression) storage
# We need k+1 rows
appendafes(task,k + 1)
# The first affine expression = gamma
putafeg(task,1, gamma)
# The remaining k expressions comprise GT*x, we add them row by row
# In more realisic scenarios it would be better to extract nonzeros and input in sparse form
subj = [1:n...]
for i in 1:k
putafefrow(task,i + 1, subj, GT[i,:])
end
# Input the affine conic constraint (gamma, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
qdom = appendquadraticconedomain(task,k + 1)
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,1, "risk")
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
optimize(task)
# Display solution summary for quick inspection of results
solutionsummary(task,MSK_STREAM_LOG)
writedata(task,"portfolio_1_basic.ptf");
# Check if the interior point solution is an optimal point
if getsolsta(task, MSK_SOL_ITR) != MSK_SOL_STA_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
# Read the results
xx = getxxslice(task,MSK_SOL_ITR, x_ofs+1,x_ofs+n+1)
expret = mu' * xx
(xx,expret)
end
end # portfolio()
let w = 59.0,
mu = [0.07197349, 0.15518171, 0.17535435, 0.0898094 , 0.42895777, 0.39291844, 0.32170722, 0.18378628],
x0 = [8.0, 5.0, 3.0, 5.0, 2.0, 9.0, 3.0, 6.0],
gamma = 36.0,
GT = [ 0.30758 0.12146 0.11341 0.11327 0.17625 0.11973 0.10435 0.10638
0. 0.25042 0.09946 0.09164 0.06692 0.08706 0.09173 0.08506
0. 0. 0.19914 0.05867 0.06453 0.07367 0.06468 0.01914
0. 0. 0. 0.20876 0.04933 0.03651 0.09381 0.07742
0. 0. 0. 0. 0.36096 0.12574 0.10157 0.0571
0. 0. 0. 0. 0. 0.21552 0.05663 0.06187
0. 0. 0. 0. 0. 0. 0.22514 0.03327
0. 0. 0. 0. 0. 0. 0. 0.2202 ],
(k,n) = size(GT)
let (xx,expret) = portfolio(mu,x0,w,gamma,GT)
println("Expected return $(expret) for gamma $(gamma)")
println("Solution vector = $(xx)")
end
end
portfolio_2_frontier.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_2_frontier.jl
#
# Description : Implements a basic portfolio optimization model.
using Mosek
using Printf
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
alphas :: Vector{Float64},
GT :: Array{Float64,2})
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
#putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
s_ofs = n
# Constraints offsets
budget_ofs = 0
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,n+1)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]")
# No short-selling - x^l = 0, x^u = inf
putvarbound(task,x_ofs+j, MSK_BK_LO, 0.0, Inf)
end
putvarname(task, s_ofs+1, "s")
putvarbound(task,s_ofs+1, MSK_BK_FR, -Inf, Inf)
# One linear constraint: total budget
appendcons(task,1)
putconname(task,1,"budget")
for j in 1:n
# Coefficients in the first row of A
putaij(task,budget_ofs+1, x_ofs+j, 1.0)
end
putconbound(task, budget_ofs+1, MSK_BK_FX, totalBudget, totalBudget)
# Input (gamma, GTx) in the AFE (affine expression) storage
# We build the following F and g for variables [x, s]:
# [0, 1] [0 ]
# F = [0, 0], g = [0.5]
# [GT,0] [0 ]
# We need k+2 rows
appendafes(task,k + 2)
# The first affine expression = alpha
putafefentry(task,1,s_ofs+1,1.0)
putafeg(task,2,0.5)
# The remaining k expressions comprise GT*x, we add them row by row
# In more realisic scenarios it would be better to extract nonzeros and input in sparse form
subj = [1:n...]
for i in 1:k
putafefrow(task,i + 2, subj, GT[i,:])
end
# Input the affine conic constraint (alpha, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
qdom = appendrquadraticconedomain(task,k + 2)
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,1, "risk")
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
expret = Float64[]
stddev = Float64[]
for alpha in alphas
putcj(task,s_ofs+1,-alpha)
optimize(task)
writedata(task,"portfolio_2_frontier-$alpha.ptf")
# Check if the interior point solution is an optimal point
if getsolsta(task, MSK_SOL_ITR) != MSK_SOL_STA_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
# Display solution summary for quick inspection of results
solutionsummary(task,MSK_STREAM_LOG)
# Read the results
r = mu' * getxxslice(task,MSK_SOL_ITR, x_ofs+1,x_ofs+n+1)
slvl = getxxslice(task,MSK_SOL_ITR, s_ofs+1,s_ofs+2)
push!(expret,r)
push!(stddev,sqrt(slvl[1]))
end
(expret,stddev)
end
end # portfolio()
let w = 1.0,
mu = [0.07197, 0.15518, 0.17535, 0.08981, 0.42896, 0.39292, 0.32171, 0.18379],
x0 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
GT = [0.30758 0.12146 0.11341 0.11327 0.17625 0.11973 0.10435 0.10638
0. 0.25042 0.09946 0.09164 0.06692 0.08706 0.09173 0.08506
0. 0. 0.19914 0.05867 0.06453 0.07367 0.06468 0.01914
0. 0. 0. 0.20876 0.04933 0.03651 0.09381 0.07742
0. 0. 0. 0. 0.36096 0.12574 0.10157 0.0571
0. 0. 0. 0. 0. 0.21552 0.05663 0.06187
0. 0. 0. 0. 0. 0. 0.22514 0.03327
0. 0. 0. 0. 0. 0. 0. 0.2202 ],
alphas = [0.0, 0.01, 0.1, 0.25, 0.30, 0.35, 0.4, 0.45, 0.5, 0.75, 1.0, 1.5, 2.0, 3.0, 10.0],
(k,n) = size(GT),
N = length(alphas)
let (expret,stddev) = portfolio(mu,x0,w,alphas,GT)
for i in 1:N
@printf("alpha = %2.2e, exp. ret. = %2.3e, std. dev. = %2.3e\n",alphas[i],expret[i],stddev[i])
end
end
end
portfolio_3_impact.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_3_impact.jl
#
# Description : Implements a basic portfolio optimization model.
##
using Mosek
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
gamma :: Float64,
GT :: Array{Float64,2},
m :: Vector{Float64})
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
c_ofs = n
z_ofs = 2*n
# Constraints offsets
budget_ofs = 0
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,3*n)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]");
putvarname(task, c_ofs+j, "c[$(j)]");
putvarname(task, z_ofs+j, "z[$(j)]");
end
# No short-selling - x^l = 0, x^u = inf
putvarboundsliceconst(task,x_ofs+1,x_ofs+n+1, MSK_BK_LO, 0.0, Inf);
putvarboundsliceconst(task,c_ofs+1,c_ofs+n+1, MSK_BK_FR, -Inf, Inf);
putvarboundsliceconst(task,z_ofs+1,z_ofs+n+1, MSK_BK_FR, -Inf, Inf);
# Linear constraint: total budget
let budget_ofs = getnumcon(task)
appendcons(task,1);
putconname(task,1,"budget");
for j in 1:n
# Coefficients in the first row of A
putaij(task,budget_ofs+1, x_ofs+j, 1.0)
putaij(task,budget_ofs+1, c_ofs+j, m[j])
end
end
putconbound(task, budget_ofs+1, MSK_BK_FX, totalBudget, totalBudget)
# - Absolute value
let zabs1_ofs = getnumcon(task),
zabs2_ofs = zabs1_ofs+n
appendcons(task,2*n)
for i in 1:n
putconname(task,zabs1_ofs + i, "zabs1[$i]")
putaij(task,zabs1_ofs + i, x_ofs + i, -1.0)
putaij(task,zabs1_ofs + i, z_ofs + i, 1.0)
putconbound(task,zabs1_ofs + i, MSK_BK_LO, -x0[i], Inf)
putconname(task,zabs2_ofs + i, "zabs2[$i]")
putaij(task,zabs2_ofs + i, x_ofs + i, 1.0)
putaij(task,zabs2_ofs + i, z_ofs + i, 1.0)
putconbound(task,zabs2_ofs + i, MSK_BK_LO, x0[i], Inf)
end
end
let qdom = appendquadraticconedomain(task,k + 1)
afei = getnumafe(task)
acci = getnumacc(task)
# Input (gamma, GTx) in the AFE (affine expression) storage
# We need k+1 rows
appendafes(task,k + 1)
# The first affine expression = gamma
putafeg(task,afei+1, gamma)
# The remaining k expressions comprise GT*x, we add them row by row
# In more realisic scenarios it would be better to extract nonzeros and input in sparse form
subj = [1:n...]
for i in 1:k
putafefrow(task,afei+i+1, subj, GT[i,:])
end
# Input the affine conic constraint (gamma, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,acci+1, "risk")
end
let dom = appendprimalpowerconedomain(task,3,[2.0,1.0])
afei = getnumafe(task)
acci = getnumacc(task)
afe0 = afei
appendafes(task,2*n+1)
putafeg(task,afe0+1,1.0)
afei += 1
for i in 1:n
putafefentry(task,afei + i, c_ofs + i, 1.0);
putafefentry(task,afei + n + i, z_ofs + i, 1.0);
end
accafes = Int64[ k for i in 1:n for k in [ afei + i, afe0+1, afei + n + i ] ]
accb = zeros(n*3)
accdom = Int64[ dom for i in 1:n ]
appendaccs(task,accdom,accafes,accb)
for i in 1:n
putaccname(task,acci+i,"market_impact[$i]")
end
end
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
optimize(task)
# Check if the interior point solution is an optimal point
if getsolsta(task, MSK_SOL_ITR) != MSK_SOL_STA_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
# Display solution summary for quick inspection of results
solutionsummary(task,MSK_STREAM_LOG)
writedata(task,"portfolio_3_impact.ptf");
# Read the results
xx = getxxslice(task,MSK_SOL_ITR, x_ofs+1,x_ofs+n+1)
expret = mu' * xx
(xx,expret)
end
end # portfolio()
let w = 1.0,
mu = [0.07197, 0.15518, 0.17535, 0.08981, 0.42896, 0.39292, 0.32171, 0.18379],
x0 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
GT = [0.30758 0.12146 0.11341 0.11327 0.17625 0.11973 0.10435 0.10638
0. 0.25042 0.09946 0.09164 0.06692 0.08706 0.09173 0.08506
0. 0. 0.19914 0.05867 0.06453 0.07367 0.06468 0.01914
0. 0. 0. 0.20876 0.04933 0.03651 0.09381 0.07742
0. 0. 0. 0. 0.36096 0.12574 0.10157 0.05710
0. 0. 0. 0. 0. 0.21552 0.05663 0.06187
0. 0. 0. 0. 0. 0. 0.22514 0.03327
0. 0. 0. 0. 0. 0. 0. 0.22020],
gamma = 0.36,
(k,n) = size(GT),
market_impact = 0.01 * ones(n)
let (xx,expret) = portfolio(mu,x0,w,gamma,GT,market_impact)
println("Solution x = $(xx)")
println("Expected return $(expret) for gamma $(gamma)")
end
end
portfolio_4_transcost.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_4_transcost.jl
#
# Description : Implements a basic portfolio optimization model
# with fixed setup costs and transaction costs
# as a mixed-integer problem.
using Mosek
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
gamma :: Float64,
GT :: Array{Float64,2},
f :: Vector{Float64},
g :: Vector{Float64})
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
y_ofs = n
z_ofs = 2*n
# Constraints offsets
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,3*n)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]");
putvarname(task, z_ofs+j, "z[$(j)]");
putvarname(task, y_ofs+j, "y[$(j)]");
# No short-selling - x^l = 0, x^u = inf
putvartype(task, y_ofs+j, MSK_VAR_TYPE_INT);
end
putvarboundsliceconst(task,x_ofs+1,x_ofs+n+1, MSK_BK_LO, 0.0, Inf);
putvarboundsliceconst(task,y_ofs+1,y_ofs+n+1, MSK_BK_RA, 0.0, 1.0);
putvarboundsliceconst(task,z_ofs+1,z_ofs+n+1, MSK_BK_FR, -Inf, Inf);
# One linear constraint: total budget
let coni = getnumcon(task)
appendcons(task,1);
putconname(task,coni+1,"budget")
consub = Int32[ coni+1 for i in 1:n ]
putaijlist(task,consub,[x_ofs+1:x_ofs+n...], ones(n))
putaijlist(task,consub,[z_ofs+1:z_ofs+n...],g)
putaijlist(task,consub,[y_ofs+1:y_ofs+n...],f)
putconbound(task, coni+1, MSK_BK_FX, totalBudget, totalBudget)
end
let coni = getnumcon(task)
appendcons(task,2*n)
for i in 1:n
putconname(task,coni+i,"zabs1[$i]")
putconname(task,coni+n+i,"zabs2[$i]")
putaij(task,coni + i, x_ofs + i, -1.0)
putaij(task,coni + i, z_ofs + i, 1.0)
putconbound(task,coni + i, MSK_BK_LO, -x0[i], Inf)
putaij(task,coni + n + i, x_ofs + i, 1.0)
putaij(task,coni + n + i, z_ofs + i, 1.0)
putconbound(task,coni + n + i, MSK_BK_LO, x0[i], Inf)
end
end
# - Switch
let coni = getnumcon(task)
appendcons(task,n)
for i in 1:n
putconname(task,coni+i, "switch[$i]")
putaij(task,coni + i, z_ofs + i, 1.0)
putaij(task,coni + i, y_ofs + i, -totalBudget)
putconbound(task,coni + i, MSK_BK_UP, -Inf, 0.0)
end
end
let afei = getnumafe(task),
acci = getnumacc(task)
# Input (gamma, GTx) in the AFE (affine expression) storage
# We need k+1 rows
appendafes(task,k + 1)
# The first affine expression = gamma
putafeg(task,afei+1, gamma)
# The remaining k expressions comprise GT*x, we add them row by row
# In more realisic scenarios it would be better to extract nonzeros and input in sparse form
subj = [1:n...]
for i in 1:k
putafefrow(task,afei+i+1, subj, GT[i,:])
end
# Input the affine conic constraint (gamma, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
qdom = appendquadraticconedomain(task,k + 1)
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,acci+1, "risk")
end
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
optimize(task)
# Check if the integer solution is an optimal point
if getsolsta(task, MSK_SOL_ITG) != MSK_SOL_STA_INTEGER_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
writedata(task,"portfolio_4_transcost.ptf")
# Display solution summary for quick inspection of results
solutionsummary(task,MSK_STREAM_LOG)
writedata(task,"portfolio_4_transcost.ptf");
# Read the results
xx = getxxslice(task,MSK_SOL_ITG, x_ofs+1,x_ofs+n+1)
expret = mu' * xx
(xx,expret)
end
end # portfolio()
let w = 1.0,
mu = [0.07197, 0.15518, 0.17535, 0.08981, 0.42896, 0.39292, 0.32171, 0.18379],
x0 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
GT = [ 0.30758 0.12146 0.11341 0.11327 0.17625 0.11973 0.10435 0.10638
0. 0.25042 0.09946 0.09164 0.06692 0.08706 0.09173 0.08506
0. 0. 0.19914 0.05867 0.06453 0.07367 0.06468 0.01914
0. 0. 0. 0.20876 0.04933 0.03651 0.09381 0.07742
0. 0. 0. 0. 0.36096 0.12574 0.10157 0.0571
0. 0. 0. 0. 0. 0.21552 0.05663 0.06187
0. 0. 0. 0. 0. 0. 0.22514 0.03327
0. 0. 0. 0. 0. 0. 0. 0.2202 ],
(k,n) = size(GT),
f = 0.01 * ones(n),
g = 0.001 * ones(n),
gamma = 0.36
let (xx,expret) = portfolio(mu,x0,w,gamma,GT,f,g)
println("Expected return $(expret) for gamma $(gamma)")
println("Solution vector = $(xx)")
end
end
portfolio_5_card.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_5_card.jl
#
# Description : Implements a basic portfolio optimization model.
using Mosek
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
gamma :: Float64,
GT :: Array{Float64,2},
K :: Int)
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
# putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
y_ofs = n
z_ofs = 2*n
# Constraints offsets
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,3*n)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]");
putvarname(task, y_ofs+j, "y[$(j)]");
putvarname(task, z_ofs+j, "z[$(j)]");
putvartype(task, y_ofs+j, MSK_VAR_TYPE_INT);
end
# No short-selling - x^l = 0, x^u = inf
putvarboundsliceconst(task,x_ofs+1,x_ofs+n+1, MSK_BK_LO, 0.0, Inf);
putvarboundsliceconst(task,y_ofs+1,y_ofs+n+1, MSK_BK_RA, 0.0, 1.0);
putvarboundsliceconst(task,z_ofs+1,z_ofs+n+1, MSK_BK_FR, -Inf, Inf);
# One linear constraint: total budget
let coni = getnumcon(task)
appendcons(task,1);
putconname(task,coni+1,"budget")
consub = Int32[ coni+1 for i in 1:n ]
putaijlist(task,consub,[x_ofs+1:x_ofs+n...], ones(n))
putconbound(task, coni+1, MSK_BK_FX, totalBudget, totalBudget)
end
let coni = getnumcon(task)
appendcons(task,1)
putarow(task,coni+1,[y_ofs+1:y_ofs+1+n...],ones(n));
putconbound(task,coni+1,MSK_BK_UP, 0.0, K)
putconname(task,coni+1,"cardinality")
end
let coni = getnumcon(task)
appendcons(task,2*n)
for i in 1:n
putconname(task,coni+i,"zabs1[$i]")
putconname(task,coni+n+i,"zabs2[$i]")
putaij(task,coni + i, x_ofs + i, -1.0)
putaij(task,coni + i, z_ofs + i, 1.0)
putconbound(task,coni + i, MSK_BK_LO, -x0[i], Inf)
putaij(task,coni + n + i, x_ofs + i, 1.0)
putaij(task,coni + n + i, z_ofs + i, 1.0)
putconbound(task,coni + n + i, MSK_BK_LO, x0[i], Inf)
end
end
# - Switch
let coni = getnumcon(task)
appendcons(task,n)
for i in 1:n
putconname(task,coni+i, "switch[$i]")
putaij(task,coni + i, z_ofs + i, 1.0)
putaij(task,coni + i, y_ofs + i, -totalBudget)
putconbound(task,coni + i, MSK_BK_UP, -Inf, 0.0)
end
end
let afei = getnumafe(task),
acci = getnumacc(task)
# Input (gamma, GTx) in the AFE (affine expression) storage
# We need k+1 rows
appendafes(task,k + 1)
# The first affine expression = gamma
putafeg(task,afei+1, gamma)
# The remaining k expressions comprise GT*x, we add them row by row
# In more realisic scenarios it would be better to extract nonzeros and input in sparse form
subj = [1:n...]
for i in 1:k
putafefrow(task,afei+i+1, subj, GT[i,:])
end
# Input the affine conic constraint (gamma, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
qdom = appendquadraticconedomain(task,k + 1)
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,acci+1, "risk")
end
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
optimize(task)
if getsolsta(task, MSK_SOL_ITG) != MSK_SOL_STA_INTEGER_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
# Display solution summary for quick inspection of results
solutionsummary(task,MSK_STREAM_LOG)
writedata(task,"portfolio_5_card-$K.ptf");
# Read the results
xx = getxxslice(task,MSK_SOL_ITG, x_ofs+1,x_ofs+n+1)
expret = mu' * xx
(xx,expret)
end
end # portfolio()
let w = 1.0,
mu = [0.07197, 0.15518, 0.17535, 0.08981, 0.42896, 0.39292, 0.32171, 0.18379],
x0 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
GT = [ 0.30758 0.12146 0.11341 0.11327 0.17625 0.11973 0.10435 0.10638
0. 0.25042 0.09946 0.09164 0.06692 0.08706 0.09173 0.08506
0. 0. 0.19914 0.05867 0.06453 0.07367 0.06468 0.01914
0. 0. 0. 0.20876 0.04933 0.03651 0.09381 0.07742
0. 0. 0. 0. 0.36096 0.12574 0.10157 0.0571
0. 0. 0. 0. 0. 0.21552 0.05663 0.06187
0. 0. 0. 0. 0. 0. 0.22514 0.03327
0. 0. 0. 0. 0. 0. 0. 0.2202 ],
(k,n) = size(GT),
gamma = 0.25
for K in 1:n
res = Float64[]
let (xx,expret) = portfolio(mu,x0,w,gamma,GT,K)
println("cardinality = $K")
println("\tExpected return $(expret) for gamma $(gamma)")
println("\tSolution vector = $(xx)")
end
end
end
portfolio_6_factor.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : portfolio_6_factor.jl
#
# Description : Implements a portfolio optimization model using factor model
#
using Mosek
using LinearAlgebra
function portfolio( mu :: Vector{Float64},
x0 :: Vector{Float64},
w :: Float64,
gammas :: Vector{Float64},
theta :: Vector{Float64},
GT :: Array{Float64,2})
(k,n) = size(GT)
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Directs the log task stream
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
totalBudget = sum(x0)+w
#Offset of variables into the API variable.
x_ofs = 0
# Constraints offsets
budget_ofs = 0
# Holding variable x of length n
# No other auxiliary variables are needed in this formulation
appendvars(task,n)
# Setting up variable x
for j in 1:n
# Optionally we can give the variables names
putvarname(task, x_ofs+j, "x[$(j)]");
# No short-selling - x^l = 0, x^u = inf
putvarbound(task,x_ofs+j, MSK_BK_LO, 0.0, Inf);
end
# One linear constraint: total budget
appendcons(task,1);
putconname(task,1,"budget");
for j in 1:n
# Coefficients in the first row of A
putaij(task,budget_ofs+1, x_ofs+j, 1.0)
end
putconbound(task, budget_ofs+1, MSK_BK_FX, totalBudget, totalBudget)
# Input (gamma, G_factor_T x, diag(sqrt(theta))*x) in the AFE (affine expression) storage
# We need k+n+1 rows and we fill them in in three parts
appendafes(task,n+k+1)
# 1. The first affine expression = gamma, will be specified later
# 2. The next k expressions comprise G_factor_T*x, we add them row by row
# transposing the matrix G_factor on the fly
subj = [1:n...]
for i in 1:k
putafefrow(task,i + 1, subj, GT[i,:])
end
# 3. The remaining n rows contain sqrt(theta) on the diagonal
for i in 1:n
putafefentry(task,k + 1 + i, subj[i], sqrt(theta[i]))
end
# Input the affine conic constraint (gamma, GT*x) \in QCone
# Add the quadratic domain of dimension k+1
qdom = appendquadraticconedomain(task,n + k + 1)
# Add the constraint
appendaccseq(task,qdom,1,nothing)
putaccname(task,1, "risk")
# Objective: maximize expected return mu^T x
putclist(task,[x_ofs+1:x_ofs+n...],mu)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
res = Tuple{Float64,Float64}[]
for gamma in gammas
putafeg(task,1, gamma)
optimize(task)
if getsolsta(task, MSK_SOL_ITR) != MSK_SOL_STA_OPTIMAL
# See https://docs.mosek.com/latest/juliaapi/accessing-solution.html about handling solution statuses.
error("Solution not optimal")
end
writedata(task,"portfolio_6_factor-$(gamma).ptf");
# Read the results
xx = getxxslice(task,MSK_SOL_ITR, x_ofs+1,x_ofs+n+1)
expret = mu' * xx
push!(res,(gamma,expret))
end
res
end
end # portfolio()
let w = 1.0,
mu = [0.07197, 0.15518, 0.17535, 0.08981, 0.42896, 0.39292, 0.32171, 0.18379],
x0 = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
gammas = [0.24, 0.28, 0.32, 0.36, 0.4, 0.44, 0.48],
# Factor exposure matrix, n x 2
B = [ 0.4256 0.1869
0.2413 0.3877
0.2235 0.3697
0.1503 0.4612
1.5325 -0.2633
1.2741 -0.2613
0.6939 0.2372
0.5425 0.2116 ],
# Factor covariance matrix, 2x2
S_F = [ 0.0620 0.0577
0.0577 0.0908 ],
# Specific risk components
theta = [0.0720, 0.0508, 0.0377, 0.0394, 0.0663, 0.0224, 0.0417, 0.0459],
S_sqrt_theta = Matrix(Diagonal(sqrt.(theta))),
P = Matrix(cholesky(S_F).L),
BP = B * P,
G = [ BP S_sqrt_theta ],
GT = Matrix(G') # 10 x 8
for (gamma,expret) in portfolio(mu,x0,w,gammas,theta,GT)
println("Expected return $expret for gamma $gamma");
end
print(P)
end
pow1.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : pow1.jl
#
# Purpose: Demonstrates how to solve the problem
#
# maximize x^0.2*y^0.8 + z^0.4 - x
# st x + y + 0.5z = 2
# x,y,z >= 0
#
using Mosek
printstream(msg::AbstractString) = print(msg)
csub = [ 4, 5, 1 ]
cval = [ 1.0, 1.0, -1.0]
asub = [ 1, 2, 3]
aval = [ 1.0, 1.0, 0.5]
numvar = 5
numcon = 1
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
appendcons(task,numcon)
appendvars(task,numvar)
# Set up the linear part of the problem
putclist(task,csub, cval)
putarow(task,1, asub, aval)
putconbound(task,1, MSK_BK_FX, 2.0, 2.0)
putvarboundsliceconst(task,1, numvar+1,MSK_BK_FR,-Inf,Inf)
# Input the cones
pc1 = appendprimalpowerconedomain(task,3,[0.2, 0.8])
pc2 = appendprimalpowerconedomain(task,3,[4.0, 6.0])
appendafes(task,6)
putafefentrylist(task,
[1, 2, 3, 4, 6], # Rows
[1, 2, 4, 3, 5], #Columns
[1.0, 1.0, 1.0, 1.0, 1.0])
putafeg(task,5,1.0)
# Append the two conic constraints
appendacc(task,
pc1, # Domain
[1, 2, 3], # Rows from F
nothing)
appendacc(task,
pc2, # Domain
[4, 5, 6], # Rows from F
nothing)
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Optimize the task
optimize(task)
writedata(task,"pow1.ptf")
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Output a solution
xx = getxx(task,MSK_SOL_ITR)
println("Optimal solution: $(xx[1:3])")
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility.")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
qcqo1.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : qcqo1.jl
#
# Purpose: Demonstrates how to solve small quadratic and quadratically
# constrained optimization problem using the MOSEK Python API.
##
using Mosek
using Printf
# Since the actual value of Infinity is ignores, we define it solely
# for symbolic purposes:
# Set up and input bounds and linear coefficients
bkc = [ MSK_BK_LO ]
blc = [ 1.0 ]
buc = [ Inf ]
bkx = [ MSK_BK_LO
MSK_BK_LO
MSK_BK_LO ]
blx = [ 0.0, 0.0, 0.0 ]
bux = [ Inf, Inf, Inf ]
c = [ 0.0, -1.0, 0.0 ]
asub = [ 1 ,2, 3 ]
aval = [ 1.0, 1.0, 1.0 ]
numvar = length(bkx)
numcon = length(bkc)
# Create a task
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
#Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
#Optionally add a constant term to the objective.
putcfix(task,0.0)
# Set the linear term c_j in the objective.
putclist(task,[1:numvar;],c)
# Set the bounds on variable j
# blx[j] <= x_j <= bux[j]
putvarboundslice(task,1,numvar+1,bkx,blx,bux)
# Input column j of A
putarow(task,1,asub,aval)
putconbound(task,1,bkc[1],blc[1],buc[1])
# Set up and input quadratic objective
qsubi = [ 1, 2, 3, 3 ]
qsubj = [ 1, 2, 1, 3 ]
qval = [ 2.0, 0.2, -1.0, 2.0 ]
putqobj(task,qsubi,qsubj,qval)
# The lower triangular part of the Q^0
# matrix in the first constraint is specified.
# This corresponds to adding the term
# - x0^2 - x1^2 - 0.1 x2^2 + 0.2 x0 x2
qsubi = [ 1, 2, 3, 3 ]
qsubj = [ 1, 2, 3, 1 ]
qval = [ -2.0, -2.0, -0.2, 0.2 ]
# put Q^0 in constraint with index 0.
putqconk(task,1, qsubi,qsubj, qval)
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
# Optimize the task
optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Output a solution
xx = getxx(task,MSK_SOL_ITR)
@printf("Optimal solution: %s\n", xx')
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
qo1.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : qo1.jl
#
# Purpose: Demonstrate how to solve a quadratic
# optimization problem using the MOSEK Python API.
##
using Mosek
using Printf, SparseArrays
# Define a stream printer to grab output from MOSEK
bkc = [ MSK_BK_LO ]
blc = [ 1.0 ]
buc = [ Inf ]
bkx = [ MSK_BK_LO, MSK_BK_LO, MSK_BK_LO ]
blx = [ 0.0, 0.0, 0.0 ]
bux = [ Inf, Inf, Inf ]
numvar = length(bkx)
numcon = length(bkc)
c = [ 0.0, -1.0, 0.0 ]
A = sparse( [ 1, 1, 1 ],
[ 1, 2, 3 ],
[ 1.0, 1.0, 1.0 ],
numcon, numvar )
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
# Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Set the linear term c_j in the objective.
putclist(task,[1:numvar;],c)
# Set the bounds on variable j
# blx[j] <= x_j <= bux[j]
putvarboundslice(task,1,numvar+1,bkx,blx,bux)
putacolslice(task,1,numvar+1,
A.colptr[1:numvar],A.colptr[2:numvar+1],
A.rowval,A.nzval)
# Set up and input quadratic objective
qsubi = [ 1, 2, 3, 3 ]
qsubj = [ 1, 2, 1, 3 ]
qval = [ 2.0, 0.2, -1.0, 2.0 ]
putqobj(task,qsubi,qsubj,qval)
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
# Optimize
r = optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
xx = getxx(task,MSK_SOL_ITR)
println("Optimal solution:")
println(xx)
elseif solsta in [ MSK_SOL_STA_DUAL_INFEAS_CER,
MSK_SOL_STA_PRIM_INFEAS_CER ]
println("Primal or dual infeasibility certificate found.\n")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
@printf("Other solution status (%d)\n",solsta)
end
end
reoptimization.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : reoptimization.jl
#
# Purpose : Demonstrates how to solve a linear
# optimization problem using the MOSEK API
# and modify and re-optimize the problem.
using Mosek
let numcon = 3,
numvar = 3,
c = [1.5, 2.5, 3.0 ],
bkc = [ MSK_BK_UP,
MSK_BK_UP,
MSK_BK_UP
],
blc = [ -Inf,
-Inf,
-Inf ],
buc = [ 100000.0,
50000.0,
60000.0
],
bkx = [ MSK_BK_LO,
MSK_BK_LO,
MSK_BK_LO ],
blx = [ 0.0, 0.0, 0.0 ],
bux = [ +Inf,
+Inf,
+Inf ],
aptrb = Int64[ 1,4,7 ],
aptre = Int64[ 4,7,10 ],
asub = Int32[ 1, 2, 3,
1, 2, 3,
1, 2, 3 ],
aval = [ 2.0, 3.0, 2.0,
4.0, 2.0, 3.0,
3.0, 3.0, 2.0 ]
maketask() do task # 126
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Append the constraints.
appendcons(task,numcon)
# Append the variables.
appendvars(task,numvar)
# Put C.
for j in 1:numvar
putcj(task,j, c[j])
end
# Put constraint bounds
for i in 1:numcon
putconbound(task,i, bkc[i], blc[i], buc[i])
end
# Put variable bounds.
for j in 1:numvar
putvarbound(task,j, bkx[j], blx[j], bux[j])
end
# Put A.
putacolslice(task,1,numvar+1,aptrb,aptre,asub,aval)
# A maximization problem
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
# Solve the problem
optimize(task)
xx = getxx(task,MSK_SOL_BAS) # Request the basic solution.
println("x = $xx")
#***************** Make a change to the A matrix *****************
putaij(task,1,1, 3.0)
optimize(task);
xx = getxx(task,MSK_SOL_BAS)
println("x = $xx")
#**************** Add a new variable *****************************
# Get index of new variable.
varidx = getnumvar(task)+1
# Append a new variable x_3 to the problem
appendvars(task,1)
numvar += 1
# Set bounds on new varaible
putvarbound(task,
varidx,
MSK_BK_LO,
0.0,
Inf);
# Change objective
putcj(task,varidx, 1.0)
# Put new values in the A matrix
let acolsub = Int32[1,3],
acolval = Float64[4.0, 1.0]
putacol(task,varidx, # column index
acolsub,
acolval)
end
# Change optimizer to simplex free and reoptimize
putintparam(task,MSK_IPAR_OPTIMIZER, MSK_OPTIMIZER_FREE_SIMPLEX)
optimize(task)
xx = getxx(task,MSK_SOL_BAS) # Request the basic solution.
println("x = $xx")
#********************* Add a new constraint **************************
# Get index of new constraint.
conidx = getnumcon(task)+1
# Append a new constraint
appendcons(task,1)
numcon += 1;
# Set bounds on new constraint
putconbound(task,
conidx,
MSK_BK_UP,
-Inf,
30000.0)
# Put new values in the A matrix
let arowsub = [1, 2, 3, 4 ],
arowval = [1.0, 2.0, 1.0, 1.0]
putarow(task,conidx, # row index
arowsub,
arowval)
end
optimize(task)
xx = getxx(task,MSK_SOL_BAS)
println("x = $xx")
#********************* Change constraint bounds *******************
let newbkc = [MSK_BK_UP,
MSK_BK_UP,
MSK_BK_UP,
MSK_BK_UP],
newblc = [ -Inf,
-Inf,
-Inf,
-Inf],
newbuc = [ 80000.0, 40000.0, 50000.0, 22000.0 ]
putconboundslice(task,1, numcon+1, newbkc, newblc, newbuc)
end
optimize(task)
xx = getxx(task,MSK_SOL_BAS) # Request the basic solution.
println("x = $xx")
end
end
response.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : response.jl
#
# Purpose : This example demonstrates proper response handling
# for problems solved with the interior-point optimizers.
#
##
using Mosek
cqo1_ptf = "
Task ''
Objective obj
Minimize + x4 + x5 + x6
Constraints
c1 [1] + x1 + x2 + 2 x3
k1 [QUAD(3)]
@ac1: + x4
@ac2: + x1
@ac3: + x2
k2 [RQUAD(3)]
@ac4: + x5
@ac5: + x6
@ac6: + x3
Variables
x4
x1 [0;+inf]
x2 [0;+inf]
x5
x6
x3 [0;+inf]
"
try
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
if length(ARGS) < 1
readptfstring(task,cqo1_ptf)
else
readdata(task,ARGS[1])
end
# (Optional) uncomment to see what happens when solution status is unknown
# putintparam(task,MSK_IPAR_INTPNT_MAX_ITERATIONS, 1)
# Optimize
trmcode = optimize(task)
solutionsummary(task,MSK_STREAM_LOG)
# We expect solution status OPTIMAL
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Optimal solution. Fetch and print it.
println("An optimal interior-point solution is located.")
numvar = getnumvar(task)
xx = getxx(task,MSK_SOL_ITR)
println("x = $xx")
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Dual infeasibility certificate found.")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal infeasibility certificate found.")
elseif solsta == MSK_SOL_STA_UNKNOWN
# The solutions status is unknown. The termination code
# indicates why the optimizer terminated prematurely.
println("The solution status is unknown.")
(symname, desc) = getcodedesc(trmcode)
println(" Termination code: $symname $desc")
else
println("An unexpected solution status $solsta is obtained.")
end
end
catch e
println("En error occurred: $(e.rcode), $(e.msg)")
end
sdo1.jl
##
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : sdo1.jl
#
# Purpose: Demonstrates how to solve a small mixed semidefinite and conic quadratic
# optimization problem using the MOSEK Julia API.
##
using Mosek
using Printf, SparseArrays
printstream(msg::String) = print(msg)
# Bound keys for constraints
bkc = [ MSK_BK_FX
MSK_BK_FX]
# Bound values for constraints
blc = [1.0, 0.5]
buc = [1.0, 0.5]
A = sparse( [1,2,2],[1,2,3],[1.0, 1.0, 1.0])
conesub = [1, 2, 3]
barci = [1, 2, 2, 3, 3]
barcj = [1, 1, 2, 2, 3]
barcval = [2.0, 1.0, 2.0, 1.0, 2.0]
barai = Any[ [1, 2, 3],
[1, 2, 3, 2, 3, 3] ]
baraj = Any[ [1, 2, 3],
[1, 1, 1, 2, 2, 3] ]
baraval = Any[ [1.0, 1.0, 1.0],
[1.0, 1.0, 1.0, 1.0, 1.0, 1.0] ]
numvar = 3
numcon = length(bkc)
barvardim = [3]
# Create a task object and attach log stream printer
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,printstream)
# Append 'numvar' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Append 'numcon' empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
# Append matrix variables of sizes in 'BARVARDIM'.
# The variables will initially be fixed at zero.
appendbarvars(task,barvardim)
# Set the linear term c_0 in the objective.
putcj(task, 1, 1.0)
# Set the bounds on variable j
# blx[j] <= x_j <= bux[j]
putvarboundsliceconst(task,1,numvar+1,
MSK_BK_FR,
-Inf,
+Inf)
# Set the bounds on constraints.
# blc[i] <= constraint_i <= buc[i]
putconboundslice(task,1,numcon+1, bkc,blc,buc)
# Append the conic quadratic cone
let afei = getnumafe(task)+1,
dom = appendquadraticconedomain(task,3)
appendafes(task,3)
putafefentrylist(task,[1,2,3],
[1,2,3],
[1.0,1.0,1.0])
appendaccseq(task,dom,afei,nothing)
end
# Input row i of A
putacolslice(task,1,numvar+1,
A.colptr[1:numvar], A.colptr[2:numvar+1],
A.rowval,A.nzval)
symc = appendsparsesymmat(task,barvardim[1],
barci,
barcj,
barcval)
syma0 = appendsparsesymmat(task,barvardim[1],
barai[1],
baraj[1],
baraval[1])
syma1 = appendsparsesymmat(task,barvardim[1],
barai[2],
baraj[2],
baraval[2])
putbarcj(task,1, [symc], [1.0])
putbaraij(task,1, 1, [syma0], [1.0])
putbaraij(task,2, 1, [syma1], [1.0])
# Input the objective sense (minimize/maximize)
putobjsense(task,MSK_OBJECTIVE_SENSE_MINIMIZE)
# Solve the problem and print summary
optimize(task)
writedata(task,"sdo1.ptf")
solutionsummary(task,MSK_STREAM_MSG)
# Get status information about the solution
prosta = getprosta(task,MSK_SOL_ITR)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Output a solution
xx = getxx(task,MSK_SOL_ITR)
barx = getbarxj(task,MSK_SOL_ITR, 1)
@printf("Optimal solution: \n xx = %s\n barx = %s\n", xx',barx')
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility.\n")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("Unknown solution status")
else
println("Other solution status")
end
end
sdo2.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : sdo2.jl
#
# Purpose : Solves the semidefinite problem with two symmetric variables:
#
# min <C1,X1> + <C2,X2>
# st. <A1,X1> + <A2,X2> = b
# (X2)_{1,2} <= k
#
# where X1, X2 are symmetric positive semidefinite,
#
# C1, C2, A1, A2 are assumed to be constant symmetric matrices,
# and b, k are constants.
#
using Mosek
# Input data
let numcon = 2, # Number of constraints.
numbarvar = 2,
dimbarvar = Int32[3, 4], # Dimension of semidefinite variables
# Objective coefficients concatenated
Cj = Int32[ 1, 1, 2, 2, 2, 2 ], # Which symmetric variable (j)
Ck = Int32[ 1, 3, 1, 2, 2, 3 ], # Which entry (k,l)->v
Cl = Int32[ 1, 3, 1, 1, 2, 3 ],
Cv = [ 1.0, 6.0, 1.0, -3.0, 2.0, 1.0 ],
# Equality constraints coefficients concatenated
Ai = Int32[ 1, 1, 1, 1, 1, 1 ], # Which constraint (i = 0)
Aj = Int32[ 1, 1, 1, 2, 2, 2 ], # Which symmetric variable (j)
Ak = Int32[ 1, 3, 3, 2, 2, 4 ], # Which entry (k,l)->v
Al = Int32[ 1, 1, 3, 1, 2, 4 ],
Av = [ 1.0, 1.0, 2.0, 1.0, -1.0, -3.0 ],
# The second constraint - one-term inequality
A2i = Int32[ 2 ], # Which constraint (i = 1)
A2j = Int32[ 2 ], # Which symmetric variable (j = 1)
A2k = Int32[ 2 ], # Which entry A(1,0) = A(0,1) = 0.5
A2l = Int32[ 1 ],
A2v = [ 0.5 ],
bkc = [ MSK_BK_FX,
MSK_BK_UP ],
blc = [ 23.0, 0.0 ],
buc = [ 23.0, -3.0 ]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Append numcon empty constraints.
# The constraints will initially have no bounds.
appendcons(task,numcon)
# Append numbarvar semidefinite variables.
appendbarvars(task,dimbarvar)
# Set objective (6 nonzeros).
putbarcblocktriplet(task,Cj, Ck, Cl, Cv)
# Set the equality constraint (6 nonzeros).
putbarablocktriplet(task,Ai, Aj, Ak, Al, Av)
# Set the inequality constraint (1 nonzero).
putbarablocktriplet(task,A2i, A2j, A2k, A2l, A2v)
# Set constraint bounds
putconboundslice(task,1, 3, bkc, blc, buc)
# Run optimizer
optimize(task)
solutionsummary(task,MSK_STREAM_MSG)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
# Retrieve the soution for all symmetric variables
println("Solution (lower triangular part vectorized):")
for i in 1:numbarvar
barx = getbarxj(task,MSK_SOL_ITR, i)
println("X$i: $barx")
end
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER || if solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility certificate found.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("The status of the solution could not be determined.")
else
println("Other solution status.")
end
end
end
end
sdo_lmi.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : sdo_lmi.jl
#
# Purpose : To solve a problem with an LMI and an affine conic constrained problem with a PSD term
#
# minimize Tr [1, 0; 0, 1]*X + x(1) + x(2) + 1
#
# subject to Tr [0, 1; 1, 0]*X - x(1) - x(2) >= 0
# x(1) [0, 1; 1, 3] + x(2) [3, 1; 1, 0] - [1, 0; 0, 1] >> 0
# X >> 0
#
using Mosek
let numafe = 4, # Number of affine expressions.
numvar = 2, # Number of scalar variables
dimbarvar = [2], # Dimension of semidefinite cone
lenbarvar = [2 * (2 + 1) / 2], # Number of scalar SD variables
barc_j = [1, 1],
barc_k = [1, 2],
barc_l = [1, 2],
barc_v = [1.0, 1.0],
afeidx = [1, 1, 2, 3, 3, 4],
varidx = [1, 2, 2, 1, 2, 1],
f_val = Float64[-1, -1, 3, sqrt(2.0), sqrt(2.0), 3],
g = Float64[0, -1, 0, -1],
barf_i = [1, 1],
barf_j = [1, 1],
barf_k = [1, 2],
barf_l = [1, 1],
barf_v = [0.0, 1.0]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
# Append 'NUMAFE' empty affine expressions.
appendafes(task,numafe)
# Append 'NUMVAR' variables.
# The variables will initially be fixed at zero (x=0).
appendvars(task,numvar)
# Append 'NUMBARVAR' semidefinite variables.
appendbarvars(task,dimbarvar)
# Optionally add a constant term to the objective.
putcfix(task,1.0)
# Set the linear term c_j in the objective.
putcj(task,1, 1.0)
putcj(task,2, 1.0)
for j in 1:numvar
putvarbound(task,j, MSK_BK_FR, -0.0, 0.0)
end
# Set the linear term barc_j in the objective.
putbarcblocktriplet(task, barc_j, barc_k, barc_l, barc_v)
# Set up the affine conic constraints
# Construct the affine expressions
# F matrix
putafefentrylist(task,afeidx, varidx, f_val)
# g vector
putafegslice(task,1, 5, g)
# barF block triplets
putafebarfblocktriplet(task,barf_i, barf_j, barf_k, barf_l, barf_v)
# Append R+ domain and the corresponding ACC
appendacc(task,appendrplusdomain(task,1), [1], nothing)
# Append SVEC_PSD domain and the corresponding ACC
appendacc(task,appendsvecpsdconedomain(task,3), [2,3,4], nothing)
# Run optimizer
optimize(task)
# Print a summary containing information
# about the solution for debugging purposes
solutionsummary(task,MSK_STREAM_MSG)
solsta = getsolsta(task,MSK_SOL_ITR)
if solsta == MSK_SOL_STA_OPTIMAL
xx = getxx(task,MSK_SOL_ITR)
barx = getbarxj(task,MSK_SOL_ITR,1); # Request the interior solution.
println("Optimal primal solution, x = $xx, barx = $barx")
elseif solsta == MSK_SOL_STA_PRIM_INFEAS_CER || solsta == MSK_SOL_STA_DUAL_INFEAS_CER
println("Primal or dual infeasibility certificate found.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("The status of the solution could not be determined.")
else
println("Other solution status.")
end
end
end
sensitivity.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : sensitivity.jl
#
# Purpose : To demonstrate how to perform sensitivity
# analysis from the API on a small problem:
#
# minimize
#
# obj: +1 x11 + 2 x12 + 5 x23 + 2 x24 + 1 x31 + 2 x33 + 1 x34
# st
# c1: + x11 + x12 <= 400
# c2: + x23 + x24 <= 1200
# c3: + x31 + x33 + x34 <= 1000
# c4: + x11 + x31 = 800
# c5: + x12 = 100
# c6: + x23 + x33 = 500
# c7: + x24 + x34 = 500
#
# The example uses basis type sensitivity analysis.
using Mosek
prob_ptf = "Task ''
Objective
Minimize + 1 x11 + 2 x12 + 5 x23 + 2 x24 + 1 x31 + 2 x33 + 1 x34
Constraints
c1 [-inf;400 ] + x11 + x12
c2 [-inf;1200] + x23 + x24
c3 [-inf;1000] + x31 + x33 + x34
c4 [800] + x11 + x31
c5 [100] + x12
c6 [500] + x23 + x33
c7 [500] + x24 + x34
Variables
x11 [0;+inf]
x12 [0;+inf]
x23 [0;+inf]
x24 [0;+inf]
x31 [0;+inf]
x33 [0;+inf]
x34 [0;+inf]
"
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
readptfstring(task,prob_ptf)
# A maximization problem
optimize(task)
checkmem(task,@__FILE__,@__LINE__)
# Analyze upper bound on c1 and the equality constraint on c4
let subi = [1, 4],
marki = [MSK_MARK_UP, MSK_MARK_UP],
# Analyze lower bound on the variables x12 and x31
subj = [2, 5],
markj = [MSK_MARK_LO,MSK_MARK_LO],
(leftpricei,
rightpricei,
leftrangei,
rightrangei,
leftpricej,
rightpricej,
leftrangej,
rightrangej) = primalsensitivity(task,
subi,
marki,
subj,
markj)
checkmem(task,@__FILE__,@__LINE__)
println("Results from sensitivity analysis on bounds: ")
println("For constraints:")
for i in 1:2
println("leftprice = $(leftpricei[i]) | rightprice = $(rightpricei[i]) | leftrange = $(leftrangei[i]) | rightrange = $(rightrangei[i])");
end
println("For variables:")
for i in 1:2
println("leftprice = $(leftpricej[i]) rightprice = $(rightpricej[i]) leftrange = $(leftrangej[i]) rightrange = $(rightrangej[i])")
end
end
let subc = Int32[2, 5],
(leftprice,
rightprice,
leftrange,
rightrange) = dualsensitivity(task,subc)
println("Results from sensitivity analysis on objective coefficients:")
for i in 1:2
println("leftprice = $(leftprice[i]) rightprice = $(rightprice[i]) leftrange = $(leftrange[i]) rightrange = $( rightrange[i])")
end
end
end
simple.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : simple.jl
#
# Purpose : Demonstrates a very simple example using MOSEK by
# reading a problem file, solving the problem and
# writing the problem+solution to a file.
using Mosek
if length(ARGS) < 1
println("Syntax: simple FILENAME [ OUTFILE ]")
else
let filename = ARGS[1],
outfile = if length(ARGS) > 1 ARGS[2] else Nothing end
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# We assume that a problem file was given as the first command
# line argument (received in `args')
readdata(task,filename)
# Solve the problem
optimize(task)
# Print a summary of the solution
solutionsummary(task,MSK_STREAM_LOG)
# If an output file was specified, save problem to file
if outfile != Nothing
# If using OPF format, these parameters will specify what to include in output
putintparam(task,MSK_IPAR_OPF_WRITE_SOLUTIONS, MSK_ON)
putintparam(task,MSK_IPAR_OPF_WRITE_PROBLEM, MSK_ON)
putintparam(task,MSK_IPAR_OPF_WRITE_HINTS, MSK_OFF)
putintparam(task,MSK_IPAR_OPF_WRITE_PARAMETERS, MSK_OFF)
putintparam(task,MSK_IPAR_PTF_WRITE_SOLUTIONS, MSK_ON)
writedata(task,outfile)
end
end
end
end
solutionquality.jl
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : solutionquality.jl
#
# Purpose : To demonstrate how to examine the quality of a solution.
using Mosek
if length(ARGS) < 1
println("Usage: solutionquality FILENAME")
else
filename = ARGS[1]
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
#putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# We assume that a problem file was given as the first command
# line argument (received in `args')
readdata(task,filename)
# Solve the problem
optimize(task)
# System.Out.Println (a summary of the solution
solutionsummary(task,MSK_STREAM_LOG)
solsta = getsolsta(task,MSK_SOL_BAS)
(pobj,
pviolcon,
pviolvar,
pviolbarvar,
pviolcones,
pviolitg,
dobj,
dviolcon,
dviolvar,
dviolbarvar,
dviolcones) = getsolutioninfo(task,MSK_SOL_BAS)
if solsta == MSK_SOL_STA_OPTIMAL
abs_obj_gap = abs(dobj - pobj)
rel_obj_gap = abs_obj_gap / (1.0 + min(abs(pobj), abs(dobj)))
max_primal_viol = max(pviolcon, pviolvar)
max_primal_viol = max(max_primal_viol, pviolbarvar)
max_primal_viol = max(max_primal_viol, pviolcones)
max_dual_viol = max(dviolcon, dviolvar)
max_dual_viol = max(max_dual_viol, dviolbarvar)
max_dual_viol = max(max_dual_viol, dviolcones)
# Assume the application needs the solution to be within
# 1e-6 ofoptimality in an absolute sense. Another approach
# would be looking at the relative objective gap
println("Customized solution information.")
println(" Absolute objective gap: $abs_obj_gap")
println(" Relative objective gap: $rel_obj_gap")
println(" Max primal violation : $max_primal_viol")
println(" Max dual violation : $max_dual_viol")
accepted = true
if rel_obj_gap > 1e-6
println("Warning: The relative objective gap is LARGE.")
accepted = false
end
# We will accept a primal infeasibility of 1e-8 and
# dual infeasibility of 1e-6. These number should chosen problem
# dependent.
if max_primal_viol > 1e-8
println("Warning: Primal violation is too LARGE")
accepted = false
end
if max_dual_viol > 1e-6
println("Warning: Dual violation is too LARGE.")
accepted = false
end
if accepted
numvar = getnumvar(task)
println("Optimal primal solution")
xx = getxxslice(task,MSK_SOL_BAS,1,numvar+1)
println(" xx = $xx")
else
# print etailed information about the solution
analyzesolution(task,MSK_STREAM_LOG, MSK_SOL_BAS)
end
elseif solsta == MSK_SOL_STA_DUAL_INFEAS_CER || solsta == MSK_SOL_STA_PRIM_INFEAS_CER
println("Primal or dual infeasibility certificate found.")
elseif solsta == MSK_SOL_STA_UNKNOWN
println("The status of the solution is unknown.")
else
println("Other solution status")
end
end
end
solvebasis.jl
##
#
# File: solvebasis.jl
#
# Purpose : To demonstrate the usage of
# MSK_solvewithbasis on the problem:
#
# maximize x0 + x1
# st.
# x0 + 2.0 x1 <= 2
# x0 + x1 <= 6
# x0 >= 0, x1>= 0
#
# The problem has the slack variables
# xc0, xc1 on the constraints
# and the variabels x0 and x1.
#
# maximize x0 + x1
# st.
# x0 + 2.0 x1 -xc1 = 2
# x0 + x1 -xc2 = 6
# x0 >= 0, x1>= 0,
# xc1 <= 0 , xc2 <= 0
##
using Mosek
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putobjname(task,"solvebasis")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
numcon = 2
numvar = 2
c = Float64[1.0, 1.0]
ptrb = Int64[1,3]
ptre = Int64[3,4]
asub = Int32[1,2,
1,2]
aval = Float64[1.0, 1.0,
2.0, 1.0]
bkc = Boundkey[MSK_BK_UP
MSK_BK_UP]
blc = Float64[-Inf
-Inf]
buc = Float64[2.0
6.0]
bkx = Boundkey[MSK_BK_LO
MSK_BK_LO]
blx = Float64[0.0
0.0]
bux = Float64[Inf
Inf]
w1 = Float64[2.0, 6.0]
w2 = Float64[1.0, 0.0]
inputdata(task,
Int32(numcon), Int32(numvar),
c,
0.0,
ptrb,
ptre,
asub,
aval,
bkc,
blc,
buc,
bkx,
blx,
bux)
putobjsense(task,MSK_OBJECTIVE_SENSE_MAXIMIZE)
r = optimize(task)
if r != MSK_RES_OK
println("Mosek warning: $r")
end
basis = initbasissolve(task)
#List basis variables corresponding to columns of B
varsub = Int32[1, 2]
for i in 1:numcon
if basis[varsub[i]] < numcon
println("Basis variable no $i is xc$(basis[i])")
else
println("Basis variable no $i is x$(basis[i]-numcon)")
# solve Bx = w1
# varsub contains index of non-zeros in b.
# On return b contains the solution x and
# varsub the index of the non-zeros in x.
nz = 2
nz = solvewithbasis(task,false, nz, varsub, w1)
println("nz = $nz")
println("Solution to Bx = $w1")
for i in 1:nz
if basis[varsub[i]] < numcon
println("xc $(basis[varsub[i]]) = $(w1[varsub[i]])")
else
println("x$(basis[varsub[i]] - numcon) = $(w1[varsub[i]])")
end
end
# Solve B^Tx = w2
nz = 1
varsub[1] = 1
nz = solvewithbasis(task,true,nz,varsub,w2)
println(nz)
println("Solution to B^Tx = $(w2)")
for i in 1:nz
if basis[varsub[i]] < numcon
println("xc$(basis[varsub[i]]) = $(w2[varsub[i]])")
else
println("x$(basis[varsub[i]] - numcon) = $(w2[varsub[i]])")
end
end
end
end
end
solvelinear.jl
#
# Copyright : Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File : solvelinear.jl
##
# Purpose : To demonstrate the usage of MSK_solvewithbasis
# when solving the linear system:
#
# 1.0 x1 = b1
# -1.0 x0 + 1.0 x1 = b2
#
# with two different right hand sides
#
# b = (1.0, -2.0)
#
# and
#
# b = (7.0, 0.0)
using Mosek
function setup(task :: Mosek.Task,
aval :: Vector{Float64},
asub :: Vector{Int32},
ptrb :: Vector{Int64},
ptre :: Vector{Int64},
numvar :: Int32)
appendvars(task,numvar)
appendcons(task,numvar)
putacolslice(task,1,numvar+1,ptrb,ptre,asub,aval)
putconboundsliceconst(task,1,numvar+1,MSK_BK_FX,0.0,0.0)
putvarboundsliceconst(task,1,numvar+1,MSK_BK_FR,-Inf,Inf)
# Define a basic solution by specifying status keys for variables
# & constraints.
deletesolution(task,MSK_SOL_BAS)
putskcslice(task,MSK_SOL_BAS, 1, numvar+1, fill(MSK_SK_FIX,numvar))
putskxslice(task,MSK_SOL_BAS, 1, numvar+1, fill(MSK_SK_BAS,numvar))
return initbasissolve(task)
end
let numcon = Int32(2),
numvar = Int32(2),
aval = [ -1.0 ,
1.0, 1.0 ],
asub = Int32[ 2,
1, 2 ],
ptrb = Int64[1, 2],
ptre = Int64[2, 4]
# bsub = new int[numvar];
# b = new double[numvar];
# basis = new int[numvar];
maketask() do task
# Use remote server: putoptserverhost(task,"http://solve.mosek.com:30080")
putstreamfunc(task,MSK_STREAM_LOG,msg -> print(msg))
# Put A matrix and factor A. Call this function only once for a
# given task.
basis = setup(task,
aval,
asub,
ptrb,
ptre,
numvar)
# now solve rhs
let b = Float64[1,-2],
bsub = Int32[1,2],
nz = solvewithbasis(task,false, 2, bsub, b)
println("Solution to Bx = b:")
# Print solution and show correspondents to original variables in the problem
for i in 1:nz
if basis[bsub[i]] <= numcon
println("This should never happen")
else
println("x $(basis[bsub[i]] - numcon) = $(b[bsub[i]])")
end
end
end
let b = Float64[7,0],
bsub = Int32[1,0],
nz = solvewithbasis(task,false,1, bsub, b)
println("Solution to Bx = b:")
# Print solution and show correspondents to original variables in the problem
for i in 1:nz
if (basis[bsub[i]] <= numcon)
println("This should never happen")
else
println("x $(basis[bsub[i]] - numcon) = $(b[bsub[i]])")
end
end
end
end
end
sparsecholesky.jl
#
# Copyright: Copyright (c) MOSEK ApS, Denmark. All rights reserved.
#
# File: sparsecholesky.jl
#
# Purpose: Demonstrate the sparse Cholesky factorization.
#
using Mosek
function printsparse(n :: Int32,
perm :: Vector{Int32},
diag :: Vector{Float64},
lnzc :: Vector{Int32},
lptrc :: Vector{Int64},
lensubnval :: Int64,
lsubc :: Vector{Int32},
lvalc :: Vector{Float64})
println("P = $(perm)")
println("diag(D) = $(diag)")
l = zeros(Float64,(n,n))
for j in 1:n
for i in lptrc[j]:lptrc[j]+lnzc[j]-1
l[lsubc[i],j] = lvalc[i]
end
end
println("L = $l")
end
function main()
# Example from the manual
# Observe that anzc, aptrc, asubc and avalc only specify the lower triangular part.
n = Int32(4)
anzc = Int32[4, 1, 1, 1]
asubc = Int32[1, 2, 3, 4,
2,
3,
4]
aptrc = Int64[1, 5, 6, 7]
avalc = Float64[4.0, 1.0, 1.0, 1.0,
1.0,
1.0,
1.0]
b = Float64[13.0, 3.0, 4.0, 5.0]
(perm,
diag,
lnzc,
lptrc,
lensubnval,
lsubc,
lvalc) = computesparsecholesky(Int32(0), #Mosek chooses number of threads
Int32(1), #Apply reordering heuristic
1.0e-14, #Singularity tolerance
anzc, aptrc, asubc, avalc)
printsparse(n, perm, diag, lnzc, lptrc, lensubnval, lsubc, lvalc)
# Permuted b is stored as x.
x = b[perm]
# Compute inv(L)*x.
@show x
sparsetriangularsolvedense(MSK_TRANSPOSE_NO, lnzc, lptrc, lsubc, lvalc, x)
# Compute inv(L^T)*x.
sparsetriangularsolvedense(MSK_TRANSPOSE_YES, lnzc, lptrc, lsubc, lvalc, x)
println("\nSolution A x = b, x = $([ x[j] for i in 1:n for j in 1:n if perm[j] == i ])")
n = Int32(3)
anzc = Int32[3, 2, 1]
asub = Int32[0, 1, 2, 1, 2, 2]
aptr = Int64[0, 3, 5, ]
avalc = Float64[1.0, 1.0, 1.0, 1.0, 1.0, 1.0]
(perm,
diag,
lnzc,
lptrc,
lensubnval,
lsubc,
lvalc) = computesparsecholesky(0, #Mosek chooses number of threads
1, #Apply reordering heuristic
1.0e-14, #Singularity tolerance
anzc, aptrc, asubc, avalc)
printsparse(n, perm, diag, lnzc, lptrc, lensubnval, lsubc, lvalc)
end
main()
