Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Wrapper for the Mixed Integer Linear Minimization Oracle #135

Merged
merged 86 commits into from
Oct 19, 2023
Merged
Show file tree
Hide file tree
Changes from 9 commits
Commits
Show all changes
86 commits
Select commit Hold shift + click to select a range
7a4bf49
Supertype for the Mixed Integer Linear Oracles.
Jul 25, 2023
19dde3e
BLMO wrapper and functions that have to be implemented.
Oct 13, 2023
b00426d
BLMO implemented for solvers supporting MathOptInterface.
Oct 13, 2023
5f09349
Rework build_LMO for the BLMO wrapper.
Oct 13, 2023
59bb4c9
Add the new files.
Oct 13, 2023
6f2dbd7
Started on adapting function for BLMO.
Oct 13, 2023
ccea3d9
Bounded Linear Minimzation Oracle as subtype of Linear Minimization O…
Oct 16, 2023
8191716
Adapted build_LMO for BLMO.
Oct 16, 2023
cbc8166
BLMO for all solvers supporting MOI.
Oct 16, 2023
fa434fc
Add BLMO wrapper and MOI BLMO to source.
Oct 16, 2023
ed5aa04
Adapting file ultilities.jl.
Oct 16, 2023
8174e2d
Boscia.compute_extreme_point is the same function as FrankWolfe.compu…
Oct 16, 2023
cea0ff9
Adapt time_tracking_lmo.jl.
Oct 16, 2023
c553a74
Adapt problem.jl.
Oct 16, 2023
c5df452
Adapt problem.jl II.
Oct 16, 2023
a619cd2
Adapt node.jl
Oct 16, 2023
7b757bb
Is not needed any longer.
Oct 16, 2023
ba01559
Convert between FW.MathOptLMO and MathOptBLMO.
Oct 16, 2023
e47d6d3
Convert MathOptBLMO into MathOptLMO.
Oct 16, 2023
a85593e
Get list of all variables, binaries and integers respectivily.
Oct 16, 2023
bae27ad
Deleted infeasible_pairwise.jl.
Oct 16, 2023
4c2482e
Adapt heuristics.jl.
Oct 16, 2023
dd90fbd
Adapt callback.jl
Oct 16, 2023
75d6bb5
Adapt strong branching.
Oct 16, 2023
a1b04d8
Adapt interface.jl.
Oct 16, 2023
f6295df
Rename files for more clearity.
Oct 16, 2023
f3918d8
Separate file for the build_LMO function.
Oct 16, 2023
a93213c
Fix syntax issues.
Oct 16, 2023
c5c03dc
Rather write Base.convert.
Oct 16, 2023
5441367
Have a solve function for MathOptLMO. Converts internally and calls o…
Oct 16, 2023
32e2f45
Minor change.
Oct 16, 2023
b2c0d4a
Better explanation of blmo
dhendryc Oct 16, 2023
0d8d0bc
Fix incremental compilation warning.
Oct 17, 2023
223794e
Signuature fix.
Oct 17, 2023
edf9b1c
minor change.
Oct 17, 2023
1733515
Redirect to FrankWolfe.compute_extreme_point for MathOptLMO.
Oct 17, 2023
4130c16
Syntax fix.
Oct 17, 2023
0c6413d
Merge changes
Oct 17, 2023
7090af6
Mix up between lower and upper bounds.
Oct 17, 2023
db6e3c3
Syntax fix.
Oct 17, 2023
c567780
Syntax fix.
Oct 17, 2023
f2a386e
Export Base.copy.
Oct 17, 2023
9d7c507
Syntax fix.
Oct 17, 2023
064d74c
Update src/MOI_bounded_oracle.jl
matbesancon Oct 17, 2023
9064c2c
No extra copy function needed.
Oct 17, 2023
d7fbc36
Merge changes
Oct 17, 2023
60d1b44
Changes to the Strong Branching struct.
Oct 17, 2023
701140a
Syntax fix.
Oct 17, 2023
2682cd8
Change to MathOptBLMO.
Oct 17, 2023
1b07904
Syntax change.
Oct 17, 2023
b56ae3c
Signature change.
Oct 17, 2023
edda014
Delete unnecessary code.
Oct 18, 2023
0deade1
Clean up.
Oct 18, 2023
5748d5e
Some rewriting.
Oct 18, 2023
f514281
Explicitly hand down all the arguments.
Oct 18, 2023
cda0fc3
Syntax fix.
Oct 18, 2023
23f760c
Update src/lmo_wrapper.jl
matbesancon Oct 18, 2023
3bf1ad4
Update src/lmo_wrapper.jl
matbesancon Oct 18, 2023
e310223
LMO over the cube.
Oct 18, 2023
297178b
Approx planted point with LMO over cube.
Oct 18, 2023
70674ee
Split up get_bound.
Oct 18, 2023
2f58f3c
Minor change.
Oct 18, 2023
098e7a7
Fomratting change.
Oct 18, 2023
f6fe5f9
Minor fix.
Oct 18, 2023
13d8e93
minor fix.
Oct 18, 2023
89fa456
Merge changes
Oct 18, 2023
e32de1b
Minor fix.
Oct 19, 2023
b79f45a
Move the cube LMO into the source files.
Oct 19, 2023
9057ee7
Cube BLMO now part of source, sees all dependencies.
Oct 19, 2023
005db86
Minor change.
Oct 19, 2023
512d5cb
Merge changes.
Oct 19, 2023
5e293ce
Account for numerical inaccuracies.
Oct 19, 2023
f6f1fc2
Renaming.2
Oct 19, 2023
7321fde
Renaming of fields.
Oct 19, 2023
b6d9fa5
Beware of numerical rounding issues.
Oct 19, 2023
32bd9f1
Both nodes can be pruned if one of them actually reaches the incumben…
Oct 19, 2023
b039daf
Update src/MOI_bounded_oracle.jl
matbesancon Oct 19, 2023
b91d306
Update src/cube_blmo.jl
matbesancon Oct 19, 2023
1a5ead0
Update src/MOI_bounded_oracle.jl
matbesancon Oct 19, 2023
231319e
Update src/cube_blmo.jl
matbesancon Oct 19, 2023
670bce8
cleanup
matbesancon Oct 19, 2023
854dc7b
readd file
matbesancon Oct 19, 2023
15c96b1
Merge changes.
Oct 19, 2023
6e0f597
Fix merge conflicts.
Oct 19, 2023
ca50fc3
Minor change.
Oct 19, 2023
1a9ef28
Hand down kwargs...
Oct 19, 2023
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion src/Boscia.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@ using Dates
const MOI = MathOptInterface
const MOIU = MOI.Utilities

