PyProcar: A Python library for electronic structure pre/post-processing

TL;DR

PyProcar is a versatile Python library that enables detailed visualization and analysis of electronic structures, including band structures, Fermi surfaces, and atomic effects, compatible with multiple electronic structure codes.

Contribution

It introduces a flexible, command-line compatible tool for pre/post-processing electronic structure data, supporting multiple formats and advanced features like band unfolding and property projection.

Findings

Supports band structure and Fermi surface visualization with property coloring.

Compatible with VASP, ABINIT, and adaptable to other formats.

Enables analysis of atomic effects, spin textures, and band unfolding.

Abstract

The PyProcar Python package plots the band structure and the Fermi surface as a function of site and/or s,p,d,f - projected wavefunctions obtained for each $k$-point in the Brillouin zone and band in an electronic structure calculation. This can be performed on top of any electronic structure code, as long as the band and projection information is written in the PROCAR format, as done by the VASP and ABINIT codes. PyProcar can be easily modified to read other formats as well. This package is particularly suitable for understanding atomic effects into the band structure, Fermi surface, spin texture, etc. PyProcar can be conveniently used in a command line mode, where each one of the parameters define a plot property. In the case of Fermi-surfaces, the package is able to plot the surface with colors depending on other properties such as the electron velocity or spin projection. The mesh used to calculate the property does not need to be the same as the one used to obtain the Fermi surface. A file with a specific property evaluated for each $k$-point in a $k-$mesh and for each band can be used to project other properties such as electron-phonon mean path, Fermi velocity, electron effective mass, etc. Another existing feature refers to the band unfolding of supercell calculations into predefined unit cells.

Equations

$$|⟨Y_{lm}^{\alpha }⟩{\varphi }_{nk}{|}^2$$

$$\frac{1}{2}\sum _{\alpha ,\beta =1}^2{\sigma }_{\alpha ,\beta }^i⟨{\psi }_{n,\mathbf{k}}^{\alpha }⟩Y_{lm}^{\alpha }⟨Y_{lm}^{\beta }⟩{\psi }_{n,\mathbf{k}}^{\beta }$$

$${\varphi }_{n\mathbf{k}}$$

Full text

††journal: Computer Physics Communications

PROGRAM SUMMARY

Program Title: PyProcar

Licensing provisions: GPLv3

Programming language: Python

*Supplementary material:

Nature of problem:* To automate, simplify and serialize the analysis of band structure and Fermi surface, especially for high throughput calculations.

Solution method: Implementation of a Python library able to handle, combine, parse, extract, plot and even repair data from density functional calculations. PyProcar uses color maps on the band structures or Fermi surfaces to give a simple representation of the relevant characteristics of the electronic structure. *Additional comments:

*Features: PyProcar can produce high-quality figures of band structures and Fermi surfaces (2D and 3D), projection of atomic orbitals, atoms, and/or spin components.

Restrictions: Only the VASP package is currently fully supported, the latest version of Abinit is partially supported (it will be fully supported in the the Abinit versions 9.x). The PROCAR file format can easily be implemented within any DFT code.

1 Introduction

Density Functional Theory (DFT) is one of the most widely used methodologies for electronic structure calculations in materials science [1, 2, 3]. Its very good quality/computational-cost ratio, together with the emergence of high performance computing, has reshaped the field of computational materials research. Decades efforts in developing DFT for simulation programs have led to highly efficient DFT codes which include, among others, Abinit [4], VASP [5, 6], Siesta [7] or Quantum Espresso [8, 9], which are capable of exploring remarkable material properties. Thorough calculations with these DFT codes have already been done for a wide range of compounds and the accumulated data can be publicly accessed from databases such as Materials Project [10], AFLOW [11], Materiae [12]. This new approach to materials has had a large impact into condensed matter physics developing new paradigms such as the quantum materials [13]. However, countless possibilities of exploring novel materials still exist and thus the necessity of efficient and reliable DFT pre- and post-processing tools.

As computational capabilities increase, the size and complexity of the systems under study increase as well, hindering the analysis and the capability of abstraction on key aspects of these systems. As physicists, we encountered this problem and realized that we require a tool designed to gain insight from our calculations but spending less time on coding. In particular, we have three requirements for such a tool: (i) the generation of a graphical representation of the quantity of interest -surface states, spin-texture, orbital projection, etc.- should take only a single line of input. This approach also enables the scripting of similar analysis, a must for high throughput calculations. (ii) We also realized the need of producing high-quality graphics, ideally using a vector format, suitable for further post-processing. (iii) Finally, there are several issues preventing a straightforward analysis of electronic structure results, e.g., the size of output data files generated by DFT codes can be extremely large, thus making the post-processing of data very slow by excessive memory usage. In such a case, the post-processing tool, not the researcher, should be capable of taking care of large data files in a smart way and enabling the user to extract essential information while ignoring less interesting information from a large data file.

