r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1# Copyright (C) 2010, Jesper Friis

2# (see accompanying license files for details).

3"""Definition of the Spacegroup class.

5This module only depends on NumPy and the space group database.

6"""

8import os

9import warnings

10from functools import total_ordering

11from typing import Union

13import numpy as np

15__all__ = ['Spacegroup']

18class SpacegroupError(Exception):

19 """Base exception for the spacegroup module."""

20 pass

23class SpacegroupNotFoundError(SpacegroupError):

24 """Raised when given space group cannot be found in data base."""

25 pass

28class SpacegroupValueError(SpacegroupError):

29 """Raised when arguments have invalid value."""

30 pass

33# Type alias

34_SPACEGROUP = Union[int, str, 'Spacegroup']

37@total_ordering

38class Spacegroup:

39 """A space group class.

41 The instances of Spacegroup describes the symmetry operations for

42 the given space group.

44 Example:

46 >>> from ase.spacegroup import Spacegroup

47 >>>

48 >>> sg = Spacegroup(225)

49 >>> print('Space group', sg.no, sg.symbol)

50 Space group 225 F m -3 m

51 >>> sg.scaled_primitive_cell

52 array([[ 0. , 0.5, 0.5],

53 [ 0.5, 0. , 0.5],

54 [ 0.5, 0.5, 0. ]])

55 >>> sites, kinds = sg.equivalent_sites([[0,0,0]])

56 >>> sites

57 array([[ 0. , 0. , 0. ],

58 [ 0. , 0.5, 0.5],

59 [ 0.5, 0. , 0.5],

60 [ 0.5, 0.5, 0. ]])

61 """

62 no = property(

63 lambda self: self._no,

64 doc='Space group number in International Tables of Crystallography.')

65 symbol = property(

66 lambda self: self._symbol,

67 doc='Hermann-Mauguin (or international) symbol for the space group.')

68 setting = property(lambda self: self._setting,

69 doc='Space group setting. Either one or two.')

70 lattice = property(lambda self: self._symbol[0],

71 doc="""Lattice type:

73 P primitive

74 I body centering, h+k+l=2n

75 F face centering, h,k,l all odd or even

76 A,B,C single face centering, k+l=2n, h+l=2n, h+k=2n

77 R rhombohedral centering, -h+k+l=3n (obverse); h-k+l=3n (reverse)

78 """)

79 centrosymmetric = property(lambda self: self._centrosymmetric,

80 doc='Whether a center of symmetry exists.')

81 scaled_primitive_cell = property(

82 lambda self: self._scaled_primitive_cell,

83 doc='Primitive cell in scaled coordinates as a matrix with the '

84 'primitive vectors along the rows.')

85 reciprocal_cell = property(

86 lambda self: self._reciprocal_cell,

87 doc='Tree Miller indices that span all kinematically non-forbidden '

88 'reflections as a matrix with the Miller indices along the rows.')

89 nsubtrans = property(lambda self: len(self._subtrans),

90 doc='Number of cell-subtranslation vectors.')

92 def _get_nsymop(self):

93 """Returns total number of symmetry operations."""

94 if self.centrosymmetric:

95 return 2 * len(self._rotations) * len(self._subtrans)

96 else:

97 return len(self._rotations) * len(self._subtrans)

99 nsymop = property(_get_nsymop, doc='Total number of symmetry operations.')

100 subtrans = property(

101 lambda self: self._subtrans,

102 doc='Translations vectors belonging to cell-sub-translations.')

103 rotations = property(

104 lambda self: self._rotations,

105 doc='Symmetry rotation matrices. The invertions are not included '

106 'for centrosymmetrical crystals.')

107 translations = property(

108 lambda self: self._translations,

109 doc='Symmetry translations. The invertions are not included '

110 'for centrosymmetrical crystals.')

112 def __init__(self, spacegroup: _SPACEGROUP, setting=1, datafile=None):

