Efficient wide-range calculation of free energies in solids and liquids using reversible-scaling molecular dynamics

Abstract

We elaborate a novel and efficient method to obtain multiphase Helmholtz free energies from molecular dynamics (MD) simulations over wide ranges of volume and temperature in materials that can be described by temperature-independent ion forces, with both higher accuracy and order-of-magnitude cost savings compared to direct thermodynamicintegration techniques. Our method leverages and significantly extends the technique of reversible scaling molecular dynamics (RSMD) proposed by de Koning et al. [Phys. Rev. Lett. 83, 3973 (1999)], which allows a free-energy difference in a given phase at constant volume to be calculated as a function of temperature from a single MD simulation. In mechanically stable solid phases, our approach carefully combines quasiharmonic lattice dynamics at low temperatures with an accurate and fully isolated RSMD simulation of the anharmonic vibrational free energy at high temperatures to produce a seamless free energy from zero temperature to above melt along constant-volume isochores. In the liquid, we combine a unique calculation of the free energy along a high-temperature reference isotherm with isochoric RSMD simulations from that temperature to below melt. In metastable solid phases that are mechanically unstable at low temperature, we use two-phase MD melt simulations together with the liquid free energy to obtain the solid free energy along the solidus melt line and then perform isochoric RSMD simulations to temperatures above and below that point. While our free-energy method is general, we have specifically adapted it here to the case of metals in which the ion forces are well described by model generalized pseudopotential theory (MGPT) multi-ion interatomic potentials, and additive electron-thermal free-energy contributions can be 2 included. Then using refined Ta6.8x MGPT potentials, we have converged total free energies and their components to very high and unprecedented sub−milli-Rydberg (mRy) numerical accuracy in the stable-bcc, liquid, and metastable-fcc phases of tantalum for volumes ranging from up to 26% expansion to nearly two-fold compression and for temperatures to 25,000 K. In turn, we have successfully used the free energies so obtained to calculate physically accurate thermodynamic properties and gain new insight into their behavior, including sensitive thermodynamic derivatives, bcc and fcc melt curves, and a multiphase equation of state for Ta over the same temperature range and for pressures as high as 600 GPa. We show that the anharmonic free-energy component in the bcc solid, although only 1-5 mRy in magnitude for Ta, can have a significant (15-20%) effect on thermal expansivity, the Grüneisen parameter, and melt temperatures. We further show that the electron-thermal free-energy component can similarly impact the specific heat and thermal expansivity in both the solid and the liquid, while only minimally affecting (to ≤ 3%) the bcc and fcc melt curves.

Accurate and precise ab initio anharmonic free-energy calculations for metallic crystals: Application to hcp Fe at high temperature and pressure

“…The integration path starts from T = 0 K, where the system behaves harmonically (hence A ah 0 = 0), and proceeds to increasing temperatures while recording averages needed to evaluate the FE at the temperature of interest. The method has been applied to obtain the anharmonic FE of several systems [58][59][60] . An appealing feature of a T -path is that integration to each point along the path yields the FE at that state.…”

Polymorphism and melt in high-pressure tantalum. II. Orthorhombic phases

“…The three theoretical melt curves displayed in Fig. 1 derive from the accurate MGPT free energies of Moriarty and Haskins [14], from the QMD, Z-melt simulations of Burakovsky et al [12], and from the DFT-based thermodynamic-perturbation-theory (TPT) calculations of Taioli et al [25]. The MGPT and QMD melt curves closely agree with each other up to pressures of ~ 250 GPa, as well as with the initial low-pressure melting slope obtained from isobaric expansion data [5].…”

“…MGPT calculations considered here we use the same high-quality Ta6.8x multi-ion potentials employed in Paper I and in Ref. [14]. We find, however, that the simplified Z-…”

“…The long-standing controversy in bcc transition metals between "flat" static melt curves obtained from early laser-heated diamond-anvil-cell (DAC) measurements [1][2][3] and "steep" dynamic melt curves inferred from dynamic isobaric [4][5][6][7] and shock [8][9][10][11] measurements has renewed research interest in the high-temperature (T), high-pressure (P) phase diagrams of these materials, with tantalum (Ta) emerging as an important specific example that has attracted considerable recent attention [12][13][14][15][16][17][18][19][20][21]. In 2010, Burakovsky et al [12] used small-cell (< 150 atoms) quantum molecular dynamics (QMD) melt simulations, based on first-principles density functional theory (DFT) [22] and Z-melt methodology [23], to predict a hexagonal omega (hex-w) phase in Ta at high temperature above 70 GPa.…”

“…In 2010, Burakovsky et al [12] used small-cell (< 150 atoms) quantum molecular dynamics (QMD) melt simulations, based on first-principles density functional theory (DFT) [22] and Z-melt methodology [23], to predict a hexagonal omega (hex-w) phase in Ta at high temperature above 70 GPa. Shortly thereafter, in Paper I of this series [13], we examined possible high-T,P polymorphism in five cubic and hexagonal phases of Ta with complementary DFT-based MGPT (model generalized pseudopotential theory) multi-ion interatomic potentials [13,14], which allow accurate treatment of much larger system sizes in MGPT-MD simulations. We found significant melt size effects for all non-bcc phases studied, requiring a minimum of ~ 500 atoms in the solid to simulate a reliable bulk melt curve for these phases.…”