Skip to content

Refactor some constraint into GenericConstraint #590

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 16 commits into from
Apr 18, 2024
86 changes: 33 additions & 53 deletions src/constraints/GreaterThanConstraint.jl
Original file line number Diff line number Diff line change
@@ -1,69 +1,49 @@
mutable struct GreaterThanConstraint <: Constraint
lhs::AbstractExpr
rhs::AbstractExpr
size::Tuple{Int,Int}
dual::Union{Value,Nothing}

function GreaterThanConstraint(lhs::AbstractExpr, rhs::AbstractExpr)
if sign(lhs) == ComplexSign() || sign(rhs) == ComplexSign()
error(
"Cannot create inequality constraint between expressions of sign $(sign(lhs)) and $(sign(rhs))",
)
end
if lhs.size == rhs.size || lhs.size == (1, 1)
sz = rhs.size
if lhs.size == (1, 1) && rhs.size != (1, 1)
lhs = lhs * ones(rhs.size)
end
elseif rhs.size == (1, 1)
sz = lhs.size
if rhs.size == (1, 1) && lhs.size != (1, 1)
rhs = rhs * ones(lhs.size)
end
else
error(
"Cannot create inequality constraint between expressions of size $(lhs.size) and $(rhs.size)",
)
end
return new(lhs, rhs, sz, nothing)
end
function set_with_size(::Type{MOI.Nonnegatives}, sz::Tuple{Int,Int})
return MOI.Nonnegatives(prod(sz))
end

head(io::IO, ::GreaterThanConstraint) = print(io, "≥")
head(io::IO, ::MOI.Nonnegatives) = print(io, "≥")

function is_feasible(f, ::MOI.Nonnegatives, tol)
return all(f .>= tol)
end

function vexity(c::GreaterThanConstraint)
vex = -vexity(c.lhs) + (vexity(c.rhs))
if vex == ConcaveVexity()
function vexity(vex, ::MOI.Nonnegatives)
if vex == ConvexVexity()
return NotDcp()
end
return vex
end

function _add_constraint!(
context::Context{T},
c::GreaterThanConstraint,
) where {T}
f = conic_form!(context, c.lhs - c.rhs)
if f isa AbstractVector
if !all(f .>= -CONSTANT_CONSTRAINT_TOL[])
@warn "Constant constraint is violated"
context.detected_infeasible_during_formulation[] = true
function _promote_size(lhs::AbstractExpr, rhs::AbstractExpr)
if sign(lhs) == ComplexSign() || sign(rhs) == ComplexSign()
error(
"Cannot create inequality constraint between expressions of sign $(sign(lhs)) and $(sign(rhs))",
)
end
if lhs.size == rhs.size || lhs.size == (1, 1)
sz = rhs.size
if lhs.size == (1, 1) && rhs.size != (1, 1)
lhs = lhs * ones(rhs.size)
end
elseif rhs.size == (1, 1)
sz = lhs.size
if rhs.size == (1, 1) && lhs.size != (1, 1)
rhs = rhs * ones(lhs.size)
end
return
else
error(
"Cannot create inequality constraint between expressions of size $(lhs.size) and $(rhs.size)",
)
end
set = MOI.Nonnegatives(MOI.output_dimension(f))
context.constr_to_moi_inds[c] = MOI_add_constraint(context.model, f, set)
return
return lhs, rhs
end

Base.:>=(lhs::AbstractExpr, rhs::AbstractExpr) = GreaterThanConstraint(lhs, rhs)
function Base.:>=(lhs::AbstractExpr, rhs::AbstractExpr)
lhs, rhs = _promote_size(lhs, rhs)
return GenericConstraint{MOI.Nonnegatives}(lhs - rhs)
end

Base.:>=(lhs::AbstractExpr, rhs::Value) = >=(lhs, constant(rhs))

Base.:>=(lhs::Value, rhs::AbstractExpr) = >=(constant(lhs), rhs)

function populate_dual!(model::MOI.ModelLike, c::GreaterThanConstraint, indices)
ret = MOI.get(model, MOI.ConstraintDual(), indices)
c.dual = output(reshape(ret, c.size))
return
end
64 changes: 11 additions & 53 deletions src/constraints/LessThanConstraint.jl
Original file line number Diff line number Diff line change
@@ -1,67 +1,25 @@
mutable struct LessThanConstraint <: Constraint
lhs::AbstractExpr
rhs::AbstractExpr
size::Tuple{Int,Int}
dual::Union{Value,Nothing}

function LessThanConstraint(lhs::AbstractExpr, rhs::AbstractExpr)
if sign(lhs) == ComplexSign() || sign(rhs) == ComplexSign()
error(
"Cannot create inequality constraint between expressions of sign $(sign(lhs)) and $(sign(rhs))",
)
end
if lhs.size == rhs.size || lhs.size == (1, 1)
sz = rhs.size
if lhs.size == (1, 1) && rhs.size != (1, 1)
lhs = lhs * ones(rhs.size)
end
elseif rhs.size == (1, 1)
sz = lhs.size
if rhs.size == (1, 1) && lhs.size != (1, 1)
rhs = rhs * ones(lhs.size)
end
else
error(
"Cannot create inequality constraint between expressions of size $(lhs.size) and $(rhs.size)",
)
end
return new(lhs, rhs, sz, nothing)
end
function set_with_size(::Type{MOI.Nonpositives}, sz::Tuple{Int,Int})
return MOI.Nonpositives(prod(sz))
end

head(io::IO, ::LessThanConstraint) = print(io, "≤")
head(io::IO, ::MOI.Nonpositives) = print(io, "≤")

function is_feasible(f, ::MOI.Nonpositives, tol)
return all(f .<= -tol)
end

function vexity(c::LessThanConstraint)
vex = vexity(c.lhs) + (-vexity(c.rhs))
function vexity(vex, ::MOI.Nonpositives)
if vex == ConcaveVexity()
return NotDcp()
end
return vex
end

function _add_constraint!(context::Context{T}, lt::LessThanConstraint) where {T}
f = conic_form!(context, lt.lhs - lt.rhs)
if f isa AbstractVector
# a trivial constraint without variables like `5 <= 0`
if !all(f .<= CONSTANT_CONSTRAINT_TOL[])
@warn "Constant constraint is violated"
context.detected_infeasible_during_formulation[] = true
end
return
end
set = MOI.Nonpositives(MOI.output_dimension(f))
context.constr_to_moi_inds[lt] = MOI_add_constraint(context.model, f, set)
return
function Base.:<=(lhs::AbstractExpr, rhs::AbstractExpr)
lhs, rhs = _promote_size(lhs, rhs)
return GenericConstraint{MOI.Nonpositives}(lhs - rhs)
end

Base.:<=(lhs::AbstractExpr, rhs::AbstractExpr) = LessThanConstraint(lhs, rhs)

Base.:<=(lhs::AbstractExpr, rhs::Value) = <=(lhs, constant(rhs))

Base.:<=(lhs::Value, rhs::AbstractExpr) = <=(constant(lhs), rhs)

function populate_dual!(model::MOI.ModelLike, c::LessThanConstraint, indices)
ret = MOI.get(model, MOI.ConstraintDual(), indices)
c.dual = output(reshape(ret, c.size))
return
end
66 changes: 20 additions & 46 deletions src/constraints/PositiveSemidefiniteConeConstraint.jl
Original file line number Diff line number Diff line change
@@ -1,66 +1,40 @@
mutable struct PositiveSemidefiniteConeConstraint <: Constraint
child::AbstractExpr
size::Tuple{Int,Int}
dual::Union{Value,Nothing}

function PositiveSemidefiniteConeConstraint(child::AbstractExpr)
if child.size[1] != child.size[2]
error("Positive semidefinite expressions must be square")
end
return new(child, child.size, nothing)
function set_with_size(
::Type{MOI.PositiveSemidefiniteConeSquare},
sz::Tuple{Int,Int},
)
if sz[1] != sz[2]
error("Positive semidefinite expressions must be square")
end
return MOI.PositiveSemidefiniteConeSquare(sz[1])
end

head(io::IO, ::PositiveSemidefiniteConeConstraint) = print(io, "sdp")
head(io::IO, ::MOI.PositiveSemidefiniteConeSquare) = print(io, "sdp")

AbstractTrees.children(c::PositiveSemidefiniteConeConstraint) = (c.child,)

function vexity(c::PositiveSemidefiniteConeConstraint)
if !(vexity(c.child) in (AffineVexity(), ConstVexity()))
function vexity(vex, ::MOI.PositiveSemidefiniteConeSquare)
if !(vex in (AffineVexity(), ConstVexity()))
return NotDcp()
end
return AffineVexity()
end

function _add_constraint!(
context::Context,
c::PositiveSemidefiniteConeConstraint,
)
if vexity(c.child) == ConstVexity()
x = evaluate(c.child)
tol = CONSTANT_CONSTRAINT_TOL[]
if !(x ≈ transpose(x))
@warn "constant SDP constraint is violated"
context.detected_infeasible_during_formulation[] = true
elseif evaluate(LinearAlgebra.eigmin(c.child)) < -tol
@warn "constant SDP constraint is violated"
context.detected_infeasible_during_formulation[] = true
end
return
function is_feasible(x, ::MOI.PositiveSemidefiniteConeSquare, tol)
if !(x ≈ transpose(x))
@warn "constant SDP constraint is violated"
return false
elseif LinearAlgebra.eigmin(x) < -tol
@warn "constant SDP constraint is violated"
return false
end
f = conic_form!(context, c.child)
set = MOI.PositiveSemidefiniteConeSquare(c.size[1])
context.constr_to_moi_inds[c] = MOI_add_constraint(context.model, f, set)
return
end

function populate_dual!(
model::MOI.ModelLike,
c::PositiveSemidefiniteConeConstraint,
indices,
)
dual = MOI.get(model, MOI.ConstraintDual(), indices)
c.dual = output(reshape(dual, c.size))
return
return true
end

function LinearAlgebra.isposdef(x::AbstractExpr)
if iscomplex(x)
return PositiveSemidefiniteConeConstraint(
return GenericConstraint{MOI.PositiveSemidefiniteConeSquare}(
[real(x) -imag(x); imag(x) real(x)],
)
end
return PositiveSemidefiniteConeConstraint(x)
return GenericConstraint{MOI.PositiveSemidefiniteConeSquare}(x)
end

⪰(x::AbstractExpr, y::AbstractExpr) = isposdef(x - y)
Expand Down
41 changes: 40 additions & 1 deletion src/expressions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -26,11 +26,50 @@
#
#############################################################################

const Value = Union{Number,AbstractArray}

abstract type AbstractExpr end

abstract type Constraint end

const Value = Union{Number,AbstractArray}
mutable struct GenericConstraint{S<:MOI.AbstractSet} <: Constraint
child::AbstractExpr
set::S
dual::Union{Value,Nothing}
function GenericConstraint(child, set::MOI.AbstractSet)
return new{typeof(set)}(child, set)
end
function GenericConstraint{S}(child) where {S<:MOI.AbstractSet}
return GenericConstraint(child, set_with_size(S, size(child)))
end
end

head(io::IO, c::GenericConstraint) = head(io, c.set)

AbstractTrees.children(c::GenericConstraint) = (c.child,)

function vexity(c::GenericConstraint)
return vexity(vexity(c.child), c.set)
end

function _add_constraint!(context::Context, c::GenericConstraint)
if vexity(c.child) == ConstVexity()
x = evaluate(c.child)
if !is_feasible(x, c.set, CONSTANT_CONSTRAINT_TOL[])
context.detected_infeasible_during_formulation[] = true
end
return
end
f = conic_form!(context, c.child)
context.constr_to_moi_inds[c] = MOI_add_constraint(context.model, f, c.set)
return
end

function populate_dual!(model::MOI.ModelLike, c::GenericConstraint, indices)
ret = MOI.get(model, MOI.ConstraintDual(), indices)
c.dual = output(reshape(ret, c.child.size))
return
end

# We commandeer `==` to create a constraint.
# Therefore we define `isequal` to still have a notion of equality
Expand Down
24 changes: 15 additions & 9 deletions test/test_constraints.jl
Original file line number Diff line number Diff line change
Expand Up @@ -183,30 +183,36 @@ function test_LessThanConstraint_dual_maximize()
return
end

### constraints/PositiveSemidefiniteConeConstraint
### constraints/GenericConstraint{MOI.PositiveSemidefiniteConeSquare}

function test_PositiveSemidefiniteConeConstraint()
function test_GenericConstraint_PositiveSemidefiniteConeSquare()
@test_throws(
ErrorException("Positive semidefinite expressions must be square"),
Convex.PositiveSemidefiniteConeConstraint(Variable(2, 3)),
Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}(
Variable(2, 3),
),
)
X = Variable(2, 2)
c = Convex.PositiveSemidefiniteConeConstraint(X)
c = Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}(X)
p = minimize(tr(X), [c, X >= [1 2; 3 4]])
solve!(p, SCS.Optimizer; silent_solver = true)
@test isapprox(X.value, [2.25 3; 3 4]; atol = 1e-3)
y = (c.dual + c.dual') / 2
@test isapprox(y[1], 1; atol = 1e-3)
@test (0 ⪯ X) isa Convex.PositiveSemidefiniteConeConstraint
@test (-X ⪯ 0) isa Convex.PositiveSemidefiniteConeConstraint
@test (-X ⪯ constant(0)) isa Convex.PositiveSemidefiniteConeConstraint
@test (constant(0) ⪯ X) isa Convex.PositiveSemidefiniteConeConstraint
@test (0 ⪯ X) isa
Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}
@test (-X ⪯ 0) isa
Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}
@test (-X ⪯ constant(0)) isa
Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}
@test (constant(0) ⪯ X) isa
Convex.GenericConstraint{MOI.PositiveSemidefiniteConeSquare}
@test_throws(ErrorException("Set PSD not understood"), X in :PSD)
@test vexity(X ⪯ square(Variable())) == Convex.NotDcp()
return
end

function test_PositiveSemidefiniteConeConstraint_violated()
function test_GenericConstraint_PositiveSemidefiniteConeSquare_violated()
X = constant([1 2; 3 4])
p = satisfy([X ⪰ 0])
@test_logs (:warn,) (:warn,) solve!(p, SCS.Optimizer)
Expand Down
Loading