outsource intersection to set modules

docs/src/lib/concrete_binary_operations/intersection.md

@@ -11,15 +11,13 @@ CurrentModule = LazySets
 
 ```@docs
 intersection(::AbstractSingleton, ::LazySet)
-intersection(::Line2D, ::Line2D)
 intersection(::AbstractHyperrectangle, ::AbstractHyperrectangle)
 intersection(::Interval, ::HalfSpace)
 intersection(::Interval, ::Hyperplane)
 intersection(::Interval, ::LazySet)
 intersection(::AbstractHPolygon, ::AbstractHPolygon)
 intersection(::AbstractPolyhedron{N}, ::AbstractPolyhedron{N}) where {N}
 intersection(::Union{VPolytope, VPolygon}, ::Union{VPolytope, VPolygon})
-intersection(::VPolygon, ::VPolygon; ::Bool=true)
 intersection(::UnionSet, ::LazySet)
 intersection(::UnionSetArray, ::LazySet)
 intersection(::Universe, ::LazySet)

docs/src/lib/sets/Line2D.md

@@ -23,6 +23,7 @@ rand(::Type{Line2D})
 project(::AbstractVector, ::Line2D)
 σ(::AbstractVector, ::Line2D)
 translate(::Line2D, ::AbstractVector)
+intersection(::Line2D, ::Line2D)
 ```
 
 ```@meta

docs/src/lib/sets/VPolygon.md

@@ -33,6 +33,7 @@ permute(::VPolygon, ::AbstractVector{Int})
 translate(::VPolygon, ::AbstractVector)
 translate!(::VPolygon, ::AbstractVector)
 convex_hull(::VPolygon, ::VPolygon)
+intersection(::VPolygon, ::VPolygon; ::Bool=true)
 minkowski_sum(::VPolygon, ::VPolygon)
 ```
 

src/ConcreteOperations/intersection.jl

@@ -43,47 +43,6 @@ function intersection(S1::AbstractSingleton, S2::AbstractSingleton)
     return _isapprox(element(S1), element(S2)) ? S1 : EmptySet{N}(dim(S1))
 end
 
-"""
-    intersection(L1::Line2D, L2::Line2D)
-
-Compute the intersection of two two-dimensional lines.
-
-### Input
-
-- `L1` -- line
-- `L2` -- line
-
-### Output
-
-Three outcomes are possible:
-
-- If the lines are identical, the result is the first line.
-- If the lines are parallel and not identical, the result is the empty set.
-- Otherwise the result is the set with the unique intersection point.
-
-### Algorithm
-
-We first check whether the lines are parallel.
-If not, we use [Cramer's rule](https://en.wikipedia.org/wiki/Cramer%27s_rule)
-to compute the intersection point.
-
-### Examples
-
-The line ``y = x`` intersected with the line ``y = -x + 1`` respectively with
-itself:
-
-```jldoctest
-julia> intersection(Line2D([-1.0, 1], 0.0), Line2D([1.0, 1], 1.0))
-Singleton{Float64, Vector{Float64}}([0.5, 0.5])
-
-julia> intersection(Line2D([1.0, 1], 1.0), Line2D([1.0, 1], 1.0))
-Line2D{Float64, Vector{Float64}}([1.0, 1.0], 1.0)
-```
-"""
-function intersection(L1::Line2D, L2::Line2D)
-    return _intersection_line2d(L1, L2)
-end
-
 # this method can also be called with `HalfSpace` arguments
 function _intersection_line2d(L1, L2)
     det = right_turn(L1.a, L2.a)
