# Set Interfaces

This section of the manual describes the interfaces for different set types. Every set that fits the description of an interface should also implement it. This helps in several ways:

- avoid code duplicates,
- provide functions for many sets at once,
- allow changes in the source code without changing the API.

The interface functions are outlined in the interface documentation. For implementations of the interfaces see the corresponding sub-pages linked in the respective sections.

The naming convention is such that all interface names (with the exception of the main abstract type `LazySet`

) should be preceded by `Abstract`

.

The following diagram shows the interface hierarchy.

- Set Interfaces
- General sets (LazySet)
- Convex sets (ConvexSet)
- Centrally symmetric sets (AbstractCentrallySymmetric)
- Polyhedra (AbstractPolyhedron)
- Polytopes (AbstractPolytope)
- Polygons (AbstractPolygon)
- Polygons in constraint representation (AbstractHPolygon)
- Centrally symmetric polytopes (AbstractCentrallySymmetricPolytope)
- Zonotopes (AbstractZonotope)
- Hyperrectangles (AbstractHyperrectangle)
- Singletons (AbstractSingleton)
- Affine maps (AbstractAffineMap)
- Star sets (AbstractStar)
- Polynomial zonotope sets (AbstractPolynomialZonotope)

## General sets (LazySet)

Every set in this library is a subtype of the abstract type `LazySet`

.

`LazySets.LazySet`

— Type`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(S::LazySet)`

– the ambient dimension of`S`

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

`σ(d::AbstractVector, S::LazySet)`

– the support vector of`S`

in a given direction`d`

The method

`ρ(d::AbstractVector, S::LazySet)`

– the support function of`S`

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)
18-element Vector{Any}:
AbstractAffineMap
AbstractPolynomialZonotope
Bloating
CachedMinkowskiSumArray
CartesianProduct
CartesianProductArray
Complement
ConvexSet
Intersection
IntersectionArray
LazySets.AbstractStar
MinkowskiSum
MinkowskiSumArray
Polygon
QuadraticMap
Rectification
UnionSet
UnionSetArray
```

If we only consider *concrete* subtypes, then:

```
julia> concrete_subtypes = subtypes(LazySet, true);
julia> length(concrete_subtypes)
54
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
LazySets.AbstractStar
Line
Line2D
LineSegment
LinearMap
MinkowskiSum
MinkowskiSumArray
Polygon
QuadraticMap
Rectification
ResetMap
RotatedHyperrectangle
SimpleSparsePolynomialZonotope
Singleton
SparsePolynomialZonotope
Star
SymmetricIntervalHull
Translation
UnionSet
UnionSetArray
Universe
VPolygon
VPolytope
ZeroSet
Zonotope
```

### Plotting

Plotting via the `Plots`

package is available for one- or two-dimensional sets. The default algorithm is to plot an outer approximation using the support function (1D) respectively the support vector (2D). This means that (1) plotting will fail if these functionalities are not available (e.g., for lazy `Intersection`

s) and (2) that plots of non-convex sets can be misleading. The implementation below internally relies on the function `plot_recipe`

. For some set types (e.g., `Intersection`

), the default implementation is overridden.

`RecipesBase.apply_recipe`

— Method`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
```

`RecipesBase.apply_recipe`

— Method```
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
```

`LazySets.plot_vlist`

— Method`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 `ε`

.

For three-dimensional sets, we support `Makie`

:

`LazySets.plot3d`

— Function```
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)
```

`LazySets.plot3d!`

— Function```
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.

### Globally defined set functions

`LazySets.isconvextype`

— Method`isconvextype(X::Type{<:LazySet})`

Check whether the given `LazySet`

type is convex.

**Input**

`X`

– subtype of`LazySet`

**Output**

`true`

if the given set type is guaranteed to be convex by using only type information, and `false`

otherwise.

**Notes**

Since this operation only acts on types (not on values), it can return false negatives, i.e., there may be instances where the set is convex, even though the answer of this function is `false`

. The examples below illustrate this point.

**Examples**

A ball in the infinity norm is always convex, hence we get:

```
julia> isconvextype(BallInf)
true
```

For instance, 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 type of its arguments. Consider for instance the lazy intersection. The intersection of two convex sets is always convex, hence we get:

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

`LazySets.low`

— Method`low(X::LazySet)`

Return a vector with the lowest coordinates of the set in each canonical direction.

**Input**

`X`

– set

**Output**

A vector with the lower coordinate of the set in each dimension.

**Notes**

See also `low(X::LazySet, i::Int)`

.

The result is the lowermost corner of the box approximation, so it is not necessarily contained in `X`

.

`LazySets.high`

— Method`high(X::LazySet)`

Return a vector with the highest coordinate of the set in each canonical direction.

**Input**

`X`

– set

**Output**

A vector with the highest coordinate of the set in each dimension.

**Notes**

See also `high(X::LazySet, i::Int)`

.

The result is the uppermost corner of the box approximation, so it is not necessarily contained in `X`

.

`Base.extrema`

— Method`extrema(X::LazySet, i::Int)`

Return the lower and higher coordinate of a set in a given dimension.

**Input**

`X`

– set`i`

– dimension of interest

**Output**

The lower and higher coordinate of the set in the given dimension.

**Notes**

The result is equivalent to `(low(X, i), high(X, i))`

, but sometimes it can be computed more efficiently.

**Algorithm**

The bounds are computed with `low`

and `high`

.

`Base.extrema`

— Method`extrema(X::LazySet)`

Return two vectors with the lowest and highest coordinate of `X`

in each canonical direction.

**Input**

`X`

– set

**Output**

Two vectors with the lowest and highest coordinates of `X`

in each dimension.

**Notes**

See also `extrema(X::LazySet, i::Int)`

.

The result is equivalent to `(low(X), high(X))`

, but sometimes it can be computed more efficiently.

The resulting points are the lowermost and uppermost corners of the box approximation, so they are not necessarily contained in `X`

.

**Algorithm**

The bounds are computed with `low`

and `high`

by default.

`LazySets.convex_hull`

— Method`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.

`LazySets.triangulate`

— 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.

`LazySets.basetype`

— Function`basetype(T::Type{<:LazySet})`

Return the base type of the given set type (i.e., without type parameters).

**Input**

`T`

– set type

**Output**

The base type of `T`

.

`basetype(S::LazySet)`

Return the base type of the given set (i.e., without type parameters).

**Input**

`S`

– set

**Output**

The base type of `S`

.

**Examples**

```
julia> Z = rand(Zonotope);
julia> basetype(Z)
Zonotope
julia> basetype(Z + Z)
MinkowskiSum
julia> basetype(LinearMap(rand(2, 2), Z + Z))
LinearMap
```

`LazySets.isboundedtype`

— Method`isboundedtype(T::Type{<:LazySet})`

Check whether a set type only represents bounded sets.

**Input**

`T`

– set type

**Output**

`true`

if the set type only represents bounded sets. 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).

**Notes**

By default this function returns `false`

. All set types that can determine boundedness should override this behavior.

`LazySets.isbounded`

— Method`isbounded(S::LazySet)`

Check whether a set is bounded.

**Input**

`S`

– 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.

`LazySets._isbounded_unit_dimensions`

— Method`_isbounded_unit_dimensions(S::LazySet)`

Check whether a set is bounded in each unit dimension.

**Input**

`S`

– 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.

`LazySets.is_polyhedral`

— Method`is_polyhedral(S::LazySet)`

Trait for polyhedral sets.

**Input**

`S`

– set

**Output**

`true`

only if the set behaves like an `AbstractPolyhedron`

.

**Notes**

The answer is conservative, i.e., may sometimes be `false`

even if the set is polyhedral.

`LinearAlgebra.norm`

— Function`norm(S::LazySet, [p]::Real=Inf)`

Return the norm of a set. It is the norm of the enclosing ball (of the given $p$-norm) of minimal volume that is centered in the origin.

**Input**

`S`

– set`p`

– (optional, default:`Inf`

) norm

**Output**

A real number representing the norm.

`IntervalArithmetic.radius`

— Function`radius(S::LazySet, [p]::Real=Inf)`

Return the radius of a set. It is the radius of the enclosing ball (of the given $p$-norm) of minimal volume with the same center.

**Input**

`S`

– set`p`

– (optional, default:`Inf`

) norm

**Output**

A real number representing the radius.

`LazySets.diameter`

— Function`diameter(S::LazySet, [p]::Real=Inf)`

Return the diameter of a set. It is the maximum distance between any two elements of the set, or, equivalently, the diameter of the enclosing ball (of the given $p$-norm) of minimal volume with the same center.

**Input**

`S`

– set`p`

– (optional, default:`Inf`

) norm

**Output**

A real number representing the diameter.

`Base.isempty`

— Method```
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 or not. 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.

`LazySets.affine_map`

— Method`affine_map(M::AbstractMatrix, X::LazySet, v::AbstractVector; kwargs...)`

Compute the concrete affine map $M·X + v$.

**Input**

`M`

– linear map`X`

– set`v`

– translation vector

**Output**

A set representing the affine map $M·X + v$.

**Algorithm**

The implementation applies the functions `linear_map`

and `translate`

.

`LazySets.exponential_map`

— Method`exponential_map(M::AbstractMatrix, X::LazySet)`

Compute the concrete exponential map of `M`

and `X`

, i.e., `exp(M) * X`

.

**Input**

`M`

– matrix`X`

– set

**Output**

A set representing the exponential map of `M`

and `X`

.

**Algorithm**

The implementation applies the functions `exp`

and `linear_map`

.

`LazySets.an_element`

— Method`an_element(S::LazySet)`

Return some element of a set.

**Input**

`S`

– set

**Output**

An element of a set.

**Algorithm**

An element of the set is obtained by evaluating its support vector along direction $[1, 0, …, 0]$. This may fail for unbounded sets.

`LazySets.tosimplehrep`

— Method`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)`

.

`LazySets.reflect`

— Method`reflect(P::LazySet)`

Concrete reflection of a set `P`

, resulting in the reflected set `-P`

.

**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 | a ∈ A, b ∈ B\} = A ⊕ (-B)\]

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

.

`LazySets.is_interior_point`

— Method```
is_interior_point(d::AbstractVector{N}, X::LazySet{N};
p=N(Inf), ε=_rtol(N)) where {N}
```

Check whether the point `d`

is contained in the interior of the set `X`

.

**Input**

`d`

– point`X`

– set`p`

– (optional; default:`N(Inf)`

) norm of the ball used to apply the error tolerance`ε`

– (optional; default:`_rtol(N)`

) error tolerance of check

**Output**

Boolean which indicates if the point `d`

is contained in `X`

.

**Algorithm**

The implementation checks if a `Ballp`

of norm `p`

with center `d`

and radius `ε`

is contained in the set `X`

. This is a numerical check for `d ∈ interior(X)`

with error tolerance `ε`

.

`LazySets.isoperationtype`

— Method`isoperationtype(X::Type{<:LazySet})`

Check whether the given set type is an operation or not.

**Input**

`X`

– set type

**Output**

`true`

if the given set type is a set operation and `false`

otherwise.

**Notes**

This fallback implementation returns an error that `isoperationtype`

is not implemented. Subtypes of `LazySet`

should dispatch on this function as required.

See also `isoperation(X<:LazySet)`

.

**Examples**

```
julia> isoperationtype(BallInf)
false
julia> isoperationtype(LinearMap)
true
```

`LazySets.isoperation`

— Method`isoperation(X::LazySet)`

Check whether a set is an instance of a set operation or not.

**Input**

`X`

– set

**Output**

`true`

if `X`

is an instance of a set-based operation and `false`

otherwise.

**Notes**

This fallback 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
```

`LazySets.isequivalent`

— Method`isequivalent(X::LazySet, Y::LazySet)`

Check whether two sets are equal in the mathematical sense, i.e., equivalent.

**Input**

`X`

– set`Y`

– set

**Output**

`true`

iff `X`

is equivalent to `Y`

(up to some precision).

**Algorithm**

First we check `X ≈ Y`

, which returns `true`

if and only if `X`

and `Y`

