fully debugged the symmetry code

pyxtal/__init__.py

@@ -810,7 +810,7 @@ def _subgroup_by_splitter(self, splitter, eps=0.05, mut_lat=False):
                         rot = g1s[0].affine_matrix[:3,:3].T
                     else:
                         #for special wyc, needs to get better treatment
-                        op = wp1.get_euclidean_operation(self.lattice.matrix, id)
+                        op = wp1.get_euclidean_generator(self.lattice.matrix, id)
                         rot = op.affine_matrix[:3, :3].T 
 
                     # xyz in new lattice

pyxtal/symmetry.py

@@ -50,7 +50,7 @@
 t_subgroup = loadfn(rf("pyxtal",'database/t_subgroup.json'))
 k_subgroup = loadfn(rf("pyxtal",'database/k_subgroup.json'))
 wyc_sets = loadfn(rf("pyxtal",'database/wyckoff_sets.json'))
-
+hex_cell = np.array([[1, -0.5, 0], [0, np.sqrt(3) / 2, 0], [0, 0, 1]])
 # --------------------------- Group class -----------------------------
 class Group:
     """
@@ -742,6 +742,7 @@ def get_dof(self):
     def __str__(self, supress=False):
         if self.dim not in list(range(4)):
             return "invalid crystal dimension. Must be a number between 0 and 3."
+        if not hasattr(self, "site_symm"): self.get_site_symmetry()   
         s = "Wyckoff position " + str(self.multiplicity) + self.letter + " in "
         if self.dim == 3:
             s += "space "
@@ -753,10 +754,7 @@ def __str__(self, supress=False):
             s += "Point group " + self.symbol
         if self.dim != 0:
             s += "group " + str(self.number)
-        #s += " with site symmetry " + self.site_symm
-        s += " with site symmetry " + ss_string_from_ops(
-                self.symmetry_m[0], self.number, dim=self.dim
-        )
+        s += " with site symmetry " + self.site_symm
 
         if not supress:
             for op in self.ops:
@@ -1143,7 +1141,31 @@ def __len__(self):
         return self.multiplicity
 
     def get_site_symmetry(self):
-        self.site_symm = ss_string_from_ops(self.symmetry[0], self.number, dim=self.dim)
+        if self.euclidean:
+            ops = self.get_euclidean_symmetries()
+        else:
+            ops = self.symmetry[0]
+        self.site_symm = ss_string_from_ops(ops, self.number, dim=self.dim)
+
+
+    def get_euclidean_symmetries(self):
+        """
+        return the symmetry operation object at the Euclidean space 
+
+        Returns:
+            list of pymatgen SymmOp object
+        """
+        ops = []
+        for op in self.symmetry[0]:
+            hat = SymmOp.from_rotation_and_translation(hex_cell, [0, 0, 0])
+            ops.append(hat * op * hat.inverse)
+        return ops
+
+    def print_ops(self, ops=None):
+        if ops is None:
+            ops = self.ops
+        for op in ops:
+            print(op.as_xyz_string())
 
     def gen_pos(self):
         """
@@ -1262,7 +1284,7 @@ def merge(self, pt, lattice, tol, orientations=None):
             else:
                 return pt, wp, valid_ori
 
-    def get_euclidean_operation(self, cell, idx=0):
+    def get_euclidean_generator(self, cell, idx=0):
         """
         return the symmetry operation object at the Euclidean space 
 
@@ -1280,6 +1302,7 @@ def get_euclidean_operation(self, cell, idx=0):
 
         return op
 
+
     def set_euclidean(self):
         convert = False
         if self.dim == 3:

pyxtal/test_all.py

@@ -142,6 +142,17 @@ def test_molecular_nodiag(self):
 
 
 class TestWP(unittest.TestCase):
+    def test_wp_site_symm(self):
+        data = [(143, 1, '3 . .'), 
+                (230, 6, '. 3 2'), 
+                (160, 1, '. . m'),
+                (160, 2, '3 m .')]
+        for d in data:
+            (sg, i, symbol) = d
+            wp = Group(sg)[i]
+            wp.get_site_symmetry()
+            self.assertTrue(wp.site_symm == symbol)
+
     def test_wp_label(self):
         symbol = str(wp1.multiplicity) + wp1.letter
         self.assertTrue(symbol == "8b")
@@ -206,7 +217,7 @@ def check_error(spg, pt, cell):
                     op0 = wp[i]
                     p1 = op0.operate(pt)
 
-                    op1 = wp.get_euclidean_operation(cell.matrix, i)
+                    op1 = wp.get_euclidean_generator(cell.matrix, i)
                     if wp.euclidean:
                         p2 = np.dot(op1.operate(p0), cell.inv_matrix)
                     else:

pyxtal/wyckoff_site.py

@@ -428,7 +428,7 @@ def _get_coords_and_species(self, absolute=False, PBC=False, first=False, unitce
             # Rotate the molecule (Euclidean metric)
             #op2_m = self.wp.generators_m[point_index]
 
-            op2_m = self.wp.get_euclidean_operation(self.lattice.matrix, point_index)
+            op2_m = self.wp.get_euclidean_generator(self.lattice.matrix, point_index)
             rot = op2_m.affine_matrix[:3, :3].T
             #if self.diag and self.wp.index > 0:
             #    tau = op2.translation_vector
@@ -569,7 +569,7 @@ def get_mol_object(self, id=0):
             #op_m = self.wp.generators_m[id]
             #rot = op_m.affine_matrix[0:3][:, 0:3].T
             #tau = op_m.affine_matrix[0:3][:, 3]
-            op0 = self.wp.get_euclidean_operation(self.lattice.matrix, id)
+            op0 = self.wp.get_euclidean_generator(self.lattice.matrix, id)
             rot = op0.rotation_matrix.T
             tmp = np.dot(coord0, rot) 
             # Add absolute center to molecule