This paper focuses on the presentation of PyProcar, a robust, open-source Python library used for pre- and post-processing of the electronic structure data coming from DFT calculations. PyProcar is capable of performing a multitude of tasks including plotting plain and spin/atom/orbital projected band structures and Fermi surfaces- both in 2D and 3D, Fermi velocity plots, unfolding bands of a super cell, comparing band structures from multiple DFT calculations and generating a $k$-path for a given crystal structure. Our code has matured throughout the last decade and currently has over forty users who are members of the PyProcar forum [14]. To retrieve a reliable number of users of an Open Source code is a difficult task due to inconsistencies in keeping track of downloads or installations, nevertheless, according to Google BigQuery [15] PyProcar has been downloaded almost 9000 times and PyPI Download Stats [16] counts more than 100 installations of PyProcar per month based on the download analytics from the Python Packaging Index (PyPI) [17]. Even though, PyProcar is written to parse the PROCAR file from the VASP package (and from a fork of Abinit, which will be included in the next main production version), due to its object-oriented approach, PyProcar can be easily adapted to handle the output from other DFT codes. In Section 2 we will briefly explain some basic aspects of PyProcar. Then in Section 3 we present some examples illustrating the capabilities of the PyProcar code for analysis of electronic structure.

2 Library Overview

2.1 Electronic structure projection

The projection of the Kohn-Sham states over atomic orbitals and spins can give a large amount of information about the calculated system. However, such amount of information needs to be post-processed to extract the physical insight. In the VASP code, this information is written into the PROCAR file, which arranges the information of projections in blocks, as seen in Fig. 1 This block is repeated if calculations are spin polarized and for noncollinear spin calculations, three more blocks are present which correspond to $S_{x}$, $S_{y}$ and $S_{z}$ spin directions. There also exist a special format of PROCAR file that includes the phases of the wave function, as described in Fig. 2, which is an important quantity required for band unfolding. The direction of the atomic orbitals ($p_{x}$, $p_{y}$, $p_{z}$, etc.) is defined with respect to the Cartesian coordinates $x$, $y$ and $z$.

The site projected wave function in the PROCAR file is calculated by projecting the Kohn-Sham wave functions onto spherical harmonics that are non-zero within spheres of a Wigner-Seitz radius around each ion by:

[TABLE]

where, $Y^{\alpha}_{lm}$ are the spherical harmonics centered at ion index $\alpha$ with angular moment $l$ and magnetic quantum number $m$, and $\phi_{n\mathbf{k}}$ are the Kohn-Sham wave functions. In general, for a non-collinear electronic structure calculation the same equation is generalized to:

[TABLE]

where $\sigma^{i}$ are the Pauli matrices with $i=x,y,z$ and the spinor wavefunction $\phi_{n\mathbf{k}}$ is now defined as

[TABLE]

This projection is performed for every $k$-point used in the DFT calculation, for every energy band and every atom of the crystal. Files can be stored as the absolute value of the projections, as it is usually needed for the band analysis or by also including the complex nature of the projections, which are used for the electron band folding. The raw data written on the PROCAR file captures most of the material details, essential to investigate interesting materials properties. But on the downside, these files are often too big to be analyzed without additional computational support. For instance, a spin non-polarized calculation involving 10 atoms, with 100 bands and 100 $k$-points gives almost $10^{6}$ entries in the PROCAR file, and the file size grows with the number of atoms as $\mathcal{O}(n^{2})$. Therefore, it is essential to reduce the PROCAR file size in an efficient way without losing the key information stored in the file. The PyProcar parser can be used to manipulate, reduce, filter, and comprehend the information stored in the PROCAR file.

2.2 PyProcar Design

PyProcar is an object-oriented code, and despite its name, only the file-related classes (mostly the parsing) are related to the VASP’s PROCAR and OUTCAR. The OUTCAR gives a summary of the calculation including information about the electronic steps, eigenvalues, Fermi energy, forces on the atoms, etc. PyProcar code parses the OUTCAR file and stores the Fermi energy, which is later used to shift the Fermi level to zero in band structure plots. The reciprocal lattice vectors are also parsed from OUTCAR. PyProcar is currently capable of parsing the outputs from a recent version of Abinit and given the code’s flexibility, the extension to other codes is trivial. Once the data is stored in memory, PyProcar can be used regardless of the employed DFT code. PyProcar consists of several classes and functions, as mentioned below. The functions provide a high-level interaction with the user, allowing the relevant task to get done with minimum instructions. For convenience, and to allow an easy scripting, all these functions are fully decoupled. Also, the manipulation of the data can be saved/exported as a new file on the disk.

