Rewrite IntervalBox to contain an SVector of Intervals (#152)

* Arithmetic operations and functions on IntervalBox return IntervalBox

* Add some tests of functions on IntervalBoxes

* More tests...

* Rewrite IntervalBox to contain an SVector

* Modify arithmetic operations for new IntervalBox

* Fix display of IntervalBox

* Add missing functions and fix tests

* Rename arithmetic_functions -> arithmetic

* Fix broadcasting

* Fix tests. So far no iteration for IntervalBox

* Add setindex

* Add propagate_inbounds

* Add iteration for IntervalBox

* Add tests for getindex and setindex

* Typo

* Fix broadcasting and add tests

* Fix broadcasting, take 2

src/IntervalArithmetic.jl

@@ -29,6 +29,11 @@ import Base:
     isinteger, setdiff,
     parse
 
+import Base:  # for IntervalBox
+    broadcast, dot, length,
+    getindex, setindex,
+    start, next, done, eltype
+
 export
     AbstractInterval, Interval,
     interval,

src/display.jl

@@ -245,10 +245,10 @@ function representation(X::IntervalBox, format=nothing)
     end
 
     if display_params.format == :full
-        return string("IntervalBox(", join(X, ", "), ")")
+        return string("IntervalBox(", join(X.v, ", "), ")")
 
     else
-        return join(X, " × ")
+        return join(X.v, " × ")
     end
 
 end

src/multidim/arithmetic.jl

@@ -0,0 +1,25 @@
+# This file is part of the IntervalArithmetic.jl package; MIT licensed
+
++(a::IntervalBox, b::IntervalBox) = IntervalBox( a.v .+ b.v )
++(a::IntervalBox, b::Real) = IntervalBox( a.v .+ b )
++(a::Real, b::IntervalBox) = IntervalBox( a .+ b.v )
+
+-(a::IntervalBox, b::IntervalBox) = IntervalBox( a.v .- b.v )
+-(a::IntervalBox, b::Real) = IntervalBox( a.v .- b )
+-(a::Real, b::IntervalBox) = IntervalBox( a .- b.v )
+-(a::IntervalBox) = IntervalBox( .- a.v )
+
+*(a::IntervalBox, b::Real) = IntervalBox( a.v .* b )
+*(a::Real, b::IntervalBox) = IntervalBox( a .* b.v )
+
+/(a::IntervalBox, b::Real) = IntervalBox( a.v ./ b )
+
+
+# broadcasting:
+
+# wrap decides whether to wrap the result in an IntervalBox or not, based on the return type
+wrap(v::SVector{N,T} where {N,T<:Interval}) = IntervalBox(v)
+wrap(v) = v
+
+Base.broadcast(f, X::IntervalBox) = wrap(f.(X.v))
+Base.broadcast(f, X::IntervalBox, Y::IntervalBox) = wrap(f.(X.v, Y.v))

src/multidim/intervalbox.jl

@@ -1,55 +1,76 @@
 # This file is part of the IntervalArithmetic.jl package; MIT licensed
 
 """An `IntervalBox` is an `N`-dimensional rectangular box, given
-by a Cartesian product of `N` `Interval`s.
+by a Cartesian product of a vector of `N` `Interval`s.
 """
-struct IntervalBox{N,T} <: StaticVector{N, Interval{T}}
-    data::NTuple{N,Interval{T}}
+struct IntervalBox{N,T}
+    v::SVector{N, Interval{T}}
 end
 
-# StaticArrays.Size{N,T}(::Type{IntervalBox{N,T}}) = StaticArrays.Size(N) # @pure not needed, I think...
-Base.@propagate_inbounds Base.getindex(a::IntervalBox, i::Int) = a.data[i]
+# IntervalBox(x::Interval) = IntervalBox( SVector(x) )  # single interval treated as tuple with one element
+
+IntervalBox(x::Interval...) = IntervalBox(SVector(x))
+IntervalBox(x::Tuple{T}) where {T<:Interval} = IntervalBox(SVector(x))
+
+Base.@propagate_inbounds Base.getindex(X::IntervalBox, i) = X.v[i]
+
+setindex(X::IntervalBox, y, i) = IntervalBox( setindex(X.v, y, i) )
+
+# iteration:
+
+
+start(X::IntervalBox{N,T}) where {N,T} = 1
+
+next(X::IntervalBox{N,T}, state) where {N,T} = (X[state], state+1)
+
+done(X::IntervalBox{N,T}, state) where {N,T} = state > N
+
+eltype(::Type{IntervalBox{N,T}}) where {N,T} = Interval{T} # Note that this is defined for the type
+
+# length(X::IntervalBox{N,T}) where {N,T} = N
 
 
-IntervalBox(x::Interval) = IntervalBox( (x,) )  # single interval treated as tuple with one element
 
 
 ## arithmetic operations
 # Note that standard arithmetic operations are implemented automatically by FixedSizeArrays.jl
 
-mid(X::IntervalBox) = mid.(X)
+mid(X::IntervalBox) = mid.(X.v)
 
 
 ## set operations
 
 # TODO: Update to use generator
 ⊆(X::IntervalBox{N,T}, Y::IntervalBox{N,T}) where {N,T} =
-    all(X .⊆ Y)
+    all(X.v .⊆ Y.v)
 
 ∩(X::IntervalBox{N,T}, Y::IntervalBox{N,T}) where {N,T} =
-    IntervalBox(X .∩ Y)
+    IntervalBox(X.v .∩ Y.v)
 ∪(X::IntervalBox{N,T}, Y::IntervalBox{N,T}) where {N,T} =
-    IntervalBox(X .∪ Y)
+    IntervalBox(X.v .∪ Y.v)
 
 #=
 On Julia 0.6 can now write
 ∩{N,T}(X::IntervalBox{N,T}, Y::IntervalBox{N,T}) = IntervalBox(NTuple{N, Interval{Float64}}( (X[i] ∩ Y[i]) for i in 1:N))
 =#
 
 
-isempty(X::IntervalBox) = any(isempty, X)
+isempty(X::IntervalBox) = any(isempty, X.v)
 
-diam(X::IntervalBox) = maximum(diam.(X))
+diam(X::IntervalBox) = maximum(diam.(X.v))
 
-emptyinterval(X::IntervalBox{N,T}) where {N,T} = IntervalBox(emptyinterval.(X))
+emptyinterval(X::IntervalBox{N,T}) where {N,T} = IntervalBox(emptyinterval.(X.v))
 
 
 import Base.×
 ×(a::Interval...) = IntervalBox(a...)
-×(a::Interval, b::IntervalBox) = IntervalBox(a, b...)
-×(a::IntervalBox, b::Interval) = IntervalBox(a..., b)
-×(a::IntervalBox, b::IntervalBox) = IntervalBox(a..., b...)
+×(a::Interval, b::IntervalBox) = IntervalBox(a, b.v...)
+×(a::IntervalBox, b::Interval) = IntervalBox(a.v..., b)
+×(a::IntervalBox, b::IntervalBox) = IntervalBox(a.v..., b.v...)
 
 IntervalBox(x::Interval, ::Type{Val{n}}) where {n} = IntervalBox(SVector(ntuple(i->x, Val{n})))
 
 IntervalBox(x::Interval, n::Int) = IntervalBox(x, Val{n})
+
+dot(x::IntervalBox, y::IntervalBox) = dot(x.v, y.v)
+length(x::IntervalBox) = length(x.v)

src/multidim/multidim.jl

@@ -1,2 +1,3 @@
 include("intervalbox.jl")
 include("setdiff.jl")
+include("arithmetic.jl")

src/multidim/setdiff.jl

@@ -43,7 +43,7 @@ Algorithm: Start from the total overlap (in all directions);
 expand each direction in turn.
 """
 function setdiff(A::IntervalBox{N,T}, B::IntervalBox{N,T}) where {N,T}
-    X = [labelled_setdiff(a,b) for (a, b) in zip(A, B)]
+    X = [labelled_setdiff(a, b) for (a, b) in zip(A.v, B.v)]
     # ordered such that the first in each is the excluded interval
 
     first = [ i[1] for i in X ]
@@ -61,7 +61,7 @@ function setdiff(A::IntervalBox{N,T}, B::IntervalBox{N,T}) where {N,T}
         for which in X[dimension][2:end]
             excluded[dimension] = which[1]
             push!(result_list,
-                IntervalBox(excluded[1:dimension]..., A[dimension+1:end]...))
+                IntervalBox(excluded[1:dimension]..., A[dimension+1:N]...))
         end
     end
 

test/multidim_tests/multidim.jl

@@ -1,14 +1,27 @@
 using IntervalArithmetic
 using Base.Test
+using StaticArrays
 
 
 @testset "Operations on boxes" begin
     A = IntervalBox(1..2, 3..4)
     B = IntervalBox(0..2, 3..6)
 
-    @test 2*A == IntervalBox(2..4, 6..8)
+    @test 2*A == A*2 == IntervalBox(2..4, 6..8)
+    @test typeof(2*A) == IntervalBox{2, Float64}
     @test A + B == IntervalBox(1..4, 6..10)
+    @test 2 + A == IntervalBox(3..4,5..6)
+    @test A + 2 == IntervalBox(3..4,5..6)
+    @test -A == IntervalBox((-2)..(-1), (-4)..(-3))
+    @test 2 - A == IntervalBox(0..1, (-2)..(-1))
+    @test B - 2 == IntervalBox((-2)..0, 1..4)
     @test dot(A, B) == @interval(9, 28)
+    @test A .* B == IntervalBox(0..4, 9..24)
+    @test A ./ A == IntervalBox((0.5)..2, (0.75)..(4/3))
+    @test 1 ./ B == IntervalBox((0.5)..Inf, (1/6)..(1/3))
+    @test B ./ 1 == B
+    @test A .^ 2 == IntervalBox(1..4, 9..16)
+    @test B .^ 0.5 == IntervalBox(@interval(0,sqrt(2)), @interval(sqrt(3),sqrt(6)))
 
     @test A ⊆ B
     @test A ∩ B == A
@@ -32,7 +45,17 @@ using Base.Test
     @test isa(Y, IntervalBox)
     @test length(Y) == 1
     @test Y == IntervalBox( (Interval(1., 2.),) )
+    @test typeof(Y) == IntervalBox{1, Float64}
+end
+
+@testset "Functions on boxes" begin
+    A = IntervalBox(1..2, 3..4)
 
+    @test exp.(A) == IntervalBox(exp(A[1]), exp(A[2]))
+    @test typeof(exp.(A)) == IntervalBox{2,Float64}
+    @test log.(A) == IntervalBox(log(A[1]), log(A[2]))
+    @test sqrt.(A) == IntervalBox(sqrt(A[1]), sqrt(A[2]))
+    @test inv.(A) == IntervalBox(inv(A[1]), inv(A[2]))
 end
 
 # @testset "@intervalbox tests" begin
@@ -121,3 +144,28 @@ end
     @test IntervalBox(1..2, 3) == IntervalBox(1..2, Val{3})
 
 end
+
+@testset "getindex and setindex" begin
+    X = IntervalBox(3..4, 5..6)
+    @test X[1] == 3..4
+    @test X[2] == 5..6
+    @test_throws BoundsError X[3]
+
+    @test setindex(X, 5..5, 2) == IntervalBox(3..4, 5..5)
+    @test_throws BoundsError setindex(X, 5..5, 3)
+end
+
+@testset "Iteration" begin
+    X = IntervalBox(3..4, 5..6)
+    Y = collect(X)
+    @test Y == [3..4, 5..6]
+    @test eltype(Y) == Interval{Float64}
+end
+
+@testset "Broadcasting" begin
+    X = IntervalBox(3..4, 5..6)
+
+    @test sin.(X) == IntervalBox(sin(X[1]), sin(X[2]))
+    @test mid.(X) == SVector(mid(X[1]), mid(X[2]))
+    @test diam.(X) == SVector(diam(X[1]), diam(X[2]))
+end