113 """Returns a new Spacegroup instance.

115 Parameters:

117 spacegroup : int | string | Spacegroup instance

118 The space group number in International Tables of

119 Crystallography or its Hermann-Mauguin symbol. E.g.

120 spacegroup=225 and spacegroup='F m -3 m' are equivalent.

121 setting : 1 | 2

122 Some space groups have more than one setting. `setting`

123 determines Which of these should be used.

124 datafile : None | string

125 Path to database file. If `None`, the the default database

126 will be used.

127 """

128 if isinstance(spacegroup, Spacegroup):

129 for k, v in spacegroup.__dict__.items():

130 setattr(self, k, v)

131 return

132 if not datafile:

133 datafile = get_datafile()

134 with open(datafile, 'r') as fd:

137 def __repr__(self):

138 return 'Spacegroup(%d, setting=%d)' % (self.no, self.setting)

140 def todict(self):

141 return {'number': self.no, 'setting': self.setting}

143 def __str__(self):

144 """Return a string representation of the space group data in

145 the same format as found the database."""

146 retval = []

147 # no, symbol

148 retval.append('%-3d %s\n' % (self.no, self.symbol))

149 # setting

150 retval.append(' setting %d\n' % (self.setting))

151 # centrosymmetric

152 retval.append(' centrosymmetric %d\n' % (self.centrosymmetric))

153 # primitive vectors

154 retval.append(' primitive vectors\n')

155 for i in range(3):

156 retval.append(' ')

157 for j in range(3):

158 retval.append(' %13.10f' % (self.scaled_primitive_cell[i, j]))

159 retval.append('\n')

160 # primitive reciprocal vectors

161 retval.append(' reciprocal vectors\n')

162 for i in range(3):

163 retval.append(' ')

164 for j in range(3):

165 retval.append(' %3d' % (self.reciprocal_cell[i, j]))

166 retval.append('\n')

167 # sublattice

168 retval.append(' %d subtranslations\n' % self.nsubtrans)

169 for i in range(self.nsubtrans):

170 retval.append(' ')

171 for j in range(3):

172 retval.append(' %13.10f' % (self.subtrans[i, j]))

173 retval.append('\n')

174 # symmetry operations

175 nrot = len(self.rotations)

176 retval.append(' %d symmetry operations (rot+trans)\n' % nrot)

177 for i in range(nrot):

178 retval.append(' ')

179 for j in range(3):

180 retval.append(' ')

181 for k in range(3):

182 retval.append(' %2d' % (self.rotations[i, j, k]))

183 retval.append(' ')

184 for j in range(3):

185 retval.append(' %13.10f' % self.translations[i, j])

186 retval.append('\n')

187 retval.append('\n')

188 return ''.join(retval)

190 def __eq__(self, other):

191 return self.no == other.no and self.setting == other.setting

193 def __ne__(self, other):

194 return not self.__eq__(other)

196 def __lt__(self, other):

197 return self.no < other.no or (self.no == other.no

198 and self.setting < other.setting)

200 def __index__(self):

201 return self.no

203 __int__ = __index__

205 def get_symop(self):

206 """Returns all symmetry operations (including inversions and

207 subtranslations) as a sequence of (rotation, translation)

208 tuples."""

209 symop = []

210 parities = [1]

211 if self.centrosymmetric:

212 parities.append(-1)

213 for parity in parities:

214 for subtrans in self.subtrans:

215 for rot, trans in zip(self.rotations, self.translations):

216 newtrans = np.mod(trans + subtrans, 1)

217 symop.append((parity * rot, newtrans))

218 return symop

220 def get_op(self):

221 """Returns all symmetry operations (including inversions and

222 subtranslations), but unlike get_symop(), they are returned as

223 two ndarrays."""

224 if self.centrosymmetric:

225 rot = np.tile(np.vstack((self.rotations, -self.rotations)),

226 (self.nsubtrans, 1, 1))

227 trans = np.tile(np.vstack((self.translations, -self.translations)),

228 (self.nsubtrans, 1))

