Merge pull request #29 from JuliaOpt/fix8

Fix #8

Author: blegat

REQUIRE

@@ -1,5 +1,6 @@
 julia 0.6
-MultivariatePolynomials 0.1.0
+MultivariatePolynomials 0.1.1
+SemialgebraicSets 0.0.2
 MultivariateMoments
 JuMP
 PolyJuMP 0.1.0

src/constraint.jl

@@ -19,11 +19,8 @@ function addpolyconstraint!(m::JuMP.Model, p, s::ZeroPoly, domain::FullSpace)
     JuMP.addVectorizedConstraint(m, constraints)
 end
 
-function addpolyconstraint!(m::JuMP.Model, p, s::ZeroPoly, domain::AlgebraicSet)
-    if !isempty(domain.p)
-        warn("Equality on algebraic set has not been implemented yet, ignoring the domain")
-    end
-    addpolyconstraint!(m, p, FullSpace())
+function addpolyconstraint!(m::JuMP.Model, p, s::ZeroPoly, domain::AbstractAlgebraicSet)
+    addpolyconstraint!(m, rem(p, ideal(domain)), s, FullSpace())
 end
 
 function addpolyconstraint!(m::JuMP.Model, p, s::ZeroPoly, domain::BasicSemialgebraicSet)
@@ -42,22 +39,15 @@ function addpolyconstraint!(m::JuMP.Model, P::Matrix{PT}, ::PSDCone, domain::Abs
     addpolyconstraint!(m, p, NonNegPoly(), domain)
 end
 
-function addpolyconstraint!(m::JuMP.Model, p, ::Union{NonNegPoly, SOSCone}, domain::FullSpace)
-    # FIXME If p is a MatPolynomial, p.x will not be correct
-    Z = getmonomialsforcertificate(monomials(p))
-    slack = createpoly(m, Poly{true}(Z), :Cont)
-    q = p - slack
+function addpolyconstraint!(m::JuMP.Model, p, ::NonNegPolySubCones, domain::AbstractAlgebraicSet)
+    r = rem(p, ideal(domain))
+    X = getmonomialsforcertificate(monomials(r))
+    slack = createpoly(m, Poly{true}(X), :Cont)
+    q = r - slack
     lincons = addpolyconstraint!(m, q, ZeroPoly(), domain)
     SOSConstraint(slack, lincons, monomials(q))
 end
 
-function addpolyconstraint!(m::JuMP.Model, p, s::NonNegPolySubCones, domain::AlgebraicSet)
-    if !isempty(equalities(domain))
-        warn("Equality on algebraic set has not been implemented yet, ignoring the domain")
-    end
-    addpolyconstraint!(m, p, s, FullSpace())
-end
-
 function addpolyconstraint!(m::JuMP.Model, p, set::NonNegPolySubCones, domain::BasicSemialgebraicSet)
     mindeg, maxdeg = extdegree(p)
     for q in domain.p

test/BPT12e399.jl

@@ -0,0 +1,31 @@
+# Adapted from:
+# Example 3.99 of
+# Blekherman, G.; Parrilo, P. & Thomas, R.
+# Semidefinite Optimization and Convex Algebraic Geometry
+# Society for Industrial and Applied Mathematics, 2012
+
+@testset "[BPT12] Example 3.99 with $solver" for solver in sdp_solvers
+    @polyvar x y
+
+    m = SOSModel(solver = solver)
+
+    @variable m α
+
+    @constraint(m, x^2 + α*y <= 10, domain = @set x^2 + y^2 == 1)
+
+    @objective m Max α
+
+    status = solve(m)
+
+    @test status == :Optimal
+
+    @test getvalue(α) ≈ 6 atol=1e-6
+
+    @objective m Min α
+
+    status = solve(m)
+
+    @test status == :Optimal
+
+    @test getvalue(α) ≈ -6 atol=1e-6
+end

test/domain.jl

@@ -1,2 +1,3 @@
 include("sosdemo7.jl")
 include("sosdemo8.jl")
+include("BPT12e399.jl")

test/solvers.jl

@@ -10,7 +10,7 @@ function try_import(name::Symbol)
 end
 
 mos = try_import(:Mosek)
-csd = false && try_import(:CSDP)
+csd = try_import(:CSDP)
 scs = try_import(:SCS)
 
 isscs(solver) = contains(string(typeof(solver)),"SCSSolver")