import Base: convert
import MathOptSetDistances as MOD

include("time_tracking_lmo.jl")
Expand All @@ -23,7 +24,6 @@ include("node.jl")
include("custom_bonobo.jl")
include("callbacks.jl")
include("problem.jl")
include("infeasible_pairwise.jl")
include("heuristics.jl")
include("strong_branching.jl")
include("utilities.jl")
Expand Down
347 changes: 305 additions & 42 deletions src/MOI_bounded_oracle.jl
Original file line number Diff line number Diff line change
Expand Up @@ -17,13 +17,46 @@ function MathOptBLMO(lmo::FrankWolfe.MathOptLMO)
return MathOptBLMO(lmo.o, lmo.use_modfify)
end

"""
Convert object of Type MathOptLMO into MathOptBLMO and viceversa.
"""
function convert(::Type{MathOptBLMO}, lmo::FrankWolfe.MathOptLMO)
return MathOptBLMO(lmo.o, lmo.use_modify)
end
function convert(::Type{FrankWolfe.MathOptLMO}, nlmo::MathOptBLMO)
return FrankWolfe.MathOptLMO(blmo.o, blmo.use_modfify)
end


################## Necessary to implement ####################
"""
compute_extreme_point

Is implemented in the FrankWolfe package in file "moi_oracle.jl".
"""

"""
Get list of variables indices.
If the problem has n variables, they are expected to contiguous and ordered from 1 to n.
"""
function get_list_of_variables(blmo::MathOptBLMO)
v_indices = MOI.get(blmo.o, MOI.ListOfVariableIndices())
if v_indices != MOI.VariableIndex.(1:n)
error("Variables are expected to be contiguous and ordered from 1 to N")
end
return v_indices
end

"""
Get list of binary and integer variables, respectively.
"""
function get_binary_variables(blmo::MathOptBLMO)
return MOI.get(blmo.o, MOI.ListOfConstraintIndices{MOI.VariableIndex,MOI.ZeroOne}())
end
function get_integer_variables(blmo::MathOptBLMO)
return MOI.get(blmo.o, MOI.ListOfConstraintIndices{MOI.VariableIndex,MOI.Integer}())
end

