diff --git a/README.md b/README.md index 8b4f784..069d098 100644 --- a/README.md +++ b/README.md @@ -51,6 +51,7 @@ Each algorithm supports only a subset of the attributes. Consult the algorithm's docstring for details on which attributes it supports, and how it uses them in the solution process. + * `MOO.ObjectiveAbsoluteTolerance(index::Int)` * `MOO.ObjectivePriority(index::Int)` * `MOO.ObjectiveRelativeTolerance(index::Int)` * `MOO.ObjectiveWeight(index::Int)` diff --git a/src/MOO.jl b/src/MOO.jl index ff9c1a7..1f9e602 100644 --- a/src/MOO.jl +++ b/src/MOO.jl @@ -194,6 +194,21 @@ end default(::ObjectiveRelativeTolerance) = 0.0 +""" + ObjectiveAbsoluteTolerance(index::Int) <: AbstractAlgorithmAttribute -> Float64 + +Assign a `Float64` tolerance to objective number `index`. This is most commonly +used to constrain an objective to a range in absolute terms to the optimal +objective value of that objective. + +Defaults to `0.0`. +""" +struct ObjectiveAbsoluteTolerance <: AbstractAlgorithmAttribute + index::Int +end + +default(::ObjectiveAbsoluteTolerance) = 0.0 + ### RawOptimizerAttribute function MOI.supports(model::Optimizer, attr::MOI.RawOptimizerAttribute) diff --git a/src/algorithms/EpsilonConstraint.jl b/src/algorithms/EpsilonConstraint.jl index 5c53428..0adf851 100644 --- a/src/algorithms/EpsilonConstraint.jl +++ b/src/algorithms/EpsilonConstraint.jl @@ -11,17 +11,25 @@ bi-objective programs. ## Supported optimizer attributes - * `MOO.SolutionLimit()`: with a slight abuse of notation, `EpsilonConstraint` - divides the width of the first-objective's domain in objective space by - `SolutionLimit` to obtain the epsilon to use when iterating. Thus, there - can be at most `SolutionLimit` solutions returned, but there may be fewer. - If no value is set, the default is `100`, instead of the typical - `default(::SolutionLimit)`. + * `MOO.ObjectiveAbsoluteTolerance(1)`: if set, `EpsilonConstraint` uses this + tolerance as the epsilon by which it partitions the first-objective's space. + If the objective is a pure integer program, set the tolerance to `1` to + enumerate all non-dominated solutions. Note that you can set only the + tolerance for the first objective index; the tolerances for other objective + indices are ignored. + + * `MOO.SolutionLimit()`: if `MOO.ObjectiveAbsoluteTolerance(1)` is not set + then, with a slight abuse of notation, `EpsilonConstraint` divides the width + of the first-objective's domain in objective space by `SolutionLimit` to + obtain the epsilon to use when iterating. Thus, there can be at most + `SolutionLimit` solutions returned, but there may be fewer. If no value is + set, the default is `100`, instead of the typical `default(::SolutionLimit)`. """ mutable struct EpsilonConstraint <: AbstractAlgorithm solution_limit::Union{Nothing,Int} + atol::Union{Nothing,Float64} - EpsilonConstraint() = new(nothing) + EpsilonConstraint() = new(nothing, nothing) end default(::EpsilonConstraint, ::SolutionLimit) = 100 @@ -37,6 +45,24 @@ function MOI.get(alg::EpsilonConstraint, attr::SolutionLimit) return something(alg.solution_limit, default(alg, attr)) end +MOI.supports(::EpsilonConstraint, ::ObjectiveAbsoluteTolerance) = true + +function MOI.set( + alg::EpsilonConstraint, + attr::ObjectiveAbsoluteTolerance, + value, +) + if attr.index == 1 + alg.atol = value + end + return +end + +function MOI.get(alg::EpsilonConstraint, attr::ObjectiveAbsoluteTolerance) + @assert attr.index == 1 + return something(alg.atol, default(alg, attr)) +end + function optimize_multiobjective!( algorithm::EpsilonConstraint, model::Optimizer, @@ -55,8 +81,11 @@ function optimize_multiobjective!( a, b = solution_1[1].y[1], solution_2[1].y[1] left, right = min(a, b), max(a, b) # Compute the epsilon that we will be incrementing by each iteration - n_points = MOI.get(algorithm, SolutionLimit()) - ε = abs(right - left) / (n_points - 1) + ε = MOI.get(algorithm, ObjectiveAbsoluteTolerance(1)) + if iszero(ε) + n_points = MOI.get(algorithm, SolutionLimit()) + ε = abs(right - left) / (n_points - 1) + end solutions = SolutionPoint[] f1, f2 = MOI.Utilities.eachscalar(model.f) MOI.set(model.inner, MOI.ObjectiveFunction{typeof(f2)}(), f2) @@ -66,10 +95,10 @@ function optimize_multiobjective!( MOI.LessThan{Float64}, MOI.GreaterThan{Float64}, ) - ci = MOI.add_constraint(model, f1, SetType(left + ε)) + ci = MOI.add_constraint(model, f1, SetType(left)) variables = MOI.get(model.inner, MOI.ListOfVariableIndices()) - for i in 1:n_points - rhs = left + (i - 1) * ε + rhs = left + while rhs <= right + ε / 2 MOI.set(model, MOI.ConstraintSet(), ci, SetType(rhs)) MOI.optimize!(model.inner) if MOI.get(model.inner, MOI.TerminationStatus()) != MOI.OPTIMAL @@ -83,6 +112,7 @@ function optimize_multiobjective!( if isempty(solutions) || !(Y ≈ solutions[end].y) push!(solutions, SolutionPoint(X, Y)) end + rhs += ε end MOI.delete(model, ci) return MOI.OPTIMAL, unique(solutions) diff --git a/test/algorithms/EpsilonConstraint.jl b/test/algorithms/EpsilonConstraint.jl index 256bd41..cc1b679 100644 --- a/test/algorithms/EpsilonConstraint.jl +++ b/test/algorithms/EpsilonConstraint.jl @@ -66,6 +66,50 @@ function test_biobjective_knapsack() return end +function test_biobjective_knapsack_atol() + p1 = [77, 94, 71, 63, 96, 82, 85, 75, 72, 91, 99, 63, 84, 87, 79, 94, 90] + p2 = [65, 90, 90, 77, 95, 84, 70, 94, 66, 92, 74, 97, 60, 60, 65, 97, 93] + w = [80, 87, 68, 72, 66, 77, 99, 85, 70, 93, 98, 72, 100, 89, 67, 86, 91] + model = MOO.Optimizer(HiGHS.Optimizer) + MOI.set(model, MOO.Algorithm(), MOO.EpsilonConstraint()) + MOI.set(model, MOO.ObjectiveAbsoluteTolerance(1), 1.0) + MOI.set(model, MOI.Silent(), true) + x = MOI.add_variables(model, length(w)) + MOI.add_constraint.(model, x, MOI.ZeroOne()) + MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE) + f = MOI.Utilities.operate( + vcat, + Float64, + [sum(1.0 * p[i] * x[i] for i in 1:length(w)) for p in [p1, p2]]..., + ) + MOI.set(model, MOI.ObjectiveFunction{typeof(f)}(), f) + MOI.add_constraint( + model, + sum(1.0 * w[i] * x[i] for i in 1:length(w)), + MOI.LessThan(900.0), + ) + MOI.optimize!(model) + results = Dict( + [955, 906] => [2, 3, 5, 6, 9, 10, 11, 14, 15, 16, 17], + [949, 915] => [1, 2, 5, 6, 8, 9, 10, 11, 15, 16, 17], + [948, 939] => [1, 2, 3, 5, 6, 8, 10, 11, 15, 16, 17], + [943, 940] => [2, 3, 5, 6, 8, 9, 10, 11, 15, 16, 17], + [936, 942] => [1, 2, 3, 5, 6, 10, 11, 12, 15, 16, 17], + [935, 947] => [2, 5, 6, 8, 9, 10, 11, 12, 15, 16, 17], + [934, 971] => [2, 3, 5, 6, 8, 10, 11, 12, 15, 16, 17], + [927, 972] => [2, 3, 5, 6, 8, 9, 10, 11, 12, 16, 17], + [918, 983] => [2, 3, 4, 5, 6, 8, 10, 11, 12, 16, 17], + ) + @test MOI.get(model, MOI.ResultCount()) == 9 + for i in 1:MOI.get(model, MOI.ResultCount()) + x_sol = MOI.get(model, MOI.VariablePrimal(i), x) + X = findall(elt -> elt > 0.9, x_sol) + Y = MOI.get(model, MOI.ObjectiveValue(i)) + @test results[round.(Int, Y)] == X + end + return +end + function test_biobjective_knapsack_min() p1 = [77, 94, 71, 63, 96, 82, 85, 75, 72, 91, 99, 63, 84, 87, 79, 94, 90] p2 = [65, 90, 90, 77, 95, 84, 70, 94, 66, 92, 74, 97, 60, 60, 65, 97, 93]