Primary data structures¶
Here are the primary data structures available in crystals
. You will work directly with these structures,
maybe even either creating them directly.
Crystal¶
- class crystals.Crystal(unitcell: Iterable[Union[Atom, AtomicStructure]], lattice_vectors: Union[Sequence[Sequence[Sequence[Sequence[Sequence[Any]]]]], _SupportsArray[dtype], Sequence[_SupportsArray[dtype]], Sequence[Sequence[_SupportsArray[dtype]]], Sequence[Sequence[Sequence[_SupportsArray[dtype]]]], Sequence[Sequence[Sequence[Sequence[_SupportsArray[dtype]]]]], bool, int, float, complex, str, bytes, Sequence[Union[bool, int, float, complex, str, bytes]], Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]], Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]], Sequence[Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]]]], source: Optional[str] = None, **kwargs)¶
Bases:
AtomicStructure
,Lattice
The
Crystal
class is a set-like container that represent crystalline structures. In addition to constructing theCrystal
object yourself, other constructors are also available (and preferred):Crystal.from_cif()
: create an instance from a CIF file;Crystal.from_pdb()
: create an instance from a Protein Data Bank entry;Crystal.from_database()
: create an instance from the internal database of CIF files;Crystal.from_cod()
: create an instance from a Crystallography Open Database entry;Crystal.from_mp()
: create an instance from the Materials Project database;Crystal.from_pwscf()
: create an instance from the output of the PWSCF program;Crystal.from_ase()
: create an instance from anase.Atoms
instance;Crystal.from_poscar()
: create an instance from VASP POSCAR files.
- Parameters
unitcell (iterable of
Atom
orAtomicStructure
) – Unit cell atoms or substructures. It is assumed that the atoms are in fractional coordinates.lattice_vectors (iterable of array_like) – Lattice vectors. If
lattice_vectors
is provided as a 3x3 array, it is assumed that each lattice vector is a row.source (str or None, optional) – Provenance, e.g. filename. Only used for bookkeeping.
- asymmetric_cell() Set[Atom] ¶
Calculates the asymmetric cell that generates the crystal unit cell.
New in version 1.2.0.
- classmethod from_ase(atoms, **kwargs)¶
Returns a Crystal object created from an ASE Atoms object. Keyword arguments are passed to the class constructor.
- Parameters
atoms (ase.Atoms) – Atoms group.
- classmethod from_cif(path: PathLike, **kwargs)¶
Returns a Crystal object created from a CIF 1.0, 1.1 or 2.0 file. Keyword arguments are passed to the Crystal constructor.
- Parameters
path (path-like) – File path
- classmethod from_cod(num: int, revision: Optional[int] = None, download_dir: Optional[PathLike] = None, overwrite: bool = False, **kwargs)¶
Returns a Crystal object built from the Crystallography Open Database. Keyword arguments are passed to the class constructor.
- Parameters
num (int) – COD identification number.
revision (int or None, optional) – Revision number. If None (default), the latest revision is used.
download_dir (path-like object, optional) – Directory where to save the CIF file. Default is a local folder in the current directory
overwrite (bool, optional) – Whether or not to overwrite files in cache if they exist. If no revision number is provided, files will always be overwritten.
- classmethod from_database(name: str, **kwargs)¶
Returns a Crystal object create from the internal CIF database. Keyword arguments are passed to the class constructor.
- Parameters
name (path-like) – Name of the database entry. Available items can be retrieved from Crystal.builtins
- classmethod from_mp(query: str, api_key: Optional[str] = None, download_dir: Optional[PathLike] = None, overwrite: bool = False, **kwargs)¶
Returns a Crystal object built from the Materials Project. Keyword arguments are passed to the class constructor.
- Parameters
query (str) – The query can be a Materials Project material id (e.g., “mp-1234”), a formula, e.g. (“Fe2O3”), or a chemical system (“-” separated list of elemments, e.g., “Li-Fe-O”).
api_key (str or None, optional) – An API key for accessing the Materials Project REST interface. Please obtain your API key at https://www.materialsproject.org/dashboard. If None (default),
crystals
will look for your API key in the MATERIALS_PROJECT_API_KEY environment variable.download_dir (path-like object or None, optional) – Directory where to save the CIF file. This is used for caching.
overwrite (bool, optional) – Whether or not to overwrite files in cache if they exist. If True, a new file will be downloaded, possibly overwriting previously-downloaded file.
- classmethod from_pdb(ID: str, download_dir: Optional[PathLike] = None, overwrite: bool = False, **kwargs)¶
Returns a Crystal object created from a Protein DataBank entry. Keyword arguments are passed to the class constructor.
- Parameters
ID (str) – Protein DataBank identification. The correct .pdb file will be downloaded, cached and parsed.
download_dir (path-like object, optional) – Directory where to save the PDB file.
overwrite (bool, optional) – Whether or not to overwrite files in cache if they exist. If no revision number is provided, files will always be overwritten.
- classmethod from_poscar(path: PathLike, **kwargs)¶
Returns a Crystal object created from a VASP’s POSCAR file. Keyword arguments are passed to the class constructor.
New in version 1.4.0.
- Parameters
path (path-like) – File path
- classmethod from_pwscf(path: PathLike, **kwargs)¶
Returns a Crystal object created from an output file of PWSCF. Keyword arguments are passed to the class constructor.
- Parameters
path (path-like) – File path
- groupby(by: str, symprec: float = 0.01, angle_tolerance: float = -1.0) Dict[Any, Iterable[Atom]] ¶
Group unit cell atoms by some measure of site symmetry, for example Wyckoff letters or crystallographic orbits.
New in version 1.6.0.
- Parameters
by (str) – Measure of site symmetry by which to group atoms. Can be one of “wyckoffs”, “crystallographic_orbits”, or “equivalent_atoms”. See the documentation of spglib for a description of these quantities.
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
angle_tolerance (float, optional) – Symmetry-search tolerance in degrees. If the value is negative (default), an internally optimized routine is used to judge symmetry.
- Returns
groups (dict) – Dictionary where each key is a distinct site symmetry (e.g. a distinct orbit), and values are iterable of
Atom
which share the same site symmetry.- Raises
RuntimeError – If symmetry-determination has not succeeded.
ValueError – If by is not one of the supported values.
Notes
Note that crystals generated from the Protein Data Bank are often incomplete; in such cases the site symmetry information will be incorrect.
- ideal(symprec: float = 0.01) Crystal ¶
Returns a Crystal object with an idealized unit cell.
- Parameters
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
- Returns
ideal (Crystal) – Crystal with idealized cell.
:raises RuntimeError : If an ideal cell could not be found.:
Notes
Optional atomic properties (e.g magnetic moment) might be lost in the symmetrization.
- indexed_by(lattice: Union[Lattice, ndarray]) Crystal ¶
Return a crystal structure, indexed by another lattice/crystal structure.
- primitive(symprec: float = 0.01) Crystal ¶
Returns a Crystal object in the primitive unit cell.
- Parameters
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
- Returns
primitive (Crystal) – Crystal with primitive cell. Even if the crystal already has a primitive cell, a new crystal is returned.
:raises RuntimeError : If primitive cell could not be found.:
Notes
Optional atomic properties (e.g magnetic moment) might be lost in the reduction.
- reciprocal_symmetry_operations(symprec: float = 0.01) Iterable[ndarray] ¶
Get the symmetry operations that the reciprocal unit cell respects. These symmetry operations are expressed in reciprocal fractional coordinates.
- Parameters
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
- Returns
sym_ops (iterable of array_like, shapes (4,4)) – Iterable of affine transforms, where
m[:3,:3]
is the rotation part, whilem[:3,-1]
is the translation.
:raises RuntimeError : if symmetry-determination has not succeeded.:
See also
Crystal.symmetry_operations
symmetry operations in lattice basis
- supercell(n1: int, n2: int, n3: int) Supercell ¶
Create a supercell from this crystal, i.e. an atomic structure where the crystal unit cell is duplicated along lattice vectors.
- Parameters
n1, n2, n3 (int) – Repeat along the a1, a2, and a3 lattice vectors. For example,
(1, 1, 1)
represents the trivial supercell.- Returns
cell (Supercell) – Iterable of crystals.Atom objects following the supercell dimensions.
- symmetry(symprec: float = 0.01, angle_tolerance: float = -1.0) Dict[str, Any] ¶
Returns a dictionary containing space-group information. This information is computed from the crystal unit cell.
- Parameters
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
angle_tolerance (float, optional) – Symmetry-search tolerance in degrees. If the value is negative (default), an internally optimized routine is used to judge symmetry.
- Returns
info (dict) –
Dictionary of space-group information. The following keys are available:
'international_symbol'
: International Tables of Crystallography space-group symbol (short);'international_full'
: International Tables of Crystallography space-group full symbol;'hall_symbol'
: Hall symbol;'hm_symbol'
: Hermann-Mauguin symbol;'centering'
: Centering-type (“P”, “F”, etc.);'pointgroup'
: International Tables of Crystallography point-group;'international_number'
: International Tables of Crystallography space-group number (between 1 and 230);'hall_number'
: Hall number (between 1 and 531).
:raises RuntimeError : if symmetry-determination has not succeeded.:
Notes
Note that crystals generated from the Protein Data Bank are often incomplete; in such cases the space-group information will be incorrect.
- symmetry_operations(symprec: float = 0.01) Iterable[ndarray] ¶
Get the symmetry operations that the crystal unit cell respects. These symmetry operations are expressed in fractional coordinates.
- Parameters
symprec (float, optional) – Symmetry-search distance tolerance in Cartesian coordinates [Angstroms].
- Returns
sym_ops (iterable of array_like, shapes (4,4)) – Iterable of affine transforms, where
m[:3,:3]
is the rotation part, whilem[:3,-1]
is the translation.
:raises RuntimeError : if symmetry-determination has not succeeded.:
See also
Crystal.reciprocal_symmetry_operations
symmetry operations in reciprocal basis
- to_ase(**kwargs)¶
Convert a into an
ase.Atoms
object. Keyword arguments are passed toase.Atoms
constructor.Note that some information may be lost in the translation. However, we guarantee that reading a structure from a file, and then writing back to the same format is idempotent.
- Returns
atoms (ase.Atoms) – Group of atoms ready for ASE’s routines.
:raises ImportError : If ASE is not installed:
See also
Crystal.to_cif
write a structure to a .cif file.
Crystal.to_xyz
write a structure to a .xyz file.
- to_cif(filename: PathLike)¶
Convert this
Crystal
instance to a CIF file.Note that some information may be lost in the translation. However, we guarantee that reading a structure from a file, and then writing back to the same format is idempotent.
- Parameters
filename (path-like) – Path to a file. If the file already exists, it will be overwritten.
See also
Crystal.to_xyz
write a structure to a .xyz file.
Crystal.to_ase
convert a structure into an
ase.Atoms
object.
- to_poscar(filename: PathLike, **kwargs)¶
Convert this
Crystal
instance to a POSCAR file. Keyword arguments are passed towriters.write_poscar()
.Note that some information may be lost in the translation. However, we guarantee that reading a structure from a file, and then writing back to the same format is idempotent.
- Parameters
filename (path-like) – Path to a file. If the file already exists, it will be overwritten.
kwargs – Keyword arguments are passed to
writers.write_poscar()
.
- to_xyz(filename: PathLike)¶
Convert this
Crystal
instance to a XYZ file.Note that some information may be lost in the translation. However, we guarantee that reading a structure from a file, and then writing back to the same format is idempotent.
- Parameters
filename (path-like) – Path to a file. If the file already exists, it will be overwritten.
See also
Crystal.to_cif
write a structure to a .cif file.
Crystal.to_ase
convert a structure into an
ase.Atoms
object.
- property centering: CenteringType¶
Centering type of this crystals.
- property hall_number: int¶
Hall number (between 1 and 531).
- property hall_symbol: str¶
Hall symbol.
- property hm_symbol: str¶
Hermann-Mauguin symbol.
- property international_full: str¶
International Tables of Crystallography space-group full symbol.
- property international_number: int¶
International Tables of Crystallography space-group number (between 1 and 230).
- property international_symbol: str¶
International Tables of Crystallography space-group short symbol.
- property pointgroup: str¶
International Tables of Crystallography point-group.
Supercell¶
- class crystals.Supercell(dimensions: Tuple[int, int, int], **kwargs)¶
Bases:
Crystal
The
Supercell
class is a set-like container that represents a supercell of crystalline structures. It has all the same attributes as aCrystal
, with the addition ofSupercell.scaled_lattice_vectors
andSupercell.dimensions
.It is strongly recommended that you do not instantiate a
Supercell
by hand, but rather create aCrystal
object and use theCrystal.supercell()
method.To iterate over all atoms in the supercell, use this object as an iterable.
- property scaled_lattice_vectors: Tuple[ndarray, ndarray, ndarray]¶
3-tuple of lattice vectors scaled by the dimensions of the supercell. If the supercell is (2x2x2) unit cells, then the scaled lattice vectors are all 2x of the unit cell lattice vectors.
New in version 1.5.0.
Atom¶
To deal with atoms with coordinates, take a look at the Atom
class:
- class crystals.Atom(element: Union[str, int, Element], coords: Union[Sequence[Sequence[Sequence[Sequence[Sequence[Any]]]]], _SupportsArray[dtype], Sequence[_SupportsArray[dtype]], Sequence[Sequence[_SupportsArray[dtype]]], Sequence[Sequence[Sequence[_SupportsArray[dtype]]]], Sequence[Sequence[Sequence[Sequence[_SupportsArray[dtype]]]]], bool, int, float, complex, str, bytes, Sequence[Union[bool, int, float, complex, str, bytes]], Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]], Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]], Sequence[Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]]]], lattice: Optional[Union[Lattice, Sequence[Sequence[Sequence[Sequence[Sequence[Any]]]]], _SupportsArray[dtype], Sequence[_SupportsArray[dtype]], Sequence[Sequence[_SupportsArray[dtype]]], Sequence[Sequence[Sequence[_SupportsArray[dtype]]]], Sequence[Sequence[Sequence[Sequence[_SupportsArray[dtype]]]]], bool, int, float, complex, str, bytes, Sequence[Union[bool, int, float, complex, str, bytes]], Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]], Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]], Sequence[Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]]]]] = None, displacement: Optional[Union[Sequence[Sequence[Sequence[Sequence[Sequence[Any]]]]], _SupportsArray[dtype], Sequence[_SupportsArray[dtype]], Sequence[Sequence[_SupportsArray[dtype]]], Sequence[Sequence[Sequence[_SupportsArray[dtype]]]], Sequence[Sequence[Sequence[Sequence[_SupportsArray[dtype]]]]], bool, int, float, complex, str, bytes, Sequence[Union[bool, int, float, complex, str, bytes]], Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]], Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]], Sequence[Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]]]]] = None, magmom: Optional[float] = None, occupancy: float = 1.0, tag: Optional[Any] = None, electronic_structure: Optional[ElectronicStructure] = None, **kwargs)¶
Bases:
Element
Container object for atomic data.
- Parameters
element (str or int) – Chemical element symbol or atomic number.
coords (array-like, shape (3,)) – Coordinates of the atom in fractional form.
lattice (Lattice or array-like, shape (3,3)) – Lattice on which the atom is positioned.
displacement (array-like or None, optional) – Atomic maximum displacement [Angs].
magmom (float, optional) – Magnetic moment. If None (default), the ground-state magnetic moment is used.
occupancy (float, optional) – Fractional occupancy. If None (default), occupancy is set to 1.0.
tag (int or None, optional) – Tag an atom with a unique identifier. Useful to keep track of atom order, for example in PWSCF output files. This is mostly for internal use.
electronic_structure (ElectronicStructure or None, optional) – Electronic orbital structure for this atom. If None (default), the ground state for this element will be used.
- classmethod from_ase(atom) Atom ¶
Returns an Atom instance from an ASE atom
- Parameters
atom (ase.Atom)
- transform(*matrices: Union[Sequence[Sequence[Sequence[Sequence[Sequence[Any]]]]], _SupportsArray[dtype], Sequence[_SupportsArray[dtype]], Sequence[Sequence[_SupportsArray[dtype]]], Sequence[Sequence[Sequence[_SupportsArray[dtype]]]], Sequence[Sequence[Sequence[Sequence[_SupportsArray[dtype]]]]], bool, int, float, complex, str, bytes, Sequence[Union[bool, int, float, complex, str, bytes]], Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]], Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]], Sequence[Sequence[Sequence[Sequence[Union[bool, int, float, complex, str, bytes]]]]]]) Atom ¶
Return an Atom with fractional coordinates transformed according to symmetry operators.
- Parameters
matrices (ndarrays, shape {(3,3), (4,4)}) – Transformation matrices.
- Returns
atm (Atom) – Transformed atom. The original atom is left unchanged.
- property coords_cartesian: ndarray¶
Real-space position of the atom on the lattice, in Angstroms.
- Returns
pos (~numpy.ndarray, shape (3,)) – Atomic position
:raises RuntimeError : if this atom is not place on a lattice:
Element¶
To access atomic data only, take a look at the Element
class:
- class crystals.Element(element: Union[str, int, Element], *args, **kwargs)¶
Bases:
object
Class representing an abtract chemical element, but no particular atom. This class gives access to elemental properties, like atomic number, atomic mass, full element name, etc.
- Parameters
element (str, int, or
Element
) – Elemental symbol (e.g. “He”), element name (e.g. “Helium”), atomic number, or another Element instance.
:raises ValueError : if the element is not valid.:
- property atomic_number: int¶
Atomic number
- property element_full: str¶
Full element name, e.g. “Hydrogen”
- property magnetic_moment_ground: float¶
Ground state magnetic moment.
- property mass: float¶
Atomic mass [u]
- property name: str¶
Full element name, e.g. “Hydrogen”
- property symbol: str¶
Elemental symbol, e.g. “He”
Electronic Orbitals¶
If you need to specify the orbital structure of atoms by hand, you will need to know about ElectronicStructure
class:
- class crystals.ElectronicStructure(shells: Dict[Union[str, Orbital], int])¶
Description of the atomic orbital structure.
- Parameters
shells (dict[Orbital,int]) – Dictionary containing the number of electrons in each Orbital, e.g. {“1s”: 2}.
:raises ValueError : if the electronic structure is not representable:
Examples
Electronic structures can be specified by hand:
>>> ElectronicStructure({"1s":2, "2s":2, "2p":2}) < ElectronicStructure: 1s²2s²2p² >
A shortcut exists for atomic ground states:
>>> ElectronicStructure.ground_state("Ti") < ElectronicStructure: 1s²2s²2p⁶3s²3p⁶4s²3d² >
Notes
Shells are allowed to be filled our of order deliberately, given that unusual electronic structures can arise from ultrafast photoexcitation.
Electronic orbitals are described by the Orbital
class:
- class crystals.Orbital(value)¶
Enumeration of electronic orbitals, used to described atomic orbital structure.
We note that Orbital instances are ordered according to the Madelung rule.
Examples
>>> Orbital("1s") <Orbital.one_s: '1s'> >>> Orbital.four_p == Orbital("4p") True
Enumerations¶
To represent lattice systems and centering types, the following enumerations are used:
- class crystals.LatticeSystem(value)¶
Lattice system enumeration.
Equivalent to Crystal family except that the hexagonal crystal family is split between the rhombohedral system and hexagonal system.
- class crystals.CenteringType(value)¶
Enumeration of possible centering types. Together with the lattice system, these centering types defined all 14 Bravais lattices in 3D.
The possible centering types are:
'P'
: Primitive'I'
: Body-centered'F'
: Face-centered'C'
: Base-centered'R'
: Rhombohedral in hexagonal setting.