copy_to incorrectly maps scalar constraints

Question

copy_to incorrectly maps scalar constraints

Closed this issue 4 years ago · 3 comments

using Clp
const MOI = Clp.MathOptInterface
model = MOI.Utilities.Model{Float64}()
x = MOI.add_variable(model)
con = [
    MOI.add_constraint(
        model, 
        MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(1.0, x)], 0.0),
        MOI.EqualTo(1.0),
    )
    for i = 1:2
]
clp = Clp.Optimizer()
index_map = MOI.copy_to(clp, model)
index_map[con[1]] === index_map[con[2]]  # True! But should be false.

Original issue text

I created a simple model in JuMP v0.21.5 (Julia Version 1.5.2):

max 10𝑋1+15𝑋2
Subject to
6𝑋1+20𝑋2≤2400.0
20𝑋1+20𝑋2≤4000.0
4𝑋1+9𝑋2≤900.0
𝑋1≤100.0
𝑋2≤100.0
𝑋1≥0.0
𝑋2≥0.0

After optimizing (Solver: Clp 0.8.1) I do get the correct objective value as well as the correct variable values(X1=100, X2=55.55555). But for the following it seems as there is a problem. In earlier versions of Julia this code used to work.

JuMP.has_duals(LP)
Out: true

JuMP.dual(CapacityRestriction[1]), JuMP.dual(CapacityRestriction[2]), JuMP.dual(CapacityRestriction[3])
Out: (0.0, 0.0, 0.0)

JuMP.shadow_price(CapacityRestriction[1]), JuMP.shadow_price(CapacityRestriction[2]), JuMP.shadow_price(CapacityRestriction[3])
Out: (-0.0, -0.0, -0.0)