229 trans += np.repeat(self.subtrans, 2 * len(self.rotations), axis=0)

230 trans = np.mod(trans, 1)

231 else:

232 rot = np.tile(self.rotations, (self.nsubtrans, 1, 1))

233 trans = np.tile(self.translations, (self.nsubtrans, 1))

234 trans += np.repeat(self.subtrans, len(self.rotations), axis=0)

235 trans = np.mod(trans, 1)

236 return rot, trans

238 def get_rotations(self):

239 """Return all rotations, including inversions for

240 centrosymmetric crystals."""

241 if self.centrosymmetric:

242 return np.vstack((self.rotations, -self.rotations))

243 else:

244 return self.rotations

246 def equivalent_reflections(self, hkl):

247 """Return all equivalent reflections to the list of Miller indices

248 in hkl.

250 Example:

252 >>> from ase.spacegroup import Spacegroup

253 >>> sg = Spacegroup(225) # fcc

254 >>> sg.equivalent_reflections([[0, 0, 2]])

255 array([[ 0, 0, -2],

256 [ 0, -2, 0],

257 [-2, 0, 0],

258 [ 2, 0, 0],

259 [ 0, 2, 0],

260 [ 0, 0, 2]])

261 """

262 hkl = np.array(hkl, dtype='int', ndmin=2)

263 rot = self.get_rotations()

264 n, nrot = len(hkl), len(rot)

265 R = rot.transpose(0, 2, 1).reshape((3 * nrot, 3)).T

266 refl = np.dot(hkl, R).reshape((n * nrot, 3))

267 ind = np.lexsort(refl.T)

268 refl = refl[ind]

269 diff = np.diff(refl, axis=0)

273 def equivalent_lattice_points(self, uvw):

274 """Return all lattice points equivalent to any of the lattice points

275 in `uvw` with respect to rotations only.

277 Only equivalent lattice points that conserves the distance to

278 origo are included in the output (making this a kind of real

279 space version of the equivalent_reflections() method).

281 Example:

283 >>> from ase.spacegroup import Spacegroup

284 >>> sg = Spacegroup(225) # fcc

285 >>> sg.equivalent_lattice_points([[0, 0, 2]])

286 array([[ 0, 0, -2],

287 [ 0, -2, 0],

288 [-2, 0, 0],

289 [ 2, 0, 0],

290 [ 0, 2, 0],

291 [ 0, 0, 2]])

293 """

294 uvw = np.array(uvw, ndmin=2)

295 rot = self.get_rotations()

296 n, nrot = len(uvw), len(rot)

297 directions = np.dot(uvw, rot).reshape((n * nrot, 3))

298 ind = np.lexsort(directions.T)

299 directions = directions[ind]

300 diff = np.diff(directions, axis=0)

304 def symmetry_normalised_reflections(self, hkl):

305 """Returns an array of same size as *hkl*, containing the

306 corresponding symmetry-equivalent reflections of lowest

307 indices.

309 Example:

311 >>> from ase.spacegroup import Spacegroup

312 >>> sg = Spacegroup(225) # fcc

313 >>> sg.symmetry_normalised_reflections([[2, 0, 0], [0, 2, 0]])

314 array([[ 0, 0, -2],

315 [ 0, 0, -2]])

316 """

317 hkl = np.array(hkl, dtype=int, ndmin=2)

318 normalised = np.empty(hkl.shape, int)

319 R = self.get_rotations().transpose(0, 2, 1)

320 for i, g in enumerate(hkl):

321 gsym = np.dot(R, g)

322 j = np.lexsort(gsym.T)[0]

323 normalised[i, :] = gsym[j]

324 return normalised

326 def unique_reflections(self, hkl):

327 """Returns a subset *hkl* containing only the symmetry-unique

328 reflections.

330 Example:

332 >>> from ase.spacegroup import Spacegroup

333 >>> sg = Spacegroup(225) # fcc

334 >>> sg.unique_reflections([[ 2, 0, 0],

335 ... [ 0, -2, 0],

336 ... [ 2, 2, 0],

337 ... [ 0, -2, -2]])

338 array([[2, 0, 0],

339 [2, 2, 0]])

340 """

