Named Constraints Error with Gurobi

[pulsipher]

I am trying to implement the FAC example from the JuMP paper, but using POI to update the M value. However, this is incurring some strange errors. (I am using JuMP v1.1 and POI v0.4)

using JuMP, ParametricOptInterface, Gurobi
const POI = ParametricOptInterface

function create_model(G, F, M_start)
    # Initialize the model object
    model = Model(() -> POI.Optimizer(Gurobi.Optimizer()))

    # Add the variables
    @variable(model, d)
    @variable(model, 0 ≤ y[1:F, 1:2] ≤ 1)
    @variable(model, z[0:G, 0:G, 1:F], Bin)
    @variable(model, 0 ≤ s[0:G, 0:G, 1:F])
    @variable(model, r[0:G, 0:G, 1:F, 1:2])
    @variable(model, Mp in POI.Parameter(M_start))
    
    # Set the objective
    @objective(model, Min, d)
    
    # Add the constraints
    @constraint(model, [i ∈ 0:G, j ∈ 0:G], sum(z[i,j,f] for f ∈ 1:F) == 1)
    @constraint(model, con[i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], s[i,j,f] == d + Mp*(1 - z[i,j,f]))
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,1] == i/G - y[f,1])
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,2] == j/G - y[f,2])
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,1]^2 + r[i,j,f,2]^2 ≤ s[i,j,f]^2)

    # Return the model
    return model
end

m = create_model(4, 4, 4.0)
optimize!(m) 

ERROR: The index MathOptInterface.ConstraintIndex{MathOptInterface.ScalarQuadraticFunction{Float64}, MathOptInterface.EqualTo{Float64}}(1) is invalid. Note that an index becomes invalid after it has been deleted.
Stacktrace:
  [1] _info
    @ C:\Users\bbgui\.julia\packages\Gurobi\cnAYk\src\MOI_wrapper\MOI_wrapper.jl:2332 [inlined]
  [2] set
    @ C:\Users\bbgui\.julia\packages\Gurobi\cnAYk\src\MOI_wrapper\MOI_wrapper.jl:2487 [inlined]
  [3] set
    @ C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\ParametricOptInterface.jl:446 [inlined]
  [4] set
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Bridges\bridge_optimizer.jl:1420 [inlined]
  [5] _pass_attribute(dest::MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, src::MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}, index_map::MathOptInterface.Utilities.IndexMap, cis_src::Vector{MathOptInterface.ConstraintIndex{MathOptInterface.ScalarQuadraticFunction{Float64}, MathOptInterface.EqualTo{Float64}}}, attr::MathOptInterface.ConstraintName)
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\copy.jl:134
  [6] pass_attributes(dest::MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, src::MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}, index_map::MathOptInterface.Utilities.IndexMap, cis_src::Vector{MathOptInterface.ConstraintIndex{MathOptInterface.ScalarQuadraticFunction{Float64}, MathOptInterface.EqualTo{Float64}}})
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\copy.jl:119
  [7] _pass_constraints(dest::MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, src::MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}, index_map::MathOptInterface.Utilities.IndexMap, variable_constraints_not_added::Vector{Any})
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\copy.jl:313
  [8] default_copy_to(dest::MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, src::MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}})
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\copy.jl:481
  [9] #copy_to#7
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Bridges\bridge_optimizer.jl:421 [inlined]
 [10] copy_to
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Bridges\bridge_optimizer.jl:421 [inlined]
 [11] optimize!
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\MathOptInterface.jl:80 [inlined]
 [12] optimize!(m::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}})
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\cachingoptimizer.jl:310
 [13] optimize!(model::Model; ignore_optimize_hook::Bool, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ JuMP C:\Users\bbgui\.julia\packages\JuMP\0C6kd\src\optimizer_interface.jl:161
 [14] optimize!(model::Model)
    @ JuMP C:\Users\bbgui\.julia\packages\JuMP\0C6kd\src\optimizer_interface.jl:143
 [15] top-level scope
    @ REPL[13]:1

Strangely, if I remove the constraint name a different error is incurred:

using JuMP, ParametricOptInterface, Gurobi
const POI = ParametricOptInterface

function create_model(G, F, M_start)
    # Initialize the model object
    model = Model(() -> POI.Optimizer(Gurobi.Optimizer()))

    # Add the variables
    @variable(model, d)
    @variable(model, 0 ≤ y[1:F, 1:2] ≤ 1)
    @variable(model, z[0:G, 0:G, 1:F], Bin)
    @variable(model, 0 ≤ s[0:G, 0:G, 1:F])
    @variable(model, r[0:G, 0:G, 1:F, 1:2])
    @variable(model, Mp in POI.Parameter(M_start))
    
    # Set the objective
    @objective(model, Min, d)
    
    # Add the constraints
    @constraint(model, [i ∈ 0:G, j ∈ 0:G], sum(z[i,j,f] for f ∈ 1:F) == 1)
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], s[i,j,f] == d + Mp*(1 - z[i,j,f]))
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,1] == i/G - y[f,1])
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,2] == j/G - y[f,2])
    @constraint(model, [i ∈ 0:G, j ∈ 0:G, f ∈ 1:F], r[i,j,f,1]^2 + r[i,j,f,2]^2 ≤ s[i,j,f]^2)

    # Return the model
    return model
