Fixes for complex numbers (#75)

src/constraint.jl

@@ -67,9 +67,9 @@ bridges(::Type{<:MOI.AbstractFunction}, ::Type{<:MOI.AbstractSet}) = []
 """
     bridges(S::Type{<:MOI.AbstractSet})
 
-Return a list of bridges that may be needed to bridge constrained variables
-in `S` but not the bridges that may be needed by constraints added by the
-bridges.
+Return a list of bridges that may be needed to bridge variables constrained
+in `S` on creation but not the bridges that may be needed by constraints added
+by the bridges.
 """
 bridges(S::Type{<:MOI.AbstractSet}) = bridges(MOI.VectorOfVariables, S)
 
@@ -89,6 +89,11 @@ function bridges(::Type{<:MOI.AbstractVectorFunction},
 end
 
 """
+    bridgeable(c::JuMP.AbstractConstraint, S::Type{<:MOI.AbstractSet})
+
+Wrap the constraint `c` in `JuMP.BridgeableConstraint`s that may be needed to
+bridge variables constrained in `S` on creation.
+
     bridgeable(c::JuMP.AbstractConstraint, F::Type{<:MOI.AbstractFunction},
                S::Type{<:MOI.AbstractSet})
 
@@ -121,14 +126,20 @@ function bridgeable(c::JuMP.AbstractConstraint,
     return c
 end
 
+_coef_type(::Type{<:MOI.AbstractFunction}) = Float64
+_coef_type(::Type{<:MOI.Utilities.TypedLike{T}}) where {T} = T
+
 function bridgeable(c::JuMP.AbstractConstraint,
                     F::Type{<:MOI.AbstractFunction},
                     S::Type{<:MOI.AbstractSet})
     bridge_types = bridges(F, S)
     for bridge_type in bridge_types
         c = BridgeableConstraint(c, bridge_type)
         concrete_bridge_type = MOIB.Constraint.concrete_bridge_type(
-            bridge_type{Float64}, F, S)
+            bridge_type{_coef_type(F)},
+            F,
+            S,
+        )
         for (ST,) in MOIB.added_constrained_variable_types(concrete_bridge_type)
             c = bridgeable(c, ST)
         end
@@ -141,7 +152,8 @@ end
 
 ### @constraint macro ###
 
-non_constant(a::Vector{<:Number}) = convert(Vector{AffExpr}, a)
+non_constant(a::Vector{<:Real}) = convert(Vector{AffExpr}, a)
+non_constant(a::Vector{<:Complex}) = convert(Vector{GenericAffExpr{ComplexF64,VariableRef}}, a)
 non_constant(a::Vector{<:JuMP.AbstractJuMPScalar}) = a
 non_constant_coefficients(p) = non_constant(coefficients(p))
 

test/Mock/zero_polynomial.jl

@@ -1,10 +1,15 @@
-config = MOI.Test.Config()
+@testset "Model" begin
+    optimize!(mock) = MOIU.mock_optimize!(mock, [1.0, 1.0, 0.0],
+        (MOI.VectorAffineFunction{Float64}, MOI.Zeros) => [[-1.0], [-1.0]])
+    mock = bridged_mock(optimize!)
+    Tests.real_zero_polynomial_test(mock, MOI.Test.Config())
+end
 
 @testset "Model" begin
     optimize!(mock) = MOIU.mock_optimize!(mock, [1.0, 1.0, 0.0],
         (MOI.VectorAffineFunction{Float64}, MOI.Zeros) => [[-1.0], [-1.0]])
     mock = bridged_mock(optimize!)
-    Tests.zero_polynomial_test(mock, config)
+    Tests.complex_zero_polynomial_test(mock, MOI.Test.Config())
 end
 
 # The VectorOfVariables-in-ZeroPolynomialSet is bridged by VectorFunctionizeBridge
@@ -15,5 +20,5 @@ end
         (MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) => [-1.0],
         (MOI.VectorAffineFunction{Float64}, MOI.Zeros) => [[-1.0], [-1.0]])
     mock = bridged_mock(optimize!, model = NoFreeVariable{Float64}())
-    Tests.zero_polynomial_test(mock, config)
+    Tests.real_zero_polynomial_test(mock, MOI.Test.Config())
 end

test/Tests/zero_polynomial.jl

@@ -4,8 +4,11 @@ using MultivariateMoments
 using PolyJuMP
 using DynamicPolynomials
 
-function zero_polynomial_test(optimizer::MOI.AbstractOptimizer,
-                              config::MOI.Test.Config)
+function _zero_polynomial_test(
+    optimizer::MOI.AbstractOptimizer,
+    config::MOI.Test.Config,
+    complex::Bool,
+)
     atol = config.atol
     rtol = config.rtol
 
@@ -16,8 +19,10 @@ function zero_polynomial_test(optimizer::MOI.AbstractOptimizer,
     @variable(model, γ)
 
     @polyvar x y
-    cref = @constraint(model, (α - β) * x * y == 0)
-    cγ = @constraint(model, γ * x * y in PolyJuMP.ZeroPoly())
+    z = complex ? 0 + 0im : 0
+    cref = @constraint(model, (α - β) * x * y == z)
+    p = γ * x * y
+    cγ = @constraint(model, p in PolyJuMP.ZeroPoly())
 
     @objective(model, Max, α + γ)
     optimize!(model)
@@ -28,9 +33,11 @@ function zero_polynomial_test(optimizer::MOI.AbstractOptimizer,
     @test primal_status(model) == MOI.FEASIBLE_POINT
     @test value(α) ≈ 1.0 atol=atol rtol=rtol
     @test value(β) ≈ 1.0 atol=atol rtol=rtol
+    @test value(γ) ≈ 0.0 atol=atol rtol=rtol
     @test value(UpperBoundRef(β)) ≈ 1.0 atol=atol rtol=rtol
     @test value(cref) isa MultivariatePolynomials.AbstractPolynomial{Float64}
     @test value(cref) ≈ 0.0 * x * y atol=atol rtol=rtol
+    @test value(cγ) ≈ 0.0 * x * y atol=atol rtol=rtol
 
     @test dual_status(model) == MOI.FEASIBLE_POINT
     @test dual(UpperBoundRef(β)) ≈ -1.0 atol=atol rtol=rtol
@@ -52,4 +59,18 @@ function zero_polynomial_test(optimizer::MOI.AbstractOptimizer,
     end
 end
 
-linear_tests["zero_polynomial"] = zero_polynomial_test
+function real_zero_polynomial_test(
+    optimizer::MOI.AbstractOptimizer,
+    config::MOI.Test.Config,
+)
+    return _zero_polynomial_test(optimizer, config, false)
+end
+function complex_zero_polynomial_test(
+    optimizer::MOI.AbstractOptimizer,
+    config::MOI.Test.Config,
+)
+    return _zero_polynomial_test(optimizer, config, false)
+end
+
+linear_tests["zero_polynomial"] = real_zero_polynomial_test
+linear_tests["zero_polynomial"] = complex_zero_polynomial_test

test/constraint.jl

@@ -19,11 +19,15 @@ function _isequal(x::AbstractArray, y::AbstractArray)
     size(x) == size(y) && all(_isequal.(x, y))
 end
 
+_con_constant(a::Real) = convert(GenericAffExpr{Float64,VariableRef}, a)
+_con_constant(a::Complex) = convert(GenericAffExpr{ComplexF64,VariableRef}, a)
+_con_constant(a) = a
+
 # `MOI.Utilities.Model` canonicalizes the constraints so we need to
 # canonicalize them as well for the printing tests.
 function _canon(model, p::MP.APL)
     return MP.polynomial(map(MP.terms(p)) do t
-        coef = MP.coefficient(t)
+        coef = _con_constant(MP.coefficient(t))
         moi = JuMP.moi_function(coef)
         jump = JuMP.jump_function(model, MOI.Utilities.canonical(moi))
         return MP.term(jump, MP.monomial(t))
@@ -32,9 +36,9 @@ end
 _canon(model, p::Matrix) = _canon.(model, p)
 
 function _test_constraint(m, cref, S::Type, jump_set::PolyJuMP.PolynomialSet,
-                 p, ineqs, eqs, basis=MB.MonomialBasis, kwargs=[])
+                 p, ineqs, eqs, basis=MB.MonomialBasis, kwargs=[]; T=Float64)
     @test cref isa JuMP.ConstraintRef{
-        Model, <:MOI.ConstraintIndex{MOI.VectorAffineFunction{Float64},
+        Model, <:MOI.ConstraintIndex{MOI.VectorAffineFunction{T},
                                      <:S}}
     c = JuMP.constraint_object(cref)
     set = JuMP.moi_set(c)
@@ -63,7 +67,7 @@ function _test_constraint(m, cref, S::Type, jump_set::PolyJuMP.PolynomialSet,
     # == 0 is not defined either
     # c.p and p can be matrices
     @test _isequal(JuMP.reshape_vector(JuMP.jump_function(c),
-                                       JuMP.shape(c)), p)
+                                       JuMP.shape(c)), p_canon)
     if isempty(ineqs)
         if isempty(eqs)
             @test isa(set.domain, FullSpace)
@@ -132,6 +136,10 @@ function test_NonNeg(var)
     dom = @set x^2 + y^2 == 1 && x^3 + x*y^2 + y >= 1
     S = DummyPolyModule.NonNeg
     jump_set = DummyPolyModule.DummyNonNeg()
+    _test_constraint(m, @constraint(m, x >= y),
+            S, jump_set, x - y, [], [])
+    _test_constraint(m, @constraint(m, im * x >= y),
+            S, jump_set, im * x - y, [], [], T=ComplexF64)
     _test_constraint(m, @constraint(m, p >= q + 1, domain = @set y >= 1 && dom),
             S, jump_set, p - q - 1, [y - 1, x^3 + x*y^2 + y - 1],
             [x^2 + y^2 - 1])
@@ -164,12 +172,18 @@ function test_ZeroPolynomialSet(var)
         S = PolyJuMP.ZeroPolynomialSet{FullSpace}
         _test_constraint(m, @constraint(m, p == q),
                 S, jump_set, p - q, [], [])
+        @test PolyJuMP.ZeroPolynomialBridge in m.bridge_types
         _test_constraint(m, @constraint(m, p - q in PolyJuMP.ZeroPoly()),
                 S, jump_set, p - q, [], [])
+        _test_constraint(m, @constraint(m, x == y),
+                S, jump_set, x - y, [], [])
+        _test_constraint(m, @constraint(m, x == im * y),
+                S, jump_set, x - im * y, [], [], T=ComplexF64)
     end
     S = PolyJuMP.ZeroPolynomialSet
     _test_constraint(m, @constraint(m, p == q, domain = @set x == 1 && f(x, y)),
             S, jump_set, p - q, [], [x + y - 2, x - 1])
+    @test PolyJuMP.ZeroPolynomialInAlgebraicSetBridge in m.bridge_types
     _test_constraint(m, @constraint(m, p - q in PolyJuMP.ZeroPoly(), domain = @set x == 1 && f(x, y)),
             S, jump_set, p - q, [], [x + y - 2, x - 1])
     S = PolyJuMP.PlusMinusSet

test/testpolymodule.jl

@@ -25,6 +25,31 @@ Base.copy(set::NonNeg) = set
 struct DummyNonNeg <: PolyJuMP.PolynomialSet end
 
 JuMP.reshape_set(::NonNeg, ::PolyJuMP.PolynomialShape) = DummyNonNeg()
+
+struct DummyNonNegBridge{T,F} <: MOIB.Constraint.AbstractBridge end
+function PolyJuMP.bridges(
+    ::Type{<:MOI.AbstractVectorFunction},
+    ::Type{<:NonNeg},
+)
+    return [DummyNonNegBridge]
+end
+function MOI.Bridges.added_constrained_variable_types(::Type{<:DummyNonNegBridge})
+    return Tuple{Type}[]
+end
+function MOI.Bridges.added_constraint_types(::Type{DummyNonNegBridge{T,F}}) where {T, F}
+    return Tuple{Type,Type}[(F, MOI.EqualTo{T})]
+end
+function MOI.Bridges.Constraint.concrete_bridge_type(
+    ::Type{<:DummyNonNegBridge{T}},
+    F::Type{<:MOI.AbstractVectorFunction},
+    ::Type{<:NonNeg}) where {T}
+    # This tests that `T` matches the coefficient type of `F`
+    # as otherwise it would error.
+    G = MOI.Utilities.promote_operation(-, T, F, MOI.VectorOfVariables)
+    H = MOI.Utilities.scalar_type(G)
+    return DummyNonNegBridge{T,H}
+end
+
 function JuMP.moi_set(cone::DummyNonNeg,
                       monos::AbstractVector{<:AbstractMonomial};
                       domain::AbstractSemialgebraicSet=FullSpace(),
@@ -38,7 +63,11 @@ function JuMP.build_constraint(_error::Function, p::AbstractPolynomialLike,
     coefs = PolyJuMP.non_constant_coefficients(p)
     monos = monomials(p)
     set = JuMP.moi_set(s, monos; kwargs...)
-    return JuMP.VectorConstraint(coefs, set, PolyJuMP.PolynomialShape(monos))
+    return PolyJuMP.bridgeable(
+        JuMP.VectorConstraint(coefs, set, PolyJuMP.PolynomialShape(monos)),
+        JuMP.moi_function_type(typeof(coefs)),
+        typeof(set),
+    )
 end
 
 struct MatrixPolynomialShape{MT <: AbstractMonomial,