"""
Get the index of the integer variable the bound is working on.
"""
Expand Down Expand Up @@ -120,6 +153,104 @@ function is_integer_constraint(blmo::MathOptBLMO, idx::Int)
return false, -1
end

"""
Is a given point v linear feasible for the model?
"""
function is_linear_feasible(blmo::MathOptBLMO, v::AbstractVector)
return is_linear_feasible(blmo.o, v)
end
function is_linear_feasible(o::MOI.ModelLike, v::AbstractVector)
valvar(f) = v[f.value]
for (F, S) in MOI.get(o, MOI.ListOfConstraintTypesPresent())
isfeasible = is_linear_feasible_subroutine(o, F, S, valvar)
if !isfeasible
return false
end
end
# satisfies all constraints
return true
end
# function barrier for performance
function is_linear_feasible_subroutine(o::MOI.ModelLike, ::Type{F}, ::Type{S}, valvar) where {F,S}
if S == MOI.ZeroOne || S <: MOI.Indicator || S == MOI.Integer
return true
end
cons_list = MOI.get(o, MOI.ListOfConstraintIndices{F,S}())
for c_idx in cons_list
func = MOI.get(o, MOI.ConstraintFunction(), c_idx)
val = MOIU.eval_variables(valvar, func)
set = MOI.get(o, MOI.ConstraintSet(), c_idx)
# @debug("Constraint: $(F)-$(S) $(func) = $(val) in $(set)")
dist = MOD.distance_to_set(MOD.DefaultDistance(), val, set)
scip_tol = 1e-6
if o isa SCIP.Optimizer
scip_tol = MOI.get(o, MOI.RawOptimizerAttribute("numerics/feastol"))
end
if dist > 5000.0 * scip_tol
@debug("Constraint: $(F)-$(S) $(func) = $(val) in $(set)")
@debug("Distance to set: $(dist)")
return false
end
end
return true
end

"""
Define a copy function for strong branching
"""
function Base.copy(blmo::MathOptBLMO)
return MathOptBLMO(blmo.o, blmo.use_modify)
end

"""
Read global bounds from the problem
"""
function build_global_bounds(blmo::MathOptBLMO)
global_bounds = Boscia.IntegerBounds()
for idx in integer_variables
for ST in (MOI.LessThan{Float64}, MOI.GreaterThan{Float64})
cidx = MOI.ConstraintIndex{MOI.VariableIndex,ST}(idx)
# Variable constraints to not have to be explicitly given, see Buchheim example
if MOI.is_valid(blmo.o, cidx)
s = MOI.get(blmo.o, MOI.ConstraintSet(), cidx)
push!(global_bounds, (idx, s))
end
end
cidx = MOI.ConstraintIndex{MOI.VariableIndex,MOI.Interval{Float64}}(idx)
if MOI.is_valid(blmo.o, cidx)
x = MOI.VariableIndex(idx)
s = MOI.get(blmo.o, MOI.ConstraintSet(), cidx)
MOI.delete(blmo.o, cidx)
MOI.add_constraint(blmo.o, x, MOI.GreaterThan(s.lower))
MOI.add_constraint(blmo.o, x, MOI.LessThan(s.upper))
push!(global_bounds, (idx, MOI.GreaterThan(s.lower)))
push!(global_bounds, (idx, MOI.LessThan(s.upper)))
end
@assert !MOI.is_valid(lmo.o, cidx)
end
return global_bounds
end

"""
Add explicit bounds for binary variables.
"""
function explicit_bounds_binary_var(blmo::MathOptBLMO, global_bounds::IntegerBounds, binary_variables)
# adding binary bounds explicitly
for idx in binary_variables
cidx = MOI.ConstraintIndex{MOI.VariableIndex,MOI.LessThan{Float64}}(idx)
if !MOI.is_valid(blmo.o, cidx)
MOI.add_constraint(blmo.o, MOI.VariableIndex(idx), MOI.LessThan(1.0))
end
@assert MOI.is_valid(blmo.o, cidx)
cidx = MOI.ConstraintIndex{MOI.VariableIndex,MOI.GreaterThan{Float64}}(idx)
if !MOI.is_valid(blmo.o, cidx)
MOI.add_constraint(blmo.o, MOI.VariableIndex(idx), MOI.GreaterThan(0.0))
end
global_bounds[idx, :greaterthan] = MOI.GreaterThan(0.0)
global_bounds[idx, :lessthan] = MOI.LessThan(1.0)
end
end