341 hkl = np.array(hkl, dtype=int, ndmin=2)

342 hklnorm = self.symmetry_normalised_reflections(hkl)

343 perm = np.lexsort(hklnorm.T)

344 iperm = perm.argsort()

350 def equivalent_sites(self,

351 scaled_positions,

352 onduplicates='error',

353 symprec=1e-3,

354 occupancies=None):

355 """Returns the scaled positions and all their equivalent sites.

357 Parameters:

359 scaled_positions: list | array

360 List of non-equivalent sites given in unit cell coordinates.

362 occupancies: list | array, optional (default=None)

363 List of occupancies corresponding to the respective sites.

365 onduplicates : 'keep' | 'replace' | 'warn' | 'error'

366 Action if `scaled_positions` contain symmetry-equivalent

367 positions of full occupancy:

369 'keep'

371 'replace'

372 replace

373 'warn'

374 like 'keep', but issue an UserWarning

375 'error'

376 raises a SpacegroupValueError

378 symprec: float

379 Minimum "distance" betweed two sites in scaled coordinates

380 before they are counted as the same site.

382 Returns:

384 sites: array

385 A NumPy array of equivalent sites.

386 kinds: list

387 A list of integer indices specifying which input site is

388 equivalent to the corresponding returned site.

390 Example:

392 >>> from ase.spacegroup import Spacegroup

393 >>> sg = Spacegroup(225) # fcc

394 >>> sites, kinds = sg.equivalent_sites([[0, 0, 0], [0.5, 0.0, 0.0]])

395 >>> sites

396 array([[ 0. , 0. , 0. ],

397 [ 0. , 0.5, 0.5],

398 [ 0.5, 0. , 0.5],

399 [ 0.5, 0.5, 0. ],

400 [ 0.5, 0. , 0. ],

401 [ 0. , 0.5, 0. ],

402 [ 0. , 0. , 0.5],

403 [ 0.5, 0.5, 0.5]])

404 >>> kinds

405 [0, 0, 0, 0, 1, 1, 1, 1]

406 """

407 kinds = []

408 sites = []

410 scaled = np.array(scaled_positions, ndmin=2)

412 for kind, pos in enumerate(scaled):

413 for rot, trans in self.get_symop():

414 site = np.mod(np.dot(rot, pos) + trans, 1.)

415 if not sites:

416 sites.append(site)

417 kinds.append(kind)

418 continue

419 t = site - sites

421 (abs(t) < symprec) | (abs(abs(t) - 1.0) < symprec), axis=1)

424 for ind in inds:

425 # then we would just add the same thing again -> skip

426 if kinds[ind] == kind:

427 pass

428 elif onduplicates == 'keep':

429 pass

430 elif onduplicates == 'replace':

431 kinds[ind] = kind

432 elif onduplicates == 'warn':

433 warnings.warn('scaled_positions %d and %d '

434 'are equivalent' %

435 (kinds[ind], kind))

436 elif onduplicates == 'error':

437 raise SpacegroupValueError(

438 'scaled_positions %d and %d are equivalent' %

439 (kinds[ind], kind))

440 else:

441 raise SpacegroupValueError(

442 'Argument "onduplicates" must be one of: '

443 '"keep", "replace", "warn" or "error".')

444 else:

445 sites.append(site)

446 kinds.append(kind)

448 return np.array(sites), kinds

450 def symmetry_normalised_sites(self,

451 scaled_positions,

452 map_to_unitcell=True):

453 """Returns an array of same size as *scaled_positions*,

454 containing the corresponding symmetry-equivalent sites of

455 lowest indices.

457 If *map_to_unitcell* is true, the returned positions are all

458 mapped into the unit cell, i.e. lattice translations are

459 included as symmetry operator.

461 Example:

463 >>> from ase.spacegroup import Spacegroup

464 >>> sg = Spacegroup(225) # fcc

465 >>> sg.symmetry_normalised_sites([[0.0, 0.5, 0.5], [1.0, 1.0, 0.0]])

466 array([[ 0., 0., 0.],

467 [ 0., 0., 0.]])

468 """

