vasppy package¶
Submodules¶
vasppy.atom module¶
vasppy.calculation module¶
-
class
vasppy.calculation.
Calculation
(title, energy, stoichiometry)[source]¶ Bases:
object
class describing a single VASP calculation
-
vasppy.calculation.
delta_E
(reactants, products, check_balance=True)[source]¶ Calculate the change in energy for reactants –> products.
Parameters: - reactants (list(vasppy.Calculation) – A list of vasppy.Calculation objects. The initial state.
- products (list(vasppy.Calculation) – A list of vasppy.Calculation objects. The final state.
- (bool (check_balance) – optional): Check that the reaction stoichiometry is balanced. Default: True.
Returns: (float) The change in energy.
-
vasppy.calculation.
delta_stoichiometry
(reactants, products)[source]¶ Calculate the change in stoichiometry for reactants –> products.
Parameters: - reactants (list(vasppy.Calculation) – A list of vasppy.Calculation objects. The initial state.
- products (list(vasppy.Calculation) – A list of vasppy.Calculation objects. The final state.
Returns: The change in stoichiometry.
Return type: (Counter)
-
vasppy.calculation.
energy_string_to_float
(string)[source]¶ Convert a string of a calculation energy, e.g. ‘-1.2345 eV’ to a float.
Parameters: string (str) – The string to convert. - Return
- (float)
-
vasppy.calculation.
import_calculations_from_file
(filename)[source]¶ Construct a list of Calculation objects by reading a YAML file. Each YAML document should include ‘title’, ‘stoichiometry’, and ‘energy’ fields. e.g.:
title: my calculation stoichiometry: - A: 1 - B: 2 energy: -0.1234 eV
Separate calculations should be distinct YAML documents, separated by —
Parameters: filename (str) – Name of the YAML file to read. Returns: A dictionary of Calculation objects. For each Calculation object, the ‘title’ field from the YAML input is used as the dictionary key. Return type: (dict(vasppy.Calculation))
vasppy.cell module¶
-
class
vasppy.cell.
Cell
(matrix)[source]¶ Bases:
object
-
angles
()[source]¶ The cell angles (in degrees).
Parameters: None – Returns: The cell angles. Return type: (list(alpha,beta,gamma))
-
cartesian_to_fractional_coordinates
(coordinates)[source]¶ Convert a set of Cartesian coordinates to fractional coordinates in the cell.
Parameters: coordinates (np.array(dim(N,3))) – The set of Cartesian coordinates. Returns: The corresponding set of fractional coordinates. Return type: (np.array(dim(N,3)))
-
dr
(r1, r2, cutoff=None)[source]¶ Calculate the distance between two fractional coordinates in the cell.
Parameters: - r1 (np.array) – fractional coordinates for position 1.
- r2 (np.array) – fractional coordinates for position 2.
- (optional (cutoff) – Bool): If set, returns None for distances greater than the cutoff. Default None (unset).
Returns: the distance between r1 and r2.
Return type: (float)
-
fractional_to_cartesian_coordinates
(coordinates)[source]¶ Convert a set of fractional coordinates in the cell to Cartesian coordinates.
Parameters: coordinates (np.array(dim(N,3))) – The set of fractional coordinates. Returns: The corresponding set of Cartesian coordinates. Return type: (np.array(dim(N,3)))
-
inside_cell
(r)[source]¶ Given a fractional-coordinate, if this lies outside the cell return the equivalent point inside the cell.
Parameters: r (np.array) – Fractional coordinates of a point (this may be outside the cell boundaries). Returns: Fractional coordinates of an equivalent point, inside the cell boundaries. Return type: (np.array)
-
lengths
()[source]¶ The cell lengths.
Parameters: None – Returns: The cell lengths. Return type: (np.array(a,b,c))
-
minimum_image
(r1, r2)[source]¶ Find the minimum image vector from point r1 to point r2.
Parameters: - r1 (np.array) – fractional coordinates of point r1.
- r2 (np.array) – fractional coordinates of point r2.
Returns: the fractional coordinate vector from r1 to the nearest image of r2.
Return type: (np.array)
-
minimum_image_dr
(r1, r2, cutoff=None)[source]¶ Calculate the shortest distance between two points in the cell, accounting for periodic boundary conditions.
Parameters: - r1 (np.array) – fractional coordinates of point r1.
- r2 (np.array) – fractional coordinates of point r2.
- ( (cutoff) – obj: float, optional): if set, return zero if the minimum distance is greater than cutoff. Defaults to None.
Returns: The distance between r1 and r2.
Return type: (float)
-
nearest_image
(origin, point)[source]¶ Find the fractional_coordinates of the nearest periodic image to a point of origin.
Parameters: - origin (np.array) – fractional coordinates of the point of origin.
- point (np.array) – fractional coordinates of the other point.
Returns: the fractional coordinates of the nearest image of point to origin.
Return type: (np.array)
-
-
vasppy.cell.
angle
(x, y)[source]¶ Calculate the angle between two vectors, in degrees.
Parameters: - x (np.array) – one vector.
- y (np.array) – the other vector.
Returns: the angle between x and y in degrees.
Return type: (float)
-
vasppy.cell.
rotation_matrix
(axis, theta)[source]¶ Return the 3D rotation matrix associated with counterclockwise rotation about the given axis by theta radians.
Parameters: - axis (np.array) – length 3 numpy array defining the axis of rotation.
- theta (float) – rotation angle in radians.
Returns: the corredponding rotation matrix.
Return type: (np.array)
vasppy.configuration module¶
-
class
vasppy.configuration.
Configuration
(cell, atoms)[source]¶ Bases:
object
A Configuration object stores a single structure.
vasppy.doscar module¶
-
class
vasppy.doscar.
Doscar
(filename, ispin=2, lmax=2, lorbit=11, spin_orbit_coupling=False, read_pdos=True, species=None)[source]¶ Bases:
object
Contains all the data in a VASP DOSCAR file, and methods for manipulating this.
-
number_of_channels
¶
-
number_of_header_lines
= 6¶
-
pdos_select
(atoms=None, spin=None, l=None, m=None)[source]¶ Returns a subset of the projected density of states array.
-
vasppy.grid module¶
vasppy.kpoints module¶
vasppy.optics module¶
functions for working with optical properties from vasprun.xml
-
vasppy.optics.
absorption_coefficient
(dielectric)[source]¶ Calculate the optical absorption coefficient from an input set of pymatgen vasprun dielectric constant data.
Parameters: dielectric (list) – A list containing the dielectric response function in the pymatgen vasprun format.
element 0: list of energieselement 1: real dielectric tensors, in[xx, yy, zz, xy, xz, yz]
format.element 2: imaginary dielectric tensors, in[xx, yy, zz, xy, xz, yz]
format.Returns: absorption coefficient using eV as frequency units. Return type: (np.array) Notes
The absorption coefficient is calculated as
\[\alpha = \frac{2\sqrt{2} \pi}{\lambda} \sqrt{-\epsilon_1+\sqrt{\epsilon_1^2+\epsilon_2^2}}\]
-
vasppy.optics.
matrix_eigvals
(matrix)[source]¶ Calculate the eigenvalues of a matrix.
Parameters: matrix (np.array) – The matrix to diagonalise. Returns: Array of the matrix eigenvalues. Return type: (np.array)
-
vasppy.optics.
parse_dielectric_data
(data)[source]¶ Convert a set of 2D vasprun formatted dielectric data to the eigenvalues of each corresponding 3x3 symmetric numpy matrices.
Parameters: data (list) – length N list of dielectric data. Each entry should be a list of [xx, yy, zz, xy, xz, yz ]
dielectric tensor elements.Returns: - a Nx3 numpy array. Each row contains the eigenvalues
- for the corresponding row in data.
Return type: (np.array)
-
vasppy.optics.
to_matrix
(xx, yy, zz, xy, yz, xz)[source]¶ Convert a list of matrix components to a symmetric 3x3 matrix. Inputs should be in the order xx, yy, zz, xy, yz, xz.
Parameters: - xx (float) – xx component of the matrix.
- yy (float) – yy component of the matrix.
- zz (float) – zz component of the matrix.
- xy (float) – xy component of the matrix.
- yz (float) – yz component of the matrix.
- xz (float) – xz component of the matrix.
Returns: The matrix, as a 3x3 numpy array.
Return type: (np.array)
vasppy.outcar module¶
-
vasppy.outcar.
fermi_energy_from_outcar
(filename='OUTCAR')[source]¶ Finds and returns the fermi energy. Args: -filename: the name of the outcar file to be read
Returns: The fermi energy as found in the OUTCAR Return type: (Float)
-
vasppy.outcar.
final_energy_from_outcar
(filename='OUTCAR')[source]¶ Finds and returns the energy from a VASP OUTCAR file, by searching for the last energy(sigma->0) entry.
Parameters: filename (Str, optional) – OUTCAR filename. Defaults to ‘OUTCAR’. Returns: The last energy read from the OUTCAR file. Return type: (Float)
-
vasppy.outcar.
potcar_eatom_list_from_outcar
(filename='OUTCAR')[source]¶ Returns a list of EATOM values for the pseudopotentials used.
Parameters: filename (Str, optional) – OUTCAR filename. Defaults to ‘OUTCAR’. Returns: A list of EATOM values, in the order they appear in the OUTCAR. Return type: (List(Float))
-
vasppy.outcar.
reciprocal_lattice_from_outcar
(filename)[source]¶ Finds and returns the reciprocal lattice vectors, if more than one set present, it just returns the last one. :param filename: The name of the outcar file to be read :type filename: Str
Returns: The reciprocal lattice vectors. Return type: List(Float)
vasppy.pimaim module¶
vasppy.polyhedron module¶
vasppy.poscar module¶
vasppy.procar module¶
-
class
vasppy.procar.
Procar
(spin=1)[source]¶ Bases:
object
-
effective_mass_calc
(k_point_indices, band_index, reciprocal_lattice, spin=1, printing=False)[source]¶
-
-
vasppy.procar.
area_of_a_triangle_in_cartesian_space
(a, b, c)[source]¶ Returns the area of a triangle defined by three points in Cartesian space.
Parameters: - a (np.array) – Cartesian coordinates of point A.
- b (np.array) – Cartesian coordinates of point B.
- c (np.array) – Cartesian coordinates of point C.
Returns: the area of the triangle.
Return type: (float)
-
vasppy.procar.
least_squares_effective_mass
(cartesian_k_points, eigenvalues)[source]¶ Calculate the effective mass using a least squares quadratic fit.
Parameters: - cartesian_k_points (np.array) – Cartesian reciprocal coordinates for the k-points
- eigenvalues (np.array) – Energy eigenvalues at each k-point to be used in the fit.
Returns: The fitted effective mass
Return type: (float)
Notes
If the k-points do not sit on a straight line a ValueError will be raised.
-
vasppy.procar.
points_are_in_a_straight_line
(points, tolerance=1e-07)[source]¶ Check whether a set of points fall on a straight line. Calculates the areas of triangles formed by triplets of the points. Returns False is any of these areas are larger than the tolerance.
Parameters: - points (list(np.array)) – list of Cartesian coordinates for each point.
- (optional (tolerance) – float): the maximum triangle size for these points to be considered colinear. Default is 1e-7.
Returns: True if all points fall on a straight line (within the allowed tolerance).
Return type: (bool)
-
vasppy.procar.
two_point_effective_mass
(cartesian_k_points, eigenvalues)[source]¶ Calculate the effective mass given eigenvalues at two k-points. Reimplemented from Aron Walsh’s original effective mass Fortran code.
Parameters: - cartesian_k_points (np.array) – 2D numpy array containing the k-points in (reciprocal) Cartesian coordinates.
- eigenvalues (np.array) – numpy array containing the eigenvalues at each k-point.
Returns: The effective mass
Return type: (float)
vasppy.rdf module¶
vasppy.summary module¶
-
class
vasppy.summary.
Summary
(directory='.')[source]¶ Bases:
object
TODO Document Summary class
-
functional
¶ String description of the calculation functional.
- Recognises:
- PBE
- PBEsol
- PBE-based hybrids:
- PBE0 (alpha=0.25, no screening)
- HSE06 (alpha=0.25, mu=0.2)
- generic hybrids (alpha=?, no screening)
- generic screened hybrids (alpha=?, mu=?)
Returns: String describing the calculation functional. Return type: (Str)
-
parse_vasprun
()[source]¶ Read in vasprun.xml as a pymatgen Vasprun object.
Parameters: None – Returns: None - None:
- If the vasprun.xml is not well formed this method will catch the ParseError and set self.vasprun = None.
-
stoich
¶
-
supported_flags
= {'cbm': 'Vasprun conduction band minimum', 'converged': 'converged', 'description': 'Description', 'directory': 'directory', 'eatom': 'POTCAR EATOM values', 'ediffg': 'ediffg', 'encut': 'encut', 'energy': 'Energy', 'functional': 'functional', 'ibrion': 'ibrion', 'k-points': 'k-points', 'lreal': 'LREAL', 'md5': 'md5', 'nelect': 'NELECT', 'notes': 'Notes', 'plus_u': 'Dudarev +U parameters', 'potcar': 'POTCAR', 'status': 'Status', 'stoichiometry': 'Stoichiometry', 'title': 'Title', 'track': 'tracking for files', 'type': 'Type', 'vbm': 'Vasprun valence band maximum', 'version': 'VASP executable version'}¶
-
-
vasppy.summary.
find_vasp_calculations
()[source]¶ Returns a list of all subdirectories that contain either a vasprun.xml file or a compressed vasprun.xml.gz file.
Parameters: None – Returns: list of all VASP calculation subdirectories. Return type: (List)
-
vasppy.summary.
potcar_spec
(filename)[source]¶ Returns a dictionary specifying the pseudopotentials contained in a POTCAR file.
Parameters: filename (Str) – The name of the POTCAR file to process. Returns: - A dictionary of pseudopotential filename: dataset pairs, e.g.
- { ‘Fe_pv’: ‘PBE_54’, ‘O’, ‘PBE_54’ }
Return type: (Dict)