JuMP semidefinite programming - formulating via constraints versus variables

I have the following two model definitions in JuMP. They are formally equivalent. Will one of them be more efficient than the other?

Version 1:

a = [1.0 0.0; 0.0 1.0]
b = [1.0 0.0; 0.0 -1.0]
@variable(model, x)
@constraint(model, a + x * b >= 0, PSDCone())

Version 2:

a = [1.0 0.0; 0.0 1.0]
b = [1.0 0.0; 0.0 -1.0]
@variable(model, x)
@variable(model, X[1:2, 1:2], PSD)
for i in 1:2, j in 1:2
    @constraint(model, a[i,j] + x * b[i,j] == X[i,j])

After defining the model in either manner, I will run

@objective(model, Max, x)

The performance of the code optimize!(model) is what I’m concerned about.

P.S. Of course the above problem is trivial. In practice I’ll deal with matrices of size N much bigger than 2, but still just one variable x or at most a few variables.

In the documentation for Semidefinite constraints you can read

Where possible, prefer constructing a matrix of Semidefinite variables using the @variable macro, rather than adding a constraint like @constraint(model, X >= 0, PSDCone()) . In some solvers, adding the constraint via @constraint is less efficient, and can result in additional intermediate variables and constraints being added to the model.

But in my problem, declaring a semidefinite matrix variable (version 2) means that I then have to add N^2 constraints for the value of each matrix component. Could this become a source of inefficiency?

Will one of them be more efficient than the other?

It likely depends on the solver. But since version 1 is the easiest to read, you should probably just stick with that unless it is too slow.

In more detail, a + x * b >= 0, PSDCone() is a VectorAffineFunction-in-PositiveSemidefiniteConeSquare. If solvers support that it will be efficient, otherwise we’d rewrite to version 2 anyway.

The documentation is really about the difference between @constraint(model, X >= 0, PSDCone()) and @variable(model, X[1:n, 1:n], PSD), i.e., a variable instead of an affine function.

1 Like

For large matrix sizes, I get much better accuracy if I use Dualization.jl (at least with the Mosek solver), probably due to how my problem is formulated, following suggestions in another thread.

1 Like