Approximations

Module Approximations.jl – polygonal approximation of convex sets through support vectors.

decompose(S::LazySet{N};
          [set_type]::Type{<:Union{HPolygon, Hyperrectangle, Interval}}=Hyperrectangle,
          [ε]::Real=Inf,
          [blocks]::AbstractVector{Int}=default_block_structure(S, set_type),
          [block_types]::Dict{Type{<:LazySet}, AbstractVector{<:AbstractVector{Int}}}(),
          [directions]::Union{Type{<:AbstractDirections}, Nothing}=nothing
         )::CartesianProductArray where {N<:Real}

Decompose a high-dimensional set into a Cartesian product of overapproximations of the projections over the specified subspaces.

Input

• S – set

• set_type – (optional, default: Hyperrectangle) type of set approximation for each subspace

• ε – (optional, default: Inf) error bound for polytopic approximation

• blocks – (optional, default: [2, …, 2] or [1, …, 1] if set_type is an interval) block structure - a vector with the size of each block

• block_types – (optional, default: Interval for 1D and Hyperrectangle for mD blocks) a mapping from set types to blocks

• directions – (optional, default: nothing) template direction type, or nothing

Output

A CartesianProductArray containing the low-dimensional approximated projections.

Algorithm

For each block a specific project method is called, dispatched on the set_type argument.

Notes

If directions is different from nothing, the template directions are used together with blocks. Otherwise, if block_types is given, the options set_type and blocks are ignored.

Examples

The decompose function supports different options, such as: supplying different dimensions for the decomposition, defining the target set of the decomposition, or specifying the degree of accuracy of the target decomposition. You can also choose to make the approximations in low dimensions using template directions. These options are exemplified below.

Different dimensions

By default, decompose returns a Cartesian product of 2D Hyperrectangle sets. For example:

julia> import LazySets.Approximations:decompose

julia> S = Ball2(zeros(4), 1.);

julia> array(decompose(S))
2-element Array{LazySets.LazySet{Float64},1}:
 LazySets.Hyperrectangle{Float64}([0.0, 0.0], [1.0, 1.0])
 LazySets.Hyperrectangle{Float64}([0.0, 0.0], [1.0, 1.0])

Other block sizes can be specified using the blocks option, which refers to each block size of the partition:

julia> array(decompose(S, blocks=[1, 3]))
2-element Array{LazySets.LazySet{Float64},1}:
 LazySets.Hyperrectangle{Float64}([0.0], [1.0])
 LazySets.Hyperrectangle{Float64}([0.0, 0.0, 0.0], [1.0, 1.0, 1.0])

julia> array(decompose(S, blocks=[4]))
1-element Array{LazySets.LazySet{Float64},1}:
 LazySets.Hyperrectangle{Float64}([0.0, 0.0, 0.0, 0.0], [1.0, 1.0, 1.0, 1.0])

Different set types

We can also decompose using polygons in constraint representation, through the set_type optional argument:

julia> all([ai isa HPolygon for ai in array(decompose(S, set_type=HPolygon))])
true

For decomposition into 1D subspaces, we can use Interval:

julia> all([ai isa Interval for ai in array(decompose(S, set_type=Interval))])
true

However, if you need to specify different set types for different blocks, the interface presented so far does not apply. In the paragraph Advanced different set types input we explain the input block_types, that can be used precisely for that purpose.

Refining the decomposition I: $ε$-close approximation

The $ε$ option can be used to refine, that is obtain a more accurate decomposition in those blocks where HPolygon types are used, and it relies on the iterative refinement algorithm provided in the Approximations module.

To illustrate this, consider the unit 4D ball in the 2-norm. Using smaller $ε$ implies a better precision, thus more constraints in each 2D decomposition:

julia> S = Ball2(zeros(4), 1.);

julia> d(ε, bi) = array(decompose(S, set_type=HPolygon, ε=ε))[bi]
d (generic function with 1 method)

julia> [length(constraints_list(d(ε, 1))) for ε in [Inf, 0.1, 0.01]]
3-element Array{Int64,1}:
  4
  8
 32

Refining the decomposition II: template polyhedra

Another way to refine the decomposition is using template polyhedra. The idea is to specify a set of template directions, and on each block, compute the polytopic overapproximation obtained by evaluating the support function of the given input set over the template directions.

For example, octagonal 2D approximations of the ball S are obtained with:

julia> B = decompose(S, directions=OctDirections);

julia> length(B.array) == 2 && all(dim(bi) == 2 for bi in B.array)
true

See template_directions.jl for the available template directions. Note that, in contrast to the polygonal $ε$-close approximation, this method can be applied for blocks of any size.

julia> B = decompose(S, directions=OctDirections, blocks=[4]);

julia> length(B.array) == 1 && dim(B.array[1]) == 4
true

Advanced different set types input

We can define different set types for different blocks, using the optional block_types input argument. It is a dictionary where the keys correspond to set types, and the values correspond to the blocks, namely the initial and final block indices should be given.

For example:

julia> S = Ball2(zeros(3), 1.);

julia> array(decompose(S, block_types=Dict(Interval=>[1:1], Hyperrectangle=>[2:3])))
2-element Array{LazySets.LazySet{Float64},1}:
 LazySets.Interval{Float64,IntervalArithmetic.Interval{Float64}}([-1, 1])
 LazySets.Hyperrectangle{Float64}([0.0, 0.0], [1.0, 1.0])

We can additionally pass ε, which is automatically used for each HPolygon type block.

julia> S = Ball2(zeros(8), 1.);

julia> bt = Dict(Interval=>[1:1], Hyperrectangle=>[2:4], HPolygon=>[5:6, 7:8]);

julia> [typeof(ai) for ai in array(decompose(S, block_types=bt, ε=0.01))]
4-element Array{DataType,1}:
 LazySets.Interval{Float64,IntervalArithmetic.Interval{Float64}}
 LazySets.Hyperrectangle{Float64}
 LazySets.HPolygon{Float64}
 LazySets.HPolygon{Float64}

default_block_structure(S::LazySet, set_type::Type{<:LazySet})::AbstractVector{Int}

Compute the default block structure.

Input

• S – set

• set_type – target set type

Output

A vector representing the block structure, such that:

• If the target set_type is an interval, the default is blocks of size 1.

• Otherwise, the default is blocks of size 2. Depending on the dimension, the last block has size 1 or 2.

project(S::LazySet{N},
        block::AbstractVector{Int},
        set_type::Type{<:LazySet},
        [n]::Int=dim(S),
        [ε]::Real=Inf
       )::LazySet{N} where {N<:Real}

Default implementation for projecting a high-dimensional set to a given set type with possible 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

• ε – (optional, default: Inf) ignored

Output

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

Algorithm

1. Project the set S with M⋅S, where M is the identity matrix in the block coordinates and zero otherwise.

2. Overapproximate the projected lazy set using overapproximate.

project(S::LazySet{N},
        block::AbstractVector{Int},
        set_type::Type{<:HPolygon},
        [n]::Int=dim(S),
        [ε]::Real=Inf
       )::HPolygon where {N<:Real}

Project a high-dimensional set to a two-dimensional polygon with a certified error bound.

Input

• S – set

• block – block structure - a vector with the dimensions of interest

• set_type – HPolygon - used for dispatch

• n – (optional, default: dim(S)) ambient dimension of the set S

• ε – (optional, default: Inf) error bound for polytopic approximation

Output

A HPolygon representing the epsilon-close approximation of the box approximation of the projection of S.

Notes

block must have length 2.

Algorithm

If ε < Inf, the algorithm proceeds as follows:

1. Project the set S with M⋅S, where M is the identity matrix in the block coordinates and zero otherwise.

2. Overapproximate the set with the given error bound ε.

If ε == Inf, the algorithm uses a box approximation.

project(S::LazySet{N},
        block::AbstractVector{Int},
        set_type::Type{<:Hyperrectangle},
        [n]::Int=dim(S),
        [ε]::Real=Inf
       )::Hyperrectangle where {N<:Real}

Project a high-dimensional set to a low-dimensional hyperrectangle.

Input

• S – set

• block – block structure - a vector with the dimensions of interest

• set_type – Hyperrectangle - used for dispatch

• n – (optional, default: dim(S)) ambient dimension of the set S

• ε – (optional, default: Inf) - used for dispatch, ignored

Output

The box approximation of the projection of S.

project(S::LazySet{N},
        block::AbstractVector{Int},
        directions::AbstractDirections{N},
        n::Int
       )::HPolytope where {N<:Real}

Project a high-dimensional set to a low-dimensional set using template directions.

Input

• S – set

• block – block structure - a vector with the dimensions of interest

• directions – template directions

• n – (optional, default: dim(S)) ambient dimension of the set S

Output

The template direction approximation of the projection of S.

overapproximate(X::S, ::Type{S}) where {S <: LazySet}

Overapproximating a set of type S with type S is a no-op.

Input

• X – set

• Type{S} – set type

Output

The input set.

overapproximate(S::LazySet{N},
                ::Type{<:HPolygon},
                [ε]::Real=Inf)::HPolygon where {N<:Real}

Return an approximation of a given 2D convex set. If no error tolerance is given, or is Inf, the result is a box-shaped polygon. Otherwise the result is an ε-close approximation as a polygon.

Input

• S – convex set, assumed to be two-dimensional

• HPolygon for dispatch

• ε – (optional, default: Inf) error bound

Output

A polygon in constraint representation.

overapproximate(S::LazySet, ε::Real)::HPolygon

Alias for overapproximate(S, HPolygon, ε).

overapproximate(S::LazySet, Type{<:Hyperrectangle})::Hyperrectangle

Return an approximation of a given set as a hyperrectangle.

Input

• S – set

• Hyperrectangle for dispatch

Output

A hyperrectangle.

overapproximate(S::LazySet)::Hyperrectangle

Alias for overapproximate(S, Hyperrectangle).

overapproximate(S::ConvexHull{N, Zonotope{N}, Zonotope{N}},
                 ::Type{<:Zonotope})::Zonotope where {N<:Real}

Overapproximate the convex hull of two zonotopes.

Input

• S – convex hull of two zonotopes of the same order

• Zonotope for dispatch

Algorithm

This function implements the method proposed in Reachability of Uncertain Linear Systems Using Zonotopes, A. Girard, HSCC 2005. The convex hull of two zonotopes of the same order, that we write $Z_j = ⟨c^{(j)}, g^{(j)}_1, …, g^{(j)}_p⟩$ for $j = 1, 2$, can be overapproximated as follows:

It should be noted that the output zonotope is not necessarily the minimal enclosing zonotope, which is in general expensive in high dimensions. This is further investigated in: Zonotopes as bounding volumes, L. J. Guibas et al, Proc. of Symposium on Discrete Algorithms, pp. 803-812.

overapproximate(X::LazySet, dir::AbstractDirections)::HPolytope

Overapproximating a set with template directions.

Input

• X – set

• dir – direction representation

Output

A HPolytope overapproximating the set X with the directions from dir.

overapproximate(::LazySet{N}, ::Type{LazySets.Interval}) where {N<:Real}

Return the overapproximation of a real unidimensional set with an interval.

Input

• S – one-dimensional set

• Interval for dispatch

Output

An interval.

ballinf_approximation(S::LazySet{N})::BallInf{N} where {N<:Real}

Overapproximate a convex set by a tight ball in the infinity norm.

Input

• S – convex set

Output

A tight ball in the infinity norm.

Algorithm

The center and radius of the box are obtained by evaluating the support function of the given convex set along the canonical directions.

box_approximation(S::LazySet)::Hyperrectangle

Overapproximate a convex set by a tight hyperrectangle.

Input

• S – convex set

Output

A tight hyperrectangle.

Algorithm

The center of the hyperrectangle is obtained by averaging the support function of the given set in the canonical directions, and the lengths of the sides can be recovered from the distance among support functions in the same directions.

interval_hull

Alias for box_approximation.

box_approximation_symmetric(S::LazySet{N})::Hyperrectangle{N} where {N<:Real}

Overapproximate a convex set by a tight hyperrectangle centered in the origin.

Input

• S – convex set

Output

A tight hyperrectangle centered in the origin.

Algorithm

The center of the box is the origin, and the radius is obtained by computing the maximum value of the support function evaluated at the canonical directions.

symmetric_interval_hull

Alias for box_approximation_symmetric.

box_approximation_helper(S::LazySet)

Common code of box_approximation and box_approximation_symmetric.

Input

• S – convex set

Output

A tuple containing the data that is needed to construct a tightly overapproximating hyperrectangle.

• c – center

• r – radius

Algorithm

The center of the hyperrectangle is obtained by averaging the support function of the given convex set in the canonical directions. The lengths of the sides can be recovered from the distance among support functions in the same directions.

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

Return the norm of a convex 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 – convex set

• p – (optional, default: Inf) norm

Output

A real number representing the norm.

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

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

Input

• S – convex set

• p – (optional, default: Inf) norm

Output

A real number representing the radius.

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

Return the diameter of a convex 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 – convex set

• p – (optional, default: Inf) norm

Output

A real number representing the diameter.

LocalApproximation{N<:Real}

Type that represents a local approximation in 2D.

Fields

• p1 – first inner point

• d1 – first direction

• p2 – second inner point

• d2 – second direction

• q – intersection of the lines l1 ⟂ d1 at p1 and l2 ⟂ d2 at p2

• refinable – states if this approximation is refinable

• err – error upper bound

Notes

The criteria for refinable are determined in the method new_approx.

PolygonalOverapproximation{N<:Real}

Type that represents the polygonal approximation of a convex set.

Fields

• S – convex set

• approx_list – vector of local approximations

new_approx(S::LazySet, p1::Vector{N}, d1::Vector{N}, p2::Vector{N},
           d2::Vector{N}) where {N<:Real}

Input

• S – convex set

• p1 – first inner point

• d1 – first direction

• p2 – second inner point

• d2 – second direction

Output

A local approximation of S in the given directions.

addapproximation!(Ω::PolygonalOverapproximation,
    p1::Vector{N}, d1::Vector{N}, p2::Vector{N}, d2::Vector{N}) where {N<:Real}

Input

• Ω – polygonal overapproximation of a convex set

• p1 – first inner point

• d1 – first direction

• p2 – second inner point

• d2 – second direction

Output

The list of local approximations in Ω of the set Ω.S is updated in-place and the new approximation is returned by this function.

refine(Ω::PolygonalOverapproximation, i::Int)::Tuple{LocalApproximation, LocalApproximation}

Refine a given local approximation of the polygonal approximation of a convex set, by splitting along the normal direction to the approximation.

Input

• Ω – polygonal overapproximation of a convex set

• i – integer index for the local approximation to be refined

Output

The tuple consisting of the refined right and left local approximations.

tohrep(Ω::PolygonalOverapproximation{N})::AbstractHPolygon where {N<:Real}

Convert a polygonal overapproximation into a concrete polygon.

Input

• Ω – polygonal overapproximation of a convex set

Output

A polygon in constraint representation.

approximate(S::LazySet{N},
            ε::N)::PolygonalOverapproximation{N} where {N<:Real}

Return an ε-close approximation of the given 2D convex set (in terms of Hausdorff distance) as an inner and an outer approximation composed by sorted local Approximation2D.

Input

• S – 2D convex set

• ε – error bound

Output

An ε-close approximation of the given 2D convex set.

AbstractDirections{N}

Abstract type for template direction representations.

Notes

All subtypes should implement the standard iterator methods from Base and the function dim(d<:AbstractDirections)::Int.

BoxDirections{N} <: AbstractDirections{N}

Box direction representation.

Fields

• n – dimension

OctDirections{N} <: AbstractDirections{N}

Octagon direction representation.

Fields

• n – dimension

Notes

Octagon directions consist of all vectors that are zero almost everywhere except in two dimensions $i$, $j$ (possibly $i = j$) where it is $±1$.

BoxDiagDirections{N} <: AbstractDirections{N}

Box-diagonal direction representation.

Fields

• n – dimension

Notes

Box-diagonal directions can be seen as the union of diagonal directions (all entries are ±1) and box directions (one entry is ±1, all other entries are 0). The iterator first enumerates all diagonal directions, and then all box directions.