Skip to content

Geometry & Molecular Operations

mol_ops

Referenced by:

SEA dataclass 

SEA(label: str, subset: array, axis: array)

Symmetry equivalent atoms (SEA).

SEAs are atoms that can be swapped with no distinguishable change in the molecule.

Parameters:

  • label (str or None) –

    Rotor type of the SEA set (e.g. Single Atom, Linear, Spherical, Regular Polygon, Oblate Symmetric Top).

  • subset (array) –

    Atom indices in the molecule that belong to this SEA set, shape (N,).

  • axis (array or None) –

    Candidate rotational symmetry axis, shape (3,).

Returned by:

transform 

transform(positions: array, M: array) -> array

Transform coordinates of molecule by matrix M and return new positions.

Parameters:

  • positions (array) –

    Molecule positions.

  • M (array) –

    Transformation matrix (e.g. rotation, reflection, etc.), shape (3,3)

Returns:

  • Atoms –

    Molecule with transformed atom coordinates

Source code in minimalsym/core/mol_ops.py 
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
@njit
def transform(positions: "np.array", M: "np.array") -> "np.array":
    """
    Transform coordinates of molecule by matrix M and return new positions.

    Parameters
    ----------
    positions: np.array
        Molecule positions.
    M: np.array
        Transformation matrix (e.g. rotation, reflection, etc.), shape (3,3)

    Returns
    -------
    ase.Atoms
        Molecule with transformed atom coordinates
    """
    return np.dot(positions, np.transpose(M))

distance_matrix 

distance_matrix(positions)

Calculate the interatomic distance matrix as all pairwise distances between atoms.

Parameters:

  • positions –

    Numpy array of positions of Molecule.

Returns:

  • array –

    Interatomic distance matrix, shape(len(mol),len(mol))

Source code in minimalsym/core/mol_ops.py 
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
@njit
def distance_matrix(positions):
    """
    Calculate the interatomic distance matrix as all pairwise distances between atoms.

    Parameters
    ----------
    positions: np.array
        Numpy array of positions of Molecule.

    Returns
    -------
    np.array
        Interatomic distance matrix, shape(len(mol),len(mol))
    """
    dm = np.zeros((len(positions), len(positions)))
    for i in range(len(positions)):
        for j in range(i, len(positions)):
            dx = positions[i,0] - positions[j,0]
            dy = positions[i,1] - positions[j,1]
            dz = positions[i,2] - positions[j,2]
            dm[i,j] = np.sqrt(dx*dx + dy*dy + dz*dz)
            dm[j,i] = dm[i,j]
    return dm

find_SEAs 

find_SEAs(mol)

Find sets of symmetry equivalent atoms. Permutations of the distance matrix reveal which atoms form symmetry equivalent sets.

Parameters:

  • mol –

    Molecule object.

Returns:

  • List[SEA] –

    List of symmetry equivalent atom sets

Source code in minimalsym/core/mol_ops.py 
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
def find_SEAs(mol):
    """
    Find sets of symmetry equivalent atoms.
    Permutations of the distance matrix reveal which atoms form symmetry equivalent sets.

    Parameters
    ----------
    mol: ase.Atoms
        Molecule object.

    Returns
    -------
    List[SEA]
        List of symmetry equivalent atom sets
    """
    sea_id = _jit_find_SEAs(len(mol), mol.positions, mol.info["geom_tol"])
    n_seas = int(sea_id.max()) + 1
    SEAs = []
    for k in range(n_seas):
        subset = np.where(sea_id == k)[0].astype(np.int64)
        SEAs.append(SEA("", subset, np.zeros(3)))
    return SEAs

get_SEAs_from_atom_map 

get_SEAs_from_atom_map(atom_map)

Find sets of symmetry equivalent atoms based on a Symtext.atom_map. Form symmetry equivalent sets from who each atom maps to.

Parameters:

  • atom_map –

    Symtext.atom_map

Returns:

  • List[SEA] –

    List of symmetry equivalent atom sets

