From e2e6ac82aff9806542cf7f83d848207e1968fcc4 Mon Sep 17 00:00:00 2001
From: Maja Warlich <maja@warlich.name>
Date: Fri, 3 Sep 2021 12:16:31 +0200
Subject: [PATCH] No-slip flip & corners fixed.

---
 lbmpy/sparse/mapping.py           |  70 +++----
 lbmpy_tests/test_sparse_lbm.ipynb | 297 +++++++-----------------------
 2 files changed, 95 insertions(+), 272 deletions(-)

diff --git a/lbmpy/sparse/mapping.py b/lbmpy/sparse/mapping.py
index d50992d..0d639cc 100644
--- a/lbmpy/sparse/mapping.py
+++ b/lbmpy/sparse/mapping.py
@@ -192,15 +192,20 @@ class SparseLbPeriodicityMapper:
         self.density_flag = mapping.density_flag
         self.ubb_flag = mapping.ubb_flag
     
-    def cell_idx(self, coordinate: Tuple[int, ...]) -> np.uint32:
-        """Maps from coordinates (x,y,z) or (x,y) tuple to the list index. Raises ValueError if coordinate not found."""
-        coordinate = np.array(coordinate, dtype=self.mapping._coordinate_arr.dtype)
-        left = np.searchsorted(self.mapping._coordinate_arr, coordinate, sorter=self.mapping._sorter, side='left')
-        right = np.searchsorted(self.mapping._coordinate_arr, coordinate, sorter=self.mapping._sorter, side='right')
-        if left + 1 != right:
-            raise IndexError("Coordinate not found")
+    def get_fluid_border_coord(self):
+        fluid_border_coord = []
+        for coord in self.mapping.fluid_coordinates:
+            print(coord)
+    
+    def get_read_idx(self, read, write_idx, direction_idx):  
+        if self.flag_arr[tuple(read)] & self.no_slip_flag: # Read cell is no-slip: flip PDF!
+            #read from write cell, inverse direction
+            return pdf_index(write_idx, inverse_idx(self.method.stencil, direction_idx), len(self.mapping))
+            #periodic_idx_array.append([direction_idx, inverse_idx(self.method.stencil, direction_idx), write, write]) #nur zu debug Zwecken
         else:
-            return self.mapping._sorter[left]
+            return pdf_index(self.mapping.cell_idx(tuple(read)), direction_idx, len(self.mapping))
+            #periodic_idx_array.append([direction_idx, write, read]) #nur zu debug Zwecken
+        
         
     def create_index_arr(self): # erstellt index arrays für ALLE fluid Zellen, die sich am Rand der domain befinden.
         # ein index array für alle Werte, die innerhalb des Blocks verschickt werden
@@ -211,8 +216,10 @@ class SparseLbPeriodicityMapper:
         inner_idx_array = []
         write = [0,0]
         fluid_boundary_mask = self.fluid_flag | self.ubb_flag | self.density_flag
-        for direction_idx, direction in enumerate(stencil):
-            if all(d_i == 0 for d_i in direction):#(direction_idx != 5 and direction_idx != 1) or direction_idx < 5: # direction (0,0) irrelevant
+        for i in range(1, 4, 1):
+            print(i)
+        for direction_idx, direction in enumerate(self.method.stencil):
+            if all(d_i == 0 for d_i in direction):# direction (0,0) irrelevant
                 continue
             print("\n New direction:", direction, ", ", direction_idx)
             for pos in range(0,2): # einmal für x, einmal für y Richtung ... 
@@ -223,7 +230,6 @@ class SparseLbPeriodicityMapper:
                     print("(periodic:)")
                     periodic_idx_array = []
                     coord = int((self.domain_size[pos]-1)*(1-direction[pos])/2)
-                    print("first")
                     start = 1 if direction[sop] == 1 else 0
                     end = self.domain_size[sop]-1 if direction[sop] == -1 else self.domain_size[sop]
                     for i in range(start, end):
@@ -233,23 +239,16 @@ class SparseLbPeriodicityMapper:
                         if not (self.flag_arr[tuple(write)] & self.fluid_flag):
                             continue
                         read = [(write_i - dir_i)%ds_i for write_i, dir_i, ds_i in zip(write, direction, self.domain_size)]
-                        write_idx = pdf_index(self.cell_idx(tuple(write)), direction_idx, len(self.mapping))
-                        read_idx = pdf_index(self.cell_idx(tuple(read)), direction_idx, len(self.mapping))
-                        if self.flag_arr[tuple(read)] & self.no_slip_flag: # Read cell is no-slip: flip PDF!
-                            #read from write cell, inverse direction
-                            read_idx = pdf_index(write_idx, inverse_idx(stencil, direction_idx), len(self.mapping))
+                        write_idx = pdf_index(self.mapping.cell_idx(tuple(write)), direction_idx, len(self.mapping))
+                        read_idx = self.get_read_idx(read, write_idx, direction_idx)
                         periodic_idx_array.append([direction_idx, write_idx, read_idx])
                         # "Die Zelle "write" bekommt ihren neuen Wert der jeweiligen direction von der Zelle "read"
-                        print("write:", write, "read:", read) 
-                        #periodic_idx_array.append([direction_idx, write, read]) #nur zu debug Zwecken
                     result.append(tuple(periodic_idx_array))
                     # inner: wird zwischen benachbarten Zellen *im gleichen Block* geschickt
                     print("(inner:)")
                     pos_bound = int((self.domain_size[pos]-1)*(1+direction[pos])/2)
                     pos_mid = pos_bound+direction[pos]*(-self.domain_size[pos]+1)
                     sop_position = [int((self.domain_size[sop]-1)*(1+direction[sop])/2)] if direction[sop] != 0 else [0, self.domain_size[sop]-1]
-                    #print("pos_bound:", pos_bound, "pos_mid", pos_mid, "sop_position:", sop_position)
-                    print("second")
                     for b in sop_position:
                         for i in range(pos_bound, pos_mid, -direction[pos]):
                             write = [0,0]
@@ -258,22 +257,15 @@ class SparseLbPeriodicityMapper:
                             if not (self.flag_arr[tuple(write)] & self.fluid_flag):
                                 continue
                             read = [write_i - dir_i for write_i, dir_i in zip(write, direction)]
-                            write_idx = pdf_index(self.cell_idx(tuple(write)), direction_idx, len(self.mapping))
-                            read_idx = pdf_index(self.cell_idx(tuple(read)), direction_idx, len(self.mapping))
-                            if self.flag_arr[tuple(read)] & self.no_slip_flag: # Read cell is no-slip: flip PDF!
-                                #read from write cell, inverse direction
-                                read_idx = pdf_index(write_idx, inverse_idx(stencil, direction_idx), len(self.mapping))
+                            write_idx = pdf_index(self.mapping.cell_idx(tuple(write)), direction_idx, len(self.mapping))
+                            read_idx = self.get_read_idx(read, write_idx, direction_idx)
                             inner_idx_array.append([direction_idx, write_idx, read_idx])
-                            print("write:", write, "read:", read)
-                            #inner_idx_array.append([direction_idx, write, read]) #for debug
                 if direction[pos] == 0: #spricht directions 1, 2, 3 und 4 an
                     # inner: wird zwischen benachbarte Zellen *im gleichen Block* geschickt
                     print("(inner:)")
-                    print("third")
                     pos_low = 1
                     pos_high = self.domain_size[pos]-1
                     sop_position = int((self.domain_size[sop]-1)*(1+direction[sop])/2)
-                    #print("pos_low:", pos_low, "pos_high:", pos_high, "sop_position:", sop_position)
                     for i in range(pos_low, pos_high):
                         write = [0,0]
                         write[pos] = i
@@ -281,29 +273,19 @@ class SparseLbPeriodicityMapper:
                         if not (self.flag_arr[tuple(write)] & self.fluid_flag):
                             continue
                         read = [write_i - dir_i for write_i, dir_i in zip(write, direction)]
-                        write_idx = pdf_index(self.cell_idx(tuple(write)), direction_idx, len(self.mapping))
-                        read_idx = pdf_index(self.cell_idx(tuple(read)), direction_idx, len(self.mapping))
-                        if self.flag_arr[tuple(read)] & self.no_slip_flag: # Read cell is no-slip: flip PDF!
-                            #read from write cell, inverse direction
-                            read_idx = pdf_index(write_idx, inverse_idx(stencil, direction_idx), len(self.mapping))
-                        inner_idx_array.append([direction_idx, write_idx, read_idx])
-                        print("write:", write, "read:", read)
-                        #inner_idx_array.append([direction_idx, write, read]) #for debug
+                        write_idx = pdf_index(self.mapping.cell_idx(tuple(write)), direction_idx, len(self.mapping))
+                        read_idx = self.get_read_idx(read, write_idx, direction_idx)
+                        inner_idx_array.append([direction_idx, write, read]) #for debug
             #Four corners: extra periodic_idx_array for each direction 5, 6, 7, 8
             if (direction[0]*direction[1] != 0):
                 write = [int((self.domain_size[0]-1)*(1-direction[0])/2),int((self.domain_size[1]-1)*(1-direction[1])/2)]
                 if not (self.flag_arr[tuple(write)] & self.fluid_flag):
                     continue
                 read = [(write_i - dir_i)%ds_i for write_i, dir_i, ds_i in zip(write, direction, self.domain_size)]
-                write_idx = pdf_index(self.cell_idx(tuple(write)), direction_idx, len(self.mapping))
-                read_idx = pdf_index(self.cell_idx(tuple(read)), direction_idx, len(self.mapping))
-                if self.flag_arr[tuple(read)] & self.no_slip_flag: # Read cell is no-slip: flip PDF!
-                    #read from write cell, inverse direction
-                    read_idx = pdf_index(write_idx, inverse_idx(stencil, direction_idx), len(self.mapping))
+                write_idx = pdf_index(self.mapping.cell_idx(tuple(write)), direction_idx, len(self.mapping))
+                read_idx = self.get_read_idx(read, write_idx, direction_idx)
                 periodic_idx_array.append([direction_idx, write_idx, read_idx])
