exporting common intersection preprocess part

src/Approximations/overapproximate.jl

@@ -1,4 +1,4 @@
-import LazySets: block_to_dimension_indices, combine_resulted_cpa
+import LazySets: block_to_dimension_indices, combine_resulted_cpa, get_constrained_subset
 """
     overapproximate(X::S, ::Type{S}) where {S<:LazySet}
 
@@ -539,6 +539,8 @@ variables, or "blocks") which are constrained in `Y`
 Hence we first search for constrained blocks and then take the intersection of a lower-dimensional Cartesian product of these blocks
 with the projection of `Y` onto the variables of these blocks. (This projection is syntactic and exact.)
 The result is a `CartesianProductArray` with the same block structure as in `X`.
+However, when we decompose back our set we overapproximate during projection operation,
+therefore we need to specify overapproximation strategy (it is Hyperrectangle by default)
 """
 function overapproximate(cap::Intersection{N,
                                             <:CartesianProductArray{N},
@@ -547,7 +549,7 @@ function overapproximate(cap::Intersection{N,
 
     X, Y = cap.X, cap.Y
 
-    low_set, vars, block_structure = get_constrained_subset(X, Y)
+    low_set, vars, block_structure, blocks = get_constrained_subset(X, Y)
 
     low_intersection = intersection(low_set, project(Y, vars))
 
@@ -563,7 +565,7 @@ end
 # symmetric method
 function overapproximate(cap::Intersection{N,
                                             <:AbstractPolyhedron{N},
-                                            <:CartesianProductArray{N, T}},
-                            ::Type{CartesianProductArray}, oa, kwargs...) where {N, T<:LazySet{N}}
-    overapproximate(Intersection(cap.Y, cap.X), oa, kwargs)
+                                            <:CartesianProductArray{N}},
+                            ::Type{CartesianProductArray}, oa) where {N}
+    overapproximate(Intersection(cap.Y, cap.X), oa)
 end

src/CartesianProduct.jl

@@ -629,7 +629,7 @@ This vector represents the mapping "second block from dimension 2 to dimension 4
 These blocks contain the dimensions specified in `[2, 4, 8]`.
 """
 function block_to_dimension_indices(cpa::CartesianProductArray{N}, vars::Vector{Int}) where {N}
-    result = fill(length(array(cpa)), (-1, -1))
+    result = fill((-1, -1), length(array(cpa)))
     non_empty_length = 0
     start_index, end_index = 1, 0
     v_i = 1

src/LazySets.jl

@@ -16,7 +16,6 @@ export Approximations
 # Auxiliary functions
 # ===================
 include("helper_functions.jl")
-include("intersection_helper.jl")
 include("comparisons.jl")
 include("macros.jl")
 
@@ -91,24 +90,25 @@ include("Rectification.jl")
 # =============================
 include("convert.jl")
 
+# =====================
+# Approximations module
+# =====================
+include("Approximations/Approximations.jl")
+
 # ===========================
 # Concrete operations on sets
 # ===========================
 include("concrete_intersection.jl")
 include("is_intersection_empty.jl")
 include("is_subset.jl")
 include("difference.jl")
+include("intersection_helper.jl")
 
 # =======
 # Aliases
 # =======
 include("aliases.jl")
 
-# =====================
-# Approximations module
-# =====================
-include("Approximations/Approximations.jl")
-
 # ==========================
 # Parallel algorithms module
 # ==========================

src/intersection_helper.jl

@@ -1,4 +1,5 @@
-function get_constrained_blocks(X::CartesianProductArray{N}, Y::HPolyhedron{N})
+# wrapper for block_to_dimension_indices
+function get_constrained_blocks(X::CartesianProductArray{N}, Y::HPolyhedron{N}) where {N<:Real}
     if isbounded(Y)
         blocks, non_empty_length = block_to_dimension_indices(X)
     else
@@ -8,11 +9,30 @@ function get_constrained_blocks(X::CartesianProductArray{N}, Y::HPolyhedron{N})
     return blocks, non_empty_length
 end
 
+"""
+    get_constrained_subset(X::CartesianProductArray{N, S}, Y::HPolyhedron{N}) where {N<:Real, S<:LazySet{N}}
+
+Preprocess step for intersection between cartesian product array and H-polyhedron.
+Returns low-dimensional HPolytope subset of X in constrained dimensions,
+constrained variables and variables in corresponding blocks, original block structure
+of low-dimensional subset and list of constrained blocks.
+
+### Input
+
+- `X` -- Cartesian product array of convex sets
+- `Y` -- H-polyhedron
+
+### Output
+
+A tuples of low-dimensional set, list of constrained dimensions, original block
+structure of low-dimensional set and corresponding blocks indicies
+"""
 function get_constrained_subset(X::CartesianProductArray{N, S}, Y::HPolyhedron{N}) where {N<:Real, S<:LazySet{N}}
 
     blocks, non_empty_length = get_constrained_blocks(X, Y)
 
-    array = Vector{S}(undef, non_empty_length)
+    array = Vector{S}()
+    sizehint!(array, non_empty_length)
     vars = Vector{Int}()
     block_structure = Vector{UnitRange{Int}}()
 
@@ -21,14 +41,14 @@ function get_constrained_subset(X::CartesianProductArray{N, S}, Y::HPolyhedron{N
         start_index, end_index = blocks[i]
         block_end = last_var + end_index - start_index
         if (start_index != -1)
-            array[i] = X.array[i]
-            append!(vars, start_index:end_index)
+            push!(array, X.array[i])
+            append!(vars, start_index : end_index)
             push!(block_structure, last_var:(block_end))
         end
         last_var += block_end + 1
     end
 
     low_set = HPolytope(constraints_list(CartesianProductArray(array)));
 
-    return low_set, vars, block_structure
+    return low_set, vars, block_structure, blocks
 end

src/is_intersection_empty.jl

@@ -1,3 +1,4 @@
+import LazySets.Approximations: project
 export is_intersection_empty,
        isdisjoint
 
@@ -1385,7 +1386,7 @@ function is_intersection_empty(X::CartesianProductArray{N},
                                Y::HPolyhedron{N}) where {N<:Real}
     approx_low_set, vars, block_structure = get_constrained_subset(X, Y)
 
-    return isdisjoint(approx_low_set, Approximations.project(Y, vars))
+    return isdisjoint(approx_low_set, project(Y, vars))
 end
 
 # symmetric method