Source code for ase.spacegroup.xtal

from __future__ import print_function
# Copyright (C) 2010, Jesper Friis
# (see accompanying license files for details).

"""
A module for ASE for simple creation of crystalline structures from
knowledge of the space group.

"""

import numpy as np

import ase
from ase.symbols import string2symbols
from ase.spacegroup import Spacegroup
from ase.geometry import cellpar_to_cell
from ase.utils import basestring

__all__ = ['crystal']


[docs]def crystal(symbols=None, basis=None, occupancies=None, spacegroup=1, setting=1, cell=None, cellpar=None, ab_normal=(0, 0, 1), a_direction=None, size=(1, 1, 1), onduplicates='warn', symprec=0.001, pbc=True, primitive_cell=False, **kwargs): """Create an Atoms instance for a conventional unit cell of a space group. Parameters: symbols : str | sequence of str | sequence of Atom | Atoms Element symbols of the unique sites. Can either be a string formula or a sequence of element symbols. E.g. ('Na', 'Cl') and 'NaCl' are equivalent. Can also be given as a sequence of Atom objects or an Atoms object. basis : list of scaled coordinates Positions of the unique sites corresponding to symbols given either as scaled positions or through an atoms instance. Not needed if *symbols* is a sequence of Atom objects or an Atoms object. occupancies : list of site occupancies Occupancies of the unique sites. Defaults to 1.0 and thus no mixed occupancies are considered if not explicitly asked for. If occupancies are given, the most dominant species will yield the atomic number. spacegroup : int | string | Spacegroup instance Space group given either as its number in International Tables or as its Hermann-Mauguin symbol. setting : 1 | 2 Space group setting. cell : 3x3 matrix Unit cell vectors. cellpar : [a, b, c, alpha, beta, gamma] Cell parameters with angles in degree. Is not used when `cell` is given. ab_normal : vector Is used to define the orientation of the unit cell relative to the Cartesian system when `cell` is not given. It is the normal vector of the plane spanned by a and b. a_direction : vector Defines the orientation of the unit cell a vector. a will be parallel to the projection of `a_direction` onto the a-b plane. size : 3 positive integers How many times the conventional unit cell should be repeated in each direction. onduplicates : 'keep' | 'replace' | 'warn' | 'error' Action if `basis` contain symmetry-equivalent positions: 'keep' - ignore additional symmetry-equivalent positions 'replace' - replace 'warn' - like 'keep', but issue an UserWarning 'error' - raises a SpacegroupValueError symprec : float Minimum "distance" betweed two sites in scaled coordinates before they are counted as the same site. pbc : one or three bools Periodic boundary conditions flags. Examples: True, False, 0, 1, (1, 1, 0), (True, False, False). Default is True. primitive_cell : bool Wheter to return the primitive instead of the conventional unit cell. Keyword arguments: All additional keyword arguments are passed on to the Atoms constructor. Currently, probably the most useful additional keyword arguments are `info`, `constraint` and `calculator`. Examples: Two diamond unit cells (space group number 227) >>> diamond = crystal('C', [(0,0,0)], spacegroup=227, ... cellpar=[3.57, 3.57, 3.57, 90, 90, 90], size=(2,1,1)) >>> ase.view(diamond) # doctest: +SKIP A CoSb3 skutterudite unit cell containing 32 atoms >>> skutterudite = crystal(('Co', 'Sb'), ... basis=[(0.25,0.25,0.25), (0.0, 0.335, 0.158)], ... spacegroup=204, cellpar=[9.04, 9.04, 9.04, 90, 90, 90]) >>> len(skutterudite) 32 """ sg = Spacegroup(spacegroup, setting) if (not isinstance(symbols, basestring) and hasattr(symbols, '__getitem__') and len(symbols) > 0 and isinstance(symbols[0], ase.Atom)): symbols = ase.Atoms(symbols) if isinstance(symbols, ase.Atoms): basis = symbols symbols = basis.get_chemical_symbols() if isinstance(basis, ase.Atoms): basis_coords = basis.get_scaled_positions() if cell is None and cellpar is None: cell = basis.cell if symbols is None: symbols = basis.get_chemical_symbols() else: basis_coords = np.array(basis, dtype=float, copy=False, ndmin=2) if occupancies is not None: # find the identical sites (doubt this will ever become performance # critical -- and if so, maybe use the neighborlist module) ident_sites = np.linalg.norm(basis_coords[:, None]-basis_coords[None,:], axis=-1) < symprec # for convenience, only look at the upper triangle (i.e. set the remainder to False) np.fill_diagonal(ident_sites, False) # this is already the reduced form... unique = [i for i in range(len(basis_coords))] occupancies_dict = {} for i, mask in enumerate(ident_sites): if i not in unique: continue occ = {symbols[i] : occupancies[i]} for close in np.nonzero(mask): for j in close: occ.update({symbols[j] : occupancies[j]}) occupancies_dict[i] = occ.copy() sites, kinds = sg.equivalent_sites(basis_coords, onduplicates=onduplicates, symprec=symprec) symbols = parse_symbols(symbols) if occupancies is None: symbols = [symbols[i] for i in kinds] else: # make sure that we put the dominant species there symbols = [sorted(occupancies_dict[i].items(), key=lambda x : x[1])[-1][0] for i in kinds] if cell is None: cell = cellpar_to_cell(cellpar, ab_normal, a_direction) info = dict(spacegroup=sg) if primitive_cell: info['unit_cell'] = 'primitive' else: info['unit_cell'] = 'conventional' if 'info' in kwargs: info.update(kwargs['info']) if occupancies is not None: info['occupancy'] = occupancies_dict kwargs['info'] = info atoms = ase.Atoms(symbols, scaled_positions=sites, cell=cell, # use tags to identify sites, and in particular the occupancy tags=kinds, pbc=pbc, **kwargs) if isinstance(basis, ase.Atoms): for name in basis.arrays: if not atoms.has(name): array = basis.get_array(name) atoms.new_array(name, [array[i] for i in kinds], dtype=array.dtype, shape=array.shape[1:]) if primitive_cell: from ase.build import cut prim_cell = sg.scaled_primitive_cell atoms = cut(atoms, a=prim_cell[0], b=prim_cell[1], c=prim_cell[2]) if size != (1, 1, 1): atoms = atoms.repeat(size) return atoms
def parse_symbols(symbols): """Return `sumbols` as a sequence of element symbols.""" if isinstance(symbols, basestring): symbols = string2symbols(symbols) return symbols