##################### Optional to implement ################

Expand Down Expand Up @@ -211,48 +342,6 @@ function get_BLMO_solve_data(blmo::MathOptBLMO)
return opt_times, numberofnodes, simplex_iterations
end

"""
Is a given point v linear feasible for the model?
"""
function is_linear_feasible(blmo::MathOptBLMO)
return is_linear_feasible(blmo.o)
end
function is_linear_feasible(o::MOI.ModelLike, v::AbstractVector)
valvar(f) = v[f.value]
for (F, S) in MOI.get(o, MOI.ListOfConstraintTypesPresent())
isfeasible = is_linear_feasible_subroutine(o, F, S, valvar)
if !isfeasible
return false
end
end
# satisfies all constraints
return true
end
# function barrier for performance
function is_linear_feasible_subroutine(o::MOI.ModelLike, ::Type{F}, ::Type{S}, valvar) where {F,S}
if S == MOI.ZeroOne || S <: MOI.Indicator || S == MOI.Integer
return true
end
cons_list = MOI.get(o, MOI.ListOfConstraintIndices{F,S}())
for c_idx in cons_list
func = MOI.get(o, MOI.ConstraintFunction(), c_idx)
val = MOIU.eval_variables(valvar, func)
set = MOI.get(o, MOI.ConstraintSet(), c_idx)
# @debug("Constraint: $(F)-$(S) $(func) = $(val) in $(set)")
dist = MOD.distance_to_set(MOD.DefaultDistance(), val, set)
scip_tol = 1e-6
if o isa SCIP.Optimizer
scip_tol = MOI.get(o, MOI.RawOptimizerAttribute("numerics/feastol"))
end
if dist > 5000.0 * scip_tol
@debug("Constraint: $(F)-$(S) $(func) = $(val) in $(set)")
@debug("Distance to set: $(dist)")
return false
end
end
return true
end

"""
Is a given point v indicator feasible, i.e. meets the indicator constraints? If applicable.
"""
Expand Down Expand Up @@ -305,3 +394,177 @@ end
function get_tol(o::MOI.AbstractOptimizer)
return 1e-06
end

"""
Find best solution from the solving process.
"""
function find_best_solution(f::Function, blmo::MathOptBLMO, vars, domain_oracle)
return find_best_solution(f, blmo.o, vars, domain_oracle)
end

"""
List of all variable pointers. Depends on how you save your variables internally.

Is used in `find_best_solution`.
"""
function get_variables_pointers(blmo::BoundedLinearMinimizationOracle, tree)
return [MOI.VariableIndex(var) for var in 1:tree.root.problem.nvars]
end

"""
Deal with infeasible vertex if necessary, e.g. check what caused it etc.
"""
function check_infeasible_vertex(blmo::MathOptBLMO, tree)
node = tree.nodes[tree.root.current_node_id[]]
node_bounds = node.local_bounds
for list in (node_bounds.lower_bounds, node_bounds.upper_bounds)
for (idx, set) in list
c_idx = MOI.ConstraintIndex{MOI.VariableIndex, typeof(set)}(idx)
@assert MOI.is_valid(state.tlmo.blmo.o, c_idx)
set2 = MOI.get(state.tlmo.blmo.o, MOI.ConstraintSet(), c_idx)
if !(set == set2)
MOI.set(tlmo.blmo.o, MOI.ConstraintSet(), c_idx, set)
set3 = MOI.get(tlmo.blmo.o, MOI.ConstraintSet(), c_idx)
@assert (set3 == set) "$((idx, set3, set))"
end
end
end
end

