Repository: JuliaReach/ReachabilityAnalysis.jl (julia)

Issue:
Add "discrete time" approximation model
see ecs. (14) in the [BFFPSV18 paper](https://arxiv.org/pdf/1801.09526.pdf)

Pull request description:
Closes #208.

Patch:
diff --git a/src/Continuous/discretization.jl b/src/Continuous/discretization.jl
index e37de0db08..068be4952b 100644
--- a/src/Continuous/discretization.jl
+++ b/src/Continuous/discretization.jl
@@ -32,8 +32,15 @@ function Backward(; exp::Symbol=:base, setops::Symbol=:lazy, sih::Symbol=:concre
     Backward(Val(exp), Val(setops), Val(sih))
 end
 
-struct Discrete <: AbstractApproximationModel
-#
+# no bloating or "discrete time" approximation model ref Eqs. (14) in [[BFFPSV18]]
+struct NoBloating{EM, SO} <: AbstractApproximationModel
+    exp::EM
+    setops::SO
+end
+
+# convenience constructor using symbols
+function NoBloating(; exp::Symbol=:base, setops::Symbol=:lazy)
+    NoBloating(Val(exp), Val(setops))
 end
 
 struct CorrectionHull{EM} <: AbstractApproximationModel
@@ -178,3 +185,37 @@ function discretize(ivp::IVP{<:CLCS, <:LazySet}, δ, alg::CorrectionHull)
     ivp_discr = ConstrainedLinearDiscreteSystem(Φ, X)
     return InitialValueProblem(ivp_discr, Ω0)
 end
+
+# ============================================================
+# NoBloating Approximation
+# ============================================================
+
+# homogeneous case
+function discretize(ivp::IVP{<:CLCS, <:LazySet}, δ, alg::NoBloating)
+    A = state_matrix(ivp)
+    X0 = initial_state(ivp)
+
+    Φ = _exp(A, δ, alg.exp)
+
+    Ω0 = copy(X0)
+    X = stateset(ivp)
+    Sdiscr = ConstrainedLinearDiscreteSystem(Φ, X)
+    return InitialValueProblem(Sdiscr, Ω0)
+end
+
+# inhomogeneous case
+function discretize(ivp::IVP{<:CLCCS, <:LazySet}, δ, alg::NoBloating)
+    A = state_matrix(ivp)
+    X0 = initial_state(ivp)
+
+    Φ = _exp(A, δ, alg.exp)
+    M = Φ₁(A, δ, alg.exp)
+    U = next_set(inputset(ivp), 1) # inputset(ivp)
+    Ud = M * U # TODO assert alg.setops  / alg.inputsetops are lazy
+    In = IdentityMultiple(one(eltype(A)), size(A, 1))
+
+    Ω0 = copy(X0)
+    X = stateset(ivp)
+    Sdiscr = ConstrainedLinearControlDiscreteSystem(Φ, In, X, Ud)
+    return InitialValueProblem(Sdiscr, Ω0)
+end
diff --git a/src/Initialization/exports.jl b/src/Initialization/exports.jl
index 811eda4604..7d46561020 100644
--- a/src/Initialization/exports.jl
+++ b/src/Initialization/exports.jl
@@ -18,6 +18,7 @@ export
     Backward,
     Discrete,
     CorrectionHull,
+    NoBloating,
 
 # Flowpipes
     flowpipe,


Test patch:
diff --git a/test/algorithms/GLGM06.jl b/test/algorithms/GLGM06.jl
index a071e148aa..32d957902a 100644
--- a/test/algorithms/GLGM06.jl
+++ b/test/algorithms/GLGM06.jl
@@ -1,16 +1,16 @@
 @testset "GLGM06 algorithm" begin
 
     # one-dimensional
-    prob, tspan = exponential_1d()
-    sol = solve(prob, tspan=tspan, GLGM06(δ=0.01))
+    prob, _ = exponential_1d()
+    sol = solve(prob, T=5.0, GLGM06(δ=0.01))
     @test isa(sol.alg, GLGM06)
     @test setrep(sol) <: Zonotope
     @test setrep(sol) == Zonotope{Float64,Array{Float64,1},Array{Float64,2}}
     @test dim(sol) == 1
 
     # higher-dimensional homogeneous
-    prob, tspan = linear5D_homog()
-    sol = solve(prob, tspan=tspan, GLGM06(δ=0.01))
+    prob, _ = linear5D_homog()
+    sol = solve(prob, T=5.0, GLGM06(δ=0.01))
     @test dim(sol) == 5
 
     # static option TODO
@@ -21,4 +21,15 @@
     # TODO higher-dimensional with input
     #prob, tspan = linear5D()
     #sol = solve(prob, tspan=tspan, GLGM06(δ=0.01))
+
+    # use approx model for "discrete-time" reachability
+    sol = solve(prob, T=5.0, GLGM06(δ=0.01, approx_model=NoBloating()))
+
+    # higher-dimensional inhomogeneous
+    prob, _ = linear5D_homog()
+    sol = solve(prob, T=5.0, GLGM06(δ=0.01))
+    @test dim(sol) == 5
+
+    # use approx model for "discrete-time" reachability
+    sol = solve(prob, T=5.0, GLGM06(δ=0.01, approx_model=NoBloating()))
 end