jump-dev / MathOptInterface.jl

A data structure for mathematical optimization problems
http://jump.dev/MathOptInterface.jl/
Other
392 stars 87 forks source link

Revisit PR#2508 #2513

Closed odow closed 3 months ago

odow commented 3 months ago

See https://github.com/jump-dev/MathOptInterface.jl/actions/runs/9491359380/job/26156620575

Likely caused by https://github.com/jump-dev/MathOptInterface.jl/pull/2508

Found working on https://github.com/jump-dev/MathOptInterface.jl/pull/2495

odow commented 3 months ago

2508 is wrong, because they are not the same duals:

import CSDP
import MathOptInterface as MOI
function main()
    model = MOI.instantiate(CSDP.Optimizer, with_bridge_type = Float64)
    MOI.Bridges.remove_bridge(model, MOI.Bridges.Variable.ZerosBridge{Float64})
    config = MOI.Test.Config(; rtol = 1e-3, atol = 1e-3)
    MOI.Test.test_conic_PositiveSemidefiniteConeTriangle_VectorOfVariables(model, config)
    return model
end
model = main()
julia> model = main()
CSDP 6.2.0
Iter:  0 Ap: 0.00e+00 Pobj:  0.0000000e+00 Ad: 0.00e+00 Dobj:  0.0000000e+00 
Iter:  1 Ap: 1.00e+00 Pobj: -1.3091666e+01 Ad: 9.00e-01 Dobj: -1.4400002e+00 
Iter:  2 Ap: 1.00e+00 Pobj: -6.5001984e+00 Ad: 9.15e-01 Dobj: -8.2898589e-01 
Iter:  3 Ap: 9.50e-01 Pobj: -2.1617676e+00 Ad: 8.66e-01 Dobj: -1.4438101e+00 
Iter:  4 Ap: 1.00e+00 Pobj: -2.0165778e+00 Ad: 9.38e-01 Dobj: -1.9577837e+00 
Iter:  5 Ap: 1.00e+00 Pobj: -2.0005894e+00 Ad: 9.66e-01 Dobj: -1.9980351e+00 
Iter:  6 Ap: 1.00e+00 Pobj: -2.0000291e+00 Ad: 9.82e-01 Dobj: -1.9999373e+00 
Iter:  7 Ap: 1.00e+00 Pobj: -2.0000041e+00 Ad: 1.00e+00 Dobj: -1.9999974e+00 
Iter:  8 Ap: 1.00e+00 Pobj: -2.0000009e+00 Ad: 1.00e+00 Dobj: -2.0000001e+00 
Iter:  9 Ap: 1.00e+00 Pobj: -2.0000001e+00 Ad: 1.00e+00 Dobj: -2.0000000e+00 
Iter: 10 Ap: 9.70e-01 Pobj: -2.0000000e+00 Ad: 9.70e-01 Dobj: -2.0000000e+00 
Success: SDP solved
Primal objective value: -2.0000000e+00 
Dual objective value: -2.0000000e+00 
Relative primal infeasibility: 2.72e-16 
Relative dual infeasibility: 2.12e-09 
Real Relative Gap: -2.72e-10 
XZ Relative Gap: 6.89e-09 
DIMACS error measures: 2.72e-16 0.00e+00 4.34e-09 0.00e+00 -2.72e-10 6.89e-09
MOI.get(model, a, bridge.slack_in_set) = [1.0000000029485168, -1.0000000016370894, 1.0000000029485168]
MOI.get(model, a, bridge.equality) = [1.000000000703211, -2.0000000032741787, 1.000000000703211]
Test Failed at /Users/oscar/.julia/dev/MathOptInterface/src/Test/test_conic.jl:4857
  Expression: ≈(MOI.get(model, MOI.ConstraintDual(), cX), cXv, atol = atol, rtol = rtol)
   Evaluated: [1.000000000703211, -2.0000000032741787, 1.000000000703211] ≈ [1.0, -1.0, 1.0] (atol=0.001, rtol=0.001)
odow commented 3 months ago
julia> function main2()
           model = MOI.instantiate(CSDP.Optimizer, with_bridge_type = Float64)
           MOI.set(model, MOI.Silent(), true)
           x, _ = MOI.add_constrained_variables(model, MOI.Nonnegatives(2))
           g = MOI.Utilities.operate(vcat, Float64, 1.0 * x[1], 1.0, 1.0 * x[2])
           ci = MOI.add_constraint(model, g, MOI.PositiveSemidefiniteConeTriangle(2))
           MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE)
           f = 1.0 * x[1] + 1.0 * x[2]
           MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f)
           MOI.optimize!(model)
           dual = MOI.get(model, MOI.ConstraintDual(), ci)
           return model, dual
       end
main2 (generic function with 1 method)

julia> model, dual = main2();
MOI.get(model, a, bridge.slack_in_set) = [1.0000000200504315, -1.000000020050426, 1.0000000200504315]
MOI.get(model, a, bridge.equality) = [1.0000000100252102, -2.000000040100852, 1.0000000100252102]

julia> dual
3-element Vector{Float64}:
  1.0000000100252102
 -2.000000040100852
  1.0000000100252102

julia> model
MOIB.LazyBridgeOptimizer{MOIU.CachingOptimizer{CSDP.Optimizer, MOIU.UniversalFallback{MOIU.Model{Float64}}}}
├ Variable bridges: none
├ Constraint bridges:
│ ├ MOIB.Constraint.ScalarizeBridge{Float64, MOI.ScalarAffineFunction{Float64}, MOI.EqualTo{Float64}}
│ └ MOIB.Constraint.VectorSlackBridge{Float64, MOI.VectorAffineFunction{Float64}, MOI.PositiveSemidefiniteConeTriangle}
├ Objective bridges: none
└ model: MOIU.CachingOptimizer
  ├ state: ATTACHED_OPTIMIZER
  ├ mode: AUTOMATIC
  ├ model_cache: MOIU.UniversalFallback{MOIU.Model{Float64}}
  │ ├ ObjectiveSense: MIN_SENSE
  │ ├ ObjectiveFunctionType: MOI.ScalarAffineFunction{Float64}
  │ ├ NumberOfVariables: 5
  │ └ NumberOfConstraints: 5
  │   ├ MOI.ScalarAffineFunction{Float64} in MOI.EqualTo{Float64}: 3
  │   ├ MOI.VectorOfVariables in MOI.Nonnegatives: 1
  │   └ MOI.VectorOfVariables in MOI.PositiveSemidefiniteConeTriangle: 1
  └ optimizer: CSDP.Optimizer
    ├ ObjectiveSense: unknown
    ├ ObjectiveFunctionType: unknown
    ├ NumberOfVariables: unknown
    └ NumberOfConstraints: unknown
odow commented 3 months ago

I've reverted #2508 for now, but @blegat probably wants to take another look at a proper fix.

blegat commented 3 months ago

It makes sense, because the constraint is VectorAffineFunction-in-S and we transform to VectorAffineFunction-in-Zeros, if the scalar product between S and Zeros is not the same, we should scale the dual.

odow commented 3 months ago

How do we find the arbitrary scaling though?