updated linear map to use correct rows and updated tests

src/Approximations/overapproximate.jl

@@ -520,55 +520,51 @@ end
     overapproximate(lm::LinearMap{N, <:CartesianProductArray{N,
                             <:LazySet{N}}}, ::Type{<:CartesianProductArray}, oa=Hyperrectangle)  where {N}
 
-Overapproximating a lazy linear map of cartesian product array with template directions for each block.
+Decompose a lazy linear map of a cartesian product array.
 
 ### Input
 
 - `lm`   -- lazy linear map of cartesian product array
-- Type{<:CartesianProductArray} -- type to specify decomposed overapproximation
-- `oa` -- (concrete) direction representation
+- `CartesianProductArray` -- type for dispatch
+- `oa` -- approximation option for decomposition
 
 ### Output
 
-An `CartesianProductArray` with overapproximation of the each block with the directions from `dir`.
+A `CartesianProductArray` representing the decomposed linear map.
 """
 function overapproximate(lm::LinearMap{N, <:CartesianProductArray{N,
                             <:LazySet{N}}}, ::Type{<:CartesianProductArray}, oa=Hyperrectangle) where {N}
 
     @assert size(lm.M, 2) == dim(lm.X) "matrix needs to be commensurate with the cartesian product"
 
-    cp = lm.X
+    cpa = lm.X
     M = lm.M
 
     col_end_ind = 0
 
-    return decomposed_overapproximation(cp, M, col_end_ind, oa)
-end
-
-#template directions
-function decomposed_overapproximation(cpa::CartesianProductArray{N,<:LazySet{N}},
-                                      M::AbstractMatrix{N},
-                                      col_end_ind::Int, oa) where {N}
     array = allocate_result(oa, N)
     sizehint!(array,length(cpa.array))
 
-    for i in 1:size(M, 2)
-        ms = overapproximate_row_blocks(cpa, M, i)
+    i, n = 1, 0
+    for bi in cpa.array
+        n += dim(bi)
+        ms = blocks_linear_map(cpa, M, i, n)
         push!(array, overapproximate(ms, oa))
+        i += n
     end
 
     result = CartesianProductArray(array)
     return result
 end
 
-#overapproximation of linear map to each block
-function overapproximate_row_blocks(cpa::CartesianProductArray{N,<:LazySet{N}},
-                                    M::AbstractMatrix{N}, i::Int) where {N}
+
+function blocks_linear_map(cpa::CartesianProductArray{N,<:LazySet{N}},
+                                    M::AbstractMatrix{N}, row_start_ind::Int, row_end_ind::Int) where {N}
     col_start_ind, col_end_ind = 1, 0
     h_min_sum = MinkowskiSumArray()
     for bi in cpa.array
         col_end_ind += dim(bi)
-        push!(h_min_sum.array, LinearMap(M[i : i, col_start_ind : col_end_ind], bi))
+        push!(h_min_sum.array, LinearMap(M[row_start_ind : row_end_ind, col_start_ind : col_end_ind], bi))
         col_start_ind = col_end_ind + 1
     end
 

test/unit_overapproximate.jl

@@ -53,6 +53,30 @@ for N in [Float64, Rational{Int}, Float32]
     M = Diagonal(N[2, 2])
     OA = overapproximate(M*H, Hyperrectangle)
     @test OA isa Hyperrectangle && OA.center == N[0, 0] && OA.radius == N[1, 2]
+
+    #overapproximation of Minkowski sum of linear maps for each block in the row block
+    i1 = Interval(N[0, 1])
+    h = Hyperrectangle(low=N[3, 4], high=N[5, 7])
+    M = N[1 2 3; 4 5 6; 7 8 9]
+    cpa = CartesianProductArray([i1, h])
+    lm = M * cpa
+
+    oa = Approximations.overapproximate(lm, Hyperrectangle)
+    oa_box = Approximations.overapproximate(lm, Approximations.BoxDirections)
+    d_oa_d_hp = Approximations.overapproximate(lm, CartesianProductArray)
+    d_oa_d_box = Approximations.overapproximate(lm, CartesianProductArray, Approximations.BoxDirections)
+    oa_d_hp = Approximations.overapproximate(d_oa_d_hp)
+    oa_d_box = Approximations.overapproximate(d_oa_d_box, Approximations.BoxDirections)
+
+    @test oa == oa_d_hp
+    @test oa_box == oa_d_box
+
+    for (oax, set_type) in [(d_oa_d_hp, Hyperrectangle), (d_oa_d_box, HPolytope)]
+        @test oax isa CartesianProductArray
+        arr = oax.array
+        @test length(arr) == 2 && dim(arr[1]) == 1 && dim(arr[2]) == 2
+        @test all(X -> X isa set_type, arr)
+    end
 end
 
 # tests that do not work with Rational{Int}
@@ -130,21 +154,4 @@ for N in [Float64, Float32]
     @test Y_polygon ⊆ Y_zonotope
     @test !(Y_zonotope ⊆ Y_polygon)
 
-    #Decomposed approximation of lazy linear map of CartesianProductArray
-    i1 = Interval([0, 1])
-    h = Hyperrectangle(low=[3., 3.], high=[5., 5.])
-    M = [1. 2. 3.; 4. 5. 6.; 7. 8. 9.]
-    cpa = CartesianProductArray([i1, h])
-    lm = M * cpa
-
-    oa = overapproximate(lm, Hyperrectangle)
-    oa_box = overapproximate(lm, Approximations.BoxDirections)
-    d_oa_d_hp = overapproximate(lm, CartesianProductArray)
-    d_oa_d_box = overapproximate(lm, CartesianProductArray, Approximations.BoxDirections)
-    oa_d_hp = overapproximate(d_oa_d_hp)
-    oa_d_box = overapproximate(d_oa_d_box, Approximations.BoxDirections)
-
-    @test oa == oa_d_hp
-    @test oa_box == oa_d_box
-
 end