@@ -675,42 +634,6 @@ function intersection(P1::Union{VPolygon,VPolytope},
     return VPolytope(Pint)
 end
 
-"""
-    intersection(P1::VPolygon, P2::VPolygon; apply_convex_hull::Bool=true)
-
-Compute the intersection of two polygons in vertex representation.
-
-### Input
-
-- `P1` -- polygon in vertex representation
-- `P2` -- polygon in vertex representation
-- `apply_convex_hull` -- (default, optional: `true`) if `false`, skip the
-                         computation of the convex hull of the resulting polygon
-
-### Output
-
-A `VPolygon`, or an `EmptySet` if the intersection is empty.
-
-### Algorithm
-
-This function applies the [Sutherland–Hodgman polygon clipping
-algorithm](https://en.wikipedia.org/wiki/Sutherland%E2%80%93Hodgman_algorithm).
-The implementation is based on the one found in
-[rosetta code](http://www.rosettacode.org/wiki/Sutherland-Hodgman_polygon_clipping#Julia).
-"""
-function intersection(P1::VPolygon, P2::VPolygon; apply_convex_hull::Bool=true)
-    v1 = vertices_list(P1)
-    v2 = vertices_list(P2)
-    v12 = _intersection_vrep_2d(v1, v2)
-
-    if isempty(v12)
-        N = promote_type(eltype(P1), eltype(P2))
-        return EmptySet{N}(2)
-    else
-        return VPolygon(v12; apply_convex_hull=apply_convex_hull)
-    end
-end
-
 """
     intersection(cup::UnionSet, X::LazySet)
 
@@ -830,7 +753,6 @@ The set `X`.
 end
 
 # disambiguations
-intersection(U::Universe, ::Universe) = U
 @commutative intersection(U::Universe, P::AbstractPolyhedron) = P
 @commutative intersection(U::Universe, S::AbstractSingleton) = S
 @commutative intersection(U::Universe, X::Interval) = X

src/Sets/Line2D/Line2DModule.jl

@@ -3,8 +3,9 @@ module Line2DModule
 using Reexport, Requires
 
 using ..LazySets: AbstractPolyhedron, AbstractLinearMapAlgorithm,
-                  _constraints_list_hyperplane, _linear_map_hrep,
-                  _non_element_halfspace, _σ_hyperplane_halfspace
+                  _constraints_list_hyperplane, _intersection_line2d,
+                  _linear_map_hrep, _non_element_halfspace,
+                  _σ_hyperplane_halfspace
 using LinearAlgebra: dot
 using Random: AbstractRNG, GLOBAL_RNG
 using ReachabilityBase.Arrays: nonzero_indices
@@ -14,7 +15,7 @@ using ReachabilityBase.Require: require
 
 @reexport import ..API: an_element, constraints_list, dim, isbounded, isempty,
                         isoperationtype, isuniversal, rand, ∈, project, σ,
-                        translate
+                        translate, intersection
 @reexport import ..LazySets: constrained_dimensions, _linear_map_hrep_helper
 @reexport using ..API
 
@@ -23,7 +24,6 @@ export Line2D
 include("Line2D.jl")
 
 include("an_element.jl")
-include("constrained_dimensions.jl")
 include("constraints_list.jl")
 include("dim.jl")
 include("isbounded.jl")
@@ -36,6 +36,9 @@ include("linear_map.jl")
 include("project.jl")
 include("support_vector.jl")
 include("translate.jl")
+include("intersection.jl")
+
+include("constrained_dimensions.jl")
 
 include("init.jl")
 

src/Sets/Line2D/intersection.jl

@@ -0,0 +1,40 @@
+"""
+    intersection(L1::Line2D, L2::Line2D)
+
+Compute the intersection of two two-dimensional lines.
+
+### Input
+
+- `L1` -- line
+- `L2` -- line
+
+### Output
+
+Three outcomes are possible:
+
+- If the lines are identical, the result is the first line.
+- If the lines are parallel and not identical, the result is the empty set.
+- Otherwise the result is the set with the unique intersection point.
+
+### Algorithm
+
+We first check whether the lines are parallel.
+If not, we use [Cramer's rule](https://en.wikipedia.org/wiki/Cramer%27s_rule)
+to compute the intersection point.
+
+### Examples
+
+The line ``y = x`` intersected with the line ``y = -x + 1`` respectively with
+itself:
+
+```jldoctest
+julia> intersection(Line2D([-1.0, 1], 0.0), Line2D([1.0, 1], 1.0))
+Singleton{Float64, Vector{Float64}}([0.5, 0.5])
+
+julia> intersection(Line2D([1.0, 1], 1.0), Line2D([1.0, 1], 1.0))
+Line2D{Float64, Vector{Float64}}([1.0, 1.0], 1.0)
+```
+"""
+function intersection(L1::Line2D, L2::Line2D)
+    return _intersection_line2d(L1, L2)
+end

src/Sets/Universe/UniverseModule.jl

@@ -12,7 +12,7 @@ using ReachabilityBase.Require: require
                         diameter, dim, isbounded, isboundedtype, isempty,
                         isoperationtype, isuniversal, norm, radius, rand,
                         reflect, ∈, permute, project, scale, scale!, ρ, σ,
-                        translate, translate!, cartesian_product
+                        translate, translate!, cartesian_product, intersection
 @reexport import ..LazySets: constrained_dimensions, linear_map_inverse,
                              tosimplehrep
 @reexport using ..API
@@ -44,6 +44,7 @@ include("support_function.jl")
 include("support_vector.jl")
 include("translate.jl")
 include("cartesian_product.jl")
+include("intersection.jl")
 
 include("constrained_dimensions.jl")
 include("linear_map_inverse.jl")

src/Sets/Universe/intersection.jl

@@ -0,0 +1 @@
+intersection(U::Universe, ::Universe) = U

src/Sets/VPolygon/VPolygonModule.jl

@@ -3,8 +3,8 @@ module VPolygonModule
 using Reexport, Requires
 
 using ..LazySets: AbstractPolygon, LazySet, AbstractHPolygon, halfspace_left,
-                  is_right_turn, _area_vlist, _linear_map_vrep,
-                  _minkowski_sum_vrep_2d
+                  is_right_turn, _area_vlist, _intersection_vrep_2d,
+                  _linear_map_vrep, _minkowski_sum_vrep_2d
 using ..HPolygonModule: HPolygon
 using LinearAlgebra: dot
 using Random: AbstractRNG, GLOBAL_RNG, shuffle
@@ -37,6 +37,7 @@ include("project.jl")
 include("support_vector.jl")
 include("translate.jl")
 include("convex_hull.jl")
+include("intersection.jl")
 include("minkowski_sum.jl")
 
 include("remove_redundant_vertices.jl")

src/Sets/VPolygon/intersection.jl

@@ -0,0 +1,37 @@
+"""
+    intersection(P1::VPolygon, P2::VPolygon; apply_convex_hull::Bool=true)
+
+Compute the intersection of two polygons in vertex representation.
+
+### Input
+
+- `P1` -- polygon in vertex representation
+- `P2` -- polygon in vertex representation
+- `apply_convex_hull` -- (default, optional: `true`) if `false`, skip the
+                         computation of the convex hull of the resulting polygon
+
+### Output
+
+A `VPolygon`, or an `EmptySet` if the intersection is empty.
+
+### Algorithm
+
+This function applies the [Sutherland–Hodgman polygon clipping
+algorithm](https://en.wikipedia.org/wiki/Sutherland%E2%80%93Hodgman_algorithm).
+The implementation is based on the one found in
+[rosetta code](http://www.rosettacode.org/wiki/Sutherland-Hodgman_polygon_clipping#Julia).
+"""
+function intersection(P1::VPolygon, P2::VPolygon; apply_convex_hull::Bool=true)
+    v1 = vertices_list(P1)
+    v2 = vertices_list(P2)
+    v12 = _intersection_vrep_2d(v1, v2)
+
+    if isempty(v12)
+        require(@__MODULE__, :LazySets; fun_name="intersection")
+
+        N = promote_type(eltype(P1), eltype(P2))
+        return EmptySet{N}(2)
+    else
+        return VPolygon(v12; apply_convex_hull=apply_convex_hull)
+    end
+end