469 scaled = np.array(scaled_positions, ndmin=2)

470 normalised = np.empty(scaled.shape, float)

471 rot, trans = self.get_op()

472 for i, pos in enumerate(scaled):

473 sympos = np.dot(rot, pos) + trans

474 if map_to_unitcell:

475 # Must be done twice, see the scaled_positions.py test

476 sympos %= 1.0

477 sympos %= 1.0

478 j = np.lexsort(sympos.T)[0]

479 normalised[i, :] = sympos[j]

480 return normalised

482 def unique_sites(self,

483 scaled_positions,

484 symprec=1e-3,

486 map_to_unitcell=True):

487 """Returns a subset of *scaled_positions* containing only the

488 symmetry-unique positions. If *output_mask* is True, a boolean

489 array masking the subset is also returned.

491 If *map_to_unitcell* is true, all sites are first mapped into

492 the unit cell making e.g. [0, 0, 0] and [1, 0, 0] equivalent.

494 Example:

496 >>> from ase.spacegroup import Spacegroup

497 >>> sg = Spacegroup(225) # fcc

498 >>> sg.unique_sites([[0.0, 0.0, 0.0],

499 ... [0.5, 0.5, 0.0],

500 ... [1.0, 0.0, 0.0],

501 ... [0.5, 0.0, 0.0]])

502 array([[ 0. , 0. , 0. ],

503 [ 0.5, 0. , 0. ]])

504 """

505 scaled = np.array(scaled_positions, ndmin=2)

506 symnorm = self.symmetry_normalised_sites(scaled, map_to_unitcell)

507 perm = np.lexsort(symnorm.T)

508 iperm = perm.argsort()

509 xmask = np.abs(np.diff(symnorm[perm], axis=0)).max(axis=1) > symprec

514 else:

517 def tag_sites(self, scaled_positions, symprec=1e-3):

518 """Returns an integer array of the same length as *scaled_positions*,

519 tagging all equivalent atoms with the same index.

521 Example:

523 >>> from ase.spacegroup import Spacegroup

524 >>> sg = Spacegroup(225) # fcc

525 >>> sg.tag_sites([[0.0, 0.0, 0.0],

526 ... [0.5, 0.5, 0.0],

527 ... [1.0, 0.0, 0.0],

528 ... [0.5, 0.0, 0.0]])

529 array([0, 0, 0, 1])

530 """

531 scaled = np.array(scaled_positions, ndmin=2)

532 scaled %= 1.0

533 scaled %= 1.0

534 tags = -np.ones((len(scaled), ), dtype=int)

535 mask = np.ones((len(scaled), ), dtype=bool)

536 rot, trans = self.get_op()

537 i = 0

540 sympos = np.dot(rot, pos) + trans

541 # Must be done twice, see the scaled_positions.py test

542 sympos %= 1.0

543 sympos %= 1.0

544 m = ~np.all(np.any(np.abs(scaled[np.newaxis, :, :] -

545 sympos[:, np.newaxis, :]) > symprec,

546 axis=2),

547 axis=0)

548 assert not np.any((~mask) & m)

549 tags[m] = i

551 i += 1

552 return tags

555def get_datafile():

556 """Return default path to datafile."""

557 return os.path.join(os.path.dirname(__file__), 'spacegroup.dat')

560def format_symbol(symbol):

561 """Returns well formatted Hermann-Mauguin symbol as extected by

562 the database, by correcting the case and adding missing or

563 removing dublicated spaces."""

564 fixed = []

565 s = symbol.strip()

566 s = s[0].upper() + s[1:].lower()

567 for c in s:

568 if c.isalpha():

569 if len(fixed) and fixed[-1] == '/':

