Passed
Push — master ( f4ca5c...e2877c )
by Daniel
07:55
created

amd.periodicset.PeriodicSet.ndim()   A

Complexity

Conditions 1

Size

Total Lines 3
Code Lines 3

Duplication

Lines 0
Ratio 0 %

Importance

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