From 39263a4eea9e6f5a41c49dcaa79e5910ee576614 Mon Sep 17 00:00:00 2001 From: rafabench_psr Date: Fri, 11 Dec 2020 14:20:07 -0300 Subject: [PATCH 1/2] Add farkas certificate for variable bounds --- src/MOI/MOI_wrapper.jl | 40 ++++++ src/api.jl | 2 +- test/MathOptInterface/MOI_Wrapper.jl | 191 +++++++++++++++++++++++++-- 3 files changed, 224 insertions(+), 9 deletions(-) diff --git a/src/MOI/MOI_wrapper.jl b/src/MOI/MOI_wrapper.jl index 1be315fe..0a125b5e 100644 --- a/src/MOI/MOI_wrapper.jl +++ b/src/MOI/MOI_wrapper.jl @@ -2567,6 +2567,29 @@ function _dual_multiplier(model::Optimizer) return MOI.get(model, MOI.ObjectiveSense()) == MOI.MIN_SENSE ? 1.0 : -1.0 end +""" + _farkas_variable_dual(model::Optimizer, col::Cint) + +Return a Farkas dual associated with the variable bounds of `col`. + +Compute the Farkas dual as: + + ā * x = λ' * A * x <= λ' * b = -β + sum(āᵢ * Uᵢ | āᵢ < 0) + sum(āᵢ * Lᵢ | āᵢ > 0) + +The Farkas dual of the variable is ā, and it applies to the upper bound if ā < 0, +and it applies to the lower bound if ā > 0. +""" +function _farkas_variable_dual(model::Optimizer, col::Int64) + nvars = length(model.variable_info) + nrows = length(model.affine_constraint_info) + mstart = Array{Cint}(undef,2) + mrwind = Array{Cint}(undef,nrows) + dmatval = Array{Float64}(undef,nrows) + ncoeffs = Array{Cint}(undef,1) + getcols(model.inner, mstart, mrwind, dmatval, nrows, ncoeffs, col, col) + return dmatval' * model.cached_solution.linear_dual +end + function MOI.get( model::Optimizer, attr::MOI.ConstraintDual, c::MOI.ConstraintIndex{MOI.SingleVariable, MOI.LessThan{Float64}} @@ -2574,6 +2597,10 @@ function MOI.get( _throw_if_optimize_in_progress(model, attr) MOI.check_result_index_bounds(model, attr) column = _info(model, c).column + if model.cached_solution.has_dual_certificate + dual = -_farkas_variable_dual(model, column) + return min(dual, 0.0) + end reduced_cost = model.cached_solution.variable_dual[column] sense = MOI.get(model, MOI.ObjectiveSense()) # The following is a heuristic for determining whether the reduced cost @@ -2601,6 +2628,10 @@ function MOI.get( _throw_if_optimize_in_progress(model, attr) MOI.check_result_index_bounds(model, attr) column = _info(model, c).column + if model.cached_solution.has_dual_certificate + dual = -_farkas_variable_dual(model, column) + return max(dual,0.0) + end reduced_cost = model.cached_solution.variable_dual[column] sense = MOI.get(model, MOI.ObjectiveSense()) # The following is a heuristic for determining whether the reduced cost @@ -2628,6 +2659,9 @@ function MOI.get( _throw_if_optimize_in_progress(model, attr) MOI.check_result_index_bounds(model, attr) column = _info(model, c).column + if model.cached_solution.has_dual_certificate + return -_farkas_variable_dual(model, column) + end reduced_cost = model.cached_solution.variable_dual[column] return _dual_multiplier(model) * reduced_cost end @@ -2639,6 +2673,9 @@ function MOI.get( _throw_if_optimize_in_progress(model, attr) MOI.check_result_index_bounds(model, attr) column = _info(model, c).column + if model.cached_solution.has_dual_certificate + return -_farkas_variable_dual(model, column) + end reduced_cost = model.cached_solution.variable_dual[column] return _dual_multiplier(model) * reduced_cost end @@ -2650,6 +2687,9 @@ function MOI.get( _throw_if_optimize_in_progress(model, attr) MOI.check_result_index_bounds(model, attr) row = _info(model, c).row + if model.cached_solution.has_dual_certificate + return model.cached_solution.linear_dual[row] + end return _dual_multiplier(model) * model.cached_solution.linear_dual[row] end diff --git a/src/api.jl b/src/api.jl index c185944d..3c83a74b 100644 --- a/src/api.jl +++ b/src/api.jl @@ -1546,7 +1546,7 @@ Returns the nonzeros in the constraint matrix for the columns in a given range """ function getcols(prob::XpressProblem, _mstart, _mrwind, _dmatval, maxcoeffs, ncoeffs, first::Integer, last::Integer) - @checked Lib.XPRSgetcols(prob, _mstart, _mrwind, _dmatval, maxcoeffs, ncoeffs, Cint(first), Cint(last)) + @checked Lib.XPRSgetcols(prob, _mstart, _mrwind, _dmatval, maxcoeffs, ncoeffs, Cint(first-1), Cint(last-1)) end # # Disable 64Bit versions do to reliability issues. diff --git a/test/MathOptInterface/MOI_Wrapper.jl b/test/MathOptInterface/MOI_Wrapper.jl index 0ccffa02..37673246 100644 --- a/test/MathOptInterface/MOI_Wrapper.jl +++ b/test/MathOptInterface/MOI_Wrapper.jl @@ -56,7 +56,7 @@ end "solve_qcp_edge_cases", "solve_qp_edge_cases", - # Issue #105: Farkas dual for variable bound not implemented. + # These tests require extra parameters to obtain certificates. "solve_farkas_equalto_upper", "solve_farkas_equalto_lower", "solve_farkas_lessthan", @@ -67,6 +67,14 @@ end "solve_farkas_variable_lessthan_max", ], ) + MOIT.solve_farkas_equalto_upper(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_equalto_lower(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_lessthan(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_greaterthan(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_interval_upper(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_interval_lower(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_variable_lessthan(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) + MOIT.solve_farkas_variable_lessthan_max(BRIDGED_CERTIFICATE_OPTIMIZER,CONFIG) MOIT.solve_qcp_edge_cases(BRIDGED_OPTIMIZER, CONFIG_LOW_TOL) MOIT.solve_qp_edge_cases(BRIDGED_OPTIMIZER, CONFIG_LOW_TOL) # MOIT.delete_soc_variables(OPTIMIZER, CONFIG_LOW_TOL) @@ -81,18 +89,13 @@ end infeas_certificates = false ), [ # These tests require extra parameters to obtain certificates. - "linear8a", "linear8b", "linear8c", - # TODO: This requires an infeasiblity certificate for a variable bound. - "linear12", + "linear8a", "linear8b", "linear8c","linear12", ] ) MOIT.linear8atest(BRIDGED_CERTIFICATE_OPTIMIZER, CONFIG) MOIT.linear8btest(BRIDGED_CERTIFICATE_OPTIMIZER, CONFIG) MOIT.linear8ctest(BRIDGED_CERTIFICATE_OPTIMIZER, CONFIG) - - MOIT.linear12test( - BRIDGED_OPTIMIZER, MOIT.TestConfig(infeas_certificates = false) - ) + MOIT.linear12test(BRIDGED_CERTIFICATE_OPTIMIZER, CONFIG) end @testset "Quadratic tests" begin @@ -280,3 +283,175 @@ end @test MOI.get(model, Xpress.ConstraintConflictStatus(), c2) == false end end + +@testset "test_farkas_dual_min" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.SingleVariable}(), + MOI.SingleVariable(x[1]), + ) + clb = MOI.add_constraint.( + model, MOI.SingleVariable.(x), MOI.GreaterThan(0.0) + ) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] > 1e-6 + @test clb_dual[2] > 1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ -2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ -c_dual atol = 1e-6 +end + +@testset "test_farkas_dual_min_interval" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.SingleVariable}(), + MOI.SingleVariable(x[1]), + ) + clb = MOI.add_constraint.( + model, MOI.SingleVariable.(x), MOI.Interval(0.0, 10.0) + ) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] > 1e-6 + @test clb_dual[2] > 1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ -2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ -c_dual atol = 1e-6 +end + +@testset "test_farkas_dual_min_equalto" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.SingleVariable}(), + MOI.SingleVariable(x[1]), + ) + clb = MOI.add_constraint.(model, MOI.SingleVariable.(x), MOI.EqualTo(0.0)) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] > 1e-6 + @test clb_dual[2] > 1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ -2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ -c_dual atol = 1e-6 +end + +@testset "test_farkas_dual_min_ii" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), + MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(-1.0, x[1])], 0.0), + ) + clb = MOI.add_constraint.( + model, MOI.SingleVariable.(x), MOI.LessThan(0.0) + ) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([-2.0, -1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] < -1e-6 + @test clb_dual[2] < -1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ 2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ c_dual atol = 1e-6 +end + +@testset "test_farkas_dual_max" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.SingleVariable}(), + MOI.SingleVariable(x[1]), + ) + clb = MOI.add_constraint.( + model, MOI.SingleVariable.(x), MOI.GreaterThan(0.0) + ) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([2.0, 1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] > 1e-6 + @test clb_dual[2] > 1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ -2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ -c_dual atol = 1e-6 +end + +@testset "test_farkas_dual_max_ii" begin + model = Xpress.Optimizer(OUTPUTLOG = 0, PRESOLVE = 0) + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) + MOI.set( + model, + MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(), + MOI.ScalarAffineFunction([MOI.ScalarAffineTerm(-1.0, x[1])], 0.0), + ) + clb = MOI.add_constraint.( + model, MOI.SingleVariable.(x), MOI.LessThan(0.0) + ) + c = MOI.add_constraint( + model, + MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.([-2.0, -1.0], x), 0.0), + MOI.LessThan(-1.0), + ) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.INFEASIBLE + @test MOI.get(model, MOI.DualStatus()) == MOI.INFEASIBILITY_CERTIFICATE + clb_dual = MOI.get.(model, MOI.ConstraintDual(), clb) + c_dual = MOI.get(model, MOI.ConstraintDual(), c) + @test clb_dual[1] < -1e-6 + @test clb_dual[2] < -1e-6 + @test c_dual[1] < -1e-6 + @test clb_dual[1] ≈ 2 * c_dual atol = 1e-6 + @test clb_dual[2] ≈ c_dual atol = 1e-6 +end \ No newline at end of file From 5ad939ec8492502d93887eae6bf3113a05c98cdb Mon Sep 17 00:00:00 2001 From: rafabench_psr Date: Fri, 11 Dec 2020 16:20:25 -0300 Subject: [PATCH 2/2] fix coefficients vector size --- src/MOI/MOI_wrapper.jl | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/src/MOI/MOI_wrapper.jl b/src/MOI/MOI_wrapper.jl index 0a125b5e..0e9cdfbe 100644 --- a/src/MOI/MOI_wrapper.jl +++ b/src/MOI/MOI_wrapper.jl @@ -2582,12 +2582,14 @@ and it applies to the lower bound if ā > 0. function _farkas_variable_dual(model::Optimizer, col::Int64) nvars = length(model.variable_info) nrows = length(model.affine_constraint_info) - mstart = Array{Cint}(undef,2) - mrwind = Array{Cint}(undef,nrows) - dmatval = Array{Float64}(undef,nrows) ncoeffs = Array{Cint}(undef,1) + getcols(model.inner, C_NULL, C_NULL, C_NULL, nrows, ncoeffs, col, col) + ncoeffs_ = ncoeffs[1] + mstart = Array{Cint}(undef,2) + mrwind = Array{Cint}(undef,ncoeffs_) + dmatval = Array{Float64}(undef,ncoeffs_) getcols(model.inner, mstart, mrwind, dmatval, nrows, ncoeffs, col, col) - return dmatval' * model.cached_solution.linear_dual + return sum(v * model.cached_solution.linear_dual[i + 1] for (i, v) in zip(mrwind, dmatval)) end function MOI.get(