The low level work is carried by classes, each handling one specific task or sub-task, they are fairly complex and no interaction with the user is expected unless a new feature is needed. By design, the classes are as loosely connected as possible. Each class uses the logging module to provide a user-defined level of verbosity, which is very useful for debugging.

Fig. 3 displays an overview of the PyProcar library. PyProcar can be used to generate files required for DFT calculations such as a suitable KPOINTS file for both self-consistent and non self-consistent DFT calculations. The structures can be generated manually, or from one of many databases publicly available such as Materials Project [10], AFLOW [11], etc. Once the DFT calculation is complete, PyProcar uses the VASP generated outputs, PROCAR and OUTCAR, for further post-processing.

The core classes within PyProcar consist of Procarparser, Utilsprocar and Procarselect which parses the PROCAR data and stores them in organized arrays for later calculations or analysis. These core classes are explained briefly below.

Utilsprocar: This class contains modules to parse the OUTCAR file from a DFT calculation. It reads and stores the Fermi energy and the reciprocal lattice vectors. The Fermi energy is used to shift the Fermi level to zero in band structure calculations. The availability of the reciprocal lattice vectors enable Utilsprocar to convert the $k$-mesh grid between direct and Cartesian coordinates. The Fermi energy and reciprocal lattice vector could also be provided manually through command line inputs. 2. 2.

ProcarParser: As the name suggests, this class contains modules required to parse the PROCAR file. It reads information related to bands, $k$-points and orbitals, and saves the data in the memory. It is designed to be somewhat resilient to errors derived from the fixed-format of the PROCAR file (i.e., missing blank spaces, use of *** characters as an index, etc.) 3. 3.

ProcarSelect After the PROCAR file is parsed and stored in memory, this class manipulates the required orbital(s), atom(s) and or spin information separately.

2.3 Installation

The latest stable version of PyProcar, version 3.8.5 at the time of writing this paper, can be installed using the Python Packaging Index (pip) using the following command:

⬇

pip install pyprocar

The project’s GitHub repository is located at https://github.com/romerogroup/pyprocar. An easy to follow documentation with examples can be found at https://romerogroup.github.io/pyprocar/. PyProcar is supported by both Python 2.x and 3.x.

3 Examples: MgB2, BiSb, TaSb, and SrVO3

3.1 The MgB2, BiSb, TaSb and SrVO3 crystal structures

In this section, we choose four example systems, namely MgB2, BiSb, TaSb, and SrVO3, to illustrate the capabilities of the PyProcar code. Fig. 4 shows the crystal structure of these systems. The DFT calculations were performed using the Vienna Ab Initio Simulation Package (VASP) version 5.4.4. Similar calculations can be performed using the forthcoming version of Abinit, version 9.x, where the keyword prtprocar needs to be added in the input.

For the strongly correlated cubic perovskite SrVO3 (crystal symmetry $Pm\bar{3}m$ and lattice constant 3.84 Å), the plain and projected band structures were evaluated within the generalized gradient approximation (GGA) using PBE functional [18]. An 8 $\times$ 8 $\times$ 8 Monkhorst-Pack [19] $k$-mesh and an energy cutoff of 600 eV were required for electron wavefunction convergence. We considered 10 valence electrons of Sr (3s23p64s2), 5 valence electrons of V (3d34s2) and 6 valence electrons of O (2s22p4) in the PAW pseudopotential.

The superconducting material MgB2 of the crystal symmetry $P6/mmm$ and lattice parameters 3.07 Å and 3.53 Å for $a$ and $c$, respectively, was used to showcase the Fermi surface plotting and band unfolding capabilities of PyProcar. For the former, we used a 10 $\times$ 10 $\times$ 10 Monkhort-Pack grid with an energy cut-off of 700 eV within the GGA with the PBE exchange-correlation function. The latter was done with a 5 $\times$ 5 $\times$ 5 Monkhort-Pack grid with an energy cut-off of 500 eV within the GGA with the PBEsol exchange-correlation function [20]. Both of these calculations were done with 2 valence electrons of Mg (3s2) and 3 valence electrons of B (2s22p1) in the PAW pseudopotential.

To demonstrate the 2D spin texture plotting capability of PyProcar we use the Rashba semiconductor BiSb monolayer which belongs to the space group $P3m1$ with lattice parameters $a=b=$ 4.26 Å. The calculation, which included spin-orbit coupling (SOC) was performed with a $k$-mesh grid of 10 $\times$ 10 $\times$ 1 and an energy cut-off of 650 eV with PBE using 15 valence electrons of Bi (5d106s26p3) and 5 valence electrons of Sb (5s25p3) in the PAW pseudopotential.

