fix parallel reach plot bug

Author: JianqiangDing

example/benchmark_boundary_analysis_vs_naive_partition_cmp.py

@@ -1,19 +1,26 @@
 import numpy as np
 from pybdr.geometry import Interval, Zonotope, Geometry
 from pybdr.geometry.operation import cvt2, partition, boundary
-from pybdr.algorithm import ASB2008CDC
+from pybdr.algorithm import ASB2008CDCParallel
 from pybdr.dynamic_system import NonLinSys
 from pybdr.model import *
 from pybdr.util.visualization import plot, plot_cmp
 from pybdr.util.functional import performance_counter, performance_counter_start
 
 if __name__ == '__main__':
+    # tuples_list = [(1, 'a', 'x', 'u'), (2, 'b', 'y', 'v'), (3, 'c', 'z', 'w')]
+    #
+    # lists_list = [list(t) for t in zip(*tuples_list)]
+    #
+    # print(lists_list)
+    #
+    # exit(False)
     time_start = performance_counter_start()
     # init dynamic system
-    system = NonLinSys(Model(lotka_volterra_2d, [2, 1]))
+    # system = NonLinSys(Model(lotka_volterra_2d, [2, 1]))
 
     # settings for the computation
-    options = ASB2008CDC.Options()
+    options = ASB2008CDCParallel.Options()
     options.t_end = 2.2
     options.step = 0.005
     options.tensor_order = 3
@@ -26,29 +33,37 @@
     Zonotope.REDUCE_METHOD = Zonotope.REDUCE_METHOD.GIRARD
     Zonotope.ORDER = 50
 
-    z = Interval.identity(2) * 0.5 + 3
+    init_set = Interval.identity(2) * 0.5 + 3
 
-    options.r0 = [cvt2(z, Geometry.TYPE.ZONOTOPE)]
-    _, tp_whole, _, _ = ASB2008CDC.reach(system, options)
+    _, tp_whole = ASB2008CDCParallel.reach(lotka_volterra_2d, [2, 1], options,
+                                           cvt2(init_set, Geometry.TYPE.ZONOTOPE))
 
+    # NAIVE PARTITION
     # --------------------------------------------------------
     # options.r0 = partition(z, 1, Geometry.TYPE.ZONOTOPE)
+    # xs = partition(init_set, 1, Geometry.TYPE.ZONOTOPE)
     # 4
-    # ASB2008CDC cost: 43.344963666000005s
+    # ASB2008CDCParallel.reach_parallel cost: 20.126129667s
     # --------------------------------------------------------
-    # options.r0 = partition(z, 0.5, Geometry.TYPE.ZONOTOPE)
+    # xs = partition(init_set, 0.5, Geometry.TYPE.ZONOTOPE)
     # 9
-    # ASB2008CDC cost: 3868912500001s
+    # ASB2008CDCParallel.reach_parallel cost: 23.938516459000002s
     # --------------------------------------------------------
-    options.r0 = partition(z, 0.2, Geometry.TYPE.ZONOTOPE)
+    xs = partition(init_set, 0.2, Geometry.TYPE.ZONOTOPE)
     # 36
-    # ASB2008CDC cost: 317.59988937500003s
+    # ASB2008CDCParallel.reach_parallel cost: 65.447113125s
     # --------------------------------------------------------
-    # options.r0 = partition(z, 0.1, Geometry.TYPE.ZONOTOPE)
 
-    print(len(options.r0))
+    #  BOUNDAYR ANALYSIS
+    # --------------------------------------------------------
+    xs = boundary(init_set, 1, Geometry.TYPE.ZONOTOPE)
+    # 8
+    # ASB2008CDCParallel.reach_parallel cost: 22.185758250000003s
+
+    print(len(xs))
+
+    _, tp_part_00 = ASB2008CDCParallel.reach_parallel(lotka_volterra_2d, [2, 1], options, xs)
 
-    _, tp_part_00, _, _ = ASB2008CDC.reach(system, options)
-    performance_counter(time_start, "ASB2008CDC")
+    performance_counter(time_start, "ASB2008CDCParallel.reach_parallel")
 
     plot_cmp([tp_whole, tp_part_00], [0, 1], cs=["#FF5722", "#303F9F"])

pybdr/algorithm/asb2008cdc_parallel.py