have the same type and approximately the same values (checked with `LazySets._isapprox`

). If that fails, we check 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
```

`LazySets.surface`

— Method`surface(X::LazySet)`

Compute the surface area of a set.

**Input**

`X`

– set

**Output**

A real number representing the surface area of `X`

.

`LazySets.area`

— Method`area(X::LazySet{N}) where {N}`

Compute the area of a two-dimensional polytopic set using the Shoelace formula.

**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.

`LazySets.concretize`

— Method`concretize(X::LazySet)`

Construct a concrete representation of a (possibly lazy) set.

**Input**

`X`

– set

**Output**

A concrete representation of `X`

(as far as possible).

**Notes**

Since not every lazy set has a concrete set representation in this library, the result may be partially lazy.

`LazySets.complement`

— Method`complement(X::LazySet)`

Return the complement of a polyhedral set.

**Input**

`X`

– polyhedral set

**Output**

A `UnionSetArray`

of half-spaces, i.e., the output is the union of the linear constraints which are obtained by complementing each constraint of `X`

.

**Algorithm**

The principle used in this implementation is that for any pair of sets $(X, Y)$ we have that $(X ∩ Y)^C = X^C ∪ Y^C$. In particular, we can apply this rule for each constraint that defines a polyhedral set. Hence the concrete complement can be represented as the set union of the complement of each constraint.

`Polyhedra.polyhedron`

— Method`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.

`LazySets.project`

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

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

**Input**

`S`

– set`block`

– block structure - a vector with the dimensions of interest`nothing`

– (default:`nothing`

)`n`

– (optional, default:`dim(S)`

) ambient dimension of the set`S`

**Output**

A set representing the projection of the set `S`

to block `block`

.

**Algorithm**

We apply the function `linear_map`

.

`LazySets.project`

— Method```
project(S::LazySet, block::AbstractVector{Int}, set_type::Type{TS},
[n]::Int=dim(S); [kwargs...]) where {TS<:LazySet}
```

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

**Input**

`S`

– set`block`

– block structure - a vector with the dimensions of interest`set_type`

– target set type`n`

– (optional, default:`dim(S)`

) ambient dimension of the set`S`

**Output**

A set of type `set_type`

representing an overapproximation of the projection of `S`

.

**Algorithm**

- Project the set
`S`

with`M⋅S`

, where`M`

is the identity matrix in the block

coordinates and zero otherwise.

- Overapproximate the projected set using
`overapproximate`

and`set_type`

.

`LazySets.project`

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

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

**Input**

`S`

– 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(S)`

) ambient dimension of the set`S`

**Output**

A set representing the epsilon-close approximation of the projection of `S`

.

**Notes**

Currently we only support `HPolygon`

as set type, which implies that the set must be two-dimensional.

**Algorithm**

- Project the set
`S`

with`M⋅S`

, where`M`

is the identity matrix in the block

coordinates and zero otherwise.

- Overapproximate the projected set with the given error bound
`ε`

.

`LazySets.project`

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

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

**Input**

`S`

– set`block`

– block structure - a vector with the dimensions of interest`ε`

– error bound for approximation`n`

– (optional, default:`dim(S)`

) ambient dimension of the set`S`

**Output**

A set representing the epsilon-close approximation of the projection of `S`

.

**Algorithm**

- Project the set
`S`

with`M⋅S`

, where`M`

is the identity matrix in the block

coordinates and zero otherwise.

- Overapproximate the projected set with the given error bound
`ε`

.

The target set type is chosen automatically.

`ReachabilityBase.Arrays.rectify`

— Function`rectify(X::LazySet, [concrete_intersection]::Bool=false)`

Concrete rectification of a set.

**Input**

`X`

– set`concrete_intersection`

– (optional, default:`false`

) flag to compute concrete intersections for intermediate results

**Output**

A set corresponding to the rectification of `X`

, which is in general a union of linear maps of intersections.

**Algorithm**

For each dimension in which `X`

is both positive and negative, we split `X`

into these two parts. Additionally we project the negative part to zero.

`SparseArrays.permute`

— Function`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`

.

`Base.rationalize`

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

Approximate a set of floating-point numbers as a set whose entries are rationals of the given integer type.

**Input**

`T`

– (optional, default:`Int`

) integer type to represent the rationals`X`

– set which has 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 of `X`

where each numerical component is of type `Rational{T}`

.

`LazySets.singleton_list`

— Method`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 `Singleton`

s.

**Notes**

This function relies on `vertices_list`

, which raises an error if the set is not polytopic (e.g., unbounded).

`LazySets.constraints`

— Method`constraints(X::LazySet)`

Construct an iterator over the constraints of a polyhedral set.

**Input**

`X`

– polyhedral set

**Output**

An iterator over the constraints of `X`

.

`LazySets.vertices`

— Method`vertices(X::LazySet)`

Construct an iterator over the vertices of a polytopic set.

**Input**

`X`

– polytopic set

**Output**

An iterator over the vertices of `X`

.

`MiniQhull.delaunay`

— Function`delaunay(X::LazySet)`

Compute the Delaunay triangulation of the given polytopic set.

**Input**

`X`

– polytopic set

**Output**

A union of polytopes in vertex representation.

**Notes**

This function requires that you have properly installed the package MiniQhull.jl, including the library Qhull.

The method works in arbitrary dimension and the requirement is that the list of vertices of `X`

can be obtained.

`LazySets.chebyshev_center_radius`

— Method```
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`

.

`LazySets.plot_recipe`

— Method`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).

The following methods are also defined for `LazySet`

but cannot be documented due to a bug in the documentation package.

`LazySets.low`

— Method`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.

`LazySets.high`

— Method`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.

`LazySets.an_element`

— Method```
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.

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

Return some element of a universe.

**Input**

`U`

– universe

**Output**

The origin.

### Support function and support vector

Every `LazySet`

type must define a function `σ`

to compute the support vector. The support function, `ρ`

, can optionally be defined; otherwise, a fallback definition based on `σ`

is used.

`LazySets.σ`

— Function`σ`

Function to compute the support vector σ.

`LazySets.support_vector`

— Function`support_vector`

Alias for the support vector σ.

`LazySets.ρ`

— Method`ρ(d::AbstractVector, S::LazySet)`

Evaluate the support function of a set in a given direction.

**Input**

`d`

– direction`S`

– set

**Output**

The evaluation of the support function of the set `S`

for the direction `d`

.

`LazySets.support_function`

— Function`support_function`

Alias for the support function ρ.

### Set functions that override Base functions

`Base.:==`

— Method`==(X::LazySet, Y::LazySet)`

Check whether two sets use exactly the same set representation.

**Input**

`X`

– set`Y`

– set

**Output**

`true`

iff`X`

is equal to`Y`

.

**Notes**

The check is purely syntactic and the sets need to have the same base type. For instance, `X::VPolytope == Y::HPolytope`

returns `false`

even if `X`

and `Y`

represent the same polytope. However `X::HPolytope{Int64} == Y::HPolytope{Float64}`

is a valid comparison.

**Algorithm**

We recursively compare the fields of `X`

and `Y`

until a mismatch is found.

**Examples**

```
julia> HalfSpace([1], 1) == HalfSpace([1], 1)
true
julia> HalfSpace([1], 1) == HalfSpace([1.0], 1.0)
true
julia> Ball1([0.0], 1.0) == Ball2([0.0], 1.0)
false
```

`Base.:≈`

— Method`≈(X::LazySet, Y::LazySet)`

Check whether two sets of the same type are approximately equal.

**Input**

`X`

– set`Y`

– set of the same base type as`X`

**Output**

`true`

iff`X`

is equal to`Y`

.

**Notes**

The check is purely syntactic and the sets need to have the same base type. For instance, `X::VPolytope ≈ Y::HPolytope`

returns `false`

even if `X`

and `Y`

represent the same polytope. However `X::HPolytope{Int64} ≈ Y::HPolytope{Float64}`

is a valid comparison.

**Algorithm**

We recursively compare the fields of `X`

and `Y`

until a mismatch is found.

**Examples**

```
julia> HalfSpace([1], 1) ≈ HalfSpace([1], 1)
true
julia> HalfSpace([1], 1) ≈ HalfSpace([1.00000001], 0.99999999)
true
julia> Ball1([0.0], 1.0) ≈ Ball2([0.0], 1.0)
false
```

`Base.copy`

— Method`copy(S::LazySet)`

Return a copy of a set by copying its values recursively.

**Input**

`S`

– set

**Output**

A copy of `S`

.

**Notes**

This function computes a `copy`

of each field in `S`

. See the documentation of `?copy`

for further details.

`Base.eltype`

— Function`eltype(::Type{<:LazySet{N}}) where {N}`

Return the numeric type (`N`

) of the given set type.

**Input**

`T`

– set type

**Output**

The numeric type of `T`

.

`eltype(::LazySet{N}) where {N}`

Return the numeric type (`N`

) of the given set.

**Input**

`X`

– set

**Output**

The numeric type of `X`

.

### Aliases for set types

`LazySets.CompactSet`

— Type`CompactSet`

An alias for compact set types.

**Notes**

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

`LazySets.NonCompactSet`

