Документация Engee

Dualization

Цель данного руководства — объяснить использование Dualization.jl для повышения производительности некоторых моделей конической оптимизации. Важно отметить два момента:

  1. JuMP переформулирует задачи в соответствии с входными требованиями

решателя, из-за чего размер задачи может увеличиваться вследствие добавления ослабляющих переменных и ограничений.

  1. Решение двойственной формы конической модели может быть более эффективным, чем

прямой.

Пакет Dualization.jl устраняет эти проблемы, позволяя решать двойственную задачу вместо прямой, изменив всего одну строку кода.

В этом руководстве используются следующие пакеты.

using JuMP
import Dualization
import SCS

Вводная информация

Решатели конической оптимизации обычно принимают одну из двух входных формулировок.

Первая — это стандартная коническая форма:

в которой имеется набор линейных ограничений равенства , а переменные принадлежат конусу .

Вторая — это геометрическая коническая форма:

в которой аффинная функция принадлежит конусу , а переменные свободны.

Преобразование между этими двумя формами не представляет никакой сложности. Например, чтобы перейти от геометрической конической формы к стандартной, мы вводим ослабляющие переменные :

А чтобы перейти от стандартной конической формы к геометрической, можно переписать ограничение равенства как функцию, принадлежащую конусу {0}:

С теоретической точки зрения эти две формулировки эквивалентны, и если реализовать модель в стандартной конической форме и передать ее решателю геометрической конической формы (или наоборот), то JuMP автоматически переформулирует задачу в правильную форму.

Однако на практике переформулировки вызывают проблемы, так как дополнительные ослабляющие переменные и ограничения могут значительно усложнить задачу и, следовательно, затруднить ее решение.

Также следует отметить, что многие задачи содержат комбинацию конических ограничений и переменных, поэтому их нельзя однозначно отнести ни к одной из двух формулировок. В таких случаях JuMP переформулирует только те переменные и ограничения, которые необходимы для преобразования задачи в требуемую форму.

Прямая и двойственная формулировки

Двойственность играет большую роль в конической оптимизации. Подробное описание конической двойственности см. в разделе Duality.

Полезно будет отметить, что если прямая задача имеет стандартную коническую форму, то двойственная задача имеет геометрическую коническую форму, и наоборот. Более того, прямая и двойственная задачи могут содержать разное количество переменных и ограничений, причем то, где их будет меньше, зависит от задачи. Поэтому вместо переформулирования задачи из одной формы в другую может быть эффективнее решить двойственную задачу вместо прямой.

Для демонстрации мы используем вариацию примера из раздела Maximum cut via SDP.

Прямая формулировка (в стандартной конической форме):

model_primal = Model()
@variable(model_primal, X[1:2, 1:2], PSD)
@objective(model_primal, Max, sum([1 -1; -1 1] .* X))
@constraint(model_primal, primal_c[i = 1:2], 1 - X[i, i] == 0)
print(model_primal)
Max X[1,1] - 2 X[1,2] + X[2,2]
Subject to
 primal_c[1] : -X[1,1] = -1
 primal_c[2] : -X[2,2] = -1
 [X[1,1], X[1,2], X[2,2]] ∈ MathOptInterface.PositiveSemidefiniteConeTriangle(2)

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

Двойственная по отношению к model_primal задача:

model_dual = Model()
@variable(model_dual, y[1:2])
@objective(model_dual, Min, sum(y))
@constraint(model_dual, dual_c, [y[1]-1 1; 1 y[2]-1] in PSDCone())
print(model_dual)
Min y[1] + y[2]
Subject to
 dual_c : [y[1] - 1  1
  1         y[2] - 1] ∈ PSDCone()

В этой задаче две скалярные переменные решения и ограничение в виде положительно полуопределенной матрицы 2x2.

Если вы раньше не сталкивались с конической двойственностью, попробуйте вывести двойственную задачу на основе описания в разделе Duality. При этом нужно знать, что двойственный конус по отношению к PSDCone — это PSDCone.

При решении model_primal с помощью SCS.Optimizer SCS сообщает о трех переменных (variables n: 3), пяти строках в матрице ограничений (constraints m: 5) и пяти ненулевых элементах в матрице (nnz(A): 5):

set_optimizer(model_primal, SCS.Optimizer)
optimize!(model_primal)
@assert is_solved_and_feasible(model_primal; dual = true)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 3, constraints m: 5
cones: 	  z: primal zero / dual free vars: 2
	  s: psd vars: 3, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 5, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.65e+01  1.60e-01  5.09e+01 -2.91e+01  1.00e-01  2.18e-03
    50| 1.74e-08  2.70e-10  4.88e-08 -4.00e+00  1.00e-01  2.28e-03
