Delay between call to optimize! and COSMO solving problem

[junglegobs]

I noticed that for large(ish) problems there is a significant delay between calling optimize!(m) and COSMO actually running. See example below:

using COSMO, JuMP, Cbc

I = 1:10_000 # Gets significantly slower at 10_000
m = Model(COSMO.Optimizer)
x = @variable(m, [i=I], lower_bound=0.0)
con1 = @constraint(m, [i=I], x[i]*rand() >= rand())
obj = @objective(m, Min, sum(rand()*x[i] for i in I))
@time optimize!(m)

m = Model(Cbc.Optimizer)
x = @variable(m, [i=I], lower_bound=0.0)
con1 = @constraint(m, [i=I], x[i]*rand() >= rand())
obj = @objective(m, Min, sum(rand()*x[i] for i in I))
@time optimize!(m)

REPL output for COSMO:

Runtime: 0.097s (96.76ms)
42.261334 seconds (1.95 M allocations: 24.543 GiB, 5.02% gc time)

REPL output for Cbc:

Total time (CPU seconds):       0.01   (Wallclock seconds):       0.01
0.065991 seconds (191.12 k allocations: 16.796 MiB, 25.94% gc time)

This output is from the third time running the above script. This delay strikes me as excessive, especially when compared to Cbc. Am I doing something wrong? Is this a known issue?

[migarstka]

This must be an issue with the way the problem is transformed from JuMP/MathOptInterface to COSMO, so most likely in MOIWrapper.jl.

Thanks for the MWE, I will profile the code and see if there is an obvious bottleneck.

Edit: Does the problem still exist if you would assemble the constraint in matrix-vector-form and just create one constraint via @constraint ?

[junglegobs]

Yes. Reworked example:

using COSMO, JuMP, Cbc, LinearAlgebra

n = 10_000 # Gets significantly slower at 10_000
A = diagm(rand(n))
b = rand(n)
c = rand(n)

m = Model(COSMO.Optimizer)
@variable(m, x[i=1:n], lower_bound=0.0)
con1 = @constraint(m, A*x .>= b)
obj = @objective(m, Min, c'*x)

@time optimize!(m)

m = Model(Cbc.Optimizer)
@variable(m, x[i=1:n], lower_bound=0.0)
con1 = @constraint(m, A*x .>= b)
obj = @objective(m, Min, c'*x)

@time optimize!(m)

Outputs:

# COSMO:
Runtime: 0.1s (100.27ms)
 44.785131 seconds (1.95 M allocations: 24.543 GiB, 3.22% gc time)

# Cbc
Total time (CPU seconds):       0.01   (Wallclock seconds):       0.01
  0.053005 seconds (191.12 k allocations: 16.796 MiB)

[migarstka]

The issue seems to be that the way the problem is created in JuMP results in 2n scalar constraints which have to be transformed back into COSMOs internal constraint type, i.e. 1 vector constraint for all the inequality constraints.

A quick fix for now would be to pass the problem with just 2 vector constraints:

using COSMO, JuMP, Cbc, LinearAlgebra
using MathOptInterface
n = 10_000 # Gets significantly slower at 10_000
A = diagm(rand(n))
b = rand(n)
c = rand(n)

m = Model(COSMO.Optimizer)
@variable(m, x[i=1:n]);
con2 = @constraint(m, x in MathOptInterface.Nonnegatives(n)); # x >= 0
con1 = @constraint(m, A*x - b in MathOptInterface.Nonnegatives(n)); # Ax - b >= 0 <-> Ax >= b
obj = @objective(m, Min, c'*x);
@time optimize!(m)

[junglegobs]

Thanks for the speedy reply and quick fix! Out of interest is this a particularly difficult / annoying issue to solve? (Just curious if this is likely to be addressed soon).

Edit: also because it makes much more sense for me to just use Gurobi for now instead of spending half a day converting my problem into linear algebraic form. I would ideally like to use COSMO however, so that my code can be easily replicated.

[migarstka]

I'll look into it. Shouldn't take very long to fix.

[migarstka]

99% of the time is spent on merging all the scalar constraints xi >=0 into one vector constraint (here).

@odow @blegat Is there any performance penalty of removing support for all the scalar sets in the solver and let bridges do the conversion instead?

[blegat]

We had similar performance issues with SCS that we fixed during 2020. See jump-dev/SCS.jl#181 and jump-dev/SCS.jl#186 and PR mentioned there. Having a quick look, https://github.com/oxfordcontrol/COSMO.jl/blob/master/src/MOIWrapper.jl#L456 will allocate a lot for each constraint and https://github.com/oxfordcontrol/COSMO.jl/blob/master/src/MOIWrapper.jl#L177 is slow as it is type unstable with F and S, you should use a function barrier. We created a conic form in MatrixOptInterface that SCS will use soon. We are trying to make as many solvers as possible use it to make sure they all use the same copy that is precisely benchmarked and optimized for every case. Could COSMO use it if we add support for a quadratic objective?

[migarstka]

I don't see why not. COSMO supports all the constraints that SCS supports + quadratic objective + intervall constraints

[migarstka]

I mostly fixed the slow constraint set merging issue: #129

[migarstka]

Fixed in new release v0.7.9