PyProcar can also be used to investigate the band degeneracy and Dirac/Weyl points in topological materials. We use the material TaSb belonging to the space group $P\bar{6}m2$ with lattice parameters a = b = 3.58 Å, c = 3.81 Å to demonstrate this feature. An energy cut-off of 600 eV and $k$-mesh grid of 10 $\times$ 10 $\times$ 10 along with PBE and SOC were used for the calculation considering 5 valence electrons of Ta (5d36s2) and 5 valence electrons of Sb (5s25p3) in the PAW pseudopotential.

3.2 Generating a $k$-path

In order to plot a band structure, one must define a set of $k$-points following a desired $k$-path in momentum space. A common option is to employ the AFLOW software framework [11], the SeeK-path module in Materials Cloud [21] or the Bilbao Crystallographic Server [22, 23, 24] to generate the $k$-path. However, PyProcar’s $k$-path generation utility enables a the user to automatically generate a suitable and sufficient $k$-path given the crystal structure, typically read from the POSCAR file. The $k$-path generation utility within PyProcar is based on the Python library seekpath developed by Hinuma et al. [25]. The red line in Fig. 4 (c) depicts such a $k$-path of a hexagonal lattice generated using the algorithms of seekpath.

General format:

⬇

pyprocar.kpath(infile,grid_size,with_time_reversal,recipe,threshhold,symprec,angle_tolerence)

Usage:

⬇

pyprocar.kpath(‘POSCAR’,40,True,‘hpkot’,1e-07,1e-05,-1.0)

More details regarding these parameters can be found in Ref. [26].

3.3 Repair

This utility is used to repair the ill-formatting of the PROCAR file due to the erroneous file handling in Fortran, particularly in a VASP calculation. This prevents issues arising from the lack of white space between a number and a negative sign, for instance 0.000000-0.5000000. Typically, pyprocar.repair() is recommended to be applied before using any other utility.

Usage:

⬇

pyprocar.repair(‘PROCAR’,‘PROCAR-repaired’)

3.4 $k$-mesh generator

This utility can be used to generate a 2D $k$-mesh centered at a given $k$-point and in a given $k$-plane. This is particularly useful in computing 2D spin-textures and plotting 2D Fermi surfaces. For example, the following command creates a 2D $k_{x}-k_{y}$ mesh centered at the $\Gamma$ point ($k_{z}=0$) ranging from coordinates (-0.5, -0.5, 0.0) to (0.5, 0.5, 0.0) with a grid size of 0.02:

General format:

⬇

pyprocar.generate2dkmesh(x1,y1,x2,y2,grid_size)

Usage:

⬇

pyprocar.generate2dkmesh(-0.5,-0.5,0.5,0.5,0.02)

3.5 Band structure

PyProcar goes beyond the conventional plain band structure to plot the projected bands that carry even more information, which will be described shortly. The projected bands are color coded in an informative manner to portray fine details. PyProcar is capable of labeling the $k$-path names automatically, however, the user can manually input them as desired. This will be showcased in the next section.

No projection

This is the most basic type of band structure. No projection information is contained here. See Fig. 5 (a) for a plain band structure of SrVO3. In order to use the plain mode one sets mode=‘plain’. elimit sets the energy window limits. outcar specifies the OUTCAR file. For Abinit calculations, abinit_output is used instead. color lets the user use any color available in the matplotlib [27] package.

Usage:

⬇

pyprocar.bandsplot(‘PROCAR-repaired’,outcar=‘OUTCAR’,elimit=[-2,2],mode=‘plain’,color=‘blue’) 2. 2.

Spin projection

For collinear spin polarized and noncollinear spin calculations of DFT codes, PyProcar is able to plot the bands of each spin channel or direction separately. An example for a collinear spin polarized calculation is given in Fig. 5 (b) where blue corresponds to spin down channel and red to spin up channel. For this case setting spin=‘0’ plots the total density from both spin channels and spin=‘1’ plots the spin channels separately.

Usage:

⬇

pyprocar.bandsplot(‘PROCAR_repaired’,outcar=‘OUTCAR’,elimit=[-5,5],kticks=[0,39,79,119,159],knames=[‘G’,‘X’,‘M’,‘G’,‘R’],cmap=‘seismic’,

mode=‘parametric’,spin=‘1’)

For noncollinear spin calculations, spin=1,2,3 corresponds to spins oriented in Sx, Sy and Sz directions respectively. For parametric plots such as spin, atom and orbitals, the user should set mode=‘parametric’. knames and kticks corresponds to the labels and the number of grid points between the high symmetry points in the $k$-path used for the band structure calculation. At the end of this section we explain how to retrieve this automatically. LaTeXentries, such r‘\Gamma’ also can be used. cmap refers to the matplotlib color map used for the parametric plotting and can be modified by using the same color maps used in matplotlib. 3. 3.

