Coverage for /builds/hweiske/ase/ase/spacegroup/symmetrize.py: 98.78%
82 statements
« prev ^ index » next coverage.py v7.2.7, created at 2024-04-22 11:22 +0000
1"""
2Provides utility functions for FixSymmetry class
3"""
4import numpy as np
6from ase.utils import atoms_to_spglib_cell
8__all__ = ['refine_symmetry', 'check_symmetry']
11def print_symmetry(symprec, dataset):
12 print("ase.spacegroup.symmetrize: prec", symprec,
13 "got symmetry group number", dataset["number"],
14 ", international (Hermann-Mauguin)", dataset["international"],
15 ", Hall ", dataset["hall"])
18def refine_symmetry(atoms, symprec=0.01, verbose=False):
19 """
20 Refine symmetry of an Atoms object
22 Parameters
23 ----------
24 atoms - input Atoms object
25 symprec - symmetry precicion
26 verbose - if True, print out symmetry information before and after
28 Returns
29 -------
31 spglib dataset
33 """
34 import spglib
36 # test orig config with desired tol
37 dataset = check_symmetry(atoms, symprec, verbose=verbose)
39 # set actual cell to symmetrized cell vectors by copying
40 # transformed and rotated standard cell
41 std_cell = dataset['std_lattice']
42 trans_std_cell = dataset['transformation_matrix'].T @ std_cell
43 rot_trans_std_cell = trans_std_cell @ dataset['std_rotation_matrix']
44 atoms.set_cell(rot_trans_std_cell, True)
46 # get new dataset and primitive cell
47 dataset = check_symmetry(atoms, symprec=symprec, verbose=verbose)
48 res = spglib.find_primitive(atoms_to_spglib_cell(atoms), symprec=symprec)
49 prim_cell, prim_scaled_pos, prim_types = res
51 # calculate offset between standard cell and actual cell
52 std_cell = dataset['std_lattice']
53 rot_std_cell = std_cell @ dataset['std_rotation_matrix']
54 rot_std_pos = dataset['std_positions'] @ rot_std_cell
55 pos = atoms.get_positions()
56 dp0 = (pos[list(dataset['mapping_to_primitive']).index(0)] - rot_std_pos[
57 list(dataset['std_mapping_to_primitive']).index(0)])
59 # create aligned set of standard cell positions to figure out mapping
60 rot_prim_cell = prim_cell @ dataset['std_rotation_matrix']
61 inv_rot_prim_cell = np.linalg.inv(rot_prim_cell)
62 aligned_std_pos = rot_std_pos + dp0
64 # find ideal positions from position of corresponding std cell atom +
65 # integer_vec . primitive cell vectors
66 # here we are assuming that primitive vectors returned by find_primitive
67 # are compatible with std_lattice returned by get_symmetry_dataset
68 mapping_to_primitive = list(dataset['mapping_to_primitive'])
69 std_mapping_to_primitive = list(dataset['std_mapping_to_primitive'])
70 pos = atoms.get_positions()
71 for i_at in range(len(atoms)):
72 std_i_at = std_mapping_to_primitive.index(mapping_to_primitive[i_at])
73 dp = aligned_std_pos[std_i_at] - pos[i_at]
74 dp_s = dp @ inv_rot_prim_cell
75 pos[i_at] = (aligned_std_pos[std_i_at] - np.round(dp_s) @ rot_prim_cell)
76 atoms.set_positions(pos)
78 # test final config with tight tol
79 return check_symmetry(atoms, symprec=1e-4, verbose=verbose)
82def check_symmetry(atoms, symprec=1.0e-6, verbose=False):
83 """
84 Check symmetry of `atoms` with precision `symprec` using `spglib`
86 Prints a summary and returns result of `spglib.get_symmetry_dataset()`
87 """
88 import spglib
89 dataset = spglib.get_symmetry_dataset(atoms_to_spglib_cell(atoms),
90 symprec=symprec)
91 if verbose:
92 print_symmetry(symprec, dataset)
93 return dataset
96def is_subgroup(sup_data, sub_data, tol=1e-10):
97 """
98 Test if spglib dataset `sub_data` is a subgroup of dataset `sup_data`
99 """
100 for rot1, trns1 in zip(sub_data['rotations'], sub_data['translations']):
101 for rot2, trns2 in zip(sup_data['rotations'], sup_data['translations']):
102 if np.all(rot1 == rot2) and np.linalg.norm(trns1 - trns2) < tol:
103 break
104 else:
105 return False
106 return True
109def prep_symmetry(atoms, symprec=1.0e-6, verbose=False):
110 """
111 Prepare `at` for symmetry-preserving minimisation at precision `symprec`
113 Returns a tuple `(rotations, translations, symm_map)`
114 """
115 import spglib
117 dataset = spglib.get_symmetry_dataset(atoms_to_spglib_cell(atoms),
118 symprec=symprec)
119 if verbose:
120 print_symmetry(symprec, dataset)
121 rotations = dataset['rotations'].copy()
122 translations = dataset['translations'].copy()
123 symm_map = []
124 scaled_pos = atoms.get_scaled_positions()
125 for (rot, trans) in zip(rotations, translations):
126 this_op_map = [-1] * len(atoms)
127 for i_at in range(len(atoms)):
128 new_p = rot @ scaled_pos[i_at, :] + trans
129 dp = scaled_pos - new_p
130 dp -= np.round(dp)
131 i_at_map = np.argmin(np.linalg.norm(dp, axis=1))
132 this_op_map[i_at] = i_at_map
133 symm_map.append(this_op_map)
134 return (rotations, translations, symm_map)
137def symmetrize_rank1(lattice, inv_lattice, forces, rot, trans, symm_map):
138 """
139 Return symmetrized forces
141 lattice vectors expected as row vectors (same as ASE get_cell() convention),
142 inv_lattice is its matrix inverse (reciprocal().T)
143 """
144 scaled_symmetrized_forces_T = np.zeros(forces.T.shape)
146 scaled_forces_T = np.dot(inv_lattice.T, forces.T)
147 for (r, t, this_op_map) in zip(rot, trans, symm_map):
148 transformed_forces_T = np.dot(r, scaled_forces_T)
149 scaled_symmetrized_forces_T[:, this_op_map] += transformed_forces_T
150 scaled_symmetrized_forces_T /= len(rot)
151 symmetrized_forces = (lattice.T @ scaled_symmetrized_forces_T).T
153 return symmetrized_forces
156def symmetrize_rank2(lattice, lattice_inv, stress_3_3, rot):
157 """
158 Return symmetrized stress
160 lattice vectors expected as row vectors (same as ASE get_cell() convention),
161 inv_lattice is its matrix inverse (reciprocal().T)
162 """
163 scaled_stress = np.dot(np.dot(lattice, stress_3_3), lattice.T)
165 symmetrized_scaled_stress = np.zeros((3, 3))
166 for r in rot:
167 symmetrized_scaled_stress += np.dot(np.dot(r.T, scaled_stress), r)
168 symmetrized_scaled_stress /= len(rot)
170 sym = np.dot(np.dot(lattice_inv, symmetrized_scaled_stress), lattice_inv.T)
171 return sym