vibe-sketcherDocumentation
Source code for vibesketcher.document.spacegroups
"""Space groups from Hall symbols, in pure Python (#199).
The operations of every setting are generated from its Hall symbol (S. R.
Hall, Acta Cryst. A37, 517 (1981); the explicit-origin conventions of
R. W. Grosse-Kunstleve, Acta Cryst. A55, 383 (1999)): the lattice symbol adds
its centring translations (and inversion when it starts with ``-``), each
matrix symbol adds a generator, an optional ``(vx vy vz)`` shifts the origin
in twelfths, and the group is closed under multiplication modulo lattice
translations. The table of settings is ``_space_group_table`` (from Gemmi's
table, which follows the International Tables, Vol. A);
``scripts/generate_space_group_table.py`` checks every setting against Gemmi.
Operations act on fractional column vectors: x' = W x + t, with t reduced to
[0, 1). This module has no third-party dependency, so it can be checked by an
interpreter that only has Gemmi.
"""
from __future__ import annotations
import re
from dataclasses import dataclass
from fractions import Fraction
Matrix = tuple[tuple[int, int, int], tuple[int, int, int], tuple[int, int, int]]
Vector = tuple[Fraction, Fraction, Fraction]
_I: Matrix = ((1, 0, 0), (0, 1, 0), (0, 0, 1))
_H = Fraction(1, 2)
_T1, _T2 = Fraction(1, 3), Fraction(2, 3)
_Q = Fraction(1, 4)
_ZERO: Vector = (Fraction(0), Fraction(0), Fraction(0))
_CENTRING: dict[str, tuple[Vector, ...]] = {
"P": (),
"A": ((Fraction(0), _H, _H),),
"B": ((_H, Fraction(0), _H),),
"C": ((_H, _H, Fraction(0)),),
"I": ((_H, _H, _H),),
"R": ((_T2, _T1, _T1), (_T1, _T2, _T2)),
"S": ((_T1, _T1, _T2), (_T2, _T2, _T1)),
"T": ((_T1, _T2, _T1), (_T2, _T1, _T2)),
"F": ((Fraction(0), _H, _H), (_H, Fraction(0), _H), (_H, _H, Fraction(0))),
}
# Proper rotations about the principal axes (Hall 1981, table 3).
_ROTATIONS: dict[tuple[str, int], Matrix] = {
("z", 2): ((-1, 0, 0), (0, -1, 0), (0, 0, 1)),
("z", 3): ((0, -1, 0), (1, -1, 0), (0, 0, 1)),
("z", 4): ((0, -1, 0), (1, 0, 0), (0, 0, 1)),
("z", 6): ((1, -1, 0), (1, 0, 0), (0, 0, 1)),
("x", 2): ((1, 0, 0), (0, -1, 0), (0, 0, -1)),
("x", 3): ((1, 0, 0), (0, 0, -1), (0, 1, -1)),
("x", 4): ((1, 0, 0), (0, 0, -1), (0, 1, 0)),
("x", 6): ((1, 0, 0), (0, 1, -1), (0, 1, 0)),
("y", 2): ((-1, 0, 0), (0, 1, 0), (0, 0, -1)),
("y", 3): ((-1, 0, 1), (0, 1, 0), (-1, 0, 0)),
("y", 4): ((0, 0, 1), (0, 1, 0), (-1, 0, 0)),
("y", 6): ((0, 0, 1), (0, 1, 0), (-1, 0, 1)),
}
# Two-fold axes along face diagonals, relative to the preceding axis.
_PRIMED: dict[tuple[str, str], Matrix] = {
("z", "'"): ((0, -1, 0), (-1, 0, 0), (0, 0, -1)),
("z", '"'): ((0, 1, 0), (1, 0, 0), (0, 0, -1)),
("x", "'"): ((-1, 0, 0), (0, 0, -1), (0, -1, 0)),
("x", '"'): ((-1, 0, 0), (0, 0, 1), (0, 1, 0)),
("y", "'"): ((0, 0, -1), (0, -1, 0), (-1, 0, 0)),
("y", '"'): ((0, 0, 1), (0, -1, 0), (1, 0, 0)),
}
_BODY_DIAGONAL: Matrix = ((0, 0, 1), (1, 0, 0), (0, 1, 0))
_TRANSLATIONS: dict[str, Vector] = {
"a": (_H, Fraction(0), Fraction(0)),
"b": (Fraction(0), _H, Fraction(0)),
"c": (Fraction(0), Fraction(0), _H),
"n": (_H, _H, _H),
"u": (_Q, Fraction(0), Fraction(0)),
"v": (Fraction(0), _Q, Fraction(0)),
"w": (Fraction(0), Fraction(0), _Q),
"d": (_Q, _Q, _Q),
}
_AXIS_VECTOR = {"x": (1, 0, 0), "y": (0, 1, 0), "z": (0, 0, 1)}
[docs]
@dataclass(frozen=True)
class Operation:
"""x' = rotation @ x + translation (fractional coordinates, translation in [0, 1))."""
rotation: Matrix
translation: Vector
[docs]
def apply(self, point: tuple[float, float, float]) -> tuple[float, float, float]:
r, t = self.rotation, self.translation
return (
r[0][0] * point[0] + r[0][1] * point[1] + r[0][2] * point[2] + float(t[0]),
r[1][0] * point[0] + r[1][1] * point[1] + r[1][2] * point[2] + float(t[1]),
r[2][0] * point[0] + r[2][1] * point[1] + r[2][2] * point[2] + float(t[2]),
)
[docs]
def triplet(self) -> str:
"""The coordinate triplet, e.g. ``-y,x-y,z+1/3``."""
parts = []
for row, shift in zip(self.rotation, self.translation, strict=True):
text = ""
for coefficient, name in zip(row, "xyz", strict=True):
if coefficient:
sign = "-" if coefficient < 0 else ("+" if text else "")
text += sign + ("" if abs(coefficient) == 1 else str(abs(coefficient))) + name
if shift:
text += f"+{shift}" if text else str(shift)
parts.append(text or "0")
return ",".join(parts)
def _reduce(t: tuple[Fraction, ...]) -> Vector:
return (t[0] % 1, t[1] % 1, t[2] % 1)
def _mul(a: Matrix, b: Matrix) -> Matrix:
return tuple( # type: ignore[return-value]
tuple(sum(a[i][k] * b[k][j] for k in range(3)) for j in range(3)) for i in range(3)
)
def _mat_vec(a: Matrix, v: Vector) -> Vector:
return tuple(sum((a[i][k] * v[k] for k in range(3)), Fraction(0)) for i in range(3)) # type: ignore[return-value]
[docs]
def compose(first: Operation, second: Operation) -> Operation:
"""``second`` after ``first``: x -> W2 (W1 x + t1) + t2."""
rotation = _mul(second.rotation, first.rotation)
shifted = _mat_vec(second.rotation, first.translation)
return Operation(
rotation, _reduce(tuple(s + t for s, t in zip(shifted, second.translation, strict=True)))
)
_TOKEN = re.compile(r"^(-?)([12346])([xyz'\"*]?)([abcnuvwd12345]*)$")
def _generator(
token: str, index: int, previous: tuple[str, int] | None
) -> tuple[Operation, tuple[str, int]]:
match = _TOKEN.match(token)
if match is None:
raise ValueError(f"cannot read the Hall matrix symbol {token!r}")
improper, order_text, axis, extras = match.groups()
order = int(order_text)
if not axis:
if order == 1 or index == 0:
axis = "z"
elif index == 1 and order == 2:
if previous is None:
raise ValueError("a second two-fold axis needs a preceding one")
axis = "x" if previous[1] in (2, 4) else "'"
elif index == 2 and order == 3:
axis = "*"
else:
raise ValueError(f"the Hall symbol leaves the axis of {token!r} undetermined")
reference = previous[0] if previous and previous[0] in "xyz" else "z"
if order == 1:
rotation = _I
elif axis in "xyz":
rotation = _ROTATIONS[(axis, order)]
elif axis in "'\"":
if order != 2:
raise ValueError(f"{token!r}: only two-fold axes are primed")
rotation = _PRIMED[(reference, axis)]
else:
if order != 3:
raise ValueError(f"{token!r}: only three-fold axes lie along the body diagonal")
rotation = _BODY_DIAGONAL
translation = [Fraction(0), Fraction(0), Fraction(0)]
for symbol in extras:
if symbol.isdigit():
if axis not in "xyz":
raise ValueError(f"{token!r}: screw components need a principal axis")
step = Fraction(int(symbol), order)
for i, component in enumerate(_AXIS_VECTOR[axis]):
translation[i] += step * component
else:
for i in range(3):
translation[i] += _TRANSLATIONS[symbol][i]
if improper:
rotation = tuple(tuple(-x for x in row) for row in rotation) # type: ignore[assignment]
used_axis = axis if axis in "xyz" else reference
return Operation(rotation, _reduce(tuple(translation))), (used_axis, order)
[docs]
def operations_from_hall(hall: str) -> list[Operation]:
"""All operations of the setting the Hall symbol names (centring included)."""
text = hall.strip()
shift: Vector = _ZERO
origin = re.search(r"\(([^)]*)\)\s*$", text)
if origin:
values = origin.group(1).replace(",", " ").split()
if len(values) != 3:
raise ValueError(f"cannot read the origin shift in {hall!r}")
shift = tuple(Fraction(int(v), 12) for v in values) # type: ignore[assignment]
text = text[: origin.start()].strip()
tokens = text.split()
if not tokens:
raise ValueError("an empty Hall symbol")
lattice = tokens[0]
centro = lattice.startswith("-")
letter = lattice.lstrip("-").upper()
if letter not in _CENTRING:
raise ValueError(f"unknown lattice symbol in {hall!r}")
generators = [Operation(_I, v) for v in _CENTRING[letter]]
if centro:
generators.append(Operation(((-1, 0, 0), (0, -1, 0), (0, 0, -1)), _ZERO))
previous: tuple[str, int] | None = None
for index, token in enumerate(tokens[1:]):
op, previous = _generator(token, index, previous)
generators.append(op)
group = {Operation(_I, _ZERO)}
frontier = list(group)
while frontier:
fresh = []
for a in frontier:
for g in generators:
for product in (compose(a, g), compose(g, a)):
if product not in group:
group.add(product)
fresh.append(product)
if len(group) > 192:
raise ValueError(f"the Hall symbol {hall!r} does not close to a space group")
frontier = fresh
if any(shift):
# Change of origin by v: W' = W, t' = t + (I - W) v (Grosse-Kunstleve 1999).
moved = set()
for op in group:
wv = _mat_vec(op.rotation, shift)
moved_t = tuple(t + v - w for t, v, w in zip(op.translation, shift, wv, strict=True))
moved.add(Operation(op.rotation, _reduce(moved_t)))
group = moved
return sorted(group, key=lambda op: (op.rotation != _I, op.key()))
[docs]
@dataclass(frozen=True)
class Setting:
number: int
hm: str
extension: str
qualifier: str
hall: str
standard: bool
[docs]
@property
def symbol(self) -> str:
return f"{self.hm}:{self.extension}" if self.extension else self.hm
[docs]
def crystal_system(number: int) -> str:
if not 1 <= number <= 230:
raise ValueError("space-group numbers run from 1 to 230")
for upper, name in ((2, "triclinic"), (15, "monoclinic"), (74, "orthorhombic"),
(142, "tetragonal"), (167, "trigonal"), (194, "hexagonal")): # fmt: skip
if number <= upper:
return name
return "cubic"
[docs]
def settings() -> list[Setting]:
from vibesketcher.document._space_group_table import SETTINGS
return [Setting(*row) for row in SETTINGS]
def _compact(text: str) -> str:
return text.replace(" ", "").replace("_", "").lower()
[docs]
def find(symbol: str | int) -> Setting:
"""A setting by number, Hermann-Mauguin symbol or Hall symbol.
``"R -3 c"``, ``"R-3c"``, ``"R -3 c:H"``, ``"167"``, ``167`` and ``"-R 3 2\\"c"``
all name the same setting; ``":R"`` selects rhombohedral axes, ``":1"``/``":2"``
an origin choice.
"""
table = settings()
if isinstance(symbol, int) or str(symbol).strip().isdigit():
number = int(symbol)
numbered = [s for s in table if s.number == number]
if not numbered:
raise ValueError(f"no space group number {number}")
return _preferred(numbered)
text = str(symbol).strip()
name, _, extension = text.partition(":")
wanted = _compact(name)
candidates = [
s for s in table
if _compact(s.hm) == wanted or _compact(s.hall) == _compact(text)
or _compact(s.hm.replace(" 1", "")) == wanted
] # fmt: skip
if extension:
candidates = [s for s in candidates if s.extension.lower() == extension.strip().lower()]
if not candidates:
raise ValueError(f"unknown space group {symbol!r}")
return _preferred(candidates)
def _preferred(candidates: list[Setting]) -> Setting:
"""Without an extension: hexagonal axes, origin choice 2, then the reference setting.
Origin choice 2 (inversion centre at the origin) is the Bilbao Crystallographic
Server's default; give ``:1`` explicitly for origin choice 1.
"""
return sorted(candidates, key=lambda s: (s.extension not in ("", "H", "2"), not s.standard))[0]