General sets (LazySet)

LazySet{N}

Abstract type for the set types in LazySets.

Notes

LazySet types should be parameterized with a type N, typically N<:Real, for using different numeric types.

Every concrete LazySet must define the following method:

• dim(X::LazySet) – the ambient dimension of X

While not strictly required, it is useful to define the following method:

• σ(d::AbstractVector, X::LazySet) – the support vector of X in a given direction d

The method

• ρ(d::AbstractVector, X::LazySet) – the support function of X in a given direction d

is optional because there is a fallback implementation relying on σ. However, for potentially unbounded sets (which includes most lazy set types) this fallback cannot be used and an explicit method must be implemented.

The subtypes of LazySet (including abstract interfaces):

julia> subtypes(LazySet, false)
17-element Vector{Any}:
 AbstractAffineMap
 AbstractPolynomialZonotope
 Bloating
 CachedMinkowskiSumArray
 CartesianProduct
 CartesianProductArray
 Complement
 ConvexSet
 Intersection
 IntersectionArray
 MinkowskiSum
 MinkowskiSumArray
 Polygon
 QuadraticMap
 Rectification
 UnionSet
 UnionSetArray

If we only consider concrete subtypes, then:

julia> concrete_subtypes = subtypes(LazySet, true);

julia> length(concrete_subtypes)
53

julia> println.(concrete_subtypes);
AffineMap
Ball1
Ball2
BallInf
Ballp
Bloating
CachedMinkowskiSumArray
CartesianProduct
CartesianProductArray
Complement
ConvexHull
ConvexHullArray
DensePolynomialZonotope
Ellipsoid
EmptySet
ExponentialMap
ExponentialProjectionMap
HParallelotope
HPolygon
HPolygonOpt
HPolyhedron
HPolytope
HalfSpace
Hyperplane
Hyperrectangle
Intersection
IntersectionArray
Interval
InverseLinearMap
Line
Line2D
LineSegment
LinearMap
MinkowskiSum
MinkowskiSumArray
Polygon
QuadraticMap
Rectification
ResetMap
SimpleSparsePolynomialZonotope
Singleton
SparsePolynomialZonotope
Star
SymmetricIntervalHull
Tetrahedron
Translation
UnionSet
UnionSetArray
Universe
VPolygon
VPolytope
ZeroSet
Zonotope

dim(X::LazySet)

Compute the ambient dimension of a set.

Input

• X – set

Output

The ambient dimension of the set.

plot_lazyset(X::LazySet{N}, [ε]::Real=N(PLOT_PRECISION); ...) where {N}

Plot a set.

Input

• X – set
• ε – (optional, default: PLOT_PRECISION) approximation error bound

Notes

This recipe just defines the default plotting options and then calls the function plot_recipe, which then implements the set-specific plotting.

The argument ε is ignored by some set types, e.g., for polyhedra (subtypes of AbstractPolyhedron).

Examples

julia> B = Ball2(ones(2), 0.1);

julia> plot(B, 1e-3)  # default accuracy value (explicitly given for clarity here)

julia> plot(B, 1e-2)  # faster but less accurate than the previous call

plot_list(list::AbstractVector{VN}, [ε]::Real=N(PLOT_PRECISION),
          [Nφ]::Int=PLOT_POLAR_DIRECTIONS; [same_recipe]=false; ...)
    where {N, VN<:LazySet{N}}

Plot a list of sets.

Input

• list – list of sets (1D or 2D)
• ε – (optional, default: PLOT_PRECISION) approximation error bound
• Nφ – (optional, default: PLOT_POLAR_DIRECTIONS) number of polar directions (used to plot lazy intersections)
• same_recipe – (optional, default: false) switch for faster plotting but without individual plot recipes (see notes below)

Notes

For each set in the list we apply an individual plot recipe.

The option same_recipe provides access to a faster plotting scheme where all sets in the list are first converted to polytopes and then plotted in one single run. This, however, is not suitable when plotting flat sets (line segments, singletons) because then the polytope plot recipe does not deliver good results. Hence by default we do not use this option. For plotting a large number of (non-flat) polytopes, we highly advise activating this option.

Examples

julia> B1 = BallInf(zeros(2), 0.4);

julia> B2 = BallInf(ones(2), 0.4);

julia> plot([B1, B2])

Some of the sets in the list may not be plotted precisely but rather overapproximated first. The second argument ε controls the accuracy of this overapproximation.

julia> Bs = [BallInf(zeros(2), 0.4), Ball2(ones(2), 0.4)];

julia> plot(Bs, 1e-3)  # default accuracy value (explicitly given for clarity)

julia> plot(Bs, 1e-2)  # faster but less accurate than the previous call

plot_vlist(X::S, ε::Real) where {S<:LazySet}

Return a list of vertices used for plotting a two-dimensional set.

Input

• X – two-dimensional set
• ε – precision parameter

Output

A list of vertices of a polygon P. For convex X, P usually satisfies that the Hausdorff distance to X is less than ε.

plot3d(S::LazySet; [backend]=default_polyhedra_backend(S), [alpha]=1.0,
       [color]=:blue, [colormap]=:viridis, [colorrange]=nothing,
       [interpolate]=false, [linewidth]=1, [overdraw]=false, [shading]=true,
       [transparency]=true, [visible]=true)

Plot a three-dimensional set using Makie.

Input

