Merge pull request jump-dev#1677 from JuliaOpt/bl/vecquad

✨ Add support for vector of quadratic constraints

src/aff_expr.jl

@@ -295,18 +295,21 @@ function jump_function(model::AbstractModel, f::MOI.VectorAffineFunction)
 end
 
 """
-    _fillvaf!(terms, offset::Int, oi::Int, aff::AffExpr)
+    _fill_vaf!(terms::Vector{<:MOI.VectorAffineTerm}, offset::Int, oi::Int,
+               aff::AbstractJuMPScalar)
 
-Fills the vectors terms at indices starting at `offset+1` with the terms of `aff`.
-The output index for all terms is `oi`.
+Fills the vectors terms at indices starting at `offset+1` with the affine terms
+of `aff`. The output index for all terms is `oi`. Return the index of the last
+term added.
 """
-function _fillvaf!(terms, offset::Int, oi::Int, aff::AffExpr)
+function _fill_vaf!(terms::Vector{<:MOI.VectorAffineTerm}, offset::Int, oi::Int,
+                    aff::AbstractJuMPScalar)
     i = 1
     for (coef, var) in linear_terms(aff)
         terms[offset+i] = MOI.VectorAffineTerm(Int64(oi), MOI.ScalarAffineTerm(coef, index(var)))
         i += 1
     end
-    offset + length(linear_terms(aff))
+    return offset + length(linear_terms(aff))
 end
 
 function MOI.VectorAffineFunction(affs::Vector{AffExpr})
@@ -316,7 +319,7 @@ function MOI.VectorAffineFunction(affs::Vector{AffExpr})
     offset = 0
     for (i, aff) in enumerate(affs)
         constant[i] = aff.constant
-        offset = _fillvaf!(terms, offset, i, aff)
+        offset = _fill_vaf!(terms, offset, i, aff)
     end
     MOI.VectorAffineFunction(terms, constant)
 end

src/quad_expr.jl

@@ -198,13 +198,21 @@ function Base.convert(::Type{GenericQuadExpr{C, V}}, v::Union{Real,AbstractVaria
 end
 GenericQuadExpr{C, V}() where {C, V} = zero(GenericQuadExpr{C, V})
 
+"""
+    moi_quadratic_term(t::Tuple)
+
+Return the MOI.ScalarQuadraticTerm for the quadratic term `t`, element of the
+[`quadterms`](@ref) iterator. Note that the `JuMP.VariableRef`s are transformed
+into `MOI.VariableIndex`s hence the owner model information is lost.
+"""
+function moi_quadratic_term(t::Tuple)
+    return MOI.ScalarQuadraticTerm(t[2] == t[3] ? 2t[1] : t[1], index(t[2]),
+                                   index(t[3]))
+end
 function MOI.ScalarQuadraticFunction(q::QuadExpr)
     assert_isfinite(q)
-    qterms = MOI.ScalarQuadraticTerm{Float64}[MOI.ScalarQuadraticTerm(
-                                               t[2] == t[3] ? 2t[1] : t[1],
-                                               index(t[2]),
-                                               index(t[3]))
-                                               for t in quadterms(q)]
+    qterms = MOI.ScalarQuadraticTerm{Float64}[moi_quadratic_term(t)
+                                              for t in quadterms(q)]
     moi_aff = MOI.ScalarAffineFunction(q.aff)
     return MOI.ScalarQuadraticFunction(moi_aff.terms,
                                        qterms, moi_aff.constant)
@@ -238,6 +246,46 @@ end
 function jump_function(model::AbstractModel, aff::MOI.ScalarQuadraticFunction)
     return QuadExpr(model, aff)
 end
+function jump_function(model::AbstractModel, f::MOI.VectorQuadraticFunction)
+    return QuadExpr[QuadExpr(model, f) for f in MOIU.eachscalar(f)]
+end
+
+"""
+    _fill_vqf!(terms::Vector{<:MOI.VectorQuadraticTerm}, offset::Int, oi::Int,
+               quad::AbstractJuMPScalar)
+
+Fills the vectors terms at indices starting at `offset+1` with the quadratic
+terms of `quad`. The output index for all terms is `oi`. Return the index of the
+last term added.
+"""
+function _fill_vqf!(terms::Vector{<:MOI.VectorQuadraticTerm}, offset::Int,
+                    oi::Int, aff::AbstractJuMPScalar)
+    i = 1
+    for term in quadterms(aff)
+        terms[offset + i] = MOI.VectorQuadraticTerm(Int64(oi),
+                                                    moi_quadratic_term(term))
+        i += 1
+    end
+    return offset + length(quadterms(aff))
+end
+
+function MOI.VectorQuadraticFunction(quads::Vector{QuadExpr})
+    num_quad_terms = sum(quad -> length(quadterms(quad)), quads)
+    quad_terms = Vector{MOI.VectorQuadraticTerm{Float64}}(undef, num_quad_terms)
+    num_aff_terms = sum(quad -> length(linear_terms(quad)), quads)
+    lin_terms = Vector{MOI.VectorAffineTerm{Float64}}(undef, num_aff_terms)
+    constants = Vector{Float64}(undef, length(quads))
+    quad_offset = 0
+    aff_offset = 0
+    for (i, quad) in enumerate(quads)
+        quad_offset = _fill_vqf!(quad_terms, quad_offset, i, quad)
+        aff_offset = _fill_vaf!(lin_terms, aff_offset, i, quad)
+        constants[i] = constant(quad)
+    end
+    MOI.VectorQuadraticFunction(lin_terms, quad_terms, constants)
+end
+moi_function(a::Vector{<:GenericQuadExpr}) = MOI.VectorQuadraticFunction(a)
+
 
 # Copy a quadratic expression to a new model by converting all the
 # variables to the new model's variables

test/constraint.jl

@@ -163,25 +163,26 @@ function constraints_test(ModelType::Type{<:JuMP.AbstractModel})
     end
 
     @testset "QuadExpr constraints" begin
-        m = ModelType()
-        @variable(m, x)
-        @variable(m, y)
+        model = ModelType()
+        @variable(model, x)
+        @variable(model, y)
 
-        cref = @constraint(m, x^2 + x <= 1)
+        cref = @constraint(model, x^2 + x <= 1)
         c = JuMP.constraint_object(cref)
         @test JuMP.isequal_canonical(c.func, x^2 + x)
         @test c.set == MOI.LessThan(1.0)
 
-        cref = @constraint(m, y*x - 1.0 == 0.0)
+        cref = @constraint(model, y*x - 1.0 == 0.0)
         c = JuMP.constraint_object(cref)
         @test JuMP.isequal_canonical(c.func, x*y)
         @test c.set == MOI.EqualTo(1.0)
 
-        # TODO: VectorQuadraticFunctions
-        # cref = @constraint(m, [x^2 - 1] in MOI.SecondOrderCone(1))
-        # c = JuMP.constraint_object(cref)
-        # @test JuMP.isequal_canonical(c.func, -1 + x^2)
-        # @test c.set == MOI.SecondOrderCone(1)
+        cref = @constraint(model, [ 2x - 4x*y + 3x^2 - 1,
+                                   -3y + 2x*y - 2x^2 + 1] in SecondOrderCone())
+        c = JuMP.constraint_object(cref)
+        @test JuMP.isequal_canonical(c.func[1], -1 + 3x^2 - 4x*y + 2x)
+        @test JuMP.isequal_canonical(c.func[2],  1 - 2x^2 + 2x*y - 3y)
+        @test c.set == MOI.SecondOrderCone(2)
     end
 
     @testset "Syntax error" begin