MosekMathProgModel numvar field is broken

[joehuchette]

The field is not updated when variables are added internally (for instance, when auxiliary variables are created in addquadconstr!); this bit me for pretty much the whole weekend when I couldn't find out why multiple quadratic constraints were breaking :(

It seems that there is no reason to carry around this state when it's liable to break like this and a call to getnumvar(m.task) is so easy. I'm not sure if there was a reason to do it this way that I'm missing @ulfworsoe @wflu, but if not I'll go through the MathProgBase interface and change it.

[ulfworsoe]

Yes, it is a hack, and there probably is a bug. I will see if I can fix it. The reason it is that SOCP a constraint can have the form xz-y'y>0, while MOSEK accepts the form 1/2 xz-y'y>0. I have to add an extra variable x=2w and add the constraint as 1/2wz-y'y>0. The numvar should keep track of the number of variable the user has added, since I don't want to report the internal variable values. Of course this breaks if the user adds variables after a cone.
I don't think it was possible to add variable on the go when I wrote this. I'll see if I can work around it, but ultimately I may have build up the whole model in the Julia struct and commit it to MOSEK when optimize!() is called.

[joehuchette]

One fix would be to change the MathProgBase interface, as I don't think that this will be the last time that auxiliary variables will bite us. One simple solution would be to have addconstr! and similar functions return integers indicating the number of auxiliary variables introduced when adding the constraint. Then at least the changes are transparent to the user (e.g. a JuMP model could pad the number of variables appropriately).

cc @mlubin @IainNZ

[mlubin]

If/when solvers need to add auxiliary variables is a solver-specific detail that I don't think users of the interface should have to deal with, and it will add unneeded complexity to JuMP if we have to pad variables for this. The second option, building up the model in a Julia struct and then committing it to MOSEK, seems better in my opinion.

[joehuchette]

I agree that it'd be nice to keep this complexity out of JuMP, but it seems like a lot of work at the solver level. I think it's also possible to "quarantine" the auxiliary variables until after all user variables are added, and then we only need to carry around an integer equal to the number of auxiliary variables added so far. E.g. in the norm stuff I did, I need to know the solver-level indices for aux variables I'm introducing, but if I just add these constraints after all user variables are added, I just need to know the offset.

[ulfworsoe]

I don't think the interface should expose non-user-generated variables. Should I propose one change in the interface, I would say that addconstr! and addvar! might handles that give access to their solution values. I think this was discussed before somewhere.

[ulfworsoe]

I think I have fixed this now. I have added a mapping from user indexes to native indexes for both variables and constraints. MosekMathProgModel now contains varmap and conmap. Auxiliary vars/cons do not appear in varmap and conmap.

[joehuchette]

Hmm this seems to break my code using the interface. The only bit I dropped down to the Mosek level was in a Mosek.getnumvar(m.internalModel.task) call, where I should really be using numvar(m.internalModel). It only seems to be broken when using SOCP constraints, SDP and LP stuff seem fine.

[mlubin]

Mapping between user indices and native indices seems like a reasonable approach to me.

[joehuchette]

Yeah, it's probably the way to go, especially when the number of variables isn't extremely large.

[joehuchette]

Minimal working example: adding two SOCP constraints and two new variables is enough to break it. JuMP model:

using JuMP, Mosek

m = Model()
@defMatrixVar(m,x[1])
@defMatrixVar(m,y[1])
addConstraint(m, norm(x) <= 1)
addConstraint(m, norm(y) <= 1)
setObjective(m, :Max, x[1]+y[1])
solve(m) # :Unbounded

Corresponding MathProgBase interface calls:

using Mosek, MathProgBase
m = MathProgBase.model(Mosek.MosekSolver())

f = [1.0,1.0]
A = spzeros(0,2)
lb = [-Inf,-Inf]
ub = [Inf,Inf]

MathProgBase.loadproblem!(m, A, lb, ub, f, Float64[], Float64[], :Max)

MathProgBase.addvar!(m, 0.0, Inf, 0.0)
MathProgBase.addconstr!(m, [3], [1.0], 1.0, 1.0)
MathProgBase.addquadconstr!(m, [], [], [1,3], [1,3], [1.0,-1.0], '<', 0.0)

MathProgBase.addvar!(m, 0.0, Inf, 0.0)
MathProgBase.addconstr!(m, [4], [1.0], 1.0, 1.0)
MathProgBase.addquadconstr!(m, [], [], [2,4], [2,4], [1.0,-1.0], '<', 0.0)

MathProgBase.optimize!(m)
stat=MathProgBase.status(m) # = :Unbounded

[ulfworsoe]

Oops. I forgot to map the variable indexes in addconstr!. A fix has been committed.

[joehuchette]

Fixes that example, but I'm still having problems with a maximum ellipse, which is giving incorrect answers now. It's a bit more complicated, but I'll try to narrow down the MPB calls.

[ulfworsoe]

I have implemented the conic interface for MOSEK and fixed at least one addquadconstr!() bug. This may have fixed the bug - but I have not been able to verify this (the maximum ellipse examples does not run with the public version of JuMP).