------------------------------------------------------------------
status:  solved
timings: total: 2.28e-03s = setup: 1.44e-04s + solve: 2.14e-03s
	 lin-sys: 1.36e-05s, cones: 1.88e-04s, accel: 3.87e-06s
------------------------------------------------------------------
objective = -4.000000
------------------------------------------------------------------

(В матрице ограничений пять строк по той причине, что SCS ожидает задачи в геометрической конической форме и поэтому в JuMP переменное ограничение X, PSD было переформулировано в аффинное ограничение X .+ 0 in PSDCone().)

Получаем следующее решение:

value.(X)
2×2 Matrix{Float64}:
  1.0  -1.0
 -1.0   1.0
dual.(primal_c)
2-element Vector{Float64}:
 1.9999999997299098
 1.9999999997299098
objective_value(model_primal)
3.9999999506359964

При решении model_dual с помощью SCS.Optimizer SCS сообщает о двух переменных (variables n: 2), трех строках в матрице ограничений (constraints m: 3) и двух ненулевых элементах в матрице (nnz(A): 2):

set_optimizer(model_dual, SCS.Optimizer)
optimize!(model_dual)
@assert is_solved_and_feasible(model_dual; dual = true)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 2, constraints m: 3
cones: 	  s: psd vars: 3, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 2, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.23e+01  1.00e+00  2.73e+01 -9.03e+00  1.00e-01  2.94e-04
    50| 1.13e-07  1.05e-09  3.23e-07  4.00e+00  1.00e-01  3.80e-04
------------------------------------------------------------------
status:  solved
timings: total: 3.81e-04s = setup: 6.18e-05s + solve: 3.19e-04s
	 lin-sys: 8.08e-06s, cones: 5.23e-05s, accel: 2.94e-06s
------------------------------------------------------------------
objective = 4.000000
------------------------------------------------------------------

Решение получается таким:

dual.(dual_c)
2×2 Matrix{Float64}:
  1.0  -1.0
 -1.0   1.0
value.(y)
2-element Vector{Float64}:
 2.0000001592726786
 2.0000001592726786
objective_value(model_dual)
4.000000318545357

Данная задача настолько мала, что сравнивать время ее решения не имеет смысла, но в целом следует ожидать, что форма model_dual решится быстрее, чем model_primal, так как она содержит меньше переменных и ограничений. Разница особенно заметна при решении масштабных задач оптимизации.

dual_optimizer

Вывод конической двойственной задачи вручную сложен, и высока вероятность ошибок. Пакет Dualization.jl предоставляет мета-решатель Dualization.dual_optimizer, который заключает любой совместимый с MathOptInterface решатель в интерфейс, автоматически формулирующий и решающий двойственную по отношению к входной задаче.

Для демонстрации мы используем Dualization.dual_optimizer для решения model_primal:

set_optimizer(model_primal, Dualization.dual_optimizer(SCS.Optimizer))
optimize!(model_primal)
@assert is_solved_and_feasible(model_primal; dual = true)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 2, constraints m: 3
cones: 	  s: psd vars: 3, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 2, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.23e+01  1.00e+00  2.73e+01 -9.03e+00  1.00e-01  3.55e-04
    50| 1.13e-07  1.05e-09  3.23e-07  4.00e+00  1.00e-01  4.35e-04
------------------------------------------------------------------
status:  solved
timings: total: 4.36e-04s = setup: 5.39e-05s + solve: 3.82e-04s
	 lin-sys: 8.55e-06s, cones: 6.48e-05s, accel: 3.15e-06s
------------------------------------------------------------------
objective = 4.000000
------------------------------------------------------------------

Производительность такая же, как при решении model_dual, а в X возвращается правильное решение:

value.(X)
2×2 Matrix{Float64}:
  1.0  -1.0
 -1.0   1.0
dual.(primal_c)
2-element Vector{Float64}:
 2.0000001592726786
 2.0000001592726786

Более того, если мы используем dual_optimizer применительно к model_dual, то получим ту же производительность, что и при решении model_primal:

set_optimizer(model_dual, Dualization.dual_optimizer(SCS.Optimizer))
optimize!(model_dual)
@assert is_solved_and_feasible(model_dual; dual = true)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 3, constraints m: 5
cones: 	  z: primal zero / dual free vars: 2
	  s: psd vars: 3, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 5, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.65e+01  1.60e-01  5.09e+01 -2.91e+01  1.00e-01  3.46e-04
    50| 1.74e-08  2.70e-10  4.88e-08 -4.00e+00  1.00e-01  4.24e-04
