Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Remove symmetry elements for SaveINS #38605

Open
wants to merge 16 commits into
base: release-next
Choose a base branch
from
73 changes: 71 additions & 2 deletions Framework/PythonInterface/plugins/algorithms/SaveINS.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,41 @@
# Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
# SPDX - License - Identifier: GPL - 3.0 +
from mantid.api import AlgorithmFactory, FileProperty, FileAction, WorkspaceProperty, PythonAlgorithm
from mantid.kernel import Direction
from mantid.kernel import Direction, V3D
from mantid.geometry import SymmetryOperationFactory, SpaceGroupFactory
from os import path, makedirs
import numpy as np


class SaveINS(PythonAlgorithm):
LATT_TYPE_MAP = {type: itype + 1 for itype, type in enumerate(["P", "I", "R", "F", "A", "B", "C"])}
IDENTIY_OP = SymmetryOperationFactory.createSymOp("x,y,z")
INVERSION_OP = SymmetryOperationFactory.createSymOp("-x,-y,-z")
ROTATION_OPS = {1: [IDENTIY_OP, INVERSION_OP], -1: [IDENTIY_OP]}
CENTERING_OPS = {
1: [IDENTIY_OP],
2: [
SymmetryOperationFactory.createSymOp("x+1/2,y+1/2,z+1/2"),
],
3: [
SymmetryOperationFactory.createSymOp("x+1/3,y+2/3,z+2/3"),
SymmetryOperationFactory.createSymOp("x+2/3,y+1/3,z+1/3"),
],
4: [
SymmetryOperationFactory.createSymOp("x,y+1/2,z+1/2"),
SymmetryOperationFactory.createSymOp("x+1/2,y,z+1/2"),
SymmetryOperationFactory.createSymOp("x+1/2,y+1/2,z"),
],
5: [
SymmetryOperationFactory.createSymOp("x,y+1/2,z+1/2"),
],
6: [
SymmetryOperationFactory.createSymOp("x+1/2,y,z+1/2"),
],
7: [
SymmetryOperationFactory.createSymOp("x+1/2,y+1/2,z"),
],
}
DUMMY_WAVELENGTH = 1.0

def category(self):
Expand Down Expand Up @@ -113,7 +140,7 @@ def PyExec(self):
f_handle.write(f"LATT {latt_type}\n")

# print sym operations
for sym_str in spgr.getSymmetryOperationStrings():
for sym_str in self._get_shelx_symmetry_operators(spgr, latt_type):
f_handle.write(f"SYMM {sym_str}\n")

# print atom info
Expand All @@ -140,5 +167,47 @@ def PyExec(self):
f_handle.write("HKLF 2\n") # tells SHELX the columns saved in the reflection file
f_handle.write("END")

def _symmetry_operation_key(self, W1, w1, W2=np.eye(3), w2=np.zeros(3)):
"""Symmetry element key (9 element tuple) for comparison"""
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Perhaps a bit more explanation as to why this works as a unique identifier?

W = W1 @ W2
w = W1 @ w2 + w1
w[w < 0] += 1
w[w > 1] -= 1
return tuple(np.round(W, 0).astype(int).flatten().tolist() + np.round(w, 3).tolist())

def _symmetry_matrix_vector(self, symop):
"""Symmetry rotation matrix and translation vector"""
W = np.linalg.inv(np.column_stack([symop.transformHKL(V3D(*vec)) for vec in np.eye(3)]))
w = np.array(symop.transformCoordinates(V3D(0, 0, 0)))
return W, w

def _get_shelx_symmetry_operators(self, spgr, latt_type):
"""Get SHELX SYMM records https://cci.lbl.gov/cctbx/shelx.html"""
indentity = self.IDENTIY_OP.getIdentifier()
inversion = self.INVERSION_OP.getIdentifier()
latt_numb = abs(latt_type)
latt_sign = 1 if latt_type > 0 else -1
sym_ops = spgr.getSymmetryOperations()
sym_ops_list = []
sym_ops_set = set()
for sym_op in sym_ops:
# space group symmetry operator
W1, w1 = self._symmetry_matrix_vector(sym_op)
sym_key = sym_op.getIdentifier()
if sym_key != indentity and sym_key != inversion:
S1 = self._symmetry_operation_key(W1, w1)
if S1 not in sym_ops_set:
sym_ops_list.append(sym_key)
# identity(/inversion) rotation symmetry operator
for rot in self.ROTATION_OPS[latt_sign]:
Copy link
Contributor

@RichardWaiteSTFC RichardWaiteSTFC Jan 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Seems like this is a repeat of earlier just with sym_op instead of rot (is that correct?)

for rot in self.ROTATION_OPS[latt_sign]:
W2, _ = self._symmetry_matrix_vector(rot)
# lattice centering translation symmetry operator
for cent in self.CENTERING_OPS[latt_numb]:
_, w2 = self._symmetry_matrix_vector(cent)
# equivalent symmetry operator generated from rotation and translation
S3 = self._symmetry_operation_key(np.eye(3), np.zeros(3), W2, w2)

Could it be refactored into a separate function?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I re-organized the code a bit so it is easier to follow.

W2, _ = self._symmetry_matrix_vector(rot)
# lattice centering translation symmetry operator
for cent in self.CENTERING_OPS[latt_numb]:
_, w2 = self._symmetry_matrix_vector(cent)
# equivalent symmetry operator generated from rotation and translation
S3 = self._symmetry_operation_key(W1, w1, W2, w2)
sym_ops_set.add(S3)
return sym_ops_list


AlgorithmFactory.subscribe(SaveINS)
Original file line number Diff line number Diff line change
Expand Up @@ -25,9 +25,6 @@ def setUpClass(cls):
"ZERR 4 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000\n",
"LATT 1\n",
"SYMM -x+1/2,y+1/2,-z+1/2\n",
"SYMM -x,-y,-z\n",
"SYMM x+1/2,-y+1/2,z+1/2\n",
"SYMM x,y,z\n",
"NEUT\n",
]
cls.file_end = ["UNIT 48 36 12 8 4\n", "MERG 0\n", "HKLF 2\n", "END"]
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
- Fix :ref:`SaveINS<algm-SaveINS-v1>` that saved all symmetry records to file. Only the minimum are needed that can be generated by translation/rotation corresponding to the lattice type.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this needs to move a level deeper, into the Bugfixes folder

Loading