570 fixed.append(c)

571 else:

572 fixed.append(' ' + c + ' ')

573 elif c.isspace():

574 fixed.append(' ')

575 elif c.isdigit():

576 fixed.append(c)

577 elif c == '-':

578 fixed.append(' ' + c)

579 elif c == '/':

580 fixed.append(c)

581 s = ''.join(fixed).strip()

582 return ' '.join(s.split())

585# Functions for parsing the database. They are moved outside the

586# Spacegroup class in order to make it easier to later implement

587# caching to avoid reading the database each time a new Spacegroup

588# instance is created.

591def _skip_to_blank(f, spacegroup, setting):

592 """Read lines from f until a blank line is encountered."""

593 while True:

595 if not line:

596 raise SpacegroupNotFoundError(

598 % (spacegroup, setting))

599 if not line.strip():

600 break

603def _skip_to_nonblank(f, spacegroup, setting):

604 """Read lines from f until a nonblank line not starting with a

605 hash (#) is encountered and returns this and the next line."""

606 while True:

608 if not line1:

609 raise SpacegroupNotFoundError(

610 'invalid spacegroup %s, setting %i not found in data base' %

611 (spacegroup, setting))

612 line1.strip()

613 if line1 and not line1.startswith('#'):

615 break

616 return line1, line2

619def _read_datafile_entry(spg, no, symbol, setting, f):

620 """Read space group data from f to spg."""

622 floats = {'0.0': 0.0, '1.0': 1.0, '0': 0.0, '1': 1.0, '-1': -1.0}

623 for n, d in [(1, 2), (1, 3), (2, 3), (1, 4), (3, 4), (1, 6), (5, 6)]:

624 floats['{0}/{1}'.format(n, d)] = n / d

625 floats['-{0}/{1}'.format(n, d)] = -n / d

627 spg._no = no

628 spg._symbol = symbol.strip()

629 spg._setting = setting

631 # primitive vectors

633 spg._scaled_primitive_cell = np.array(

634 [[float(floats.get(s, s)) for s in f.readline().split()]

635 for i in range(3)],

636 dtype=float)

637 # primitive reciprocal vectors

639 spg._reciprocal_cell = np.array([[int(i) for i in f.readline().split()]

640 for i in range(3)],

641 dtype=int)

642 # subtranslations

644 spg._subtrans = np.array(

645 [[float(floats.get(t, t)) for t in f.readline().split()]

646 for i in range(spg._nsubtrans)],

647 dtype=float)

648 # symmetry operations

650 symop = np.array([[float(floats.get(s, s)) for s in f.readline().split()]

651 for i in range(nsym)],

652 dtype=float)

653 spg._nsymop = nsym

654 spg._rotations = np.array(symop[:, :9].reshape((nsym, 3, 3)), dtype=int)

655 spg._translations = symop[:, 9:]

659 if isinstance(spacegroup, int):

660 pass

661 elif isinstance(spacegroup, str):

662 spacegroup = ' '.join(spacegroup.strip().split())

663 compact_spacegroup = ''.join(spacegroup.split())

664 else:

665 raise SpacegroupValueError('`spacegroup` must be of type int or str')

666 while True:

667 line1, line2 = _skip_to_nonblank(f, spacegroup, setting)

668 _no, _symbol = line1.strip().split(None, 1)

669 _symbol = format_symbol(_symbol)

670 compact_symbol = ''.join(_symbol.split())

671 _setting = int(line2.strip().split()[1])

672 _no = int(_no)

674 condition = (

675 (isinstance(spacegroup, int) and _no == spacegroup

676 and _setting == setting)

677 or (isinstance(spacegroup, str)

678 and compact_symbol == compact_spacegroup) and

679 (setting is None or _setting == setting))

681 if condition:

682 _read_datafile_entry(spg, _no, _symbol, _setting, f)

683 break

684 else:

685 _skip_to_blank(f, spacegroup, setting)

688def parse_sitesym_element(element):

