Passed
Push — master ( bc96b3...eb28f3 )
by Daniel
06:34
created

amd.periodicset.PeriodicSet.density()   A

Complexity

Conditions 2

Size

Total Lines 7
Code Lines 7

Duplication

Lines 0
Ratio 0 %

Importance

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