Merge pull request #245 from JuliaReach/schillic/241

Fix some piracies

src/IntervalMatrices.jl

@@ -38,29 +38,17 @@ if vIA >= v"0.22"
     import Base: intersect
     export ±, midpoint_radius
 
-    function ±(x::Number, y::Number)
-        return x + interval(-y, y)
-    end
-
     function Base.:(==)(A::AbstractMatrix{<:Interval}, B::AbstractMatrix{<:Interval})
         return size(A) == size(B) && all(map((a, b) -> isequal_interval(a, b), A, B))
     end
-
-    function intersect(a::Interval{T}, b::Interval{T}) where {T}
-        lo = max(inf(a), inf(b))
-        hi = min(sup(a), sup(b))
-        if lo > hi
-            return emptyinterval(T)
-        end
-        return interval(lo, hi)
-    end
-    intersect(a::Interval{T}, b::Interval{S}) where {T,S} = intersect(promote(a, b)...)
-
-    Base.in(x::Number, y::Interval) = inf(y) <= x <= sup(y)
 else
     import IntervalArithmetic: ±, midpoint_radius
 
     issubset_interval(x::Interval, y::Interval) = issubset(x, y)
+
+    in_interval(x::Number, y::Interval) = inf(y) <= x <= sup(y)
+
+    intersect_interval(a::Interval, b::Interval) = intersect(a, b)
 end
 
 # ================================

src/exponential.jl

@@ -96,7 +96,7 @@ function _exp_remainder(A::IntervalMatrix{T}, t, p; n=checksquare(A)) where {T}
     end
     M = exp(C * t)
     Y = M - Q
-    Γ = IntervalMatrix(fill(zero(T) ± one(T), (n, n)))
+    Γ = IntervalMatrix(fill(interval(-one(T), one(T)), (n, n)))
     E = Γ * Y
     return E
 end
@@ -111,7 +111,7 @@ function _exp_remainder_series(A::IntervalMatrix{T}, t, p; n=checksquare(A)) whe
     @assert c < 1 "the remainder of the matrix exponential could not be " *
                   "computed because a convergence condition is not satisfied: $c ≥ 1 " *
                   "but it should be smaller than 1; try choosing a larger order"
-    Γ = IntervalMatrix(fill(zero(T) ± one(T), (n, n)))
+    Γ = IntervalMatrix(fill(interval(-one(T), one(T)), (n, n)))
     return Γ * ((nA * t)^(p + 1) * (1 / factorial(p + 1) * 1 / (1 - c)))
 end
 

src/matrix.jl

@@ -18,10 +18,10 @@ parametrized in the number field, the interval type, and the matrix type.
 ### Examples
 
 ```jldoctest
-julia> A = IntervalMatrix([-0.9±0.1 0±0; 0±0 -0.9±0.1])
+julia> A = IntervalMatrix([-1 .. -0.8 0 .. 0; 0 .. 0 -1 .. -0.8])
 2×2 IntervalMatrix{Float64, Interval{Float64}, Matrix{Interval{Float64}}}:
- [-1.00001, -0.7999999]   [0.0, 0.0]
-  [0.0, 0.0]             [-1.00001, -0.7999999]
+ [-1.0, -0.7999999]   [0.0, 0.0]
+  [0.0, 0.0]         [-1.0, -0.7999999]
 ```
 
 An interval matrix proportional to the identity matrix can be built using the
@@ -176,7 +176,7 @@ function ±(C::MT, S::MT) where {T,MT<:AbstractMatrix{T}}
                                               "radii matrix should match, but they are $(size(C)) " *
                                               "and $(size(S)) respectively"))
 
-    return IntervalMatrix(map((x, y) -> x ± y, C, S))
+    return IntervalMatrix(map((x, y) -> interval(x - y, x + y), C, S))
 end
 
 for op in (:Adjoint, :Bidiagonal, :Diagonal, :Hermitian,

src/operations/setops.jl

@@ -24,7 +24,7 @@ function ∈(M::AbstractMatrix, A::AbstractIntervalMatrix)
     m, n = size(A)
     @inbounds for j in 1:n
         for i in 1:m
-            if M[i, j] ∉ A[i, j]
+            if !in_interval(M[i, j], A[i, j])
                 return false
             end
         end
@@ -78,7 +78,7 @@ function ∩(A::IntervalMatrix, B::IntervalMatrix)
     @assert size(A) == size(B) "incompatible matrix sizes (A: $(size(A)), B: " *
                                "$(size(B)))"
 
-    return IntervalMatrix(map((x, y) -> x ∩ y, A, B))
+    return IntervalMatrix(map((x, y) -> intersect_interval(x, y), A, B))
 end
 
 """