[joehuchette]

I am testing this on the sdp2 branch of JuMP, and I'm having issues with adding a constraint of the form addConstraint(mod, t^2 <= Z[1,1] * Z[2,2]):

linearidx = Int32[]
linearval = Float64[]
quadrowidx = Int32[11,12,13]
quadcolidx = Int32[11,12,13]
quadval = [1.0,1.0,-1.0]
ptrb = [1,3,5]
ptre = [3,5,7]
subj = Int32[0,14,0,15,0,16]
MOSEK error 1203: The index value -1 occurring in argument '(asub+aptrb[0])[0]' is too small.
ERROR: MosekError(1203,"The index value -1 occurring in argument '(asub+aptrb[0])[0]' is too small.")
 in putarowslice at /Users/huchette/.julia/v0.4/Mosek/src/msk_functions.jl:2721
 in putarowslice at /Users/huchette/.julia/v0.4/Mosek/src/msk_functions.jl:2699
 in addquadconstr! at /Users/huchette/.julia/v0.4/Mosek/src/MosekLowLevelQCQO.jl:135
 in addquadconstr! at /Users/huchette/.julia/v0.4/Mosek/src/MosekLowLevelQCQO.jl:31
 in addnorm2 at /Users/huchette/.julia/v0.4/JuMP/src/sdp_solve.jl:234
 in addPrimalConstraint at /Users/huchette/.julia/v0.4/JuMP/src/sdp_solve.jl:57
 in solveSDP at /Users/huchette/.julia/v0.4/JuMP/src/sdp_solve.jl:407
 in solve at /Users/huchette/.julia/v0.4/JuMP/src/solvers.jl:15
 in include at ./boot.jl:242
 in include_from_node1 at loading.jl:128
 in process_options at ./client.jl:293
 in _start at ./client.jl:375
 in _start at /usr/local/julia/usr/lib/julia/sys.dylib
while loading /Users/huchette/.julia/v0.4/JuMP/examples/maxellipse.jl, in expression starting on line 86

It looks like x (which populates the first row of subj here) isn't getting filled with meaningful indices.

[ulfworsoe]

Ah. I see the problem: Z is an SDP variable, and I have completely forgotten to add support for SDP in addquadconstr!().

[mlubin]

This seems like an awkward way to combine SOC and SDP constraints.

[joehuchette]

Z isn't an SDP variable itself, though its components are used in SDP constraints: https://github.com/JuliaOpt/JuMP.jl/blob/sdp2/examples/maxellipse.jl#L76-L78

[ulfworsoe]

Miles: Yes, it is a bit awkward, but it just requires me to split the A matrix into SDP and non-SDP parts. I think a more natural way would be to extend the conic interface with something like

addconicconstr!(m,rowidxs,colidxs,vals,b,conetype)
addconicvar!(m,length,conetype)

[ulfworsoe]

I have fixed a bunch of problems, including handling linear and quadratic constraints separetely. Also:

• Handling of SDP variables in conic-quadratic constraints should also work, but I wasn't able to run all tests in the sdp2 branch of JuMP (some issue with getValue([X])).
• It fails some non-linear tests in the current JuMP. I'll file that as a separate issue and look at it as soon as possible.
  I'll close this issue as I think the original problem has been fixed.

[joehuchette]

@ulfworsoe I'm still seeing the same Mosek error as before, for both the maxellipse.jl and robust_uncertainty.jl examples.

[mlubin]

I think a more natural way would be to extend the conic interface with something like ...

@joehuchette, what do you think about folding the SDP interface into addconicvar! and addconicconstr! as @ulfworsoe suggested?

[joehuchette]

Yes, that seems reasonable to me.

[madeleineudell]

I think I'm seeing an error related to this: on the warmstart branch of Convex, trying

using Convex, Mosek

set_default_solver(MosekSolver())

x = Variable()
y = Variable()

p = minimize(x+y, x>=0, y>=0)
solve!(p)

y.value = 4
fix!(y)
solve!(p)

gives the error

MOSEK error 1241: Too small maximum number of variables 0 is specified. Currently, the number of variables is 3.
ERROR: MosekError(1241,"Too small maximum number of variables 0 is specified. Currently, the number of variables is 3.")
 in putmaxnumvar at /Users/madeleine/.julia/v0.3/Mosek/src/msk_functions.jl:2992
 in loadconicproblem! at /Users/madeleine/.julia/v0.3/Mosek/src/MosekConicInterface.jl:235
 in load_problem! at /Users/madeleine/.julia/v0.3/Convex/src/solution.jl:81
 in solve! at /Users/madeleine/.julia/v0.3/Convex/src/solution.jl:27

(If this problem is unrelated I'll open a separate issue.)

[ulfworsoe]

I think this is a separate issue. It appears that MOSEK has changed behaviour so putmaxnumvar cannot remove variables.

[ulfworsoe]

I have replace the putmax... with resizetask. Could you check if it fixes the problem?

[madeleineudell]

Great! That fixes it.

[ulfworsoe]

Closing