add dual relative entropy cone and bridge

docs/src/manual/standard_form.md

@@ -75,13 +75,14 @@
 | [`RotatedSecondOrderCone(d)`](@ref) | ``\{ (t,u,x) \in \mathbb{R}^{d} : 2tu \ge \lVert x \rVert_2^2, t \ge 0,u \ge 0 \}``                                                        |
 | [`ExponentialCone()`](@ref)         | ``\{ (x,y,z) \in \mathbb{R}^3 : y \exp (x/y) \le z, y > 0 \}``                                                                             |
 | [`DualExponentialCone()`](@ref)     | ``\{ (u,v,w) \in \mathbb{R}^3 : -u \exp (v/u) \le \exp(1) w, u < 0 \}``                                                                    |
-| [`GeometricMeanCone(d)`](@ref)      | ``\{ (t,x) \in \mathbb{R}^{1+n} : x \ge 0, t \le \sqrt[n]{x_1 x_2 \cdots x_n} \}`` where ``n`` is ``d - 1``                                |
-| [`DualGeometricMeanCone(d)`](@ref)  | ``\{ (u,v) \in \mathbb{R}^{1+n} : v \ge 0, 0 \ge u \ge -n \sqrt[n]{\prod_i v_i} \}``, where ``n`` is ``d - 1``                             |
+| [`GeometricMeanCone(d)`](@ref)      | ``\{ (t,x) \in \mathbb{R}^{1+n} : x \ge 0, t \le \sqrt[n]{\prod_i x_i} \}`` where ``n`` is ``d - 1``                                |
+| [`DualGeometricMeanCone(d)`](@ref)  | ``\{ (t,x) \in \mathbb{R}^{1+n} : x \ge 0, 0 \ge t \ge -n \sqrt[n]{\prod_i x_i} \}`` where ``n`` is ``d - 1``                             |
 | [`PowerCone(α)`](@ref)              | ``\{ (x,y,z) \in \mathbb{R}^3 : x^{\alpha} y^{1-\alpha} \ge \|z\|, x \ge 0,y \ge 0 \}``                                                    |
 | [`DualPowerCone(α)`](@ref)          | ``\{ (u,v,w) \in \mathbb{R}^3 : \left(\frac{u}{\alpha}\right)^{\alpha}\left(\frac{v}{1-\alpha}\right)^{1-\alpha} \ge \|w\|, u,v \ge 0 \}`` |
 | [`NormOneCone(d)`](@ref)            | ``\{ (t,x) \in \mathbb{R}^{d} : t \ge \sum_i \lvert x_i \rvert \}``                                                                        |
 | [`NormInfinityCone(d)`](@ref)       | ``\{ (t,x) \in \mathbb{R}^{d} : t \ge \max_i \lvert x_i \rvert \}``                                                                        |
-| [`RelativeEntropyCone(d)`](@ref)    | ``\{ (u, v, w) \in \mathbb{R}^{d} : u \ge \sum_i w_i \log (\frac{w_i}{v_i}), v_i \ge 0, w_i \ge 0 \}``                                     |
+| [`RelativeEntropyCone(d)`](@ref)    | ``\{ (u, v, w) \in \mathbb{R}^{d} : u \ge \sum_i w_i \log (\frac{w_i}{v_i}), v_i > 0, w_i > 0 \}``                                     |
+| [`DualRelativeEntropyCone(d)`](@ref)| ``\{ (u, v, w) \in \mathbb{R}^{d} : \forall i, w_i \ge u(\log(\frac{u}{v_i}) - 1), v_i > 0, u > 0 \}``                                   |
 | [`HyperRectangle(l, u)`](@ref)      | ``\{ x \in \bar{\mathbb{R}}^d: x_i \in [l_i, u_i] \forall i=1,\ldots,d \}``                                                                |
 | [`NormCone(p, d)`](@ref)            | ``\{ (t,x) \in \mathbb{R}^{d} : t \ge \left(\sum\limits_i \lvert x_i \rvert^p\right)^{\frac{1}{p}} \}``                                    |
 | [`VectorNonlinearOracle`](@ref)     | ``\{x \in \mathbb{R}^{dimension}: l \le f(x) \le u \}``                                                                                    |