689 """Parses one element from a single site symmetry in the form used

690 by the International Tables.

692 Examples:

694 >>> parse_sitesym_element("x")

695 ([(0, 1)], 0.0)

696 >>> parse_sitesym_element("-1/2-y")

697 ([(1, -1)], -0.5)

698 >>> parse_sitesym_element("z+0.25")

699 ([(2, 1)], 0.25)

700 >>> parse_sitesym_element("x-z+0.5")

701 ([(0, 1), (2, -1)], 0.5)

705 Parameters

706 ----------

708 element: str

709 Site symmetry like "x" or "-y+1/4" or "0.5+z".

712 Returns

713 -------

715 list[tuple[int, int]]

716 Rotation information in the form '(index, sign)' where index is

717 0 for "x", 1 for "y" and 2 for "z" and sign is '1' for a positive

718 entry and '-1' for a negative entry. E.g. "x" is '(0, 1)' and

719 "-z" is (2, -1).

721 float

722 Translation information in fractional space. E.g. "-1/4" is

723 '-0.25' and "1/2" is '0.5' and "0.75" is '0.75'.

726 """

727 element = element.lower()

728 is_positive = True

729 is_frac = False

730 sng_trans = None

731 fst_trans = []

732 snd_trans = []

733 rot = []

735 for char in element:

736 if char == "+":

737 is_positive = True

738 elif char == "-":

739 is_positive = False

740 elif char == "/":

741 is_frac = True

742 elif char in "xyz":

743 rot.append((ord(char) - ord("x"), 1 if is_positive else -1))

744 elif char.isdigit() or char == ".":

745 if sng_trans is None:

746 sng_trans = 1.0 if is_positive else -1.0

747 if is_frac:

748 snd_trans.append(char)

749 else:

750 fst_trans.append(char)

752 trans = 0.0 if not fst_trans else (sng_trans * float("".join(fst_trans)))

753 if is_frac:

754 trans /= float("".join(snd_trans))

756 return rot, trans

759def parse_sitesym_single(sym, out_rot, out_trans, sep=",",

760 force_positive_translation=False):

761 """Parses a single site symmetry in the form used by International

762 Tables and overwrites 'out_rot' and 'out_trans' with data.

764 Parameters

765 ----------

767 sym: str

768 Site symmetry in the form used by International Tables

769 (e.g. "x,y,z", "y-1/2,x,-z").

771 out_rot: np.array

772 A 3x3-integer array representing rotations (changes are made inplace).

774 out_rot: np.array

775 A 3-float array representing translations (changes are made inplace).

777 sep: str

778 String separator ("," in "x,y,z").

780 force_positive_translation: bool

781 Forces fractional translations to be between 0 and 1 (otherwise

782 negative values might be accepted). Defaults to 'False'.

785 Returns

786 -------

788 Nothing is returned: 'out_rot' and 'out_trans' are changed inplace.

791 """

792 out_rot[:] = 0.0

793 out_trans[:] = 0.0

795 for i, element in enumerate(sym.split(sep)):

796 e_rot_list, e_trans = parse_sitesym_element(element)

797 for rot_idx, rot_sgn in e_rot_list:

798 out_rot[i][rot_idx] = rot_sgn

799 out_trans[i] = \

800 (e_trans % 1.0) if force_positive_translation else e_trans

803def parse_sitesym(symlist, sep=',', force_positive_translation=False):

