#1209 - Concrete projection of Hyperplane & Line2D

docs/src/lib/sets/HalfSpace.md

@@ -27,7 +27,7 @@ tosimplehrep(::AbstractVector{LC}) where {N, LC<:LinearConstraint{N}}
 remove_redundant_constraints
 remove_redundant_constraints!
 complement(::HalfSpace)
-project(::HalfSpace, ::AbstractVector{Int})
+project(::HalfSpace{N}, ::AbstractVector{Int}) where {N}
 ```
 Inherited from [`LazySet`](@ref):
 * [`norm`](@ref norm(::LazySet, ::Real))

src/Sets/HalfSpace.jl

@@ -756,10 +756,10 @@ julia> project(H, [1])
 Universe{Float64}(1)
 ```
 """
-function project(H::HalfSpace, block::AbstractVector{Int})
+function project(H::HalfSpace{N}, block::AbstractVector{Int}) where {N}
     if constrained_dimensions(H) ⊆ block
         return HalfSpace(H.a[block], H.b)
     else
-        return Universe(length(block))
+        return Universe{N}(length(block))
     end
 end

src/Sets/Hyperplane.jl

@@ -478,6 +478,14 @@ function translate(hp::Hyperplane, v::AbstractVector; share::Bool=false)
     return Hyperplane(a, b)
 end
 
+function project(hp::Hyperplane{N}, block::AbstractVector{Int}) where {N}
+    if constrained_dimensions(hp) ⊆ block
+        return Hyperplane(hp.a[block], hp.b)
+    else
+        return Universe{N}(length(block))
+    end
+end
+
 # ============================================
 # Functionality that requires ModelingToolkit
 # ============================================

src/Sets/Line2D.jl

@@ -368,3 +368,29 @@ function translate(L::Line2D, v::AbstractVector; share::Bool=false)
     b = L.b + dot(L.a, v)
     return Line2D(a, b)
 end
+
+function project(L::Line2D{N}, block::AbstractVector{Int}) where {N}
+    m = length(block)
+    if m == 2
+        @assert ispermutation(block, 1:2) "invalid dimensions $block for projection"
+        return L  # no projection
+    elseif m == 1
+        # projection to dimension i
+        cdims = constrained_dimensions(L)
+        if length(cdims) == 1
+            @inbounds if cdims[1] == block[1]
+                # L: aᵢxᵢ = b where aᵢ ≠ 0
+                return Singleton([L.b / L.a[cdims[1]]])
+            else
+                # L: aⱼxⱼ = b where i ≠ j
+                return Universe{N}(1)
+            end
+        else
+            # L is constrained in both dimensions
+            @assert length(cdims) == 2
+            return Universe{N}(1)
+        end
+    else
+        throw(ArgumentError("cannot project a two-dimensional line to $m dimensions"))
+    end
+end

test/unit_HalfSpace.jl

@@ -117,6 +117,11 @@ for N in [Float64, Rational{Int}, Float32]
         end
     end
 
+    # projection
+    H = HalfSpace(N[1, -1], N(0))  # x <= y
+    @test project(H, [1]) == project(H, [2]) == Universe{N}(1)
+    @test project(H, [1, 2]) == H
+
     # conversion of the normal vector
     hs_sev = HalfSpace(SingleEntryVector(2, 3, N(1)), N(1))
     hs_vec = convert(HalfSpace{N, Vector{N}}, hs_sev)

test/unit_Hyperplane.jl

@@ -82,6 +82,11 @@ for N in [Float64, Rational{Int}, Float32]
         end
     end
 
+    # projection
+    H = Hyperplane(N[1, -1], N(0))  # x = y
+    @test project(H, [1]) == project(H, [2]) == Universe{N}(1)
+    @test project(H, [1, 2]) == H
+
     @test_throws ArgumentError linear_map(M, H, algorithm="inv")
     M = N[2 2; 0 1] # invertible matrix
 

test/unit_Line2D.jl

@@ -91,6 +91,15 @@ for N in [Float64, Rational{Int}, Float32]
     @test lm isa Line2D{N, Vector{N}}
     @test lm.a ≈ N[1/2, -2] && lm.b ≈ N(0)
 
+    # projection
+    L = Line2D(N[1, -1], N(0))  # x = y
+    @test project(L, [1]) == project(L, [2]) == Universe{N}(1)
+    @test project(L, [1, 2]) == L
+    L = Line2D(N[2, 0], N(4))  # x = 2
+    @test project(L, [1]) == Singleton(N[2])
+    @test project(L, [2]) == Universe{N}(1)
+    @test project(L, [1, 2]) == L
+
     # translation
     @test translate(l1, N[1, 2]) == Line2D(a1, N(3))
 end