docs/src/reference/standard_form.md

@@ -94,6 +94,7 @@ DualExponentialCone
 PowerCone
 DualPowerCone
 RelativeEntropyCone
+DualRelativeEntropyCone
 NormSpectralCone
 NormNuclearCone
 SOS1

src/Bridges/Constraint/Constraint.jl

@@ -40,6 +40,7 @@ function add_all_bridges(model, ::Type{T}) where {T}
     MOI.Bridges.add_bridge(model, CountDistinctToMILPBridge{T})
     MOI.Bridges.add_bridge(model, CountGreaterThanToMILPBridge{T})
     MOI.Bridges.add_bridge(model, DualGeoMeanBridge{T})
+    MOI.Bridges.add_bridge(model, DualRelativeEntropyBridge{T})
     # * ExponentialConeToScalarNonlinearFunctionBridge{T}
     #      This bridge is not added by default because it starts with a convex
     #      conic constraint and adds a nonlinear constraint that local NLP

src/Bridges/Constraint/bridges/DualRelativeEntropyBridge.jl

@@ -0,0 +1,172 @@
+# Copyright (c) 2017: Miles Lubin and contributors
+# Copyright (c) 2017: Google Inc.
+#
+# Use of this source code is governed by an MIT-style license that can be found
+# in the LICENSE.md file or at https://opensource.org/licenses/MIT.
+
+"""
+    DualRelativeEntropyBridge{T,F,G} <: Bridges.Constraint.AbstractBridge
+
+`DualRelativeEntropyBridge` implements the following reformulation that converts a
+[`MOI.DualRelativeEntropyCone`](@ref) into an [`MOI.DualExponentialCone`](@ref):
+
+  * ``w_i \\ge u \\left(\\log \\left(\\frac{u}{v_i}\\right) - 1\\right)`` into
+    ``(-u, w_i, v_i) \\in DualExponentialCone``.
+
+## Source node
+
+`DualRelativeEntropyBridge` supports:
+
+  * `G` in [`MOI.DualRelativeEntropyCone`](@ref)
+
+## Target nodes
+
+`DualRelativeEntropyBridge` creates:
+
+  * `F` in [`MOI.DualExponentialCone`](@ref)
+"""
+struct DualRelativeEntropyBridge{T,F,G} <: AbstractBridge
+    exp_indices::Vector{MOI.ConstraintIndex{F,MOI.DualExponentialCone}}
+end
+
+const DualRelativeEntropy{T,OT<:MOI.ModelLike} =
+    SingleBridgeOptimizer{DualRelativeEntropyBridge{T},OT}
+
+function bridge_constraint(
+    ::Type{DualRelativeEntropyBridge{T,F,G}},
+    model::MOI.ModelLike,
+    f::G,
+    s::MOI.DualRelativeEntropyCone,
+) where {T,F,G}
+    f_scalars = MOI.Utilities.eachscalar(f)
+    d = MOI.dimension(s)
+    v_dim = div(d - 1, 2)
+    exp_indices = MOI.ConstraintIndex{F,MOI.DualExponentialCone}[]
+    for i in 1:v_dim
+        fi = MOI.Utilities.operate(
+            vcat,
+            T,
+            MOI.Utilities.operate(-, T, f_scalars[1]), # -u
+            f_scalars[1+i+v_dim],                      # w_i
+            f_scalars[1+i],                            # v_i
+        )
+        push!(
+            exp_indices,
+            MOI.add_constraint(model, fi, MOI.DualExponentialCone()),
+        )
+    end
+    return DualRelativeEntropyBridge{T,F,G}(exp_indices)
+end
+
+function MOI.supports_constraint(
+    ::Type{DualRelativeEntropyBridge{T}},
+    ::Type{<:MOI.AbstractVectorFunction},
+    ::Type{MOI.DualRelativeEntropyCone},
+) where {T}
+    return true
+end
+
+function MOI.Bridges.added_constrained_variable_types(
+    ::Type{<:DualRelativeEntropyBridge},
+)
+    return Tuple{Type}[]
+end
+
+function MOI.Bridges.added_constraint_types(
+    ::Type{DualRelativeEntropyBridge{T,F,G}},
+) where {T,F,G}
+    return Tuple{Type,Type}[(F, MOI.DualExponentialCone)]
+end
+
+function concrete_bridge_type(
+    ::Type{<:DualRelativeEntropyBridge{T}},
+    G::Type{<:MOI.AbstractVectorFunction},
+    ::Type{MOI.DualRelativeEntropyCone},
+) where {T}
+    S = MOI.Utilities.scalar_type(G)
+    U = MOI.Utilities.promote_operation(-, T, S, S)
+    F = MOI.Utilities.promote_operation(vcat, T, U, S)
+    return DualRelativeEntropyBridge{T,F,G}
+end
+
+function MOI.get(
+    bridge::DualRelativeEntropyBridge{T,F},
+    ::MOI.NumberOfConstraints{F,MOI.DualExponentialCone},
+)::Int64 where {T,F}
+    return length(bridge.exp_indices)
+end
+
+function MOI.get(
+    bridge::DualRelativeEntropyBridge{T,F},
+    ::MOI.ListOfConstraintIndices{F,MOI.DualExponentialCone},
+) where {T,F}
+    return copy(bridge.exp_indices)
+end
+
+function MOI.delete(model::MOI.ModelLike, bridge::DualRelativeEntropyBridge)
+    for exp_index_i in bridge.exp_indices
+        MOI.delete(model, exp_index_i)
+    end
+    return
+end
+
+function MOI.get(
+    model::MOI.ModelLike,
+    ::MOI.ConstraintFunction,
+    bridge::DualRelativeEntropyBridge{T,F,G},
+) where {T,F,G}
+    d = length(bridge.exp_indices)
+    func = MOI.Utilities.zero_with_output_dimension(F, 1 + 2 * d)
+    exp_1 = MOI.Utilities.eachscalar(
+        MOI.get(model, MOI.ConstraintFunction(), bridge.exp_indices[begin]),
+    )
+    MOI.Utilities.operate_output_index!(-, T, 1, func, exp_1[1])
+    for i in eachindex(bridge.exp_indices)
+        exp_i = MOI.Utilities.eachscalar(
+            MOI.get(model, MOI.ConstraintFunction(), bridge.exp_indices[i]),
+        )
+        MOI.Utilities.operate_output_index!(+, T, 1 + i, func, exp_i[3])
+        MOI.Utilities.operate_output_index!(+, T, 1 + d + i, func, exp_i[2])
+    end
+    return MOI.Utilities.convert_approx(G, func)
+end
+
+function MOI.get(
+    ::MOI.ModelLike,
+    ::MOI.ConstraintSet,
+    bridge::DualRelativeEntropyBridge,
+)
+    return MOI.DualRelativeEntropyCone(1 + 2 * length(bridge.exp_indices))
+end
+
+function MOI.get(
+    model::MOI.ModelLike,
+    attr::MOI.ConstraintPrimal,
+    bridge::DualRelativeEntropyBridge{T},
+) where {T}
+    d = length(bridge.exp_indices)
+    primal = zeros(T, 1 + 2d)
+    primal[1] = -MOI.get(model, attr, bridge.exp_indices[begin])[1]
+    for i in eachindex(bridge.exp_indices)
+        primal_i = MOI.get(model, attr, bridge.exp_indices[i])
+        primal[1+i] = primal_i[3]
+        primal[1+i+d] = primal_i[2]
+    end
+    return primal
+end
+
+function MOI.get(
+    model::MOI.ModelLike,
+    attr::MOI.ConstraintDual,
+    bridge::DualRelativeEntropyBridge{T},
+) where {T}
+    d = length(bridge.exp_indices)
+    dual = zeros(T, 1 + 2d)
+    for (i, ci) in enumerate(bridge.exp_indices)
+        dual_i = MOI.get(model, attr, ci)
+        dual[1] -= dual_i[1]
+        dual[1+i] = dual_i[3]
+        dual[1+i+d] = dual_i[2]
+    end
+    return dual
+end