Atom projection

The projection of atoms onto bands can provide information such as which atoms contribute to the electronic states near the Fermi level. PyProcar counts each row of ions in the PROCAR file, starting from zero. In the example of a five atom SrVO3, the indexes of atoms for Sr, V and the three O atoms would be 1,2 and 3,4,5 respectively. It is also possible to include more than one type of atom. See Fig. 5 (c) for an example of V atom projected bands in SrVO3.

Usage:

⬇

pyprocar.bandsplot(‘PROCAR_repaired’,outcar=‘OUTCAR’,elimit=[-5,5],kticks=[0,39,79,119,159],knames=[‘G’,‘X’,‘M’,‘G’,‘R’],cmap=‘seismic’, mode=‘parametric’,atoms=[1]) 4. 4.

Orbital projection

The projection of atomic orbitals onto bands is also useful to identify the contribution of orbitals to bands. For instance, to identify correlated $d$ or $f$ orbitals in a strongly correlated material near the Fermi level. It is possible to include more than one type of orbital projection. Fig. 5 (d) displays an orbital projected band structure of SrVO3. The mapping of the index of orbitals to be used in orbitals is as follows (this is the same order from the PROCAR file, see Fig. 2).

[TABLE]

Usage: To project all five $d$-orbitals:

⬇

pyprocar.bandsplot(‘PROCAR_repaired’,outcar=‘OUTCAR’,elimit=[-5,5],kticks=[0,39,79,119,159],knames=[‘G’,‘X’,‘M’,‘G’,‘R’],cmap=‘seismic’,mode=‘parametric’,orbitals=[4,5,6,7,8])

If a KPOINTS file from a VASP calculation is present, PyProcar automatically retrieves kticks and knames and labels the band structure plots accordingly for any of the above cases and also for the bandscompare() function in section 3.10. For example:

⬇

pyprocar.bandsplot(‘PROCAR_repaired’,outcar=‘OUTCAR’,mode=’plain’,kpointsfile=‘KPOINTS’)

One or many of the above can be combined together to allow the user to probe into more specific queries such as a collinear spin projection of a certain orbital of a certain atom.

Different modes of band structures are useful for obtaining information for different cases. The four modes available within PyProcar are plain, scatter, parametric and atomic. The plain bands contain no projection information. The scatter mode creates a scatter plot of points. The parametric mode interpolates between points to create bands which are also projectable. Finally, the atomic mode is useful to plot energy levels for atoms.

These projected band structures play a vital role in revealing the physics of a system. For instance, the V4+ projected bands displayed in Fig. 5 (c) correctly captures the contribution from the V-$3d$ electrons with O-$2p$ states from -6 eV to -1 eV . Furthermore, in Fig. 5 (d) where $e_{g}$ states are projected, one can clearly notice the octahedral crystal field splitting in SrVO3 where $t_{2g}$ states exist below $e_{g}$ states, as expected.

3.6 3D Fermi surface

PyProcar’s 3D Fermi surface utility is able to generate Fermi surface plots projected over spin, atoms and orbitals or a combination of one or many of each. This utility is also capable of projecting external properties that are provided on a mesh grid in momentum space. This feature is useful when one wants to project properties that are not provided in a PROCAR file such as Fermi velocity, electron-phonon coupling and electron effective mass. We divide this section into three sub sections, plain Fermi surface, projection of properties from PROCAR and projection of properties from external file.

Plain Fermi surface

The ProcarParser class provides eigenvalues on a momentum space mesh grid. In all of the 3D Fermi surface functions in PyProcar the eigenvalues are interpolated using Fourier transform. This type of interpolation is suitable for eigenvalues as they are periodic in nature. The points in the mesh grid is generated in reduced space and some of the points sampled by VASP might not be in the first Brillouin Zone (BZ). The points outside of the first BZ can be returned to the first BZ using the reciprocal lattice vectors. As transforming the points to the first BZ will cause the distortion in the shape of the iso-surface, it is better to add the points missing from the first BZ then remove the points lying in the second BZ. In order to identify these points, we first create the boundaries of the first BZ by creating the Wigner-Seitz cell from the reciprocal lattice vectors. The BZ is a convex hull created by the points located at the corner of the BZ. To figure out if a point is inside the first BZ or outside the following procedure is performed on each point. A point is added to the collection of the points on the BZ corners and a new convex hull is created, if the convex hull is the same as the original BZ, the point is located inside, and if the BZ is different, the point is located outside of the first BZ. Since this operation can be slow in a scripting language like Python, we use data parallelism. This task is offloaded to a number of worker processes, which can be defined by the user. After obtaining enough points we use Lewiner marching cubes algorithm provided by scikit-image package[28]. This function provides vertices and faces required to create the Fermi surface. The last step that is needed is to transform the coordinates of the vertices from the reduced coordinates to the Cartesian coordinates. To visualize the Fermi surface, we have provided four different plotting packages, mayavi[29], matplotlib[30], plotly[31] and ipyvolume [32] which can be chosen by the user. The figures shown here are produced using mayavi.

