Correct optimization criteria to improve mesh quality; jit-ify some slow stuff

Author: FranzBangar

.vscode/launch.json

@@ -6,7 +6,7 @@
     "configurations": [
         {
             "name": "Python: classy example",
-            "type": "python",
+            "type": "debugpy",
             "request": "launch",
             "program": "${file}",
             "cwd": "${fileDirname}",

classy_blocks.code-workspace

@@ -6,13 +6,6 @@
     ],
     "settings": {
         "python.languageServer": "Pylance",
-        "cSpell.words": [
-            "astype",
-            "linspace",
-            "Pylance",
-            "scipy"
-        ],
-        "editor.autoClosingBrackets": "never",
         "editor.defaultFormatter": "charliermarsh.ruff",
         "ruff.organizeImports": true
     },

src/classy_blocks/assemble/assembler.py

@@ -1,7 +1,5 @@
 from typing import Tuple, get_args
 
-import numpy as np
-
 from classy_blocks.assemble.depot import Depot
 from classy_blocks.assemble.dump import AssembledDump
 from classy_blocks.assemble.settings import Settings
@@ -15,25 +13,26 @@
 from classy_blocks.lists.patch_list import PatchList
 from classy_blocks.lists.vertex_list import VertexList
 from classy_blocks.lookup.point_registry import HexPointRegistry
-from classy_blocks.optimize.grid import HexGrid
 from classy_blocks.util import constants
 
 
 class MeshAssembler:
-    def __init__(self, depot: Depot, settings: Settings):
+    def __init__(self, depot: Depot, settings: Settings, merge_tol=constants.TOL):
         self.depot = depot
         self.settings = settings
+        self.merge_tol = merge_tol
 
         # once the mesh is assembled, adding new stuff to depot will break things;
         # better (and faster) is to cache status quo
         self._operations = self.depot.operations
-        self._grid = HexGrid.from_elements(self._operations)
+
+        self._points = HexPointRegistry.from_operations(self._operations, self.merge_tol)
 
     def _create_blocks(self, vertex_list: VertexList) -> BlockList:
         block_list = BlockList()
 
         for iop, operation in enumerate(self._operations):
-            op_indexes = self._grid.addressing[iop]
+            op_indexes = self._points.cell_addressing[iop]
             op_vertices = [vertex_list.vertices[i] for i in op_indexes]
 
             # duplicate vertices on slave patches
@@ -111,14 +110,11 @@ def _create_patches(self, block_list: BlockList) -> Tuple[PatchList, FaceList]:
         return patch_list, face_list
 
     def _update_neighbours(self, block_list: BlockList) -> None:
-        points = np.array([[v.position for v in b.vertices] for b in block_list.blocks])
-
-        navigator = HexPointRegistry(HexPointRegistry.flatten(points, 8 * len(block_list.blocks)), constants.TOL)
-        block_list.update_neighbours(navigator)
+        block_list.update_neighbours(self._points)
 
     def assemble(self) -> AssembledDump:
         # Create reused/indexes vertices from operations' points
-        vertex_list = VertexList([Vertex(pos, i) for i, pos in enumerate(self._grid.points)])
+        vertex_list = VertexList([Vertex(pos, i) for i, pos in enumerate(self._points.unique_points)])
         # Create blocks from vertices; when there's a slave patch specified in an operation,
         # duplicate vertices for that patch
         block_list = self._create_blocks(vertex_list)

src/classy_blocks/items/wires/axis.py

@@ -27,6 +27,7 @@ def add_neighbour(self, other: "Axis") -> None:
             for nei_wire in other.wires:
                 if this_wire.is_coincident(nei_wire):
                     self.neighbours.add(other)
+                    break
 
     def add_inline(self, other: "Axis") -> None:
         """Adds an axis that comes before/after this one"""
@@ -35,6 +36,7 @@ def add_inline(self, other: "Axis") -> None:
             for this_wire in self.wires:
                 for nei_wire in other.wires:
                     this_wire.add_inline(nei_wire)
+                    break
 
     def is_aligned(self, other: "Axis") -> bool:
         """Returns True if wires of the other axis are aligned

src/classy_blocks/util/functions.py

@@ -229,7 +229,8 @@ def divide_arc(center: NPPointType, point_1: NPPointType, point_2: NPPointType,
     return result
 
 
-def arc_mid(center: PointType, point_1: PointType, point_2: PointType) -> PointType:
+@jit(nopython=True, cache=True)
+def arc_mid(center: NPPointType, point_1: NPPointType, point_2: NPPointType) -> PointType:
     """Returns the midpoint of the specified arc in 3D space"""
     return divide_arc(center, point_1, point_2, 1)[0]