The third shadow price is not 0. Why doesn`t it show the correct solution.
Thanks for your help.

Answer 1 · 2020-10-26T19:06:25.000Z

-0.0 is equivalent to 0.0:

julia> -0.0 == 0.0
true

However, you should never compare floating point numbers with strict equality. Always use isapprox.

I'm going to close because this isn't a bug in JuMP. If you have follow-up questions, please post on the community forum. Before doing so, please read PSA: make it easier to help you, and provide a reproducible example.

Answer 2 · 2020-10-26T19:08:49.000Z

~~Or is your point that the dual and shadow price should not be 0? Please provide reproducible code. I don't know which constraints are CapacityRestriction.~~

Maybe there is something weird going on:

julia> b = [2400, 4000, 900]
3-element Array{Int64,1}:
 2400
 4000
  900

julia> A = [6 20; 20 20; 4 9]
3×2 Array{Int64,2}:
  6  20
 20  20
  4   9

julia> using Clp, JuMP

julia> model = Model(Clp.Optimizer)
A JuMP Model
Feasibility problem with:
Variables: 0
Model mode: AUTOMATIC
CachingOptimizer state: EMPTY_OPTIMIZER
Solver name: Clp

julia> @variable(model, 0 <= x[1:2] <= 100)
2-element Array{VariableRef,1}:
 x[1]
 x[2]

julia> @objective(model, Max, 10 * x[1] + 15 * x[2])
10 x[1] + 15 x[2]

julia> @constraint(model, con, A * x .<= b)
3-element Array{ConstraintRef{Model,MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.LessThan{Float64}},ScalarShape},1}:
 6 x[1] + 20 x[2] ≤ 2400.0
 20 x[1] + 20 x[2] ≤ 4000.0
 4 x[1] + 9 x[2] ≤ 900.0

julia> optimize!(model)
Coin0506I Presolve 2 (-1) rows, 2 (0) columns and 4 (-2) elements
Clp0006I 0  Obj -0 Dual inf 37.499998 (2)
Clp0006I 1  Obj 1833.3333
Clp0000I Optimal - objective value 1833.3333
Coin0511I After Postsolve, objective 1833.3333, infeasibilities - dual 0 (0), primal 0 (0)
Clp0032I Optimal objective 1833.333333 - 1 iterations time 0.002, Presolve 0.00

julia> dual.(con)
3-element Array{Float64,1}:
 -1.6666666666666665
 -1.6666666666666665
 -1.6666666666666665

Compare with pure MOI:

julia> b = [2400.0, 4000.0, 900.0]
3-element Array{Float64,1}:
 2400.0
 4000.0
  900.0

julia> A = [6.0 20.0; 20.0 20.0; 4.0 9.0]
3×2 Array{Float64,2}:
  6.0  20.0
 20.0  20.0
  4.0   9.0

julia> model = MOI.Utilities.Model{Float64}()
MOIU.Model{Float64}

julia> x = MOI.add_variables(model, 2)
2-element Array{MathOptInterface.VariableIndex,1}:
 MathOptInterface.VariableIndex(1)
 MathOptInterface.VariableIndex(2)

julia> MOI.add_constraint.(model, MOI.SingleVariable.(x), MOI.Interval(0.0, 100.0))
2-element Array{MathOptInterface.ConstraintIndex{MathOptInterface.SingleVariable,MathOptInterface.Interval{Float64}},1}:
 MathOptInterface.ConstraintIndex{MathOptInterface.SingleVariable,MathOptInterface.Interval{Float64}}(1)
 MathOptInterface.ConstraintIndex{MathOptInterface.SingleVariable,MathOptInterface.Interval{Float64}}(2)

julia> MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE)
MAX_SENSE::OptimizationSense = 1

julia> MOI.set(
           model, 
           MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(),
           MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([10.0, 15.0], x), 0.0),
       )
MathOptInterface.ScalarAffineFunction{Float64}(MathOptInterface.ScalarAffineTerm{Float64}[MathOptInterface.ScalarAffineTerm{Float64}(10.0, MathOptInterface.VariableIndex(1)), MathOptInterface.ScalarAffineTerm{Float64}(15.0, MathOptInterface.VariableIndex(2))], 0.0)

julia> con = [
           MOI.add_constraint(
               model, 
               MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(A[i, :], x), 0.0),
               MOI.LessThan(bi),
           )
           for (i, bi) in enumerate(b)
       ]
3-element Array{MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.LessThan{Float64}},1}:
 MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.LessThan{Float64}}(1)
 MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.LessThan{Float64}}(2)
 MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.LessThan{Float64}}(3)

julia> clp = Clp.Optimizer()
Clp.Optimizer

julia> MOI.copy_to(clp, model)
MathOptInterface.Utilities.IndexMap with 7 entries:
  VariableIndex(1)                 => VariableIndex(1)
  VariableIndex(2)                 => VariableIndex(2)
  ConstraintIndex{ScalarAffineFun… => ConstraintIndex{ScalarAffineFunction{Float64},LessThan{Float64…
  ConstraintIndex{ScalarAffineFun… => ConstraintIndex{ScalarAffineFunction{Float64},LessThan{Float64…
  ConstraintIndex{ScalarAffineFun… => ConstraintIndex{ScalarAffineFunction{Float64},LessThan{Float64…
  ConstraintIndex{SingleVariable,… => ConstraintIndex{SingleVariable,Interval{Float64}}(2)
  ConstraintIndex{SingleVariable,… => ConstraintIndex{SingleVariable,Interval{Float64}}(1)

julia> MOI.optimize!(clp)
Coin0506I Presolve 2 (-1) rows, 2 (0) columns and 4 (-2) elements
Clp0006I 0  Obj -0 Dual inf 37.499998 (2)
Clp0006I 1  Obj 1833.3333
Clp0000I Optimal - objective value 1833.3333
Coin0511I After Postsolve, objective 1833.3333, infeasibilities - dual 0 (0), primal 0 (0)
Clp0032I Optimal objective 1833.333333 - 1 iterations time 0.002, Presolve 0.00

julia> MOI.get.(clp, MOI.ConstraintDual(), con)
3-element Array{Float64,1}:
  0.0
  0.0
 -1.6666666666666665

Answer 3 · 2020-10-27T11:08:11.000Z

I edited the package source code as you suggested and now everything works fine. Thanks for your help!