LGG09{N, AM, TN<:AbstractDirections} <: AbstractContinuousPost

Implementation of Girard - Le Guernic algorithm for reachability analysis using support functions.


  • δ – step-size of the discretization
  • approx_model – (optional, default: Forward) approximation model; see Notes below for possible options
  • template – (alias: dirs) struct that holds the directions (either lazily or concretely) for each support function evaluation defining the template
  • vars – (optional, default: all variables are computed) an integer or a vector of integers specifying the variables of interest to automatically construct a template using canonical directions; requires that n (or dim) is specified as well
  • static – (optional, default: false) if true, use statically sized arrays
  • threaded – (optional, default: false) if true, use multi-threading to compute different template directions in parallel
  • sparse – (optional, default: false) if true, convert the matrix exponential obtained after discretization to a sparse matrix
  • cache – (optional, default: true) if true, use a cache for intermediate computations in the set recurrence loop
  • vars – (optional, default: missing) used to specify the variables instead of passing the template


The type fields are:

  • N – number type of the step-size
  • AM – type of the approximation model
  • TN – type of the abstract directions that define the template

The flag threaded=true specifies that the support functions along different directions should be computed in parallel. See the section on Multi-threading for details on how to setup the number of threads.


This is an implementation of the algorithm from [LGG09].

These methods are described at length in the dissertation [LG09].

reach_homog_krylov_LGG09!(out, Ω₀::LazySet, Aᵀδ::AbstractMatrix,
                          ℓ::AbstractVector, NSTEPS;
                          hermitian=false, m=min(30, size(Aᵀδ, 1)), tol=1e-7)


We compute the sequence:

\[ ρ(ℓ, Ω₀), ρ(ℓ, Φ Ω₀), ρ(ℓ, Φ^2 Ω₀), ρ(ℓ, Φ^3 Ω₀), ...\]

Using Krylov subspace approximations to compute the action of Φ := exp(Aδ) over the direction ℓ.

The method is (see [1]):

out[1] <- ρ(ℓ, Ω₀)

out[2] <- ρ(ℓ, Φ Ω₀) = ρ(Φᵀ ℓ, Ω₀)

out[3] <- ρ(ℓ, Φ^2 Ω₀) = ρ((Φᵀ)^2 ℓ, Ω₀)

out[4] <- ρ(ℓ, Φ^3 Ω₀) = ρ((Φᵀ)^3 ℓ, Ω₀)

and so on.


[1] Reach Set Approximation through Decomposition with Low-dimensional Sets and High-dimensional Matrices. Sergiy Bogomolov, Marcelo Forets, Goran Frehse, Frédéric Viry, Andreas Podelski and Christian Schilling (2018) HSCC'18 Proceedings of the 21st International Conference on Hybrid Systems: Computation and Control: 41–50.