dfttools

dfttools.formatters

This submodule contains routines presenting data (unit cell) in various text formats.

dfttools.formatters.openmx_input(cell, populations, l=None, r=None, tolerance=1e-10, indent=4)[source]

Generates OpenMX minimal input file with atomic structure.

Args:

cell (UnitCell): input unit cell;

populations (dict): a dict with initial electronic populations data;

Kwargs:

l (UnitCell): left lead;

r (UnitCell): right lead;

tolerance (float): tolerance for checking whether left-center-right unit cells can be stacked;

indent (int): size of indent;

Returns:

String with OpenMX input file formatted data.
dfttools.formatters.pyscf_cell(cell, **kwargs)[source]

Constructs a unit cell object in pyscf.

Args:

cell (UnitCell): a unit cell object to convert from;

Kwargs are passed to pyscf.pbc.gto.M.

Returns:

A Pyscf Cell object.
dfttools.formatters.qe_input(cell=None, relax_triggers=0, parameters={}, inline_parameters={}, pseudopotentials={}, indent=4)[source]

Generates Quantum Espresso input file.

Kwargs:

cell (UnitCell): a unit cell with atomic coordinates;

relax_triggers (array,int): array with triggers for relaxation written asadditional columns in the input file;

parameters (dict): parameters for the input file;

inline_parameters (dict): a dict of inline parameters such as crystal_b, etc;

pseudopotentials (dict): a dict of pseudopotential file names;

indent (int): size of indent;

Returns:

String contating Quantum Espresso input file contents.
dfttools.formatters.siesta_input(cell, indent=4)[source]

Generates Siesta minimal input file with atomic structure.

Args:

cell (UnitCell): input unit cell;

Kwargs:

indent (int): size of indent;

Returns:

String with Siesta input file contents.
dfttools.formatters.xsf_grid(grid, cell, npl=6)[source]

Generates an xcrysden file with the data on the grid.

Args:

grid (Grid): data on the grid;

cell (UnitCell): structural data;

Kwargs:

npl (int): numbers per line in the grid section;

Returns:

A string contating XSF-formatted data.
dfttools.formatters.xsf_structure(*cells)[source]

Generates an xcrysden file with the structure.

Args:

cells (list): unit cells with atomic coordinates;

Returns:

A string contating XSF-formatted data.

dfttools.presentation

This submodule contains data visualization routines.

dfttools.presentation.matplotlib2svgwrite(fig, svg, insert, size, **kwargs)[source]

Saves a matplotlib image to an existing svgwrite object.

Args:

fig (matplotlib.figure.Figure): a figure to save;

svg (svgwrite.Drawing): an svg drawing to save to;

insert (tuple): a tuple of ints defining destination to insert a drawing;

size (tuple): size of the inserted image;

Kwargs:

The kwargs are passed to fig.savefig used to print the plot.
dfttools.presentation.matplotlib_bands(cell, axes, show_fermi=True, energy_range=None, energy_units='eV', energy_units_name=None, coordinate_units=None, coordinate_units_name=None, threshold=0.01, weights=None, weights_color=None, weights_size=None, optimize_visible=False, edge_names=[], mark_points=None, project=None, **kwargs)[source]

Plots basic band structure using pyplot.

Args:

cell (UnitCell): cell with the band structure;

axes (matplotlib.axes.Axes): axes to plot on;

Kwargs:

show_fermi (bool): shows the Fermi level if specified;

energy_range (array): 2 floats defining plot energy range. The units of energy are defined by the units keyword;

energy_units (str, float): either a field from numericalunits package or a float with energy units;

energy_units_name (str): a string used for the units. Used only if the energy_units keyword is a float;

coordinate_units (str, float): either a field from numericalunits package or a float with coordinate units or None;

coordinate_units_name (str): a string used for the coordinate units. Used only if the coordinate_units keyword is a float;

threshold (float): threshold for determining edges of k point path;

weights, weights_color (array): a 2D array with weights on the band structure which will be converted to color according to current colormap;

weights_size (array): a 2D array with weights on the band structure which will be converted to line thickness;

optimize_visible (bool): draw only visible lines;

