I have to solve a basic parametric MILP.:
\begin{equation} \begin{split} \mbox{maximize} & \quad \alpha \\ \mbox{s.t.} & \quad \alpha \leq - \lambda_j(\psi)(r_{*j} - z_j(x)), \,\, j=1,\dots,d \\ & \quad x\in \mathbb X \subset\mathbb R^{n}, \ \alpha \ge 0 \end{split} \end{equation}
The parameter is a weight value \lambda_j(\psi) of dimension d from 2 to 6.
The number of weights to consider varies between 100 and 10000.
I am repeating the optimization 20 times.
To avoid to rebuild the model from scratch, I modify it at each iteration, in deleting and adding d constraints concerned by a weight.
I implemented it, and I solved instances respectively with GLPK, HiGHS, Gurobi, and CPLEX. The results collected are valid (for the 4 MIP solvers) but the time consumed are drastically differents… in favor of GLPK.
Ok, my comparison is unfair; the preprocessing is not disabled for the others solvers…but the elapsed times are very amazing:
first run (with 1000 weights) to compile:
GLPK.Optimizer | time(s)= 2.8956220149993896
with 1000 weights:
GLPK.Optimizer | time(s)= 0.2676060199737549
GLPK.Optimizer | time(s)= 0.260451078414917
GLPK.Optimizer | time(s)= 0.2592289447784424
HiGHS.Optimizer | time(s)= 15.928501844406128
HiGHS.Optimizer | time(s)= 15.65535306930542
HiGHS.Optimizer | time(s)= 15.687047004699707
Gurobi.Optimizer | time(s)= 1.675018072128296
Gurobi.Optimizer | time(s)= 1.4513440132141113
Gurobi.Optimizer | time(s)= 1.4410531520843506
CPLEX.Optimizer | time(s)= 105.86109781265259
CPLEX.Optimizer | time(s)= 118.24663805961609
CPLEX.Optimizer | time(s)= 116.49558091163635
with 10000 weights:
GLPK.Optimizer | time(s)= 12.330639839172363
HiGHS.Optimizer | time(s)= 1553.121855020523
Gurobi.Optimizer | time(s)= 154.71333193778992
too long with CPLEX
I am wondering about the cause of these differences; is it
- a consequence of by my approach for modifying the model?
- inhered to the MIP solver activity?
- something internal to the packages bridging MOI and the corresponding solver?
Any opinion/explanation/advice?
Here the code (extracted from my implementation and simplified to illustrate this) used to report here the elapsed times:
using JuMP, GLPK, HiGHS, Gurobi, CPLEX
function L( solver::DataType,
p::Matrix{Int64}, w::Vector{Int64}, c::Int64,
rp::Vector{Int64}, λ_ψ::Vector{Vector{Float64}}
)
L = (Float64)[] # a series of L(S,r*,ψ) values
d,n = size(p) # number of objectives and number of variables
# define a JuMP model
model = Model(solver)
set_silent(model)
# setup X according (1)
@variable(model, x[1:n], Bin) # the n binary variables of 01UKP
@constraint(model, sum(w[i] * x[i] for i in 1:n) ≤ c) # the constraint of 01UKP
# declare the objective functions
@expression(model, z[k=1:d], sum(p[k,j] * x[j] for j in 1:n)) # the d objectives of 01UKP
# ------------------------------------------------------------------------------
# Concerning the Chebychev part of the model
@variable(model, α ≥ 0)
@objective(model, Max, α)
for i in 1:length(λ_ψ)
# add the d constraints for a given λ_ψ to the model
@constraint(model, con[k=1:d], α ≤ - λ_ψ[i][k] * (rp[k] - z[k]) )
optimize!(model)
@assert is_solved_and_feasible(model) "Error: optimal solution not found"
push!(L,objective_value(model))
for k=1:d
delete(model, con[k])
end
unregister(model, :con)
end
return L
end
p = [ 13 10 3 16 12 11 1 9 19 13 ; 1 10 3 13 12 19 16 13 11 9 ]
w = [ 4, 4, 3, 5, 5, 3, 2, 3, 5, 4 ]
c = 19
rp = [40,40]
λ_ψ = fill([1.0238499551062132, 4.66018735258353],1000)
solver = GLPK.Optimizer
#solver = HiGHS.Optimizer
#solver = Gurobi.Optimizer
#solver = CPLEX.Optimizer
start = time()
for _ in 1:20
L(solver, p, w, c, rp, λ_ψ) # => 25.596248877655327 for this weight
end
t_elapsed = time() - start
println(solver, " | time(s)= ", t_elapsed)
My configuration is the following:
- macBook pro M2
- julia 1.10.5
and packages are:
- JuMP v1.26.0
- GLPK v1.2.1
- HiGHS v1.18.1
- Gurobi v1.7.4
- CPLEX v1.1.1