• S – set
• backend – (optional, default: default_polyhedra_backend(S)) backend for polyhedral computations
• alpha – (optional, default: 1.0) float in [0,1]; the alpha or transparency value
• color – (optional, default: :blue) Symbol or Colorant; the color of the main plot element (markers, lines, etc.), which can be a color symbol/string like :red
• colormap – (optional, default: :viridis) the color map of the main plot; use available_gradients() to see which gradients are available, which can also be used as [:red, :black]
• colorrange – (optional, default: nothing, which falls back to Makie.Automatic()) a tuple (min, max) where min and max specify the data range to be used for indexing the colormap
• interpolate – (optional, default: false) a boolean for heatmap and images; toggles color interpolation between nearby pixels
• linewidth – (optional, default: 1) a number that specifies the width of the line in line and linesegments plots
• overdraw – (optional, default: false)
• shading – (optional, default: true) a boolean that toggles shading (for meshes)
• transparency – (optional, default: true) if true, the set is transparent, otherwise it is displayed as a solid object
• visible – (optional, default: true) a boolean that toggles visibility of the plot

For a complete list of attributes and usage see Makie's documentation.

Notes

This plot recipe works by computing the list of constraints of S and converting to a polytope in H-representation. Then, this polytope is transformed with Polyhedra.Mesh and plotted using the mesh function.

If the function constraints_list is not applicable to your set S, try overapproximation first; e.g. via

julia> Sapprox = overapproximate(S, SphericalDirections(10))

julia> using Polyhedra, GLMakie

julia> plot3d(Sapprox)

The number 10 above corresponds to the number of directions considered; for better resolution use higher values (but it will take longer).

For efficiency consider using the CDDLib backend, as in

julia> using CDDLib

julia> plot3d(Sapprox, backend=CDDLib.Library())

Examples

The functionality requires both Polyhedra and a Makie backend. After loading LazySets, do using Polyhedra, GLMakie (or another Makie backend).

julia> using LazySets, Polyhedra, GLMakie

julia> plot3d(10 * rand(Hyperrectangle, dim=3))

julia> plot3d!(10 * rand(Hyperrectangle, dim=3), color=:red)

plot3d!(S::LazySet; backend=default_polyhedra_backend(S), [alpha]=1.0,
       [color]=:blue, [colormap]=:viridis, [colorrange]=nothing,
       [interpolate]=false, [linewidth]=1, [overdraw]=false, [shading]=true,
       [transparency]=true, [visible]=true)

Plot a three-dimensional set using Makie.

Input

See plot3d for the description of the inputs. For a complete list of attributes and usage see Makie's documentation.

Notes

See the documentation of plot3d for examples.

○(c, a)

Convenience constructor of Ellipsoids or Ball2s depending on the type of a.

Input

• c – center
• a – additional parameter (either a shape matrix for Ellipsoid or a radius for Ball2)

Output

A Ellipsoids or Ball2s depending on the type of a.

Notes

"○" can be typed by \bigcirc<tab>.

Algorithm

The default implementation returns false. All set types that can determine convexity should override this behavior.

Examples

A ball in the infinity norm is always convex:

julia> isconvextype(BallInf)
true

The union (UnionSet) of two sets may in general be either convex or not. Since convexity cannot be decided by just using type information, isconvextype returns false.

julia> isconvextype(UnionSet)
false

However, the type parameters of set operations allow to decide convexity in some cases by falling back to the convexity information of the argument types. Consider the lazy intersection. The intersection of two convex sets is always convex:

julia> isconvextype(Intersection{Float64, BallInf{Float64}, BallInf{Float64}})
true

Algorithm

The default implementation applies low in each dimension.

Algorithm

The default implementation applies high in each dimension.

Algorithm

The default implementation computes the extrema via low and high.

Algorithm

The default implementation computes the extrema via low and high.

convex_hull(X::LazySet; kwargs...)

Compute the convex hull of a polytopic set.

Input

• X – polytopic set

Output

The set X itself if its type indicates that it is convex, or a new set with the list of the vertices describing the convex hull.

Algorithm

For non-convex sets, this method relies on the vertices_list method.

triangulate(X::LazySet)

Triangulate a three-dimensional polyhedral set.

Input

• X – three-dimensional polyhedral set

Output

A tuple (p, c) where p is a matrix, with each column containing a point, and c is a list of 3-tuples containing the indices of the points in each triangle.

Notes

Note that some sets may still represent an unbounded set even though their type actually does not (example: HPolytope, because the construction with non-bounding linear constraints is allowed).

Algorithm

The default implementation returns false. All set types that can determine boundedness should override this behavior.

isbounded(X::LazySet; [algorithm]="support_function")

Check whether a set is bounded.

Input

• X – set
• algorithm – (optional, default: "support_function") algorithm choice, possible options are "support_function" and "stiemke"

Output

true iff the set is bounded.

Algorithm

See the documentation of _isbounded_unit_dimensions or _isbounded_stiemke for details.

_isbounded_unit_dimensions(X::LazySet)

Check whether a set is bounded in each unit dimension.

Input

• X – set

Output

true iff the set is bounded in each unit dimension.

Algorithm

This function asks for upper and lower bounds in each ambient dimension.

Algorithm

The default implementation returns false. All set types that can determine the polyhedral property should override this behavior.

isfeasible(constraints::AbstractVector{<:HalfSpace}, [witness]::Bool=false;
           [solver]=nothing)

Check for feasibility of a list of linear constraints.

Input

• constraints – list of linear constraints
• witness – (optional; default: false) flag for witness production
• solver – (optional; default: nothing) LP solver

Output

true if the linear constraints are feasible, and false otherwise.