edge_names (list): the edges names to appear on the band structure;

mark_points (list): marks specific points on the band structure, the first number in each list element is interpreted as k-point while the second number is band number;

project (array): projects k-points along specified direction instead of unfolding the entire bands path. If coordinate_units specified the direction is expressed in the unit cell vectors, otherwise cartesian basis is used;

The rest of kwargs are passed to matplotlib.collections.LineCollection.

Returns:

A plotted LineCollection.
dfttools.presentation.matplotlib_bands_density(cell, axes, energies, show_fermi=True, energy_range=None, units='eV', units_name=None, weights=None, on_top_of=None, use_fill=False, orientation='landscape', gaussian_spread=None, method='default', **kwargs)[source]

Plots density of bands (density of states).

The cell values are considered to be band energies.

Args:

cell (Grid,UnitCell): a unit cell with the band structure, possibly on the grid;

axes (matplotlib.axes.Axes): axes to plot on;

energies (int,array): energies to calculate density at. The integer value has the meaning of number of points to cover the range energy_range. Otherwise the units of energy are defined by the units keyword;

Kwargs:

show_fermi (bool): shows the Fermi level if specified;

energy_range (array): 2 floats defining plot energy range. The units of energy are defined by the units keyword;

units (str, float): either a field from numericalunits package or a float with energy units;

units_name (str): a string used for the units. Used only if the units keyword is a float;

weights (array): a 2D array with weights on the band structure;

on_top_of (array): a 2D array with weights on the band structure to plot on top of;

use_fill (bool): fill the area below plot;

orientation (str): either ‘portrait’ or ‘landscape’ - orientation of the plot;

gaussian_spread (float): the gaussian spread for the density of states. This value is used only if the provided cell is not a Grid;

method (bool): method to calculate density: ‘default’, ‘gaussian’ or ‘optimal’;

The rest of kwargs are passed to pyplot plotting functions.

Returns:

A plotted Line2D or a PolyCollection, depending on use_fill.
dfttools.presentation.matplotlib_scalar(grid, axes, origin, plane, units='angstrom', units_name=None, show_cell=False, normalize=True, ppu=None, isolines=None, window=None, margins=0.1, scale_bar=None, scale_bar_location=1, **kwargs)[source]

Plots scalar values on the grid using imshow.

Args:

grid (Grid): a 3D grid to be plotted;

axes (matplotlib.axes.Axes): axes to plot on;

origin (array): origin of the 2D slice to be plotted in the units of grid;

plane (str, int): the plotting plane: either ‘x’,’y’ or ‘z’ or a correspondint int.

Kwargs:

units (str, float): either a field from numericalunits package or a float with energy units;

units_name (str): a string used for the units. Used only if the units keyword is a float;

show_cell (bool): if True then projected unit cell boundaries are shown on the final image;

normalize (bool): normalize data before plotting such that the minimum is set at zero and the maximum is equal to one;

ppu (float): points per unit for the raster image;

isolines (array): plot isolines at the specified levels;

window (array): 4 values representing a window to plot the data: minimum and maximum ‘x’ coordinate and minimum and maximum ‘y’ coordinate;

margins (float): adds margins to the grid where the data is interpolated;

scale_bar (int): adds a scal bar to the image at the specified location;

scale_bar_location (int): location of the scale bar;

The rest of kwargs are passed to pyplot.imshow or pyplot.contour.

Returns:

A matplotlib.image.AxesImage plotted.
dfttools.presentation.svgwrite_unit_cell(cell, svg, camera=None, camera_top=None, insert=(0, 0), size=(600, 600), circle_size=0.4, circle_opacity=None, margin=6, show_cell=False, show_atoms=True, show_bonds=True, show_legend=True, show_numbers=False, fadeout_strength=0.8, bg=(255, 255, 255), bond_ratio=1, hook_atomic_color=None, coordinates='right', invisible=None, title=None)[source]

Creates an svg drawing of a unit cell.

Args:

cell (UnitCell): the cell to be visualized;

svg (str, svgwrite.Drawing): either file name to save the drawing to or an svgwrite.Drawing object to draw with.

