Incorrect SDPA file output

[kibaekkim]

I think the SDPA file is not correctly written. Here is an example:

using JuMP

m = Model()
@variable(m, x[1:3,1:3], PSD)
@variable(m, 12.34 <= y[1:2] <= 123.4)
@objective(m, Min, x[1,1] + 4*x[1,2] + 6*x[1,3] + 2*x[2,1] + 9*x[2,2] + 3*x[3,1] + 7*x[3,3])
@constraint(m, x[1,1] + 2*x[1,3] + 3*x[2,2] + 14*x[2,3] + 5*x[3,3] + y[1] == 11)
@constraint(m, 4*x[1,2] + 16*x[1,3] + 6*x[2,2] + 4*x[3,3] + y[2] == 19)

JuMP.write_to_file(m, "ex.dat-s", format=MOI.FileFormats.FORMAT_SDPA)

The file ex.dat-s should have

8
9
-1 -1 -1 -1 -1 -1 -1 -1 3
1.0 6.0 9.0 9.0 0.0 7.0 0.0 0.0
0 1 1 1 -12.34
7 1 1 1 1.0
0 2 1 1 -12.34
8 2 1 1 1.0
0 3 1 1 123.4
7 3 1 1 -1.0
0 4 1 1 123.4
8 4 1 1 -1.0
0 5 1 1 -11.0
1 5 1 1 1.0
3 5 1 1 3.0
4 5 1 1 2.0
5 5 1 1 14.0
6 5 1 1 5.0
7 5 1 1 1.0
0 6 1 1 11.0
1 6 1 1 -1.0
3 6 1 1 -3.0
4 6 1 1 -2.0
5 6 1 1 -14.0
6 6 1 1 -5.0
7 6 1 1 -1.0
0 7 1 1 -19.0
2 7 1 1 4.0
3 7 1 1 6.0
4 7 1 1 16.0
6 7 1 1 4.0
8 7 1 1 1.0
0 8 1 1 19.0
2 8 1 1 -4.0
3 8 1 1 -6.0
4 8 1 1 -16.0
6 8 1 1 -4.0
8 8 1 1 -1.0
1 9 1 1 1.0
2 9 1 2 0.5
3 9 2 2 1.0
4 9 1 3 0.5
5 9 2 3 0.5
6 9 3 3 1.0

There are some (not all) findings:

• The objective coefficients are wrong.
• The lower bounds of y variable are negated.

[odow]

PSD enforces x to be symmetric. So your objective is correct. The lower bounds are also correct because they are represented as [x - l] in NonNegatives(1). The upper bounds are represented as [u - x] in Nonnegatives(1).

You can check the file is read/written correctly by round-tripping it.

using JuMP
using SCS

m = Model(SCS.Optimizer)
@variable(m, x[1:3,1:3], PSD)
@variable(m, 12.34 <= y[1:2] <= 123.4)
@objective(m, Min, x[1,1] + 4*x[1,2] + 6*x[1,3] + 2*x[2,1] + 9*x[2,2] + 3*x[3,1] + 7*x[3,3])
@constraint(m, x[1,1] + 2*x[1,3] + 3*x[2,2] + 14*x[2,3] + 5*x[3,3] + y[1] == 11)
@constraint(m, 4*x[1,2] + 16*x[1,3] + 6*x[2,2] + 4*x[3,3] + y[2] == 19)
write_to_file(m, "mod.sdpa")
optimize!(m)

m2 = read_from_file("mod.sdpa")
set_optimizer(m2, SCS.Optimizer)
optimize!(m2)

isapprox(objective_value(m), objective_value(m2), atol=1e-6) # true

[kibaekkim]

Interesting.. You are right. I did not see that there are symmetric elements in the objective function.

I tested the output file with the solvers available in NEOS, as well as scipsdp (with Mosek). They all resulted in infeasibility. And, I still see that. I just tested the output file with DSDP solver. It returned the same result: infeasibility.

But, in fact, I have not tried reading again and solving it in JuMP, except that I just verified your script... There are a couple possibilities:

1. The SDPA reader in NEOS and scip-sdp are buggy.
2. The SDPA format written in JuMP is not fully compatible to others, although it can still be read by JuMP again.

Any idea?

[odow]

Probably 2. Can you make an SDPA file using a different modeling language/solver?

[kibaekkim]

Good idea. This is truss1.dat-s, an example from DSDP5.8 package.

6
7
2 2 2 2 2 2 1
-1.0 -0.0 -2.0 -0.0 -0.0 -0.0
0 7 1 1 -1.0
1 1 2 2 -1.0
1 2 2 2 -1.0
1 3 2 2 -1.0
1 4 2 2 -1.0
1 5 2 2 -1.0
1 6 2 2 -1.0
2 2 1 2 -1.000000999999999918
2 5 1 2 -5.0e-01
2 6 1 2 3.240558000000000158e-07
3 2 1 2 -7.137334999999999900e-08
3 5 1 2 4.999998999999999416e-01
3 6 1 2 1.0
4 3 1 2 -5.000000999999999474e-01
4 4 1 2 -9.999993999999998717e-01
4 6 1 2 -3.240558000000000158e-07
5 3 1 2 -5.000002000000000058e-01
5 4 1 2 2.486843000000000082e-07
5 6 1 2 -1.0
6 1 1 1 -1.0
6 2 1 1 -1.0
6 3 1 1 -1.0
6 4 1 1 -1.0
6 5 1 1 -1.0
6 6 1 1 -1.0
6 7 1 1 1.0

DSDP found the objective value of 9, whereas SCS called from JuMP claimed Infeasible.

[odow]

We likely have implemented the spec wrong then, but done so consistently for reading and writing.

[odow]

I've never really looked into SDPA. Here's another bug: https://github.com/JuliaOpt/MathOptInterface.jl/issues/1073

@blegat did you have issues? Are there different SDPA variants? e.g., between SDPA and Mosek?

[blegat]

After https://github.com/JuliaOpt/MathOptInterface.jl/pull/1075, I get:

using JuMP
model = JuMP.read_from_file("truss1.dat-s", format=MOI.FileFormats.FORMAT_SDPA)
using SCS
set_optimizer(model, SCS.Optimizer)
optimize!(model)
termination_status(model) # OPTIMAL::TerminationStatusCode = 1
objective_value(model) # -2.2499990764908504

So it seems to be improved but the objective value is not yet the 9 found by DSDP. The model read is

julia> for (i, x) in enumerate(all_variables(model))
       set_name(x, "x$i")
       end

julia> print(model)
Min -x1 - 2 x3
Subject to
 [-x6, 0, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [-x6, -2.000002 x2 - 1.427467e-7 x3, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [-x6, -1.0000002 x4 - 1.0000004 x5, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [-x6, -1.9999987999999997 x4 + 4.973686e-7 x5, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [-x6, -x2 + 0.9999997999999999 x3, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [-x6, 6.481116e-7 x2 + 2 x3 - 6.481116e-7 x4 - 2 x5, -x1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)
 [x6 + 1] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(1)

@kibaekkim Would you expect something different ? EDIT It seems the off-diagonal entries are scaled by two but they shouldn't. For instance, -2.000002 x2 should be -1.000001 x2.

[kibaekkim]

@blegat Unfortunately, I was not able to find the formulation for the problem instance. I took the instance file from DSDP package. Maybe you can try this https://github.com/JuliaOpt/MathOptInterface.jl/issues/1071#issuecomment-614979535?

[blegat]

After https://github.com/JuliaOpt/MathOptInterface.jl/pull/1076, I get -9 :tada: