Merge pull request #286 from JuliaReach/schillic/exp_remainder

Fix `_exp_remainder`

Author: schillic

src/exponential.jl

@@ -92,10 +92,12 @@ function _exp_remainder(A::IntervalMatrix{T}, t, p; n=checksquare(A)) where {T}
         Cⁱ *= C
     end
     M = exp(C * t)
-    Y = M - Q
-    Γ = IntervalMatrix(fill(interval(-one(T), one(T)), (n, n)))
-    E = Γ * Y
-    return E
+    E = Matrix{Interval{T}}(undef, n, n)
+    @inbounds for i in eachindex(M)
+        y = M[i] - Q[i]
+        E[i] = (y >= zero(T)) ? interval(-y, y) : interval(y, -y)
+    end
+    return IntervalMatrix(E)
 end
 
 # Estimates the sum of the series in the matrix exponential. See

test/exponential.jl

@@ -5,6 +5,7 @@ using IntervalMatrices: TaylorOverapproximation,
                         horner,
                         scale_and_square,
                         _truncated_exponential_series,
+                        _exp_remainder,
                         _exp_remainder_series
 
 @testset "Interval matrix exponential" begin
@@ -44,6 +45,14 @@ using IntervalMatrices: TaylorOverapproximation,
     @test underapp1 == exp(M; alg=TaylorUnderapproximation(4))
 end
 
+@testset "Interval matrix exponential remainder" begin
+    M = IntervalMatrix([interval(1.5, 2) interval(0); interval(0) interval(0.5, 1)])
+    E = _exp_remainder(M, 1, 1)
+    e1 = exp(2) - 3
+    e2 = exp(1) - 2
+    @test E ⩵ IntervalMatrix([interval(-e1, e1) interval(0); interval(0) interval(-e2, e2)])
+end
+
 @testset "Interval matrix correction terms" begin
     m = IntervalMatrix([interval(-1.1, 0.9) interval(-4.1, -3.9);
                         interval(3.9, 4.1) interval(-1.1, 0.9)])