Coverage for /builds/kinetik161/ase/ase/spacegroup/xtal.py: 93.42%
76 statements
« prev ^ index » next coverage.py v7.2.7, created at 2023-12-10 11:04 +0000
« prev ^ index » next coverage.py v7.2.7, created at 2023-12-10 11:04 +0000
1# Copyright (C) 2010, Jesper Friis
2# (see accompanying license files for details).
4"""
5A module for ASE for simple creation of crystalline structures from
6knowledge of the space group.
8"""
10from typing import Any, Dict
12import numpy as np
13from scipy import spatial
15import ase
16from ase.geometry import cellpar_to_cell
17from ase.spacegroup import Spacegroup
18from ase.symbols import string2symbols
20__all__ = ['crystal']
23def crystal(symbols=None, basis=None, occupancies=None, spacegroup=1, setting=1,
24 cell=None, cellpar=None,
25 ab_normal=(0, 0, 1), a_direction=None, size=(1, 1, 1),
26 onduplicates='warn', symprec=0.001,
27 pbc=True, primitive_cell=False, **kwargs) -> ase.Atoms:
28 """Create an Atoms instance for a conventional unit cell of a
29 space group.
31 Parameters:
33 symbols : str | sequence of str | sequence of Atom | Atoms
34 Element symbols of the unique sites. Can either be a string
35 formula or a sequence of element symbols. E.g. ('Na', 'Cl')
36 and 'NaCl' are equivalent. Can also be given as a sequence of
37 Atom objects or an Atoms object.
38 basis : list of scaled coordinates
39 Positions of the unique sites corresponding to symbols given
40 either as scaled positions or through an atoms instance. Not
41 needed if *symbols* is a sequence of Atom objects or an Atoms
42 object.
43 occupancies : list of site occupancies
44 Occupancies of the unique sites. Defaults to 1.0 and thus no mixed
45 occupancies are considered if not explicitly asked for. If occupancies
46 are given, the most dominant species will yield the atomic number.
47 The occupancies in the atoms.info['occupancy'] dictionary will have
48 integers keys converted to strings. The conversion is done in order
49 to avoid unexpected conversions when using the JSON serializer.
50 spacegroup : int | string | Spacegroup instance
51 Space group given either as its number in International Tables
52 or as its Hermann-Mauguin symbol.
53 setting : 1 | 2
54 Space group setting.
55 cell : 3x3 matrix
56 Unit cell vectors.
57 cellpar : [a, b, c, alpha, beta, gamma]
58 Cell parameters with angles in degree. Is not used when `cell`
59 is given.
60 ab_normal : vector
61 Is used to define the orientation of the unit cell relative
62 to the Cartesian system when `cell` is not given. It is the
63 normal vector of the plane spanned by a and b.
64 a_direction : vector
65 Defines the orientation of the unit cell a vector. a will be
66 parallel to the projection of `a_direction` onto the a-b plane.
67 size : 3 positive integers
68 How many times the conventional unit cell should be repeated
69 in each direction.
70 onduplicates : 'keep' | 'replace' | 'warn' | 'error'
71 Action if `basis` contain symmetry-equivalent positions:
72 'keep' - ignore additional symmetry-equivalent positions
73 'replace' - replace
74 'warn' - like 'keep', but issue an UserWarning
75 'error' - raises a SpacegroupValueError
76 symprec : float
77 Minimum "distance" betweed two sites in scaled coordinates
78 before they are counted as the same site.
79 pbc : one or three bools
80 Periodic boundary conditions flags. Examples: True,
81 False, 0, 1, (1, 1, 0), (True, False, False). Default
82 is True.
83 primitive_cell : bool
84 Whether to return the primitive instead of the conventional
85 unit cell.
87 Keyword arguments:
89 All additional keyword arguments are passed on to the Atoms
90 constructor. Currently, probably the most useful additional
91 keyword arguments are `info`, `constraint` and `calculator`.
93 Examples:
95 Two diamond unit cells (space group number 227)
97 >>> diamond = crystal('C', [(0,0,0)], spacegroup=227,
98 ... cellpar=[3.57, 3.57, 3.57, 90, 90, 90], size=(2,1,1))
99 >>> ase.view(diamond) # doctest: +SKIP
101 A CoSb3 skutterudite unit cell containing 32 atoms
103 >>> skutterudite = crystal(('Co', 'Sb'),
104 ... basis=[(0.25,0.25,0.25), (0.0, 0.335, 0.158)],
105 ... spacegroup=204, cellpar=[9.04, 9.04, 9.04, 90, 90, 90])
106 >>> len(skutterudite)
107 32
108 """
109 sg = Spacegroup(spacegroup, setting)
110 if (
111 not isinstance(symbols, str) and
112 hasattr(symbols, '__getitem__') and
113 len(symbols) > 0 and
114 isinstance(symbols[0], ase.Atom)):
115 symbols = ase.Atoms(symbols)
116 if isinstance(symbols, ase.Atoms):
117 basis = symbols
118 symbols = basis.get_chemical_symbols()
119 if isinstance(basis, ase.Atoms):
120 basis_coords = basis.get_scaled_positions()
121 if cell is None and cellpar is None:
122 cell = basis.cell
123 if symbols is None:
124 symbols = basis.get_chemical_symbols()
125 else:
126 basis_coords = np.array(basis, dtype=float, copy=False, ndmin=2)
128 if occupancies is not None:
129 occupancies_dict = {}
131 for index, coord in enumerate(basis_coords):
132 # Compute all distances and get indices of nearest atoms
133 dist = spatial.distance.cdist(coord.reshape(1, 3), basis_coords)
134 indices_dist = np.flatnonzero(dist < symprec)
136 occ = {symbols[index]: occupancies[index]}
138 # Check nearest and update occupancy
139 for index_dist in indices_dist:
140 if index == index_dist:
141 continue
142 else:
143 occ.update({symbols[index_dist]: occupancies[index_dist]})
145 occupancies_dict[str(index)] = occ.copy()
147 sites, kinds = sg.equivalent_sites(basis_coords,
148 onduplicates=onduplicates,
149 symprec=symprec)
151 # this is needed to handle deuterium masses
152 masses = None
153 if 'masses' in kwargs:
154 masses = kwargs['masses'][kinds]
155 del kwargs['masses']
157 symbols = parse_symbols(symbols)
159 if occupancies is None:
160 symbols = [symbols[i] for i in kinds]
161 else:
162 # make sure that we put the dominant species there
163 symbols = [sorted(occupancies_dict[str(i)].items(),
164 key=lambda x: x[1])[-1][0] for i in kinds]
166 if cell is None:
167 cell = cellpar_to_cell(cellpar, ab_normal, a_direction)
169 info: Dict[str, Any] = {}
170 info['spacegroup'] = sg
171 if primitive_cell:
172 info['unit_cell'] = 'primitive'
173 else:
174 info['unit_cell'] = 'conventional'
176 if 'info' in kwargs:
177 info.update(kwargs['info'])
179 if occupancies is not None:
180 info['occupancy'] = occupancies_dict
182 kwargs['info'] = info
184 atoms = ase.Atoms(symbols,
185 scaled_positions=sites,
186 cell=cell,
187 pbc=pbc,
188 masses=masses,
189 **kwargs)
191 if isinstance(basis, ase.Atoms):
192 for name in basis.arrays:
193 if not atoms.has(name):
194 array = basis.get_array(name)
195 atoms.new_array(name, [array[i] for i in kinds],
196 dtype=array.dtype, shape=array.shape[1:])
198 if kinds:
199 atoms.new_array('spacegroup_kinds', np.asarray(kinds, dtype=int))
201 if primitive_cell:
202 from ase.build import cut
203 prim_cell = sg.scaled_primitive_cell
205 # Preserve calculator if present:
206 calc = atoms.calc
207 atoms = cut(atoms, a=prim_cell[0], b=prim_cell[1], c=prim_cell[2])
208 atoms.calc = calc
210 if size != (1, 1, 1):
211 atoms = atoms.repeat(size)
212 return atoms
215def parse_symbols(symbols):
216 """Return `sumbols` as a sequence of element symbols."""
217 if isinstance(symbols, str):
218 symbols = string2symbols(symbols)
219 return symbols