end

m = create_model(2, 2, 4.0)
optimize!(m) 

Gurobi Optimizer version 9.5.0 build v9.5.0rc5 (win64)
Thread count: 8 physical cores, 16 logical processors, using up to 16 threads
Optimize a model with 63 rows, 77 columns and 144 nonzeros
Model fingerprint: 0x6c4035d3
Model has 18 quadratic constraints
Variable types: 59 continuous, 18 integer (18 binary)
Coefficient statistics:
  Matrix range     [1e+00, 4e+00]
  QMatrix range    [1e+00, 1e+00]
  Objective range  [1e+00, 1e+00]
  Bounds range     [1e+00, 1e+00]
  RHS range        [5e-01, 4e+00]
Presolve removed 21 rows and 21 columns
Presolve time: 0.00s
Presolved: 42 rows, 56 columns, 102 nonzeros
Presolved model has 18 quadratic constraint(s)
Variable types: 47 continuous, 9 integer (9 binary)

Root relaxation: objective 0.000000e+00, 0 iterations, 0.00 seconds (0.00 work units)

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

H    0     0                       0.7071069    0.00000   100%     -    0s
*    0     0               0       0.7071069    0.00000   100%     -    0s
H    0     0                       0.7071068    0.00000   100%     -    0s
*    0     0               0       0.7071068    0.00000   100%     -    0s
H    0     0                       0.7071068    0.00000   100%     -    0s
     0     0    0.00000    0    6    0.70711    0.00000   100%     -    0s
H    0     0                       0.5590172    0.00000   100%     -    0s
     0     0    0.00000    0    8    0.55902    0.00000   100%     -    0s
     0     0    0.00000    0    8    0.55902    0.00000   100%     -    0s
H    0     0                       0.5590168    0.00000   100%     -    0s
     0     2    0.00000    0    8    0.55902    0.00000   100%     -    0s
H    3     8                       0.5590168    0.11881  78.7%  11.3    0s
*   35     4               4       0.5590168    0.54668  2.21%   5.0    0s
*   36     4               4       0.5590165    0.54668  2.21%   4.9    0s
*   39     0               4       0.5590161    0.54668  2.21%   4.7    0s

Cutting planes:
  Gomory: 1
  Implied bound: 2
  MIR: 11

Explored 45 nodes (252 simplex iterations) in 0.07 seconds (0.02 work units)
Thread count was 16 (of 16 available processors)

Solution count 8: 0.559016 0.559017 0.559017 ... 0.707107

Optimal solution found (tolerance 1.00e-04)
Warning: max constraint violation (1.9584e-06) exceeds tolerance
Best objective 5.590161137321e-01, best bound 5.590161137321e-01, gap 0.0000%

User-callback calls 473, time in user-callback 0.00 sec
ERROR: Gurobi Error 10005: Unable to retrieve attribute 'Pi'
Stacktrace:
  [1] _check_ret
    @ C:\Users\bbgui\.julia\packages\Gurobi\cnAYk\src\MOI_wrapper\MOI_wrapper.jl:317 [inlined]
  [2] get(model::Gurobi.Optimizer, attr::MathOptInterface.ConstraintDual, c::MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.EqualTo{Float64}})
    @ Gurobi C:\Users\bbgui\.julia\packages\Gurobi\cnAYk\src\MOI_wrapper\MOI_wrapper.jl:3099
  [3] calculate_parameters_in_ci!(param_dual_cum_sum::Vector{Float64}, optimizer::Gurobi.Optimizer, param_array::Vector{MathOptInterface.ScalarAffineTerm{Float64}}, ci::MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.EqualTo{Float64}})
    @ ParametricOptInterface C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\duals.jl:106
  [4] update_duals_with_quadratic_constraint_cache!(param_dual_cum_sum::Vector{Float64}, model::ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}, quadratic_constraint_cache_pc_inner::MathOptInterface.Utilities.DoubleDicts.DoubleDictInner{MathOptInterface.ScalarQuadraticFunction{Float64}, MathOptInterface.EqualTo{Float64}, Vector{MathOptInterface.ScalarAffineTerm{Float64}}})
    @ ParametricOptInterface C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\duals.jl:90
  [5] update_duals_with_quadratic_constraint_cache!(param_dual_cum_sum::Vector{Float64}, model::ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer})
    @ ParametricOptInterface C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\duals.jl:69
  [6] calculate_dual_of_parameters(model::ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer})
    @ ParametricOptInterface C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\duals.jl:8
  [7] optimize!
    @ C:\Users\bbgui\.julia\packages\ParametricOptInterface\3MNVD\src\ParametricOptInterface.jl:1280 [inlined]
  [8] optimize!
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Bridges\bridge_optimizer.jl:348 [inlined]
  [9] optimize!
    @ C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\MathOptInterface.jl:81 [inlined]
 [10] optimize!(m::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{ParametricOptInterface.Optimizer{Float64, Gurobi.Optimizer}}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}})
    @ MathOptInterface.Utilities C:\Users\bbgui\.julia\packages\MathOptInterface\ElZ45\src\Utilities\cachingoptimizer.jl:310
 [11] optimize!(model::Model; ignore_optimize_hook::Bool, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ JuMP C:\Users\bbgui\.julia\packages\JuMP\0C6kd\src\optimizer_interface.jl:161
 [12] optimize!(model::Model)
    @ JuMP C:\Users\bbgui\.julia\packages\JuMP\0C6kd\src\optimizer_interface.jl:143
 [13] top-level scope
    @ REPL[20]:1

I could be doing something wrong, but I suspect this is a bug.

[guilhermebodin]

Hi @pulsipher I dont have Gurobi on my computer right now. I tried it on Xpress and it worked. I am wondering if this is a problem with POI itself or some integration with solvers. I am checking with some colleagues that might have Gurobi.

[rafabench]

Hi @pulsipher. The second error occurs because the ParametricOptInterface tries to evaluate the dual of the constraints and it gives an error because it's a MIP. You should change the line:

model = Model(() -> POI.Optimizer(Gurobi.Optimizer()))

To

model = Model(() -> POI.Optimizer(Gurobi.Optimizer(),evaluate_duals=false))

It worked for Gurobi and CPLEX in my computer.
The first error we're still looking into it.

[guilhermebodin]

Some reference to the docs about the POI options. https://jump.dev/ParametricOptInterface.jl/dev/reference/#ParametricOptInterface.Optimizer

[pulsipher]

Thanks for all the help guys, this helps me get the code working!

[guilhermebodin]

Both cases worked on Xpress on my computer for this environment

[4076af6c] JuMP v1.1.1
[b8f27783] MathOptInterface v1.5.0
[9e70acf3] Xpress v0.15.1

I am closing the issue.

[pulsipher]

So this seems to be a bug when using POI and Gurobi jointly with named constraints then, it doesn't happen with POI and other solvers and nor does it occur with just Gurobi.

[guilhermebodin]

Well, we can re-open, my main problem to debug this is that I don't have any Gurobi licenses, I can wait until @rafabench has a little more time to dig into the specifics of this integration