"""
Behavior for strong branching.
"""
function Bonobo.get_branching_variable(
tree::Bonobo.BnBTree,
branching::PartialStrongBranching{MathOptBLMO},
node::Bonobo.AbstractNode,
)
xrel = Bonobo.get_relaxed_values(tree, node)
max_lowerbound = -Inf
max_idx = -1
# copy problem and remove integer constraints
filtered_src = MOI.Utilities.ModelFilter(tree.root.problem.lmo.lmo.o) do item
if item isa Tuple
(_, S) = item
if S <: Union{MOI.Indicator,MOI.Integer,MOI.ZeroOne}
return false
end
end
return !(item isa MOI.ConstraintIndex{<:Any,<:Union{MOI.ZeroOne,MOI.Integer,MOI.Indicator}})
end
index_map = MOI.copy_to(branching.optimizer, filtered_src)
# sanity check, otherwise the functions need permuted indices
for (v1, v2) in index_map
if v1 isa MOI.VariableIndex
@assert v1 == v2
end
end
relaxed_lmo = FrankWolfe.MathOptLMO(branching.optimizer)
@assert !isempty(node.active_set)
active_set = copy(node.active_set)
empty!(active_set)
num_frac = 0
for idx in Bonobo.get_branching_indices(tree.root)
if !isapprox(xrel[idx], round(xrel[idx]), atol=tree.options.atol, rtol=tree.options.rtol)
# left node: x_i <= floor(̂x_i)
fxi = floor(xrel[idx])
# create LMO
boundsLeft = copy(node.local_bounds)
if haskey(boundsLeft.upper_bounds, idx)
delete!(boundsLeft.upper_bounds, idx)
end
push!(boundsLeft.upper_bounds, (idx => MOI.LessThan(fxi)))
build_LMO(
relaxed_lmo,
tree.root.problem.integer_variable_bounds,
boundsLeft,
Bonobo.get_branching_indices(tree.root),
)
MOI.optimize!(relaxed_lmo.o)
#MOI.set(relaxed_lmo.o, MOI.Silent(), false)
if MOI.get(relaxed_lmo.o, MOI.TerminationStatus()) == MOI.OPTIMAL
empty!(active_set)
for (λ, v) in node.active_set
if v[idx] <= xrel[idx]
push!(active_set, ((λ, v)))
end
end
@assert !isempty(active_set)
FrankWolfe.active_set_renormalize!(active_set)
_, _, primal_relaxed, dual_gap_relaxed, _ =
FrankWolfe.blended_pairwise_conditional_gradient(
tree.root.problem.f,
tree.root.problem.g,
relaxed_lmo,
active_set,
verbose=false,
epsilon=branching.solving_epsilon,
max_iteration=branching.max_iteration,
)
left_relaxed = primal_relaxed - dual_gap_relaxed
else
@debug "Left non-optimal status $(MOI.get(relaxed_lmo.o, MOI.TerminationStatus()))"
left_relaxed = Inf
end
#right node: x_i >= floor(̂x_i)
cxi = ceil(xrel[idx])
boundsRight = copy(node.local_bounds)
if haskey(boundsRight.lower_bounds, idx)
delete!(boundsRight.lower_bounds, idx)
end
push!(boundsRight.lower_bounds, (idx => MOI.GreaterThan(cxi)))
build_LMO(
relaxed_lmo,
tree.root.problem.integer_variable_bounds,
boundsRight,
Bonobo.get_branching_indices(tree.root),
)
MOI.optimize!(relaxed_lmo.o)
if MOI.get(relaxed_lmo.o, MOI.TerminationStatus()) == MOI.OPTIMAL
empty!(active_set)
for (λ, v) in node.active_set
if v[idx] >= xrel[idx]
push!(active_set, (λ, v))
end
end
if isempty(active_set)
@show xrel[idx]
@show length(active_set)
@info [active_set.atoms[idx] for idx in eachindex(active_set)]
error("Empty active set, unreachable")
end
FrankWolfe.active_set_renormalize!(active_set)
_, _, primal_relaxed, dual_gap_relaxed, _ =
FrankWolfe.blended_pairwise_conditional_gradient(
tree.root.problem.f,
tree.root.problem.g,
relaxed_lmo,
active_set,
verbose=false,
epsilon=branching.solving_epsilon,
max_iteration=branching.max_iteration,
)
right_relaxed = primal_relaxed - dual_gap_relaxed
else
@debug "Right non-optimal status $(MOI.get(relaxed_lmo.o, MOI.TerminationStatus()))"
right_relaxed = Inf
end
# lowest lower bound on the two branches
lowerbound_increase = min(left_relaxed, right_relaxed)
if lowerbound_increase > max_lowerbound
max_lowerbound = lowerbound_increase
max_idx = idx
end
num_frac += 1
end
end
@debug "strong branching: index $max_idx, lower bound $max_lowerbound"
if max_idx <= 0 && num_frac != 0
error("Infeasible node! Please check constraints! node lb: $(node.lb)")
max_idx = -1
end
if max_idx <= 0
max_idx = -1
end
return max_idx
end
Loading