-                #periodic_idx_array.append([direction_idx, write, read]) #nur zu debug Zwecken
                 result.append(tuple(periodic_idx_array))
-                print("(Ecke) write:", write, "read:", read)
             
         # result enthält *mehrere* index arrays
         #result = list(dict.fromkeys(result)) # entferne doppelte index_arrays: speziell Ecken der Domain
diff --git a/lbmpy_tests/test_sparse_lbm.ipynb b/lbmpy_tests/test_sparse_lbm.ipynb
index 612764e..6b5769d 100644
--- a/lbmpy_tests/test_sparse_lbm.ipynb
+++ b/lbmpy_tests/test_sparse_lbm.ipynb
@@ -142,7 +142,7 @@
    "outputs": [
     {
      "data": {
-      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1gAAAFpCAYAAACBJomJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3df6xlZ3kv9u+Tia00hiub6x+4HoOtaBRw6MXQkePI/5hLktpuFEMUJFsqWK6jCRG+BYmrypcrhbT3HxQFUrhQW87FwlaJES04jLiTOI5L5SDxw2Mz+CcWE4fgYaYeDC0GkRs6M0//OGvo7uGc2TPea82ZM/vzkZb2Xut913rf4y1+PH7e91nV3QEAAGBxP7fREwAAADhdCLAAAABGIsACAAAYiQALAABgJAIsAACAkQiwAAAARrJwgFVVF1fVF6rq6ap6sqrevUafq6vqB1W1Zzj+cNFxAQAAjqqqu6rqYFU9sU77a6rqS1X1T1X1r1e1XVNVz1TV3qq6beb6K6rqgar65vB5zrx5jJHBOpTkvd392iRXJnlXVV22Rr+/7e7Lh+N/HGFcAACAoz6R5JpjtH8/yX+X5E9mL1bVliQfS3JtksuS3DgTz9yW5MHu3pbkweH8mBYOsLr7QHc/Onz/YZKnk1y06HMBAACOV3c/lJUgar32g939cJL/Z1XTFUn2dvez3f2TJJ9Kcv3Qdn2Su4fvdyd5y7x5jLoHq6ouSfKGJF9Zo/nXqurrVfWXVfUrY44LAADwEl2U5LmZ8335/xJGF3T3gWQlsZTk/HkP+/mxZlVVL0vymSTv6e4XVzU/muTV3f2jqrouyV8k2bbOc3Yk2ZEkZ5111n/5mte8ZqwpkuTxg89v9BQAgCXzX5x/wUZP4bTzyCOPvNDd5230PI7Xf/Wms/p73z+80DMeeeyfnkzyn2Yu3dnddy700BW1xrV+qQ8bJcCqqjOyElx9srs/u7p9NuDq7l1V9T9X1bnd/cIafe9McmeSbN++vXfv3j3GFBlc+uEPbvQUAIAls/vd793oKZx2quofNnoOJ+KF7x/OV+7futAzzrjw7/5Td28faUqz9iW5eOZ8a5L9w/fnq+rC7j5QVRcmOTjvYWNUEawkH0/ydHd/aJ0+rxz6paquGMb93qJjAwAALOjhJNuq6tKqOjPJDUl2Dm07k9w0fL8pyefmPWyMDNZVSd6e5PGq2jNce1+SVyVJd9+R5HeT/EFVHUryj0lu6O6XnHYDAAA2k87hPjLpCFV1b5Krk5xbVfuSvD/JGclKTFJVr0yyO8k/S3Kkqt6T5LLufrGqbk1yf5ItSe7q7ieHx34gyaer6pYk307ytnnzWDjA6u4vZu11i7N9Pprko4uOBQAAbD6d5MhL39Z0fGN03zin/f/MyvK/tdp2Jdm1xvXvJXnzicxjtCIXAAAA6zmSaTNYp4pRy7QDAAAsMxksAABgUp3O4SUpwSDAAgAAJjf1HqxThQALAACYVCc5LMACAAAYx7JksBS5AAAAGIkMFgAAMKlOFLkAAAAYy3K8BUuABQAATKzTilwAAACMopPDyxFfKXIBAAAwFhksAABgUh17sAAAAEZSOZza6EmcFAIsAABgUp3kiD1YAAAAnAgZLAAAYHKWCAIAAIygI8ACAAAYzZEWYAEAACxsmTJYilwAAACMRAYLAACYVKdyeElyOwIsAABgcvZgAQAAjGCZ9mAJsAAAgIlVDvdyLBFcjr8SAADgJJDBAgAAJtVJjixJbkeABQAATM4eLAAAgBF024MFAADACZLBAgAAJnfEEkEAAIDFrbwHazkWzy3HXwkAAGyglT1YixxzR6i6q6oOVtUT67RXVX2kqvZW1WNV9cbh+i9X1Z6Z48Wqes/Q9kdV9Z2ZtuvmzUMGCwAAmNRJKtP+iSQfTXLPOu3XJtk2HL+a5PYkv9rdzyS5PEmqakuS7yS5b+a+P+3uPzneSchgAQAAm153P5Tk+8focn2Se3rFl5OcXVUXrurz5iR/193/8FLnIcACAAAmd7hroSPJuVW1e+bYcYJTuCjJczPn+4Zrs25Icu+qa7cOSwrvqqpz5g1iiSAAADCpTo1R5OKF7t6+wP1rlTHsnzZWnZnkt5P8m5n225P8u6Hfv0vywST/7bEGEWABAACTO7LxLxrel+TimfOtSfbPnF+b5NHufv7ohdnvVfVnST4/b5AN/ysBAIDT29Ey7YscI9iZ5B1DNcErk/yguw/MtN+YVcsDV+3RemuSNSsUzpLBAgAANr2qujfJ1VnZq7UvyfuTnJEk3X1Hkl1JrkuyN8mPk9w8c+8vJvmNJL+/6rF/XFWXZyVG/NYa7T9DgAUAAEyq89NCFdON0X3jnPZO8q512n6c5J+vcf3tJzoPARYAADC5k/AerFOCAAsAAJhUd3J444tcnBTL8VcCAACcBDJYAADAxCpH1nwN1elHgAUAAEyqY4ngcauqi6vqC1X1dFU9WVXvXqNPVdVHqmpvVT1WVW9cdFwAAGDzOAXeg3VSjJHBOpTkvd39aFW9PMkjVfVAdz810+faJNuG41eT3D58AgAAp7lO5cjEZdpPFQuHgt19oLsfHb7/MMnTSS5a1e36JPf0ii8nOXvVW5EBAAA2vVH3YFXVJUnekOQrq5ouSvLczPm+4dqBNZ6xI8mOJHnVq1415vQAAIANspmW+S1itL+yql6W5DNJ3tPdL65uXuOWXus53X1nd2/v7u3nnXfeWNMDAAA2SCc50j+30LFZjJLBqqozshJcfbK7P7tGl31JLp4535pk/xhjAwAAp7rK4SUp0z5GFcFK8vEkT3f3h9bptjPJO4Zqglcm+UF3/8zyQAAA4PQjg3Virkry9iSPV9We4dr7krwqSbr7jiS7klyXZG+SHye5eYRxAQAATikLB1jd/cWsvcdqtk8nedeiYwEAAJvTsiwRHLWKIAAAwGrdtamW+S1CgAUAAEzu8JIEWMvxVwIAAJwEMlgAAMCkOskRe7AAAADGUEuzRFCABQAATGrlPVgyWAAAAKM4vCTlH5bjrwQAADgJZLAAAIBJdcoSQQAAgLEcWZLFcwIsAABgUt3JYRksAACAcSzLEsHlyNMBAACcBDJYAADApFaKXCxHbkeABQAATO5wlmOJoAALAACYVMceLAAAAE6QDBYAADCx5dmDtRx/JQAAsKGOpBY65qmqu6rqYFU9sU57VdVHqmpvVT1WVW+caftWVT1eVXuqavfM9VdU1QNV9c3h85x58xBgAQAAkzr6ouFFjuPwiSTXHKP92iTbhmNHkttXtb+puy/v7u0z125L8mB3b0vy4HB+TAIsAABgckf65xY65unuh5J8/xhdrk9yT6/4cpKzq+rCOY+9Psndw/e7k7xl3jwEWAAAwGZwblXtnjl2nOD9FyV5buZ833AtWSl0+NdV9ciq517Q3QeSZPg8f94gilwAAACTWnnR8MJl2l9YtXzvRK01gR4+r+ru/VV1fpIHquobQ0bshMlgAQAAk5u6yMVx2Jfk4pnzrUn2J0l3H/08mOS+JFcMfZ4/uoxw+Dw4bxABFgAAMKmjLxpe5BjBziTvGKoJXpnkB919oKrOqqqXJ0lVnZXkN5M8MXPPTcP3m5J8bt4glggCAACbXlXdm+TqrOzV2pfk/UnOSJLuviPJriTXJdmb5MdJbh5uvSDJfVWVrMRHf97dfzW0fSDJp6vqliTfTvK2efMQYAEAAJOb+kXD3X3jnPZO8q41rj+b5PXr3PO9JG8+kXkIsAAAgGmNt8zvlCfAAgAAJtXJWIUqTnkCLAAAYHLLksFSRRAAAGAkMlgAAMCkjpZpXwYCLAAAYHICLAAAgBF0VBEEAAAYzbJUEVTkAgAAYCQyWAAAwLTaHiwAAIBRqCIIAAAwomUJsOzBAgAAGIkMFgAAMCll2gEAAEbUAiwAAIBxLMt7sARYAADApHqJyrQrcgEAADCSUQKsqrqrqg5W1RPrtF9dVT+oqj3D8YdjjAsAAGwO3bXQsVmMtUTwE0k+muSeY/T52+7+rZHGAwAANg1VBE9Idz9UVZeM8SwAAOD0s5myUIs4mXuwfq2qvl5Vf1lVv7Jep6raUVW7q2r3d7/73ZM4PQAAYAqdlSIXixybxckKsB5N8urufn2Sf5/kL9br2N13dvf27t5+3nnnnaTpAQAALO6kBFjd/WJ3/2j4vivJGVV17skYGwAA2GC9Uqp9kWOzOCnvwaqqVyZ5vru7qq7ISmD3vZMxNgAAsPG8aPgEVNW9Sa5Ocm5V7Uvy/iRnJEl335Hkd5P8QVUdSvKPSW7o3kxxKAAA8FJ1lqfIxVhVBG+c0/7RrJRxBwAAOG2dlCWCAADAMttclQAXIcACAAAmtywbhARYAADA5OzBAgAAGMFKqfXlCLBO1ouGAQAATnsCLAAAYHJHuhY65qmqu6rqYFU9sU57VdVHqmpvVT1WVW8crl9cVV+oqqer6smqevfMPX9UVd+pqj3Dcd28eQiwAACAya0sE3zpx3H4RJJrjtF+bZJtw7Ejye3D9UNJ3tvdr01yZZJ3VdVlM/f9aXdfPhy75k3CHiwAAGByU+/B6u6HquqSY3S5Psk93d1JvlxVZ1fVhd19IMmB4Rk/rKqnk1yU5KmXMg8ZLAAAYFKdSvdiR5Jzq2r3zLHjBKdxUZLnZs73Ddd+agjQ3pDkKzOXbx2WFN5VVefMG0SABQAAbAYvdPf2mePOE7x/rRTaTxcfVtXLknwmyXu6+8Xh8u1JfinJ5VnJcn1w3iCWCAIAAJM7Bd4zvC/JxTPnW5PsT5KqOiMrwdUnu/uzRzt09/NHv1fVnyX5/LxBZLAAAIBpDe/BWnCJ4KJ2JnnHUE3wyiQ/6O4DVVVJPp7k6e7+0OwNVXXhzOlbk6xZoXCWDBYAADC9iVNYVXVvkquzsldrX5L3JzkjSbr7jiS7klyXZG+SHye5ebj1qiRvT/J4Ve0Zrr1vqBj4x1V1+TD7byX5/XnzEGABAACbXnffOKe9k7xrjetfzNr7s9Ldbz/ReQiwAACAyU1dpv1UIcACAAAmd5wvC970BFgAAMCkOjJYAAAA4+gkSxJgKdMOAAAwEhksAABgcvZgAQAAjEWABQAAMIZS5AIAAGA0S5LBUuQCAABgJDJYAADAtNp7sAAAAMazJEsEBVgAAMBJsBwZLHuwAAAARiKDBQAATM8SQQAAgJEIsAAAAEbQSVQRBAAAGEcvSQZLkQsAAICRyGABAADTW5IMlgALAACYnj1YAAAA4ygZLAAAgBF0lmaJoCIXAAAAI5HBAgAAJlb2YAEAAIxmSZYICrAAAIDpLUmAZQ8WAADASGSwAACA6S1JBkuABQAATKujyAUAAMBYluVFw6Pswaqqu6rqYFU9sU57VdVHqmpvVT1WVW8cY1wAAGCT6AWPORaJSarqmqp6Zmi7beb6K6rqgar65vB5zrx5jFXk4hNJrjlG+7VJtg3HjiS3jzQuAABA8hJjkqrakuRjQ/tlSW6sqsuGe25L8mB3b0vy4HB+TKMEWN39UJLvH6PL9Unu6RVfTnJ2VV04xtgAAAALxCRXJNnb3c9290+SfGroe/Seu4fvdyd5y7x5nKwy7RcleW7mfN9wDQAAWALVix0jWC8mOVasckF3H0iS4fP8eYOcrCIXa5UMWfMfU1XtyErKLq961aumnBMAAHCyLF5F8Nyq2j1zfmd333kC968Xkxx3rHI8TlaAtS/JxTPnW5PsX6vj8A/pziTZvn37ktQaAQCA09hxFqqY44Xu3r7A/evFJGeucz1Jnq+qC7v7wLCc8OC8QU7WEsGdSd4xVO64MskPjqbaAAAAToL1YpKHk2yrqkur6swkNwx9j95z0/D9piSfmzfIKBmsqro3ydVZSdvtS/L+JGckSXffkWRXkuuS7E3y4yQ3jzEuAACwSUy8Nu2lxiTdfaiqbk1yf5ItSe7q7ieHx34gyaer6pYk307ytnnzGCXA6u4b57R3kneNMRYAALD5TP2i4UViku7elZUAbPX17yV584nM42TtwQIAAJbZklRXOFl7sAAAAE57MlgAAMD0liSDJcACAAAmNeLLgk95AiwAAGB6i79oeFMQYAEAANNbkgyWIhcAAAAjkcECAAAmZw8WAADAWARYAAAAI1iiKoL2YAEAAIxEBgsAAJjekmSwBFgAAMD0BFgAAADjsAcLAACAEyLAAgAAGIklggAAwPSWZImgAAsAAJjWEr0HS4AFAABMT4AFAAAwkiUJsBS5AAAAGIkMFgAAMKmKPVgAAADjEWABAACMYImqCNqDBQAAMBIZLAAAYHpLksESYAEAANMTYAEAAIxjWfZgCbAAAIDpLUmApcgFAADASGSwAACAaXWWJoMlwAIAACa3LHuwLBEEAACm1wsec1TVNVX1TFXtrarb1mg/p6ruq6rHquqrVfW64fovV9WemePFqnrP0PZHVfWdmbbr5s1DBgsAAJjclBmsqtqS5GNJfiPJviQPV9XO7n5qptv7kuzp7rdW1WuG/m/u7meSXD7znO8kuW/mvj/t7j853rnIYAEAAJvdFUn2dvez3f2TJJ9Kcv2qPpcleTBJuvsbSS6pqgtW9Xlzkr/r7n94qRMRYAEAANNbfInguVW1e+bYMfP0i5I8N3O+b7g26+tJfidJquqKJK9OsnVVnxuS3Lvq2q3DssK7quqceX+mAAsAAJjWosHVSoD1QndvnznunBmh1hl11geSnFNVe5L8qyRfS3Lopw+oOjPJbyf5X2fuuT3JL2VlCeGBJB+c96fagwUAAEyqsnYENKJ9SS6eOd+aZP9sh+5+McnNSVJVleTvh+Ooa5M82t3Pz9zz0+9V9WdJPj9vIjJYAADAZvdwkm1VdemQibohyc7ZDlV19tCWJL+X5KEh6DrqxqxaHlhVF86cvjXJE/MmIoMFAABMb8Iqgt19qKpuTXJ/ki1J7uruJ6vqnUP7HUlem+Seqjqc5Kkktxy9v6p+MSsVCH9/1aP/uKouH2b/rTXaf4YACwAAmNzULxru7l1Jdq26dsfM9y8l2bbOvT9O8s/XuP72E52HAAsAAJjexAHWqUKABQAATG9JAixFLgAAAEYigwUAAEyrp9+DdaoQYAEAANNbkgBrlCWCVXVNVT1TVXur6rY12q+uqh9U1Z7h+MMxxgUAADaH6sWOzWLhDFZVbUnysazUjd+X5OGq2tndT63q+rfd/VuLjgcAAGxCmyhIWsQYGawrkuzt7me7+ydJPpXk+hGeCwAAsKmMEWBdlOS5mfN9w7XVfq2qvl5Vf1lVv7Lew6pqR1Xtrqrd3/3ud0eYHgAAsNGWZYngGAFWrXFt9T+CR5O8urtfn+TfJ/mL9R7W3Xd29/bu3n7eeeeNMD0AAGBD9QjHJjFGgLUvycUz51uT7J/t0N0vdvePhu+7kpxRVeeOMDYAALAZCLCO28NJtlXVpVV1ZpIbkuyc7VBVr6yqGr5fMYz7vRHGBgAAOGUsXEWwuw9V1a1J7k+yJcld3f1kVb1zaL8jye8m+YOqOpTkH5Pc0N2bKA4FAABeqsrm2ke1iFFeNDws+9u16todM98/muSjY4wFAABsQgIsAACAcdSSLGATYAEAANPaZIUqFjFGkQsAAAAigwUAAJwEilwAAACMRYAFAAAwDhksAACAsSxJgKXIBQAAwEhksAAAgGm1JYIAAADjEWABAAAsrrI8GSx7sAAAAEYigwUAAEyvlyOFJcACAAAmtyxLBAVYAADAtDqKXAAAAIyljmz0DE4ORS4AAABGIoMFAABMb0mWCMpgAQAAk6te7Jj7/KprquqZqtpbVbet0X5OVd1XVY9V1Ver6nUzbd+qqserak9V7Z65/oqqeqCqvjl8njNvHgIsAABgWp2VMu2LHMdQVVuSfCzJtUkuS3JjVV22qtv7kuzp7n+R5B1JPryq/U3dfXl3b5+5dluSB7t7W5IHh/NjEmABAACTmziDdUWSvd39bHf/JMmnkly/qs9lWQmS0t3fSHJJVV0w57nXJ7l7+H53krfMm4gACwAA2AzOrardM8eOmbaLkjw3c75vuDbr60l+J0mq6ookr06ydWjrJH9dVY+seu4F3X0gSYbP8+dNUpELAABgeosXuXhh1fK9WXUcI34gyYerak+Sx5N8Lcmhoe2q7t5fVecneaCqvtHdD72USQqwAACASVWOr1DFAvYluXjmfGuS/bMduvvFJDcnSVVVkr8fjnT3/uHzYFXdl5Ulhw8leb6qLuzuA1V1YZKD8yZiiSAAADCtRQtczClykeThJNuq6tKqOjPJDUl2znaoqrOHtiT5vSQPdfeLVXVWVb186HNWkt9M8sTQb2eSm4bvNyX53LyJyGABAACbWncfqqpbk9yfZEuSu7r7yap659B+R5LXJrmnqg4neSrJLcPtFyS5byWplZ9P8ufd/VdD2weSfLqqbkny7SRvmzcXARYAADC5iZcIprt3Jdm16todM9+/lGTbGvc9m+T16zzze0nefCLzEGABAADTmzjAOlUIsAAAgMlNncE6VQiwAACAaXWSI8sRYakiCAAAMBIZLAAAYHrLkcASYAEAANOzBwsAAGAs818WfFoQYAEAAJNblgyWIhcAAAAjkcECAACm1VHkAgAAYAyVpOzBAgAAGMmRjZ7AyWEPFgAAwEhksAAAgMlZIggAADAGRS4AAADG0l40DAAAMBYvGgYAAOCEjBJgVdU1VfVMVe2tqtvWaK+q+sjQ/lhVvXGMcQEAgE2ie7Fjk1h4iWBVbUnysSS/kWRfkoeramd3PzXT7dok24bjV5PcPnwCAACnu07Ke7CO2xVJ9nb3s939kySfSnL9qj7XJ7mnV3w5ydlVdeEIYwMAAJvBkmSwxgiwLkry3Mz5vuHaifZJklTVjqraXVW7v/vd744wPQAAgJNjjACr1ri2OsQ8nj4rF7vv7O7t3b39vPPOW3hyAADAKaAXPDaJMcq070ty8cz51iT7X0IfAADgNFWbaJnfIsbIYD2cZFtVXVpVZya5IcnOVX12JnnHUE3wyiQ/6O4DI4wNAABsBkuyB2vhDFZ3H6qqW5Pcn2RLkru6+8mqeufQfkeSXUmuS7I3yY+T3LzouAAAwCbRSZakiuAYSwTT3buyEkTNXrtj5nsnedcYYwEAAJyqRgmwAAAA1lPppdmDJcACAACmJ8ACAAAYiQALAABgBEtU5GKMMu0AAABEgAUAAJwE1b3QMff5VddU1TNVtbeqbluj/Zyquq+qHquqr1bV64brF1fVF6rq6ap6sqrePXPPH1XVd6pqz3BcN28elggCAADTm3APVlVtSfKxJL+RZF+Sh6tqZ3c/NdPtfUn2dPdbq+o1Q/83JzmU5L3d/WhVvTzJI1X1wMy9f9rdf3K8c5HBAgAAJtYrAdYix7FdkWRvdz/b3T9J8qkk16/qc1mSB5Oku7+R5JKquqC7D3T3o8P1HyZ5OslFL/UvFWABAACbwblVtXvm2DHTdlGS52bO9+Vng6SvJ/mdJKmqK5K8OsnW2Q5VdUmSNyT5yszlW4dlhXdV1TnzJinAAgAAptUZI4P1QndvnznunBmh1hl11geSnFNVe5L8qyRfy8rywJUHVL0syWeSvKe7Xxwu357kl5JcnuRAkg/O+1PtwQIAAKY3bZn2fUkunjnfmmT/bIchaLo5Saqqkvz9cKSqzshKcPXJ7v7szD3PH/1eVX+W5PPzJiKDBQAATG7iKoIPJ9lWVZdW1ZlJbkiy8/83ftXZQ1uS/F6Sh7r7xSHY+niSp7v7Q6vuuXDm9K1Jnpg3ERksAABgehNWEezuQ1V1a5L7k2xJcld3P1lV7xza70jy2iT3VNXhJE8luWW4/aokb0/y+LB8MEne1927kvxxVV2eleWG30ry+/PmIsACAAA2vSEg2rXq2h0z37+UZNsa930xa+/hSne//UTnIcACAACm1UmOTJfBOpUIsAAAgIkd17usTgsCLAAAYHoCLAAAgJEsSYClTDsAAMBIZLAAAIBpKXIBAAAwlk76yEZP4qQQYAEAANOzBwsAAIATIYMFAABMyx4sAACAES3JEkEBFgAAMD0BFgAAwBh6aQIsRS4AAABGIoMFAABMq5Mc8R4sAACAcSzJEkEBFgAAMD0BFgAAwBh6ad6DpcgFAADASGSwAACAaXXSrcgFAADAOJZkiaAACwAAmN6SFLmwBwsAAGAkMlgAAMC0ur1oGAAAYDRLskRQgAUAAEyuZbAAAADG0EuTwVLkAgAAYCQyWAAAwLQ63oMFAAAwmrYHCwAAYGGdpJckg7XQHqyqekVVPVBV3xw+z1mn37eq6vGq2lNVuxcZEwAA2GS6VzJYixxzVNU1VfVMVe2tqtvWaD+nqu6rqseq6qtV9bp59x5vvDNr0SIXtyV5sLu3JXlwOF/Pm7r78u7evuCYAAAAP1VVW5J8LMm1SS5LcmNVXbaq2/uS7Onuf5HkHUk+fBz3nki8k2TxAOv6JHcP3+9O8pYFnwcAAJyG+kgvdMxxRZK93f1sd/8kyaeyEqvMuiwrQVK6+xtJLqmqC+bce8LxzqIB1gXdfWCY5IEk56/Tr5P8dVU9UlU7FhwTAADYbKZdInhRkudmzvcN12Z9PcnvJElVXZHk1Um2zrn3eOOdn5pb5KKq/ibJK9do+rfz7p1xVXfvr6rzkzxQVd/o7ofWGW9HkqNB2I+q6pkTGGdM5yZ5YYPG5sT5vTYfv9nm4zfbXPxem89p+5vVe/71Rk9hKhv5m716g8Z9SX6Y/+v+v+n/7dwFH/MLq+o53Nnddw7fa43+q9NeH0jy4arak+TxJF9Lcug47z1ucwOs7v719dqq6vmqurC7D1TVhUkOrvOM/cPnwaq6LytpuDUDrOEf0p1rtZ1MVbXbfrHNw++1+fjNNh+/2ebi99p8/Gabj9/s+HX3NRMPsS/JxTPnW5PsXzWHF5PcnCRVVUn+fjh+8Rj3Hle8M2vRJYI7k9w0fL8pyedWd6iqs6rq5Ue/J/nNJE8sOC4AAMBRDyfZVlWXVtWZSW7ISqzyU1V19tCWJL+X5KEh6DrWvXPjndUWfQ/WB5J8uqpuSfLtJG8bJv+fJ/kP3X1dkguS3LcSJObnk/x5d//VguMCAIX/Td0AAAN7SURBVAAkSbr7UFXdmuT+JFuS3NXdT1bVO4f2O5K8Nsk9VXU4yVNJbjnWvcOj14x3jqW6l+OFXyeqqnbMrOnkFOf32nz8ZpuP32xz8XttPn6zzcdvxloEWAAAACNZdA8WAAAAAwHWOqrqbVX1ZFUdqSrVYU5hVXVNVT1TVXurau7btdlYVXVXVR2sKsVuNoGquriqvlBVTw//nfjujZ4Tx1ZVv1BVX62qrw+/2f+w0XNivqraUlVfq6rPb/RcmK+qvlVVj1fVnlVlw0GAdQxPZOVFZGuWk+fUUFVbknwsybVZeTv3jVV12cbOijk+kWTqUq2M51CS93b3a5NcmeRd/jN2yvunJP+yu1+f5PIk11TVlRs8J+Z7d5KnN3oSnJA3dfflyrSzmgBrHd39dHdv1EuOOX5XJNnb3c9290+SfCrJ9Rs8J45heMn49zd6Hhyf7j7Q3Y8O33+Ylf8DeNGx72Ij9YofDadnDIcN16ewqtqa5L9O8h82ei7A4gRYbHYXJXlu5nxf/J8/mERVXZLkDUm+srEzYZ5hudmerLwQ84Hu9pud2v6nJP99kiMbPRGOWyf566p6pKp2bPRkOLUs+h6sTa2q/ibJK9do+rfdPfclYpwSao1r/k0tjKyqXpbkM0neM7yUkVNYdx9OcnlVnZ2Vd1G+rrvtezwFVdVvJTnY3Y9U1dUbPR+O21Xdvb+qzk/yQFV9Y1ihAcsdYHX3r2/0HFjYviQXz5xvTbJ/g+YCp6WqOiMrwdUnu/uzGz0fjl93/99V9X9kZd+jAOvUdFWS366q65L8QpJ/VlX/S3f/Nxs8L46hu/cPnwer6r6sbFkQYJHEEkE2v4eTbKuqS6vqzCQ3JNm5wXOC00ZVVZKPJ3m6uz+00fNhvqo6b8hcpar+syS/nuQbGzsr1tPd/6a7t3b3JVn537D/XXB1aquqs6rq5Ue/J/nN+BcYzBBgraOq3lpV+5L8WpL/WFX3b/Sc+FndfSjJrUnuz8rm+09395MbOyuOparuTfKlJL9cVfuq6paNnhPHdFWStyf5l0M54j3Dv2nn1HVhki9U1WNZ+ZdQD3S30t8wnguSfLGqvp7kq0n+Y3f/1QbPiVNIdduuAgAAMAYZLAAAgJEIsAAAAEYiwAIAABiJAAsAAGAkAiwAAICRCLAAAABGIsACAAAYiQALAABgJP8vV+h5pWq2hLcAAAAASUVORK5CYII=\n",
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0IAAAFpCAYAAACxqpioAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAarUlEQVR4nO3dbaxlV3kf8P8TMxQwjkjrATt+wVRyoQmqDRk5RpYq85YYB8VtRSQjBVIUaZqUtFAhRZBI0H7rhwglBIQ1CpSgUFAKIbXAYEwCAqTgMOMMxmagOIiWiZ3YhmLjQKGe+/TDPaaX4b7NnLPv2efu30/amvOyzlkLjgz+z7PWs6u7AwAAMCU/tuwFAAAA7DVBCAAAmBxBCAAAmBxBCAAAmBxBCAAAmBxBCAAAmJy5g1BVXVJVn6iqE1V1d1W9ZpMx11bVQ1V1fHa9cd55AQAAHlNV/2GWR+6qqvdW1RO2G/+4Bcz5aJLXdfcdVXVekmNVdVt3f/G0cZ/u7pcuYD4AAIAfqKqLkvz7JD/V3d+tqj9OcmOSd231mbkrQt19X3ffMXv87SQnklw07/cCAACcgccleWJVPS7Jk5Lcu93ghZ4RqqrLkjwnye2bvP28qvp8VX2kqn56kfMCAADT1d1/k+R3kvyvJPcleai7P7bdZxaxNS5JUlVPTvKBJK/t7odPe/uOJE/v7keq6vokf5rk8i2+53CSw0ly7rnn/syznvWsRS2RJPm/dy17BZyh/3Hnk5a9BIBR+Sf/7DvLXgJn6sCzl72CfefYsWMPdvfBZa9jt37++ef2N755aq7vOHbn9+5O8n82vHSku48kSVX9RJIbkjwjybeS/Leq+uXu/qOtvq+6e64FzSY+kORDSW7t7jfvYvzXkhzq7ge3G3fo0KE+evTo3Ovj/1v7203zJyP28z955bKXADAqt957fNlL4Az92AVfWfYS9p2qOtbdh5a9jt36mSue0LffevFc33Hgwr/e8j9zVf1Skuu6+1dnz1+Z5Oru/rdbfd8iusZVknckObFVCKqqC2bjUlVXzeb9xrxzAwAAZH1L3NVV9aRZ7nhh1nsXbGkRW+OuSfKKJF+oqsf+iua3klyaJN19U5KXJfn1qno0yXeT3NiLKEUBAAAroHOq14b79u7bq+r9WT+S82iSv0pyZLvPzB2EuvszSWqHMW9N8tZ55wIAAFZPJ1nLsHWQ7n5TkjftdvzCmiUAAABsZS3DVYTOxkLbZwMAAKwCFSEAAGBQnc6pkbUIEIQAAIDBDX1G6EwJQgAAwKA6ySlBCAAAmJqxVYQ0SwAAACZHRQgAABhUJ5olAAAA0zOuuwgJQgAAwMA6rVkCAAAwMZ2cGlcO0iwBAACYHhUhAABgUB1nhAAAgMmpnEotexE/RBACAAAG1UnWnBECAABYLhUhAABgcLbGAQAAk9IRhAAAgAlaa0EIAACYkDFWhDRLAAAAJkdFCAAAGFSncmpkNRhBCAAAGJwzQgAAwKSM8YyQIAQAAAyscqrHtTVuXKsBAADYAypCAADAoDrJ2shqMIIQAAAwOGeEAACASel2RggAAGDpVIQAAIDBrdkaBwAATMn6fYTGtRlNEAIAAAY2vjNCghAAADCoMbbPHtdqAAAA9oCKEAAAMLhTrVkCAAAwIZ3SLAEAAJietZE1SxjXagAAgH3nsfbZ81w7qapnVtXxDdfDVfXarcarCAEAACuvu7+c5MokqapzkvxNkg9uNV4QAgAABtWpvW6W8MIkf93d/3OrAYIQAAAwuAXcR+j8qjq64fmR7j6yxdgbk7x3uy8ThAAAgEF1J6fmb5bwYHcf2mlQVT0+yS8mecN24zRLAAAA9pOXJLmju/9uu0EqQgAAwMAqa9mzM0Ivzw7b4hJBCAAAGFhnIVvjdlRVT0ry4iT/Zqexc6+mqi6pqk9U1YmquruqXrPJmKqqt1TVPVV1Z1U9d955AQCA1TH0fYSSpLu/093/qLsf2mnsIipCjyZ5XXffUVXnJTlWVbd19xc3jHlJkstn188mefvsTwAAYJ/rVNb2tn32juauCHX3fd19x+zxt5OcSHLRacNuSPLuXvfZJE+pqgvnnRsAAOBsLPSMUFVdluQ5SW4/7a2Lknx9w/OTs9fu2+Q7Dic5nCSXXnrpIpcHAAAsyW63t+2Vha2mqp6c5ANJXtvdD5/+9iYf6c2+p7uPdPeh7j508ODBRS0PAABYkk6y1j8217VoC6kIVdWBrIeg93T3n2wy5GSSSzY8vzjJvYuYGwAAGLvKqb1rn70ri+gaV0nekeREd795i2E3J3nlrHvc1Uke6u4f2RYHAADsP/u1InRNklck+UJVHZ+99ltJLk2S7r4pyS1Jrk9yT5LvJHnVAuYFAAA4K3MHoe7+TDY/A7RxTCd59bxzAQAAq2lsW+MW2jUOAADgdN01yPa2eQhCAADA4E6NLAiNazUAAAB7QEUIAAAYVCdZc0YIAACYlhrd1jhBCAAAGNT6fYRUhAAAgIk5NbL2BONaDQAAwB5QEQIAAAbVKVvjAACA6Vkb2WY0QQgAABhUd3JKRQgAAJiasW2NG1d9CgAAYA+oCAEAAINab5YwrhqMIAQAAAzuVMa1NU4QAgAABtVxRggAAGDpVIQAAICBOSMEAABM0JozQgAAwJS4oSoAADBJY9saN67VAAAA7AEVIQAAYFDrN1S1NQ4AAJgYzRIAAIBJcUNVAACAEVARAgAABje2rnGCEAAAMKzWLAEAAJiYzviaJYyrPgUAAOxLa7Oq0Nleu1FVT6mq91fVl6rqRFU9b6uxKkIAAMB+8XtJPtrdL6uqxyd50lYDBSEAAGBQe9E+u6p+PMk/T/Kvk6S7v5/k+1uNF4QAAIDB7UGzhH+c5IEk/6WqrkhyLMlruvvvNxvsjBAAADCoznzng2Yh6vyqOrrhOnzaNI9L8twkb+/u5yT5+ySv32pNKkIAAMDgFtA17sHuPrTN+yeTnOzu22fP359tgpCKEAAAsPK6+2+TfL2qnjl76YVJvrjVeBUhAABgWL0nZ4SS5N8lec+sY9xXk7xqq4GCEAAAMKi96BqXJN19PMl22+d+QBACAAAGt0cVoV1zRggAAJgcFSEAAGBQj7XPHhNBCAAAGFwLQgAAwNQs4D5CCyUIAQAAg+q9a5+9a5olAAAAk7OQIFRV76yq+6vqri3ev7aqHqqq47PrjYuYFwAAWA3dNde1aIvaGveuJG9N8u5txny6u1+6oPkAAICVsU+7xnX3p6rqskV8FwAAsP+MrWvcXp4Rel5Vfb6qPlJVP73VoKo6XFVHq+roAw88sIfLAwAAhtBZb5Ywz7VoexWE7kjy9O6+IsnvJ/nTrQZ295HuPtTdhw4ePLhHywMAAKZkT4JQdz/c3Y/MHt+S5EBVnb8XcwMAAEvW6y2057kWbU/uI1RVFyT5u+7uqroq6wHsG3sxNwAAsHz78oaqVfXeJNcmOb+qTiZ5U5IDSdLdNyV5WZJfr6pHk3w3yY3dQ+Q6AABgbDrja5awqK5xL9/h/bdmvb02AADA0u3J1jgAAGDK9ul9hAAAALYztoMxghAAADC4fXlGCAAAYCvrLbDHFYT26oaqAAAAo6EiBAAADE6zBAAAYHI0SwAAACZnbGeEBCEAAGBQnRpdENIsAQAAmBwVIQAAYHAjOyIkCAEAAAMb4X2EBCEAAGB4IysJOSMEAABMjooQAAAwOFvjAACAyXFDVQAAYFI6KkIAAMDUdJKRBSHNEgAAgMlREQIAAAbnjBAAADA9ghAAADAtpVkCAAAwQXtQEaqqryX5dpJTSR7t7kNbjRWEAACA/eT53f3gToMEIQAAYFg9vvsIaZ8NAAAMr+e8kvOr6uiG6/AWs3ysqo5t8f4PqAgBAAB7YO6K0IPbnfmZuaa7762qpya5raq+1N2f2mygihAAALAvdPe9sz/vT/LBJFdtNVYQAgAAhjf/1rhtVdW5VXXeY4+T/FySu7Yab2scAAAwvOHbZz8tyQerKlnPOf+1uz+61WBBCAAAGFYnGbhrXHd/NckVux0vCAEAAIPrPbih6plwRggAAJgcFSEAAGB4I6sICUIAAMDwBj4jdKYEIQAAYHClIgQAAEzKLu8FtJc0SwAAACZHRQgAABhYOSMEAABM0Mi2xglCAADA8EYWhJwRAgAAJkdFCAAAGN7IKkKCEAAAMKyOZgkAAMD0jO2Gqgs5I1RV76yq+6vqri3er6p6S1XdU1V3VtVzFzEvAACwInrOa8EW1SzhXUmu2+b9lyS5fHYdTvL2Bc0LAABwxhYShLr7U0m+uc2QG5K8u9d9NslTqurCRcwNAABwpvbqjNBFSb6+4fnJ2Wv37dH8AADAEo3tjNBeBaHNWkRs+l9FVR3O+va5XHrppUOuCQAA2Csj6xq3VzdUPZnkkg3PL05y72YDu/tIdx/q7kMHDx7ck8UBAAADmrdRwoibJezk5iSvnHWPuzrJQ91tWxwAALAUC9kaV1XvTXJtkvOr6mSSNyU5kCTdfVOSW5Jcn+SeJN9J8qpFzAsAAKyI/XhGqLtfvsP7neTVi5gLAABYPVNtlgAAAEzZyILQXp0RAgAAGA0VIQAAYHgjqwgJQgAAwKCqnRECAACmaGQ3VBWEAACA4Y2sIqRZAgAAMDkqQgAAwOCcEQIAAKZHEAIAACZlhF3jnBECAAAmR0UIAAAY3sgqQoIQAAAwPEEIAACYGmeEAAAAlkwQAgAA9oWqOqeq/qqqPrTTWFvjAACA4e3N1rjXJDmR5Md3GqgiBAAADGt2H6F5rp1U1cVJfiHJH+xmSSpCAADA8OavCJ1fVUc3PD/S3Uc2PP/dJL+Z5LzdfJkgBAAADG/+IPRgdx/a7I2qemmS+7v7WFVdu5svszUOAABYddck+cWq+lqS9yV5QVX90XYfEIQAAIBBVYY9I9Tdb+jui7v7siQ3Jvnz7v7l7T5jaxwAADC8kd1QVRACAACGtcvObwuZqvuTST650zhb4wAAgMlREQIAAIZnaxwAADA5ghAAADA1e3VGaLcEIQAAYHgjC0KaJQAAAJOjIgQAAAyrM7qKkCAEAAAMzhkhAABgegQhAABgasZWEdIsAQAAmBwVIQAAYHgjqwgJQgAAwLB0jQMAAKamZteYOCMEAABMjooQAAAwPFvjAACAqRlb+2xBCAAAGJ4gBAAATM7IgpBmCQAAwOSoCAEAAMNqZ4QAAIApGlkQWsjWuKq6rqq+XFX3VNXrN3n/2qp6qKqOz643LmJeAABgNVTPdy3a3BWhqjonyduSvDjJySSfq6qbu/uLpw39dHe/dN75AACAFbQPK0JXJbmnu7/a3d9P8r4kNyzgewEAAAaxiCB0UZKvb3h+cvba6Z5XVZ+vqo9U1U9v9WVVdbiqjlbV0QceeGABywMAAJZtbFvjFhGEapPXTl/qHUme3t1XJPn9JH+61Zd195HuPtTdhw4ePLiA5QEAAEvVC7gWbBFB6GSSSzY8vzjJvRsHdPfD3f3I7PEtSQ5U1fkLmBsAAFgF+zAIfS7J5VX1jKp6fJIbk9y8cUBVXVBVNXt81WzebyxgbgAAgDM2d9e47n60qn4jya1Jzknyzu6+u6p+bfb+TUleluTXq+rRJN9NcmN3j6xvBAAAMITKPr2h6my72y2nvXbThsdvTfLWRcwFAACsoP0YhAAAALZTI9sQJggBAADDGqjhwTwW0SwBAABgpagIAQAAg9uXzRIAAAC2JQgBAABToyIEAABMz4BBqKqekORTSf5B1jPO+7v7Tdt9RhACAABW3feSvKC7H6mqA0k+U1Uf6e7PbvUBQQgAABhWD7s1rrs7ySOzpwdm17Yzap8NAAAMr+e8dlBV51TV8ST3J7mtu2/fbrwgBAAADKqyXhGa50pyflUd3XAd3jhHd5/q7iuTXJzkqqp69nZrsjUOAABYBQ9296GdBnX3t6rqk0muS3LXVuNUhAAAgOF1z3dto6oOVtVTZo+fmORFSb603WdUhAAAgMENfB+hC5P8YVWdk/Vizx9394e2+4AgBAAADGuXDQ/O+uu770zynDP5jCAEAAAMrtaWvYIf5owQAAAwOSpCAADA8IY9I3TGBCEAAGBwAzdLOGOCEAAAMKzOji2w95ogBAAADG5sFSHNEgAAgMlREQIAAIY3soqQIAQAAAyqMr6tcYIQAAAwrO7RNUtwRggAAJgcFSEAAGBwtsYBAADTIwgBAABToyIEAABMSydZG1cS0iwBAACYHBUhAABgeOMqCAlCAADA8JwRAgAApmdkN1QVhAAAgMGNrSKkWQIAADA5KkIAAMCwOpolAAAA01JJyhkhAABgctaWvYAf5owQAAAwOSpCAADA4GyNAwAApkWzBAAAYHraDVUBAIDpcUNVAACAJVtIEKqq66rqy1V1T1W9fpP3q6reMnv/zqp67iLmBQAAVkT3fNeCzb01rqrOSfK2JC9OcjLJ56rq5u7+4oZhL0ly+ez62SRvn/0JAADsd53UPryP0FVJ7unur3b395O8L8kNp425Icm7e91nkzylqi5cwNwAAMAqGFlFaBFB6KIkX9/w/OTstTMdkySpqsNVdbSqjj7wwAMLWB4AAMAPW0QQqk1eOz2y7WbM+ovdR7r7UHcfOnjw4NyLAwAARqDnvBZsEe2zTya5ZMPzi5PcexZjAACAfapGdh+hRVSEPpfk8qp6RlU9PsmNSW4+bczNSV456x53dZKHuvu+BcwNAACsgpGdEZq7ItTdj1bVbyS5Nck5Sd7Z3XdX1a/N3r8pyS1Jrk9yT5LvJHnVvPMCAAAropOMrGvcIrbGpbtvyXrY2fjaTRsed5JXL2IuAACA01XVJUneneSCrMeuI939e1uNX0gQAgAA2Eql9+KM0KNJXtfdd1TVeUmOVdVtp93f9AcEIQAAYHgDB6FZD4L7Zo+/XVUnsn7LHkEIAABYkvmD0PlVdXTD8yPdfWSzgVV1WZLnJLl9qy8ThAAAgGEtplnCg919aKdBVfXkJB9I8trufnircYtonw0AALB0VXUg6yHoPd39J9uNVRECAAAGN3SzhKqqJO9IcqK737zTeBUhAABgeMPfUPWaJK9I8oKqOj67rt9qsIoQAAAwsF2HmbOfofszSWq341WEAACAyVERAgAAhtUZvCJ0pgQhAABgePO3z14oQQgAABjc0F3jzpQgBAAADG9kQUizBAAAYHJUhAAAgGF1krVxVYQEIQAAYGDD30foTAlCAADA8AQhAABgckYWhDRLAAAAJkdFCAAAGJZmCQAAwPR00mvLXsQPEYQAAIDhOSMEAACwXCpCAADAsJwRAgAAJmlkW+MEIQAAYHiCEAAAMC09uiCkWQIAADA5KkIAAMCwOsma+wgBAABTM7KtcYIQAAAwPEEIAACYlh7dfYQ0SwAAACZHRQgAABhWJ92aJQAAAFMzsq1xghAAADC8kTVLcEYIAACYHBUhAABgWN1uqAoAAEzQyLbGCUIAAMDgWkUIAACYlh5dRUizBAAAYHJUhAAAgGF13EcIAACYoHZGCAAAmJBO0iOrCM11Rqiq/mFV3VZVX5n9+RNbjPtaVX2hqo5X1dF55gQAAFZM93pFaJ5rB1X1zqq6v6ru2s2S5m2W8Pokf9bdlyf5s9nzrTy/u6/s7kNzzgkAAHC6dyW5breD5w1CNyT5w9njP0zyL+b8PgAAYB/qtZ7r2vH7uz+V5Ju7Xc+8Qehp3X3fbOL7kjx1q3Ul+VhVHauqw3POCQAArJqBt8adqR2bJVTVx5NcsMlbv30G81zT3fdW1VOT3FZVX5olts3mO5zksbD0SFV9+QzmWaTzkzy4pLk5c/v497pn2QsYyj7+zfYtv9lq2be/1zkXLnsFg9m3v1lSy17AUJb5mz19SfOelW/nf9/68X7/+XN+zRNO6zdwpLuPnO2XVc9xh9dZSLm2u++rqguTfLK7n7nDZ/5jkke6+3fOeuI9UFVHnWdaHX6v1eM3Wz1+s9Xi91o9frPV4zcbn6q6LMmHuvvZO42dd2vczUl+Zfb4V5L8900Wc25VnffY4yQ/l2RXnRwAAACGMG8Q+s9JXlxVX0ny4tnzVNVPVtUtszFPS/KZqvp8kr9M8uHu/uic8wIAAPxAVb03yV8keWZVnayqX91u/Fw3VO3ubyR54Sav35vk+tnjrya5Yp55luSs9xuyFH6v1eM3Wz1+s9Xi91o9frPV4zcbke5++ZmMn+uMEAAAwCqad2scAADAyhGEtlBVv1RVd1fVWlXpBjJiVXVdVX25qu6pqtcvez1sr6reWVX3V5WmKSugqi6pqk9U1YnZ/ya+ZtlrYntV9YSq+suq+vzsN/tPy14TO6uqc6rqr6rqQ8teCzurqq9V1Req6vhp7ZxZIYLQ1u5K8q+SbHq/I8ahqs5J8rYkL0nyU0leXlU/tdxVsYN3Jblu2Ytg1x5N8rru/qdJrk7yav+Mjd73krygu69IcmWS66rq6iWviZ29JsmJZS+CM/L87r5S++zVJQhtobtPdPeybubK7l2V5J7u/mp3fz/J+5LcsOQ1sY3ZzZS/uex1sDvdfV933zF7/O2s/4vaRctdFdvpdY/Mnh6YXQ4Ej1hVXZzkF5L8wbLXAlMiCLHqLkry9Q3PT8a/pMEgZjepe06S25e7EnYy22Z1PMn9SW7rbr/ZuP1ukt9MsrbshbBrneRjVXWsqg4vezGcnbnaZ6+6qvp4kgs2eeu3u/tHbg7LKNUmr/mbT1iwqnpykg8keW13P7zs9bC97j6V5MqqekqSD1bVs7vbubwRqqqXJrm/u49V1bXLXg+7dk1331tVT01yW1V9abbjgRUy6SDU3S9a9hqY28kkl2x4fnGSe5e0FtiXqupA1kPQe7r7T5a9Hnavu79VVZ/M+rk8QWicrknyi1V1fZInJPnxqvqj7v7lJa+LbczumZnuvr+qPpj1rfqC0IqxNY5V97kkl1fVM6rq8UluTHLzktcE+0ZVVZJ3JDnR3W9e9nrYWVUdnFWCUlVPTPKiJF9a7qrYSne/obsv7u7Lsv7/YX8uBI1bVZ1bVec99jjJz8VfNKwkQWgLVfUvq+pkkucl+XBV3brsNfGjuvvRJL+R5NasH+L+4+6+e7mrYjtV9d4kf5HkmVV1sqp+ddlrYlvXJHlFkhfM2sQen/3NNeN1YZJPVNWdWf/Lotu6W0tmWJynJflMVX0+yV8m+XB3f3TJa+IsVLfjFAAAwLSoCAEAAJMjCAEAAJMjCAEAAJMjCAEAAJMjCAEAAJMjCAEAAJMjCAEAAJMjCAEAAJPz/wC3LLIpRluZ6AAAAABJRU5ErkJggg==\n",
       "text/plain": [
        "<Figure size 1152x432 with 2 Axes>"
       ]
@@ -171,13 +171,13 @@
     "    flag_arr[-1, :] = 0   \n",
     "    flag_arr[:, 0] = flags[noslip]\n",
     "    flag_arr[:, -1] = flags[noslip]\n",
-    "#else:\n",
+    "else:\n",
     "#    flag_arr[0, :] = 0\n",
     "#    flag_arr[-1, :] = 0   \n",
     "#    flag_arr[:, 0] = 0\n",
     "#    flag_arr[:, -1] = 0\n",
-    "#    flag_arr[0, :] = flags[noslip]\n",
-    "#    flag_arr[-1, :] = flags[density]\n",
+    "    flag_arr[0, :] = flags[noslip]\n",
+    "    flag_arr[-1, :] = flags[noslip]\n",
     "#    flag_arr[:, -1] = flags[noslip]\n",
     "#    flag_arr[:, 0] = flags[noslip]\n",
     "\n",
@@ -227,11 +227,11 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "flag: [[1 1 1]\n",
+      "flag: [[8 8 8]\n",
       " [1 1 1]\n",
       " [1 1 1]\n",
       " [1 1 1]\n",
-      " [1 1 1]]\n",
+      " [8 8 8]]\n",
       "dir: 0\n",
       "dir: 1\n",
       "dir: 2\n",
@@ -241,7 +241,7 @@
       "dir: 6\n",
       "dir: 7\n",
       "dir: 8\n",
-      "number of fluid: 15 counter: 12\n"
+      "number of fluid: 9 counter: 6\n"
      ]
     }
    ],
@@ -300,288 +300,129 @@
      "output_type": "stream",
      "text": [
       "domain_size: (5, 3)\n",
+      "1\n",
+      "2\n",
+      "3\n",
       "\n",
       " New direction: (0, 1) ,  1\n",
       "pos is  0\n",
       "(inner:)\n",
-      "third\n",
-      "write: [1, 2] read: [1, 1]\n",
-      "write: [2, 2] read: [2, 1]\n",
-      "write: [3, 2] read: [3, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 0] read: [0, 2]\n",
-      "write: [1, 0] read: [1, 2]\n",
-      "write: [2, 0] read: [2, 2]\n",
-      "write: [3, 0] read: [3, 2]\n",
-      "write: [4, 0] read: [4, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 2] read: [0, 1]\n",
-      "write: [0, 1] read: [0, 0]\n",
-      "write: [4, 2] read: [4, 1]\n",
-      "write: [4, 1] read: [4, 0]\n",
       "\n",
       " New direction: (0, -1) ,  2\n",
       "pos is  0\n",
       "(inner:)\n",
-      "third\n",
-      "write: [1, 0] read: [1, 1]\n",
-      "write: [2, 0] read: [2, 1]\n",
-      "write: [3, 0] read: [3, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 2] read: [0, 0]\n",
-      "write: [1, 2] read: [1, 0]\n",
-      "write: [2, 2] read: [2, 0]\n",
-      "write: [3, 2] read: [3, 0]\n",
-      "write: [4, 2] read: [4, 0]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 0] read: [0, 1]\n",
-      "write: [0, 1] read: [0, 2]\n",
-      "write: [4, 0] read: [4, 1]\n",
-      "write: [4, 1] read: [4, 2]\n",
       "\n",
       " New direction: (-1, 0) ,  3\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [4, 0] read: [0, 0]\n",
-      "write: [4, 1] read: [0, 1]\n",
-      "write: [4, 2] read: [0, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 0] read: [1, 0]\n",
-      "write: [1, 0] read: [2, 0]\n",
-      "write: [2, 0] read: [3, 0]\n",
-      "write: [3, 0] read: [4, 0]\n",
-      "write: [0, 2] read: [1, 2]\n",
-      "write: [1, 2] read: [2, 2]\n",
-      "write: [2, 2] read: [3, 2]\n",
-      "write: [3, 2] read: [4, 2]\n",
       "pos is  1\n",
       "(inner:)\n",
-      "third\n",
-      "write: [0, 1] read: [1, 1]\n",
       "\n",
       " New direction: (1, 0) ,  4\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 0] read: [4, 0]\n",
-      "write: [0, 1] read: [4, 1]\n",
-      "write: [0, 2] read: [4, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [4, 0] read: [3, 0]\n",
-      "write: [3, 0] read: [2, 0]\n",
-      "write: [2, 0] read: [1, 0]\n",
-      "write: [1, 0] read: [0, 0]\n",
-      "write: [4, 2] read: [3, 2]\n",
-      "write: [3, 2] read: [2, 2]\n",
-      "write: [2, 2] read: [1, 2]\n",
-      "write: [1, 2] read: [0, 2]\n",
       "pos is  1\n",
       "(inner:)\n",
-      "third\n",
-      "write: [4, 1] read: [3, 1]\n",
       "\n",
       " New direction: (-1, 1) ,  5\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [4, 1] read: [0, 0]\n",
-      "write: [4, 2] read: [0, 1]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 2] read: [1, 1]\n",
-      "write: [1, 2] read: [2, 1]\n",
-      "write: [2, 2] read: [3, 1]\n",
-      "write: [3, 2] read: [4, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 0] read: [1, 2]\n",
-      "write: [1, 0] read: [2, 2]\n",
-      "write: [2, 0] read: [3, 2]\n",
-      "write: [3, 0] read: [4, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 2] read: [1, 1]\n",
-      "write: [0, 1] read: [1, 0]\n",
-      "(Ecke) write: [4, 0] read: [0, 2]\n",
       "\n",
       " New direction: (1, 1) ,  6\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 1] read: [4, 0]\n",
-      "write: [0, 2] read: [4, 1]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [4, 2] read: [3, 1]\n",
-      "write: [3, 2] read: [2, 1]\n",
-      "write: [2, 2] read: [1, 1]\n",
-      "write: [1, 2] read: [0, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [1, 0] read: [0, 2]\n",
-      "write: [2, 0] read: [1, 2]\n",
-      "write: [3, 0] read: [2, 2]\n",
-      "write: [4, 0] read: [3, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [4, 2] read: [3, 1]\n",
-      "write: [4, 1] read: [3, 0]\n",
-      "(Ecke) write: [0, 0] read: [4, 2]\n",
       "\n",
       " New direction: (-1, -1) ,  7\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [4, 0] read: [0, 1]\n",
-      "write: [4, 1] read: [0, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 0] read: [1, 1]\n",
-      "write: [1, 0] read: [2, 1]\n",
-      "write: [2, 0] read: [3, 1]\n",
-      "write: [3, 0] read: [4, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 2] read: [1, 0]\n",
-      "write: [1, 2] read: [2, 0]\n",
-      "write: [2, 2] read: [3, 0]\n",
-      "write: [3, 2] read: [4, 0]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [0, 0] read: [1, 1]\n",
-      "write: [0, 1] read: [1, 2]\n",
-      "(Ecke) write: [4, 2] read: [0, 0]\n",
       "\n",
       " New direction: (1, -1) ,  8\n",
       "pos is  0\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [0, 0] read: [4, 1]\n",
-      "write: [0, 1] read: [4, 2]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [4, 0] read: [3, 1]\n",
-      "write: [3, 0] read: [2, 1]\n",
-      "write: [2, 0] read: [1, 1]\n",
-      "write: [1, 0] read: [0, 1]\n",
       "pos is  1\n",
       "(periodic:)\n",
-      "first\n",
-      "write: [1, 2] read: [0, 0]\n",
-      "write: [2, 2] read: [1, 0]\n",
-      "write: [3, 2] read: [2, 0]\n",
-      "write: [4, 2] read: [3, 0]\n",
       "(inner:)\n",
-      "second\n",
-      "write: [4, 0] read: [3, 1]\n",
-      "write: [4, 1] read: [3, 2]\n",
-      "(Ecke) write: [0, 2] read: [4, 0]\n",
-      "[[1, 15, 17], [1, 18, 20], [1, 21, 23], [1, 24, 26], [1, 27, 29]]\n",
-      "[[2, 32, 30], [2, 35, 33], [2, 38, 36], [2, 41, 39], [2, 44, 42]]\n",
-      "[[3, 57, 45], [3, 58, 46], [3, 59, 47]]\n",
-      "[[4, 60, 72], [4, 61, 73], [4, 62, 74]]\n",
-      "[[5, 88, 75], [5, 89, 76]]\n",
-      "[[5, 75, 80], [5, 78, 83], [5, 81, 86], [5, 84, 89]]\n",
-      "[[5, 75, 80], [5, 78, 83], [5, 81, 86], [5, 84, 89], [5, 87, 77]]\n",
-      "[[6, 91, 102], [6, 92, 103]]\n",
-      "[[6, 93, 92], [6, 96, 95], [6, 99, 98], [6, 102, 101]]\n",
-      "[[6, 93, 92], [6, 96, 95], [6, 99, 98], [6, 102, 101], [6, 90, 104]]\n",
-      "[[7, 117, 106], [7, 118, 107]]\n",
-      "[[7, 107, 108], [7, 110, 111], [7, 113, 114], [7, 116, 117]]\n",
-      "[[7, 107, 108], [7, 110, 111], [7, 113, 114], [7, 116, 117], [7, 119, 105]]\n",
-      "[[8, 120, 133], [8, 121, 134]]\n",
-      "[[8, 125, 120], [8, 128, 123], [8, 131, 126], [8, 134, 129]]\n",
-      "[[8, 125, 120], [8, 128, 123], [8, 131, 126], [8, 134, 129], [8, 122, 132]]\n",
-      "[[1, 20, 19], [1, 23, 22], [1, 26, 25], [1, 17, 16], [1, 16, 15], [1, 29, 28], [1, 28, 27], [2, 33, 34], [2, 36, 37], [2, 39, 40], [2, 30, 31], [2, 31, 32], [2, 42, 43], [2, 43, 44], [3, 45, 48], [3, 48, 51], [3, 51, 54], [3, 54, 57], [3, 47, 50], [3, 50, 53], [3, 53, 56], [3, 56, 59], [3, 46, 49], [4, 72, 69], [4, 69, 66], [4, 66, 63], [4, 63, 60], [4, 74, 71], [4, 71, 68], [4, 68, 65], [4, 65, 62], [4, 73, 70], [5, 77, 79], [5, 80, 82], [5, 83, 85], [5, 86, 88], [5, 77, 79], [5, 76, 78], [6, 104, 100], [6, 101, 97], [6, 98, 94], [6, 95, 91], [6, 104, 100], [6, 103, 99], [7, 105, 109], [7, 108, 112], [7, 111, 115], [7, 114, 118], [7, 105, 109], [7, 106, 110], [8, 132, 130], [8, 129, 127], [8, 126, 124], [8, 123, 121], [8, 132, 130], [8, 133, 131]]\n"
+      "[[1, 9, 11], [1, 12, 14], [1, 15, 17]]\n",
+      "[[2, 20, 18], [2, 23, 21], [2, 26, 24]]\n",
+      "[]\n",
+      "[]\n",
+      "[]\n",
+      "[[5, 45, 50], [5, 48, 53], [5, 51, 123]]\n",
+      "[]\n",
+      "[[6, 54, 117], [6, 57, 56], [6, 60, 59]]\n",
+      "[]\n",
+      "[[7, 65, 66], [7, 68, 69], [7, 71, 125]]\n",
+      "[]\n",
+      "[[8, 74, 119], [8, 77, 72], [8, 80, 75]]\n",
+      "[[1, [1, 2], [1, 1]], [1, [2, 2], [2, 1]], [1, [3, 2], [3, 1]], [2, [1, 0], [1, 1]], [2, [2, 0], [2, 1]], [2, [3, 0], [3, 1]], [3, 27, 30], [3, 30, 33], [3, 33, 69], [3, 29, 32], [3, 32, 35], [3, 35, 71], [4, 42, 39], [4, 39, 36], [4, 36, 63], [4, 44, 41], [4, 41, 38], [4, 38, 65], [5, 47, 49], [5, 50, 52], [5, 53, 125], [6, 62, 58], [6, 59, 55], [6, 56, 119], [7, 63, 67], [7, 66, 70], [7, 69, 123], [8, 78, 76], [8, 75, 73], [8, 72, 117]]\n"
      ]
     },
     {
      "data": {
       "text/plain": [
-       "[[[1, 15, 17], [1, 18, 20], [1, 21, 23], [1, 24, 26], [1, 27, 29]],\n",
-       " [[2, 32, 30], [2, 35, 33], [2, 38, 36], [2, 41, 39], [2, 44, 42]],\n",
-       " [[3, 57, 45], [3, 58, 46], [3, 59, 47]],\n",
-       " [[4, 60, 72], [4, 61, 73], [4, 62, 74]],\n",
-       " [[5, 88, 75], [5, 89, 76]],\n",
-       " [[5, 75, 80], [5, 78, 83], [5, 81, 86], [5, 84, 89]],\n",
-       " [[5, 75, 80], [5, 78, 83], [5, 81, 86], [5, 84, 89], [5, 87, 77]],\n",
-       " [[6, 91, 102], [6, 92, 103]],\n",
-       " [[6, 93, 92], [6, 96, 95], [6, 99, 98], [6, 102, 101]],\n",
-       " [[6, 93, 92], [6, 96, 95], [6, 99, 98], [6, 102, 101], [6, 90, 104]],\n",
-       " [[7, 117, 106], [7, 118, 107]],\n",
-       " [[7, 107, 108], [7, 110, 111], [7, 113, 114], [7, 116, 117]],\n",
-       " [[7, 107, 108], [7, 110, 111], [7, 113, 114], [7, 116, 117], [7, 119, 105]],\n",
-       " [[8, 120, 133], [8, 121, 134]],\n",
-       " [[8, 125, 120], [8, 128, 123], [8, 131, 126], [8, 134, 129]],\n",
-       " [[8, 125, 120], [8, 128, 123], [8, 131, 126], [8, 134, 129], [8, 122, 132]],\n",
-       " [[1, 20, 19],\n",
-       "  [1, 23, 22],\n",
-       "  [1, 26, 25],\n",
-       "  [1, 17, 16],\n",
-       "  [1, 16, 15],\n",
-       "  [1, 29, 28],\n",
-       "  [1, 28, 27],\n",
-       "  [2, 33, 34],\n",
-       "  [2, 36, 37],\n",
-       "  [2, 39, 40],\n",
-       "  [2, 30, 31],\n",
-       "  [2, 31, 32],\n",
-       "  [2, 42, 43],\n",
-       "  [2, 43, 44],\n",
-       "  [3, 45, 48],\n",
-       "  [3, 48, 51],\n",
-       "  [3, 51, 54],\n",
-       "  [3, 54, 57],\n",
-       "  [3, 47, 50],\n",
-       "  [3, 50, 53],\n",
-       "  [3, 53, 56],\n",
-       "  [3, 56, 59],\n",
-       "  [3, 46, 49],\n",
-       "  [4, 72, 69],\n",
-       "  [4, 69, 66],\n",
-       "  [4, 66, 63],\n",
-       "  [4, 63, 60],\n",
-       "  [4, 74, 71],\n",
-       "  [4, 71, 68],\n",
-       "  [4, 68, 65],\n",
-       "  [4, 65, 62],\n",
-       "  [4, 73, 70],\n",
-       "  [5, 77, 79],\n",
-       "  [5, 80, 82],\n",
-       "  [5, 83, 85],\n",
-       "  [5, 86, 88],\n",
-       "  [5, 77, 79],\n",
-       "  [5, 76, 78],\n",
-       "  [6, 104, 100],\n",
-       "  [6, 101, 97],\n",
-       "  [6, 98, 94],\n",
-       "  [6, 95, 91],\n",
-       "  [6, 104, 100],\n",
-       "  [6, 103, 99],\n",
-       "  [7, 105, 109],\n",
-       "  [7, 108, 112],\n",
-       "  [7, 111, 115],\n",
-       "  [7, 114, 118],\n",
-       "  [7, 105, 109],\n",
-       "  [7, 106, 110],\n",
-       "  [8, 132, 130],\n",
-       "  [8, 129, 127],\n",
-       "  [8, 126, 124],\n",
-       "  [8, 123, 121],\n",
-       "  [8, 132, 130],\n",
-       "  [8, 133, 131]]]"
+       "[[[1, 9, 11], [1, 12, 14], [1, 15, 17]],\n",
+       " [[2, 20, 18], [2, 23, 21], [2, 26, 24]],\n",
+       " [],\n",
+       " [],\n",
+       " [],\n",
+       " [[5, 45, 50], [5, 48, 53], [5, 51, 123]],\n",
+       " [],\n",
+       " [[6, 54, 117], [6, 57, 56], [6, 60, 59]],\n",
+       " [],\n",
+       " [[7, 65, 66], [7, 68, 69], [7, 71, 125]],\n",
+       " [],\n",
+       " [[8, 74, 119], [8, 77, 72], [8, 80, 75]],\n",
+       " [[1, [1, 2], [1, 1]],\n",
+       "  [1, [2, 2], [2, 1]],\n",
+       "  [1, [3, 2], [3, 1]],\n",
+       "  [2, [1, 0], [1, 1]],\n",
+       "  [2, [2, 0], [2, 1]],\n",
+       "  [2, [3, 0], [3, 1]],\n",
+       "  [3, 27, 30],\n",
+       "  [3, 30, 33],\n",
+       "  [3, 33, 69],\n",
+       "  [3, 29, 32],\n",
+       "  [3, 32, 35],\n",
+       "  [3, 35, 71],\n",
+       "  [4, 42, 39],\n",
+       "  [4, 39, 36],\n",
+       "  [4, 36, 63],\n",
+       "  [4, 44, 41],\n",
+       "  [4, 41, 38],\n",
+       "  [4, 38, 65],\n",
+       "  [5, 47, 49],\n",
+       "  [5, 50, 52],\n",
+       "  [5, 53, 125],\n",
+       "  [6, 62, 58],\n",
+       "  [6, 59, 55],\n",
+       "  [6, 56, 119],\n",
+       "  [7, 63, 67],\n",
+       "  [7, 66, 70],\n",
+       "  [7, 69, 123],\n",
+       "  [8, 78, 76],\n",
+       "  [8, 75, 73],\n",
+       "  [8, 72, 117]]]"
       ]
      },
      "execution_count": 10,
@@ -692,12 +533,12 @@
      "output_type": "stream",
      "text": [
       "domain: (5, 3)\n",
-      "pdf size in bytes: 1080\n",
-      "pdf size: 135\n",
+      "pdf size in bytes: 648\n",
+      "pdf size: 81\n",
       "index array size in bytes: 96\n",
       "density index array size in bytes: 0\n",
       "ubb index array size in bytes: 0\n",
-      "sum: 1176\n"
+      "sum: 744\n"
      ]
     }
    ],
@@ -856,7 +697,7 @@
     {
      "data": {
       "text/plain": [
-       "<matplotlib.quiver.Quiver at 0x7f52e08f38e0>"
+       "<matplotlib.quiver.Quiver at 0x7f943ffa12e0>"
       ]
      },
      "execution_count": 20,
-- 
GitLab