trihexagonal tiling

democapsid/src/democapsid.py

@@ -16,21 +16,37 @@
 SQRT5 = np.sqrt(5)
 PHI = (1 + SQRT5) / 2
 
-# ⬢ lattice unit
-HEXAGON = np.array([
+R6 = 1
+r6 = R6 * (SQRT3 / 2)
+
+HEXAGON_1 = np.array([
         [0, 1], 
-        [SQRT3 / 2, 0.5], 
-        [SQRT3 / 2, -0.5], 
+        [r6, 0.5], 
+        [r6, -0.5], 
         [0, -1], 
-        [-(SQRT3 / 2), -0.5], 
-        [-(SQRT3 / 2), 0.5]
+        [-r6, -0.5], 
+        [-r6, 0.5]
     ]
 )
-# tilers
-TILERS = [
-    lambda coor: HEXAGON + coor
+HEXAGON_2 = [
+    [R6 / 2, r6],
+    [R6, 0],
+    [R6 / 2, -r6],
+    [-R6 / 2, -r6],
+    [-R6, 0],
+    [-R6 / 2, r6]
 ]
-# icosahedron symmetry config
+LATTICE = {
+    "hex": (
+        lambda R: (R * (SQRT3 / 2)) * np.array([[2, 0], [1, SQRT3]]),
+        [lambda coor: HEXAGON + coor]
+    ),
+    "trihex": (
+        lambda R: np.array([[2 * R, 0], [R, R * SQRT3]]),
+        [lambda coor: HEXAGON_2 + coor]
+    )
+}
+
 ICO_CONFIG = (
     (),
     (),
@@ -55,7 +71,6 @@
     )
 )
 
-
 def iter_ring(elements):
     yield from ((elements[i], elements[(i + 1) % len(elements)]) for i in range(len(elements))) if len(elements) > 1 else ()
 
@@ -99,13 +114,6 @@ def in_triangle(p, q1, q2, q3):
     )
 
 
-def update_coor_to_id(coor_to_id, coors):
-    for coor in map(tuple, coors):
-        identifier = coor_to_id.get(coor, len(coor_to_id))
-        coor_to_id[coor] = identifier
-        yield identifier
-
-
 def intersection(p1, q1, p2, q2):
     # http://paulbourke.net/geometry/pointlineplane/edge_intersection.py
 
@@ -338,25 +346,21 @@ def main(argv):
     args = parse_args(argv[1:])
     h, k, H, K, s, c = args.h, args.k, args.H, args.K, args.symmetry, args.sphericity
 
-    # hexagonal tile unit inradius
-    r = SQRT3 / 2
-    # hexagonal tile lattice basis
-    basis = np.array([[2 * r, 0], [r, SQRT3 * r]])
+    # tile
+    tile = "trihex"
+
+    # lattice basis
+    basis = LATTICE[tile][0](1)
     v1, v2, v3 = (*basis, basis[1] @ rmat(np.pi / 3)) 
+
     # Caspar-Klug vectors
-    # self.q1 = np.array([0, 0]), self.C1, self.C4, self.C1 + self.C4
-    # self.q2 = np.array([0, 0]), self.C2, self.C1, self.C2 + self.C1
-    # self.q3 = np.array([0, 0]), self.C3, self.C2, self.C3 + self.C2
     ckv = [
         [h, k] @ basis,
         [H, K] @ np.stack([v2, v3]),
         [-h - k, h] @ basis,
         [k, -h] @ np.stack([v1, v3])
     ]
     # Caspar-Klug triangles
-    # self.C1 = h * self.a1 + k * self.a2         # \vec{C}_{T}
-    # self.C2 = H * self.A1 + K * self.A2         # \vec{C}_{Q}
-    # self.C3 = (-h - k) * self.a1 + h * self.a2  # \vec{C}^{120^{\circ}}_{T}
     ckt = [
         [],
         [np.array([0, 0]), ckv[0], ckv[3]],
@@ -379,8 +383,8 @@ def main(argv):
         mesh = []
         for coor in lattice_coordinates:
             # process tile subunits
-            for tiler in TILERS:
-                path = list(iter_ring(tiler(coor @ basis)))
+            for calc_tile in LATTICE[tile][1]:
+                path = list(iter_ring(calc_tile(coor @ basis)))
                 vertices = []
                 # iterate polygon edges
                 for src, tar in path:
@@ -391,15 +395,15 @@ def main(argv):
                         # add point that at the intersetion of the polygon and triangle edges
                         if (x := intersection(src, tar, *edge)).any():
                             vertices.append((np.append(x, 1), 1))
-                # keep edges if they occuron the tile polygon path
-                edges = [
-                    (s1, t1) 
-                    for s1, t1 in iter_ring(list(range(len(vertices))))
-                    if any(on_same_line(vertices[s1][0], vertices[t1][0], np.append(s2, 1), np.append(t2, 1)) for s2, t2 in path)
-                ]
-                # if there are only two edges and they point to each other, then only keep one
-                edges = [edges[0]] if len(edges) == 2 and edges[0] == edges[1][::-1] else edges        
-                edges and mesh.append(([ele[0] for ele in vertices], edges))
+            # keep edges if they occur on the tile polygon path
+            edges = [
+                (s1, t1) 
+                for s1, t1 in iter_ring(list(range(len(vertices))))
+                if any(on_same_line(vertices[s1][0], vertices[t1][0], np.append(s2, 1), np.append(t2, 1)) for s2, t2 in path)
+            ]
+            # if there are only two edges and they point to each other, then only keep one
+            edges = [edges[0]] if len(edges) == 2 and edges[0] == edges[1][::-1] else edges        
+            edges and mesh.append(([ele[0] for ele in vertices], edges))
         mesh and meshes.append(mesh)
 
     meshes3d = [[], [], [], []]