Algorithm

This implementation solves the corresponding feasibility linear program.

Algorithm

The default implementation handles the case p == Inf using extrema. Otherwise it checks whether X is polytopic, in which case it iterates over all vertices.

Algorithm

The default implementation handles the case p == Inf using ballinf_approximation.

Algorithm

The default implementation applies the function radius and doubles the result.

isempty(P::LazySet{N}, witness::Bool=false;
        [use_polyhedra_interface]::Bool=false, [solver]=nothing,
        [backend]=nothing) where {N}

Check whether a polyhedral set is empty.

Input

• P – polyhedral set
• witness – (optional, default: false) compute a witness if activated
• use_polyhedra_interface – (optional, default: false) if true, we use the Polyhedra interface for the emptiness test
• solver – (optional, default: nothing) LP-solver backend; uses default_lp_solver(N) if not provided
• backend – (optional, default: nothing) backend for polyhedral computations in Polyhedra; uses default_polyhedra_backend(P) if not provided

Output

• If witness option is deactivated: true iff $P = ∅$
• If witness option is activated:
  • (true, []) iff $P = ∅$
  • (false, v) iff $P ≠ ∅$ and $v ∈ P$

Notes

The default value of the backend is set internally and depends on whether the use_polyhedra_interface option is set. If the option is set, we use default_polyhedra_backend(P).

Witness production is not supported if use_polyhedra_interface is true.

Algorithm

The algorithm sets up a feasibility LP for the constraints of P. If use_polyhedra_interface is true, we call Polyhedra.isempty. Otherwise, we set up the LP internally.

linear_map(M::AbstractMatrix, P::LazySet; kwargs...)

Concrete linear map of a polyhedral set.

Input

• M – matrix
• P – polyhedral set

Output

A set representing the concrete linear map.

linear_map(a::Number, X::LazySet; kwargs...)

Alias for scale(a, X; kwargs...).

Algorithm

The default implementation applies the functions linear_map and translate.

Algorithm

The default implementation applies the functions exp and linear_map.

Algorithm

The default implementation computes a support vector along direction $[1, 0, …, 0]$. This may fail for unbounded sets.

tosimplehrep(S::LazySet)

Return the simple constraint representation $Ax ≤ b$ of a polyhedral set from its list of linear constraints.

Input

• S – polyhedral set

Output

The tuple (A, b) where A is the matrix of normal directions and b is the vector of offsets.

Algorithm

This fallback implementation relies on constraints_list(S).

reflect(P::LazySet)

Concrete reflection of a set P, resulting in the reflected set -P.

Input

• P – polyhedral set

Output

The set -P, which is either of type HPolytope if P is a polytope (i.e., bounded) or of type HPolyhedron otherwise.

Algorithm

This function requires that the list of constraints of the set P is available, i.e., that it can be written as $P = \{z ∈ ℝⁿ: ⋂ sᵢᵀz ≤ rᵢ, i = 1, ..., N\}.$

This function can be used to implement the alternative definition of the Minkowski difference

\[ A ⊖ B = \{a − b \mid a ∈ A, b ∈ B\} = A ⊕ (-B)\]

by calling minkowski_sum(A, reflect(B)).

Algorithm

The default implementation determines v ∈ interior(X) with error tolerance ε by checking whether a Ballp of norm p with center v and radius ε is contained in X.

Algorithm

The default implementation checks whether the set type of the input is an operation type using isoperationtype(::Type{<:LazySet}).

Examples

julia> B = BallInf([0.0, 0.0], 1.0);

julia> isoperation(B)
false

julia> isoperation(B ⊕ B)
true

Algorithm

The default implementation first checks X ≈ Y, which returns true if and only if X and Y have the same base type and approximately the same values. If that fails, the implementation checks the double inclusion X ⊆ Y && Y ⊆ X.

Examples

julia> X = BallInf([0.1, 0.2], 0.3);

julia> Y = convert(HPolytope, X);

julia> X == Y
false

julia> isequivalent(X, Y)
true

Algorithm

The default implementation applies area for two-dimensional sets.

area(X::LazySet)

Compute the area of a two-dimensional polytopic set.

Input

• X – two-dimensional polytopic set

Output

A number representing the area of X.

Notes

This algorithm is applicable to any polytopic set X whose list of vertices can be computed via vertices_list.

Algorithm

Let m be the number of vertices of X. We consider the following instances:

• m = 0, 1, 2: the output is zero.
• m = 3: the triangle case is solved using the Shoelace formula with 3 points.
• m = 4: the quadrilateral case is solved by the factored version of the Shoelace formula with 4 points.

Otherwise, the general Shoelace formula is used; for details see the Wikipedia page.

Algorithm

The default implementation returns X. All relevant lazy set types should override this behavior, typically by recursively calling concretize on the set arguments.

Algorithm

The default implementation assumes that X is polyhedral and returns a UnionSetArray of HalfSpaces, i.e., the output is the union of the linear constraints which are obtained by complementing each constraint of X. For any pair of sets $(X, Y)$ we have the identity $(X ∩ Y)^C = X^C ∪ Y^C$. We can apply this identity for each constraint that defines a polyhedral set.

polyhedron(P::LazySet; [backend]=default_polyhedra_backend(P))

Compute a set representation from Polyhedra.jl.

Input

• P – polyhedral set
• backend – (optional, default: call default_polyhedra_backend(P)) the polyhedral computations backend

Output

A set representation in the Polyhedra library.

Notes

For further information on the supported backends see Polyhedra's documentation.

