Added Base.length for EmptySet and SingleTon

[mossr]

We're using LazySets in our new textbook with Mykel Kochenderfer and would like to use length on EmptySet and Singleton. It is already defined for AbstractArraySet and the code for the other two types is straightforward.

[mossr]

For context, the Julia standard library defines length for empty arrays

julia> length([])
0

and for "singletons"

julia> length(123)
1

[mforets]

initially i didn't make sense of length of a LazySet, since the library is mostly about "continuous" sets.

but there is a well defined notion of length for array sets such as finite unions or other structs which have an array field. then i find it natural that it also works for an empty set and for singletons as you have proposed (even if they don't have such array field).

[schillic]

My intuition for length of a set is not about the number of elements, because we usually deal with dense sets anyway. With your proposal, we would have length(E::EmptySet) == 0 but length(E ∪ E) == 2. My thoughts are more that length is just convenient to have for iterable sets (which we call array sets); you could say this is part of the Julia contract of iterables. I would thus argue that length should return 1 for all non-array sets, just as length is 1 for all non-array numbers, including the number 0. Could you explain the reasoning for EmptySet?

Another note is that, if length exists, I would also expect that getindex(X::SomeNonArraySetType, 1) = X would work then, but this is maybe not necessary.

[mforets]

CI-PR / Julia 1.6 - ubuntu-latest - x64 - pull_request (pull_request) Failing after 360m

not sure why are we testing in such an old version, i suggest to just test in the latest julia stable release

[mforets]

this is in a weird state -- one approval (me) but some questioning by @schillic

let's have another round of discussion on our next dev call (friday 3rd) -- details on https://hackmd.io/@juliareach/dev for anyone who's keen to join

@mossr it's delightful to know that you're using this library in a textbook project! say hi to Mykel for me 👋

[schillic]

@mossr: To move forward, could you explain your use case and in particular why you chose length(::EmptySet) = 0?

[mossr]

Of course. We are implementing coverage metrics for a validation textbook and in the implementation of star discrepancy (eqs. (3-4) of Dang et al., pg. 193 [1]) we need the number of points of P that are inside of b⁻ and b⁺ as Pb⁻ and Pb⁺, respectively. We take an intersection and use length to achieve this, where there exists cases where no points are inside the box, thus length of zero is used.

using LazySets
function star_discrepancy(points, a, b, lengths)
    points_norm = [(point .- a) ./ (b .- a) for point in points]
    P = UnionSetArray(Singleton.(points_norm))
    ranges = [range(0, 1, length)[1:end-1] for length in lengths]
    steps = [Float64(r.step) for r in ranges]
    ℬ = Hyperrectangle(low=zeros(length(a)), high=ones(length(a)))
    
    lbs = []
    ubs = []
    for grid_point in Iterators.product(ranges...)
        b⁻ = Hyperrectangle(low=zeros(length(a)), high=[grid_point...])
        b⁺ = Hyperrectangle(low=zeros(length(a)), high=grid_point .+ steps)
        Pb⁻, Pb⁺ = length(intersection(P, b⁻)), length(intersection(P, b⁺)) # <-- use of length here
        push!(lbs, max(abs(Pb⁻ / length(points) - volume(b⁻) / volume(ℬ)),
                       abs(Pb⁺ / length(points) - volume(b⁺) / volume(ℬ))))
        push!(ubs, max(Pb⁺ / length(points) - volume(b⁻) / volume(ℬ),
                       volume(b⁺) / volume(ℬ) - Pb⁻ / length(points)))
    end
    return maximum(lbs), maximum(ubs)
end

[1] Dang, Thao, and Tarik Nahhal. "Coverage-guided test generation for continuous and hybrid systems." Formal Methods in System Design 34 (2009): 183-213.

[schillic]

Okay, I understand now, thanks. The implementation of intersection with a UnionSetArray generally returns a UnionSetArray, but:

1. It tries to be smart and simplifies the union with just one set in it to the set itself. We could remove this special case, but it would not help because we still have the other case below.
2. We do not want to have UnionSetArrays with no sets inside (because that comes with many other problems). Hence we simplify an empty union to an EmptySet.

LazySets.jl/src/ConcreteOperations/intersection.jl

Lines 848 to 867 in 97c91fc

	@commutative function intersection(cup::UnionSetArray, X::LazySet)
	return _intersection_usa(cup, X)
	end
	function _intersection_usa(cup::UnionSetArray, X::LazySet)
	sets = [intersection(Y, X) for Y in cup]
	l = length(sets)
	filter!(!isempty, sets)
	if length(sets) > 1
	if length(sets) < l
	sets = [X for X in sets] # re-allocate set-specific array
	end
	return UnionSetArray(sets)
	elseif length(sets) == 1
	return sets[1]
	else
	N = promote_type(eltype(cup), eltype(X))
	return EmptySet{N}(dim(X))
	end
	end

Since you are only interested in the number of nonempty intersections, your proposal is one way to deal with these special cases.

However, I still do not think that length of an empty set should be 0.

Since you just have singletons, I actually find it cleaner if you avoid the set view and operate on the element level. Here are three proposals, in the order of suggestion:

1. Define P as an array of Vectors and then do:

Pb⁻ = length(filter(p -> p ∈ b⁻, P))
Pb⁻ = length(filter(∈(b⁻), P))  # alternative syntax

2. Define P as an array of Singletons and then do:

Pb⁻ = length(filter(p -> p ⊆ b⁻, P))

3. Keep P as a UnionSetArray of Singletons.

Pb⁻ = filter(!isempty, intersection.(array(P), Ref(b⁻)))

Would that be okay?

(Side note: You do not have to recompute volume(ℬ) in the loop.)

[mforets]

Closing as stale.

[mossr]

Thank you @schillic for putting those options together. We discussed how we want to proceed and will use the syntax Pb⁻ = length(filter(p -> p ∈ b⁻, P)) as suggested in option 1. Also, @mforets, Mykel says hi back :)