Add support for `setdiff(1..3, 2..4)`

[hyrodium]

[cstjean]

setdiff([1,2,3], 2..4) would be pretty neat, too.

[daanhb]

I'm just pointing out that DomainSets does this:

julia> using DomainSets

julia> setdiff(1..3, 2..4)
1..2

However, DomainSets also does the following:

julia> setdiff(1..4, 2..3)
(1..2) ∪ (3..4)

which IntervalSets can't do (and won't do because it is not type-stable). Instead, IntervalSets must throw an error in that case, like it does for unions of non-overlapping intervals.

In DomainSets we don't want to change the behaviour of IntervalSets. For unions we resorted to a new function called uniondomain to construct unions without concerns about being type stable:

julia> using DomainSets

julia> (1..2) ∪ (3..4)
ERROR: ArgumentError: Cannot construct union of disjoint sets.

julia> uniondomain(1..2, 3..4)
(1..2) ∪ (3..4)

Once the setdiff feature gets implemented in IntervalSets and throws an error, we'll change DomainSets accordingly to keep the error, while setdiffdomain will be free to violate type stability.

I think similar compatibility concerns arise when subtracting points or vectors of points:

julia> using DomainSets

julia> setdiff(1..5, 2)
(1..2 (closed–open)) ∪ (2..5 (open–closed))

julia> setdiff(1..5, [2,3,4])
(1..5) \ [2, 3, 4]

julia> 2 ∈ ans
false

Again, just pointing this out. There is no more elegant solution than this, IntervalSets should just be the most useful package for intervals it can be and we'll adapt by using setdiffdomain, uniondomain and intersectdomain in more cases.

[daanhb]

In fact starting from v0.6 we'll have:

julia> using DomainSets

julia> setdiff(1..3, 2..4)
1..2 (closed–open)

julia> 2 ∈ ans
false

which will probably just be annoying in most instances, but it is more accurate to remove the point 2. The closed interval would be

julia> closure(setdiff(1..3, 2..4))
1..2

[hyrodium]

This original issue with setdiff can be solved with complement in #123.

[dlfivefifty]

I don't see how as it is two disjoint intervals

[hyrodium]

I was thinking that the setdiff method can be defined as setdiff(i1::Interval, i2::Interval) = intersection(i1, complement(i2)).

[dlfivefifty]

But that intersection will still error

[hyrodium]

Yes, what I was trying to say is, after introducing the complement and ComplementInterval by solving #123, this issue can be solved incidentally. The implementation will be like this:

julia> using IntervalSets

julia> struct ComplementInterval{L,R,T}
           left::T
           right::T
       end

julia> complement(i::Interval{:closed,:closed,T}) where T = ComplementInterval{:open,:open,T}(i.left, i.right)
complement (generic function with 1 method)

julia> Base.setdiff(i1::Interval, i2::Interval) = intersect(i1, complement(i2))

julia> function Base.intersect(i1::Interval{:closed, :closed}, i2::ComplementInterval{:open,:open})
           isempty(i1) && return i1
           i2.right < i2.left && return i1
           i1.left ≤ i2.left ≤ i2.right ≤ i1.right && error("The returned value cannot be disjoint sets.")
           i1.left ≤ i2.left ≤ i1.right ≤ i2.right && return Interval{:closed,:open}(i1.left,i2.left)
           i2.left ≤ i1.left ≤ i2.right ≤ i1.right && return Interval{:open,:closed}(i2.right,i1.right)
       end

julia> setdiff(1..3, 2..4)
1..2 (closed–open)

julia> setdiff(1..3, 2..1)
1..3

julia> setdiff(1..3, 0..2)
2..3 (open–closed)

julia> setdiff(1..3, 2..2)
ERROR: The returned value cannot be disjoint sets.

[dlfivefifty]

Personally I think setdiff should just return a lazy type that does not error and is type stable. Then your code would be Interval(setdiff(1..3, 2..4)), that is, the user would have to explicitly dictate they want an Interval.

[hyrodium]

How should we implement the lazy version of setdiff(1..3, 2..4)? Is it related to #41? I think the behavior (i.e. whether lazy or not) of setdiff should be the same as union.

[dlfivefifty]

https://github.com/JuliaApproximation/DomainSets.jl/blob/03c614a35b05096e1bb0d55516a94dc16729f2e3/src/generic/setoperations.jl#L277

I'm not seriously proposing we do this (now)... just a more hypothetical thought that this might be the best approach.

But probably the best would be to properly support unions of multiple disjoint intervals. The only challenge there is how to handle open/closed. Probably something like this would work:

struct IntervalUnion{T} <: Domain{T}
endoints::Vector{T} # @assert issorted(endpoints). 
isopen::Vector{Bool} 
end

function in(x, d::IntervalDomain)
   for k = 1:2:length(d.endpoints)-1
          if d.isopen[k] && d.isopen[k+1]
                  if d.endpoints[k] < x < d.endpoints[k+1]
                           return true
                  end
         elseif d.isopen[k]
                   …
        else …
          end
    end
    return false
end

[hyrodium]

But probably the best would be to properly support unions of multiple disjoint intervals. The only challenge there is how to handle open/closed. Probably something like this would work:

That implementation seems great! One of my concerns is how to deal with the unbounded intervals, but that can be solved by adding more types such as:

struct LeftUnboundedIntervalUnion{T} <: Domain{T}
    endoints::Vector{T}  # @assert issorted(endpoints) && isodd(length(endpoints))
    isopen::Vector{Bool} 
end

I still think the union(a::Interval, b::Interval) and setdiff(a::Interval, b::Interval) should return Interval. If someone wants a return type IntervalUnion to avoid the error, one can just convert types before passing the function: setdiff(UnionInterval(a), UnionInterval(b)). This approach is similar to ^(2,-3) in Base, i.e. it requires converting to float to avoid errors. https://github.com/JuliaLang/julia/blob/master/base/intfuncs.jl#L252-L256

---

Btw, I think the name of the type Domain{T} is too abstract, and should be replaced with:

"""
    OrderedDomain{T}

An abstract type that represents totally ordered sets such as ``\mathbb{R}`` and ``(1,2] \cup (3,\infty)``.
Where the type parameter `T` represents the type of endpoints (or boundaries) of the set.
Note that the `T` is different from `eltype`.
"""
OrderedDomain{T}

This definition resolves the confusion as commented in https://github.com/JuliaMath/IntervalSets.jl/issues/115#issuecomment-1219874621.

To me the T in Domain{T} is dictating the "precision" of the definition of a set, in the case of the interval the precision of the endpoints. Which ApproxFun then uses to assume the precision of functions on the set.

[dlfivefifty]

This approach is similar to ^(2,-3) in Base, i.e. it requires converting to float to avoid errors.

Touché.

OrderedDomain{T}

This is not the case in DomainSets.jl

[hyrodium]

Sorry for the late reply.

This is not the case in DomainSets.jl

I was concerned that the name of IntervalSets.jl is too specific to define the abstract type Domain{T}. It seems better to define DomainSetsCore.jl which defines Domain, and OrderedDomain{T} <: Domain should be defined in IntervalSets.jl package. (x-ref: https://github.com/JuliaMath/IntervalSets.jl/pull/126#issuecomment-1332012309) (I think Domain should not have the type parameter {T} because one might want to avoid type promotions with such as (1..2) × (1.0 .. 1.5).)