Kwargs:

camera (str, array): the direction of a camera: either ‘x’,’y’ or ‘z’ or an arbitrary 3D vector;

camera_top (array): a vector pointing up;

insert (array): a top-left corner of the drawing;

size (array): size of the bounding box;

circle_size (float): size of the circles representing atoms, arbitrary units;

circle_opacity (float,array): opacity of circles;

margin (float): size of the margin in all directions;

show_cell (bool, str): if True draws the unit cell edges projected, if ‘invisible’ the unit cell is invisible;

show_atoms (bool): if True draws atoms;

show_bonds (bool): if True draws bonds;

show_legend (bool): if True draws legend;

show_numbers (bool): if True shows numbers corresponding to the atomic order in the unit cell;

fadeout_strength (float): amount of fadeout applied to more distant atoms;

bg (array): an integer array defining background color;

bond_ratio (float): scale factor to determine whether the bond is rendered;

coordinates (str): the coordinate system, either ‘left’ or ‘right’;

hook_atomic_color (function): a function accepting integer (atom ID) and a 3-element list (suggested RGB color) and returning a new color of the atom;

invisible (str,array): make specified atoms invisible. If ‘auto’ specified, creates a supercell and makes all cell replica invisible. The bonds of invisible atoms will still be present on the final image;

title (str): a title to the drawing presented in the top left corner;

Returns:

An `svgwrite.Drawing object. The object is saved if it was created inside this method.

dfttools.simple

This submodule contains commonly used shortcuts to parse the data.

dfttools.simple.get_all_parsers(*modules)[source]

Retrieves all parsers.

Kwargs:

modules (list): a list of names of parsers submodules to search at.

Returns:

A list of parsing classes.
dfttools.simple.guess_parser(f)[source]

Guesses parsers for a given data.

Args:

f (file): a file to parse.

Returns:

A list of parser candidates.
dfttools.simple.parse(f, tag, *args)[source]

Identifies and parses data.

Args:

f (file): a file to parse;

tag (str): the data tag, such as unit-cell or band-structure;

Returns:

The parsed data.
dfttools.simple.tag_method(*tags, **kwargs)[source]

A generic decorator tagging some method.

Kwargs:

take_file (bool): set to True and the File object will be passed to this method.

dfttools.types

This submodule contains key types for handling coordinate-dependent data: UnitCell, Grid and Basis.

class dfttools.types.Basis(vectors, kind='default', meta=None, units=None)[source]

A class describing a set of vectors representing a basis.

Args:

vectors (array): a 2D or a 1D array of floats representing vectors of the basis set.

Kwargs:

kind (str): a shortcut keyword for several most common basis sets:

  • ‘default’: expects vectors to be a 2D array with basis vectors in cartesian coordinates;
  • ‘orthorombic’: expects vectors to be a 1D array with dimensions of an orthorombic basis set;
  • ‘triclinic’: expects vectors to be a 1D array with 3 lengths of edges and 3 cosines of face angles.

units (str,float): optional units for the Basis. The units are stored in self.meta[‘units’] and are used only during save/load process. The string value has to correspond to one of the values in numericalunits package.

meta (dict): a metadata for this Basis.

copy()[source]

Calculates a copy.

Returns:

A deep copy of self.
edges()[source]

Computes pairs of cartesian coordinates of all edges of the basis cell.

Returns:

A list of pairs with cartesian coordinates of vertices forming edges.
faces()[source]

Computes faces and returns corresponding cartesian coordinates.

Returns:

A list of lists of coordinates defining face polygon coordinates.
static from_json(j)[source]

Restores a Basis from JSON data.

Args:

j (dict): JSON data.

Returns:

A Basis object.
generate_path(points, n, anchor=True)[source]

Generates a path given key points and the total number of points on the path.

Args:

points (array): key points of the path expressed in this basis;

n (int): the total number of points on the path.

Kwargs:

anchor (bool): force the specified points to be present in the final path. If True alters slightly the total point number;

Returns:

Path coordinates expressed in this basis.
reciprocal()[source]

Computes a reciprocal basis.

Returns:

A reciprocal basis.

Note

The multiplier is not present.

reorder_vectors(*args, **kwargs)[source]

Reorders vectors.

Args:

new (array): new mapping of vectors.

Example:

>>> basis.reorder_vectors(0, 1, 2) # does nothing
>>> basis.reorder_vectors(1, 0, 2) # swaps first and second vectors.
repeated(*args, **kwargs)[source]

Produces a new basis from a given by repeating it in all directions.

Args:

times (array): array of ints specifying how much the basis should be repeated along each of the vectors.
rotated(axis, angle, units='rad')[source]

Rotates this basis.

Args:

axis (array): axis to rotate around;

angle (float): angle to rotate;

Kwargs:

units (str): units of the angle: ‘rad’, ‘deg’ or ‘frac’.

Returns:

A rotated copy of this basis.
stack(*args, **kwargs)[source]

Stacks several basises along one of the vectors.

Args:

basises (list): basises to stack. Corresponding vectors of all basises being stacked should match.

Kwargs:

vector (str,int): a vector along which to stack, either ‘x’, ‘y’, ‘z’ or an int specifying the vector;

tolerance (float): a largest possible error in input basises’ vectors;

restrict_collinear (bool): if True will raise an exception if the vectors to stack along are not collinear

Raises:

ArgumentError: in the case of vector mismatch.

Returns:

A larger basis containing all argument cell stacked.
to_json()[source]

Prepares a JSON-compatible object representing this Basis.

Returns:

A JSON-compatible dict.
transform_from(basis, coordinates)[source]

Transforms coordinates from another basis set.

Args:

basis (Basis): a basis to transform from.

coordinates (array): an array of coordinates to be transformed.

Returns:

An array with transformed coordinates.
transform_from_cartesian(coordinates)[source]

Transforms coordinates from cartesian.

Args:

coordinates (array): an array of coordinates to be transformed.

Returns:

An array with transformed coordinates.
transform_to(basis, coordinates)[source]

Transforms coordinates to another basis set.

Args:

basis (Basis): a new basis to transform to.

coordinates (array): an array of coordinates to be transformed.

Returns:

An array with transformed coordinates.
transform_to_cartesian(coordinates)[source]

Transforms coordinates to cartesian.

Args:

coordinates (array): an array of coordinates to be transformed.

Returns:

An array with transformed coordinates.
units_aware()[source]

Checks if units for this Basis are defined.

Returns:

True if units were defined.
vertices()[source]

Computes cartesian coordinates of all vertices of the basis cell.

Returns:

Cartesian coordinates of all vertices.
volume()[source]

Computes the volume of a triclinic cell represented by the basis.

Returns:

Volume of the cell in m^3.
class dfttools.types.Grid(basis, coordinates, values, units=None)[source]

A class describing a data on a grid in a periodic environment.

Args:

basis (Basis): a crystal basis.

coordinates (array): a list of arrays of coordinates specifying grid.

values (array): a multidimensional array with data on the grid.

add(*args, **kwargs)[source]

Adds grid points from another grids to this one.

Args:

grids (arguments): grids to be merged with.

Returns:

A new grid with merged data.
apply(*args, **kwargs)[source]

Applies selection to this grid.

Inverse of Grid.discard.

Args:

selection (array): seleted grid points.

Example:

>>> selection = grid.select((0,0,0,0.5,1,1)) # Selects species in the 'left' part of the grid.
>>> grid.apply(selection) # Applies selection. Species outside the 'left' part are discarded.
as_unitCell()[source]

Converts this cell to a UnitCell.

Returns:

A new UnitCell.
cartesian()[source]

Computes cartesian coordinates.

Returns:

A multidimensional numpy array with cartesian coordinates at each grid point.
static combine_arrays(arrays)[source]

Transforms input 1D arrays of coordinates into (N+1)D mesh array where first N dimensions correspond to a particular grid point and the last dimension specifies all coordinates of this grid point.

Args:

arrays (list): a list of 1D arrays;

Returns:

A meshgrid array with coordinates.
copy()[source]

Calculates a copy.

Returns:

A copy of self.
cut(*args, **kwargs)[source]

Selects a piece of the grid and returns it as a smaller basis.

Kwargs:

piece (array): fraction of the grid to be selected. The order of coordinates in piece is x_from, y_from, ..., z_from, x_to, y_to, ..., z_to.

Returns:

A smaller grid selected.
discard(*args, **kwargs)[source]

Removes specified points from this grid.

Inverse of Grid.apply.

Args:

selection (array): points to remove.

Example:

>>> selection = grid.select((0,0,0,0.5,1,1)) # Selects points in the 'left' part of the grid.
>>> grid.discard(selection) # Discards selection. Points inside the 'left' part are removed.
explicit_coordinates()[source]

Creates an (N+1)D array with explicit coordinates at each grid point.

Returns:

An (N+1)D array with coordinates.
interpolate_to_array(points, driver=None, periodic=True, **kwargs)[source]

Interpolates values at specified points and returns an array of interpolated values. By default uses scipy.interpolate.interpn.

Args:

points (array): points to interpolate at.

Kwargs:

driver (func): interpolation driver;

periodic (bool): employs periodicity of a unit cell;

The rest of keyword arguments are passed to the driver.

Returns:

An array with values of corresponding shape.
interpolate_to_cell(points, **kwargs)[source]

Interpolates values at specified points and returns a unit cell. By default uses scipy.interpolate.interpn.

Args:

points (array): points to interpolate at.

Kwargs are passed to self.interpolate_to_array.

Returns:

A unit cell interpolated values.
interpolate_to_grid(points, **kwargs)[source]

Interpolates values at specified points and returns a grid. By default uses scipy.interpolate.interpn.

Args:

points (array): points to interpolate at.

Kwargs are passed to self.interpolate_to_array.

Returns:

A grid with interpolated values.
interpolate_to_path(points, n, anchor=True, **kwargs)[source]

Interpolates values to a path: performs path generation using Basis.generate_path and interpolates values on it.

Args:

points (array): key points of the path expressed in lattice coordinates;

n (int): the total number of points on the path.

Kwargs:

anchor (bool): force the specified points to be present in the final path. If True alters slightly the total point number;

kwargs: keywords to the ‘interpolate_to_cell’ routine.

Returns:

A unit cell interpolated values.
isolated(*args, **kwargs)[source]

Generates an isolated representation of this grid.

Symmetrically adds vacuum along all basis vectors such that resulting grid basis vectors are parallel to the initial ones.

Args:

gaps (array): size of the vacuum layer in each direction either in cartesian or in crystal units.

Kwargs:

units (str): units of the vacuum size, ‘cartesian’ or ‘crystal’

Returns:

A new isolated grid.
normalized()[source]

Moves all grid points respecting periodicity so that each coordinate becomes in the unit range 0<=x<1 in the cell basis. Sorts the data.

Returns:

A new grid with normalized data.
reorder_vectors(*args, **kwargs)[source]

Reorders vectors. Does not change output of Grid.cartesian().

Args:

new (array): new mapping of vectors.

Example:

>>> grid.reorder_vectors(0, 1, 2) # does nothing
>>> grid.reorder_vectors(1, 0, 2) # swaps first and second vectors.

Note

A call to this method does not modify the output of self.cartesian().

select(*args, **kwargs)[source]

Selects a piece of this grid.

Args:

piece (array): fraction of the grid to be selected, see examples. The order of coordinates in piece is x_from, y_from, ..., z_from, x_to, y_to, ..., z_to.

Returns:

A list of numpy arrays with bools defining whether particular grid point is selected or not.

Example:

>>> grid.select((0,0,0,1,1,1)) # select all grid points with coordinates within (0,1) range
>>> grid.select(0,0,0,1,1,1) # a simplified version of above
>>> grid.select(0,0,0,0.5,1,1) # select the 'left' part
>>> grid.select(0.5,0,0,1,1,1) # select the 'right' part
size()[source]

Retrieves the total size of points on the grid.

Returns:

Number of species in cell as an integer.
stack(*args, **kwargs)[source]

Stacks several grids along one of the vectors.

Args:

grids (list): grids to stack. Corresponding vectors of all grids being stacked should match.

Kwargs:

vector (str,int): a vector along which to stack, either ‘x’, ‘y’, ‘z’ or an int specifying the vector.

The rest of kwargs are redirected to Basis.stack.

Raises:

ArgumentError: in the case of vector mismatch.

Returns:

A bigger grid containing all argument grids stacked.
tetrahedron_density(points, resolved=False, weights=None)[source]

Convolves data to calculate density (of states). Uses the tetrahedron method from PRB 49, 16223 by E. Blochl et al. Works only in a 3D space.

Args:

points (array): values to calculate density at.

Kwargs:

resolved (bool): if True returns a spacially and index resolved density. The dimensions of the returned array are self.values.shape + points.shape.

weights (array): if specified and resolved is False convolves result with the specified weights.

Returns:

A numpy array containing density: 1D if resolved == False or a corresponding Grid if resolved == True.
static uniform(size, endpoint=False)[source]

Transform positive integers size into a meshgrid array representing a grid where grid points span uniformly zero to one intervals.

Args:

size (array): an array with positive integers;

Kwargs:

endpoint (bool): indicates whether to include x=1 into grids.

Returns:

A meshgrid array with coordinates.
class dfttools.types.UnitCell(basis, coordinates, values, c_basis=None, units=None)[source]

A class describing a crystal unit cell in a periodic environment.

Args:

basis (Basis): a crystal basis.

coordinates (array): a 2D array of coordinates of atoms (or any other instances)

values (array): an array of atoms (or any other instances) with the leading dimenstion being the same as the one of coordinates array.

Kwargs:

c_basis (str,Basis): a Basis for input coordinates or ‘cartesian’ if coordinates are passed in the cartesian basis.
add(*args, **kwargs)[source]

Adds species from another unit cells to this one.

Args:

cells (arguments): unit cells to be merged with.

Returns:

A new unit cell with merged data.
angles(*args, **kwargs)[source]

Computes angles between cell specimens.

Args:

ids (array): a set of specimen IDs to compute angles between. Several shapes are accepted:

  • nx3 array: computes n cosines of angles [n,0]-[n,1]-[n,2];
  • 1D array of length n: computes n-2 cosines of angles [n-1]-[n]-[n+1];

Returns:

A numpy array containing cosines of angles specified.

Example:

Following are the valid calls:

>>> cell.angles((0,1,2)) # angle between vectors connecting {second and first} and {second and third} species
>>> cell.angles(0,1,2) # a simplified version of above
>>> cell.angles(0,1,3,2) # two angles along path: 0-1-3 and 1-3-2
>>> cell.angles(tuple(0,1,3,2)) # same as above
>>> cell.angles((0,1,3),(1,3,2)) # same as above
apply(*args, **kwargs)[source]

Applies selection to this cell.

Inverse of UnitCell.discard.

Args:

selection (array): seleted species.

Example:

>>> selection = cell.select((0,0,0,0.5,1,1)) # Selects species in the 'left' part of the unit cell.
>>> cell.apply(selection) # Applies selection. Species outside the 'left' part are discarded.
as_grid(fill=nan)[source]

Converts this unit cell to grid.

Kwargs:

fill: default value to fill with

Returns:

A new grid with data from initial cell.
cartesian()[source]

Computes cartesian coordinates.

Returns:

A numpy array with cartesian coordinates
copy()[source]

Calculates a copy.

Returns:

A copy of self.
cut(*args, **kwargs)[source]

Selects a piece of this unit cell and returns it as a smaller unit cell.

Kwargs:

piece (array): fraction of the cell to be selected. The order of coordinates in piece is x_from, y_from, ..., z_from, x_to, y_to, ..., z_to.

Returns:

A smaller unit cell selected.
discard(*args, **kwargs)[source]

Removes specified species from cell.

Inverse of Cell.apply.

Args:

selection (array): species to remove.

Example:

>>> selection = cell.select((0,0,0,0.5,1,1)) # Selects species in the 'left' part of the unit cell.
>>> cell.discard(selection) # Discards selection. Species inside the 'left' part are removed.
distances(*args, **kwargs)[source]

Computes distances between species and specified points.

Args:

ids (array): a list of specimen IDs to compute distances between. Several shapes are accepted:

  • empty: returns a 2D matrix of all possible distances
  • nx2 array of ints: returns n distances between each pair of [n,0]-[n,1] species;
  • 1D array of ints of length n: returns n-1 distances between each pair of [n-1]-[n] species;

Returns:

A numpy array containing list of distances.
static from_json(j)[source]

Restores a UnitCell from JSON data.

Args:

j (dict): JSON data.

Returns:

A UnitCell object.
interpolate(*args, **kwargs)[source]

Interpolates values at specified points. By default uses scipy.interpolate.griddata.

Args:

points (array): points to interpolate at.

Kwargs:

driver (func): interpolation driver.

periodic (bool): employs periodicity of a unit cell.

kwargs: keywords to the driver.

Returns:

A new unit cell with interpolated data.
isolated(*args, **kwargs)[source]

Generates an isolated representation of this cell.

Symmetrically adds vacuum along all unit cell vectors such that resulting unit cell vectors are parallel to the initial ones.

Args:

gaps (array): size of the vacuum layer in each direction either in cartesian or in crystal units.

Kwargs:

units (str): units of the vacuum size, ‘cartesian’ or ‘crystal’

Returns:

A new unit cell with spacially isolated species.
isolated2(gap)[source]

Generates an isolated representation of this cell.

The resulting cell is rectangular and contains space gaps of at least “gap” size.

Args:

gap (float): size of the space gap along all self.vectors.

Returns:

A new unit cell with spacially isolated species.
normalized(sort=None)[source]

Moves all species respecting periodicity so that each coordinate becomes in the unit range 0<=x<1 in the cell basis. Sorts the data if sort provided.

Kwargs:

sort: coordinates to sort with: either ‘x’, ‘y’, ‘z’ or 0,1,2.

Returns:

A new grid with normalized data.
packed()[source]

Moves all species as close to the origin as it is possible. Does not perform translation.

Returns:

A new unit cell with packed coordinates.
reorder_vectors(*args, **kwargs)[source]

Reorders vectors.

Args:

new (array): new mapping of vectors.

Example:

>>> cell.reorder_vectors(0, 1, 2) # does nothing
>>> cell.reorder_vectors(1, 0, 2) # swaps first and second vectors.

Note

A call to this method does not modify the output of self.cartesian().

select(*args, **kwargs)[source]

Selects a piece of this cell.

Args:

piece (array): fraction of the cell to be selected, see examples. The order of coordinates in piece is x_from, y_from, ..., z_from, x_to, y_to, ..., z_to.

Returns:

A numpy array with bools defining whether particular specimen is selected or not.

Example:

>>> cell.select((0,0,0,1,1,1)) # select all species with coordinates within (0,1) range
>>> cell.select(0,0,0,1,1,1) # a simplified version of above
>>> cell.select(0,0,0,0.5,1,1) # select the 'left' part
>>> cell.select(0.5,0,0,1,1,1) # select the 'right' part
size()[source]

Retrieves the number of points or species in this unit cell.

Returns:

Number of points or species in cell.
species()[source]

Collects number of species of each kind in this cell.

Particularly useful for counting the number of atoms.

Returns:

A dictionary containing species as keys and number of atoms as values.
stack(*args, **kwargs)[source]

Stacks several cells along one of the vectors.

Args:

cells (list): cells to stack. Corresponding vectors of all cells being stacked should match.

Kwargs:

vector (str,int): a vector along which to stack, either ‘x’, ‘y’, ‘z’ or an int specifying the vector.

The rest of kwargs are redirected to Basis.stack.

Raises:

ArgumentError: in the case of vector mismatch.

Returns:

A bigger unit cell containing all argument cell stacked.
supercell(*args, **kwargs)[source]

Produces a supercell from a given unit cell.

Args:

vec (array): the supercell vectors in units of current unit cell vectors

Returns:

A new supercell.
dfttools.types.diamond_basis(a)[source]

Creates a diamond basis with a given lattice constant.

Args:

a (float): the lattice constant;

Returns:

A diamond Basis.