| Documentation | Build Status |
|---|---|
 |
|
 |
 |
JuMP extension for Set Programming : optimization with set variables and inclusion/containment constraints. This package allows the formulation of a mathematical program involving set variables and inclusion/membership constraints in addition to classical variables and constraints supported by JuMP.
- STABLE — most recently tagged version of the documentation.
- LATEST — in-development version of the documentation.
The variables can either be
- a
Polytope; - an
Ellipsoid, or a piecewise semi-ellipsoid; - a
Polyset, that is the 1-sublevel set of a polynomial of degree2d.
@variable model S Polytope(piecewise=p) # polytope defined over the pieces defined by `p`
@variable model S Ellipsoid()
@variable model S Ellipsoid(piecewise=p) # piecewise semi-ellipsoid defined over the pieces defined by `p`
@variable model S PolySet(d) # 1-sublevel set of a polynomial of degree 2d
@variable model S PolySet(d, convex=true) # Convex 1-sublevel set of a polynomial of degree 2d
@variable model S PolySet(d, symmetric=true) # 1-sublevel set of a polynomial of degree 2d symmetric around the origin
@variable model S PolySet(d, symmetric=true, point=SetProg.CenterPoint([1, 0])) # 1-sublevel set of a polynomial of degree 2d symmetric around the [1, 0]Polynomial sets require an explicit using SumOfSquares. SumOfSquares and the
polynomial packages are optional dependencies, loaded through SetProg's package
extension. Calling PolySet before loading SumOfSquares raises an error with
this instruction. Non-symmetric and piecewise ellipsoids, and polynomial gauge
printing, also require the extension. Symmetric ellipsoids and polytopes can be
modeled without it.
using SetProg, SumOfSquares
model = Model()
@variable(model, S, PolySet(degree = 4, symmetric = true, dimension = 2))The following operations are allowed:
| Operation | Description |
|---|---|
| A*S | Linear mapping |
But more operations are planned to be added:
| Operation | Description |
|---|---|
| S + x | Translation of S by x |
| S1 + S2 | Minkowski sum |
| S1 ∩ S2 | Intersection of S1 and S2 |
| S1 ∪ S2 | Union of S1 and S2 |
| polar(S) | Polar of S |
The following constraints are implemented
| Operation | Description |
|---|---|
| x ∈ S | x is contained in S |
| S1 ⊆ S2 | S1 is included in S2 |
| S1 ⊇ S2 | S1 is included in S2 |
Consider a polytope
using Polyhedra
diamond = HalfSpace([1, 1], 1) ∩ HalfSpace([-1, -1], 1) ∩ HalfSpace([1, -1], 1) ∩ HalfSpace([-1, 1], 1)
simplex = HalfSpace([1, 1], 1) ∩ HalfSpace([-1, 0], 0) ∩ HalfSpace([0, -1], 0)Pick an SDP solver (see here for a list)
using Hypatia # Optimizer
optimizer_constructor = Hypatia.OptimizerTo compute the maximal symmetric ellipsoid contained in the polytope diamond defined above (i.e. Löwner-John ellipsoid):
using SetProg
model = Model(optimizer_constructor)
@variable(model, S, Ellipsoid(symmetric=true))
@constraint(model, S ⊆ diamond)
@objective(model, Max, nth_root(volume(S)))
optimize!(model)We specify in the example that the ellipsoid is symmetric around the origin to simplify the computation as the solver does not need to look for the center so the SDP problem that need to be solved has a smaller size.
We can visualize the result with Plots as follows:
using Plots
plot(polyhedron(diamond), ratio=1)
plot!(value(S))To compute the maximal ellipsoid contained in simplex, we don't need to specify
the center but at least a point in the interior of the ellipsoid. The SDP
formulation used will then determine the center and shape of the ellipsoid
simultaneously in the same SDP. For the interior point, we take the chebyshev
center of the simplex (which can be found by solving an LP). This the center of
the sphere of maximal volume in the simplex so one might rightly guess that is is
in the interior of the maximal ellispoid contained in the simplex.
using SetProg, SumOfSquares
cheby_center, cheby_radius = chebyshevcenter(simplex, optimizer_constructor)
interior_point = SetProg.InteriorPoint(cheby_center)
model = Model(optimizer_constructor)
@variable(model, S, Ellipsoid(point=interior_point))
@constraint(model, S ⊆ simplex)
@objective(model, Max, nth_root(volume(S)))
optimize!(model)We now visualize the result:
using Plots
plot(polyhedron(simplex), ratio=1)
plot!(value(S))To compute the maximal invariant set contained in a polytope (not yet implemented):
using SetProg
model = Model(optimizer_constructor)
@variable(model, S, Polytope())
@constraint(model, S ⊆ diamond)
@constraint(model, A*S ⊆ S) # Invariance constraint
@objective(model, Max, volume(S))
optimize!(model)To compute the maximal invariant ellipsoid contained in the polytope diamond defined above:
using SetProg
model = Model(optimizer_constructor)
@variable(model, S, Ellipsoid(symmetric=true))
@constraint(model, S ⊆ diamond)
@constraint(model, A*S ⊆ S) # Invariance constraint
@objective(model, Max, nth_root(volume(S)))
optimize!(model)To compute the maximal algebraic-invariant ellipsoid (i.e. AS ⊆ ES) contained in the polytope diamond defined above:
using SetProg
model = Model(optimizer_constructor)
@variable(model, S, Ellipsoid(symmetric=true)))
@constraint(model, S ⊆ diamond)
@constraint(model, A*S ⊆ E*S) # Invariance constraint
@objective(model, Max, L1_heuristic(volume(S), ones(Polyhedra.fulldim(P))))
optimize!(model)