From a6329c2c3a228c60c58bebe14a0871cd77ca83d5 Mon Sep 17 00:00:00 2001 From: schillic Date: Sat, 20 Jul 2024 21:30:22 +0200 Subject: [PATCH] outsource intersection to set modules --- .../intersection.md | 2 - docs/src/lib/sets/Line2D.md | 1 + docs/src/lib/sets/VPolygon.md | 1 + src/ConcreteOperations/intersection.jl | 78 ------------------- src/Sets/Line2D/Line2DModule.jl | 11 ++- src/Sets/Line2D/intersection.jl | 40 ++++++++++ src/Sets/Universe/UniverseModule.jl | 3 +- src/Sets/Universe/intersection.jl | 1 + src/Sets/VPolygon/VPolygonModule.jl | 5 +- src/Sets/VPolygon/intersection.jl | 37 +++++++++ 10 files changed, 92 insertions(+), 87 deletions(-) create mode 100644 src/Sets/Line2D/intersection.jl create mode 100644 src/Sets/Universe/intersection.jl create mode 100644 src/Sets/VPolygon/intersection.jl diff --git a/docs/src/lib/concrete_binary_operations/intersection.md b/docs/src/lib/concrete_binary_operations/intersection.md index 5f72205673..04afff7a92 100644 --- a/docs/src/lib/concrete_binary_operations/intersection.md +++ b/docs/src/lib/concrete_binary_operations/intersection.md @@ -11,7 +11,6 @@ CurrentModule = LazySets ```@docs intersection(::AbstractSingleton, ::LazySet) -intersection(::Line2D, ::Line2D) intersection(::AbstractHyperrectangle, ::AbstractHyperrectangle) intersection(::Interval, ::HalfSpace) intersection(::Interval, ::Hyperplane) @@ -19,7 +18,6 @@ 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) diff --git a/docs/src/lib/sets/Line2D.md b/docs/src/lib/sets/Line2D.md index 7c87563a2d..d11bf088f7 100644 --- a/docs/src/lib/sets/Line2D.md +++ b/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 diff --git a/docs/src/lib/sets/VPolygon.md b/docs/src/lib/sets/VPolygon.md index 6f80bc0a1b..4526834f16 100644 --- a/docs/src/lib/sets/VPolygon.md +++ b/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) ``` diff --git a/src/ConcreteOperations/intersection.jl b/src/ConcreteOperations/intersection.jl index 3f2d8f4bac..accf3c8505 100644 --- a/src/ConcreteOperations/intersection.jl +++ b/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 diff --git a/src/Sets/Line2D/Line2DModule.jl b/src/Sets/Line2D/Line2DModule.jl index 78a1c33ded..d4a9610acc 100644 --- a/src/Sets/Line2D/Line2DModule.jl +++ b/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") diff --git a/src/Sets/Line2D/intersection.jl b/src/Sets/Line2D/intersection.jl new file mode 100644 index 0000000000..3661416fde --- /dev/null +++ b/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 diff --git a/src/Sets/Universe/UniverseModule.jl b/src/Sets/Universe/UniverseModule.jl index 41bcdfc190..04a5a0f26e 100644 --- a/src/Sets/Universe/UniverseModule.jl +++ b/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") diff --git a/src/Sets/Universe/intersection.jl b/src/Sets/Universe/intersection.jl new file mode 100644 index 0000000000..778b34f393 --- /dev/null +++ b/src/Sets/Universe/intersection.jl @@ -0,0 +1 @@ +intersection(U::Universe, ::Universe) = U diff --git a/src/Sets/VPolygon/VPolygonModule.jl b/src/Sets/VPolygon/VPolygonModule.jl index a555d672d0..af88e6f00a 100644 --- a/src/Sets/VPolygon/VPolygonModule.jl +++ b/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") diff --git a/src/Sets/VPolygon/intersection.jl b/src/Sets/VPolygon/intersection.jl new file mode 100644 index 0000000000..5d2633e6a9 --- /dev/null +++ b/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