From 6d91036a8956e26d1841cecd2d99cac4d3d7d90c Mon Sep 17 00:00:00 2001 From: Deborah Hendrych Date: Wed, 1 Nov 2023 21:25:49 +0100 Subject: [PATCH 01/15] Use cube as an example of how to implement the BLMO interface. --- examples/cube_blmo.jl | 162 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 162 insertions(+) create mode 100644 examples/cube_blmo.jl diff --git a/examples/cube_blmo.jl b/examples/cube_blmo.jl new file mode 100644 index 000000000..476390ab0 --- /dev/null +++ b/examples/cube_blmo.jl @@ -0,0 +1,162 @@ +## How to implement the BLMO Interface using the cube as an example + +""" + CubeBLMO + +A Bounded Linear Minimization Oracle over a cube. +""" +mutable struct CubeBLMO <: BoundedLinearMinimizationOracle + n::Int + int_vars::Vector{Int} + bounds::IntegerBounds + solving_time::Float64 +end + +CubeBLMO(n, int_vars, bounds) = CubeBLMO(n, int_vars, bounds, 0.0) + +## Necessary + +# computing an extreme point for the cube amounts to checking the sign of the gradient +function compute_extreme_point(blmo::CubeBLMO, d; kwargs...) + time_ref = Dates.now() + v = zeros(length(d)) + for i in eachindex(d) + v[i] = d[i] > 0 ? blmo.bounds[i, :greaterthan] : blmo.bounds[i, :lessthan] + end + blmo.solving_time = float(Dates.value(Dates.now() - time_ref)) + return v +end + +## + +function build_global_bounds(blmo::CubeBLMO, integer_variables) + global_bounds = IntegerBounds() + for i in 1:blmo.n + if i in integer_variables + push!(global_bounds, (i, blmo.bounds[i, :lessthan]), :lessthan) + push!(global_bounds, (i, blmo.bounds[i, :greaterthan]), :greaterthan) + end + end + return global_bounds +end + + +## Read information from problem +function get_list_of_variables(blmo::CubeBLMO) + return blmo.n, collect(1:blmo.n) +end + +# Get list of integer variables, respectively. +function get_integer_variables(blmo::CubeBLMO) + return blmo.int_vars +end + +function get_int_var(blmo::CubeBLMO, cidx) + return cidx +end + +function get_lower_bound_list(blmo::CubeBLMO) + return keys(blmo.bounds.lower_bounds) +end + +function get_upper_bound_list(blmo::CubeBLMO) + return keys(blmo.bounds.upper_bounds) +end + +function get_bound(blmo::CubeBLMO, c_idx, sense::Symbol) + @assert sense == :lessthan || sense == :greaterthan + return blmo[c_idx, sense] +end + +## Changing the bounds constraints. +function set_bound!(blmo::CubeBLMO, c_idx, value, sense::Symbol) + if sense == :greaterthan + blmo.bounds.lower_bounds[c_idx] = value + elseif sense == :lessthan + blmo.bounds.upper_bounds[c_idx] = value + else + error("Allowed values for sense are :lessthan and :greaterthan.") + end +end + +function delete_bounds!(blmo::CubeBLMO, cons_delete) + # For the cube this shouldn't happen! Otherwise we get unbounded! + if !isempty(cons_delete) + error("Trying to delete bounds of the cube!") + end +end + +function add_bound_constraint!(blmo::CubeBLMO, key, value, sense::Symbol) + # Should not be necessary + return error("Trying to add bound constraints of the cube!") +end + +## Checks + +function is_constraint_on_int_var(blmo::CubeBLMO, c_idx, int_vars) + return c_idx in int_vars +end + +function is_bound_in(blmo::CubeBLMO, c_idx, bounds) + return haskey(bounds, c_idx) +end + +function is_linear_feasible(blmo::CubeBLMO, v::AbstractVector) + for i in eachindex(v) + if !( + blmo.bounds[i, :greaterthan] ≤ v[i] + 1e-6 || + !(v[i] - 1e-6 ≤ blmo.bounds[i, :lessthan]) + ) + @debug( + "Vertex entry: $(v[i]) Lower bound: $(blmo.bounds[i, :greaterthan]) Upper bound: $(blmo.bounds[i, :lessthan]))" + ) + return false + end + end + return true +end + +function has_integer_constraint(blmo::CubeBLMO, idx) + return idx in blmo.int_vars +end + + +###################### Optional +## Safety Functions + +function build_LMO_correct(blmo::CubeBLMO, node_bounds) + for key in keys(node_bounds.lower_bounds) + if !haskey(blmo.bounds, (key, :greaterthan)) || + blmo.bounds[key, :greaterthan] != node_bounds[key, :greaterthan] + return false + end + end + for key in keys(node_bounds.upper_bounds) + if !haskey(blmo.bounds, (key, :lessthan)) || + blmo.bounds[key, :lessthan] != node_bounds[key, :lessthan] + return false + end + end + return true +end + +function check_feasibility(blmo::CubeBLMO) + for i in 1:blmo.n + if !haskey(blmo.bounds, (i, :greaterthan)) || !haskey(blmo.bounds, (i, :lessthan)) + return UNBOUNDED + end + if blmo.bounds[i, :greaterthan] > blmo.bounds[i, :lessthan] + return INFEASIBLE + end + end + return OPTIMAL +end + +function is_valid_split(tree::Bonobo.BnBTree, blmo::CubeBLMO, vidx::Int) + return blmo.bounds[vidx, :lessthan] != blmo.bounds[vidx, :greaterthan] +end + +## Logs +function get_BLMO_solve_data(blmo::CubeBLMO) + return blmo.solving_time, 0.0, 0.0 +end \ No newline at end of file From 5cd17445dfd971794bd7c4118d46f3b8e2b1a780 Mon Sep 17 00:00:00 2001 From: Deborah Hendrych Date: Wed, 1 Nov 2023 21:26:25 +0100 Subject: [PATCH 02/15] Implement some common polytopes as SimpleBoundableLMO. --- src/Boscia.jl | 2 +- src/cube_blmo.jl | 207 ------------------------------------------ src/polytope_blmos.jl | 35 +++++++ 3 files changed, 36 insertions(+), 208 deletions(-) delete mode 100644 src/cube_blmo.jl create mode 100644 src/polytope_blmos.jl diff --git a/src/Boscia.jl b/src/Boscia.jl index 3b5a14cf0..5e92fc6ca 100644 --- a/src/Boscia.jl +++ b/src/Boscia.jl @@ -27,7 +27,7 @@ include("utilities.jl") include("interface.jl") include("managed_blmo.jl") include("MOI_bounded_oracle.jl") -include("cube_blmo.jl") +include("polytope_blmos.jl") # For extensions if !isdefined(Base, :get_extension) diff --git a/src/cube_blmo.jl b/src/cube_blmo.jl deleted file mode 100644 index 9ac7f8fec..000000000 --- a/src/cube_blmo.jl +++ /dev/null @@ -1,207 +0,0 @@ - -""" - CubeBLMO - -A Bounded Linear Minimization Oracle over a cube. -""" -mutable struct CubeBLMO <: BoundedLinearMinimizationOracle - n::Int - int_vars::Vector{Int} - bounds::IntegerBounds - solving_time::Float64 -end - -CubeBLMO(n, int_vars, bounds) = CubeBLMO(n, int_vars, bounds, 0.0) - -## Necessary - -# computing an extreme point for the cube amounts to checking the sign of the gradient -function compute_extreme_point(blmo::CubeBLMO, d; kwargs...) - time_ref = Dates.now() - v = zeros(length(d)) - for i in eachindex(d) - v[i] = d[i] > 0 ? blmo.bounds[i, :greaterthan] : blmo.bounds[i, :lessthan] - end - blmo.solving_time = float(Dates.value(Dates.now() - time_ref)) - return v -end - -## - -function build_global_bounds(blmo::CubeBLMO, integer_variables) - global_bounds = IntegerBounds() - for i in 1:blmo.n - if i in integer_variables - push!(global_bounds, (i, blmo.bounds[i, :lessthan]), :lessthan) - push!(global_bounds, (i, blmo.bounds[i, :greaterthan]), :greaterthan) - end - end - return global_bounds -end - - -## Read information from problem -function get_list_of_variables(blmo::CubeBLMO) - return blmo.n, collect(1:blmo.n) -end - -# Get list of integer variables, respectively. -function get_integer_variables(blmo::CubeBLMO) - return blmo.int_vars -end - -function get_int_var(blmo::CubeBLMO, cidx) - return cidx -end - -function get_lower_bound_list(blmo::CubeBLMO) - return keys(blmo.bounds.lower_bounds) -end - -function get_upper_bound_list(blmo::CubeBLMO) - return keys(blmo.bounds.upper_bounds) -end - -function get_bound(blmo::CubeBLMO, c_idx, sense::Symbol) - @assert sense == :lessthan || sense == :greaterthan - return blmo[c_idx, sense] -end - -#function get_lower_bound(blmo::CubeBLMO, c_idx) -# return blmo.bounds[c_idx, :greaterthan] -#end - -#function get_upper_bound(blmo::CubeBLMO, c_idx) -# return blmo.bounds[c_idx, :lessthan] -#end - -## Changing the bounds constraints. -function set_bound!(blmo::CubeBLMO, c_idx, value, sense::Symbol) - if sense == :greaterthan - blmo.bounds.lower_bounds[c_idx] = value - elseif sense == :lessthan - blmo.bounds.upper_bounds[c_idx] = value - else - error("Allowed values for sense are :lessthan and :greaterthan.") - end -end - -function delete_bounds!(blmo::CubeBLMO, cons_delete) - # For the cube this shouldn't happen! Otherwise we get unbounded! - if !isempty(cons_delete) - error("Trying to delete bounds of the cube!") - end -end - -function add_bound_constraint!(blmo::CubeBLMO, key, value, sense::Symbol) - # Should not be necessary - return error("Trying to add bound constraints of the cube!") -end - -## Checks - -function is_constraint_on_int_var(blmo::CubeBLMO, c_idx, int_vars) - return c_idx in int_vars -end - -function is_bound_in(blmo::CubeBLMO, c_idx, bounds) - return haskey(bounds, c_idx) -end - -function is_linear_feasible(blmo::CubeBLMO, v::AbstractVector) - for i in eachindex(v) - if !( - blmo.bounds[i, :greaterthan] ≤ v[i] + 1e-6 || - !(v[i] - 1e-6 ≤ blmo.bounds[i, :lessthan]) - ) - @debug( - "Vertex entry: $(v[i]) Lower bound: $(blmo.bounds[i, :greaterthan]) Upper bound: $(blmo.bounds[i, :lessthan]))" - ) - return false - end - end - return true -end - -function has_integer_constraint(blmo::CubeBLMO, idx) - return idx in blmo.int_vars -end - - - -###################### Optional -## Safety Functions - -function build_LMO_correct(blmo::CubeBLMO, node_bounds) - for key in keys(node_bounds.lower_bounds) - if !haskey(blmo.bounds, (key, :greaterthan)) || - blmo.bounds[key, :greaterthan] != node_bounds[key, :greaterthan] - return false - end - end - for key in keys(node_bounds.upper_bounds) - if !haskey(blmo.bounds, (key, :lessthan)) || - blmo.bounds[key, :lessthan] != node_bounds[key, :lessthan] - return false - end - end - return true -end - -function check_feasibility(blmo::CubeBLMO) - for i in 1:blmo.n - if !haskey(blmo.bounds, (i, :greaterthan)) || !haskey(blmo.bounds, (i, :lessthan)) - return UNBOUNDED - end - if blmo.bounds[i, :greaterthan] > blmo.bounds[i, :lessthan] - return INFEASIBLE - end - end - return OPTIMAL -end - -function is_valid_split(tree::Bonobo.BnBTree, blmo::CubeBLMO, vidx::Int) - return blmo.bounds[vidx, :lessthan] != blmo.bounds[vidx, :greaterthan] -end - -## Logs -function get_BLMO_solve_data(blmo::CubeBLMO) - return blmo.solving_time, 0.0, 0.0 -end - -######################################################################## -""" - CubeSimpleBLMO{T}(lower_bounds, upper_bounds) - -Hypercube with lower and upper bounds implementing the `SimpleBoundableLMO` interface. -""" -struct CubeSimpleBLMO <: SimpleBoundableLMO - lower_bounds::Vector{Float64} - upper_bounds::Vector{Float64} - int_vars::Vector{Int} -end - -function bounded_compute_extreme_point(sblmo::CubeSimpleBLMO, d, lb, ub, int_vars; kwargs...) - v = zeros(length(d)) - for i in eachindex(d) - if i in int_vars - idx = findfirst(x -> x == i, int_vars) - v[i] = d[i] > 0 ? lb[idx] : ub[idx] - else - v[i] = d[i] > 0 ? sblmo.lower_bounds[i] : sblmo.upper_bounds[i] - end - end - return v -end - -function is_linear_feasible(sblmo::CubeSimpleBLMO, v) - for i in setdiff(eachindex(v), sblmo.int_vars) - if !(sblmo.lower_bounds[i] ≤ v[i] + 1e-6 || !(v[i] - 1e-6 ≤ blmo.upper_bounds[i])) - @debug( - "Vertex entry: $(v[i]) Lower bound: $(blmo.bounds[i, :greaterthan]) Upper bound: $(blmo.bounds[i, :lessthan]))" - ) - return false - end - end - return true -end diff --git a/src/polytope_blmos.jl b/src/polytope_blmos.jl new file mode 100644 index 000000000..2e4b36639 --- /dev/null +++ b/src/polytope_blmos.jl @@ -0,0 +1,35 @@ +""" + CubeSimpleBLMO{T}(lower_bounds, upper_bounds) + +Hypercube with lower and upper bounds implementing the `SimpleBoundableLMO` interface. +""" +struct CubeSimpleBLMO <: SimpleBoundableLMO + lower_bounds::Vector{Float64} + upper_bounds::Vector{Float64} + int_vars::Vector{Int} +end + +function bounded_compute_extreme_point(sblmo::CubeSimpleBLMO, d, lb, ub, int_vars; kwargs...) + v = zeros(length(d)) + for i in eachindex(d) + if i in int_vars + idx = findfirst(x -> x == i, int_vars) + v[i] = d[i] > 0 ? lb[idx] : ub[idx] + else + v[i] = d[i] > 0 ? sblmo.lower_bounds[i] : sblmo.upper_bounds[i] + end + end + return v +end + +function is_linear_feasible(sblmo::CubeSimpleBLMO, v) + for i in setdiff(eachindex(v), sblmo.int_vars) + if !(sblmo.lower_bounds[i] ≤ v[i] + 1e-6 || !(v[i] - 1e-6 ≤ blmo.upper_bounds[i])) + @debug( + "Vertex entry: $(v[i]) Lower bound: $(blmo.bounds[i, :greaterthan]) Upper bound: $(blmo.bounds[i, :lessthan]))" + ) + return false + end + end + return true +end From 815ea54fca73f3efe039805b755fdef9a6a07a44 Mon Sep 17 00:00:00 2001 From: Deborah Hendrych Date: Wed, 1 Nov 2023 22:05:06 +0100 Subject: [PATCH 03/15] Prob and Unit Simplex. --- src/polytope_blmos.jl | 86 +++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 86 insertions(+) diff --git a/src/polytope_blmos.jl b/src/polytope_blmos.jl index 2e4b36639..b3e53089b 100644 --- a/src/polytope_blmos.jl +++ b/src/polytope_blmos.jl @@ -9,6 +9,9 @@ struct CubeSimpleBLMO <: SimpleBoundableLMO int_vars::Vector{Int} end +""" +If the entry is positve, choose the lower bound. Else, choose the upper bound. +""" function bounded_compute_extreme_point(sblmo::CubeSimpleBLMO, d, lb, ub, int_vars; kwargs...) v = zeros(length(d)) for i in eachindex(d) @@ -33,3 +36,86 @@ function is_linear_feasible(sblmo::CubeSimpleBLMO, v) end return true end + +""" + ProbablitySimplexSimpleBLMO(N) + +Scaled Probability Simplex: ∑ x = N. +""" +struct ProbablitySimplexSimpleBLMOSimpleBLMO <: SimpleBoundableLMO + N::Float64 +end + +""" +Assign the largest possible values to the entries corresponding to the smallest entries of d. +""" +function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMOSimpleBLMO, d, lb, ub, int_vars; kwargs...) + v = zeros(length(d)) + indices = collect(1:length(d)) + perm = sortperm(d) + + for i in indices[perm] + if i in int_vars + idx = findfirst(x -> x == i, int_vars) + v[i] = min(ub[idx], sblmo.N - sum(v)) + else + v[i] = N - sum(v) + end + end + return v +end + +function is_linear_feasible(sblmo::ProbablitySimplexSimpleBLMOSimpleBLMO, v) + if sum(v .≥ 0) < length(v) + @debug "v has negative entries: $(v)" + return false + end + return isapprox(sum(v), sblmo.N, atol=1e-4, rtol=1e-2) +end + +""" + UnitSimplexSimpleBLMO(N) + +Scaled Unit Simplex: ∑ x ≤ N. +""" +struct UnitSimplexSimpleBLMOSimpleBLMO <: SimpleBoundableLMO + N::Float64 +end + +""" +For all positive entries of d, assign the corresponding lower bound. +For non-positive entries, assign largest possible value in increasing order. +""" +function bounded_compute_extreme_point(sblmo::UnitSimplexSimpleBLMOSimpleBLMO, d, lb, ub, int_vars; kwargs...) + v = zeros(length(d)) + # d[i] positive + idx_pos = findall(x-> x > 0, d) + for i in idx_pos + if i in int_vars + idx = findfirst(x -> x == i, int_vars) + v[i] = lb[idx] + else + v[i] = 0.0 + end + end + # d[i] negative + idx_neg = setdiff(eachindex(d), idx_pos) + perm = sortperm(d[idx_neg]) + for i in idx_neg[perm] + if i in int_vars + idx = findfirst(x -> x == i, int_vars) + v[i] = min(ub[idx], sblmo.N - sum(v)) + else + v[i] = N - sum(v) + end + end + return v +end + +function is_linear_feasible(sblmo::UnitSimplexSimpleBLMOSimpleBLMO, v) + if sum(v .≥ 0) < length(v) + @debug "v has negative entries: $(v)" + return false + end + return sum(v) ≤ sblmo.N + 1e-3 +end From 19559cccfed0b869d62774766c0e854f28ab9adc Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 14:12:38 +0100 Subject: [PATCH 04/15] Add source file. --- examples/approx_planted_point.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/examples/approx_planted_point.jl b/examples/approx_planted_point.jl index d027eab42..d0bb4a8b9 100644 --- a/examples/approx_planted_point.jl +++ b/examples/approx_planted_point.jl @@ -9,6 +9,8 @@ using Distributions import MathOptInterface const MOI = MathOptInterface +include("cube_blmo.jl") + n = 20 diffi = Random.rand(Bool, n) * 0.6 .+ 0.3 From a920f1e82ae4fcd1dc6ec3ae05a793a31bde6917 Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 14:26:12 +0100 Subject: [PATCH 05/15] Minor change. --- src/polytope_blmos.jl | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/polytope_blmos.jl b/src/polytope_blmos.jl index b3e53089b..386cd47b8 100644 --- a/src/polytope_blmos.jl +++ b/src/polytope_blmos.jl @@ -42,14 +42,14 @@ end Scaled Probability Simplex: ∑ x = N. """ -struct ProbablitySimplexSimpleBLMOSimpleBLMO <: SimpleBoundableLMO +struct ProbablitySimplexSimpleBLMO <: SimpleBoundableLMO N::Float64 end """ Assign the largest possible values to the entries corresponding to the smallest entries of d. """ -function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMOSimpleBLMO, d, lb, ub, int_vars; kwargs...) +function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMO, d, lb, ub, int_vars; kwargs...) v = zeros(length(d)) indices = collect(1:length(d)) perm = sortperm(d) @@ -65,7 +65,7 @@ function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMOSimpleB return v end -function is_linear_feasible(sblmo::ProbablitySimplexSimpleBLMOSimpleBLMO, v) +function is_linear_feasible(sblmo::ProbablitySimplexSimpleBLMO, v) if sum(v .≥ 0) < length(v) @debug "v has negative entries: $(v)" return false @@ -78,7 +78,7 @@ end Scaled Unit Simplex: ∑ x ≤ N. """ -struct UnitSimplexSimpleBLMOSimpleBLMO <: SimpleBoundableLMO +struct UnitSimplexSimpleBLMO <: SimpleBoundableLMO N::Float64 end @@ -86,7 +86,7 @@ end For all positive entries of d, assign the corresponding lower bound. For non-positive entries, assign largest possible value in increasing order. """ -function bounded_compute_extreme_point(sblmo::UnitSimplexSimpleBLMOSimpleBLMO, d, lb, ub, int_vars; kwargs...) +function bounded_compute_extreme_point(sblmo::UnitSimplexSimpleBLMO, d, lb, ub, int_vars; kwargs...) v = zeros(length(d)) # d[i] positive idx_pos = findall(x-> x > 0, d) @@ -112,7 +112,7 @@ function bounded_compute_extreme_point(sblmo::UnitSimplexSimpleBLMOSimpleBLMO, d return v end -function is_linear_feasible(sblmo::UnitSimplexSimpleBLMOSimpleBLMO, v) +function is_linear_feasible(sblmo::UnitSimplexSimpleBLMO, v) if sum(v .≥ 0) < length(v) @debug "v has negative entries: $(v)" return false From 3bd1304f222f3a11379d5756962af0709ea17776 Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 16:16:02 +0100 Subject: [PATCH 06/15] Fix typo. --- src/polytope_blmos.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/polytope_blmos.jl b/src/polytope_blmos.jl index 386cd47b8..b0b15485a 100644 --- a/src/polytope_blmos.jl +++ b/src/polytope_blmos.jl @@ -42,14 +42,14 @@ end Scaled Probability Simplex: ∑ x = N. """ -struct ProbablitySimplexSimpleBLMO <: SimpleBoundableLMO +struct ProbabilitySimplexSimpleBLMO <: SimpleBoundableLMO N::Float64 end """ Assign the largest possible values to the entries corresponding to the smallest entries of d. """ -function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMO, d, lb, ub, int_vars; kwargs...) +function bounded_compute_extreme_point(sblmo::ProbabilitySimplexSimpleBLMO, d, lb, ub, int_vars; kwargs...) v = zeros(length(d)) indices = collect(1:length(d)) perm = sortperm(d) @@ -65,7 +65,7 @@ function bounded_compute_extreme_point(sblmo::ProbablitySimplexSimpleBLMO, d, lb return v end -function is_linear_feasible(sblmo::ProbablitySimplexSimpleBLMO, v) +function is_linear_feasible(sblmo::ProbabilitySimplexSimpleBLMO, v) if sum(v .≥ 0) < length(v) @debug "v has negative entries: $(v)" return false From 840c6abefe36561dbbaaab0a1e28273884733c51 Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 17:07:39 +0100 Subject: [PATCH 07/15] More information in debug statement. --- src/managed_blmo.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/managed_blmo.jl b/src/managed_blmo.jl index 6c279e490..37c9f4f58 100644 --- a/src/managed_blmo.jl +++ b/src/managed_blmo.jl @@ -174,7 +174,7 @@ function is_linear_feasible(blmo::ManagedBoundedLMO, v::AbstractVector) blmo.lower_bounds[i] ≤ v[int_var] + 1e-6 || !(v[int_var] - 1e-6 ≤ blmo.upper_bounds[i]) ) @debug( - "Vertex entry: $(v[int_var]) Lower bound: $(blmo.lower_bounds[i]) Upper bound: $(blmo.upper_bounds[i]))" + "Variable: $(int_var) Vertex entry: $(v[int_var]) Lower bound: $(blmo.lower_bounds[i]) Upper bound: $(blmo.upper_bounds[i]))" ) return false end From 60a9b7599911885028d0b74e5698040da45a4b12 Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 17:16:48 +0100 Subject: [PATCH 08/15] Tests on Unit and Probability Simplex. --- src/polytope_blmos.jl | 29 +++++++++--------- test/LMO_test.jl | 70 +++++++++++++++++++++++++++++++++++++++++-- 2 files changed, 82 insertions(+), 17 deletions(-) diff --git a/src/polytope_blmos.jl b/src/polytope_blmos.jl index b0b15485a..4bda95e36 100644 --- a/src/polytope_blmos.jl +++ b/src/polytope_blmos.jl @@ -54,12 +54,15 @@ function bounded_compute_extreme_point(sblmo::ProbabilitySimplexSimpleBLMO, d, l indices = collect(1:length(d)) perm = sortperm(d) + # The lower bounds always have to be met. + v[int_vars] = lb + for i in indices[perm] if i in int_vars idx = findfirst(x -> x == i, int_vars) - v[i] = min(ub[idx], sblmo.N - sum(v)) + v[i] += min(ub[idx]-lb[idx], sblmo.N - sum(v)) else - v[i] = N - sum(v) + v[i] += N - sum(v) end end return v @@ -88,25 +91,21 @@ For non-positive entries, assign largest possible value in increasing order. """ function bounded_compute_extreme_point(sblmo::UnitSimplexSimpleBLMO, d, lb, ub, int_vars; kwargs...) v = zeros(length(d)) - # d[i] positive - idx_pos = findall(x-> x > 0, d) - for i in idx_pos - if i in int_vars - idx = findfirst(x -> x == i, int_vars) - v[i] = lb[idx] - else - v[i] = 0.0 - end + # The wloer bounds always have to be met. + v[int_vars] = lb + cont_vars = setdiff(collect(1:length(d)), int_vars) + if !isempty(cont_vars) + v[cont_vars] .= 0.0 end - # d[i] negative - idx_neg = setdiff(eachindex(d), idx_pos) + + idx_neg = findall(x-> x <= 0, d) perm = sortperm(d[idx_neg]) for i in idx_neg[perm] if i in int_vars idx = findfirst(x -> x == i, int_vars) - v[i] = min(ub[idx], sblmo.N - sum(v)) + v[i] += min(ub[idx]-lb[idx], sblmo.N - sum(v)) else - v[i] = N - sum(v) + v[i] += N - sum(v) end end return v diff --git a/test/LMO_test.jl b/test/LMO_test.jl index 35005e0fd..65ad2a982 100644 --- a/test/LMO_test.jl +++ b/test/LMO_test.jl @@ -11,8 +11,6 @@ using Dates const MOI = MathOptInterface const MOIU = MOI.Utilities -import MathOptSetDistances -const MOD = MathOptSetDistances @testset "Integer bounds" begin n = 10 integer_bound = Boscia.IntegerBounds() @@ -66,6 +64,27 @@ diffi = Random.rand(Bool, n) * 0.6 .+ 0.3 function grad!(storage, x) @. storage = x - diffi end + + lbs = zeros(n) + ubs = ones(n) + int_vars = collect(1:n) + + sblmo = Boscia.CubeSimpleBLMO(lbs, ubs, int_vars) + + x, _, result = + Boscia.solve(f, grad!, sblmo, lbs[int_vars], ubs[int_vars], int_vars, n) + + @test sum(isapprox.(x, round.(diffi), atol=1e-6, rtol=1e-2)) == n + @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) +end + +@testset "BLMO - Strong Branching" begin + function f(x) + return 0.5 * sum((x[i] - diffi[i])^2 for i in eachindex(x)) + end + function grad!(storage, x) + @. storage = x - diffi + end @testset "Partial Strong Branching" begin int_vars = collect(1:n) lbs = zeros(n) @@ -102,3 +121,50 @@ diffi = Random.rand(Bool, n) * 0.6 .+ 0.3 @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) end end + +n = 20 +x_sol = rand(1:Int(floor(n/4)), n) +N = sum(x_sol) +dir = vcat(fill(1, Int(floor(n/2))), fill(-1, Int(floor(n/2))), fill(0, mod(n,2))) +@assert length(dir) == n +diffi = x_sol + 0.3 * dir +@assert sum(diffi) == N + +@testset "Probability Simplex LMO" begin + function f(x) + return 0.5 * sum((x[i] - diffi[i])^2 for i in eachindex(x)) + end + function grad!(storage, x) + @. storage = x - diffi + end + + sblmo = Boscia.ProbabilitySimplexSimpleBLMO(N) + + x, _, result = + Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(1.0*N, n), collect(1:n), n, verbose=true) + + @test sum(isapprox.(x, x_sol, atol=1e-6, rtol=1e-2)) == n + @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) +end + +n = 20 +x_sol = rand(1:Int(floor(n/4)), n) +diffi = x_sol + 0.3*rand([-1,1], n) + +@testset "Unit Simplex LMO" begin + function f(x) + return 0.5 * sum((x[i] - diffi[i])^2 for i in eachindex(x)) + end + function grad!(storage, x) + @. storage = x - diffi + end + + N = sum(x_sol) + floor(n/2) + sblmo = Boscia.UnitSimplexSimpleBLMO(N) + + x, _, result = + Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(N, n), collect(1:n), n, verbose=true) + + @test sum(isapprox.(x, x_sol, atol=1e-6, rtol=1e-2)) == n + @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) +end From 55bc04257ff1e30cd5eba9f6ad26bd7647ff35fa Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 17:20:36 +0100 Subject: [PATCH 09/15] No logs for tests. --- test/LMO_test.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/LMO_test.jl b/test/LMO_test.jl index 65ad2a982..a337b46a1 100644 --- a/test/LMO_test.jl +++ b/test/LMO_test.jl @@ -141,7 +141,7 @@ diffi = x_sol + 0.3 * dir sblmo = Boscia.ProbabilitySimplexSimpleBLMO(N) x, _, result = - Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(1.0*N, n), collect(1:n), n, verbose=true) + Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(1.0*N, n), collect(1:n), n) @test sum(isapprox.(x, x_sol, atol=1e-6, rtol=1e-2)) == n @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) @@ -163,7 +163,7 @@ diffi = x_sol + 0.3*rand([-1,1], n) sblmo = Boscia.UnitSimplexSimpleBLMO(N) x, _, result = - Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(N, n), collect(1:n), n, verbose=true) + Boscia.solve(f, grad!, sblmo, fill(0.0, n), fill(N, n), collect(1:n), n) @test sum(isapprox.(x, x_sol, atol=1e-6, rtol=1e-2)) == n @test isapprox(f(x), f(result[:raw_solution]), atol=1e-6, rtol=1e-3) From c59ee0ab5df6f92073d11aeb893a88a6c4bab6db Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 17:35:22 +0100 Subject: [PATCH 10/15] Minor change. --- test/LMO_test.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/test/LMO_test.jl b/test/LMO_test.jl index a337b46a1..c879bd715 100644 --- a/test/LMO_test.jl +++ b/test/LMO_test.jl @@ -126,9 +126,7 @@ n = 20 x_sol = rand(1:Int(floor(n/4)), n) N = sum(x_sol) dir = vcat(fill(1, Int(floor(n/2))), fill(-1, Int(floor(n/2))), fill(0, mod(n,2))) -@assert length(dir) == n diffi = x_sol + 0.3 * dir -@assert sum(diffi) == N @testset "Probability Simplex LMO" begin function f(x) From 5ba236905e89c4cdc6afbdff6e350bad055c8c02 Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 17:41:23 +0100 Subject: [PATCH 11/15] Add Boscia identifier. --- examples/cube_blmo.jl | 40 ++++++++++++++++++++-------------------- 1 file changed, 20 insertions(+), 20 deletions(-) diff --git a/examples/cube_blmo.jl b/examples/cube_blmo.jl index 476390ab0..17f599bea 100644 --- a/examples/cube_blmo.jl +++ b/examples/cube_blmo.jl @@ -5,7 +5,7 @@ A Bounded Linear Minimization Oracle over a cube. """ -mutable struct CubeBLMO <: BoundedLinearMinimizationOracle +mutable struct CubeBLMO <: Boscia.BoundedLinearMinimizationOracle n::Int int_vars::Vector{Int} bounds::IntegerBounds @@ -17,7 +17,7 @@ CubeBLMO(n, int_vars, bounds) = CubeBLMO(n, int_vars, bounds, 0.0) ## Necessary # computing an extreme point for the cube amounts to checking the sign of the gradient -function compute_extreme_point(blmo::CubeBLMO, d; kwargs...) +function Boscia.compute_extreme_point(blmo::CubeBLMO, d; kwargs...) time_ref = Dates.now() v = zeros(length(d)) for i in eachindex(d) @@ -29,7 +29,7 @@ end ## -function build_global_bounds(blmo::CubeBLMO, integer_variables) +function Boscia.build_global_bounds(blmo::CubeBLMO, integer_variables) global_bounds = IntegerBounds() for i in 1:blmo.n if i in integer_variables @@ -42,34 +42,34 @@ end ## Read information from problem -function get_list_of_variables(blmo::CubeBLMO) +function Boscia.get_list_of_variables(blmo::CubeBLMO) return blmo.n, collect(1:blmo.n) end # Get list of integer variables, respectively. -function get_integer_variables(blmo::CubeBLMO) +function Boscia.get_integer_variables(blmo::CubeBLMO) return blmo.int_vars end -function get_int_var(blmo::CubeBLMO, cidx) +function Boscia.get_int_var(blmo::CubeBLMO, cidx) return cidx end -function get_lower_bound_list(blmo::CubeBLMO) +function Boscia.get_lower_bound_list(blmo::CubeBLMO) return keys(blmo.bounds.lower_bounds) end -function get_upper_bound_list(blmo::CubeBLMO) +function Boscia.get_upper_bound_list(blmo::CubeBLMO) return keys(blmo.bounds.upper_bounds) end -function get_bound(blmo::CubeBLMO, c_idx, sense::Symbol) +function Boscia.get_bound(blmo::CubeBLMO, c_idx, sense::Symbol) @assert sense == :lessthan || sense == :greaterthan return blmo[c_idx, sense] end ## Changing the bounds constraints. -function set_bound!(blmo::CubeBLMO, c_idx, value, sense::Symbol) +function Boscia.set_bound!(blmo::CubeBLMO, c_idx, value, sense::Symbol) if sense == :greaterthan blmo.bounds.lower_bounds[c_idx] = value elseif sense == :lessthan @@ -79,29 +79,29 @@ function set_bound!(blmo::CubeBLMO, c_idx, value, sense::Symbol) end end -function delete_bounds!(blmo::CubeBLMO, cons_delete) +function Boscia.delete_bounds!(blmo::CubeBLMO, cons_delete) # For the cube this shouldn't happen! Otherwise we get unbounded! if !isempty(cons_delete) error("Trying to delete bounds of the cube!") end end -function add_bound_constraint!(blmo::CubeBLMO, key, value, sense::Symbol) +function Boscia.add_bound_constraint!(blmo::CubeBLMO, key, value, sense::Symbol) # Should not be necessary return error("Trying to add bound constraints of the cube!") end ## Checks -function is_constraint_on_int_var(blmo::CubeBLMO, c_idx, int_vars) +function Boscia.is_constraint_on_int_var(blmo::CubeBLMO, c_idx, int_vars) return c_idx in int_vars end -function is_bound_in(blmo::CubeBLMO, c_idx, bounds) +function Boscia.is_bound_in(blmo::CubeBLMO, c_idx, bounds) return haskey(bounds, c_idx) end -function is_linear_feasible(blmo::CubeBLMO, v::AbstractVector) +function Boscia.is_linear_feasible(blmo::CubeBLMO, v::AbstractVector) for i in eachindex(v) if !( blmo.bounds[i, :greaterthan] ≤ v[i] + 1e-6 || @@ -116,7 +116,7 @@ function is_linear_feasible(blmo::CubeBLMO, v::AbstractVector) return true end -function has_integer_constraint(blmo::CubeBLMO, idx) +function Boscia.has_integer_constraint(blmo::CubeBLMO, idx) return idx in blmo.int_vars end @@ -124,7 +124,7 @@ end ###################### Optional ## Safety Functions -function build_LMO_correct(blmo::CubeBLMO, node_bounds) +function Boscia.build_LMO_correct(blmo::CubeBLMO, node_bounds) for key in keys(node_bounds.lower_bounds) if !haskey(blmo.bounds, (key, :greaterthan)) || blmo.bounds[key, :greaterthan] != node_bounds[key, :greaterthan] @@ -140,7 +140,7 @@ function build_LMO_correct(blmo::CubeBLMO, node_bounds) return true end -function check_feasibility(blmo::CubeBLMO) +function Boscia.check_feasibility(blmo::CubeBLMO) for i in 1:blmo.n if !haskey(blmo.bounds, (i, :greaterthan)) || !haskey(blmo.bounds, (i, :lessthan)) return UNBOUNDED @@ -152,11 +152,11 @@ function check_feasibility(blmo::CubeBLMO) return OPTIMAL end -function is_valid_split(tree::Bonobo.BnBTree, blmo::CubeBLMO, vidx::Int) +function Boscia.is_valid_split(tree::Bonobo.BnBTree, blmo::CubeBLMO, vidx::Int) return blmo.bounds[vidx, :lessthan] != blmo.bounds[vidx, :greaterthan] end ## Logs -function get_BLMO_solve_data(blmo::CubeBLMO) +function Boscia.get_BLMO_solve_data(blmo::CubeBLMO) return blmo.solving_time, 0.0, 0.0 end \ No newline at end of file From cec1b41b36e336dbf6c4131adcd04f2e7265539e Mon Sep 17 00:00:00 2001 From: Hendrych Date: Thu, 2 Nov 2023 18:10:09 +0100 Subject: [PATCH 12/15] Bug Fix. --- examples/approx_planted_point.jl | 5 +++-- examples/cube_blmo.jl | 13 ++++++++----- 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/examples/approx_planted_point.jl b/examples/approx_planted_point.jl index d0bb4a8b9..32da2bcfd 100644 --- a/examples/approx_planted_point.jl +++ b/examples/approx_planted_point.jl @@ -1,5 +1,6 @@ using Boscia using FrankWolfe +using Bonobo using Test using Random using SCIP @@ -49,7 +50,7 @@ diffi = Random.rand(Bool, n) * 0.6 .+ 0.3 push!(bounds, (i, 0.0), :greaterthan) push!(bounds, (i, 1.0), :lessthan) end - blmo = Boscia.CubeBLMO(n, int_vars, bounds) + blmo = CubeBLMO(n, int_vars, bounds) x, _, result = Boscia.solve(f, grad!, blmo, verbose=true) @@ -112,7 +113,7 @@ end push!(bounds, (i, 0.0), :greaterthan) push!(bounds, (i, 1.0), :lessthan) end - blmo = Boscia.CubeBLMO(n, int_vars, bounds) + blmo = CubeBLMO(n, int_vars, bounds) x, _, result = Boscia.solve(f, grad!, blmo, verbose=true) diff --git a/examples/cube_blmo.jl b/examples/cube_blmo.jl index 17f599bea..76fa6230d 100644 --- a/examples/cube_blmo.jl +++ b/examples/cube_blmo.jl @@ -1,4 +1,7 @@ ## How to implement the BLMO Interface using the cube as an example +using Boscia +using Bonobo +using Dates """ CubeBLMO @@ -8,7 +11,7 @@ A Bounded Linear Minimization Oracle over a cube. mutable struct CubeBLMO <: Boscia.BoundedLinearMinimizationOracle n::Int int_vars::Vector{Int} - bounds::IntegerBounds + bounds::Boscia.IntegerBounds solving_time::Float64 end @@ -30,7 +33,7 @@ end ## function Boscia.build_global_bounds(blmo::CubeBLMO, integer_variables) - global_bounds = IntegerBounds() + global_bounds = Boscia.IntegerBounds() for i in 1:blmo.n if i in integer_variables push!(global_bounds, (i, blmo.bounds[i, :lessthan]), :lessthan) @@ -143,13 +146,13 @@ end function Boscia.check_feasibility(blmo::CubeBLMO) for i in 1:blmo.n if !haskey(blmo.bounds, (i, :greaterthan)) || !haskey(blmo.bounds, (i, :lessthan)) - return UNBOUNDED + return Boscia.UNBOUNDED end if blmo.bounds[i, :greaterthan] > blmo.bounds[i, :lessthan] - return INFEASIBLE + return Boscia.INFEASIBLE end end - return OPTIMAL + return Boscia.OPTIMAL end function Boscia.is_valid_split(tree::Bonobo.BnBTree, blmo::CubeBLMO, vidx::Int) From b8e09270f0daf5ccbe054a7093d5cbc8257cb29f Mon Sep 17 00:00:00 2001 From: dhendryc <92737336+dhendryc@users.noreply.github.com> Date: Thu, 2 Nov 2023 19:07:34 +0100 Subject: [PATCH 13/15] Update approx_planted_point.jl --- examples/approx_planted_point.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/examples/approx_planted_point.jl b/examples/approx_planted_point.jl index 32da2bcfd..04a709ea0 100644 --- a/examples/approx_planted_point.jl +++ b/examples/approx_planted_point.jl @@ -1,6 +1,5 @@ using Boscia using FrankWolfe -using Bonobo using Test using Random using SCIP From 5ad7da0f37d2475f50f6297b4accc4a32a48e30c Mon Sep 17 00:00:00 2001 From: Deborah Hendrych Date: Fri, 3 Nov 2023 13:24:59 +0100 Subject: [PATCH 14/15] Relax tolerance. --- test/time_limit.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/time_limit.jl b/test/time_limit.jl index 429ca7b96..8f065eacb 100644 --- a/test/time_limit.jl +++ b/test/time_limit.jl @@ -63,5 +63,5 @@ time_limit = 30.0 time_taken = float(Dates.value(Dates.now() - start_time)) / 1000 @test sum(ai' * x) <= bi + 1e-6 @test f(x) <= f(result[:raw_solution]) + 1e-6 - @test time_taken <= time_limit + 10 + @test time_taken <= time_limit + 15 end From b127a3a668b8bfa1ab61fc02018d31d6ff6467e1 Mon Sep 17 00:00:00 2001 From: Deborah Hendrych Date: Fri, 3 Nov 2023 13:27:48 +0100 Subject: [PATCH 15/15] Check against the time Boscia itself measures. --- test/time_limit.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/time_limit.jl b/test/time_limit.jl index 8f065eacb..65cafbb3a 100644 --- a/test/time_limit.jl +++ b/test/time_limit.jl @@ -63,5 +63,5 @@ time_limit = 30.0 time_taken = float(Dates.value(Dates.now() - start_time)) / 1000 @test sum(ai' * x) <= bi + 1e-6 @test f(x) <= f(result[:raw_solution]) + 1e-6 - @test time_taken <= time_limit + 15 + @test result[:total_time_in_sec] <= time_limit + 5 end