Usage:

⬇

pyprocar.fermi3D(procar,outcar,bands,scale=1,mode=‘plain’,st=False,**kwargs)

The main arguments in this function are procar, outcar, bands, scale, mode and st, where procar and outcar are the names of the input PROCAR and OUTCAR files respectively, bands is an array of the bands that are desired to be plotted. Note if bands = -1, the function will try to plot all the bands provided in the PROCAR file. The $k$-mesh will be interpolated by a factor of scale in each direction. The st tag controls the spin-texture plotting, and mode determines the type of projection of colors. There are additional keyword arguments that can be accessed in the help section of this function, such as face_color, cmap, atoms, orbitals, energy, transparent, nprocess etc. Fig. 6 shows the Fermi surface of MgB2 generated using ‘plain’ mode of this function.

Surface coloring based on properties from PROCAR

Similar to the bandsplot section one can choose to project the contribution of different properties provided in the PROCAR file, such as atom, orbital and spin contributions. The projection can be represented by different color mapping schemes chosen by the user. The projection is not restricted to only one property at a time, so it can be chosen from all the provided properties. For example, one might want to see the contribution of the orbitals px, py, pz from specific atoms, this function will parse the desired contributions and projects the sum of contributions on each face. To do so, we perform an interpolation on the data provided by the DFT package and evaluate the function at the center of each face created in the previous section. To use this functionality one has to change the mode from ‘plain’ to ‘parametric’ and choose the atoms, orbitals, spin that are desired to be projected. As an example, we project the different contribution of different p orbitals of Boron atom in Fig. 7. The results show the middle cylinder is mostly comprised of px and py orbitals while the bands closer to the edges of BZ are pz orbitals. This is in agreement with the calculations performed using the SIESTA DFT package [33].

For noncolinear calculations, this function is able to plot arrows in the direction of the spinors provided in the PROCAR file. To turn this functionality on the one can set st=True to turn the spin-texture ON. The user can choose between coloring all the arrows originated from one band with a specific color, or project the contribution of that arrow in a specific Cartesian direction. We plot two examples to demonstrate this functionality in Fig. 8, the spin-texture of BiSb at 0.60 eV above the Fermi level and the spin-texture of SrVO3. To better represent the spin-texture we use the key argument transparent=True which changes the opacity of the Fermi surface to zero.

Surface coloring based on properties obtained from an external file

Similar to the previous section, this function is able to read an external file, containing information about a scalar or a vector field in BZ and project the field on the Fermi surface. This file does not need to have the same mesh grid as the PROCAR file as long as the mesh sampling is fine enough. This function performs an interpolation on the provided data and evaluates functions at the center of each face on the Fermi surface. The external file should have the following format.

{verbnobox}

[] band = ¡band number¿ ¡kx1¿ ¡ky1¿ ¡kz1¿ ¡color1¿ ¡kx2¿ ¡ky2¿ ¡kz2¿ ¡color2¿ ¡kx3¿ ¡ky3¿ ¡kz3¿ ¡color3¿ … band = ¡band number¿ …

The function matches information about the first band present in the file to the first band requested to be plotted, second band present in the file to the second band requested to be plotted, and so on. As an example, we have plotted the Fermi velocity of MgB2 calculated by numerically evaluating the energy gradient from the PROCAR file in Fig. 9. Of course, one can acquire more accurate Fermi velocity by means of the Wannier interpolation method [34], however the following plot shows close enough accuracy for an example designed to represent the capabilities of this functional mode.

3.7 Handling big files: filtering the selected data and reducing the memory requirement for post-processing

A simpler version of a PROCAR file containing only a subset of information from the original dataset can be generated with this utility. This feature is very useful when there are many bands in the PROCAR file (e.g. in heterostructures or supercells calculations) making the file size enormously large for post-processing while only bands near the Fermi-level are needed for analysis. In this case, one can filter data of selected bands near the Fermi-level. This considerably reduces the file size and makes the post-processing of data faster. For instance, in graphene/MoS2 (5:4) heterostructures [35, 36] the size of a VASP generated PROCAR file is 1.49 GB (98 atoms, 584 electrons, 804 bands in total, 63 $k$-points, and spin-orbit coupling included). However, by filtering only 30 bands above and below the Fermi-level ($i.e.$ band indexes 554–614), the file size reduces to 108 MB and data processing becomes much faster. In the same way, one could use the filter utility to filter the PROCAR file to extract information regarding particular spins, atoms, or orbitals in a relatively smaller PROCAR-new file.

The following example extracts information of bands ranging from index 50 to 70 from a PROCAR-repaired file (Fermi-level is near band #60) while ignoring all other bands located far from the Fermi-level, and stores resulting dataset in a new file named PROCAR-repaired-band50-70. Now the new PROCAR-repaired-band50-70 file can be used for further post-processing of data at relatively low memory requirements.

Usage:

⬇

pyprocar.filter(‘PROCAR-repaired’,‘PROCAR-repaired-band50-70’,bands=[50,70])

3.8 2D spin-texture

This module can be utilized to visualize the constant energy surface spin textures in a given system. This feature is particularly useful in identifying Rashba and Dresselhaus type spin-splitting effects, analyzing the topology of Fermi surface, and examining Lifshitz transitions, as demonstrated in Refs. [37, 38, 39, 40, 36, 41, 42]. To plot 2D spin texture, we require a 2D $k$-grid centered a certain special $k$-point in Brillouin zone near which we want to examine the spin-texture in $k$-space (see section 3.4 regarding generation of 2D $k$-mesh).

Usage: To plot $S_{x}$ spin component at a constant energy surface $E=E_{F}+0.60\,eV$ (spin=‘1’, ‘2’, ‘3’ for $S_{x}$, $S_{y}$, $S_{z}$, respectively)

⬇

pyprocar.fermi2D(‘PROCAR-repaired’, outcar=‘OUTCAR’, st=True, energy=0.60, noarrow=True, spin=‘1’)

For example, in Fig. 10 we plot $S_{x}$, $S_{y}$, and $S_{z}$ spin projections on the 2D spin texture of BiSb monolayer, which is a Rashba semiconductor (Rashba spin-splitting takes place the $\Gamma$ point of BZ), computed in a $k_{x}-k_{y}$ mesh centered at the $\Gamma$ point. Fig. 10(a) shows spin-texture calculated above Fermi-energy ($E_{F}$) at constant energy value $E=E_{F}+0.60\,eV$ (conduction bands), and Fig. 10(a) shows the spin-texture calculated below Fermi surface at constant energy value $E=E_{F}-0.90\,eV$ (valence bands). One can notice linear in $k$ Rashba spin-splitting effects in Fig. 10(a), and additional warping effects in Fig. 10(b) appearing due to the higher order $k^{3}$ terms in lower energy valence bands.

One could also plot spin texture using arrows instead of a heat map, as shown in Fig. 11. This can be done by setting tag: noarrow=False.

To set maximum and minimum energy values for color map, one could use vmin and vmax tags. For example, vmin=-0.5, vmax=0.5 in Fig. 10 and Fig. 11.

3.9 Identification of Weyl points

Weyl points appear in the momentum space when two spin non-degenerate electronic bands (one valence and one conduction) linearly cross or touch each other near the Fermi-level forming a topologically non-trivial gapless point in the energy-momentum space. Depending upon the dispersion of valence and conduction bands, Weyl points can be categorized into two types: (i) Type-I, when Weyl cone is not titled and the crossing bands have opposite Fermi-velocity [43, 44, 45], and (ii) type-II, when Weyl cone is titled such that Weyl points appear at the touching point of an electron and a hole pocket. In the later case, the crossing bands have unidirectional Fermi-velocity with different magnitudes. When such band crossing occurs between spin degenerate bands, thus formed gapless points are referred as Dirac points. Analogous to Weyl points, Dirac points can also be classified into two type-I and type-II categories.

We can use the PyProcar code to check the spin degeneracy of electronic bands near the crossing points and determine the type of gapless point by analyzing the dispersion of bands in the vicinity of the crossing point. Fig. 12 shows the spin projected bandstructure of TaSb, which is a topological metal hosting both type-I and type-II Weyl points. For more details regarding the topological properties of TaSb and other similar topological metals, we refer the reader to Refs. [46, 45].

3.10 Compare Bands

This module is useful to compare different bands from different materials on the same band plot (Fig. 13). The bands are plotted for the same $k$-path in order to have a meaningful comparison but they do not need to have the same number of $k$-points in each interval. The bandscompare() function contains all the parameters that are used in the bandsplot() along with an added feature of displaying a legend to help differentiate between the two different band structures. Different marker styles can be used as well.

Usage:

⬇

pyprocar.bandscompare(‘PROCAR1’,‘PROCAR2’,outcar=‘OUTCAR1’,outcar2=‘OUTCAR2’,cmap=‘seismic’,mode=‘parametric’,marker=‘*’,marker2=‘-.’,elimit=[-5,5],kpointsfile=‘KPOINTS’,legend=‘PRO1’,legend2=‘PRO2’)

3.11 Concatenating multiple calculations

Multiple PROCAR files from multiple DFT calculations can be combined with this utility. For instance, performing DFT calculations for MgB2 supercells are very computationally expensive to do in a single run. This utility is particularly useful in such cases of large systems, where one can split the bandstructure calculations along different high-symmetry directions in BZ, and then concatenate the PROCAR files for each separate $k$-paths, and finally plot the full bandstructure in a single plot. The following command concatenates the PROCAR files obtained from three separate bandstructure calculations done along $\Gamma$-K, K-M, and M-$\Gamma$ $k$-path in hexagonal Brillouin zone.

Usage:

⬇

pyprocar.cat([‘PROCAR_G-K’,‘PROCAR_K-M’,‘PROCAR_M-G’],‘PROCAR_merged’)

To concatenate PROCAR’s generated from Abinit assuming the files are all in the same directory, use the following command:

Usage:

⬇

pyprocar.mergeabinit(‘PROCAR_merged’)

3.12 Band Unfolding

Often times, we need to perform DFT calculations for a supercell geometry rather than the primitive cell. In such cases the band structure becomes quite sophisticated due to the folding of the BZ, and it is difficult to compare the band structure of supercell with that of the primitive cell. The purpose of the band unfolding scheme [47, 48, 49, 50] is to represent the bands within the primitive cell BZ. By calculating the unfolding weight function [50] and plotting the fat bands with the line width proportional to the weight, the unfolded bands can be highlighted. Here we use a $2\times 2\times 2$ MgB2 supercell as an example to show the unfolding of band structure. In the bulk structure, where the primitive cell translation symmetry is preserved, the unfolded band is exactly the same as calculated from primitive cell structure. In a structure with a B atom replaced by Al, the translation symmetry can be seen as approximated, and we can still get an unfolded band, with some band smearing out. By comparing these two, we can clearly see the shifting and smearing of the band. The ASE library [51] was employed to perform atomic and cell manipulations required for the unfolding.

Usage: First, calculate the band structure in the primitive cell BZ. The PROCAR should be produced with the phase factor included, by setting LORBIT=12 in VASP. Then the unfold module can be used to plot the unfolded band.

⬇

import numpy as np

pyprocar.unfold(

    fname=‘PROCAR’,

    poscar=‘POSCAR’,

    outcar=‘OUTCAR’,

    supercell_matrix=np.diag([2, 2, 2]),

    efermi=None,

    shift_efermi=True,

    elimit=(-5, 15),

    kticks=[0, 36, 54, 86, 110, 147, 165, 199],

    knames=[‘$\Gamma$’, ‘K’, ‘M’, ‘$\Gamma$’, ‘A’, ‘H’, ‘L’, ‘A’],

    print_kpts=False,

    show_band=True,

    savefig=‘unfolded_band.png’)

4 Conclusion

PyProcar is a user friendly, open source Python library that can be easily used for a variety of DFT pre- and post-processing calculations. We have demonstrated its capability through providing examples for each functionality through a set of four materials with unique characteristics. PyProcar’s main specialty lies in its ability to project spin, orbitals and atoms in band structures and 2D and 3D Fermi and constant energy surfaces without involving lengthy, complex syntax. The code is freely available to download at https://github.com/romerogroup/pyprocar, and via pip. An easy to follow user manual is available at https://romerogroup.github.io/pyprocar/. The PROCAR format is easy to implement in any DFT code, rendering PyProcar accessible across a wide range of DFT codes. We hope this tool will be useful to computational materials scientists in exploring state-of-the-art novel material which would in turn immensely impact the materials community.

Acknowledgements: This work used the XSEDE which is supported by National Science Foundation grant number ACI-1053575. XH and EB acknowledge the ARC project AIMED and the F.R.S-FNRS PDR project MaRePeThe (GA 19528980). The authors also acknowledge the support from the Texas Advances Computer Center (with the Stampede2 and Bridges supercomputers), the PRACE project TheDeNoMo and on the CECI facilities funded by F.R.S-FNRS (Grant No. 2.5020.1) and Tier-1 supercomputer of the Fédération Wallonie-Bruxelles funded by the Walloon Region (Grant No. 1117545). This work was supported by the DMREF-NSF 1434897, NSF OAC-1740111 and DOE DE-SC0016176 projects. We acknowledge the West Virginia University supercomputing clusters; Spruce Knob and Thorny Flat which were used for the development of the library. FM acknowledges the support from Fondecyt grants #1150806, #1191353, the Center for the Development of Nanoscience and Nanotechnology CEDENNA FB-0807 and the supercomputing infrastructure of the NLHPC (ECM-02). A special thanks goes to Dr. Guillermo Avendaño Franco for his invaluable support.