Creating a PSDCone Constraint in COSMO

Due to memory limitation issues, I am unable to use JuMP to create variables and the model I would like to use.

As such, I am trying to use COSMO directly. I am trying to create a PSDCone constraint of the form (A - diag(y)), where y is my vector of decision variables, and A is an nxn matrix.

Any help on how to do this would be greatly appreciated.

Thank you!

Hi, you can try JuMP’s direct mode to reduce memory allocation:

model = direct_model(MOI.instantiate(COSMO.Optimizer, with_bridge_type=Float64))

You can also reduct it further by using MOI directly

model = direct_model(MOI.instantiate(COSMO.Optimizer, with_bridge_type=Flaot64))

However, if you have allocation issue it’s probably because your algebra with JuMP’s expressions is allocating more than necessary.
JuMP automatically try to avoid these allocations with GitHub - jump-dev/MutableArithmetics.jl: Interface for arithmetics on mutable types in Julia, e.g., many linear algebra operations are rewritten to make this faster but it’s not done for all of them.
I would recommend trying to pin point the bottleneck with the JuMP approach and share it here so that we can see what the issue is rather than trying to use COSMO directly. It’s probably going to be more time-efficient for you than going down to a lower level interface.

Sure thing!

I am trying to follow the solution here: Out of memory when constructing large sparse SDP in JuMP - #3 by blegat

But, when I try to run the code using the solution you provided, I am still running into memory issues (100% memory on my system is used). Hence, I was investigating whether it would be possible to do this directly in COSMO and avoid creating the n^2 / 2 affine expressions all together.

Do you have a reproducible example?

I do not at the moment. But, here is what I am hoping to translate from JuMP into COSMO.

n = 4000
b = ones(n)
p = 0.01
C = sprandn(Float64, n, n, p)  
C = triu(Symmetric(sprand(n,n,p)) )
model = JuMP.Model(COSMO.Optimizer)
@variable(model, y[1:n])
@constraint(model, Symmetric(C - spdiagm(y)) >= 0, PSDCone())
@objective(model, Max, b'y)
JuMP.optimize!(model)

The issue is that as n becomes large, the line

@constraint(model, Symmetric(C - spdiagm(y)) >= 0, PSDCone())

crashes my computer. I believe this is due to the issues stated in the previous thread. As such, I am thinking that I can try to recreate this directly in COSMO, and avoid memory issues with JuMP. However, I am having difficulties figuring out how to translate the constraint above into something that COSMO can look at directly. This is because (to my understanding), COSMO requires that I provide the data matrices, and it handles the variables (see Getting Started · COSMO.jl).

Thanks!

I think your model translates to something like

using SparseArrays, LinearAlgebra, COSMO
n = 4_000
b = ones(n)
p = 0.01
C = LinearAlgebra.Symmetric(SparseArrays.sprand(n, n, p))
b1 = [C[i, j] for i in 1:n for j in 1:i]
A1 = SparseArrays.spzeros(length(b1), n)
k = 0
for i in 1:n 
    for j in 1:i
        k += 1
        if i == j
            A1[k, i] = -1.0
        end
    end
end
model = COSMO.Model()
COSMO.assemble!(
    model, 
    SparseArrays.spzeros(n, n), 
    -b, 
    [COSMO.Constraint(A1, b1, COSMO.PsdConeTriangle)],
)
results = COSMO.optimize!(model)

But it’s very large, so that’s still going to be very difficult, if not impossible, to solve. You should consider alternative ways to model/solve your problem.

Thank you so much for this portion!

I’ll definitely do some thinking to see if there are other techniques to reduce the problem size.

Thanks!

1 Like

You can use the MOI interface. Since it is using a VectorAffineFunction storing all the terms instead of an array of JuMP.AffExpr, you won’t have the issue with allocation of a quadratic number of objects. The MOI VectorAffineFunction will just allocate two vectors. The vector of terms only contain the diagonal elements