Add tests

Author: blegat

docs/src/background/duality.md

@@ -114,7 +114,7 @@ and similarly, the dual is:
 The scalar product is different from the canonical one for the sets
 [`PositiveSemidefiniteConeTriangle`](@ref), [`LogDetConeTriangle`](@ref),
 [`RootDetConeTriangle`](@ref),
-[`SetWithDotProducts`](@ref) and
+[`SetDotProducts`](@ref) and
 [`LinearCombinationInSet`](@ref).
 
 If the set ``C_i`` of the section [Duality](@ref) is one of these three cones,

docs/src/manual/standard_form.md

@@ -82,7 +82,7 @@ The vector-valued set types implemented in MathOptInterface.jl are:
 | [`RelativeEntropyCone(d)`](@ref MathOptInterface.RelativeEntropyCone) | ``\{ (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 \}`` |
 | [`HyperRectangle(l, u)`](@ref MathOptInterface.HyperRectangle)        | ``\{x \in \bar{\mathbb{R}}^d: x_i \in [l_i, u_i] \forall i=1,\ldots,d\}`` |
 | [`NormCone(p, d)`](@ref MathOptInterface.NormCone)       | ``\{ (t,x) \in \mathbb{R}^{d} : t \ge \left(\sum\limits_i \lvert x_i \rvert^p\right)^{\frac{1}{p}} \}`` |
-| [`SetWithDotProducts(s, v)`](@ref MathOptInterface.SetWithDotProducts) | The cone `s` with dot products with the fixed vectors `v`. |
+| [`SetDotProducts(s, v)`](@ref MathOptInterface.SetDotProducts) | The cone `s` with dot products with the fixed vectors `v`. |
 | [`LinearCombinationInSet(s, v)`](@ref MathOptInterface.LinearCombinationInSet) | The cone of vector `(y, x)` such that ``\sum_i y_i v_i + x`` belongs to `s`. |
 
 ## Matrix cones

docs/src/reference/standard_form.md

@@ -153,6 +153,6 @@ LogDetConeTriangle
 LogDetConeSquare
 RootDetConeTriangle
 RootDetConeSquare
-SetWithDotProducts
+SetDotProducts
 LinearCombinationInSet
 ```

src/Bridges/Constraint/Constraint.jl

@@ -108,6 +108,7 @@ function add_all_bridges(bridged_model, ::Type{T}) where {T}
     MOI.Bridges.add_bridge(bridged_model, SOS1ToMILPBridge{T})
     MOI.Bridges.add_bridge(bridged_model, SOS2ToMILPBridge{T})
     MOI.Bridges.add_bridge(bridged_model, IndicatorToMILPBridge{T})
+    MOI.Bridges.add_bridge(bridged_model, LinearCombinationBridge{T})
     return
 end
 

src/Bridges/Variable/bridges/set_dot.jl

@@ -1,41 +1,47 @@
+# 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.
+
 struct DotProductsBridge{T,S,A,V} <:
-       SetMapBridge{T,S,MOI.SetWithDotProducts{S,A,V}}
+       SetMapBridge{T,S,MOI.SetDotProducts{S,A,V}}
     variables::Vector{MOI.VariableIndex}
     constraint::MOI.ConstraintIndex{MOI.VectorOfVariables,S}
-    set::MOI.SetWithDotProducts{S,A,V}
+    set::MOI.SetDotProducts{S,A,V}
 end
 
 function supports_constrained_variable(
     ::Type{<:DotProductsBridge},
-    ::Type{<:MOI.SetWithDotProducts},
+    ::Type{<:MOI.SetDotProducts},
 )
     return true
 end
 
 function concrete_bridge_type(
     ::Type{<:DotProductsBridge{T}},
-    ::Type{MOI.SetWithDotProducts{S,A,V}},
+    ::Type{MOI.SetDotProducts{S,A,V}},
 ) where {T,S,A,V}
     return DotProductsBridge{T,S,A,V}
 end
 
 function bridge_constrained_variable(
     BT::Type{DotProductsBridge{T,S,A,V}},
     model::MOI.ModelLike,
-    set::MOI.SetWithDotProducts{S,A,V},
+    set::MOI.SetDotProducts{S,A,V},
 ) where {T,S,A,V}
     variables, constraint =
         _add_constrained_var(model, MOI.Bridges.inverse_map_set(BT, set))
     return BT(variables, constraint, set)
 end
 
-function MOI.Bridges.map_set(bridge::DotProductsBridge{T,S}, set::S) where {T,S}
+function MOI.Bridges.map_set(bridge::DotProductsBridge{T,S}, ::S) where {T,S}
     return bridge.set
 end
 
 function MOI.Bridges.inverse_map_set(
     ::Type{<:DotProductsBridge},
-    set::MOI.SetWithDotProducts,
+    set::MOI.SetDotProducts,
 )
     return set.set
 end
@@ -53,6 +59,21 @@ function MOI.Bridges.map_function(
     )
 end
 
+function MOI.Bridges.map_function(
+    bridge::DotProductsBridge{T},
+    func,
+) where {T}
+    scalars = MOI.Utilities.eachscalar(func)
+    return MOI.Utilities.vectorize([
+        MOI.Utilities.set_dot(
+            vector,
+            scalars,
+            bridge.set.set,
+        )
+        for vector in bridge.set.vectors
+    ])
+end
+
 # This returns `true` by default for `SetMapBridge`
 # but is is not supported for this bridge because `inverse_map_function`
 # is not implemented

src/Bridges/Variable/set_map.jl

@@ -149,7 +149,7 @@ function MOI.get(
     bridge::SetMapBridge,
 )
     value = MOI.get(model, attr, bridge.constraint)
-    return MOI.Bridges.map_function(typeof(bridge), value)
+    return MOI.Bridges.map_function(bridge, value)
 end
 
 function MOI.get(
@@ -224,7 +224,7 @@ function unbridged_map(
     func = MOI.VectorOfVariables(vis)
     funcs = MOI.Bridges.inverse_map_function(bridge, func)
     scalars = MOI.Utilities.eachscalar(funcs)
-    # FIXME not correct for SetWithDotProducts, it won't recover the dot product variables
+    # FIXME not correct for SetDotProducts, it won't recover the dot product variables
     return Pair{MOI.VariableIndex,F}[
         bridge.variables[i] => scalars[i] for i in eachindex(bridge.variables)
     ]

src/Test/test_conic.jl

@@ -5770,6 +5770,165 @@ function setup_test(
     return
 end
 
+"""
+The goal is to find the maximum lower bound `γ` for the polynomial `x^2 - 2x`.
+Using samples `-1` and `1`, the polynomial `x^2 - 2x - γ` evaluates at `-γ`
+and `2 - γ` respectively.
+The dot product with the gram matrix is the evaluation of `[1; x] * [1 x]` hence
+`[1; -1] * [1 -1]` and `[1; 1] * [1 1]` respectively.
+
+The polynomial version is:
+max γ
+s.t. [-γ, 2 - γ] in SetDotProducts(
+    PSD(2),
+    [[1; -1] * [1 -1], [1; 1] * [1 1]],
+)
+Its dual (moment version) is:
+min -y[1] - y[2]
+s.t. [-γ, 2 - γ] in LinearCombinationInSet(
+    PSD(2),
+    [[1; -1] * [1 -1], [1; 1] * [1 1]],
+)
+"""
+function test_conic_PositiveSemidefinite_RankOne_polynomial(
+    model::MOI.ModelLike,
+    config::Config{T},
+) where {T}
+    set = MOI.SetDotProducts(
+        MOI.PositiveSemidefiniteConeTriangle(2),
+        MOI.TriangleVectorization.([
+            MOI.PositiveSemidefiniteFactorization(T[1, -1]),
+            MOI.PositiveSemidefiniteFactorization(T[1, 1]),
+        ]),
+    )
+    @requires MOI.supports_constraint(model, MOI.VectorAffineFunction{T}, typeof(set))
+    @requires MOI.supports_incremental_interface(model)
+    @requires MOI.supports(model, MOI.ObjectiveSense())
+    @requires MOI.supports(model, MOI.ObjectiveFunction{MOI.VariableIndex}())
+    γ = MOI.add_variable(model)
+    c = MOI.add_constraint(
+        model,
+        MOI.Utilities.operate(vcat, T, T(3) - T(1) * γ, T(-1) - T(1) * γ),
+        set,
+    )
+    MOI.set(model, MOI.ObjectiveSense(), MOI.MAX_SENSE)
+    MOI.set(model, MOI.ObjectiveFunction{MOI.VariableIndex}(), γ)
+    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
+        @test ≈(MOI.get(model, MOI.ObjectiveValue()), T(-1), config)
+        if _supports(config, MOI.DualObjectiveValue)
+            @test ≈(MOI.get(model, MOI.DualObjectiveValue()), T(-1), config)
+        end
+        @test ≈(MOI.get(model, MOI.VariablePrimal(), γ), T(-1), config)
+        @test ≈(MOI.get(model, MOI.ConstraintPrimal(), c), T[4, 0], config)
+        if _supports(config, MOI.ConstraintDual)
+            @test ≈(MOI.get(model, MOI.ConstraintDual(), c), T[0, 1], config)
+        end
+    end
+    return
+end
+
+function setup_test(
+    ::typeof(test_conic_PositiveSemidefinite_RankOne_polynomial),
+    model::MOIU.MockOptimizer,
+    ::Config{T},
+) where {T<:Real}
+    A = MOI.TriangleVectorization{T,MOI.PositiveSemidefiniteFactorization{T,Vector{T}}}
+    MOIU.set_mock_optimize!(
+        model,
+        (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
+            mock,
+            T[-1],
+            (MOI.VectorAffineFunction{T}, MOI.SetDotProducts{
+                MOI.PositiveSemidefiniteConeTriangle,
+                A,
+                Vector{A},
+            }) => [T[0, 1]],
+        ),
+    )
+    return
+end
+
+"""
+The moment version of `test_conic_PositiveSemidefinite_RankOne_polynomial`
+
+We look for a measure `μ = y1 * δ_{-1} + y2 * δ_{1}` where `δ_{c}` is the Dirac
+measure centered at `c`. The objective is
+`⟨μ, x^2 - 2x⟩ = y1 * ⟨δ_{-1}, x^2 - 2x⟩ + y2 * ⟨δ_{1}, x^2 - 2x⟩ = 3y1 - y2`.
+We want `μ` to be a probability measure so `1 = ⟨μ, 1⟩ = y1 + y2`.
+"""
+function test_conic_PositiveSemidefinite_RankOne_moment(
+    model::MOI.ModelLike,
+    config::Config{T},
+) where {T}
+    set = MOI.LinearCombinationInSet(
+        MOI.PositiveSemidefiniteConeTriangle(2),
+        MOI.TriangleVectorization.([
+            MOI.PositiveSemidefiniteFactorization(T[1, -1]),
+            MOI.PositiveSemidefiniteFactorization(T[1, 1]),
+        ]),
+    )
+    @requires MOI.supports_add_constrained_variables(model, typeof(set))
+    @requires MOI.supports_incremental_interface(model)
+    @requires MOI.supports(model, MOI.ObjectiveSense())
+    @requires MOI.supports(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{T}}())
+    y, cy = MOI.add_constrained_variables(
+        model,
+        set,
+    )
+    c = MOI.add_constraint(model, T(1) * y[1] + T(1) * y[2], MOI.EqualTo(T(1)))
+    MOI.set(model, MOI.ObjectiveSense(), MOI.MIN_SENSE)
+    MOI.set(model, MOI.ObjectiveFunction{MOI.ScalarAffineFunction{T}}(), T(3) * y[1] - T(1) * y[2])
+    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
+        @test ≈(MOI.get(model, MOI.ObjectiveValue()), T(-1), config)
+        if _supports(config, MOI.DualObjectiveValue)
+            @test ≈(MOI.get(model, MOI.DualObjectiveValue()), T(-1), config)
+        end
+        @test ≈(MOI.get(model, MOI.VariablePrimal(), y), T[0, 1], config)
+        @test ≈(MOI.get(model, MOI.ConstraintPrimal(), c), T(1), config)
+        if _supports(config, MOI.ConstraintDual)
+            @test ≈(MOI.get(model, MOI.ConstraintDual(), cy), T[4, 0], config)
+            @test ≈(MOI.get(model, MOI.ConstraintDual(), c), T(-1), config)
+        end
+    end
+    return
+end
+
+function setup_test(
+    ::typeof(test_conic_PositiveSemidefinite_RankOne_moment),
+    model::MOIU.MockOptimizer,
+    ::Config{T},
+) where {T<:Real}
+    A = MOI.TriangleVectorization{T,MOI.PositiveSemidefiniteFactorization{T,Vector{T}}}
+    MOIU.set_mock_optimize!(
+        model,
+        (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
+            mock,
+            T[0, 1],
+            (MOI.ScalarAffineFunction{T}, MOI.EqualTo{T}) => [T(-1)],
+            (MOI.VectorOfVariables, MOI.LinearCombinationInSet{
+                MOI.PositiveSemidefiniteConeTriangle,
+                A,
+                Vector{A},
+            }) => [T[4, 0]],
+        ),
+    )
+    return
+end
+
 """
     _test_det_cone_helper_ellipsoid(
         model::MOI.ModelLike,

src/Utilities/set_dot.jl

@@ -86,16 +86,6 @@ function set_dot(
     return x[1] * y[1] + x[2] * y[2] + triangle_dot(x, y, set.side_dimension, 2)
 end
 
-function set_dot(
-    x::AbstractVector,
-    y::AbstractVector,
-    set::Union{MOI.SetWithDotProducts,MOI.LinearCombinationInSet},
-)
-    m = length(set.matrices)
-    return LinearAlgebra.dot(view(x, 1:m), view(y, 1:m)) +
-           set_dot(view(x, (m+1):length(x)), view(y, (m+1):length(y)), set.set)
-end
-
 """
     dot_coefficients(a::AbstractVector, set::AbstractVectorSet)
 
@@ -155,16 +145,6 @@ function dot_coefficients(a::AbstractVector, set::MOI.LogDetConeTriangle)
     return b
 end
 
-function dot_coefficients(
-    a::AbstractVector,
-    set::Union{MOI.SetWithDotProducts,MOI.LinearCombinationInSet},
-)
-    b = copy(a)
-    m = length(set.vectors)
-    b[(m+1):end] = dot_coefficients(b[(m+1):end], set.set)
-    return b
-end
-
 # For `SetDotScalingVector`, we would like to compute the dot product
 # of canonical vectors in O(1) instead of O(n)
 # See https://github.com/jump-dev/Dualization.jl/pull/135

src/sets.jl

@@ -1804,31 +1804,31 @@ function Base.getproperty(
 end
 
 """
-    SetWithDotProducts(set::MOI.AbstractSet, vectors::AbstractVector)
+    SetDotProducts(set::MOI.AbstractSet, vectors::AbstractVector)
 
 Given a set `set` of dimension `d` and `m` vectors `a_1`, ..., `a_m` given in `vectors`, this is the set:
 ``\\{ ((\\langle a_1, x \\rangle, ..., \\langle a_m, x \\rangle) \\in \\mathbb{R}^{m} : x \\in \\text{set} \\}.``
 """
-struct SetWithDotProducts{S,A,V<:AbstractVector{A}} <: AbstractVectorSet
+struct SetDotProducts{S,A,V<:AbstractVector{A}} <: AbstractVectorSet
     set::S
     vectors::V
 end
 
-function Base.:(==)(s1::SetWithDotProducts, s2::SetWithDotProducts)
+function Base.:(==)(s1::SetDotProducts, s2::SetDotProducts)
     return s1.set == s2.set && s1.vectors == s2.vectors
 end
 
-function Base.copy(s::SetWithDotProducts)
-    return SetWithDotProducts(copy(s.set), copy(s.vectors))
+function Base.copy(s::SetDotProducts)
+    return SetDotProducts(copy(s.set), copy(s.vectors))
 end
 
-dimension(s::SetWithDotProducts) = length(s.vectors)
+dimension(s::SetDotProducts) = length(s.vectors)
 
-function dual_set(s::SetWithDotProducts)
+function dual_set(s::SetDotProducts)
     return LinearCombinationInSet(s.set, s.vectors)
 end
 
-function dual_set_type(::Type{SetWithDotProducts{S,A,V}}) where {S,A,V}
+function dual_set_type(::Type{SetDotProducts{S,A,V}}) where {S,A,V}
     return LinearCombinationInSet{S,A,V}
 end
 
@@ -1843,14 +1843,22 @@ struct LinearCombinationInSet{S,A,V<:AbstractVector{A}} <: AbstractVectorSet
     vectors::V
 end
 
+function Base.:(==)(s1::LinearCombinationInSet, s2::LinearCombinationInSet)
+    return s1.set == s2.set && s1.vectors == s2.vectors
+end
+
+function Base.copy(s::LinearCombinationInSet)
+    return LinearCombinationInSet(copy(s.set), copy(s.vectors))
+end
+
 dimension(s::LinearCombinationInSet) = length(s.vectors)
 
 function dual_set(s::LinearCombinationInSet)
-    return SetWithDotProducts(s.side_dimension, s.vectors)
+    return SetDotProducts(s.side_dimension, s.vectors)
 end
 
 function dual_set_type(::Type{LinearCombinationInSet{S,A,V}}) where {S,A,V}
-    return SetWithDotProducts{S,A,V}
+    return SetDotProducts{S,A,V}
 end
 
 abstract type AbstractFactorization{T,F} <: AbstractMatrix{T} end