Source code in minimalsym/core/mol_ops.py 
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
def get_SEAs_from_atom_map(atom_map):
    """
    Find sets of symmetry equivalent atoms based on a Symtext.atom_map.
    Form symmetry equivalent sets from who each atom maps to.

    Parameters
    ----------
    atom_map: np.array
        Symtext.atom_map

    Returns
    -------
    List[SEA]
        List of symmetry equivalent atom sets
    """
    SEAs = []
    subset_list = _jit_get_SEAs_from_atom_map(atom_map)

    for subset in subset_list:
        SEAs.append(SEA("", subset, np.zeros(3)))

    return SEAs

calcmoit 

calcmoit(atoms)

Calculates the moment of inertia tensor for a list of atoms.

Parameters:

  • atoms –

    Set of atoms.

Returns:

  • array –

    Cartesian moment of inertia tensor, shape(3,3).

Source code in minimalsym/core/mol_ops.py 
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
def calcmoit(atoms):
    """
    Calculates the moment of inertia tensor for a list of atoms.

    Parameters
    ----------
    atoms: ase.Atoms
        Set of atoms.

    Returns
    -------
    np.array
        Cartesian moment of inertia tensor, shape(3,3).
    """
    atoms.translate(-atoms.get_center_of_mass())
    masses = atoms.get_masses()
    positions = atoms.positions

    return _jit_calcmoit(positions, masses)

mol_orient

mol_orient.py — Molecular orientation utilities.

Provides rotate_mol_to_symels(), which aligns paxis with z and saxis with x so that the molecule sits in the canonical frame expected by pg_to_symels().

Referenced by:

rotate_mol_to_symels 

rotate_mol_to_symels(mol, paxis, saxis)

Rotate molecule so that paxis aligns with z and saxis aligns with x.

Returns the rotated molecule and the forward/inverse rotation matrices so that computed properties can be rotated back to the original orientation.

Parameters:

  • mol (Atoms) –
  • paxis ((array, shape(3))) –
  • saxis ((array, shape(3))) –

Returns:

  • tuple(Atoms, ndarray, ndarray) –

    Rotated molecule, rotation matrix, inverse rotation matrix (shape 3×3).

Source code in minimalsym/core/mol_orient.py 
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
def rotate_mol_to_symels(mol, paxis, saxis):
    """
    Rotate molecule so that paxis aligns with z and saxis aligns with x.

    Returns the rotated molecule and the forward/inverse rotation matrices so
    that computed properties can be rotated back to the original orientation.

    Parameters
    ----------
    mol : ase.Atoms
    paxis : np.array, shape (3,)
    saxis : np.array, shape (3,)

    Returns
    -------
    tuple(ase.Atoms, np.ndarray, np.ndarray)
        Rotated molecule, rotation matrix, inverse rotation matrix (shape 3×3).
    """
    new_mol = mol.copy()
    paxis = np.asarray(paxis, dtype=np.float64)
    saxis = np.asarray(saxis, dtype=np.float64)
    positions, rmat, rmat_inv = _jit_rotate_mol_to_symels(mol.positions, paxis, saxis)
    new_mol.positions = positions
    return new_mol, rmat, rmat_inv

atom_mapping

atom_mapping.py — Atom permutation map builders.

Provides _get_atom_mapping() and _get_linear_atom_mapping(), which build the (n_atoms × n_symels) integer array that records where each atom goes under each symmetry operation.

Referenced by:

_jit_where_you_go 

_jit_where_you_go(positions, geom_tol, atom, rrep)

Return the index of the atom that atom maps to under rrep.

Source code in minimalsym/core/atom_mapping.py 
16
17
18
19
20
21
22
23
24
25
@njit
def _jit_where_you_go(positions, geom_tol, atom, rrep):
    """Return the index of the atom that *atom* maps to under *rrep*."""
    ratom = np.dot(rrep, positions[atom, :].T)
    tol2 = geom_tol*geom_tol
    for i in range(len(positions)):
        dist = positions[i,:] - ratom
        if (dist[0]*dist[0] + dist[1]*dist[1] + dist[2]*dist[2]) < tol2:
            return i
    return np.int64(-1)

_jit_get_atom_mapping 

