|
| 1 | +""" |
| 2 | + RelativeEntropyBridge{T} |
| 3 | +
|
| 4 | +The `RelativeEntropyCone` is representable with exponential cone and LP constraints, since |
| 5 | +``u \\ge \\sum_{i=1}^n w_i \\log (\\frac{w_i}{v_i})`` if and only if there exists a vector |
| 6 | +``y`` such that ``u \\ge \\sum_i y_i`` and ``y_i \\ge w_i \\log (\\frac{w_i}{v_i})`` or |
| 7 | +equivalently ``v_i \\ge w_i \\exp (\\frac{-y_i}{w_i})`` or equivalently |
| 8 | +``(-y_i, w_i, v_i) \\in ExponentialCone``, for all ``i``. |
| 9 | +""" |
| 10 | +struct RelativeEntropyBridge{T, F, G, H} <: AbstractBridge |
| 11 | + y::Vector{MOI.VariableIndex} |
| 12 | + ge_index::CI{F, MOI.GreaterThan{T}} |
| 13 | + exp_indices::Vector{CI{G, MOI.ExponentialCone}} |
| 14 | +end |
| 15 | +function bridge_constraint(::Type{RelativeEntropyBridge{T, F, G, H}}, model::MOI.ModelLike, f::H, s::MOI.RelativeEntropyCone) where {T, F, G, H} |
| 16 | + f_scalars = MOIU.eachscalar(f) |
| 17 | + d = MOI.dimension(s) |
| 18 | + v_dim = div(d - 1, 2) |
| 19 | + y = MOI.add_variables(model, v_dim) |
| 20 | + ge_index = MOIU.normalize_and_add_constraint(model, MOIU.operate(-, T, f_scalars[1], MOIU.operate(sum, T, y)), MOI.GreaterThan(zero(T)), allow_modify_function=true) |
| 21 | + w_start = 1 + v_dim |
| 22 | + exp_funcs = [MOIU.operate(vcat, T, MOIU.operate(-, T, MOI.SingleVariable(y[i])), f_scalars[w_start + i], f_scalars[1 + i]) for i in 1:v_dim] |
| 23 | + exp_indices = [MOI.add_constraint(model, exp_func_i, MOI.ExponentialCone()) for exp_func_i in exp_funcs] |
| 24 | + return RelativeEntropyBridge{T, F, G, H}(y, ge_index, exp_indices) |
| 25 | +end |
| 26 | + |
| 27 | +MOI.supports_constraint(::Type{RelativeEntropyBridge{T}}, ::Type{<:MOI.AbstractVectorFunction}, ::Type{MOI.RelativeEntropyCone}) where T = true |
| 28 | +MOIB.added_constrained_variable_types(::Type{<:RelativeEntropyBridge}) = Tuple{DataType}[] |
| 29 | +MOIB.added_constraint_types(::Type{RelativeEntropyBridge{T, F, G, H}}) where {T, F, G, H} = [(F, MOI.GreaterThan{T}), (G, MOI.ExponentialCone)] |
| 30 | +function concrete_bridge_type(::Type{<:RelativeEntropyBridge{T}}, H::Type{<:MOI.AbstractVectorFunction}, ::Type{MOI.RelativeEntropyCone}) where T |
| 31 | + S = MOIU.scalar_type(H) |
| 32 | + F = MOIU.promote_operation(-, T, S, S) |
| 33 | + G = MOIU.promote_operation(vcat, T, MOIU.promote_operation(-, T, MOI.SingleVariable), S) |
| 34 | + return RelativeEntropyBridge{T, F, G, H} |
| 35 | +end |
| 36 | + |
| 37 | +# Attributes, Bridge acting as a model |
| 38 | +MOI.get(bridge::RelativeEntropyBridge, ::MOI.NumberOfVariables) = length(bridge.y) |
| 39 | +MOI.get(bridge::RelativeEntropyBridge, ::MOI.ListOfVariableIndices) = bridge.y |
| 40 | +MOI.get(bridge::RelativeEntropyBridge{T, F, G, H}, ::MOI.NumberOfConstraints{F, MOI.GreaterThan{T}}) where {T, F, G, H} = 1 |
| 41 | +MOI.get(bridge::RelativeEntropyBridge{T, F, G, H}, ::MOI.NumberOfConstraints{G, MOI.ExponentialCone}) where {T, F, G, H} = length(bridge.y) |
| 42 | +MOI.get(bridge::RelativeEntropyBridge{T, F, G, H}, ::MOI.ListOfConstraintIndices{F, MOI.GreaterThan{T}}) where {T, F, G, H} = [bridge.ge_index] |
| 43 | +MOI.get(bridge::RelativeEntropyBridge{T, F, G, H}, ::MOI.ListOfConstraintIndices{G, MOI.ExponentialCone}) where {T, F, G, H} = bridge.exp_indices |
| 44 | + |
| 45 | +# References |
| 46 | +function MOI.delete(model::MOI.ModelLike, bridge::RelativeEntropyBridge) |
| 47 | + for exp_index_i in bridge.exp_indices |
| 48 | + MOI.delete(model, exp_index_i) |
| 49 | + end |
| 50 | + MOI.delete(model, bridge.ge_index) |
| 51 | + MOI.delete(model, bridge.y) |
| 52 | +end |
| 53 | + |
| 54 | +# Attributes, Bridge acting as a constraint |
| 55 | +function MOI.get(model::MOI.ModelLike, ::MOI.ConstraintFunction, bridge::RelativeEntropyBridge{T, F, G, H}) where {T, F, G, H} |
| 56 | + func = MOIU.zero_with_output_dimension(G, 1 + 2 * length(bridge.y)) |
| 57 | + MOIU.operate_output_index!(+, T, 1, func, MOI.get(model, MOI.ConstraintFunction(), bridge.ge_index)) |
| 58 | + w_start = 1 + length(bridge.y) |
| 59 | + for i in eachindex(bridge.y) |
| 60 | + exp_func_i = MOIU.eachscalar(MOI.get(model, MOI.ConstraintFunction(), bridge.exp_indices[i])) |
| 61 | + MOIU.operate_output_index!(-, T, 1, func, exp_func_i[1]) |
| 62 | + MOIU.operate_output_index!(+, T, 1 + i, func, exp_func_i[3]) |
| 63 | + MOIU.operate_output_index!(+, T, w_start + i, func, exp_func_i[2]) |
| 64 | + end |
| 65 | + return MOIU.convert_approx(H, MOIU.remove_variable(func, bridge.y)) |
| 66 | +end |
| 67 | +MOI.get(model::MOI.ModelLike, ::MOI.ConstraintSet, bridge::RelativeEntropyBridge) = MOI.RelativeEntropyCone(1 + 2 * length(bridge.y)) |
| 68 | +MOI.supports(::MOI.ModelLike, ::Union{MOI.ConstraintPrimalStart, MOI.ConstraintDualStart}, ::Type{<:RelativeEntropyBridge}) = true |
| 69 | +function MOI.get(model::MOI.ModelLike, attr::Union{MOI.ConstraintPrimal, MOI.ConstraintPrimalStart}, bridge::RelativeEntropyBridge{T}) where T |
| 70 | + primal = zeros(T, 1 + 2 * length(bridge.y)) |
| 71 | + primal[1] = MOI.get(model, attr, bridge.ge_index) |
| 72 | + w_start = 1 + length(bridge.y) |
| 73 | + for i in eachindex(bridge.y) |
| 74 | + exp_primal_i = MOI.get(model, attr, bridge.exp_indices[i]) |
| 75 | + primal[1] -= exp_primal_i[1] |
| 76 | + primal[1 + i] = exp_primal_i[3] |
| 77 | + primal[w_start + i] = exp_primal_i[2] |
| 78 | + end |
| 79 | + return primal |
| 80 | +end |
| 81 | +function MOI.set(model::MOI.ModelLike, attr::MOI.ConstraintPrimalStart, bridge::RelativeEntropyBridge{T}, value) where T |
| 82 | + v_dim = length(bridge.y) |
| 83 | + v_value = value[2:(v_dim + 1)] |
| 84 | + w_value = value[(v_dim + 2):end] |
| 85 | + y_value = [w_i * log(w_i / v_i) for (v_i, w_i) in zip(v_value, w_value)] |
| 86 | + MOI.set(model, attr, bridge.ge_index, value[1] - reduce(+, y_value, init=zero(T))) |
| 87 | + for i in 1:v_dim |
| 88 | + MOI.set(model, attr, bridge.exp_indices[i], [-y_value[i], w_value[i], v_value[i]]) |
| 89 | + MOI.set(model, MOI.VariablePrimalStart(), bridge.y[i], y_value[i]) |
| 90 | + end |
| 91 | + return |
| 92 | +end |
| 93 | +# Given a is dual on u - sum(y) >= 0 and (b_i, c_i, d_i) is dual on (-y_i, w_i, v_i) |
| 94 | +# in ExponentialCone, the dual on (u, v, w) in RelativeEntropyCone is (a, d_i, c_i). |
| 95 | +function MOI.get(model::MOI.ModelLike, attr::Union{MOI.ConstraintDual, MOI.ConstraintDualStart}, bridge::RelativeEntropyBridge{T}) where {T} |
| 96 | + dual = zeros(T, 1 + 2 * length(bridge.y)) |
| 97 | + dual[1] = MOI.get(model, attr, bridge.ge_index)[1] |
| 98 | + w_start = 1 + length(bridge.y) |
| 99 | + for i in eachindex(bridge.y) |
| 100 | + exp_dual_i = MOI.get(model, attr, bridge.exp_indices[i]) |
| 101 | + dual[1 + i] = exp_dual_i[3] |
| 102 | + dual[w_start + i] = exp_dual_i[2] |
| 103 | + end |
| 104 | + return dual |
| 105 | +end |
| 106 | +# Given constraint dual start of (u, v, w), constraint dual on GreaterThan is u |
| 107 | +# and on exponential cone constraint i is (r_i, w_i, v_i), but since y_i is free, |
| 108 | +# its dual is 0, so we have -r_i + value[1] == 0 hence r_i = value[1]. |
| 109 | +# Note: alternatively, we could use the Lambert W function to calculate |
| 110 | +# r_i = exp(W(-w_i / (-v_i * e))) * (-v_i * e), but this is more complicated. |
| 111 | +function MOI.set(model::MOI.ModelLike, ::MOI.ConstraintDualStart, bridge::RelativeEntropyBridge, value) |
| 112 | + MOI.set(model, MOI.ConstraintDualStart(), bridge.ge_index, value[1]) |
| 113 | + w_start = 1 + length(bridge.y) |
| 114 | + for i in eachindex(bridge.y) |
| 115 | + MOI.set(model, MOI.ConstraintDualStart(), bridge.exp_indices[i], [value[1], value[w_start + i], value[1 + i]]) |
| 116 | + end |
| 117 | + return |
| 118 | +end |
0 commit comments