src/Test/test_basic_constraint.jl

@@ -123,6 +123,7 @@ function _set(::Type{T}, ::Type{MOI.DualPowerCone}) where {T}
     return MOI.DualPowerCone(T(1//2))
 end
 _set(::Type{MOI.RelativeEntropyCone}) = MOI.RelativeEntropyCone(3)
+_set(::Type{MOI.DualRelativeEntropyCone}) = MOI.DualRelativeEntropyCone(3)
 _set(::Type{MOI.NormSpectralCone}) = MOI.NormSpectralCone(2, 3)
 _set(::Type{MOI.NormNuclearCone}) = MOI.NormNuclearCone(2, 3)
 function _set(::Type{MOI.PositiveSemidefiniteConeTriangle})
@@ -410,6 +411,7 @@ for s in [
     :PowerCone,
     :DualPowerCone,
     :RelativeEntropyCone,
+    :DualRelativeEntropyCone,
     :NormSpectralCone,
     :NormNuclearCone,
     :PositiveSemidefiniteConeSquare,

src/Test/test_conic.jl

@@ -4469,11 +4469,6 @@ function test_conic_RelativeEntropyCone(
     model::MOI.ModelLike,
     config::Config{T},
 ) where {T<:Real}
-    @requires MOI.supports_constraint(
-        model,
-        MOI.VectorAffineFunction{T},
-        MOI.RelativeEntropyCone,
-    )
     @requires MOI.supports_incremental_interface(model)
     @requires MOI.supports(model, MOI.ObjectiveFunction{MOI.VariableIndex}())
     @requires MOI.supports(model, MOI.ObjectiveSense())
@@ -4541,6 +4536,106 @@ function setup_test(
     return
 end
 
+"""
+    test_conic_DualRelativeEntropyCone(
+        model::MOI.ModelLike,
+        config::Config{T},
+    ) where {T<:Real}
+
+Test the problem:
+```
+max -dot(v, [1, 5]) - dot(w, [2, 3])
+ st  w[i] >= log(1/v[i]) - 1  (that is, (1, v, w) in DualRelativeEntropyCone(5))
+Optimal solution is:
+u = 2*log(2/1) + 3*log(3/5) ≈ -0.1461825
+```
+"""
+function test_conic_DualRelativeEntropyCone(
+    model::MOI.ModelLike,
+    config::Config{T},
+) where {T<:Real}
+    @requires MOI.supports_incremental_interface(model)
+    @requires MOI.supports(
+        model,
+        MOI.ObjectiveFunction{MOI.ScalarAffineFunction{T}}(),
+    )
+    @requires MOI.supports(model, MOI.ObjectiveSense())
+    @requires MOI.supports_constraint(
+        model,
+        MOI.VectorAffineFunction{T},
+        MOI.DualRelativeEntropyCone,
+    )
+    v = MOI.add_variables(model, 2)
+    w = MOI.add_variables(model, 2)
+    @test MOI.get(model, MOI.NumberOfVariables()) == 4
+    relentr = MOI.add_constraint(
+        model,
+        MOI.VectorAffineFunction(
+            [
+                MOI.VectorAffineTerm(2, MOI.ScalarAffineTerm(T(1), v[1])),
+                MOI.VectorAffineTerm(3, MOI.ScalarAffineTerm(T(1), v[2])),
+                MOI.VectorAffineTerm(4, MOI.ScalarAffineTerm(T(1), w[1])),
+                MOI.VectorAffineTerm(5, MOI.ScalarAffineTerm(T(1), w[2])),
+            ],
+            T[1, 0, 0, 0, 0],
+        ),
+        MOI.DualRelativeEntropyCone(5),
+    )
+    MOI.set(
+        model,
+        MOI.ObjectiveFunction{MOI.ScalarAffineFunction{T}}(),
+        MOI.ScalarAffineFunction(
+            MOI.ScalarAffineTerm.(T[-1, -5, -2, -3], [v[1], v[2], w[1], w[2]]),
+            T(0),
+        ),
+    )
+    MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE)
+    if _supports(config, MOI.optimize!)
+        @test MOI.get(model, MOI.TerminationStatus()) == MOI.OPTIMIZE_NOT_CALLED
+        MOI.optimize!(model)
+        @test MOI.get(model, MOI.TerminationStatus()) == config.optimal_status
+        @test MOI.get(model, MOI.PrimalStatus()) == MOI.FEASIBLE_POINT
+        if _supports(config, MOI.ConstraintDual)
+            @test MOI.get(model, MOI.DualStatus()) == MOI.FEASIBLE_POINT
+        end
+        u_opt = 2 * log(T(2)) + 3 * log(T(3 // 5))
+        @test ≈(MOI.get(model, MOI.ObjectiveValue()), u_opt, config)
+        @test ≈(
+            MOI.get(model, MOI.ConstraintPrimal(), relentr),
+            T[1, 2, 3//5, log(T(1//2))-1, log(T(5//3))-1],
+            config,
+        )
+        if _supports(config, MOI.ConstraintDual)
+            @test ≈(
+                MOI.get(model, MOI.ConstraintDual(), relentr),
+                T[u_opt, 1, 5, 2, 3],
+                config,
+            )
+        end
+    end
+    return
+end
+
+function setup_test(
+    ::typeof(test_conic_DualRelativeEntropyCone),
+    model::MOIU.MockOptimizer,
+    ::Config{T},
+) where {T<:Real}
+    u_opt = 2 * log(T(2 // 1)) + 3 * log(T(3 // 5))
+    MOIU.set_mock_optimize!(
+        model,
+        (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
+            mock,
+            T[2, 3//5, log(T(1//2))-1, log(T(5//3))-1],
+            (MOI.VectorAffineFunction{T}, MOI.DualRelativeEntropyCone) =>
+                [T[u_opt, 1, 5, 2, 3]],
+        ),
+    )
+    return
+end
+
+version_added(::typeof(test_conic_DualRelativeEntropyCone)) = v"1.53.0"
+
 """
     test_conic_NormSpectralCone(model::MOI.ModelLike, config::Config{T}) where {T<:Real}
 

src/Utilities/distance_to_set.jl

@@ -339,7 +339,7 @@ end
 Let `(u, v, w) = x`. If `u < 0`, return the epigraph distance `d` such that
 `(u, v, w + d)` belongs to the set.
 
-If `u >= 0` return the 2-norm of the vector `d` such that `x + d = (u, -1, z)`
+If `u >= 0` return the 2-norm of the vector `d` such that `x + d = (-1, v, z)`
 where `z` satisfies the constraints.
 """
 function distance_to_set(
@@ -494,7 +494,7 @@ end
 """
     distance_to_set(::ProjectionUpperBoundDistance, x, ::MOI.RelativeEntropyCone)
 
-Let `(u, v..., w...) = x`. If `v` and `w` are strictly positive, return the
+Let `(u, v, w) = x`. If `v` and `w` are strictly positive, return the
 epigraph distance required to increase `u` such that the constraint is satisfied.
 
 If any elements in `v` or `w` are non-positive, return the 2-norm of the vector
@@ -527,6 +527,34 @@ function distance_to_set(
     return max(sum(w[i] * log(w[i] / v[i]) for i in eachindex(w)) - u, zero(T))
 end
 
+"""
+    distance_to_set(::ProjectionUpperBoundDistance, x, ::MOI.DualRelativeEntropyCone)
+
+Let `(u, v, w) = x`. If `u` is strictly positive, return the 2-norm of the vector
+`d` such `(u, v + d, w)` satisfies the constraints.
+
+If `u` is non-positive return the 2-norm of the vector `d` such that
+`x + d = (1, z, w)` and `z` is such that the constraints are satisfied.
+"""
+function distance_to_set(
+    ::ProjectionUpperBoundDistance,
+    x::AbstractVector{T},
+    set::MOI.DualRelativeEntropyCone,
+) where {T<:Real}
+    Base.require_one_based_indexing(x)
+    _check_dimension(x, set)
+    n = div(MOI.dimension(set) - 1, 2)
+    u, v, w = x[1], @view(x[2:(n+1)]), @view(x[(n+2):end])
+    if u > 0
+        d = max.(Ref(zero(T)), u * exp.(-w ./ u .- 1) .- v)
+        return LinearAlgebra.norm(d)
+    else
+        d1 = 1 - u
+        drest = max.(Ref(zero(T)), exp.(-w .- 1) - v)
+        return LinearAlgebra.norm((d1, drest))
+    end
+end
+
 # This is the minimal L2-norm.
 function distance_to_set(
     ::ProjectionUpperBoundDistance,