_jit_get_atom_mapping(positions, geom_tol, rreps)

Build the full (n_atoms × n_symels) permutation map.

Source code in minimalsym/core/atom_mapping.py 
27
28
29
30
31
32
33
34
35
36
37
38
39
@njit
def _jit_get_atom_mapping(positions, geom_tol, rreps):
    """Build the full (n_atoms × n_symels) permutation map."""
    natoms = len(positions)
    nsymels = rreps.shape[0]
    amap = np.empty((natoms, nsymels), dtype=np.int64)
    for i in range(natoms):
        for j in range(nsymels):
            amap[i, j] = np.int64(-1)
    for atom in range(natoms):
        for s in range(nsymels):
            amap[atom, s] = _jit_where_you_go(positions, geom_tol, atom, rreps[s])
    return amap

_where_you_go 

_where_you_go(mol, atom, symel)

Find the atom index that atom maps to under symel.

Parameters:

  • mol (Atoms) –
  • atom (int) –
  • symel (Symel) –

Returns:

  • int or None –
Source code in minimalsym/core/atom_mapping.py 
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
def _where_you_go(mol, atom, symel):
    """
    Find the atom index that *atom* maps to under *symel*.

    Parameters
    ----------
    mol : ase.Atoms
    atom : int
    symel : Symel

    Returns
    -------
    int or None
    """
    w = _jit_where_you_go(mol.positions, mol.info["geom_tol"], atom, symel.rrep)
    return w if w != -1 else None

_get_atom_mapping 

_get_atom_mapping(mol, symels)

Build the (n_atoms × n_symels) atom permutation map for non-linear groups.

Parameters:

  • mol (Atoms) –
  • symels (List[Symel]) –

Returns:

  • ndarray –

    shape (n_atoms, n_symels)

Raises:

  • Exception –

    If any atom fails to map under a symmetry operation.

Source code in minimalsym/core/atom_mapping.py 
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
def _get_atom_mapping(mol, symels):
    """
    Build the (n_atoms × n_symels) atom permutation map for non-linear groups.

    Parameters
    ----------
    mol : ase.Atoms
    symels : List[Symel]

    Returns
    -------
    : np.ndarray
        shape (n_atoms, n_symels)

    Raises
    ------
    Exception
        If any atom fails to map under a symmetry operation.
    """
    rreps = np.array([s.rrep for s in symels])
    amap = _jit_get_atom_mapping(mol.positions, mol.info["geom_tol"]*1.1, rreps)
    failed = np.argwhere(amap == -1)
    if len(failed) > 0:
        atom, s = int(failed[0, 0]), int(failed[0, 1])
        raise Exception(
            f"Atom {atom} not mapped to another atom "
            f"under symel {symels[s]}\nPositions: {mol.positions}\nRreps: {rreps[s]}"
        )
    return amap

_get_linear_atom_mapping 

_get_linear_atom_mapping(mol, pg)

Build the permutation map for linear point groups (C0v / D0h).

For C0v: identity only (each atom maps to itself). For D0h: identity + inversion.

Parameters:

  • mol (Atoms) –
  • pg (PointGroup) –

Returns:

  • ndarray –
Source code in minimalsym/core/atom_mapping.py 
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
def _get_linear_atom_mapping(mol, pg):
    """
    Build the permutation map for linear point groups (C0v / D0h).

    For C0v: identity only (each atom maps to itself).
    For D0h: identity + inversion.

    Parameters
    ----------
    mol : ase.Atoms
    pg : PointGroup

    Returns
    -------
    np.ndarray
    """
    natoms = len(mol)
    amap = np.arange(natoms, dtype=int).reshape((natoms, 1))
    if pg.family == "D":
        ungerade_map = np.zeros(natoms, dtype=int)
        inversion_symel = Symel("i", None, -np.eye(3), None, None, None)
        for atom in range(natoms):
            w = _where_you_go(mol, atom, inversion_symel)
            if w is not None:
                ungerade_map[atom] = w
            else:
                raise Exception(f"Atom {atom} not mapped to another atom under symel i")
        return np.column_stack((amap, ungerade_map))
    return amap