1
|
|
|
"""Implements the :class:`PeriodicSet` class representing a periodic |
2
|
|
|
set, defined by a motif and unit cell. This models a crystal with a |
3
|
|
|
point at the center of each atom. |
4
|
|
|
|
5
|
|
|
This is the type yielded by the readers |
6
|
|
|
:class:`amd.CifReader <.io.CifReader>` and |
7
|
|
|
:class:`amd.CSDReader <.io.CSDReader>`. A :class:`PeriodicSet` can be |
8
|
|
|
passed as the first argument to :func:`amd.AMD() <.calculate.AMD>` or |
9
|
|
|
:func:`amd.PDD() <.calculate.PDD>` to calculate its invariants. |
10
|
|
|
""" |
11
|
|
|
|
12
|
|
|
from typing import Optional |
13
|
|
|
import numpy as np |
14
|
|
|
|
15
|
|
|
from .utils import ( |
16
|
|
|
cellpar_to_cell, |
17
|
|
|
cellpar_to_cell_2D, |
18
|
|
|
cell_to_cellpar, |
19
|
|
|
cell_to_cellpar_2D, |
20
|
|
|
random_cell, |
21
|
|
|
) |
22
|
|
|
|
23
|
|
|
from ase.data import atomic_masses |
24
|
|
|
|
25
|
|
|
|
26
|
|
|
class PeriodicSet: |
27
|
|
|
"""A periodic set representats a crystal by putting a point in the |
28
|
|
|
center of each atom. Mathematically it's defined by a basis (unit |
29
|
|
|
cell) and collection of points (motif) which repeats according to |
30
|
|
|
the basis. Periodic sets are defined for any dimension. |
31
|
|
|
|
32
|
|
|
:class:`PeriodicSet` s are returned by the readers in the |
33
|
|
|
:mod:`.io` module. They can be passed to |
34
|
|
|
:func:`amd.AMD() <.calculate.AMD>` or |
35
|
|
|
:func:`amd.PDD() <.calculate.PDD>` to calculate the invariants. |
36
|
|
|
|
37
|
|
|
Parameters |
38
|
|
|
---------- |
39
|
|
|
motif : :class:`numpy.ndarray` |
40
|
|
|
Cartesian (orthogonal) coordinates of the motif, shape (no |
41
|
|
|
points, dims). |
42
|
|
|
cell : :class:`numpy.ndarray` |
43
|
|
|
Cartesian (orthogonal) square array representing the unit cell, |
44
|
|
|
shape (dims, dims). Use |
45
|
|
|
:func:`amd.cellpar_to_cell <.utils.cellpar_to_cell>` to convert |
46
|
|
|
6 cell parameters to an orthogonal square matrix. |
47
|
|
|
name : str, optional |
48
|
|
|
Name of the periodic set. |
49
|
|
|
asymmetric_unit : :class:`numpy.ndarray`, optional |
50
|
|
|
Indices for the asymmetric unit, pointing to the motif. Useful |
51
|
|
|
in invariant calculations. |
52
|
|
|
wyckoff_multiplicities : :class:`numpy.ndarray`, optional |
53
|
|
|
Wyckoff multiplicities of each atom in the asymmetric unit |
54
|
|
|
(number of unique sites generated under all symmetries). Useful |
55
|
|
|
in invariant calculations. |
56
|
|
|
types : :class:`numpy.ndarray`, optional |
57
|
|
|
Array of atomic numbers of motif points. |
58
|
|
|
""" |
59
|
|
|
|
60
|
|
|
def __init__( |
61
|
|
|
self, |
62
|
|
|
motif: np.ndarray, |
63
|
|
|
cell: np.ndarray, |
64
|
|
|
name: Optional[str] = None, |
65
|
|
|
asymmetric_unit: Optional[np.ndarray] = None, |
66
|
|
|
wyckoff_multiplicities: Optional[np.ndarray] = None, |
67
|
|
|
types: Optional[np.ndarray] = None |
68
|
|
|
): |
69
|
|
|
|
70
|
|
|
self.motif = motif |
71
|
|
|
self.cell = cell |
72
|
|
|
self.name = name |
73
|
|
|
self.asymmetric_unit = asymmetric_unit |
74
|
|
|
self.wyckoff_multiplicities = wyckoff_multiplicities |
75
|
|
|
self.types = types |
76
|
|
|
|
77
|
|
|
@property |
78
|
|
|
def ndim(self): |
79
|
|
|
return self.cell.shape[0] |
80
|
|
|
|
81
|
|
|
@property |
82
|
|
|
def density(self): |
83
|
|
|
if self.types is None: |
84
|
|
|
raise ValueError('') |
85
|
|
|
cell_mass = sum(atomic_masses[n] for n in self.types) |
86
|
|
|
cell_vol = np.linalg.det(self.cell) |
87
|
|
|
return cell_mass / cell_vol |
88
|
|
|
|
89
|
|
|
def __str__(self): |
90
|
|
|
"""Returns a string representation of the format |
91
|
|
|
PeriodicSet(name, motif (m, d), t asym sites, |
92
|
|
|
abcαβγ=cellpar). |
93
|
|
|
""" |
94
|
|
|
|
95
|
|
|
if self.asymmetric_unit is None: |
96
|
|
|
n_asym_sites = len(self.motif) |
97
|
|
|
else: |
98
|
|
|
n_asym_sites = len(self.asymmetric_unit) |
99
|
|
|
|
100
|
|
|
cellpar_str = None |
101
|
|
|
if self.ndim == 1: |
102
|
|
|
cellpar_str = f', cell={self.cell[0][0]}' |
103
|
|
|
if self.ndim == 2: |
104
|
|
|
cellpar = np.round(cell_to_cellpar_2D(self.cell), 2) |
105
|
|
|
cellpar_str = ','.join(str(round(p, 2)) for p in cellpar) |
106
|
|
|
cellpar_str = f', abα={cellpar_str}' |
107
|
|
|
elif self.ndim == 3: |
108
|
|
|
cellpar = np.round(cell_to_cellpar(self.cell), 2) |
109
|
|
|
cellpar_str = ','.join(str(round(p, 2)) for p in cellpar) |
110
|
|
|
cellpar_str = f', abcαβγ={cellpar_str}' |
111
|
|
|
else: |
112
|
|
|
cellpar_str = f'' |
113
|
|
|
|
114
|
|
|
s = f'PeriodicSet(name={self.name}, ' \ |
115
|
|
|
f'motif {self.motif.shape} ({n_asym_sites} asym sites)' \ |
116
|
|
|
f'{cellpar_str})' |
117
|
|
|
|
118
|
|
|
return s |
119
|
|
|
|
120
|
|
|
def __repr__(self): |
121
|
|
|
|
122
|
|
|
if self.asymmetric_unit is None: |
123
|
|
|
n_asym_sites = len(self.motif) |
124
|
|
|
else: |
125
|
|
|
n_asym_sites = len(self.asymmetric_unit) |
126
|
|
|
|
127
|
|
|
s = f'PeriodicSet(name={self.name}, ' \ |
128
|
|
|
f'motif={self.motif} ({n_asym_sites} asym sites), ' \ |
129
|
|
|
f'cell={self.cell})' |
130
|
|
|
|
131
|
|
|
return s |
132
|
|
|
|
133
|
|
|
|
134
|
|
|
def __eq__(self, other): |
135
|
|
|
"""Used for debugging/tests. Returns True if: |
136
|
|
|
1. The unit cells are identical (element for element within tol) |
137
|
|
|
2. The motifs are the same shape, and when transformed into |
138
|
|
|
fractional coordinates every point in one motif has an identical |
139
|
|
|
(within tol) point somewhere in the other, accounting for pbc. |
140
|
|
|
""" |
141
|
|
|
|
142
|
|
|
if self.cell.shape != other.cell.shape or \ |
143
|
|
|
self.motif.shape != other.motif.shape or \ |
144
|
|
|
not np.allclose(self.cell, other.cell): |
145
|
|
|
return False |
146
|
|
|
|
147
|
|
|
fm1 = np.mod(self.motif @ np.linalg.inv(self.cell), 1) |
148
|
|
|
fm2 = np.mod(other.motif @ np.linalg.inv(other.cell), 1) |
149
|
|
|
d1 = np.abs(fm2[:, None] - fm1) |
150
|
|
|
d2 = np.abs(d1 - 1) |
151
|
|
|
diffs = np.amax(np.minimum(d1, d2), axis=-1) |
152
|
|
|
|
153
|
|
|
if not np.all((np.amin(diffs, axis=0) <= 1e-8) | (np.amin(diffs, axis=-1) <= 1e-8)): |
154
|
|
|
return False |
155
|
|
|
|
156
|
|
|
return True |
157
|
|
|
|
158
|
|
|
def __ne__(self, other): |
159
|
|
|
return not self.__eq__(other) |
160
|
|
|
|
161
|
|
|
@staticmethod |
162
|
|
|
def cubic(scale=1, dims=3): |
163
|
|
|
"""Return a :class:`PeriodicSet` representing a cubic lattice. |
164
|
|
|
""" |
165
|
|
|
|
166
|
|
|
cell = np.identity(dims) * scale |
167
|
|
|
return PeriodicSet(np.zeros((1, dims)), cell) |
168
|
|
|
|
169
|
|
|
@staticmethod |
170
|
|
|
def hexagonal(scale=1, dims=3): |
171
|
|
|
"""Dimensions 2 and 3 only. Return a :class:`PeriodicSet` |
172
|
|
|
representing a hexagonal lattice. |
173
|
|
|
""" |
174
|
|
|
|
175
|
|
|
if dims == 3: |
176
|
|
|
cell = cellpar_to_cell(scale, scale, scale, 90, 90, 120) |
177
|
|
|
elif dims == 2: |
178
|
|
|
cell = cellpar_to_cell_2D(scale, scale, 60) |
179
|
|
|
else: |
180
|
|
|
msg = 'PeriodicSet.hexagonal() only implemented for dimensions ' \ |
181
|
|
|
f'2 and 3 (passed {dims})' |
182
|
|
|
raise NotImplementedError(msg) |
183
|
|
|
|
184
|
|
|
return PeriodicSet(np.zeros((1, dims)), cell) |
185
|
|
|
|
186
|
|
|
@staticmethod |
187
|
|
|
def _random( |
188
|
|
|
n_points, |
189
|
|
|
length_bounds=(1, 2), |
190
|
|
|
angle_bounds=(60, 120), |
191
|
|
|
dims=3 |
192
|
|
|
): |
193
|
|
|
"""Dimensions 2 and 3 only. Return a :class:`PeriodicSet` with a |
194
|
|
|
chosen number of randomly placed points, in random cell with |
195
|
|
|
edges between length_bounds and angles between angle_bounds. |
196
|
|
|
""" |
197
|
|
|
|
198
|
|
|
cell = random_cell( |
199
|
|
|
length_bounds=length_bounds, |
200
|
|
|
angle_bounds=angle_bounds, |
201
|
|
|
dims=dims |
202
|
|
|
) |
203
|
|
|
frac_motif = np.random.uniform(size=(n_points, dims)) |
204
|
|
|
return PeriodicSet(frac_motif @ cell, cell) |
205
|
|
|
|