804 """Parses a sequence of site symmetries in the form used by

805 International Tables and returns corresponding rotation and

806 translation arrays.

808 Example:

810 >>> symlist = [

811 ... 'x,y,z',

812 ... '-y+1/2,x+1/2,z',

813 ... '-y,-x,-z',

814 ... 'x-1/4, y-1/4, -z'

815 ... ]

816 >>> rot, trans = parse_sitesym(symlist)

817 >>> rot

818 array([[[ 1, 0, 0],

819 [ 0, 1, 0],

820 [ 0, 0, 1]],

821 <BLANKLINE>

822 [[ 0, -1, 0],

823 [ 1, 0, 0],

824 [ 0, 0, 1]],

825 <BLANKLINE>

826 [[ 0, -1, 0],

827 [-1, 0, 0],

828 [ 0, 0, -1]],

829 <BLANKLINE>

830 [[ 1, 0, 0],

831 [ 0, 1, 0],

832 [ 0, 0, -1]]])

833 >>> trans

834 array([[ 0. , 0. , 0. ],

835 [ 0.5 , 0.5 , 0. ],

836 [ 0. , 0. , 0. ],

837 [-0.25, -0.25, 0. ]])

838 """

840 nsym = len(symlist)

841 rot = np.zeros((nsym, 3, 3), dtype='int')

842 trans = np.zeros((nsym, 3))

844 for i, sym in enumerate(symlist):

845 parse_sitesym_single(

846 sym, rot[i], trans[i], sep=sep,

847 force_positive_translation=force_positive_translation)

849 return rot, trans

852def spacegroup_from_data(no=None,

853 symbol=None,

854 setting=None,

855 centrosymmetric=None,

856 scaled_primitive_cell=None,

857 reciprocal_cell=None,

858 subtrans=None,

859 sitesym=None,

860 rotations=None,

861 translations=None,

862 datafile=None):

863 """Manually create a new space group instance. This might be

864 useful when reading crystal data with its own spacegroup

865 definitions."""

866 if no is not None and setting is not None:

867 spg = Spacegroup(no, setting, datafile)

868 elif symbol is not None:

869 spg = Spacegroup(symbol, None, datafile)

870 else:

871 raise SpacegroupValueError('either *no* and *setting* '

872 'or *symbol* must be given')

873 if not isinstance(sitesym, list):

874 raise TypeError('sitesym must be a list')

876 have_sym = False

877 if centrosymmetric is not None:

878 spg._centrosymmetric = bool(centrosymmetric)

879 if scaled_primitive_cell is not None:

880 spg._scaled_primitive_cell = np.array(scaled_primitive_cell)

881 if reciprocal_cell is not None:

882 spg._reciprocal_cell = np.array(reciprocal_cell)

883 if subtrans is not None:

884 spg._subtrans = np.atleast_2d(subtrans)

885 spg._nsubtrans = spg._subtrans.shape[0]

886 if sitesym is not None:

887 spg._rotations, spg._translations = parse_sitesym(sitesym)

888 have_sym = True

889 if rotations is not None:

890 spg._rotations = np.atleast_3d(rotations)

891 have_sym = True

892 if translations is not None:

893 spg._translations = np.atleast_2d(translations)

894 have_sym = True

895 if have_sym:

896 if spg._rotations.shape[0] != spg._translations.shape[0]:

897 raise SpacegroupValueError('inconsistent number of rotations and '

898 'translations')

899 spg._nsymop = spg._rotations.shape[0]

900 return spg

903def get_spacegroup(atoms, symprec=1e-5):

904 """Determine the spacegroup to which belongs the Atoms object.

906 This requires spglib: https://atztogo.github.io/spglib/ .

908 Parameters:

910 atoms: Atoms object

911 Types, positions and unit-cell.

912 symprec: float

913 Symmetry tolerance, i.e. distance tolerance in Cartesian

914 coordinates to find crystal symmetry.

916 The Spacegroup object is returned.

917 """

919 # Example:

920 # (We don't include the example in docstring to appease doctests

921 # when import fails)

922 # >>> from ase.build import bulk

923 # >>> atoms = bulk("Cu", "fcc", a=3.6, cubic=True)

924 # >>> sg = get_spacegroup(atoms)

925 # >>> sg

926 # Spacegroup(225, setting=1)

927 # >>> sg.no

928 # 225

930 import spglib

932 sg = spglib.get_spacegroup((atoms.get_cell(), atoms.get_scaled_positions(),

933 atoms.get_atomic_numbers()),

934 symprec=symprec)

935 if sg is None: