Different result each time problem is solved

[baggepinnen]

This issue was originally opened as jump-dev/Hypatia.jl#801 (comment)

The following non-random problem produces very different results each time the code is run. The problem does not appear to be within the solver Hypatia since if only the part after and including prob = Convex.minimize(γ) is run, the results are the same each time.

using Hypatia, Convex, LinearAlgebra
(A, B1, B2, C1, C2, D11, D12, D21, D22) = ([0.0 10.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; -10.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; -4.0 -40.0 -0.01 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 40.0 -4.0 0.0 -0.01 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 -3535.46761115006 3535.46761115006 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 -3535.46761115006 -3535.46761115006 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 -3535.46761115006 3535.46761115006 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 -3535.46761115006 -3535.46761115006 0.0 0.0; 512.0 5120.0 0.0 0.0 0.0 0.0 0.0 0.0 -2309.2856038177697 0.0; -5120.0 512.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -2309.2856038177697], [0.0 0.0; 0.0 0.0; 4.0 0.0; 0.0 4.0; 0.0 0.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 1.0; 0.0 0.0; 0.0 0.0; 0.0 0.0; 10725.385484097349 0.0; 0.0 0.0; 0.0 10725.385484097349; 0.0 0.0; 0.0 0.0], [-0.1 -1.0 2.49975 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 1.0 -0.1 0.0 2.49975 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 -46.28786141317 -6546.0921383757795 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 -46.28786141317 -6546.0921383757795 0.0 0.0; 100.0 1000.0 0.0 0.0 0.0 0.0 0.0 0.0 -448.77718277317985 0.0; -1000.0 100.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -448.77718277317985], [-1.0 -10.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 10.0 -1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0], [0.1 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.0 0.0; 10000.0 0.0; 0.0 10000.0; 0.0 0.0; 0.0 0.0], [1.0 0.0; 0.0 1.0], [0.0 0.0; 0.0 0.0])

ϵ = 1e-3

n = size(A,1)
ny1 = size(C1, 1)
nu1 = size(B1, 2)
ny2 = size(C2, 1)
nu2 = size(B2, 2)
γ  = Convex.Variable()
X1 = Convex.Variable(n,n)
Y1 = Convex.Variable(n,n)
An = Convex.Variable(n,n)
Bn = Convex.Variable(n, ny2)
Cn = Convex.Variable(nu2, n)
Dn = Convex.Variable(nu2, ny2)

M11 = A*Y1 + Y1*A' + B2*Cn + Cn'*B2'
M21 = A' + An + (B2*Dn*C2)'
M31 = (B1 + B2*Dn*D21)'
M41 = C1*Y1 + D12*Cn
M22 = X1*A + A'*X1 + Bn*C2 + C2'*Bn'
M32 = (X1*B1 + Bn*D21)'
M42 = C1 + D12*Dn*C2
M43 = D11 + D12*Dn*D21
M33 = -γ*Matrix(I(nu1))
M44 = -γ*Matrix(I(ny1))
M = [
    M11 M21' M31' M41'
    M21 M22  M32' M42'
    M31 M32  M33  M43'
    M41 M42  M43  M44
]
prob = Convex.minimize(γ)
prob.constraints += γ > 0
prob.constraints += M ⪯ -ϵ*I(size(M,1))
prob.constraints += [X1 I(n); I(n) Y1] ⪰ ϵ*I(2n)

opt = Hypatia.Optimizer(rescale=false)
Convex.solve!(prob, opt; silent_solver=false, verbose=true)

Tested on Convex v0.15.0

[odow]

Does this happen with any other solver? Does it find different optimal solutions? Or is the problem being solved different each time.

Can you provide the logs of different solves?

If it's finding equivalent optimal solutions, then I don't really regard that as a bug. I could just be different order of variables/constraints due to dictionary iteration or something.

[chriscoey]

If indeed the variables/constraints are getting ordered randomly due to e.g. dictionary iteration, I think to avoid this kind of issue in future, Convex.jl should use e.g. ordered dictionaries. So that the MOI model ends up being exactly the same each time.

[ericphanson]

I thought we did use OrderedDicts; I went to double check and saw there is at least one which isn’t:

Convex.jl/src/conic_form.jl

Line 110 in 4f52b6f

const ConicConstrToConstr = Dict{ConicConstr,Constraint}

Someone could try replacing that (and any others) with an OrderedDict and seeing if that fixes it.