PyHMA: A VASP post-processor for precise measurement of crystalline anharmonic properties using harmonically mapped averaging

Abstract

We introduce a new Python package (pyHMA) that interfaces with VASP to compute (classical) anharmonic properties of crystalline systems by post-processing data from NVT Born–Oppenheimer ab initio molecular dynamics (AIMD) simulation. It is based on the recently developed harmonically mapped averaging (HMA) method, which leverages the analytically known harmonic behavior to reformulate the direct/conventional ensemble averages in order to significantly improve precision, for a given CPU time. The package consists of two stages: reading AIMD data from vasprun.xml file(s) and then computing anharmonic properties. While the first stage is MD package-dependent, the second one is universal, given that it receives data in the required format. To demonstrate the usage of pyHMA, we compute anharmonic energy and pressure of aluminum fcc crystal at high pressure (≈115GPa) and up to 4000 K (near melting). We further compute anharmonic free energy as a function of temperature, using thermodynamic integration of the HMA anharmonic energy. Although pyHMA currently interfaces with VASP to compute HMA anharmonic energy and pressure, it is moduled in such a way to allow for interfacing with other codes (e.g., LAMMPS) by adding a new reader and can compute other HMA anharmonic properties (e.g., heat capacity) by adding a new method, once relevant data are available. Program summaryProgram title: pyHMACPC Library link to program files:http://dx.doi.org/10.17632/bzgfk52msk.1Licensing provisions: MPL-2.0Programming language: Python 3.7Nature of problem: Theormodynamic properties (e.g., energy, pressure, and heat capacity) of crystalline systems can be decomposed into: lattice (or, property at 0 K), quasiharmonic, and anharmonic contributions. Although the first two are feasible to compute using only a few single-point density functional theory (DFT) calculations, measuring anharmonic contribution requires running ab initio molecular dynamics (AIMD) simulation, which is computationally very demanding using direct ensemble averaging.Solution method: In pyHMA, we are adopting the harmonically mapped averaging (HMA) technique that provides order(s) of magnitude higher precision, in comparison to direct/conventional (Conv) averaging. The package works as a post-processor to VASP AIMD output to provide very precise (and accurate) estimate of anharmonic properties, for a given DFT model, with application to energy and pressure (at this time).Additional comments: The term anharmonicity is commonly used in literature to qualitatively describe a system with no equilibrium configuration at 0 K (i.e., imaginary frequencies); in other words, it refers to a “non-harmonic” potential-energy surface. Here, however, we define anharmonic contribution of some property X as the residual in excess of the harmonic approximation; Xah≡X−(Xlat+Xqh). Therefore, this specific definition is meaningless if the system does not have equilibrium lattice configuration at 0 K. For this reason, pyHMA checks forces on the first configuration to make sure the system has an equilibrium configuration (i.e., zero forces). In addition, when using pyHMA to measure anharmonic free energy using thermodynamic integration from 0 K (Sec. 3.3), only ground-state DFT must be used; using finite-temperature DFT (i.e., Fermi–Dirac smearing; ISMEAR=-1 and SIGMA=kBT), as often done with metals, cannot be used as the PES in this case is temperature-dependent, which is not accounted for in the integration. This contribution, however, can still be included using free-energy perturbation methods as described elsewhere [1]. On the other hand, for properties that do not require temperature integration (e.g., energy and pressure), pyHMA reads the electronic free-energy surface (F), rather than the ground-state energy (E0); hence, electronic contribution is accounted for.