Algorithm

This default implementation uses tosimplehrep, which computes the constraint representation of P. Set types preferring the vertex representation should implement their own method.

project(X::LazySet, block::AbstractVector{Int}, [::Nothing=nothing],
        [n]::Int=dim(X); [kwargs...])

Project a set to a given block by using a concrete linear map.

Input

• X – set
• block – block structure - a vector with the dimensions of interest
• nothing – (default: nothing) needed for dispatch
• n – (optional, default: dim(X)) ambient dimension of the set X

Output

A set representing the projection of X to block block.

Algorithm

We apply the function linear_map.

project(X::LazySet, block::AbstractVector{Int}, set_type::Type{T},
        [n]::Int=dim(X); [kwargs...]) where {T<:LazySet}

Project a set to a given block and set type, possibly involving an overapproximation.

Input

• X – set
• block – block structure - a vector with the dimensions of interest
• set_type – target set type
• n – (optional, default: dim(X)) ambient dimension of the set X

Output

A set of type set_type representing an overapproximation of the projection of X.

Algorithm

1. Project the set X with M⋅X, where M is the identity matrix in the block

coordinates and zero otherwise.

1. Overapproximate the projected set using overapproximate and set_type.

project(X::LazySet, block::AbstractVector{Int},
        set_type_and_precision::Pair{T, N}, [n]::Int=dim(X);
        [kwargs...]) where {T<:UnionAll, N<:Real}

Project a set to a given block and set type with a certified error bound.

Input

• X – set
• block – block structure - a vector with the dimensions of interest
• set_type_and_precision – pair (T, ε) of a target set type T and an error bound ε for approximation
• n – (optional, default: dim(X)) ambient dimension of the set X

Output

A set representing the epsilon-close approximation of the projection of X.

Notes

Currently we only support HPolygon as set type, which implies that the set must be two-dimensional.

Algorithm

1. Project the set X with M⋅X, where M is the identity matrix in the block

coordinates and zero otherwise.

1. Overapproximate the projected set with the given error bound ε.

project(X::LazySet, block::AbstractVector{Int}, ε::Real, [n]::Int=dim(X);
        [kwargs...])

Project a set to a given block and set type with an error bound.

Input

• X – set
• block – block structure - a vector with the dimensions of interest
• ε – error bound for approximation
• n – (optional, default: dim(X)) ambient dimension of the set X

Output

A set representing the epsilon-close approximation of the projection of X.

Algorithm

1. Project the set X with M⋅X, where M is the identity matrix in the block

coordinates and zero otherwise.

1. Overapproximate the projected set with the given error bound ε.

The target set type is chosen automatically.

Algorithm

For each dimension in which X is both positive and negative, we split X into these two parts and then project the negative part to zero.

permute(S::Singleton, p::AbstractVector{Int})

Permute the dimensions according to a permutation vector.

Input

• S – singleton
• p – permutation vector

Output

A new Singleton with the permuted dimensions.

permute(V::VPolygon, p::AbstractVector{Int})

Permute the dimensions according to a permutation vector.

Input

• P – polygon in vertex representation
• p – permutation vector

Output

The permuted polygon in vertex representation.

permute(X::LazySet, p::AbstractVector{Int})

Permute the dimensions of a set according to a given permutation vector.

Input

• X – set
• p – permutation vector

Output

A new set corresponding to X where the dimensions have been permuted according to p.

permute(U::Universe, p::AbstractVector{Int})

Permute the dimensions according to a permutation vector.

Input

• U – universe
• p – permutation vector

Output

The same universe.

permute(H::Hyperrectangle, p::AbstractVector{Int})

Permute the dimensions according to a permutation vector.

Input

• H – hyperrectangle
• p – permutation vector

Output

A permuted hyperrectangle.

rationalize(::Type{T}, X::LazySet, tol::Real) where {T<:Integer}

Approximate a set represented with floating-point numbers as a set represented with rationals of the given integer type.

Input

• T – (optional, default: Int) integer type to represent the rationals
• X – set represented by floating-point components
• tol – (optional, default: eps(N)) tolerance of the result; each rationalized component will differ by no more than tol with respect to the floating-point value

Output

A set of the same base type as X where each numerical component is of type Rational{T}.

Algorithm

The default implementation rationalizes each field.

singleton_list(P::LazySet)

Return the vertices of a polytopic set as a list of singletons.

Input

• P – polytopic set

Output

A list of the vertices of P as Singletons.

Notes

This function relies on vertices_list, which raises an error if the set is not polytopic (e.g., unbounded).

Algorithm

The default implementation computes all constraints via constraints_list.

Algorithm

The default implementation computes all vertices via vertices_list.

delaunay(X::LazySet)

Compute the Delaunay triangulation of the given polytopic set.

Input

• X – polytopic set
• compute_triangles_3d – (optional; default: false) flag to compute the 2D triangulation of a 3D set

Output

A union of polytopes in vertex representation.

Notes

This implementation requires the package MiniQhull.jl, which uses the library Qhull.

The method works in arbitrary dimension and the requirement is that the list of vertices of X can be obtained.

chebyshev_center_radius(P::LazySet{N};
                        [backend]=default_polyhedra_backend(P),
                        [solver]=default_lp_solver_polyhedra(N; presolve=true)
                       ) where {N}

Compute a Chebyshev center and the corresponding radius of a polytopic set.

Input

• P – polytopic set
• backend – (optional; default: default_polyhedra_backend(P)) the backend for polyhedral computations
• solver – (optional; default: default_lp_solver_polyhedra(N; presolve=true)) the LP solver passed to Polyhedra

