Add a MappedSketch with automatic smoothing

examples/shape/monocylinder.py

@@ -0,0 +1,52 @@
+import os
+
+import numpy as np
+
+import classy_blocks as cb
+from classy_blocks.construct.flat.sketches.disk import QuarterDisk
+from classy_blocks.construct.flat.sketches.mapped import MappedSketch
+from classy_blocks.construct.shapes.shape import LoftedShape
+from classy_blocks.util import functions as f
+
+# An example case with a Shape, created from a custom sketch.
+# In this case, a cylinder with a single block in the middle and four blocks
+# surrounding it.
+mesh = cb.Mesh()
+
+center_point = np.array([0, 0, 0])
+radius_point = np.array([1, 0, 0])
+
+
+class MonoCylinderSketch(MappedSketch):
+    def add_edges(self):
+        for i in (1, 2, 3, 4):
+            self.faces[i].add_edge(1, cb.Origin(center_point))
+
+
+# a cylinder with a single block in the center
+angles = np.linspace(0, 2 * np.pi, num=4, endpoint=False)
+
+inner_points = [f.rotate(radius_point * QuarterDisk.diagonal_ratio, a, [0, 0, 1], center_point) for a in angles]
+outer_points = [f.rotate(radius_point, a, [0, 0, 1], center_point) for a in angles]
+
+locations = inner_points + outer_points
+quads = [
+    (0, 1, 2, 3),
+    (0, 4, 5, 1),
+    (1, 5, 6, 2),
+    (2, 6, 7, 3),
+    (3, 7, 4, 0),
+]
+
+base = MonoCylinderSketch(locations, quads, 5)
+
+cylinder = LoftedShape(base, base.copy().translate([0, 0, 1]))
+cylinder.operations[1].chop(2, count=10)
+
+cylinder.operations[1].chop(1, count=5)
+cylinder.operations[2].chop(1, count=5)
+
+cylinder.operations[1].chop(0, count=3)
+
+mesh.add(cylinder)
+mesh.write(os.path.join("..", "case", "system", "blockMeshDict"), debug_path="debug.vtk")

src/classy_blocks/construct/flat/sketches/mapped.py

@@ -0,0 +1,106 @@
+from typing import Dict, List, Set, Tuple
+
+import numpy as np
+
+from classy_blocks.construct.flat.face import Face
+from classy_blocks.construct.flat.sketches.sketch import Sketch
+from classy_blocks.types import NPPointListType, NPPointType, PointListType
+from classy_blocks.util import functions as f
+from classy_blocks.util.constants import DTYPE
+
+QuadType = Tuple[int, int, int, int]
+
+
+class Quad:
+    """A helper class for tracking positions-faces-indexes-neighbours-whatnot"""
+
+    def __init__(self, positions: NPPointListType, indexes: Tuple[int, int, int, int]):
+        self.indexes = indexes
+        self.face = Face(np.take(positions, list(indexes), axis=0))
+
+    @property
+    def connections(self) -> List[Set[int]]:
+        return [{self.indexes[i], self.indexes[(i + 1) % 4]} for i in range(4)]
+
+
+def get_connections(quad: QuadType) -> List[Set[int]]:
+    return [{quad[i], quad[(i + 1) % 4]} for i in range(4)]
+
+
+def get_all_connections(quads) -> List[Set[int]]:
+    return f.flatten_2d_list([get_connections(quad) for quad in quads])
+
+
+def find_neighbours(quads: List[QuadType]) -> Dict[int, Set[int]]:
+    """Returns a dictionary point:[neighbour points] as defined by quads"""
+    length = int(max(np.array(quads, dtype=DTYPE).ravel())) + 1
+
+    neighbours: Dict[int, Set[int]] = {i: set() for i in range(length)}
+    connections = get_all_connections(quads)
+
+    for connection in connections:
+        connection = list(connection)
+        neighbours[connection[0]].add(connection[1])
+        neighbours[connection[1]].add(connection[0])
+
+    return neighbours
+
+
+def get_fixed_points(quads) -> Set[int]:
+    """Returns indexes of points that can be smoothed"""
+    connections = get_all_connections(quads)
+    fixed_points: Set[int] = set()
+
+    for edge in connections:
+        if connections.count(edge) == 1:
+            fixed_points.update(edge)
+
+    return fixed_points
+
+
+def smooth(positions, quads, iterations: int) -> NPPointListType:
+    neighbours = find_neighbours(quads)
+    fixed_points = get_fixed_points(quads)
+
+    for _ in range(iterations):
+        for point_index, point_neighbours in neighbours.items():
+            if point_index in fixed_points:
+                continue
+
+            nei_positions = np.take(positions, list(point_neighbours), axis=0)
+            positions[point_index] = np.average(nei_positions, axis=0)
+
+    return positions
+
+
+class MappedSketch(Sketch):
+    """A sketch that is created from predefined points.
+    The points are connected to form quads which define Faces.
+    There is only a single grid 'tier' in these sketches (faces = grid[0])."""
+
+    def __init__(self, positions: PointListType, quads: List[Tuple[int, int, int, int]], smooth_iter: int = 0):
+        positions = np.array(positions, dtype=DTYPE)
+
+        if smooth_iter > 0:
+            positions = smooth(positions, quads, smooth_iter)
+
+        self.quads = [Quad(positions, quad) for quad in quads]
+
+        self.add_edges()
+
+    def add_edges(self) -> None:
+        """An optional method that will add edges to wherever they need to be."""
+
+    @property
+    def faces(self):
+        """A 'flattened' grid"""
+        return [quad.face for quad in self.quads]
+
+    @property
+    def grid(self):
+        return [self.faces]
+
+    @property
+    def center(self) -> NPPointType:
+        """Center of this sketch"""
+        return np.average([face.center for face in self.faces], axis=0)

src/classy_blocks/construct/shapes/shape.py

@@ -93,11 +93,6 @@ def set_end_patch(self, name: str) -> None:
     def operations(self):
         return f.flatten_2d_list(self.lofts)
 
-    @property
-    @abc.abstractmethod
-    def shell(self) -> List[Loft]:
-        """Operations on the outside of the shape"""
-
     @property
     def grid(self):
         """Analogous to Sketch's grid but corresponsing operations are returned"""

tests/test_construct/test_flat/test_sketch.py

@@ -0,0 +1,43 @@
+import unittest
+
+from classy_blocks.construct.flat.sketches.mapped import find_neighbours, get_fixed_points
+
+
+class MappedSketchTests(unittest.TestCase):
+    def test_find_neighbours(self):
+        # A random blocking (quadding)
+        quads = [(1, 2, 7, 6), (2, 3, 4, 7), (7, 4, 5, 6), (0, 1, 6, 8)]
+
+        neighbours = find_neighbours(quads)
+
+        self.assertDictEqual(
+            neighbours,
+            {
+                0: {8, 1},
+                1: {0, 2, 6},
+                2: {1, 3, 7},
+                3: {2, 4},
+                4: {3, 5, 7},
+                5: {4, 6},
+                6: {8, 1, 5, 7},
+                7: {2, 4, 6},
+                8: {0, 6},
+            },
+        )
+
+    def test_fixed_points(self):
+        # Monocylinder, core is quads[0]
+        quads = quads = [
+            (0, 1, 2, 3),
+            (0, 4, 5, 1),
+            (1, 5, 6, 2),
+            (2, 6, 7, 3),
+            (3, 7, 4, 0),
+        ]
+
+        fixed_points = get_fixed_points(quads)
+
+        self.assertSetEqual(
+            fixed_points,
+            {4, 5, 6, 7},
+        )