@@ -189,6 +189,10 @@ def reach(cls, dyn: Callable, dims, opts: Options, x: Zonotope):
 
     @classmethod
     def reach_parallel(cls, dyn: Callable, dims, opts: Options, xs: [Zonotope]):
+
+        def ll_decompose(ll):
+            return [list(group) for group in zip(*ll)]
+
         # init containers for storing the results
         ri = []
 
@@ -203,6 +207,6 @@ def reach_parallel(cls, dyn: Callable, dims, opts: Options, xs: [Zonotope]):
                 except Exception as e:
                     raise e
 
-            ri = [list(group) for group in zip(*ri)]
+            ri = ll_decompose(ri)
 
-            return ri
+            return ll_decompose(ri[0]), ll_decompose(ri[1])

pybdr/algorithm/reach_linear_zono_algo3_parallel.py

@@ -161,6 +161,9 @@ def reach(cls, lin_sys: LinearSystemSimple, opts: Settings, x: Zonotope):
     @classmethod
     def reach_parallel(cls, lin_sys: LinearSystemSimple, opts: Settings, xs: [Zonotope]):
 
+        def ll_decompose(ll):
+            return [list(group) for group in zip(*ll)]
+
         with ProcessPoolExecutor() as executor:
             partial_reach = partial(cls.reach, lin_sys, opts)
 
@@ -174,6 +177,6 @@ def reach_parallel(cls, lin_sys: LinearSystemSimple, opts: Settings, xs: [Zonoto
                 except Exception as exc:
                     raise exc
 
-            ri = [list(group) for group in zip(*ri)]
+            ri = ll_decompose(ri)
 
-            return ri
+            return ri[0], ll_decompose(ri[1]), ri[2], ri[3]

pybdr/util/visualization/plot.py

@@ -139,11 +139,13 @@ def __2d_plot_cmp(collections, dims, width, height, xlim, ylim, cs, filled):
     for i in range(len(collections)):
         this_color = plt.cm.turbo(i / len(collections)) if cs is None else cs[i]
 
-        geos = (
-            [collections[i]]
-            if not isinstance(collections[i], list)
-            else list(itertools.chain.from_iterable(collections[i]))
-        )
+        if isinstance(collections[i][0], list):
+            geos = list(itertools.chain.from_iterable(collections[i]))
+        else:
+            geos = collections[i]
+
+        # geos = [collections[i]] if not isinstance(collections[i], list) else list(
+        #     itertools.chain.from_iterable(collections[i]))
 
         for geo in geos:
             if isinstance(geo, np.ndarray):

test/algorithm/test_reach_linear_zono_algo3_parallel.py

@@ -16,20 +16,22 @@ def test_reach_linear_zono_algo3_case_00():
     lin_sys = LinearSystemSimple(xa, ub)
     opts = ReachLinearZonoAlgo3Parallel.Settings()
     opts.t_end = 5
-    opts.step = 0.04
+    opts.step = 0.01
     opts.eta = 4
     # opts.x0 = cvt2(Interval([0.9, 0.9], [1.1, 1.1]), Geometry.TYPE.ZONOTOPE)
     x0 = Interval([0.9, 0.9], [1.1, 1.1])
     opts.u = cvt2(ub @ Interval(0.1, 0.3), Geometry.TYPE.ZONOTOPE)  # u must contain origin
-    x0_bounds = boundary(x0, 0.0001, Geometry.TYPE.ZONOTOPE)
+    x0_bounds = boundary(x0, 0.1, Geometry.TYPE.ZONOTOPE)
 
     print(len(x0_bounds))
 
     _, ri, _, _ = ReachLinearZonoAlgo3Parallel.reach_parallel(lin_sys, opts, x0_bounds)
 
     performance_counter(time_tag, "reach_linear_zono_algo3_parallel")
 
-    # plot(ri, [0, 1])
+    print(len(ri))
+
+    plot(ri, [0, 1])
 
 
 def test_reach_linear_zono_algo3_case_01():
@@ -49,6 +51,8 @@ def test_reach_linear_zono_algo3_case_01():
 
     x0_bounds = boundary(x0, 0.1, Geometry.TYPE.ZONOTOPE)
 
+    performance_counter(time_tag, "boundary")
+
     print(len(x0_bounds))
 
     _, ri, _, _ = ReachLinearZonoAlgo3Parallel.reach_parallel(lin_sys, opts, x0_bounds)