------------------------------------------------------------------
status:  solved
timings: total: 4.25e-04s = setup: 5.45e-05s + solve: 3.71e-04s
	 lin-sys: 1.04e-05s, cones: 5.56e-05s, accel: 3.21e-06s
------------------------------------------------------------------
objective = -4.000000
------------------------------------------------------------------
dual.(dual_c)
2×2 Matrix{Float64}:
  1.0  -1.0
 -1.0   1.0
value.(y)
2-element Vector{Float64}:
 1.9999999997299098
 1.9999999997299098

Смешанный пример

Пример из раздела Maximum cut via SDP хорошо поддается определению, так как прямая задача имеет стандартную коническую форму, а двойственная — геометрическую коническую форму. Однако на практике многие модели содержат комбинацию этих двух формулировок. Примером может служить задача из раздела The minimum distortion problem:

D = [0 1 1 1; 1 0 2 2; 1 2 0 2; 1 2 2 0]
model = Model()
@variable(model, c²)
@variable(model, Q[1:4, 1:4], PSD)
@objective(model, Min, c²)
for i in 1:4, j in (i+1):4
    @constraint(model, D[i, j]^2 <= Q[i, i] + Q[j, j] - 2 * Q[i, j])
    @constraint(model, Q[i, i] + Q[j, j] - 2 * Q[i, j] <= c² * D[i, j]^2)
end
@constraint(model, Q[1, 1] == 0)
@constraint(model, c² >= 1)
c² ≥ 1

В этой формулировке переменная Q имеет вид , но также имеются свободная переменная c², линейное ограничение равенства Q[1, 1] == 0 и ряд линейных ограничений неравенства. Вместо того чтобы пытаться вывести формулировку, которую JuMP передаст в SCS, и двойственную по отношению к ней формулировку, простейшим решением будет попытаться решить задачу с использованием dual_optimizer и без него, чтобы узнать, какая формулировка эффективнее.

set_optimizer(model, SCS.Optimizer)
optimize!(model)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 11, constraints m: 24
cones: 	  z: primal zero / dual free vars: 1
	  l: linear vars: 13
	  s: psd vars: 10, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 54, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 4.73e+00  1.00e+00  2.92e+00  1.23e+00  1.00e-01  3.44e-04
   150| 1.01e-04  3.07e-05  6.08e-05  1.33e+00  1.00e-01  1.09e-03
------------------------------------------------------------------
status:  solved
timings: total: 1.09e-03s = setup: 9.65e-05s + solve: 9.91e-04s
	 lin-sys: 1.16e-04s, cones: 4.39e-04s, accel: 1.30e-04s
------------------------------------------------------------------
objective = 1.333363
------------------------------------------------------------------
set_optimizer(model, Dualization.dual_optimizer(SCS.Optimizer))
optimize!(model)
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 14, constraints m: 24
cones: 	  z: primal zero / dual free vars: 1
	  l: linear vars: 13
	  s: psd vars: 10, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 57, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 3.71e+01  1.48e+00  2.23e+02 -1.13e+02  1.00e-01  3.86e-04
   150| 1.57e-04  2.28e-05  2.08e-04 -1.33e+00  1.00e-01  1.05e-03
------------------------------------------------------------------
status:  solved
timings: total: 1.05e-03s = setup: 9.79e-05s + solve: 9.52e-04s
	 lin-sys: 1.33e-04s, cones: 4.24e-04s, accel: 3.31e-05s
------------------------------------------------------------------
objective = -1.333460
------------------------------------------------------------------

Для этой задачи SCS сообщает, что в прямой задаче имеется variables n: 11, constraints m: 24, а в двойственной — variables n: 14, constraints m: 24. Поэтому, вероятно, следует использовать прямую формулировку, так как в ней меньше переменных и столько же ограничений.

Когда следует использовать dual_optimizer

Определенного правила касательно того, когда следует использовать dual_optimizer, не существует, так как в зависимости от конкретной задачи и выбора решателя задача может становиться больше или меньше. Однако вы можете попробовать использовать dual_optimizer, если решение задачи конической оптимизации занимает много времени или необходимо многократно решать задачи с похожей структурой, но с разными данными. В некоторых случаях решение двойственной задачи вместо прямой может давать существенный выигрыш.