Джулия Джамп: Получение всех возможных решений для мип

JuMP поддерживает возврат нескольких решений.

Документация: https://jump.dev/JuMP.jl/stable/manual/solutions/#Multiple-solutions

Рабочий процесс примерно такой:

using JuMP
model = Model()
@variable(model, x[1:10] >= 0)
# ... other constraints ...
optimize!(model)

if termination_status(model) != OPTIMAL
    error("The model was not solved correctly.")
end

an_optimal_solution = value.(x; result = 1)
optimal_objective = objective_value(model; result = 1)
for i in 2:result_count(model)
    @assert has_values(model; result = i)
    println("Solution $(i) = ", value.(x; result = i))
    obj = objective_value(model; result = i)
    println("Objective $(i) = ", obj)
    if isapprox(obj, optimal_objective; atol = 1e-8)
        print("Solution $(i) is also optimal!")
    end
end

Но вам нужен решатель, который поддерживает возврат нескольких решений и настройку правильных параметров для конкретного решателя.

См. этот пост в блоге: https://jump.dev/tutorials/2021/11/02/tutorial-multi-jdf/

Ниже приведен пример поиска всех решений для логической задачи. С такими проблемами легче справляться, поскольку пространство решений легко перечисляется (даже несмотря на то, что оно может расти экспоненциально).

Во-первых, давайте получим пакеты и определим пример проблемы:

using Random, JuMP, HiGHS, MathOptInterface

function example_knapsack()
    profit = [5, 3, 2, 7, 4]
    weight = [2, 8, 4, 2, 5]
    capacity = 10
    minprofit = 10
    model = Model(HiGHS.Optimizer)
    set_silent(model)
    @variable(model, x[1:5], Bin)
    @objective(model, FEASIBILITY_SENSE, 0)
    @constraint(model, weight' * x <= capacity)
    @constraint(model, profit' * x >= minprofit)
    return model
end

(это задача о рюкзаке из документации JuMP).

Затем мы используем рекурсию для изучения дерева всех возможных решений. Дерево не спускается по ветвям без решения (поэтому время выполнения не всегда экспоненциально):

function findallsol(model, x)
    perm = shuffle(1:length(x))
    res = Vector{Float64}[]
    _findallsol!(res, model, x, perm, 0)
    return res
end

function _findallsol!(res, model, x, perm, depth)
    n = length(x)
    depth > n && return
    optimize!(model)
    if termination_status(model) == MathOptInterface.OPTIMAL
        if depth == n
            push!(res, value.(x))
            return
        else
            idx = perm[depth+1]
            v = value(x[idx])
            newcon = @constraint(model, x[idx] == v)
            _findallsol!(res, model, x, perm, depth + 1)
            delete(model, newcon)
            newcon = @constraint(model, x[idx] == 1 - v)
            _findallsol!(res, model, x, perm, depth + 1)
            delete(model, newcon)
        end
    end
    return
end

Теперь мы можем:

julia> m = example_knapsack()
A JuMP Model
Maximization problem with:
Variables: 5
...
Names registered in the model: x

julia> res = findallsol(m, m.obj_dict[:x])
5-element Vector{Vector{Float64}}:
 [1.0, 0.0, 0.0, 1.0, 1.0]
 [0.0, 0.0, 0.0, 1.0, 1.0]
 [1.0, 0.0, 1.0, 1.0, 0.0]
 [1.0, 0.0, 0.0, 1.0, 0.0]
 [0.0, 1.0, 0.0, 1.0, 0.0]

И получаем вектор со всеми решениями.

Если рассматриваемая проблема является булевой проблемой, этот метод можно использовать как есть. В случае, если у него есть небулевые переменные, рекурсия должна будет разделить допустимое пространство каким-то ровным образом. Например, выбор переменной и разрезание ее домена пополам, а также рекурсия к каждой половине с меньшим доменом этой переменной (чтобы гарантировать завершение).

P.S. Это не оптимальный метод. Эта проблема хорошо изучена. Возможные термины для поиска: «подсчет моделей» (особенно в логической области).

(ОБНОВЛЕНИЕ: цель изменена на использование FEASIBLE)