— Type`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 or not depends on the argument(s).

### Implementations

Concrete set representations:

Lazy set operations:

- Affine map (AffineMap)
- Linear map (LinearMap)
- Exponential map (ExponentialMap)
- Exponential projection map (ExponentialProjectionMap)
- Reset map (ResetMap)
- Translation
- Bloating
- Binary Cartesian product (CartesianProduct)
- $n$-ary Cartesian product (CartesianProductArray)
- Binary convex hull (ConvexHull)
- $n$-ary convex hull (ConvexHullArray)
- Binary intersection
- $n$-ary intersection (IntersectionArray)
- Binary Minkowski sum (MinkowskiSum)
- $n$-ary Minkowski sum (MinkowskiSumArray)
- $n$-ary Minkowski sum with cache (CachedMinkowskiSumArray)
- Binary set union (UnionSet)
- $n$-ary set union (UnionSetArray)
- Complement
- Rectification

## Convex sets (ConvexSet)

Every convex set in this library implements this interface.

`LazySets.ConvexSet`

— Type`ConvexSet{N} <: LazySet{N}`

Abstract type for convex sets, i.e., sets characterized by a (possibly infinite) intersection of halfspaces, or equivalently, sets $S$ such that for any two elements $x, y ∈ S$ and $0 ≤ λ ≤ 1$ it holds that $λ·x + (1-λ)·y ∈ S$.

## Centrally symmetric sets (AbstractCentrallySymmetric)

Centrally symmetric sets such as balls of different norms are characterized by a center. Note that there is a special interface combination Centrally symmetric polytope.

`LazySets.AbstractCentrallySymmetric`

— Type`AbstractCentrallySymmetric{N} <: ConvexSet{N}`

Abstract type for centrally symmetric compact convex sets.

**Notes**

Every concrete `AbstractCentrallySymmetric`

must define the following functions:

`center(::AbstractCentrallySymmetric)`

– return the center point`center(::AbstractCentrallySymmetric, i::Int)`

– return the center point at index`i`

The subtypes of `AbstractCentrallySymmetric`

:

```
julia> subtypes(AbstractCentrallySymmetric)
3-element Vector{Any}:
Ball2
Ballp
Ellipsoid
```

This interface defines the following functions:

`LazySets.dim`

— Method`dim(S::AbstractCentrallySymmetric)`

Return the ambient dimension of a centrally symmetric set.

**Input**

`S`

– centrally symmetric set

**Output**

The ambient dimension of the set.

`LazySets.isbounded`

— Method`isbounded(S::AbstractCentrallySymmetric)`

Check whether a centrally symmetric set is bounded.

**Input**

`S`

– centrally symmetric set

**Output**

`true`

(since a set with a unique center must be bounded).

`LazySets.isuniversal`

— Method```
isuniversal(S::AbstractCentrallySymmetric{N},
[witness]::Bool=false) where {N}
```

Check whether a centrally symmetric set is universal.

**Input**

`S`

– centrally symmetric set`witness`

– (optional, default:`false`

) compute a witness if activated

**Output**

- If
`witness`

option is deactivated:`false`

- If
`witness`

option is activated:`(false, v)`

where $v ∉ S$

**Algorithm**

Centrally symmetric sets are bounded. A witness is obtained by computing the support vector in direction `d = [1, 0, …, 0]`

and adding `d`

on top.

`LazySets.an_element`

— Method`an_element(S::AbstractCentrallySymmetric)`

Return some element of a centrally symmetric set.

**Input**

`S`

– centrally symmetric set

**Output**

The center of the centrally symmetric set.

`Base.isempty`

— Method`isempty(S::AbstractCentrallySymmetric)`

Check whether a centrally symmetric set is empty.

**Input**

`S`

– centrally symmetric set

**Output**

`false`

.

`LazySets.center`

— Method`center(H::AbstractCentrallySymmetric, i::Int)`

Return the center of a centrally symmetric set along a given dimension.

**Input**

`S`

– centrally symmetric set`i`

– dimension of interest

**Output**

The center along the given dimension.

`Base.extrema`

— Method`extrema(S::AbstractCentrallySymmetric, i::Int)`

Return the lower and higher coordinate of a centrally symmetric set in a given dimension.

**Input**

`S`

– centrally symmetric set`i`

– dimension of interest

**Output**

The lower and higher coordinate of the centrally symmetric set in the given dimension.

**Notes**

The result is equivalent to `(low(S, i), high(S, i))`

.

**Algorithm**

We compute `high(S, i)`

and then compute the lowest coordinates with the help of `center(S, i)`

(which is assumed to be cheaper to obtain).

`Base.extrema`

— Method`extrema(S::AbstractCentrallySymmetric)`

Return two vectors with the lowest and highest coordinate of a centrally symmetric set.

**Input**

`S`

– centrally symmetric set

**Output**

Two vectors with the lowest and highest coordinates of `S`

.

**Notes**

The result is equivalent to `(low(S), high(S))`

.

**Algorithm**

We compute `high(S)`

and then compute the lowest coordinates with the help of `center(S)`

(which is assumed to be cheaper to obtain).

### Implementations

## Polyhedra (AbstractPolyhedron)

A polyhedron has finitely many facets (*H-representation*) and is not necessarily bounded.

`LazySets.AbstractPolyhedron`

— Type`AbstractPolyhedron{N} <: ConvexSet{N}`

Abstract type for compact convex polyhedral sets.

**Notes**

Every concrete `AbstractPolyhedron`

must define the following functions:

`constraints_list(::AbstractPolyhedron{N})`

– return a list of all facet constraints

Polyhedra are defined as the intersection of a finite number of closed half-spaces. As such, polyhedra are closed and convex but not necessarily bounded. Bounded polyhedra are called *polytopes* (see `AbstractPolytope`

).

The subtypes of `AbstractPolyhedron`

(including abstract interfaces):

```
julia> subtypes(AbstractPolyhedron)
8-element Vector{Any}:
AbstractPolytope
HPolyhedron
HalfSpace
Hyperplane
Line
Line2D
Star
Universe
```

This interface defines the following functions:

`Base.:∈`

— Method`∈(x::AbstractVector, P::AbstractPolyhedron)`

Check whether a given point is contained in a polyhedron.

**Input**

`x`

– point/vector`P`

– polyhedron

**Output**

`true`

iff $x ∈ P$.

**Algorithm**

This implementation checks if the point lies inside each defining half-space.

`LazySets.isuniversal`

— Method`isuniversal(P::AbstractPolyhedron{N}, [witness]::Bool=false) where {N}`

Check whether a polyhedron is universal.

**Input**

`P`

– polyhedron`witness`

– (optional, default:`false`

) compute a witness if activated

**Output**

- If
`witness`

option is deactivated:`true`

iff $P$ is universal - If
`witness`

option is activated:`(true, [])`

iff $P$ is universal`(false, v)`

iff $P$ is not universal and $v ∉ P$

**Algorithm**

`P`

is universal iff it has no constraints.

A witness is produced using `isuniversal(H)`

where `H`

is the first linear constraint of `P`

.

`LazySets.constrained_dimensions`

— Method`constrained_dimensions(P::AbstractPolyhedron)`

Return the indices in which a polyhedron is constrained.

**Input**

`P`

– polyhedron

**Output**

A vector of ascending indices `i`

such that the polyhedron is constrained in dimension `i`

.

**Examples**

A 2D polyhedron with constraint $x1 ≥ 0$ is constrained in dimension 1 only.

`LazySets.linear_map`

— Method```
linear_map(M::AbstractMatrix{NM},
P::AbstractPolyhedron{NP};
[algorithm]::Union{String, Nothing}=nothing,
[check_invertibility]::Bool=true,
[cond_tol]::Number=DEFAULT_COND_TOL,
[inverse]::Union{AbstractMatrix{N}, Nothing}=nothing,
[backend]=nothing,
[elimination_method]=nothing) where {NM, NP}
```

Concrete linear map of a polyhedral set.

**Input**

`M`

– matrix`P`

– polyhedral set`algorithm`

– (optional; default:`nothing`

) algorithm to be used; for the description see the Algorithm section below; possible choices are:`"inverse"`

, alias:`"inv"`

`"inverse_right"`

, alias:`"inv_right"`

`"elimination"`

, alias:`"elim"`

`"lift"`

`"vrep"`

`"vrep_chull"`

`check_invertibility`

– (optional, default:`true`

) if`true`

check whether the given matrix`M`

is invertible; set to`false`

only if you know that`M`

is invertible`cond_tol`

– (optional; default:`DEFAULT_COND_TOL`

) tolerance of matrix condition (used to check whether the matrix is invertible)`inverse`

– (optional; default:`nothing`

) matrix inverse`M⁻¹`

; use this option if you have already computed the inverse matrix of`M`

`backend`

– (optional: default:`nothing`

) polyhedra backend`elimination_method`

– (optional: default:`nothing`

) elimination method for the`"elimination"`

algorithm

**Output**

The type of the result is "as close as possible" to the the type of `P`

. Let `(m, n)`

be the size of `M`

, where `m ≠ n`

is allowed for rectangular maps.

To fix the type of the output to something different than the default value, consider post-processing the result of this function with a call to a suitable `convert`

method.

In particular, the output depends on the type of `P`

, on `m`

, and the algorithm that was used:

If the vertex-based approach was used:

- If
`P`

is a`VPolygon`

and`m = 2`

then the output is a`VPolygon`

. - If
`P`

is a`VPolytope`

then the output is a`VPolytope`

. - Otherwise the output is an
`Interval`

if`m = 1`

, a`VPolygon`

if`m = 2`

, and a`VPolytope`

in other cases.

- If
If the invertibility criterion was used:

- The types of
`HalfSpace`

,`Hyperplane`

,`Line2D`

, and subtypes of`AbstractHPolygon`

are preserved. - If
`P`

is an`AbstractPolytope`

, then the output is an`Interval`

if`m = 1`

, an`HPolygon`

if`m = 2`

, and an`HPolytope`

in other cases. - Otherwise the output is an
`HPolyhedron`

.

- The types of

**Notes**

Since the different linear-map algorithms work at the level of constraints, this method uses dispatch on two stages: once the algorithm has been defined, first the helper methods `_linear_map_hrep_helper`

(resp. `_linear_map_vrep`

) are invoked, which dispatch on the set type. Then, each helper method calls the concrete implementation of `_linear_map_hrep`

, which dispatches on the algorithm, and returns a list of constraints.

To simplify working with different algorithms and options, the types `<: AbstractLinearMapAlgorithm`

are used. These types are singleton type or types that carry only the key data for the given algorithm, such as the matrix inverse or the polyhedra backend.

New subtypes of the `AbstractPolyhedron`

interface may define their own helper methods `_linear_map_vrep`

(respectively `_linear_map_hrep_helper`

) for special handling of the constraints returned by the implementations of `_linear_map_hrep`

; otherwise the fallback implementation for `AbstractPolyhedron`

is used, which instantiates an `HPolyhedron`

.

**Algorithm**

This method mainly implements several approaches for the linear map: inverse, right inverse, transformation to vertex representation, variable elimination, and variable lifting. Depending on the properties of `M`

and `P`

, one algorithm may be preferable over the other. Details on the algorithms are given in the following subsections.

Otherwise, if the algorithm argument is not specified, a default option is chosen based on heuristics on the types and values of `M`

and `P`

:

- If the
`"inverse"`

algorithm applies, it is used. - Otherwise, if the
`"inverse_right"`

algorithm applies, it is used. - Otherwise, if the
`"lift"`

algorithm applies, it is used. - Otherwise, the
`"elimination"`

algorithm is used.

Note that the algorithms `"inverse"`

and `"inverse_right"`

do not require the external library `Polyhedra`

. However, the fallback method `"elimination"`

requires `Polyhedra`

as well as the library `CDDLib`

.

The optional keyword arguments `inverse`

and `check_invertibility`

modify the default behavior:

- If an inverse matrix is passed in
`inverse`

, the given algorithm is applied, and if none is given, either`"inverse"`

or`"inverse_right"`

is applied (in that order of preference). - If
`check_invertibility`

is set to`false`

, the given algorithm is applied, and if none is given, either`"inverse"`

or`"inverse_right"`

is applied (in that order of preference).

**Inverse**

This algorithm is invoked with the keyword argument `algorithm="inverse"`

(or `algorithm="inv"`

). The algorithm requires that `M`

is invertible, square, and dense. If you know a priori that `M`

is invertible, set the flag `check_invertibility=false`

, such that no extra checks are done. Otherwise, we check the sufficient condition that the condition number of `M`

is not too high. The threshold for the condition number can be modified from its default value, `DEFAULT_COND_TOL`

, by passing a custom `cond_tol`

.

The algorithm is described next. Assuming that the matrix $M$ is invertible (which we check via a sufficient condition,), $y = M x$ implies $x = \text{inv}(M) y$ and we can transform the polyhedron $A x ≤ b$ to the polyhedron $A \text{inv}(M) y ≤ b$.

If the dense condition on `M`

is not fullfilled, there are two suggested workarounds: either transform to a dense matrix, i.e., calling `linear_map`

with `Matrix(M)`

, or use the `"inverse_right"`

algorithm, which does not compute the inverse matrix explicitly, but uses a polyalgorithm; see the documentation of `?`

for details.

**Inverse-right**

This algorithm is invoked with the keyword argument `algorithm="inverse_right"`

(or `algorithm="inv_right"`

). This algorithm applies to square and invertible matrices `M`

. The idea is essentially the same as for the `"inverse"`

algorithm; the difference is that in `"inverse"`

the full matrix inverse is computed, and in `"inverse_right"`

only the left division on the normal vectors is used. In particular, `"inverse_right"`

is good as a workaround when `M`

is sparse (since the `inv`

function is not available for sparse matrices).

**Elimination**

This algorithm is invoked with the keyword argument `algorithm = "elimination"`

(or `algorithm = "elim"`

). The algorithm applies to any matrix `M`

(invertible or not), and any polyhedron `P`

(bounded or not).

The idea is described next. If `P : Ax <= b`

and `y = Mx`

denote the polyhedron and the linear map, respectively, we consider the vector `z = [y, x]`

, write the given equalities and the inequalities, and then eliminate the last x variables (there are `length(x)`

in total) using a call to `Polyhedra.eliminate`

to a backend library that can do variable elimination (typically `CDDLib`

with the `BlockElimination()`

algorithm). In this way we have eliminated the "old" variables `x`

and kept the "new" or transformed variables "y".

The default elimination method is block elimination. For possible options we refer to the documentation of Polyhedra, projection/elimination.

**Lift**

This algorithm is invoked with the keyword argument `algorithm="lift"`

. The algorithm applies if `M`

is rectangular of size `m × n`

with `m > n`

and full rank (i.e., of rank `n`

).

The idea is to embed the polyhedron into the `m`

-dimensional space by appending zeros, i.e. extending all constraints of `P`

to `m`

dimensions, and constraining the last `m - n`

dimensions to `0`

. The resulting matrix is extended to an invertible `m × m`

matrix, and the algorithm using the inverse of the linear map is applied. For technical details of extending `M`

to a higher-dimensional invertible matrix, see `ReachabilityBase.Arrays.extend`

.

**Vertex representation**

This algorithm is invoked with the keyword argument `algorithm="vrep"`

(or `algorithm="vrep_chull"`

). If the polyhedron is bounded, the idea is to convert it to its vertex representation and apply the linear map to each vertex.

The returned set is a polytope in vertex representation. Note that conversion of the result back to half-space representation is not computed by default, since this may be costly. If you use this algorithm and still want to convert back to half-space representation, apply `tohrep`

to the result of this method.

`LazySets.an_element`

— Method```
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.

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

Return some element of a universe.

**Input**

`U`

– universe

**Output**

The origin.

`LazySets.isbounded`

— Method`isbounded(P::AbstractPolyhedron{N}; [solver]=default_lp_solver(N)) where {N}`

Check whether a polyhedron is bounded.

**Input**

`P`

– polyhedron`solver`

– (optional, default:`default_lp_solver(N)`

) the backend used to solve the linear program

**Output**

`true`

iff the polyhedron is bounded

**Algorithm**

We first check if the polyhedron has more than `dim(P)`

constraints, which is a necessary condition for boundedness.

If so, we check boundedness via `_isbounded_stiemke`

.

`LazySets.vertices_list`

— Method`vertices_list(P::AbstractPolyhedron; check_boundedness::Bool=true)`

Return the list of vertices of a polyhedron in constraint representation.

**Input**

`P`

– polyhedron in constraint representation`check_boundedness`

– (optional, default:`true`

) if`true`

, check whether the polyhedron is bounded

**Output**

The list of vertices of `P`

, or an error if `P`

is unbounded.

**Notes**

This function throws an error if the polyhedron is unbounded. Otherwise, the polyhedron is converted to an `HPolytope`

and its list of vertices is computed.

**Examples**

```
julia> P = HPolyhedron([HalfSpace([1.0, 0.0], 1.0),
HalfSpace([0.0, 1.0], 1.0),
HalfSpace([-1.0, 0.0], 1.0),
HalfSpace([0.0, -1.0], 1.0)]);
julia> length(vertices_list(P))
4
```

`LazySets.project`

— Method```
project(P::AbstractPolyhedron{N}, block::AbstractVector{Int};
[kwargs...]) where {N}
```

Concrete projection of a polyhedral set.

**Input**

`P`

– set`block`

– block structure, a vector with the dimensions of interest

**Output**

A polyhedron representing the projection of `P`

on the dimensions specified by `block`

. If `P`

was bounded, the result is an `HPolytope`

; otherwise the result is an `HPolyhedron`

. Note that there are more specific methods for specific input types, which give a different output type; e.g., projecting a `Ball1`

results in a `Ball1`

.

**Algorithm**

- We first try to exploit the special case where each of the constraints of
`P`

and`block`

are*compatible*, which is one of the two cases described below. Let`c`

be a constraint of`P`

and let $D_c$ and $D_b$ be the set of dimensions in which`c`

resp.`block`

are constrained.- If $D_c ⊆ D_b$, then one can project the normal vector of
`c`

. - If $D_c ∩ D_b = ∅$, then the constraint becomes redundant.

- If $D_c ⊆ D_b$, then one can project the normal vector of
- In the general case, we compute the concrete linear map of the projection matrix associated to the given block structure.

**Examples**

Consider the four-dimensional cross-polytope (unit ball in the 1-norm):

`julia> P = convert(HPolytope, Ball1(zeros(4), 1.0));`

All dimensions are constrained, and computing the (trivial) projection on the whole space behaves as expected:

```
julia> constrained_dimensions(P)
4-element Vector{Int64}:
1
2
3
4
julia> project(P, [1, 2, 3, 4]) == P
true
```

Each constraint of the cross polytope is constrained in all dimensions.

Now let us take a ball in the infinity norm and remove some constraints:

```
julia> B = BallInf(zeros(4), 1.0);
julia> c = constraints_list(B)[1:2]
2-element Vector{HalfSpace{Float64, ReachabilityBase.Arrays.SingleEntryVector{Float64}}}:
HalfSpace{Float64, ReachabilityBase.Arrays.SingleEntryVector{Float64}}([1.0, 0.0, 0.0, 0.0], 1.0)
HalfSpace{Float64, ReachabilityBase.Arrays.SingleEntryVector{Float64}}([0.0, 1.0, 0.0, 0.0], 1.0)
julia> P = HPolyhedron(c);
julia> constrained_dimensions(P)
2-element Vector{Int64}:
1
2
```

Finally, we take the concrete projection onto variables `1`

and `2`

:

```
julia> project(P, [1, 2]) |> constraints_list
2-element Vector{HalfSpace{Float64, Vector{Float64}}}:
HalfSpace{Float64, Vector{Float64}}([1.0, 0.0], 1.0)
HalfSpace{Float64, Vector{Float64}}([0.0, 1.0], 1.0)
```

`LazySets._isbounded_stiemke`

— Function```
_isbounded_stiemke(constraints::AbstractVector{<:HalfSpace{N}};
solver=LazySets.default_lp_solver(N),
check_nonempty::Bool=true) where {N}
```

Check whether a list of constraints is bounded using Stiemke's theorem of alternatives.

**Input**

`constraints`

– list of constraints`backend`

– (optional, default:`default_lp_solver(N)`

) the backend used to solve the linear program`check_nonempty`

– (optional, default:`true`

) if`true`

, check the precondition to this algorithm that`P`

is non-empty

**Output**

`true`

iff the list of constraints is bounded.

**Notes**

The list of constraints represents a polyhedron.

The algorithm calls `isempty`

to check whether the polyhedron is empty. This computation can be avoided using the `check_nonempty`

flag.

**Algorithm**

The algorithm is based on Stiemke's theorem of alternatives, see, e.g., [1].

Let the polyhedron $P$ be given in constraint form $Ax ≤ b$. We assume that the polyhedron is non-empty.

Proposition 1. If $\ker(A)≠\{0\}$, then $P$ is unbounded.

Proposition 2. Assume that $ker(A)={0}$ and $P$ is non-empty. Then $P$ is bounded if and only if the following linear program admits a feasible solution: $\min∥y∥_1$ subject to $A^Ty=0$ and $y≥1$.

[1] Mangasarian, Olvi L. *Nonlinear programming.* Society for Industrial and Applied Mathematics, 1994.

Some common functions to work with linear constraints:

`LazySets.constraints_list`

— Method`constraints_list(A::AbstractMatrix{N}, b::AbstractVector)`

Convert a matrix-vector representation to a linear-constraint representation.

**Input**

`A`

– matrix`b`

– vector

**Output**

A list of linear constraints.

`LazySets.tosimplehrep`

— Method```
tosimplehrep(constraints::AbstractVector{LC};
[n]::Int=0) where {N, LC<:HalfSpace{N}}
```

Return the simple H-representation $Ax ≤ b$ from a list of linear constraints.

**Input**

`constraints`

– a list of linear constraints`n`

– (optional; default:`0`

) dimension of the constraints

**Output**

The tuple `(A, b)`

where `A`

is the matrix of normal directions and `b`

is the vector of offsets.

**Notes**

The parameter `n`

can be used to create a matrix with no constraints but a non-zero dimension.

`LazySets.remove_redundant_constraints`

— Method```
remove_redundant_constraints(constraints::AbstractVector{S};
backend=nothing) where {S<:HalfSpace}
```

Remove the redundant constraints of a given list of linear constraints.

**Input**

`constraints`

– list of constraints`backend`

– (optional, default:`nothing`

) the backend used to solve the linear program

**Output**

The list of constraints with the redundant ones removed, or an empty set if the constraints are infeasible.

**Notes**

If `backend`

is `nothing`

, it defaults to `default_lp_solver(N)`

.

**Algorithm**

See `remove_redundant_constraints!(::AbstractVector{<:HalfSpace})`

for details.

`LazySets.remove_redundant_constraints!`

— Method```
remove_redundant_constraints!(constraints::AbstractVector{S};
[backend]=nothing) where {S<:HalfSpace}
```

Remove the redundant constraints of a given list of linear constraints; the list is updated in-place.

**Input**

`constraints`

– list of constraints`backend`

– (optional, default:`nothing`

) the backend used to solve the linear program

**Output**

`true`

if the removal was successful and the list of constraints `constraints`

is modified by removing the redundant constraints, and `false`

only if the constraints are infeasible.

**Notes**

Note that the result may be `true`

even if the constraints are infeasible. For example, $x ≤ 0 && x ≥ 1$ will return `true`

without removing any constraint. To check if the constraints are infeasible, use `isempty(HPolyhedron(constraints))`

.

If `backend`

is `nothing`

, it defaults to `default_lp_solver(N)`

.

**Algorithm**

If there are `m`

constraints in `n`

dimensions, this function checks one by one if each of the `m`

constraints is implied by the remaining ones.

To check if the `k`

-th constraint is redundant, an LP is formulated using the constraints that have not yet been removed. If, at an intermediate step, it is detected that a subgroup of the constraints is infeasible, this function returns `false`

. If the calculation finished successfully, this function returns `true`

.

For details, see Fukuda's Polyhedra FAQ.

Plotting (bounded) polyhedra is available, too:

`LazySets.plot_recipe`

— Method`plot_recipe(P::AbstractPolyhedron{N}, [ε]=zero(N)) where {N}`

Convert a (bounded) polyhedron to a pair `(x, y)`

of points for plotting.

**Input**

`P`

– bounded polyhedron`ε`

– (optional, default:`0`

) ignored, used for dispatch

**Output**

A pair `(x, y)`

of points that can be plotted, where `x`

is the vector of x-coordinates and `y`

is the vector of y-coordinates.

**Algorithm**

We first assert that `P`

is bounded (i.e., that `P`

is a polytope).

One-dimensional polytopes are converted to an `Interval`

. Three-dimensional or higher-dimensional polytopes are not supported.

For two-dimensional polytopes (i.e., polygons) we compute their set of vertices using `vertices_list`

and then plot the convex hull of these vertices.

### Implementations

- Half-space (HalfSpace)
- Polyhedron in constraint representation (HPolyhedron)
- Hyperplane
- Line2D
- Line
- Universe

## Polytopes (AbstractPolytope)

A polytope is a bounded set with finitely many vertices (*V-representation*) resp. facets (*H-representation*). Note that there is a special interface combination Centrally symmetric polytope.

`LazySets.AbstractPolytope`

— Type`AbstractPolytope{N} <: AbstractPolyhedron{N}`

Abstract type for compact convex polytopic sets.

**Notes**

Every concrete `AbstractPolytope`

must define the following method:

`vertices_list(::AbstractPolytope)`

– return a list of all vertices

```
julia> subtypes(AbstractPolytope)
4-element Vector{Any}:
AbstractCentrallySymmetricPolytope
AbstractPolygon
HPolytope
VPolytope
```

A polytope is a bounded polyhedron (see `AbstractPolyhedron`

). Polytopes are compact convex sets with either of the following equivalent properties:

- They are the intersection of a finite number of closed half-spaces.
- They are the convex hull of finitely many vertices.

This interface defines the following functions:

`LazySets.isbounded`

— Method`isbounded(P::AbstractPolytope)`

Check whether a polytopic set is bounded.

**Input**

`P`

– polytopic set

**Output**

`true`

(since a polytopic set must be bounded).

`LazySets.isuniversal`

— Method`isuniversal(P::AbstractPolytope{N}, [witness]::Bool=false) where {N}`

Check whether a polytopic set is universal.

**Input**

`P`

– polytopic set`witness`

– (optional, default:`false`

) compute a witness if activated

**Output**

- If
`witness`

option is deactivated:`false`

- If
`witness`

option is activated:`(false, v)`

where $v ∉ P$ unless the list of constraints is empty (which should not happen for a normal polytope)

**Algorithm**

A witness is produced using `isuniversal(H)`

where `H`

is the first linear constraint of `P`

.

`Base.isempty`

— Method`isempty(P::AbstractPolytope)`

Check whether a polytopic set is empty.

**Input**

`P`

– polytopic set

**Output**

`true`

if the given polytopic set contains no vertices, and `false`

otherwise.

**Algorithm**

This algorithm checks whether the `vertices_list`

of `P`

is empty.

`LazySets.volume`

— Method`volume(P::AbstractPolytope; backend=default_polyhedra_backend(P))`

Compute the volume of a polytopic set.

**Input**

`P`

– polytopic set`backend`

– (optional, default:`default_polyhedra_backend(P)`

) the backend for polyhedral computations; see Polyhedra's documentation for further information

**Output**

The volume of `P`

.

**Algorithm**

The volume is computed by the `Polyhedra`

library.

### Implementations

## Polygons (AbstractPolygon)

A polygon is a two-dimensional polytope.

`LazySets.AbstractPolygon`

— Type`AbstractPolygon{N} <: AbstractPolytope{N}`

Abstract type for convex polygons (i.e., two-dimensional polytopes).

**Notes**

Every concrete `AbstractPolygon`

must define the following functions:

`tovrep(::AbstractPolygon{N})`

– transform into vertex representation`tohrep(::AbstractPolygon{N})`

– transform into constraint representation

The subtypes of `AbstractPolygon`

(including abstract interfaces):

```
julia> subtypes(AbstractPolygon)
2-element Vector{Any}:
AbstractHPolygon
VPolygon
```

This interface defines the following functions:

`LazySets.dim`

— Method`dim(P::AbstractPolygon)`

Return the ambient dimension of a convex polygon.

**Input**

`P`

– convex polygon

**Output**

The ambient dimension of the polygon, which is 2.

The following helper functions are used for sorting directions:

`LazySets.jump2pi`

— Function`jump2pi(x::N) where {N<:AbstractFloat}`

Return $x + 2π$ if $x$ is negative, otherwise return $x$.

**Input**

`x`

– real scalar

**Output**

$x + 2π$ if $x$ is negative, $x$ otherwise.

**Examples**

```
julia> using LazySets: jump2pi
julia> jump2pi(0.0)
0.0
julia> jump2pi(-0.5)
5.783185307179586
julia> jump2pi(0.5)
0.5
```

`Base.:<=`

— Method`<=(u::AbstractVector, v::AbstractVector)`

Compare two 2D vectors by their direction.

**Input**

`u`

– first 2D direction`v`

– second 2D direction

**Output**

`true`

iff $\arg(u) [2π] ≤ \arg(v) [2π]$.

**Notes**

The argument is measured in counter-clockwise fashion, with the 0 being the direction (1, 0).

**Algorithm**

The implementation checks the quadrant of each direction, and compares directions using the right-hand rule. In particular, this method does not use the arctangent.

`LazySets._leq_trig`

— Method`_leq_trig(u::AbstractVector{N}, v::AbstractVector{N}) where {N<:AbstractFloat}`

Compare two 2D vectors by their direction.

**Input**

`u`

– first 2D direction`v`

– second 2D direction

**Output**

`true`

iff $\arg(u) [2π] ≤ \arg(v) [2π]$.

**Notes**

The argument is measured in counter-clockwise fashion, with the 0 being the direction (1, 0).

**Algorithm**

The implementation uses the arctangent function with sign, `atan`

, which for two arguments implements the `atan2`

function.

`LazySets.quadrant`

— Method`quadrant(w::AbstractVector{N}) where {N}`

Compute the quadrant where the direction `w`

belongs.

**Input**

`w`

– direction

**Output**

An integer from 0 to 3, with the following convention:

```
^
1 | 0
---+-->
2 | 3
```

**Algorithm**

The idea is to encode the following logic function: $11 ↦ 0, 01 ↦ 1, 00 ↦ 2, 10 ↦ 3$, according to the convention above.

This function is inspired from AGPX's answer in: Sort points in clockwise order?

### Implementations

## Polygons in constraint representation (AbstractHPolygon)

An HPolygon is a polygon in H-representation (or constraint representation).

`LazySets.AbstractHPolygon`

— Type`AbstractHPolygon{N} <: AbstractPolygon{N}`

Abstract type for polygons in constraint representation.

**Notes**

All subtypes must satisfy the invariant that constraints are sorted counter-clockwise.

Every concrete `AbstractHPolygon`

must have the following fields:

`constraints::Vector{HalfSpace{N, AbstractVector{N}}}`

– the constraints

The subtypes of `AbstractHPolygon`

:

```
julia> subtypes(AbstractHPolygon)
2-element Vector{Any}:
HPolygon
HPolygonOpt
```

This interface defines the following functions:

`LazySets.an_element`

— Method`an_element(P::AbstractHPolygon)`

Return some element of a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation

**Output**

A vertex of the polygon in constraint representation (the first one in the order of the constraints).

`Base.:∈`

— Method`∈(x::AbstractVector, P::AbstractHPolygon)`

Check whether a given two-dimensional point is contained in a polygon in constraint representation.

**Input**

`x`

– two-dimensional point/vector`P`

– polygon in constraint representation

**Output**

`true`

iff $x ∈ P$.

**Algorithm**

This implementation checks if the point lies inside each constraint.

`Base.rand`

— Method```
rand(::Type{HPOLYGON}; [N]::Type=Float64, [dim]::Int=2,
[rng]::AbstractRNG=GLOBAL_RNG, [seed]::Union{Int, Nothing}=nothing,
[num_constraints]::Int=-1) where {HPOLYGON<:AbstractHPolygon}
```

Create a random polygon in constraint representation.

**Input**

`HPOLYGON`

– type for dispatch`N`

– (optional, default:`Float64`

) numeric type`dim`

– (optional, default: 2) dimension`rng`

– (optional, default:`GLOBAL_RNG`

) random number generator`seed`

– (optional, default:`nothing`

) seed for reseeding`num_constraints`

– (optional, default:`-1`

) number of constraints of the polygon (must be 3 or bigger; see comment below)

**Output**

A random polygon in constraint representation.

**Algorithm**

We create a random polygon in vertex representation and convert it to constraint representation. See `rand(::Type{VPolygon})`

. For non-flat polygons the number of vertices and the number of constraints are identical.

`LazySets.tohrep`

— Method`tohrep(P::HPOLYGON) where {HPOLYGON<:AbstractHPolygon}`

Build a constraint representation of the given polygon.

**Input**

`P`

– polygon in constraint representation

**Output**

The identity, i.e., the same polygon instance.

`LazySets.tovrep`

— Method`tovrep(P::AbstractHPolygon)`

Build a vertex representation of a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation

**Output**

The same polygon but in vertex representation, a `VPolygon`

.

`LazySets.addconstraint!`

— Method```
addconstraint!(P::AbstractHPolygon, constraint::HalfSpace;
[linear_search]::Bool=length(P.constraints) < 10,
[prune]::Bool=true)
```

Add a linear constraint to a polygon in constraint representation, keeping the constraints sorted by their normal directions.

**Input**

`P`

– polygon in constraint representation`constraint`

– linear constraint to add`linear_search`

– (optional, default:`length(constraints) < 10`

) flag to choose between linear and binary search`prune`

– (optional, default:`true`

) flag for removing redundant constraints in the end

`LazySets.addconstraint!`

— Method```
addconstraint!(constraints::Vector{LC}, new_constraint::HalfSpace;
[linear_search]::Bool=length(P.constraints) < 10,
[prune]::Bool=true) where {LC<:HalfSpace}
```

Add a linear constraint to a sorted vector of constrains, keeping the constraints sorted by their normal directions.

**Input**

`constraints`

– vector of linear constraints`new_constraint`

– linear constraint to add`linear_search`

– (optional, default:`length(constraints) < 10`

) flag to choose between linear and binary search`prune`

– (optional, default:`true`

) flag for removing redundant constraints in the end

**Algorithm**

If `prune`

is active, we check if the new constraint is redundant. If the constraint is not redundant, we perform the same check to the left and to the right until we find the first constraint that is not redundant.

`LinearAlgebra.normalize`

— Method`normalize(P::AbstractHPolygon{N}, p=N(2)) where {N}`

Normalize a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation`p`

– (optional, default:`2`

) norm

**Output**

A new polygon in constraint representation whose normal directions $a_i$ are normalized, i.e., such that $‖a_i‖_p = 1$ holds.

`LazySets.isredundant`

— Method`isredundant(cmid::HalfSpace, cright::HalfSpace, cleft::HalfSpace)`

Check whether a linear constraint is redundant wrt. two surrounding constraints.

**Input**

`cmid`

– linear constraint of concern`cright`

– linear constraint to the right (clockwise turn)`cleft`

– linear constraint to the left (counter-clockwise turn)

**Output**

`true`

iff the constraint is redundant.

**Algorithm**

We first check whether the angle between the surrounding constraints is < 180°, which is a necessary condition (unless the direction is identical to one of the other two constraints). If so, we next check if the angle is 0°, in which case the constraint `cmid`

is redundant unless it is strictly tighter than the other two constraints. If the angle is strictly between 0° and 180°, the constraint `cmid`

is redundant if and only if the vertex defined by the other two constraints lies inside the set defined by `cmid`

.

`LazySets.remove_redundant_constraints!`

— Method`remove_redundant_constraints!(P::AbstractHPolygon)`

Remove all redundant constraints of a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation

**Output**

The same polygon with all redundant constraints removed.

**Notes**

Since we only consider bounded polygons and a polygon needs at least three constraints to be bounded, we stop removing redundant constraints if there are three or fewer constraints left. Hence for unbounded polygons the result may be unexpected.

**Algorithm**

We go through all consecutive triples of constraints and check if the one in the middle is redundant. For this we assume that the constraints are sorted.

`LazySets.constraints_list`

— Method`constraints_list(P::AbstractHPolygon)`

Return the list of constraints defining a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation

**Output**

The list of constraints of the polygon. The implementation guarantees that the constraints are sorted counter-clockwise.

`LazySets.vertices_list`

— Method```
vertices_list(P::AbstractHPolygon{N};
apply_convex_hull::Bool=true,
check_feasibility::Bool=true) where {N}
```

Return the list of vertices of a polygon in constraint representation.

**Input**

`P`

– polygon in constraint representation`apply_convex_hull`

– (optional, default:`true`

) flag to post-process the intersection of constraints with a convex hull`check_feasibility`

– (optional, default:`true`

) flag to check whether the polygon was empty (required for correctness in case of empty polygons)

**Output**

List of vertices.

**Notes**

By construction an `AbstractHPolygon`

should not contain any redundant vertices. Still the `apply_convex_hull`

argument is activated by default to remove potential duplicate vertices. They can exist due to numeric instability.

```
julia> p = HPolygon([HalfSpace([1.0, 0.0], 1.0),
HalfSpace([0.0, 1.0], 1.0),
HalfSpace([-1.0, 0.0], -1.0),
HalfSpace([0.0, -1.0], -1.0)]);
julia> vertices_list(p, apply_convex_hull=false)
4-element Vector{Vector{Float64}}:
[1.0, 1.0]
[1.0, 1.0]
[1.0, 1.0]
[1.0, 1.0]
```

If it is known that each constraint has a "proper" distance to the next vertex, this step can be skipped.

**Algorithm**

We compute each vertex as the intersection of consecutive lines defined by the half-spaces. If `check_feasibility`

is active, we then check if the constraints of the polygon were actually feasible (i.e., they pointed in the *right* direction). For this we compute the *average* of all vertices and check membership in each constraint.

`LazySets.isbounded`

— Function`isbounded(P::AbstractHPolygon, [use_type_assumption]::Bool=true)`

Determine whether a polygon in constraint representation is bounded.

**Input**

`P`

– polygon in constraint representation`use_type_assumption`

– (optional, default:`true`

) flag for ignoring the type assumption that polygons are bounded

**Output**

`true`

if `use_type_assumption`

is activated. Otherwise, `true`

iff `P`

is bounded.

**Algorithm**

If `!use_type_assumption`

, we use `_isbounded_unit_dimensions`

.

### Implementations

- Polygon in constraint representation (HPolygon)
- Polygon in optimized constraint representation (HPolygonOpt)

## Centrally symmetric polytopes (AbstractCentrallySymmetricPolytope)

A centrally symmetric polytope is a combination of two other interfaces: Centrally symmetric sets and Polytope.

`LazySets.AbstractCentrallySymmetricPolytope`

— Type`AbstractCentrallySymmetricPolytope{N} <: AbstractPolytope{N}`

Abstract type for centrally symmetric, polytopic sets. It combines the `AbstractCentrallySymmetric`

and `AbstractPolytope`

interfaces. Such a type combination is necessary as long as Julia does not support multiple inheritance.

**Notes**

Every concrete `AbstractCentrallySymmetricPolytope`

must define the following methods:

- from
`AbstractCentrallySymmetric`

:`center(::AbstractCentrallySymmetricPolytope)`

– return the center point`center(::AbstractCentrallySymmetricPolytope, i::Int)`

– return the center point at index`i`

- from
`AbstractPolytope`

:`vertices_list(::AbstractCentrallySymmetricPolytope)`

– return a list of all vertices

The subtypes of `AbstractCentrallySymmetricPolytope`

(including abstract interfaces):

```
julia> subtypes(AbstractCentrallySymmetricPolytope)
2-element Vector{Any}:
AbstractZonotope
Ball1
```

This interface defines the following functions:

`LazySets.dim`

— Method`dim(P::AbstractCentrallySymmetricPolytope)`

Return the ambient dimension of a centrally symmetric, polytopic set.

**Input**

`P`

– centrally symmetric, polytopic set

**Output**

The ambient dimension of the polytopic set.

`LazySets.an_element`

— Method`an_element(P::AbstractCentrallySymmetricPolytope)`

Return some element of a centrally symmetric, polytopic set.

**Input**

`P`

– centrally symmetric, polytopic set

**Output**

The center of the centrally symmetric, polytopic set.

`Base.isempty`

— Method`isempty(P::AbstractCentrallySymmetricPolytope)`

Check whether a centrally symmetric, polytopic set is empty.

**Input**

`P`

– centrally symmetric, polytopic set

**Output**

`false`

.

`LazySets.isuniversal`

— Method```
isuniversal(S::AbstractCentrallySymmetricPolytope{N},
[witness]::Bool=false) where {N}
```

Check whether a centrally symmetric, polytopic set is universal.

**Input**

`S`

– centrally symmetric, polytopic set`witness`

– (optional, default:`false`

) compute a witness if activated

**Output**

- If
`witness`

option is deactivated:`false`

- If
`witness`

option is activated:`(false, v)`

where $v ∉ S$

**Algorithm**

Centrally symmetric, polytopic sets are bounded. A witness is obtained by computing the support vector in direction `d = [1, 0, …, 0]`

and adding `d`

on top.

`LazySets.center`

— Method`center(S::AbstractCentrallySymmetricPolytope, i::Int)`

Return the center of a centrally symmetric, polytopic set along a given dimension.

**Input**

`S`

– centrally symmetric, polytopic set`i`

– dimension of interest

**Output**

The center along the given dimension.

`Base.extrema`

— Method`extrema(S::AbstractCentrallySymmetricPolytope, i::Int)`

Return the lower and higher coordinate of a centrally symmetric, polytopic set in a given dimension.

**Input**

`S`

– centrally symmetric, polytopic set`i`

– dimension of interest

**Output**

The lower and higher coordinate of the centrally symmetric, polytopic set in the given dimension.

**Notes**

The result is equivalent to `(low(S, i), high(S, i))`

.

**Algorithm**

We compute `high(S, i)`

and then compute the lowest coordinates with the help of `center(S, i)`

(which is assumed to be cheaper to obtain).

`Base.extrema`

— Method`extrema(S::AbstractCentrallySymmetricPolytope)`

Return two vectors with the lowest and highest coordinate of a centrally symmetric, polytopic set.

**Input**

`S`

– centrally symmetric, polytopic set

**Output**

Two vectors with the lowest and highest coordinates of `S`

.

**Notes**

The result is equivalent to `(low(S), high(S))`

.

**Algorithm**

We compute `high(S)`

and then compute the lowest coordinates with the help of `center(S)`

(which is assumed to be cheaper to obtain).

### Implementations

## Zonotopes (AbstractZonotope)

A zonotope is a specific centrally symmetric polytope characterized by a center and a collection of generators.

`LazySets.AbstractZonotope`

— Type`AbstractZonotope{N} <: AbstractCentrallySymmetricPolytope{N}`

Abstract type for zonotopic sets.

**Notes**

Mathematically, a zonotope is defined as the set

\[Z = \left\{ c + ∑_{i=1}^p ξ_i g_i,~~ ξ_i \in [-1, 1]~~ ∀ i = 1,…, p \right\},\]

where $c \in \mathbb{R}^n$ is its center and $\{g_i\}_{i=1}^p$, $g_i \in \mathbb{R}^n$, is the set of generators. This characterization defines a zonotope as the finite Minkowski sum of line segments. Zonotopes can be equivalently described as the image of a unit infinity-norm ball in $\mathbb{R}^n$ by an affine transformation.

See `Zonotope`

for a standard implementation of this interface.

Every concrete `AbstractZonotope`

must define the following functions:

`genmat(::AbstractZonotope)`

– return the generator matrix`generators(::AbstractZonotope)`

– return an iterator over the generators

Since the functions `genmat`

and `generators`

can be defined in terms of each other, it is sufficient to only genuinely implement one of them and let the implementation of the other function call the fallback implementation `genmat_fallback`

resp. `generators_fallback`

.

The subtypes of `AbstractZonotope`

(including abstract interfaces):

```
julia> subtypes(AbstractZonotope)
5-element Vector{Any}:
AbstractHyperrectangle
HParallelotope
LineSegment
RotatedHyperrectangle
Zonotope
```

This interface defines the following functions:

`LazySets.ngens`

— Method`ngens(Z::AbstractZonotope)`

Return the number of generators of a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

An integer representing the number of generators.

`LazySets.genmat_fallback`

— Method```
genmat_fallback(Z::AbstractZonotope{N};
[gens]=generators(Z), [ngens]=nothing) where {N}
```

Fallback definition of `genmat`

for zonotopic sets.

**Input**

`Z`

– zonotopic set`gens`

– (optional; default:`generators(Z)`

) iterator over generators`ngens`

– (optional; default:`nothing`

) number of generators or`nothing`

if unknown

**Output**

A matrix where each column represents one generator of `Z`

.

**Notes**

Passing the number of generators (`ngens`

) is more efficient, since otherwise the generators have to be obtained from the iterator (`gens`

) and stored in an intermediate vector until the final result matrix can be allocated.

`LazySets.generators_fallback`

— Method`generators_fallback(Z::AbstractZonotope)`

Fallback definition of `generators`

for zonotopic sets.

**Input**

`Z`

– zonotopic set

**Output**

An iterator over the generators of `Z`

.

`LazySets.ρ`

— Method`ρ(d::AbstractVector, Z::AbstractZonotope)`

Evaluate the support function of a zonotopic set in a given direction.

**Input**

`d`

– direction`Z`

– zonotopic set

**Output**

The evaluation of the support function in the given direction.

**Algorithm**

The support value is $cᵀ d + ‖Gᵀ d‖₁$, where $c$ is the center and $G$ is the generator matrix of `Z`

.

`LazySets.σ`

— Method`σ(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.

`Base.:∈`

— Method`∈(x::AbstractVector, Z::AbstractZonotope; solver=nothing)`

Check whether a given point is contained in a zonotopic set.

**Input**

`x`

– point/vector`Z`

– zonotopic set`solver`

– (optional, default:`nothing`

) the backend used to solve the linear program

**Output**

`true`

iff $x ∈ Z$.

**Examples**

```
julia> Z = Zonotope([1.0, 0.0], [0.1 0.0; 0.0 0.1]);
julia> [1.0, 0.2] ∈ Z
false
julia> [1.0, 0.1] ∈ Z
true
```

**Notes**

If `solver == nothing`

, we fall back to `default_lp_solver(N)`

.

**Algorithm**

The membership problem is computed by stating and solving the following linear program. Let $p$ and $n$ be the number of generators and ambient dimension, respectively. We consider the minimization of $x_0$ in the $p+1$-dimensional space of elements $(x_0, ξ_1, …, ξ_p)$ constrained to $0 ≤ x_0 ≤ ∞$, $ξ_i ∈ [-1, 1]$ for all $i = 1, …, p$, and such that $x-c = Gξ$ holds. If a feasible solution exists, the optimal value $x_0 = 0$ is achieved.

`LazySets.linear_map`

— Method`linear_map(M::AbstractMatrix, Z::AbstractZonotope)`

Concrete linear map of a zonotopic set.

**Input**

`M`

– matrix`Z`

– zonotopic set

**Output**

The zonotope obtained by applying the linear map.

**Algorithm**

We apply the linear map to the center and the generators.

`LazySets.translate`

— Method`translate(Z::AbstractZonotope, v::AbstractVector)`

Translate (i.e., shift) a zonotopic set by a given vector.

**Input**

`Z`

– zonotopic set`v`

– translation vector

**Output**

A translated zonotopic set.

**Notes**

See also `translate!(Z::AbstractZonotope, v::AbstractVector)`

for the in-place version.

`LazySets.translate!`

— Method`translate!(Z::AbstractZonotope, v::AbstractVector)`

Translate (i.e., shift) a zonotopic set by a given vector in-place.

**Input**

`Z`

– zonotopic set`v`

– translation vector

**Output**

A translated zonotopic set.

**Notes**

See also `translate(Z::AbstractZonotope, v::AbstractVector)`

for the out-of-place version.

**Algorithm**

We add the translation vector to the center of the zonotopic set.

`Base.split`

— Method`split(Z::AbstractZonotope, j::Int)`

Return two zonotopes obtained by splitting the given zonotopic set.

**Input**

`Z`

– zonotopic set`j`

– index of the generator to be split

**Output**

The zonotope obtained by splitting `Z`

into two zonotopes such that their union is `Z`

and their intersection is possibly non-empty.

**Algorithm**

This function implements [Prop. 3, 1], which we state next. The zonotopic set $Z = ⟨c, g^{(1, …, p)}⟩$ is split into:

\[Z₁ = ⟨c - \frac{1}{2}g^{(j)}, (g^{(1, …,j-1)}, \frac{1}{2}g^{(j)}, g^{(j+1, …, p)})⟩ \\ Z₂ = ⟨c + \frac{1}{2}g^{(j)}, (g^{(1, …,j-1)}, \frac{1}{2}g^{(j)}, g^{(j+1, …, p)})⟩,\]

such that $Z₁ ∪ Z₂ = Z$ and $Z₁ ∩ Z₂ = Z^*$, where

\[Z^* = ⟨c, (g^{(1,…,j-1)}, g^{(j+1,…, p)})⟩.\]

[1] Althoff, M., Stursberg, O., & Buss, M. *Reachability analysis of nonlinear systems with uncertain parameters using conservative linearization*. CDC 2008.

`Base.split`

— Method```
split(Z::AbstractZonotope, gens::AbstractVector{Int},
nparts::AbstractVector{Int})
```

Split a zonotopic set along the given generators into a vector of zonotopes.

**Input**

`Z`

– zonotopic set`gens`

– vector of indices of the generators to be split`n`

– vector of integers describing the number of partitions in the corresponding generator

**Output**

The zonotopes obtained by splitting `Z`

into `2^{n_i}`

zonotopes for each generator `i`

such that their union is `Z`

and their intersection is possibly non-empty.

**Examples**

Splitting of a two-dimensional zonotopic set along its first generator:

```
julia> Z = Zonotope([1.0, 0.0], [0.1 0.0; 0.0 0.1])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.0, 0.0], [0.1 0.0; 0.0 0.1])
julia> split(Z, [1], [1])
2-element Vector{Zonotope{Float64, Vector{Float64}, Matrix{Float64}}}:
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.95, 0.0], [0.05 0.0; 0.0 0.1])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.05, 0.0], [0.05 0.0; 0.0 0.1])
```

Here, the first vector in the arguments corresponds to the zonotopic set's generator to be split, and the second vector corresponds to the exponent of `2^n`

parts that the set will be split into along the corresponding generator.

As an example, below we split a two-dimensional zonotope along both of its generators, each time into four parts.

```
julia> Z = Zonotope([1.0, 0.0], [0.1 0.0; 0.0 0.1])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.0, 0.0], [0.1 0.0; 0.0 0.1])
julia> split(Z, [1, 2], [2, 2])
16-element Vector{Zonotope{Float64, Vector{Float64}, Matrix{Float64}}}:
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.925, -0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.925, -0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.925, 0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.925, 0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.975, -0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.975, -0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.975, 0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([0.975, 0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.025, -0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.025, -0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.025, 0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.025, 0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.075, -0.075], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.075, -0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.075, 0.025], [0.025 0.0; 0.0 0.025])
Zonotope{Float64, Vector{Float64}, Matrix{Float64}}([1.075, 0.075], [0.025 0.0; 0.0 0.025])
```

`LazySets.constraints_list`

— Method`constraints_list(P::AbstractZonotope)`

Return a list of constraints defining a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

A list of constraints of the zonotopic set.

**Algorithm**

This is the (inefficient) fallback implementation for rational numbers. It first computes the vertices and then converts the corresponding polytope to constraint representation.

`LazySets.constraints_list`

— Method`constraints_list(Z::AbstractZonotope{N}) where {N<:AbstractFloat}`

Return a list of constraints defining a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

A list of constraints of the zonotopic set.

**Notes**

The main algorithm assumes that the generator matrix is full rank. The result has $2 \binom{p}{n-1}$ (with $p$ being the number of generators and $n$ being the ambient dimension) constraints, which is optimal under this assumption. If this assumption is not given, the implementation tries to work around.

**Algorithm**

We follow the algorithm presented in [1]. Three cases are not covered by that algorithm, so we handle them separately. The first case is zonotopes in one dimension. The second case is that there are fewer generators than dimensions, $p < n$, or the generator matrix is not full rank, in which case we fall back to the (slower) computation based on the vertex representation. The third case is that the zonotope is flat in some dimensions, in which case we project the zonotope to the non-flat dimensions and extend the result later.

[1] Althoff, Stursberg, Buss. *Computing Reachable Sets of Hybrid Systems Using a Combination of Zonotopes and Polytopes*. 2009.

`LazySets.vertices_list`

— Method`vertices_list(Z::AbstractZonotope; [apply_convex_hull]::Bool=true)`

Return a list of the vertices of a zonotopic set.

**Input**

`Z`

– zonotopic set`apply_convex_hull`

– (optional, default:`true`

) if`true`

, post-process the computation with the convex hull of the points

**Output**

A list of the vertices.

**Algorithm**

**Two-dimensional case**

We use a trick to speed up enumerating vertices of 2-dimensional zonotopic sets with all generators in the first quadrant or third quadrant (same sign). Namely, sort the generators by angle and add them clockwise in increasing order and counterclockwise in decreasing order. A more detailed explanation can be found here.

To avoid the cumulative sum from both directions separately, we build a 2D index matrix to sum generators for both directions in one matrix-vector product.

**General case**

If the zonotopic set has $p$ generators, each vertex is the result of summing the center with some linear combination of generators, where the combination factors are $ξ_i ∈ \{-1, 1\}$.

There are at most $2^p$ distinct vertices. Use the flag `apply_convex_hull`

to control whether a convex-hull algorithm is applied to the vertices computed by this method; otherwise, redundant vertices may be present.

`LazySets.order`

— Method`order(Z::AbstractZonotope)`

Return the order of a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

A rational number representing the order of the zonotopic set.

**Notes**

The order of a zonotopic set is defined as the quotient of its number of generators and its dimension.

`LazySets.togrep`

— Method`togrep(Z::AbstractZonotope)`

Return a generator representation of a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

The same set in generator representation. This fallback implementation returns a `Zonotope`

; however, more specific implementations may return other generator representations.

`LazySets.remove_redundant_generators`

— Method`remove_redundant_generators(Z::AbstractZonotope)`

Remove all redundant (pairwise linearly dependent) generators of a zonotopic set.

**Input**

`Z`

– zonotopic set

**Output**

A new zonotope with fewer generators, or the same zonotopic set if no generator could be removed.

**Algorithm**

By default this implementation returns the input zonotopic set. Subtypes of `AbstractZonotope`

whose generators can be removed have to define a new method.

`LazySets.reduce_order`

— Function```
reduce_order(Z::AbstractZonotope, r::Real,
[method]::AbstractReductionMethod=GIR05())
```

Reduce the order of a zonotopic set by overapproximating with a zonotope with fewer generators.

**Input**

`Z`

– zonotopic set`r`

– desired order`method`

– (optional, default:`GIR05()`

) the reduction method used

**Output**

A new zonotope with fewer generators, if possible.

**Algorithm**

The available algorithms are:

```
julia> subtypes(AbstractReductionMethod)
3-element Vector{Any}:
LazySets.ASB10
LazySets.COMB03
LazySets.GIR05
```

See the documentation of each algorithm for references. These methods split the given zonotopic set `Z`

into two zonotopes, `K`

and `L`

, where `K`

contains the most "representative" generators and `L`

contains the generators that are reduced, `Lred`

, using a box overapproximation. We follow the notation from [1]. See also [2].

- [1] Yang, X., & Scott, J. K. *A comparison of zonotope order reduction

techniques*. Automatica 2018.

- [2] Kopetzki, A. K., Schürmann, B., & Althoff, M. *Methods for order reduction

of zonotopes*. CDC 2017.

- [3] Althoff, M., Stursberg, O., & Buss, M. *Computing reachable sets of hybrid

systems using a combination of zonotopes and polytopes*. Nonlinear analysis: hybrid systems 2010.

### Order reduction methods

`LazySets.AbstractReductionMethod`

— Type`AbstractReductionMethod`

Abstract supertype for order-reduction methods of a zonotopic set.

`LazySets.ASB10`

— Type`ASB10 <: AbstractReductionMethod`

Zonotope order-reduction method from [1].

- [1] Althoff, M., Stursberg, O., & Buss, M. *Computing reachable sets of hybrid

systems using a combination of zonotopes and polytopes*. Nonlinear analysis: hybrid systems 2010.

`LazySets.COMB03`

— Type`COMB03 <: AbstractReductionMethod`

Zonotope order-reduction method from [1].

- [1] C. Combastel.
*A state bounding observer based on zonotopes*. ECC 2003.

`LazySets.GIR05`

— Type`GIR05 <: AbstractReductionMethod`

Zonotope order-reduction method from [1].

- [1] A. Girard.
*Reachability of Uncertain Linear Systems Using Zonotopes*.

HSCC 2005.

### Implementations

## Hyperrectangles (AbstractHyperrectangle)

A hyperrectangle is a special centrally symmetric polytope with axis-aligned facets.

`LazySets.AbstractHyperrectangle`

— Type`AbstractHyperrectangle{N} <: AbstractZonotope{N}`

Abstract type for hyperrectangular sets.

**Notes**

See `Hyperrectangle`

for a standard implementation of this interface.

Every concrete `AbstractHyperrectangle`

must define the following functions:

`radius_hyperrectangle(::AbstractHyperrectangle)`

– return the hyperrectangle's radius, which is a full-dimensional vector`radius_hyperrectangle(::AbstractHyperrectangle, i::Int)`

– return the hyperrectangle's radius in the`i`

-th dimension`isflat(::AbstractHyperrectangle)`

– check whether the hyperrectangle's radius is zero in some dimension

Every hyperrectangular set is also a zonotopic set; see `AbstractZonotope`

.

The subtypes of `AbstractHyperrectangle`

(including abstract interfaces):

```
julia> subtypes(AbstractHyperrectangle)
5-element Vector{Any}:
AbstractSingleton
BallInf
Hyperrectangle
Interval
SymmetricIntervalHull
```

This interface defines the following functions:

`LinearAlgebra.norm`

— Function`norm(H::AbstractHyperrectangle, [p]::Real=Inf)`

Return the norm of a hyperrectangular set.

The norm of a hyperrectangular set is defined as the norm of the enclosing ball of the given $p$-norm, of minimal volume, that is centered in the origin.

**Input**

`H`

– hyperrectangular set`p`

– (optional, default:`Inf`

) norm

**Output**

A real number representing the norm.

**Algorithm**

Recall that the norm is defined as

\[‖ X ‖ = \max_{x ∈ X} ‖ x ‖_p = max_{x ∈ \text{vertices}(X)} ‖ x ‖_p.\]

The last equality holds because the optimum of a convex function over a polytope is attained at one of its vertices.

This implementation uses the fact that the maximum is attained in the vertex $c + \text{diag}(\text{sign}(c)) r$ for any $p$-norm. Hence it suffices to take the $p$-norm of this particular vertex. This statement is proved below. Note that, in particular, there is no need to compute the $p$-norm for *each* vertex, which can be very expensive.

If $X$ is a hyperrectangle and the $n$-dimensional vectors center and radius of the hyperrectangle are denoted $c$ and $r$ respectively, then reasoning on the $2^n$ vertices we have that:

\[\max_{x ∈ \text{vertices}(X)} ‖ x ‖_p = \max_{α_1, …, α_n ∈ \{-1, 1\}} (|c_1 + α_1 r_1|^p + ... + |c_n + α_n r_n|^p)^{1/p}.\]

The function $x ↦ x^p$, $p > 0$, is monotonically increasing and thus the maximum of each term $|c_i + α_i r_i|^p$ is given by $|c_i + \text{sign}(c_i) r_i|^p$ for each $i$. Hence, $x^* := \text{argmax}_{x ∈ X} ‖ x ‖_p$ is the vertex $c + \text{diag}(\text{sign}(c)) r$.

`IntervalArithmetic.radius`

— Function`radius(H::AbstractHyperrectangle, [p]::Real=Inf)`

Return the radius of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set`p`

– (optional, default:`Inf`

) norm

**Output**

A real number representing the radius.

**Notes**

The radius is defined as the radius of the enclosing ball of the given $p$-norm of minimal volume with the same center. It is the same for all corners of a hyperrectangular set.

`LazySets.σ`

— Method`σ(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`

.

`LazySets.ρ`

— Method`ρ(d::AbstractVector, H::AbstractHyperrectangle)`

Evaluate the support function of a hyperrectangular set in a given direction.

**Input**

`d`

– direction`H`

– hyperrectangular set

**Output**

The evaluation of the support function in the given direction.

`Base.:∈`

— Method`∈(x::AbstractVector, H::AbstractHyperrectangle)`

Check whether a given point is contained in a hyperrectangular set.

**Input**

`x`

– point/vector`H`

– hyperrectangular set

**Output**

`true`

iff $x ∈ H$.

**Algorithm**

Let $H$ be an $n$-dimensional hyperrectangular set, $c_i$ and $r_i$ be the center and radius, and $x_i$ be the vector $x$ in dimension $i$, respectively. Then $x ∈ H$ iff $|c_i - x_i| ≤ r_i$ for all $i=1,…,n$.

`LazySets.vertices_list`

— Method`vertices_list(H::AbstractHyperrectangle; kwargs...)`

Return the list of vertices of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

A list of vertices. Zeros in the radius are correctly handled, i.e., the result does not contain any duplicate vertices.

**Algorithm**

First we identify the dimensions where `H`

is flat, i.e., its radius is zero. We also compute the number of vertices that we have to create.

Next we create the vertices. We do this by enumerating all vectors `v`

of length `n`

(the dimension of `H`

) with entries `-1`

/`0`

/`1`

and construct the corresponding vertex as follows:

\[ \text{vertex}(v)(i) = \begin{cases} c(i) + r(i) & v(i) = 1 \\ c(i) & v(i) = 0 \\ c(i) - r(i) & v(i) = -1. \end{cases}\]

For enumerating the vectors `v`

, we modify the current `v`

from left to right by changing entries `-1`

to `1`

, skipping entries `0`

, and stopping at the first entry `1`

(but changing it to `-1`

). This way we only need to change the vertex in those dimensions where `v`

has changed, which usually is a smaller number than `n`

.

`LazySets.constraints_list`

— Method`constraints_list(H::AbstractHyperrectangle{N}) where {N}`

Return the list of constraints of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

A list of $2n$ linear constraints, where $n$ is the dimension of `H`

.

`constraints_list(U::Universe{N}) where {N}`

Return the list of constraints defining a universe.

**Input**

`U`

– universe

**Output**

The empty list of constraints, as the universe is unconstrained.

`LazySets.high`

— Method`high(H::AbstractHyperrectangle)`

Return the higher coordinates of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

A vector with the higher coordinates of the hyperrectangular set.

`LazySets.high`

— Method`high(H::AbstractHyperrectangle, i::Int)`

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

**Input**

`H`

– hyperrectangular set`i`

– dimension of interest

**Output**

The higher coordinate of the hyperrectangular set in the given dimension.

`LazySets.low`

— Method`low(H::AbstractHyperrectangle)`

Return the lower coordinates of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

A vector with the lower coordinates of the hyperrectangular set.

`LazySets.low`

— Method`low(H::AbstractHyperrectangle, i::Int)`

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

**Input**

`H`

– hyperrectangular set`i`

– dimension of interest

**Output**

The lower coordinate of the hyperrectangular set in the given dimension.

`Base.extrema`

— Method`extrema(H::AbstractHyperrectangle)`

Return the lower and higher coordinates of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

The lower and higher coordinates of the set.

**Notes**

The result is equivalent to `(low(H), high(H))`

.

`Base.extrema`

— Method`extrema(H::AbstractHyperrectangle, i::Int)`

Return the lower and higher coordinate of a hyperrectangular set in a given dimension.

**Input**

`H`

– hyperrectangular set`i`

– dimension of interest

**Output**

The lower and higher coordinate of the set in the given dimension.

**Notes**

The result is equivalent to `(low(H, i), high(H, i))`

.

`LazySets.isflat`

— Method`isflat(H::AbstractHyperrectangle)`

Check whether a hyperrectangular set is flat, i.e., whether its radius is zero in some dimension.

**Input**

`H`

– hyperrectangular set

**Output**

`true`

iff the hyperrectangular set is flat.

**Notes**

For robustness with respect to floating-point inputs, this function relies on the result of `isapproxzero`

when applied to the radius in some dimension. Hence this function depends on the absolute zero tolerance `ABSZTOL`

.

`Base.split`

— Method```
split(H::AbstractHyperrectangle{N},
num_blocks::AbstractVector{Int}) where {N}
```

Partition a hyperrectangular set into uniform sub-hyperrectangles.

**Input**

`H`

– hyperrectangular set`num_blocks`

– number of blocks in the partition for each dimension

**Output**

A list of `Hyperrectangle`

s.

`LazySets.generators`

— Method`generators(H::AbstractHyperrectangle)`

Return an iterator over the generators of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

An iterator over the generators of `H`

.

`LazySets.genmat`

— Methodgenmat(H::AbstractHyperrectangle)

Return the generator matrix of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

A matrix where each column represents one generator of `H`

.

`LazySets.ngens`

— Method`ngens(H::AbstractHyperrectangle{N}) where {N}`

Return the number of generators of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

The number of generators.

**Algorithm**

A hyperrectangular set has one generator for each non-flat dimension.

`ReachabilityBase.Arrays.rectify`

— Method`rectify(H::AbstractHyperrectangle)`

Concrete rectification of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

The `Hyperrectangle`

that corresponds to the rectification of `H`

.

`LazySets.volume`

— Method`volume(H::AbstractHyperrectangle)`

Return the volume of a hyperrectangular set.

**Input**

`H`

– hyperrectangular set

**Output**

The volume of $H$.

**Algorithm**

The volume of the $n$-dimensional hyperrectangle $H$ with radius vector $r$ is $2ⁿ ∏ᵢ rᵢ$ where $rᵢ$ denotes the $i$-th component of $r$.

`ReachabilityBase.Arrays.distance`

— Method```
distance(x::AbstractVector, H::AbstractHyperrectangle{N};
[p]::Real=N(2)) where {N}
```

Compute the distance between a point `x`

and a hyperrectangular set `H`

with respect to the given `p`

-norm.

**Input**

`x`

– point/vector`H`

– hyperrectangular set

**Output**

A scalar representing the distance between point `x`

and hyperrectangle `H`

.

### Implementations

Concrete set representations:

Lazy set operations:

## Singletons (AbstractSingleton)

A singleton is a special hyperrectangle consisting of only one point.

`LazySets.AbstractSingleton`

— Type`AbstractSingleton{N} <: AbstractHyperrectangle{N}`

Abstract type for sets with a single value.

**Notes**

Every concrete `AbstractSingleton`

must define the following function:

`element(::AbstractSingleton)`

– return the single element

```
julia> subtypes(AbstractSingleton)
2-element Vector{Any}:
Singleton
ZeroSet
```

This interface defines the following functions:

`LazySets.σ`

— Method`σ(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.

`LazySets.ρ`

— Method`ρ(d::AbstractVector, S::AbstractSingleton)`

Evaluate the support function of a set with a single value in a given direction.

**Input**

`d`

– direction`S`

– set with a single value

**Output**

The support value in the given direction.

`Base.:∈`

— Method`∈(x::AbstractVector, S::AbstractSingleton)`

Check whether a given point is contained in a set with a single value.

**Input**

`x`

– point/vector`S`

– set with a single value

**Output**

`true`

iff $x ∈ S$.

**Notes**

This implementation performs an approximate comparison to account for imprecision in floating-point computations.

`LazySets.center`

— Method`center(S::AbstractSingleton)`

Return the center of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

The center of the set.

`LazySets.center`

— Method`center(S::AbstractSingleton, i::Int)`

Return the center of a set with a single value in a given dimension.

**Input**

`S`

– set with a single value`i`

– dimension of interest

**Output**

The `i`

-th entry of the center of the set.

`LazySets.element`

— Method`element(S::AbstractSingleton, i::Int)`

Return the i-th entry of the element of a set with a single value.

**Input**

`S`

– set with a single value`i`

– dimension of interest

**Output**

The i-th entry of the element.

`LazySets.vertices`

— Method`vertices(S::AbstractSingleton{N}) where {N}`

Construct an iterator over the vertices of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

An iterator with a single value.

`LazySets.vertices_list`

— Method`vertices_list(S::AbstractSingleton)`

Return the list of vertices of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

A list containing only a single vertex.

`LazySets.radius_hyperrectangle`

— Method`radius_hyperrectangle(S::AbstractSingleton{N}) where {N}`

Return the box radius of a set with a single value in every dimension.

**Input**

`S`

– set with a single value

**Output**

The zero vector.

`LazySets.radius_hyperrectangle`

— Method`radius_hyperrectangle(S::AbstractSingleton{N}, i::Int) where {N}`

Return the box radius of a set with a single value in a given dimension.

**Input**

`S`

– set with a single value`i`

– dimension of interest

**Output**

Zero.

`LazySets.high`

— Method`high(S::AbstractSingleton)`

Return the higher coordinates of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

A vector with the higher coordinates.

`LazySets.high`

— Method`high(S::AbstractSingleton, i::Int)`

Return the higher coordinate of a set with a single value in the given dimension.

**Input**

`S`

– set with a single value`i`

– dimension of interest

**Output**

The higher coordinate in the given dimension.

`LazySets.low`

— Method`low(S::AbstractSingleton)`

Return the lower coordinates of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

A vector with the lower coordinates.

`LazySets.low`

— Method`low(S::AbstractSingleton, i::Int)`

Return the lower coordinate of a set with a single value in the given dimension.

**Input**

`S`

– set with a single value`i`

– dimension of interest

**Output**

The lower coordinate in the given dimension.

`LazySets.generators`

— Method`generators(S::AbstractSingleton{N}) where {N}`

Return an (empty) iterator over the generators of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

An empty iterator.

`LazySets.genmat`

— Methodgenmat(S::AbstractSingleton{N}) where {N}

Return the (empty) generator matrix of a set with a single value.

**Input**

`S`

– set with a single value

**Output**

A matrix with no columns representing the generators of `S`

.

`LazySets.ngens`

— Method`ngens(S::AbstractSingleton)`

Return the number of generators of a set with a single value.

**Input**

`H`

– set with a single value

**Output**

The number of generators, which is $0$.

`LazySets.plot_recipe`

— Method`plot_recipe(S::AbstractSingleton{N}, [ε]=zero(N)) where {N}`

Convert a singleton to a pair `(x, y)`

of points for plotting.

**Input**

`S`

– singleton`ε`

– (optional, default:`0`

) ignored, used for dispatch

**Output**

A pair `(x, y)`

of one point that can be plotted.

`RecipesBase.apply_recipe`

— Method`plot_singleton(S::AbstractSingleton{N}, [ε]::Real=zero(N); ...) where {N}`

Plot a singleton.

**Input**

`S`

– singleton`ε`

– (optional, default:`0`

) ignored, used for dispatch

**Examples**

`julia> plot(Singleton([0.5, 1.0]))`

### Implementations

## Affine maps (AbstractAffineMap)

An affine map consists of a linear map and a translation.

`LazySets.AbstractAffineMap`

— Type`AbstractAffineMap{N, S<:LazySet{N}} <: LazySet{N}`

Abstract type for affine maps.

**Notes**

See `AffineMap`

for a standard implementation of this interface.

Every concrete `AbstractAffineMap`

must define the following methods:

`matrix(::AbstractAffineMap)`

– return the linear map`vector(::AbstractAffineMap)`

– return the affine translation vector`set(::AbstractAffineMap)`

– return the set that the map is applied to

The subtypes of `AbstractAffineMap`

:

```
julia> subtypes(AbstractAffineMap)
7-element Vector{Any}:
AffineMap
ExponentialMap
ExponentialProjectionMap
InverseLinearMap
LinearMap
ResetMap
Translation
```

This interface defines the following functions:

`LazySets.dim`

— Method`dim(am::AbstractAffineMap)`

Return the dimension of an affine map.

**Input**

`am`

– affine map

**Output**

The ambient dimension of an affine map.

`LazySets.σ`

— Method`σ(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.

`LazySets.ρ`

— Method`ρ(d::AbstractVector, am::AbstractAffineMap)`

Evaluate the support function of an affine map.

**Input**

`d`

– direction`am`

– affine map

**Output**

The evaluation of the support function in the given direction.

`LazySets.an_element`

— Method`an_element(am::AbstractAffineMap)`

Return some element of an affine map.

**Input**

`am`

– affine map

**Output**

An element of the affine map.

**Algorithm**

The implementation relies on the `an_element`

method of the wrapped set.

`Base.isempty`

— Method`isempty(am::AbstractAffineMap)`

Check whether an affine map is empty.

**Input**

`am`

– affine map

**Output**

`true`

iff the wrapped set is empty.

`LazySets.isbounded`

— Method`isbounded(am::AbstractAffineMap; cond_tol::Number=DEFAULT_COND_TOL)`

Check whether an affine map is bounded.

**Input**

`am`

– affine map`cond_tol`

– (optional) tolerance of matrix condition (used to check whether the matrix is invertible)

**Output**

`true`

iff the affine map is bounded.

**Algorithm**

We first check if the matrix is zero or the wrapped set is bounded. If not, we perform a sufficient check whether the matrix is invertible. If the matrix is invertible, then the map being bounded is equivalent to the wrapped set being bounded, and hence the map is unbounded. Otherwise, we check boundedness via `_isbounded_unit_dimensions`

.

`Base.:∈`

— Method`∈(x::AbstractVector, am::AbstractAffineMap)`

Check whether a given point is contained in the affine map of a convex set.

**Input**

`x`

– point/vector`am`

– affine map of a convex set

**Output**

`true`

iff $x ∈ am$.

**Algorithm**

Observe that $x ∈ M⋅S ⊕ v$ iff $M^{-1}⋅(x - v) ∈ S$. This implementation does not explicitly invert the matrix, which is why it also works for non-square matrices.

**Examples**

```
julia> am = AffineMap([2.0 0.0; 0.0 1.0], BallInf([1., 1.], 1.), [-1.0, -1.0]);
julia> [5.0, 1.0] ∈ am
false
julia> [3.0, 1.0] ∈ am
true
```

An example with a non-square matrix:

```
julia> B = BallInf(zeros(4), 1.);
julia> M = [1. 0 0 0; 0 1 0 0]/2;
julia> [0.5, 0.5] ∈ M*B
true
```

`LazySets.center`

— Method`center(am::AbstractAffineMap)`

Return the center of an affine map of a centrally-symmetric set.

**Input**

`cp`

– affine map of a centrally-symmetric set

**Output**

The center of the affine map.

**Algorithm**

The implementation relies on the `center`

method of the wrapped set.

`LazySets.vertices_list`

— Method`vertices_list(am::AbstractAffineMap; [apply_convex_hull]::Bool)`

Return the list of vertices of a (polytopic) affine map.

**Input**

`am`

– affine map of a polytopic set`apply_convex_hull`

– (optional, default:`true`

) if`true`

, apply the convex hull operation to the list of vertices transformed by the affine map

**Output**

A list of vertices.

**Algorithm**

This implementation computes all vertices of `X`

, then transforms them through the affine map, i.e., `x ↦ M*x + v`

for each vertex `x`

of `X`

. By default, the convex-hull operation is taken before returning this list. For dimensions three or higher, this operation relies on the functionality through the concrete polyhedra library `Polyhedra.jl`

.

If you are not interested in taking the convex hull of the resulting vertices under the affine map, pass `apply_convex_hull=false`

as a keyword argument.

Note that we assume that the underlying set `X`

is polytopic, either concretely or lazily, i.e., the function `vertices_list`

should be applicable.

`LazySets.constraints_list`

— Method`constraints_list(am::AbstractAffineMap)`

Return the list of constraints of a (polyhedral) affine map.

**Input**

`am`

– affine map of a polyhedral set

**Output**

The list of constraints of the affine map.

**Notes**

We assume that the underlying set `X`

is polyhedral, i.e., offers a method `constraints_list(X)`

.

**Algorithm**

This implementation uses the method to compute the list of constraints of the translation of a lazy linear map.

`LazySets.linear_map`

— Method`linear_map(M::AbstractMatrix, am::AbstractAffineMap)`

Return the linear map of a lazy affine map.

**Input**

`M`

– matrix`am`

– affine map

**Output**

A set corresponding to the linear map of the lazy affine map of a set.

### Implementations

- Affine map (AffineMap)
- Exponential map (ExponentialMap)
- Linear map (LinearMap)
- Reset map (ResetMap)
- Translation

## Star sets (AbstractStar)

`LazySets.AbstractStar`

— Type`AbstractStar{N} <: LazySet{N}`

Abstract supertype for all star set types.

**Notes**

A set $X$ is star-like (also known as generalized star) if it can be represented by a center $x₀ ∈ \mathbb{R}^n$ and $m$ vectors $v₁, …, vₘ$ forming the basis, and a predicate $P : \mathbb{R}^n → \{⊤, ⊥\}$ such that

\[ X = \{x ∈ \mathbb{R}^n : x = x₀ + \sum_{i=1}^m α_i v_i,~~\textrm{s.t. } P(α) = ⊤ \}.\]

### Implementations

## Polynomial zonotope sets (AbstractPolynomialZonotope)

`LazySets.AbstractPolynomialZonotope`

— Type`AbstractPolynomialZonotope{N} <: LazySet{N}`

Abstract type for polynomial zonotope sets.

**Notes**

Polynomial zonotopes are in general non-convex. They are always bounded.

```
julia> subtypes(AbstractPolynomialZonotope)
3-element Vector{Any}:
DensePolynomialZonotope
SimpleSparsePolynomialZonotope
SparsePolynomialZonotope
```