Output

The pair (c, r) where c is a Chebyshev center of P and r is the radius of the largest ball with center c enclosed by P.

Notes

The Chebyshev center is the center of a largest Euclidean ball enclosed by P. In general, the center of such a ball is not unique, but the radius is.

Algorithm

We call Polyhedra.chebyshevcenter.

Algorithm

The default implementation calls scale! on a copy of X.

Algorithm

The default implementation calls translate! on a copy of X.

plot_recipe(X::LazySet, [ε])

Convert a compact convex set to a pair (x, y) of points for plotting.

Input

• X – compact convex set
• ε – approximation-error bound

Output

A pair (x, y) of points that can be plotted.

Notes

We do not support three-dimensional or higher-dimensional sets at the moment.

Algorithm

One-dimensional sets are converted to an Interval.

For two-dimensional sets, we first compute a polygonal overapproximation. The second argument, ε, corresponds to the error in Hausdorff distance between the overapproximating set and X. On the other hand, if you only want to produce a fast box-overapproximation of X, pass ε=Inf.

Finally, we use the plot recipe for the constructed set (interval or polygon).

low(X::ConvexSet{N}, i::Int) where {N}

Return the lower coordinate of a convex set in a given dimension.

Input

• X – convex set
• i – dimension of interest

Output

The lower coordinate of the set in the given dimension.

high(X::ConvexSet{N}, i::Int) where {N}

Return the higher coordinate of a convex set in a given dimension.

Input

• X – convex set
• i – dimension of interest

Output

The higher coordinate of the set in the given dimension.

an_element(U::Universe{N}) where {N}

Return some element of a universe.

Input

• U – universe

Output

The origin.

an_element(P::AbstractPolyhedron{N};
           [solver]=default_lp_solver(N)) where {N}

Return some element of a polyhedron.

Input

• P – polyhedron
• solver – (optional, default: default_lp_solver(N)) LP solver

Output

An element of the polyhedron, or an error if the polyhedron is empty.

Algorithm

An element is obtained by solving a feasibility linear program.

σ(d::AbstractVector, B::Ball1)

Return the support vector of a ball in the 1-norm in the given direction.

Input

• d – direction
• B – ball in the 1-norm

Output

The support vector in the given direction.

σ(d::AbstractVector, P::VPolygon)

Return a support vector of a polygon in a given direction.

Input

• d – direction
• P – polygon in vertex representation

Output

A support vector in the given direction. If the direction has norm zero, the first vertex is returned.

Algorithm

This implementation uses a binary search algorithm when the polygon has more than 10 vertices and a brute-force search when it has 10 or fewer vertices. The brute-force search compares the projection of each vector along the given direction and runs in $O(n)$ where $n$ is the number of vertices. The binary search runs in $O(log n)$ and we follow this implementation based on an algorithm described in [1].

[1] Joseph O'Rourke, Computational Geometry in C (2nd Edition).

σ(d::AbstractVector, L::Line2D)

Return the support vector of a 2D line in a given direction.

Input

• d – direction
• L – 2D line

Output

The support vector in the given direction, which is defined the same way as for the more general Hyperplane.

σ(d::AbstractVector, X::LazySet)

Compute a support vector of a set in a given direction.

Input

• d – direction
• X – set

Output

A support vector of X in direction d.

Notes

A convenience alias support_vector is also available.

σ(d::AbstractVector, ∅::EmptySet)

Return the support vector of an empty set.

Input

• d – direction
• ∅ – empty set

Output

An error.

σ(d::AbstractVector, B::Ball2)

Return the support vector of a 2-norm ball in the given direction.

Input

• d – direction
• B – ball in the 2-norm

Output

The support vector in the given direction. If the direction has norm zero, the center is returned.

Notes

Let $c$ and $r$ be the center and radius of a ball $B$ in the 2-norm, respectively. For nonzero direction $d$ we have

\[σ(d, B) = c + r \frac{d}{‖d‖_2}.\]

This function requires computing the 2-norm of the input direction, which is performed in the given precision of the numeric datatype of both the direction and the set. Exact inputs are not supported.

σ(d::AbstractVector, B::BallInf)

Return the support vector of a ball in the infinity norm in the given direction.

Input

• d – direction
• B – ball in the infinity norm

Output

The support vector in the given direction. If the direction has norm zero, the center of the ball is returned.

σ(d::AbstractVector, L::LineSegment)

Return the support vector of a 2D line segment in a given direction.

Input

• d – direction
• L – 2D line segment

Output

The support vector in the given direction.

Algorithm

If the angle between the vector $q - p$ and $d$ is bigger than 90° and less than 270° (measured in counter-clockwise order), the result is $p$, otherwise it is $q$. If the angle is exactly 90° or 270°, or if the direction has norm zero, this implementation returns $q$.

σ(d::AbstractVector, U::Universe)

Return the support vector of a universe.

Input

• d – direction
• U – universe

Output

A vector with infinity values, except in dimensions where the direction is zero.

σ(d::AbstractVector, hs::HalfSpace)

Return the support vector of a half-space in a given direction.

Input

• d – direction
• hs – half-space

Output

The support vector in the given direction, which is only defined in the following two cases:

1. The direction has norm zero.
2. The direction is (a multiple of) the normal direction of the half-space.

In both cases the result is any point on the boundary (the defining hyperplane). Otherwise this function throws an error.

σ(d::AbstractVector, L::Line)

Return a support vector of a line in a given direction.

Input

• d – direction
• L – line

Output

A support vector in the given direction.

σ(d::AbstractVector, E::Ellipsoid)

Return the support vector of an ellipsoid in a given direction.

Input

• d – direction
• E – ellipsoid

Output

The support vector in the given direction.

Algorithm

Let $E$ be an ellipsoid of center $c$ and shape matrix $Q = BB^\mathrm{T}$. Its support vector along direction $d$ can be deduced from that of the unit Euclidean ball $\mathcal{B}_2$ using the algebraic relations for the support vector,

\[σ_{B\mathcal{B}_2 ⊕ c}(d) = c + Bσ_{\mathcal{B}_2} (B^\mathrm{T} d) = c + \dfrac{Qd}{\sqrt{d^\mathrm{T}Q d}}.\]

σ(d::AbstractVector, P::VPolytope)

Return a support vector of a polytope in vertex representation in a given direction.

Input

• d – direction
• P – polytope in vertex representation

Output

A support vector in the given direction.

Algorithm

A support vector maximizes the support function. For a polytope, the support function is always maximized in some vertex. Hence it is sufficient to check all vertices.

σ(d::AbstractVector, H::Hyperplane)

Return a support vector of a hyperplane.

Input

• d – direction
• H – hyperplane

Output

A support vector in the given direction, which is only defined in the following two cases:

1. The direction has norm zero.
2. The direction is the hyperplane's normal direction or its opposite direction.

In all cases, any point on the hyperplane is a solution. Otherwise this function throws an error.

σ(d::AbstractVector, T::Tetrahedron)

Return a support vector of a tetrahedron in a given direction.

Input

• d – direction
• T – tetrahedron

Output

A support vector in the given direction.

Algorithm

Currently falls back to the VPolytope implementation.

σ(d::AbstractVector, Z::AbstractZonotope)

Return a support vector of a zonotopic set in a given direction.

Input

• d – direction
• Z – zonotopic set

Output

A support vector in the given direction. If the direction has norm zero, the vertex with $ξ_i = 1 \ \ ∀ i = 1,…, p$ is returned.

σ(d::AbstractVector, H::AbstractHyperrectangle)

Return a support vector of a hyperrectangular set in a given direction.

Input

• d – direction
• H – hyperrectangular set

Output

A support vector in the given direction.

If the direction vector is zero in dimension $i$, the result will have the center's coordinate in that dimension. For instance, for the two-dimensional infinity-norm ball, if the direction points to the right, the result is the vector [1, 0] in the middle of the right-hand facet.

If the direction has norm zero, the result can be any point in H. The default implementation returns the center of H.

σ(d::AbstractVector, S::AbstractSingleton)

Return the support vector of a set with a single value.

Input

• d – direction
• S – set with a single value

Output

The support vector, which is the set's vector itself, irrespective of the given direction.

σ(d::AbstractVector, am::AbstractAffineMap)

Return a support vector of an affine map.

Input

• d – direction
• am – affine map

Output

A support vector in the given direction.

σ(d::AbstractVector, B::AbstractBallp)

Return the support vector of a ball in the p-norm in a given direction.

Input

• d – direction
• B – ball in the p-norm

Output

The support vector in the given direction. If the direction has norm zero, the center of the ball is returned.

Algorithm

The support vector of the unit ball in the $p$-norm along direction $d$ is:

\[σ(d, \mathcal{B}_p^n(0, 1)) = \dfrac{\tilde{v}}{‖\tilde{v}‖_q},\]

where $\tilde{v}_i = \frac{|d_i|^q}{d_i}$ if $d_i ≠ 0$ and $\tilde{v}_i = 0$ otherwise, for all $i=1,…,n$, and $q$ is the conjugate number of $p$. By the affine transformation $x = r\tilde{x} + c$, one obtains that the support vector of $\mathcal{B}_p^n(c, r)$ is

\[σ(d, \mathcal{B}_p^n(c, r)) = \dfrac{v}{‖v‖_q},\]

where $v_i = c_i + r\frac{|d_i|^q}{d_i}$ if $d_i ≠ 0$ and $v_i = 0$ otherwise, for all $i = 1, …, n$.

σ(d::AbstractVector, P::HPolygonOpt;
  [linear_search]::Bool=(length(P.constraints) < 10))

Return a support vector of an optimized polygon in a given direction.

Input

• d – direction
• P – optimized polygon in constraint representation
• linear_search – (optional, default: see below) flag for controlling whether to perform a linear search or a binary search

Output

The support vector in the given direction. The result is always one of the vertices; in particular, if the direction has norm zero, any vertex is returned.

Algorithm

Comparison of directions is performed using polar angles; see the definition of ⪯ for two-dimensional vectors.

For polygons with 10 or more constraints we use a binary search by default.

σ(d::AbstractVector, cup::UnionSet; [algorithm]="support_vector")

Return a support vector of the union of two sets in a given direction.

Input

• d – direction
• cup – union of two sets
• algorithm – (optional, default: "supportvector"): the algorithm to compute the support vector; if "supportvector", use the support vector of each argument; if "support_function" use the support function of each argument and evaluate the support vector of only one of them

Output

A support vector in the given direction.

Algorithm

The support vector of the union of two sets $X$ and $Y$ can be obtained as the vector that maximizes the support function of either $X$ or $Y$, i.e., it is sufficient to find the $\argmax(ρ(d, X), ρ(d, Y)])$ and evaluate its support vector.

The default implementation, with option algorithm="support_vector", computes the support vector of $X$ and $Y$ and then compares the support function using a dot product.

If the support function can be computed more efficiently, the alternative implementation algorithm="support_function" can be used, which evaluates the support function of each set directly and then calls only the support vector of either $X$ or $Y$.

σ(d::AbstractVector, B::Bloating)

Return the support vector of a bloated set in a given direction.

Input

• d – direction
• B – bloated set

Output

The support vector of the bloated set in the given direction.

σ(d::AbstractVector, cp::CartesianProduct)

Return a support vector of a Cartesian product.

Input

• d – direction
• cp – Cartesian product

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped sets.

σ(d::AbstractVector, cpa::CartesianProductArray)

Compute a support vector of a Cartesian product of a finite number of sets.

Input

• d – direction
• cpa – Cartesian product of a finite number of sets

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the product sets.

σ(d::AbstractVector, ch::ConvexHull)

Return a support vector of the convex hull of two sets in a given direction.

Input

• d – direction
• ch – convex hull of two sets

Output

A support vector of the convex hull in the given direction.

σ(d::AbstractVector, cha::ConvexHullArray)

Return a support vector of a convex hull of a finite number of sets in a given direction.

Input

• d – direction
• cha – convex hull of a finite number of sets

Output

A support vector in the given direction.

σ(d::AbstractVector, em::ExponentialMap;
  [backend]=get_exponential_backend())

Return a support vector of an exponential map.

Input

• d – direction
• em – exponential map
• backend – (optional; default: get_exponential_backend()) exponentiation backend

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

Notes

If $E = \exp(M)⋅X$, where $M$ is a matrix and $X$ is a set, it follows that $σ(d, E) = \exp(M)⋅σ(\exp(M)^T d, X)$ for any direction $d$.

σ(d::AbstractVector, eprojmap::ExponentialProjectionMap;
  [backend]=get_exponential_backend())

Return a support vector of a projection of an exponential map.

Input

• d – direction
• eprojmap – projection of an exponential map
• backend – (optional; default: get_exponential_backend()) exponentiation backend

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

Notes

If $S = (L⋅M⋅R)⋅X$, where $L$ and $R$ are matrices, $M$ is a matrix exponential, and $X$ is a set, it follows that $σ(d, S) = L⋅M⋅R⋅σ(R^T⋅M^T⋅L^T⋅d, X)$ for any direction $d$.

σ(d::AbstractVector, cap::Intersection)

Return a support vector of an intersection of two sets in a given direction.

Input

• d – direction
• cap – intersection of two sets

Output

A support vector in the given direction.

Algorithm

We compute the concrete intersection, which may be expensive.

σ(d::AbstractVector, ia::IntersectionArray)

Return a support vector of an intersection of a finite number of sets in a given direction.

Input

• d – direction
• ia – intersection of a finite number of sets

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the individual sets.

Algorithm

This implementation computes the concrete intersection, which can be expensive.

σ(d::AbstractVector, lm::LinearMap)

Return a support vector of the linear map.

Input

• d – direction
• lm – linear map

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

Notes

If $L = M⋅S$, where $M$ is a matrix and $S$ is a set, it follows that $σ(d, L) = M⋅σ(M^T d, S)$ for any direction $d$.

σ(d::AbstractVector, ilm::InverseLinearMap)

Return a support vector of a inverse linear map.

Input

• d – direction
• ilm – inverse linear map

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

Notes

If $L = M^{-1}⋅X$, where $M$ is a matrix and $X$ is a set, since (M^T)^{-1}=(M^{-1})^T, it follows that $σ(d, L) = M^{-1}⋅σ((M^T)^{-1} d, X)$ for any direction $d$.

σ(d::AbstractVector, ms::MinkowskiSum)

Return a support vector of a Minkowski sum of two sets.

Input

• d – direction
• ms – Minkowski sum of two sets

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the summand sets.

Algorithm

A valid support vector in direction $d$ of the Minkowski sum of two sets $X$ and $Y$ is the sum of the support vectors of $X$ and $Y$ in direction $d$.

σ(d::AbstractVector, msa::MinkowskiSumArray)

Return a support vector of a Minkowski sum of a finite number of sets in a given direction.

Input

• d – direction
• msa – Minkowski sum of a finite number of sets

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the summand sets.

σ(d::AbstractVector, cms::CachedMinkowskiSumArray)

Return a support vector of a cached Minkowski sum in a given direction.

Input

• d – direction
• cms – cached Minkowski sum

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the summand sets.

Notes

The result is cached, i.e., any further query with the same direction runs in constant time. When sets are added to the cached Minkowski sum, the query is only performed for the new sets.

σ(d::AbstractVector, rm::ResetMap)

Return a support vector of a reset map.

Input

• d – direction
• rm – reset map

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

σ(d::AbstractVector, sih::SymmetricIntervalHull)

Return a support vector of the symmetric interval hull of a set in a given direction.

Input

• d – direction
• sih – symmetric interval hull of a set

Output

A support vector of the symmetric interval hull of a set in the given direction. If the direction has norm zero, the origin is returned.

Algorithm

For each non-zero entry in d we need to either look up the bound (if it has been computed before) or compute it, in which case we store it for future queries.

σ(d::AbstractVector, tr::Translation)

Return a support vector of a translation.

Input

• d – direction
• tr – translation of a set

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

σ(d::AbstractVector, cup::UnionSetArray; [algorithm]="support_vector")

Return a support vector of the union of a finite number of sets in a given direction.

Input

• d – direction
• cup – union of a finite number of sets
• algorithm – (optional, default: "supportvector"): the algorithm to compute the support vector; if "supportvector", use the support vector of each argument; if "support_function", use the support function of each argument and evaluate the support vector of only one of them

Output

A support vector in the given direction.

Algorithm

The support vector of the union of a finite number of sets $X₁, X₂, ...$ can be obtained as the vector that maximizes the support function, i.e., it is sufficient to find the $\argmax([ρ(d, X₂), ρ(d, X₂), ...])$ and evaluate its support vector.

The default implementation, with option algorithm="support_vector", computes the support vector of all $X₁, X₂, ...$ and then compares the support function using the dot product.

If the support function can be computed more efficiently, the alternative implementation algorithm="support_function" can be used, which evaluates the support function of each set directly and then calls only the support vector of one of the $Xᵢ$.

σ(d::AbstractVector, R::Rectification)

Return a support vector of a rectification.

Input

• d – direction
• R – rectification

Output

A support vector in the given direction. If the direction has norm zero, the result depends on the wrapped set.

σ(d::AbstractVector,
  R::Rectification{N, <:AbstractHyperrectangle{N}}) where {N}

Return a support vector of the rectification of a hyperrectangular set.

Input

• d – direction
• R – rectification of a hyperrectangular set

Output

A support vector in the given direction.

Algorithm

Let $R(·)$ be the rectification of a vector respectively a set, and let $H$ be a hyperrectangle. Then $σ_{R(H)}(d) = R(σ_{H}(d))$.

σ(d::AbstractVector, R::Rectification{N, <:CartesianProduct{N}}) where {N}

Return a support vector of the rectification of a Cartesian product of two sets.

Input

• d – direction
• R – rectification of a Cartesian product of two sets

Output

A support vector in the given direction.

Notes

Note that this implementation creates new Rectification objects that do not get preserved. Hence a second support-vector query does not benefit from the computations in the first query. For this use case another implementation should be added.

Algorithm

Rectification distributes with the Cartesian product. Let $R(·)$ be the rectification of a set. We can just query a support vector for $R(X)$ and $R(Y)$ recursively: $σ_{R(X × Y)}(d) = σ_{R(X)}(d_X) × σ_{R(Y)}(d_Y)$, where $x × y$ concatenates vectors $x$ and $y$.

σ(d::AbstractVector,
  R::Rectification{N, <:CartesianProductArray{N}}) where {N}

Return a support vector of the rectification of a Cartesian product of a finite number of sets.

Input

• d – direction
• R – rectification of a Cartesian product of a finite number of sets

Output

A support vector in the given direction.

Notes

Note that this implementation creates new Rectification objects that do not get preserved. Hence a second support-vector query does not benefit from the computations in the first query. For this use case another implementation should be added.

Algorithm

Rectification distributes with the Cartesian product. Let $R(·)$ be the rectification of a set. We can just query a support vector for each subspace recursively: $σ_{R(X_1 × ⋯ × X_m)}(d) = σ_{R(X_1)}(d_{X_1}) × ⋯ × σ_{R(X_m)}(d_{X_m})$, where $x × y$ concatenates vectors $x$ and $y$.

σ(d::AbstractVector{M}, P::HPoly{N};
  solver=default_lp_solver(M, N) where {M, N}

Return a support vector of a polyhedron in constraint representation in a given direction.

Input

• d – direction
• P – polyhedron in constraint representation
• solver – (optional, default: default_lp_solver(M, N)) the backend used to solve the linear program

Output

The support vector in the given direction. If a polytope is unbounded in the given direction, we throw an error. If a polyhedron is unbounded in the given direction, the result contains ±Inf entries.

σ(d::AbstractVector, P::HPolygon;
  [linear_search]::Bool=(length(P.constraints) < 10))

Return a support vector of a polygon in a given direction.

Input

• d – direction
• P – polygon in constraint representation
• linear_search – (optional, default: see below) flag for controlling whether to perform a linear search or a binary search

Output

The support vector in the given direction. The result is always one of the vertices; in particular, if the direction has norm zero, any vertex is returned.

Algorithm

Comparison of directions is performed using polar angles; see the definition of ⪯ for two-dimensional vectors.

For polygons with 10 or more constraints we use a binary search by default.

σ(d::AbstractVector, X::Star)

Return a support vector of a star.

Input

• d – direction
• X – star

Output

A support vector in the given direction.

Algorithm

The default implementation computes a support vector via σ.

Algorithm

The default implementation recursively compares the fields of X and Y until a mismatch is found.

Examples

julia> HalfSpace([1], 1) == HalfSpace([1.0], 1.0)
true

julia> HalfSpace([1], 1) == HalfSpace([1.00000001], 0.99999999)
false

julia> BallInf([0.0], 1.0) == Hyperrectangle([0.0], [1.0])
false

Algorithm

The default implementation recursively compares the fields of X and Y until a mismatch is found.

Examples

julia> HalfSpace([1], 1) ≈ HalfSpace([1.0], 1.0)
true

julia> HalfSpace([1], 1) ≈ HalfSpace([1.00000001], 0.99999999)
true

julia> BallInf([0.0], 1.0) ≈ Hyperrectangle([0.0], [1.0])
false

CompactSet

An alias for compact set types.

Notes

Most lazy operations are not captured by this alias because whether their result is compact depends on the argument(s).

NonCompactSet

An alias for non-compact set types.

Notes

Most lazy operations are not captured by this alias because whether their result is non-compact depends on the argument(s).