Closed chriscoey closed 8 years ago
Could you provide an example script that reproduces this error ?
Here is a script (I can produce more if it is helpful):
c = [0.0,0.0,1.0,0.0,0.0,0.0,0.0]
b = [0.0,10.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0]
I = [1,12,1,21,1,2,3,7,8,10,13,16,17,19,22,2,4,7,16,2,5,7,16,2,6,10,19]
J = [1,1,2,2,3,4,4,4,4,4,4,4,4,4,4,5,5,5,5,6,6,6,6,7,7,7,7]
V = [-1.0,-1.0,-1.0,-1.0,1.0,1.0,-1.0,-0.44999999999999996,0.45,-0.4500000000000001,0.,-0.44999999999999996,0.45,-0.4500000000000001,0.,1.0,-1.0,-0.9000000000000001,-0.9000000000000001,1.0,-1.0,-0.22500000000000003,-0.22500000000000003,1.0,-1.0,-0.22500000000000003,-0.22500000000000003]
A = sparse(I, J, V, length(b), length(c))
cone_con = [(:Zero,1:1),(:NonNeg,2:2),(:NonNeg,3:6),(:SDP,7:12),(:Zero,13:15),(:SDP,16:21),(:Zero,22:24)]
cone_var = [(:Free,1:7)]
using MathProgBase, Mosek
m = MathProgBase.ConicModel(MosekSolver())
MathProgBase.loadproblem!(m, c, A, b, cone_con, cone_var)
MathProgBase.optimize!(m)
@show MathProgBase.status(m)
@show MathProgBase.getsolution(m)
@show MathProgBase.getdual(m)
Here is the output I get:
julia> MathProgBase.optimize!(m)
Computer
Platform : MACOSX/64-X86
Problem
Name :
Objective sense : min
Type : CONIC (conic optimization problem)
Constraints : 24
Cones : 0
Scalar variables : 7
Matrix variables : 2
Integer variables : 0
Optimizer started.
Conic interior-point optimizer started.
Presolve started.
Linear dependency checker started.
Linear dependency checker terminated.
Eliminator started.
Total number of eliminations : 1
Eliminator terminated.
Eliminator - tries : 1 time : 0.00
Eliminator - elim's : 1
Lin. dep. - tries : 1 time : 0.00
Lin. dep. - number : 0
Presolve terminated. Time: 0.00
Optimizer - threads : 8
Optimizer - solved problem : the primal
Optimizer - Constraints : 13
Optimizer - Cones : 0
Optimizer - Scalar variables : 8 conic : 0
Optimizer - Semi-definite variables: 2 scalarized : 12
Factor - setup time : 0.00 dense det. time : 0.00
Factor - ML order time : 0.00 GP order time : 0.00
Factor - nonzeros before factor : 59 after factor : 63
Factor - dense dim. : 0 flops : 1.01e+03
ITE PFEAS DFEAS GFEAS PRSTATUS POBJ DOBJ MU TIME
0 5.0e+00 1.0e+00 1.0e+00 0.00e+00 0.000000000e+00 0.000000000e+00 1.0e+00 0.00
1 1.3e+00 2.6e-01 2.6e-01 -1.45e-01 7.446299248e-01 1.133476135e+00 2.6e-01 0.00
2 1.8e-01 3.6e-02 3.6e-02 5.79e-01 1.297888733e+00 1.369246549e+00 3.6e-02 0.00
3 2.4e-02 4.7e-03 4.7e-03 8.91e-01 9.564969913e-01 9.704884077e-01 4.7e-03 0.00
4 3.4e-03 6.8e-04 6.8e-04 8.77e-01 8.651595034e-01 8.677722067e-01 6.8e-04 0.00
5 2.7e-04 5.4e-05 5.4e-05 9.69e-01 8.507793221e-01 8.509928746e-01 5.4e-05 0.00
6 1.6e-05 3.3e-06 3.3e-06 9.97e-01 8.496036890e-01 8.496167452e-01 3.3e-06 0.00
7 9.2e-07 1.8e-07 1.8e-07 1.00e+00 8.495321450e-01 8.495328766e-01 1.8e-07 0.00
8 8.8e-09 1.8e-09 1.8e-09 1.00e+00 8.495279655e-01 8.495279725e-01 1.8e-09 0.00
Interior-point optimizer terminated. Time: 0.00.
Optimizer terminated. Time: 0.02
0
julia> @show MathProgBase.status(m)
MathProgBase.status(m) = :Optimal
:Optimal
julia> @show MathProgBase.getsolution(m)
MathProgBase.getsolution(m) = [0.2581288811960397,0.591399084272713,0.8495279654687528,2.965654193424076,3.779037063878916,3.49217819844456e-7,3.255308073884936]
7-element Array{Float64,1}:
0.258129
0.591399
0.849528
2.96565
3.77904
3.49218e-7
3.25531
julia> @show MathProgBase.getdual(m)
ERROR: BoundsError: attempt to access 2-element Array{Array{Float64,1},1}:
[0.06662544466703094,0.04302513616541794,-0.25811904947274705,0.027784617540638323,-0.16668717064440128,0.9999999997848541]
[0.027780822124551928,0.09858424790224177,-0.1666757940939465,0.3498404168267607,-0.5914730859020215,0.9999999999546839]
at index [3]
in getdual at /Users/chris/.julia/v0.4/Mosek/src/MosekConicInterface.jl:526
for comparison, here is the output I get from SCS solver:
julia> MathProgBase.optimize!(m)
----------------------------------------------------------------------------
SCS v1.1.8 - Splitting Conic Solver
(c) Brendan O'Donoghue, Stanford University, 2012-2015
----------------------------------------------------------------------------
Lin-sys: sparse-direct, nnz in A = 25
eps = 1.00e-04, alpha = 1.80, max_iters = 20000, normalize = 1, scale = 5.00
Variables n = 7, constraints m = 24
Cones: primal zero / dual free vars: 7
linear vars: 5
sd vars: 12, sd blks: 2
Setup time: 1.91e-03s
----------------------------------------------------------------------------
Iter | pri res | dua res | rel gap | pri obj | dua obj | kap/tau | time (s)
----------------------------------------------------------------------------
0| inf inf nan -inf inf inf 5.12e-04
100| 4.04e-04 2.37e-02 3.95e-03 9.03e-01 9.14e-01 0.00e+00 1.16e-03
200| 3.17e-04 1.62e-02 4.02e-03 8.83e-01 8.94e-01 0.00e+00 1.79e-03
300| 2.56e-04 1.19e-02 3.48e-03 8.72e-01 8.81e-01 0.00e+00 2.41e-03
400| 2.10e-04 9.02e-03 3.00e-03 8.64e-01 8.73e-01 0.00e+00 3.03e-03
500| 1.74e-04 7.04e-03 2.57e-03 8.60e-01 8.67e-01 0.00e+00 3.65e-03
600| 1.45e-04 5.59e-03 2.19e-03 8.57e-01 8.63e-01 0.00e+00 4.27e-03
700| 1.21e-04 4.49e-03 1.86e-03 8.54e-01 8.59e-01 0.00e+00 4.89e-03
800| 1.01e-04 3.64e-03 1.58e-03 8.53e-01 8.57e-01 0.00e+00 5.51e-03
900| 8.46e-05 2.97e-03 1.34e-03 8.52e-01 8.56e-01 0.00e+00 6.12e-03
1000| 7.09e-05 2.44e-03 1.13e-03 8.51e-01 8.54e-01 0.00e+00 6.74e-03
1100| 5.94e-05 2.01e-03 9.58e-04 8.51e-01 8.53e-01 0.00e+00 7.36e-03
1200| 4.98e-05 1.66e-03 8.07e-04 8.50e-01 8.53e-01 0.00e+00 7.98e-03
1300| 4.17e-05 1.37e-03 6.80e-04 8.50e-01 8.52e-01 0.00e+00 8.59e-03
1400| 3.50e-05 1.14e-03 5.72e-04 8.50e-01 8.52e-01 0.00e+00 9.21e-03
1500| 2.93e-05 9.47e-04 4.80e-04 8.50e-01 8.51e-01 0.00e+00 9.82e-03
1600| 2.45e-05 7.88e-04 4.03e-04 8.50e-01 8.51e-01 0.00e+00 1.05e-02
1700| 2.05e-05 6.56e-04 3.39e-04 8.50e-01 8.51e-01 0.00e+00 1.11e-02
1800| 1.72e-05 5.47e-04 2.84e-04 8.50e-01 8.50e-01 0.00e+00 1.17e-02
1900| 1.44e-05 4.56e-04 2.38e-04 8.50e-01 8.50e-01 0.00e+00 1.23e-02
2000| 1.20e-05 3.80e-04 2.00e-04 8.50e-01 8.50e-01 0.00e+00 1.29e-02
2100| 1.01e-05 3.17e-04 1.67e-04 8.50e-01 8.50e-01 0.00e+00 1.35e-02
2200| 8.44e-06 2.65e-04 1.40e-04 8.50e-01 8.50e-01 0.00e+00 1.41e-02
2300| 7.06e-06 2.21e-04 1.17e-04 8.50e-01 8.50e-01 0.00e+00 1.48e-02
2400| 5.91e-06 1.85e-04 9.81e-05 8.50e-01 8.50e-01 0.00e+00 1.54e-02
2500| 4.94e-06 1.54e-04 8.21e-05 8.50e-01 8.50e-01 0.00e+00 1.60e-02
2600| 4.13e-06 1.29e-04 6.87e-05 8.50e-01 8.50e-01 0.00e+00 1.66e-02
2700| 3.46e-06 1.08e-04 5.75e-05 8.50e-01 8.50e-01 0.00e+00 1.72e-02
2760| 3.11e-06 9.69e-05 5.17e-05 8.50e-01 8.50e-01 0.00e+00 1.75e-02
----------------------------------------------------------------------------
Status: Solved
Timing: Solve time: 1.76e-02s
Lin-sys: nnz in L factor: 61, avg solve time: 3.53e-07s
Cones: avg projection time: 5.11e-06s
----------------------------------------------------------------------------
Error metrics:
dist(s, K) = 2.2132e-09, dist(y, K*) = 9.3131e-10, s'y/m = -9.7499e-11
|Ax + s - b|_2 / (1 + |b|_2) = 3.1057e-06
|A'y + c|_2 / (1 + |c|_2) = 9.6854e-05
|c'x + b'y| / (1 + |c'x| + |b'y|) = 5.1706e-05
----------------------------------------------------------------------------
c'x = 0.8495, -b'y = 0.8497
============================================================================
0.8495352484237209
julia> @show MathProgBase.status(m)
MathProgBase.status(m) = :Optimal
:Optimal
julia> @show MathProgBase.getsolution(m)
MathProgBase.getsolution(m) = [0.2585352976879626,0.5909978259665712,0.8495352484237209,2.950007731763664,3.7708023318365167,4.2743511306364275e-8,3.27920973115259]
7-element Array{Float64,1}:
0.258535
0.590998
0.849535
2.95001
3.7708
4.27435e-8
3.27921
julia> @show MathProgBase.getdual(m)
MathProgBase.getdual(m) = [-0.9999999998163523,0.08493642159298997,0.0,0.0,0.06370292305356898,0.0,0.06683902395033697,0.060756649730726565,-0.3656175766611347,0.02761388683088022,-0.23500428789197625,0.9999862689712925,0.0,0.0,0.0,0.027613898433779184,0.06944341154397124,-0.23500438167958973,0.34927248091608026,-0.8357844423855383,0.9999866469605179,0.0,0.0,0.0]
24-element Array{Float64,1}:
-1.0
0.0849364
0.0
0.0
0.0637029
0.0
0.066839
0.0607566
-0.365618
0.0276139
-0.235004
0.999986
0.0
0.0
0.0
0.0276139
0.0694434
-0.235004
0.349272
-0.835784
0.999987
0.0
0.0
0.0
Below is a second example where Mosek getdual fails. For this one, SCS getdual also fails (I will submit an issue to them):
c = [0.0,0.0,0.0,0.0,0.0,0.0,1.0,1.0]
b = [10.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0]
I = [1,2,10,11,13,16,19,20,22,25,1,3,10,11,13,16,19,20,22,25,1,4,10,11,13,16,19,20,22,25,1,5,10,19,1,6,10,11,13,16,19,20,22,25,1,7,10,11,13,16,19,20,22,25,8,15,9,24]
J = [1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3,4,4,4,4,5,5,5,5,5,5,5,5,5,5,6,6,6,6,6,6,6,6,6,6,7,7,8,8]
V = [1.0,1.0,-0.44999999999999996,0.45000000000000007,-0.4500000000000001,5.551115123125783e-17,-0.44999999999999996,0.45000000000000007,-0.4500000000000001,5.551115123125783e-17,1.0,1.0,-0.7681980515339464,0.31819805153394654,-0.13180194846605373,5.551115123125783e-17,-0.7681980515339464,0.31819805153394654,-0.13180194846605373,5.551115123125783e-17,1.0,1.0,-0.9000000000000001,1.102182119232618e-16,-1.3497838043956718e-32,1.232595164407831e-32,-0.9000000000000001,1.102182119232618e-16,-1.3497838043956718e-32,1.232595164407831e-32,1.0,1.0,-0.22500000000000003,-0.22500000000000003,1.0,1.0,-0.11250000000000003,0.1125,-0.11249999999999999,-1.3877787807814457e-17,-0.11250000000000003,0.1125,-0.11249999999999999,-1.3877787807814457e-17,1.0,1.0,-8.436148777472949e-34,1.3777276490407724e-17,-0.22500000000000003,-1.5407439555097887e-33,-8.436148777472949e-34,1.3777276490407724e-17,-0.22500000000000003,-1.5407439555097887e-33,1.0,-1.0,1.0,-1.0]
A = sparse(I, J, V, length(b), length(c))
cone_con = [(:NonNeg,[1]),(:NonPos,[2,3,4,5,6,7,8,9]),(:SDP,10:15),(:Zero,16:18),(:SDP,19:24),(:Zero,25:27)]
cone_var = [(:NonNeg,[1,2,3,4,5,6,7,8])]
using MathProgBase, Mosek
m = MathProgBase.ConicModel(MosekSolver())
MathProgBase.loadproblem!(m, c, A, b, cone_con, cone_var)
MathProgBase.optimize!(m)
@show MathProgBase.status(m)
@show MathProgBase.getsolution(m)
@show MathProgBase.getdual(m)
another even smaller problem on which both Mosek and SCS fail with getdual
:
c = [-0.0,-0.0,-0.0,-0.0,-0.0,-0.0,-1.0]
b = [10.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]
I = [1,2,8,9,10,11,1,3,8,9,10,11,1,4,8,9,10,11,1,5,8,1,6,8,9,10,11,1,7,8,9,10,11,8,10]
J = [1,1,1,1,1,1,2,2,2,2,2,2,3,3,3,3,3,3,4,4,4,5,5,5,5,5,5,6,6,6,6,6,6,7,7]
V = [1.0,1.0,-0.44999999999999996,0.45000000000000007,-0.4500000000000001,0.0,1.0,1.0,-0.7681980515339464,0.31819805153394654,-0.13180194846605373,0.0,1.0,1.0,-0.9000000000000001,0.0,0.0,0.0,1.0,1.0,-0.22500000000000003,1.0,1.0,-0.11250000000000003,0.1125,-0.11249999999999999,0.0,1.0,1.0,0.0,0.0,-0.22500000000000003,0.0,1.0,1.0]
A = sparse(I, J, V, length(b), length(c))
cone_con = [(:NonNeg,[1]),(:NonPos,[2,3,4,5,6,7]),(:SDP,8:10),(:Zero,11:11)]
cone_var = [(:NonNeg,[1,2,3,4,5,6]),(:Free,[7])]
using MathProgBase, Mosek
m = MathProgBase.ConicModel(MosekSolver())
MathProgBase.loadproblem!(m, c, A, b, cone_con, cone_var)
MathProgBase.optimize!(m)
@show MathProgBase.status(m)
@show MathProgBase.getsolution(m)
@show MathProgBase.getdual(m)
@ulfworsoe, any idea what might be causing this crash when accessing duals for SDPs?
I can see you already figured it out. Works in the HEAD revision, so I'll close this.
Here is an example of the error I get when I call getdual:
Unfortunately, I don't understand what's going on here at line 526: