Solid-liquid interface free energies of pure bcc metals and B2 phases.

Solid-liquid interface free energies of pure bcc metals and B2 phases S. R. Wilson, K. G. S. H. Gunawardana, and M. I. Mendelev Citation: The Journal of Chemical Physics 142, 134705 (2015); doi: 10.1063/1.4916741 View online: http://dx.doi.org/10.1063/1.4916741 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/142/13?ver=pdfcov Published by the AIP Publishing Articles you may be interested in Effect of confinement on the solid-liquid coexistence of Lennard-Jones Fluid J. Chem. Phys. 139, 174706 (2013); 10.1063/1.4827397 Step free energies at faceted solid-liquid interfaces from equilibrium molecular dynamics simulations J. Chem. Phys. 137, 214108 (2012); 10.1063/1.4769381 Solid-liquid interface free energy in binary systems: Theory and atomistic calculations for the (110) Cu–Ag interface J. Chem. Phys. 131, 054702 (2009); 10.1063/1.3197005 Atomistic simulation of solid-liquid coexistence for molecular systems: Application to triazole and benzene J. Chem. Phys. 124, 164503 (2006); 10.1063/1.2188400 Free energy of solvation of simple ions: Molecular-dynamics study of solvation of Cl − and Na + in the ice/water interface J. Chem. Phys. 123, 034706 (2005); 10.1063/1.1953578

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 137.99.31.134 On: Tue, 12 May 2015 13:44:27

THE JOURNAL OF CHEMICAL PHYSICS 142, 134705 (2015)

Solid-liquid interface free energies of pure bcc metals and B2 phases S. R. Wilson, K. G. S. H. Gunawardana, and M. I. Mendeleva) Division of Materials Sciences and Engineering, Ames Laboratory, Ames, Iowa 50011, USA

(Received 8 January 2015; accepted 23 March 2015; published online 7 April 2015) The solid-liquid interface (SLI) free energy was determined from molecular dynamics (MD) simula¯ tion for several body centered cubic (bcc) metals and B2 metallic compounds (space group: Pm3m; prototype: CsCl). In order to include a bcc metal with a low melting temperature in our study, a semi-empirical potential was developed for Na. Two additional synthetic “Na” potentials were also developed to explore the effect of liquid structure and latent heat on the SLI free energy. The obtained MD data were compared with the empirical Turnbull, Laird, and Ewing relations. All three relations are found to predict the general trend observed in the MD data for bcc metals obtained within the present study. However, only the Laird and Ewing relations are able to predict the trend obtained within the sequence of “Na” potentials. The Laird relation provides the best prediction for our MD data and other MD data for bcc metals taken from the literature. Overall, the Laird relation also agrees well with our B2 data but requires a proportionality constant that is substantially different from the bcc case. It also fails to explain a considerable difference between the SLI free energies of some B2 phases which have nearly the same melting temperature. In contrast, this difference is satisfactorily described by the Ewing relation. Moreover, the Ewing relation obtained from the bcc dataset also provides a reasonable description of the B2 data. C 2015 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4916741]

I. INTRODUCTION

The work required to form a crystal nucleus from an undercooled liquid phase impedes solidification, sometimes permitting metals to be held in liquid form for long periods of time at temperatures several hundred degrees below the freezing point. The nucleation barrier responsible for this effect is determined by a balance between two competing factors: the solid-liquid interface (SLI) free energy γ and the difference in the bulk free energies.1 The study of homogeneous nucleation has a long history, spanning at least 100 yr. Nonetheless, progress in the field has been significantly handicapped by the notorious difficulties that must be overcome to obtain accurate experimental data. In particular, experimental determination of γ is extremely difficult, and most reported values depend implicitly upon a particular theoretical interpretation of the observed nucleation rate. Because γ plays a fundamental role in nucleation, and because it is very difficult to measure, empirical guidelines that provide a quick estimate of its value from the knowledge of basic bulk material properties would be very useful. Such empirical correlations could also guide in the development of a comprehensive theory of the solid-liquid interface. The most famous empirical relation for the SLI free energy was proposed by Turnbull2 based on an analysis of available experimental data. This relation takes the form γ = CT ∆Hm ρ2/3 s ,

(1)

a)Author to whom correspondence should be addressed. Electronic mail:

[email protected] 0021-9606/2015/142(13)/134705/12/$30.00

where ∆Hm is the latent heat of melting, ρ s is the atomic density in the bulk solid phase, and CT is the Turnbull coefficient. Turnbull noticed that separate consideration of metals and nonmetals produced different values of CT (0.45 and 0.32, respectively), but the limited extent and accuracy of his experimental database did not permit him to identify any further trends. In response to severe experimental limitations, nucleation studies have increasingly turned toward molecular dynamics (MD) simulations for guidance, and several factors have recently converged to provide an excellent opportunity for MD simulation to contribute toward an enhanced understanding of the relation between the SLI free energy and bulk material properties. First, significant progress has been achieved in the development of methods to determine γ from MD simulation3–5 and their accuracy is limited only by the available computational power, which has dramatically increased over the past few decades. Second, the presence of impurities, oxidation, and similar experimental difficulties is not an issue in MD, completely bypassing major sources of experimental error. Finally, it is well known that the quality of the semi-empirical potential can influence the SLI properties that are measured in the course of classical MD simulation, and this has been considered a significant source of error in values obtained from MD simulation. However, it has been recognized that the interatomic potential may be regarded as an additional “knob” that computational scientists may twist to sample a much larger set of materials than are observed in nature. Fundamental principles from statistical mechanics and thermodynamics do not depend upon the nature of the material studied, and so should apply equally well to both

142, 134705-1

© 2015 AIP Publishing LLC


134705-2

Wilson, Gunawardana, and Mendelev

J. Chem. Phys. 142, 134705 (2015)

“real” and synthetic materials. This can be extremely useful in discovering and testing empirical relations such as Eq. (1). For example, the analysis of data obtained using a variety of different pair and embedded atom method (EAM)6 potentials for pure face centered cubic (fcc) and body centered cubic (bcc) metals enabled computational determination of the Turnbull coefficients for both the fcc metals (CTfcc = 0.524) and bcc metals (CTbcc = 0.327). Although the Turnbull relation has received the most attention, a number of alternate relations exist. One such relation was recently advanced by Laird, based on the idea that “the solid–liquid free energy is primarily entropic in origin.”8 In this approximation, the SLI free energy of fcc metals should be well-described by a hard-sphere model, which leads to the following “rule of thumb”: γ = CL k BTm ρ2/3 s ,

(2)

where k B is the Boltzmann constant, Tm is the melting temperature, and CL is a constant. An earlier approach, due to Ewing,9 is based on the hypothesis that the SLI free energy may be considered as the sum of two contributions, one for each bulk phase. Ewing proposed that the term associated with the solid phase is proportional to the latent heat ∆Hm = Tm ∆Sm , essentially equivalent to the Turnbull hypothesis. However, to this term, he added a second contribution associated with the liquid phase and arising from enhanced ordering in the liquid near the interface, which he assumed should be related to the bulk liquid structure. Ewing obtained an explicit expression for the liquid structure contribution, but his derivation only considered the effect of ordering liquid atoms along the direction normal to the interface. Because the solid phase should induce systematic variations in the liquid structure observed in the tangent plane as well as along the interface normal, a recent proposal replaced the second term in the original Ewing relation by a term proportional to the liquid excess entropy,10 2/3 γ = CT′ ρ2/3 s Tm ∆Sm + CO ρ s Tm SO ,

(3)

where both CT′ and CO are constants depending on the crystalline lattice, and SO is the liquid excess entropy (with respect to the ideal gas), which in the pair approximation is given by ∞ SO = 2π k B ρl

{g(r) ln[g(r)] − [g(r) − 1]} r 2dr.

metals. For bcc metals, only the coefficient that appears in the Turnbull relation, Eq. (1), is known. Thus, the determination of these coefficients is currently a very important task. MD simulation arguably offers the easiest way to determine these coefficients, especially if we take into account that an accurate determination of these constants requires sampling materials over a very wide range of melting temperatures, and high melting temperatures represent an additional obstacle for an experimental study. Establishing an empirical relation between the SLI free energy and bulk material properties is even more important for binary alloys, because additional computational challenges in this case can significantly slow down determination of the SLI free energy from MD simulation. One such challenge occurs in the case of solid solutions, for which SLI free energy measurements require preliminary determination of the equilibrium compositions of the bulk solid and liquid phases.11 Stoichiometric congruently melting compounds that do not require these preliminary phase diagram calculations would seem to be easier to consider, but extremely low interface kinetic coefficients in these systems can nonetheless increase the required computational times by more than an order of magnitude.12 The goal of the present study was to determine whether any of the empirical relations described above can be applied to stoichiometric congruently melting B2 compounds. The rest of the paper is organized as follows. First, we will describe the material systems used in the present study. Second, we will describe how we determined the melting temperature for these systems, which is crucial to apply the capillary fluctuation method (CFM)13 used in the present work to determine SLI free energies. Third, we will present calculations of the liquid excess entropy that appears in Eq. (3). Fourth, we will describe the MD simulations performed to determine the SLI free energy. Finally, we will test all three of these empirical relations with respect to the SLI free energy data obtained from MD simulation.

(4)

0

Here, ρl is the liquid density and g(r) is the liquid pair correlation function (PCF). It was shown in Ref. 10 that for the fcc metals, this modification to the Ewing relation provides the best agreement with MD simulation data. As mentioned above, empirical relations can be very useful to obtain a quick estimate of the value of the SLI free energy. For example, to employ Eq. (4), all we need to know is the atomic density, the latent heat, and the liquid pair correlation functions. All of these quantities can be obtained with relatively little effort, either from experiment or from MD simulation. However, the coefficients CT′ and CO depend on the crystalline structure and need to be known in advance. To our best knowledge, the coefficients that appear in the three empirical relations given above are known today only for fcc

II. MATERIAL SYSTEMS

Methods to determine the SLI free energy require using rather large simulation cell sizes and long computational times (see Sec. V), so that using semi-empirical potentials is the only option. The most popular methods to describe interatomic interaction in metallic systems are based on EAM6 and FinnisSinclair (FS)14 potentials. The total potential energy in this formalism takes the following form: U=

N −1 

N 

ϕt i t j (r i j ) +

i=1 j=i+1

N 

Φ t i (ρi ),

(5)

i=1

where t i is the elemental type of atom i, N is the number of atoms in the system, r i, j is the separation between atoms i and j, ϕt i t j (r) is the pairwise potential, and Φ t i (ρ) is the embedding energy function, and  ρi = ψ t i t j (r i j ), (6) j

where ψ t i t j (r) are density functions. If the density functions only depend on the type of atom j, i.e., ψi j (r) = ψ j (r), Eqs. (5)


134705-3


J. Chem. Phys. 142, 134705 (2015)

and (6) reduce to the embedded atom method. If the density functions depend on the types of both atoms, and ψi j (r) = ψ j i (r), Eqs. (5) and (6) reduce to the Finnis-Sinclair scheme. In the present study, we focused on alloys for which the B2 phase is either stable or at least metastable in an equimolar alloy. This requirement is satisfied by Cu–Zr and Ni–Al alloys, for which B2 is the stable phase. In addition, metastable B2 appears for a short time during solidification in Ni-Zr, before the stable B33 phase appears,15 and so B2 Ni–Zr was also included in this study. The interatomic potentials for these alloys were taken from Ref. 16 (Cu–Zr), Ref. 17 (Ni–Al 2002), Ref. 18 (Ni–Al 2009), Ref. 19 (Ni–Zr 2012), and Ref. 12 (Ni–Zr 2014). Note that these Cu–Zr and Ni–Zr potentials are of the FS type, while the Ni-Al potentials are of the EAM type. If atomic type is ignored, then the B2 lattice reduces to the bcc structure. Therefore, it is of interest to compare SLI properties of the B2 compounds and pure bcc metals. It has been well established that SLI kinetic coefficients of B2 compounds can be orders of magnitude smaller than those of pure metals (e.g., compare the data presented in Refs. 20 and 21 and in Refs. 22 and 23). However, to our best knowledge, no such comparison has been made for the SLI free energy. Among the pure elements that constitute the alloys studied in this work, only Zr solidifies into the bcc lattice and could be included for comparison. In order to broaden the proposed comparison, we therefore included several additional bcc metals in our study. The three empirical relations described in Sec. I are best compared by sampling materials across a broad range of melting temperatures, and so we included bcc metals with both high and low melting temperatures in our study. Tungsten is an obvious choice as a bcc metal with high melting temperature, while sodium is a reasonable choice for a bcc metal with low melting temperature. The EAM potential for W developed by G. J. Ackland was taken from Ref. 24. An EAM potential for Na (Na1) was developed in the present study following the same fitting procedure which was used for potential #5 in Ref. 25. This potential is available in Refs. 24 and 26. Some properties obtained using this potential are provided in Table I. For the present study, it is important that this potential

FIG. 1. Pair correlation function of liquid Na at T = 378 K. The X-ray data have been taken from Ref. 60 and were used to develop the potential Na1.

satisfactorily reproduces both the melting point data and the liquid structure (see Fig. 1). To test the Ewing relation, it is convenient to consider materials whose bulk properties are nearly identical but have different liquid structures and latent heats (see the discussion in Ref. 10). In order to develop potentials that meet this criterion, we began with the Na1 potential and followed the procedure described in Ref. 27, except that in this study, we tried to obtain less ordered liquid structure and higher latent heat. Examination of the properties provided by two additional potentials, Na2 and Na3, listed in Table I, and by the PCFs plotted in Fig. 1, shows that these two potentials are useful additions to our portfolio. To compare our results with the previous SLI free energy measurements, we included MD simulation data for three bcc metals from Ref. 7 (Fe: ABCH and MH(SA)2 potentials) and Ref. 20 (Mo and V). In the present study, we only determined the excess entropy value for these systems; all other values are taken from Refs. 7 and 20.

TABLE I. Physical properties of Na.a Property a (Å), bcc at T = 0 K Ecoh (eV/atom), bcc C11 (GPa), bcc at T = 0 K C12 (GPa), bcc at T = 0 K C44 (GPa), bcc at T = 0 K Evf (eV), bcc at T = 0 K ∆Ebcc→ fcc (eV/atom) at T = 0 K ∆Ebcc→ hcp (eV/atom) Tm (K) ∆Hm (eV/atom) ∆Vm/V(%)

Target value

Na1

Na2

Na3

4.23b

4.23 1.11 8.3 6.8 5.8 0.43 0.009 0.006 366 0.028 2.3

4.23 1.09 8.2 6.8 5.8 0.41 0.009 0.006 377 0.032 1.4

4.23 1.08 7.8 6.8 5.8 0.40 0.010 0.006 376 0.035 1.8

1.13b 8.2b 6.8b 5.8b 0.42b 0.009b 371c 0.027c

is the lattice parameter, E coh is the cohesive energy, E vf is the vacancy formation energy, C i j are elastic constants, ∆E bcc→ β is the difference in energy between the bcc and β phases, Tm is the melting temperature, and ∆H m and ∆Vm are the latent heat and change in the atomic volume upon melting, respectively. b These data have been taken from Ref. 58. c These data have been taken from Ref. 59.

a In this table, a


134705-4


III. DETERMINATION OF THE MELTING TEMPERATURE

The capillary fluctuation method used in the present study to determine SLI free energies requires prior knowledge of the melting temperature. In all cases, the melting temperature should be determined with sufficient accuracy to avoid complete solidification or melting of the domain while collecting SLI profiles (see below). We can estimate the necessary accuracy as follows. The simulations performed to determine SLI free energies typically lasted around 20 ns for bcc metals and involved domains whose cross sections were about 50 nm long in the direction normal to the SLI. Since at the start of the simulation half of the simulation cell was liquid and the other half was solid, we could allow the SLI to drift by no more than about 12.5 nm during the entire MD simulation. If we assume that a typical value of the kinetic coefficient for a bcc metal that melts around 2000 K is 0.6 m/(s K) (see Ref. 28), we find that the inaccuracy in the melting temperature should not exceed 1 K or 0.05% of the melting temperature. According to Frenkel-Wilson theory,29,30 the kinetic coefficient for different bcc metals scales approximately as 1/Tm, in which case we can conclude that for any single-component bcc system, the inaccuracy in the determination of the melting temperature should be no larger than 0.0005 Tm. As mentioned above, it is known that the SLI kinetic coefficient of B2 compounds is at least one order of magnitude smaller than that of single-component systems. Therefore, the acceptable inaccuracy in the determination of the melting temperature for B2 compounds would appear to be larger by the same factor. However, determination of the SLI free energy in these systems requires longer MD simulation times, which decreases the acceptable value of inaccuracy. The method to determine melting temperatures in this work was similar to that used in Ref. 27. First, the liquid density and the lattice parameter of the solid phase corresponding to zero pressure were determined as functions of temperature (for details, see Ref. 27). These simulations also provide an approximate value of the melting temperature because the solid phase can be overheated in MD simulation by approximately 0.2 Tm. Next, a model of the crystal phase containing ∼30 000 atoms was created, such that the domain size in the z-direction was roughly 5 times longer than in the x- and y-directions. During all subsequent MD evolutions, periodic boundary conditions were only enforced along the x- and y-directions, creating a pair of free surfaces along the z-direction. Starting from this configuration, layers ∼20 Å thick along both the top and bottom of the model were annealed for 1000 MD steps at up to 1.5 times the estimated value of the melting temperature, followed by rescaling the simulation cell size in the x- and ydirections according to the equilibrium lattice parameter at the target temperature and an equilibration of the entire model in the NVT ensemble, thus producing a model in which the crystal phase was bounded by two thin layers of liquid. Molecular dynamics simulation at constant temperature and cross section then led to melting (provided that the target temperature is above the melting temperature). To avoid any temperature gradients that could develop as a result of interface migration, a “layered thermostat” was applied. Once the thickness of

J. Chem. Phys. 142, 134705 (2015)

the crystal phase in the z-direction had been reduced to less than 50 Å, the simulation was stopped. One of the resulting models was then rescaled to simulate solidification at lower temperatures. During MD simulation, the change in the total energy as a function of time was recorded. The slope of this dependence is related to the SLI velocity as follows: V=

dE 1 , 2Sx y ∆Hm/Ω dt

(7)

where Sx y is the cross-sectional area of the simulation cell perpendicular to the z-direction and Ω is the average volume per atom in the crystal phase. The SLI velocity can be expressed as V = µ(T − Tm ),

(8)

where µ is the kinetic coefficient, so that Eq. (7) may be written as dE = A(T − Tm ), (9) dt where A is a constant. This equation could be used directly to determine the melting temperature. However, the dependence of dE/dt on T − Tm is linear only near the melting temperature. If we were to apply the least-squares method to determine the melting temperature from Eq. (9), the largest contributions would be from simulations corresponding to the largest supercooling and overheating, where the dependence of dE/dt on T − Tm may not be linear. Therefore, it is better to find the melting temperature from minimization of the quantity 2   dEi /dt −A , Ti − Tm i where the summation is over all simulation data. This ensures that all simulation data points contribute equally to the determination of the melting temperature. An example is shown in Fig. 2, from which it can be seen that the observed dependence of dE/dt on T − Tm is linear near the melting temperature. The solidification simulations used to determine melting temperatures in this work also demonstrate that the employed potentials are indeed suitable for simulations that focus on SLI properties. The study of SLIs requires elevated temperatures, and under these conditions, undesirable phase behavior can often be observed in MD simulation, such as the nucleation of an unexpected solid phase at the SLI. The fact that a potential correctly reproduces the crystal phase energies at T = 0 K (a constraint that is frequently used during potential development) cannot ensure that the correct phases will remain stable at high temperatures. For example, the formation of the bcc phase was observed during MD simulations of SLIs in hcp Mg in Ref. 31 when the potential from Ref. 32 was employed. Proper solidification demonstrates that the employed potential provides that the solid phase to be studied is at least metastable. Testing whether a given solid phase grows as expected from its equimolar melt for a given potential is even more important when considering binary alloys, since there are additional ways in which undesirable behavior can occur. For instance, the composition of the as-grown solid phase can


134705-5


J. Chem. Phys. 142, 134705 (2015)

FIG. 2. Data collected from solidification/melting simulations used in the determination of the melting temperature for Zr.

differ from that of the liquid phase, or the solid phase can be disordered in contrast to the expected ordered phase. Our MD simulations showed that all alloy potentials considered here indeed lead to proper B2 phase growth during solidification. More precisely, this means that the atomic fraction of both components of the growing crystal phase was equal to 0.5, and that the solidified structure was overall ordered. The growing phases did contain some fraction of anti-site defects. This fraction was very small in the case of Cu–Zr and Ni–Zr alloys, and it was larger for the Ni–Al alloys, which is consistent with the data published in Ref. 22. The fact that the growing crystal phase exhibits the same composition as the liquid phase from which it grows means that no partitioning occurs during solidification, and therefore, no preliminary determination of the equilibrium compositions of the solid and liquid phases is necessary in order to measure SLI free energies as described in Sec. V.

In order to determine whether Eq. (10) is a good approximation for the excess entropy, we performed an independent calculation for Zr following perturbation theory (PT), as introduced by Weeks, Chandler, and Andersen33,34 and refined by Ree et al.35 in conjunction with a hard sphere reference system. A detailed description of the method followed in the present study can be found in Refs. 36 and 37. Interestingly, the value obtained from perturbation theory almost coincides with the value of SO obtained using Eq. (4) (see Table II). A generalization of Eq. (10) for binary alloys proposed in Refs. 38 and 39 takes the following form: SO = 2π k B ρl

2  2 

L

xα x β

α=1 β=1

− [gα β (r) − 1] r 2dr.

gα β (r) ln[gα β (r)]

0

(11)

IV. DETERMINATION OF LIQUID EXCESS ENTROPY

The liquid excess entropy for pure metals was determined in this work in two ways. The first method was based on Eq. (4). While the integration domain that appears in this equation extends to infinity, in practice, PCFs can be determined only up to the half of the smallest simulation cell size L. Therefore, Eq. (4) must be truncated during numerical evaluation to account for this limitation as follows: L SO = 2π k B ρl

{g(r) ln[g(r)] − [g(r) − 1]} r 2dr.

(10)

0

We investigated how sensitively the calculation of SO using Eq. (10) depends on the upper limit L. To do so, we created a liquid Zr model containing 30 000 atoms and equilibrated this model at the melting temperature in the NVT ensemble for 40 ps. The PCF was then obtained by averaging over the next 40 ps. This calculation revealed that although noticeable oscillations in the PCF persist up to distances of at least 40 Å, the integral in Eq. (10) converges much faster, achieving an accuracy of 0.01% even for L = 20 Å, as can be seen in Fig. 3.

FIG. 3. Pair correlation function of liquid Zr at Tm used to determine the excess entropy, SO . The bottom inset shows that the oscillations can be seen up to 40 Å in our simulation. However, the top inset demonstrates that accounting for these oscillations does not considerably change the value of SO if the cutoff radius, L, in Eq. (10) is larger than 20 Å.


134705-6


J. Chem. Phys. 142, 134705 (2015)

TABLE II. Material properties calculated with different semi-empirical potentials.

System

Tm (K)

d (atom/Å3)

∆H m (eV/atom)

W Zr Na1 Na2 Na3 Fe ABCHa Fe (MH(SA)2)a Mob Vb Ni-Al 2002 Ni-Al 2009 Cu-Zr Ni-Zr 2012 Ni-Zr 2014

4107 2108 367 376 378 2358 1772 3789 3031 1528 1821 1326 1404 1369

0.0571 0.0412 0.0247 0.0251 0.0249 0.0784 0.0800 0.0593 0.0640 0.0787 0.0811 0.0568 0.0613 0.0596

0.321 0.180 0.0282 0.0317 0.0347 0.218 0.162 0.379 0.332 0.233 0.302 0.198 0.186 0.100

a The b The

Eq. (4)

∆H m k BT

0.91 0.99 0.89 0.98 1.07 1.07 1.06 1.16 1.27 1.77 1.92 1.73 1.54 0.85

−SO (meV/(K atom)) 0.315 0.249 0.286 0.256 0.231 0.319 0.301 0.316 0.322 0.317 0.277 0.306 0.296 0.337

PT −SO (meV/(K atom))

0.249

0.308 0.265 0.306 0.330 0.353

values for SO were obtained within the present study; the rest of the data have been taken from Ref. 7. values for SO were obtained within the present study; the rest of the data have been taken from Ref. 20.

The excess entropy was determined for several binary systems using perturbation theory and compared to the values of SO obtained from Eq. (11) in Table II, following the same procedure as was employed for Zr. It can be seen from this table that Eq. (11) provides a reasonable approximation to the theoretical values. However, all subsequent analyses in this work will use the liquid alloy excess entropy values obtained from perturbation theory.

V. DETERMINATION OF SLI FREE ENERGIES

SLI free energies γ(n) were estimated for all materials considered in this work within the framework of the capillary fluctuation method, which is based on capillary wave theory40–42 and was introduced in Ref. 13. Since its introduction, CFM has been applied to compute SLI free energies from MD simulation in a wide variety of material systems, including binary11,43–45 and ternary46 alloys. Although a number of works have applied CFM to disordered alloys, to date few authors have considered the special case of stoichiometric compounds, with the exceptions of a binary B1 hard-sphere system,47 and remarkably succinonitrile,48 which is an organic material that solidifies into a bcc superlattice. An estimate for γ and its dependence on crystal inclination n is obtained in CFM from measurements of the stiffness γ + γ ′′. The extraction of γ from the stiffness is accomplished by expanding γ(n) in terms of the appropriate symmetryadapted harmonics, truncating to the desired level of accuracy, and measuring stiffness values for sufficiently many different inclinations to determine all coefficients that appear in this expansion. B2 and bcc are both cubic structures, so that following Ref. 44, we may in both cases write γ(n) to first order in terms of the cubic harmonics as  ( ) 3 4 γ n x , n y , n z = γ0 1 + ε 1 n − i i 5 (  ) 17 + ε2 3 n4i + 66n2x n2y n2z − , (12) i 7

where γ0, ε 1, and ε 2 are expansion coefficients. To this order, only three stiffness values, measured for three different inclinations, are required to determine the three expansion coefficients. To increase the accuracy of our results, we instead routinely determined six stiffness values (for all pairs n, t listed in Fig. 4) and performed a least-squares fit to extract γ0, ε 1, and ε 2. The stiffness is determined in CFM by monitoring fluctuations ∆h(x) = h(x) − h0 in the interface profile h(x) from the instantaneous average interface position h0 of a nominally flat, equilibrated SLI. Long-wavelength amplitudes of the Fourier transform Ak of the interface profile h(x) are known to obey the relation

2 k BTm , (13) Ak = bW (γ + γ ′′) k 2 where brackets denote a statistical average over a sufficiently large set of independent interface profile measurements. In Eq. (13), W is the width of the domain parallel to the x-axis, b is the thickness in the y-direction (chosen as small as possible to yield an effectively two-dimensional simulation domain), k is the wavenumber, and primes indicate derivatives of γ with respect to inclination n along the x-direction. The interface normal is constrained by periodic boundary conditions parallel to the z-axis. Values of the stiffness γ + γ ′′ are obtained directly from Eq. (13). The set of wavenumbers accessible to a given stiffness measurement is determined by the simulation cell size in the xdirection, with the longest wavelength specified by the fundamental harmonic, k1 = 2π/W . Because the set of wavenumbers for which the long-wavelength approximation holds true is unknown prior to measurement, we have chosen domain sizes comparable to those used in the first application of CFM,13 i.e., roughly 250 Å in the x-direction, adjusted to satisfy periodic boundary conditions. Equation (13) is only valid in the long-wavelength, small wavenumber limit, so that only Fourier amplitudes corresponding to the lowest wavenumbers are included in the slope determination. In the present work, all stiffness values have been computed from Eq. (13) on the basis


134705-7


J. Chem. Phys. 142, 134705 (2015)

  FIG. 4. The dependence of A2k on wavelength k for Na1 and Ni-Al 2009. The strong linearity observed in  these log-log plots indicate that A2k scales as the inverse square of k, as predicted by Eq. (13).

of at least five values of k. The validity of the long-wavelength approximation is assumed to be reflected in the extent to which linearity is observed with respect to k 2 (discussed in more detail below). The statistical average on the LHS of Eq. (13) is to be taken over the distribution of thermodynamic fluctuations observed in a well-equilibrated system. As was noted above, highly accurate determination of the melting temperature prior to CFM computations is essential to obtain the required equilibration, as well as to avoid significant solidification or melting during data collection. Interfaces equilibrated at a given melting temperature were prepared as follows. A perfect crystal was constructed in a fully periodic domain roughly 250 Å in the x-direction and 400 Å in the z-direction, using a prior computation of the equilibrium lattice parameter at the melting temperature in question. The top half of this domain was then annealed in the NVT ensemble for 25 000 time steps at 1.5-2.0 times the melting temperature. Next, periodicity was broken along the z-direction, all velocities rescaled to the melting temperature, and the bottom half of the domain equilibrated for 30 000 time steps in the NVT ensemble while holding all atoms in the top half of the domain fixed. This procedure created a solid-liquid interface sandwiched between a solid-vapor interface and a liquid-vapor interface, guaranteeing zero pressure throughout the domain without the assistance of a barostat. The resulting system was then permitted to equilibrate in the NVT ensemble. CFM data collection was initiated once both (a) the thermodynamic variables of the entire system and (b) the observed mean-square profile fluctuation,  ∆h2 =

W

(h(x) − h0)2dx,

(14)

0

were observed to attain steady-state values. Due to the long run times and large domains required for this measurement, all simulations were performed in parallel in the context of the classical molecular dynamics software package LAMMPS,49 utilizing a MD time step of 2 fs. In the more general case of disordered alloys, an initial determination of the phase diagram (typically from Monte Carlo calculations in the grand canonical ensemble) must be performed to determine the equilibrium bulk solid and liquid

compositions for the alloys of interest.11 Monte Carlo methods are also used in multicomponent systems to generate the equilibrium segregation profile. In the case of stoichiometric compounds, the bulk solid and liquid compositions are known in advance, so that prior phase diagram determination is not needed. However, segregation to the SLI can still occur in the stoichiometric case, and in fact, we observe segregation in a number of our simulations. For small domain sizes, significant segregation could conceivably deplete or enrich the bulk phases on either side, resulting in a shift away from stoichiometry and an erroneous measurement. However, no significant deviations from stoichiometry in the bulk phases are observed for any of the simulations included in our computations, as can be seen in the example shown in Fig. 5. This indicates that the domain sizes we have selected are large enough, and the segregation is sufficiently weak, to permit the adjacent solid and liquid phases to behave to a good approximation as infinite reservoirs, with respect to the interface that divides them. Determination of the interface profile h(x), and the measurement of fluctuations in this profile, requires an automated method to distinguish between solid and liquid phases. The authors of Ref. 13 introduced an order parameter that measures the extent to which the local atomic configuration deviates from the known perfect lattice sites. An alternate order parameter was proposed in the context of CFM in Ref. 50, based on a known set of strong Bragg reflections in reciprocal space. Subsequent work51 with the hard-sphere system demonstrated that CFM results do not depend on this choice, as long as no intervening solid phases appear in the interface region. In this work, we closely followed the original formulation, also implemented in Refs. 7 and 52. For both B2 and bcc systems, an order parameter ϕi was computed for the ith atom in the system by minimizing the quantity ϕi =

)2 1 J ( r j − r i,(0)j , j=1 J

(15)

where r i,(0)j is the perfect lattice position, centered on atom i, closest to the jth neighbor of atom i, irrespective of the atomic species of either atom i or j. The sum is taken over the J nearest neighbors of atom i. The interface profile h(x) is then determined by binning atoms in a grid composed of cells whose


134705-8


J. Chem. Phys. 142, 134705 (2015)

FIG. 5. Concentration profiles of atomic species in the directions normal to several SLIs in the Ni-Al 2009 system.

linear dimensions are each some factor larger than the lattice parameter, calculating the average value ϕ(x, ¯ z) of ϕi in each of these bins, and then interpolating the resulting averages along the z-direction to locate a contour h(x) midway between average values of ϕ¯ that are observed in the bulk solid and liquid phases. Bin sizes in the present work were chosen to be two lattice parameters wide in the x-direction, 1.5 lattice parameters long in the z-direction, and spanned the entire thickness of the domain in the y-direction. Note that a number of authors have studied the dependence of CFM measurements on bin size (for a review, see Ref. 53); in particular, it was found in Ref. 54 that CFM values obtained were insensitive to this choice, although for other purposes, greater resolution appears optimal. In Fig. 4, we present log-log plots for Na1

and Ni–Al 2009 from which the observed dependence of A2k on k may

be extracted. In all cases, the observed slopes indicate that A2k is proportional to the inverse square of k, in agreement with Eq. (13). We found that the number of interface profile samples required in the average over Ak 2 to obtain a good approxima-

tion to the steady-state value A2k varies significantly from one material to the next. Accurate determination of the stiffness on the basis of Eq. (13) requires sufficiently dense sampling over the distribution of probable interface profiles, and so data need to be collected over at least several periods of any observed oscillation in ∆h2. Oscillation periods appear to increase with decreasing interface kinetic coefficient, and because interface kinetic coefficients for stoichiometric compounds are typically at least an order of magnitude lower than those in pure materials,23 correspondingly longer simulation times are required to determine the interface stiffness via CFM for these materials. Several of the B2 compounds considered here admit very low interface kinetic coefficient, so that due to limited resources, the uncertainty in the values presented in Table III is in these cases larger, despite simulation times in some cases in excess of 200 ns. Values obtained for the cubic harmonic expansion coefficients γ0, ε 1, and ε 2 following the procedure outlined above for five B2 systems (Cu–Zr, two Ni–Zr potentials, and two Ni–Al potentials) and five bcc systems (Zr, W, Na1, Na2, and

TABLE III. Solid-liquid interface properties calculated with different semi-empirical potentials. System

γ 0 (mJ/m2)

ε1

ε2

γ 100−γ 110 (%) 2γ 0

γ 100−γ 111 (%) 2γ 0

W Zr Na1 Na2 Na3 Fe ABCHa Fe (MH(SA)2)a Mob Vb Ni-Al 2002 Ni-Al 2009 Cu-Zr Ni-Zr 2012 Ni-Zr 2014

298 ± 5 128 ± 4 16.1 ± 0.4 17.1 ± 0.4 16.9 ± 0.5 206 ± 10 175 ± 11 249 ± 30 252 ± 10 229 ± 14 281 ± 8 201 ± 20 196 ± 12 129 ± 7

0.043 ± 0.006 0.03 ± 0.01 0.050 ± 0.006 0.045 ± 0.005 0.044 ± 0.008 0.015 0.032 0.03 ± 0.04 0.031 ± 0.009 0.00 ± 0.02 0.04 ± 0.01 0.03 ± 0.03 0.07 ± 0.02 −0.00±0.02

−0.0003 ± 0.0009 −0.002 ± 0.002 0.002 ± 0.001 0.001 ± 0.001 −0.002 ± 0.002 0.0003 0.003 0.006 ± 0.005 0.002 ± 0.002 −0.001 ± 0.003 0.002 ± 0.002 0.002 ± 0.005 −0.005 ± 0.003 0.002 ± 0.003

1.1 0.7 1.4 1.1 1.0 0.4 1.0 1.2 0.9 0.0 2.2 1.2 1.8 −0.2

1.4 1.0 1.5 1.4 1.4 0.5 1.0 0.9 1.0 0.1 2.3 1.4 2.8 −0.3

a These b These

data have been taken from Ref. 7. data have been taken from Ref. 20.


134705-9


J. Chem. Phys. 142, 134705 (2015)

FIG. 6. Comparison of the empirical relations with MD simulation data. The dashed and dotted-dashed straight lines were created using the coefficients obtained for particular empirical relations. Note that only the bcc and B2 data obtained in the present work (shown as red circles) were used to obtain these coefficients for the dashed and dotteddashed lines, respectively. The coefficients given in each plot were computed using only the bcc data obtained in the present study. If the Turnbull relation were precise, the data on the bottom right plot would fall on a horizontal line. If the Ewing relation were precise, the data on the bottom right plot would fall on a straight line.

Na3) are listed in Table III, alongside values computed for four additional bcc systems (Mo, V, and two Fe potentials) taken from the literature.7,20 While this study focused primarily on the average value of the SLI free energy, as measured by γ0, and which we will discuss extensively in Sec. VI, the CFM parameters ε 1 and ε 2 provide a measurement of the inclination dependence of the SLI free energy in the obtained dataset, which we briefly discuss here. Overall, the degree of anisotropy we observed in these systems is rather small, as can be seen by the ratios provided in the last two columns in Table III. The values reported here for bcc metals are consistent with the previous measurements, i.e., around 1%, and in all of these systems, we observe that γ100 > γ110 and γ100 > γ111, in agreement with three of the four previously published values. The same ordering is likewise observed in four of the five B2 systems studied.

The exception is Ni–Zr 2014, for which γ111 > γ110 > γ100; however, as Ni–Zr 2014 is nearly the most isotropic system included in this dataset, uncertainty in the measurements may easily influence this ordering. Figures S1 and S2 in Ref. 26 visualize this data by plotting the SLI free energy and stiffness along selected cross-sectional slices of direction space. It is important to note that in both figures, we have normalized γ(n x , n y , nz ) by γ0, and that in the case of the SLI free energy (Fig. S126), we have magnified deviations from isotropy by a factor of 5. The general trend observed in our data supports the intuitive idea that the inclination dependence is largely determined by the crystal structure, especially since the bulk properties of the bcc systems studied in this work vary widely. However, a notable spread in the magnitude of the anisotropy is clearly evident among the B2 values; some are significantly more isotropic, and others more anisotropic, than bcc systems.


134705-10


It is important to note that this variation may be strongly influenced by the relatively larger uncertainties in the B2 measurements. VI. DISCUSSION

In this work, we have computed SLI free energies for five bcc and five B2 phases. Our purpose has been to investigate empirical correlations between the SLI free energy and bulk properties such as crystal structure, liquid structure, latent heat, and melting temperature, all of which figure prominently in the Turnbull, Laird, and Ewing relations (Eqs. (1)–(3), respectively). The principal results of our investigation, which we discuss in this section, have been summarized in Figs. 6 and 7. We will first focus our discussion on the bcc data plotted in these figures, before discussing the full dataset. Of the five

J. Chem. Phys. 142, 134705 (2015)

bcc phases whose SLI free energies were computed in this work, three were modeled using interatomic potentials that reproduce properties of “real” materials (W, Zr, and Na1) and two are “synthetic” materials (Na2 and Na3) that are closely related to Na1. Figure 6 shows that overall both the Turnbull and Laird relations provide rather accurate predictions of the SLI free energy for W, Zr, and Na1, in spite of tremendous differences in melting temperatures (see Table II) and other bulk properties (for instance, the bulk moduli for Na1 and W are 7.3 GPa and 310 GPa, respectively). However, the Turnbull relation cannot account for the SLI free energies obtained for the synthetic “Na” potentials, as can be clearly seen in Fig. 7. Indeed, although the latent heat varies across this series by 23% (see Table II), we found no systematic variation in the SLI free energy (see Table III and Fig. 7). In fact, the SLI free energy computed for all three “Na” potentials is nearly the

FIG. 7. Comparison of the empirical relations with MD simulation data for the “Na” potentials. The dashed straight lines were created using the coefficients obtained for particular empirical relations. The data for “Al” potentials from Ref. 10 are also shown in the bottom right plot. If the Turnbull relation were precise, the data on the bottom right plot would fall on a horizontal line. If the Ewing relation were precise, the data on the bottom right plot would fall on a straight line.


134705-11


same, differing by at most 6%. The Laird relation is seen to provide an excellent agreement with the MD data for all “Na” potentials. The fact that the Turnbull relation fails to describe the MD data for the “Na” potentials can also be clearly seen in the bottom right plot in Fig. 7. If the Turnbull relation were precise, the data on this plot would fall on a horizontal line, which is clearly not the case. On the other hand, if the Ewing relation were precise, the data on this plot would fall on a straight line, and in fact, we do observe a clear linear trend in our data. Interestingly, the Ewing relation is also almost exact for the series of “Al” potentials considered in Ref. 10 (which were developed in the same manner as the “Na” potentials in the present study). SLI free energies for four additional bcc metals were previously reported in Refs. 7 and 20, and in Ref. 20, these values were considered in the context of the Turnbull relation. In Fig. 6, these values have been plotted alongside our own. From the top left plot in Fig. 6, it can be seen that the bcc data reported in the present study produce a much stronger linear trend than do the previously reported values, indicating much better agreement with the Turnbull relation. Moreover, the Turnbull coefficient obtained in the present study (CT = 0.39) is considerably higher than that quoted in Ref. 20 (CT = 0.29) and compares very favorably to the prediction CTbcc/CTfcc = 0.82 derived from the Spaepen-MeyerThompson (SMT) model,55–57 assuming CT(fcc) = 0.49 as reported in Ref. 10, although the SMT predictions for CTbcc and CTfcc both differ from the computed values by a factor of about 1.8. Nonetheless, the top right plot in Fig. 6 demonstrates that the combined bcc dataset can be described on the whole rather well using the Laird relation. Finally, note that the Ewing relation provides modestly better agreement with the data from Refs. 7 and 20 than does the Turnbull relation. Next, let us discuss the B2 dataset obtained in the current work. Figure 6 clearly shows a significant discrepancy between the Turnbull relation and our B2 data: when considered in isolation, the B2 data points do not fall on a straight line that passes through the origin. A significantly better linear trend is observed for the B2 dataset in the case of the Laird relation, but the slope is substantially different from that obtained using our MD data for bcc metals (CL (B2) = 0.60 vs. CL (bcc) = 0.36). If we used the Laird relation with the coefficient obtained for bcc metals to predict the SLI free energy of the Cu–Zr B2 phase, the error would be 51%. Since the B2 lattice reduces to the bcc structure if atomic type is ignored, our measurements indicate that the Laird coefficient must be sensitive to site occupancy. Moreover, the Laird relation does not seem to be able to account for the fact that the melting temperatures of Cu–Zr, Ni–Zr 2012, and Ni–Zr 2014 are similar (see Table II), but that the obtained Ni–Zr 2014 SLI free energy is about 1.5 times smaller than either the Cu–Zr or Ni–Zr 2012 SLI free energies (see Table III). Thus, although the linearity observed in the Laird relation for our B2 dataset appears much stronger than the corresponding Turnbull relation, and must therefore successfully capture fundamental contributions to the SLI free energy, it does not seem to be able to account for the entire dataset.

J. Chem. Phys. 142, 134705 (2015)

On the other hand, these issues appear to be less problematic when we take the liquid structure into consideration through the Ewing relation. First, the scatter in the B2 data about the observed trend line is significantly smaller in the Ewing relation than the corresponding scatter observed in the Turnbull relation. Second, the B2 trend line in the Ewing relation appears to coincide with the bcc trend line. And third, both the Ewing and Turnbull relations are able to predict that the Ni–Zr 2014 B2 phase should have significantly lower SLI free energy than either Cu–Zr or Ni–Zr 2012, in agreement with our measurements, from the fact that the ratio ∆Hm/k BTm is twice as small for Ni–Zr 2014 than for either Cu–Zr or Ni–Zr 2012. As a predictive relation, the common bcc/B2 trend line obtained from the Ewing relation is by no means perfect, but it is not bad: the largest relative deviation between measured values for B2 phases and values predicted by the Ewing relation based on the data for bcc metals is only 17% (for Cu–Zr). VII. CONCLUSIONS

In the present study, we determined the SLI free energy for several bcc metals and B2 compounds from MD simulation. We also developed two synthetic “Na” potentials to explore the effect of liquid structure and latent heat on the SLI free energy. The obtained MD data were compared with the empirical Turnbull, Laird, and Ewing relations. Overall, the MD data for bcc metals obtained within the present study are well described by all three relations. However, only the Laird and Ewing relations are able to provide the correct trend obtained from the sequence of “Na” potentials. The Laird relation appears to provide the best prediction for our MD data and data for bcc metals taken from the literature. Also, the Laird relation appears to predict that the bcc and B2 structures obey two distinctly different trends. However, the Laird relation has difficulty with B2 phases with similar melting temperatures but significantly different latent heats. Explicitly including a possible dependence on the liquid structure through the Ewing relation appears to alleviate this problem and also enables a reasonable prediction of our B2 data based on our bcc dataset. ACKNOWLEDGMENTS

The authors gratefully acknowledge very useful discussions with Professor X. Song. This work was supported by the U.S. Department of Energy, Office of Basic Energy Science, Division of Materials Sciences and Engineering. The research was performed at the Ames Laboratory. Ames Laboratory is operated for the U.S. Department of Energy by Iowa State University under Contract No. DE-AC02-07CH11358. 1D.

Kashchiev, Nucleation: Basic Theory with Applications (Butterworth Heinemann, Oxford, Boston, 2000). 2D. Turnbull, J. Appl. Phys. 21, 1022 (1950). 3M. Asta, C. Beckermann, A. Karma, W. Kurz, R. Napolitano, M. Plapp, G. Purdy, M. Rappaz, and R. Trivedi, Acta Mater. 57, 941 (2009). 4J. J. Hoyt, M. Asta, and A. Karma, Mater. Sci. Eng.: R: Rep. 41, 121 (2003). 5B. B. Laird and R. L. Davidchack, J. Phys. Chem. B 109, 17802 (2005). 6M. S. Daw and M. I. Baskes, Phys. Rev. B 29, 6443 (1984). 7D. Y. Sun, M. Asta, J. J. Hoyt, M. I. Mendelev, and D. J. Srolovitz, Phys. Rev. B 69, 020102 (2004). 8B. B. Laird, J. Chem. Phys. 115, 2887 (2001).


134705-12 9R.


H. Ewing, J. Cryst. Growth 11, 221 (1971).

10S. R. Wilson and M. I. Mendelev, Modell. Simul. Mater. Sci. Eng. 22, 065004

(2014). A. Becker, D. L. Olmsted, M. Asta, J. J. Hoyt, and S. M. Foiles, Phys. Rev. B 79, 054109 (2009). 12S. R. Wilson and M. I. Mendelev, Philos. Mag. 95, 224 (2015). 13J. J. Hoyt, M. Asta, and A. Karma, Phys. Rev. Lett. 86, 5530 (2001). 14M. W. Finnis and J. E. Sinclair, Philos. Mag. A 50, 45 (1984). 15D. G. Quirinale, G. E. Rustan, S. R. Wilson, M. J. Kramer, A. I. Goldman, and M. I. Mendelev, J. Phys.: Condens. Matter 27, 085004 (2015). 16M. I. Mendelev, M. J. Kramer, R. T. Ott, D. J. Sordelet, D. Yagodin, and P. Popel, Philos. Mag. 89, 967 (2009). 17Y. Mishin, M. J. Mehl, and D. A. Papaconstantopoulos, Phys. Rev. B 65, 224114 (2002). 18G. P. P. Pun and Y. Mishin, Philos. Mag. 89, 3245 (2009). 19M. I. Mendelev, M. J. Kramer, S. G. Hao, K. M. Ho, and C. Z. Wang, Philos. Mag. 92, 4454 (2012). 20J. J. Hoyt, M. Asta, and D. Y. Sun, Philos. Mag. 86, 3651 (2006). 21D. Y. Sun, M. Asta, and J. J. Hoyt, Phys. Rev. B 69, 174103 (2004). 22A. Kerrache, J. Horbach, and K. Binder, EPL 81, 58001 (2008). 23C. G. Tang and P. Harrowell, Nat. Mater. 12, 507 (2013). 24C. A. Becker, F. Tavazza, Z. T. Trautt, and R. A. B. de Macedo, Curr. Opin. Solid State Mater. Sci. 17, 277 (2013), http://www.ctcms.nist.gov/potentials. 25M. I. Mendelev, S. Han, D. J. Srolovitz, G. J. Ackland, D. Y. Sun, and M. Asta, Philos. Mag. 83, 3977 (2003). 26See supplementary material at http://dx.doi.org/10.1063/1.4916741 for a description of the Na and W potentials used in the present study, as well as to obtain figures that visualize the measured anisotropy in the SLI free energy and stiffness. 27M. I. Mendelev, M. J. Rahman, J. J. Hoyt, and M. Asta, Modell. Simul. Mater. Sci. Eng. 18, 074002 (2010). 28Y. F. Gao, Y. Yang, D. Y. Sun, M. Asta, and J. J. Hoyt, J. Cryst. Growth 312, 3238 (2010). 29J. Frenkel, Phys. Z. Sowjetunion 1, 498 (1932). 30H. A. Wilson, Philos. Mag. 50, 238 (1900). 31D. Y. Sun, M. I. Mendelev, C. A. Becker, K. Kudin, T. Haxhimali, M. Asta, J. J. Hoyt, A. Karma, and D. J. Srolovitz, Phys. Rev. B 73, 024116 (2006). 11C.

J. Chem. Phys. 142, 134705 (2015) 32X.

Y. Liu, J. B. Adams, F. Ercolessi, and J. A. Moriarty, Modell. Simul. Mater. Sci. Eng. 4, 293 (1996). 33H. C. Andersen, J. D. Weeks, and D. Chandler, Phys. Rev. A 4, 1597 (1971). 34J. D. Weeks, D. Chandler, and H. C. Andersen, J. Chem. Phys. 54, 5237 (1971). 35H. S. Kang, C. S. Lee, T. Ree, and F. H. Ree, J. Chem. Phys. 82, 414 (1985). 36V. B. Warshavsky and X. Song, J. Chem. Phys. 129, 034506 (2008). 37V. B. Warshavsky and X. Song, Phys. Rev. B 79, 014101 (2009). 38J. J. Hoyt, M. Asta, and B. Sadigh, Phys. Rev. Lett. 85, 594 (2000). 39A. Samanta, S. M. Ali, and S. K. Ghosh, Phys. Rev. Lett. 87, 245901 (2001). 40F. P. Buff, R. A. Lovett, and F. H. Stilling, Phys. Rev. Lett. 15, 621 (1965). 41J. D. Weeks, J. Chem. Phys. 67, 3106 (1977). 42D. Bedeaux and J. D. Weeks, J. Chem. Phys. 82, 972 (1985). 43M. Amini and B. B. Laird, Phys. Rev. B 78, 144112 (2008). 44M. Asta, J. J. Hoyt, and A. Karma, Phys. Rev. B 66, 100101 (2002). 45C. A. Becker, D. Olmsted, M. Asta, J. J. Hoyt, and S. M. Foiles, Phys. Rev. Lett. 98, 125701 (2007). 46A. A. Potter and J. J. Hoyt, J. Cryst. Growth 327, 227 (2011). 47R. Sibug-Aga and B. B. Laird, Phys. Rev. B 66, 144106 (2002). 48X. B. Feng and B. B. Laird, J. Chem. Phys. 124, 044707 (2006). 49S. Plimpton, J. Comput. Phys. 117, 1 (1995). 50J. R. Morris, Phys. Rev. B 66, 144104 (2002). 51R. L. Davidchack, J. R. Morris, and B. B. Laird, J. Chem. Phys. 125, 094710 (2006). 52K. A. Wu, A. Karma, J. J. Hoyt, and M. Asta, Phys. Rev. B 73, 094101 (2006). 53M. Jorge, P. Jedlovszky, and M. N. D. S. Cordeiro, J. Phys. Chem. C 114, 11169 (2010). 54R. L. C. Vink, J. Horbach, and K. Binder, J. Chem. Phys. 122, 134905 (2005). 55F. Spaepen, Acta Metall. 23, 729 (1975). 56F. Spaepen and R. B. Meyer, Scr. Metall. 10, 257 (1976). 57C. V. Thompson and F. Spaepen, Acta Metall. 31, 2021 (1983). 58S. Chantasiriwan and F. Milstein, Phys. Rev. B 58, 5996 (1998). 59A. I. Efimov, L. P. Belorukova, I. V. Vasilkova, and V. P. Chechev, Svoistva Neorganicheskih Soedinenii (Himiia, Leningrad, 1983), (in Russian). 60Y. Waseda, The Structure of Non-Crystalline Materials: Liquids and Amorphous Solids (McGraw-Hill International Book Co., New York, London, 1980).


Cascade morphology transition in bcc metals.

Fast, vacancy-free climb of prismatic dislocation loops in bcc metals.

Plastic anisotropy and dislocation trajectory in BCC metals.

Transformation of BCC and B2 High Temperature Phases to HCP and Orthorhombic Structures in the Ti-Al-Nb System. Part II: Experimental TEM Study of Microstructures.

Comparison between antigen-antibody binding energies and interfacial free energies.

Prediction of interface structures and energies via virtual screening.

Transformation of BCC and B2 High Temperature Phases to HCP and Orthorhombic Structures in the Ti-Al-Nb System. Part I: Microstructural Predictions Based on a Subgroup Relation Between Phases.

Introduction: Metals in Biology: METALS AT THE HOST-PATHOGEN INTERFACE.

MM molecular dynamics, free energies, and excited state properties.

Analytic bond-order potentials for the bcc refractory metals Nb, Ta, Mo and W.

Anisotropy and roughness of the solid-liquid interface of BCC Fe.

Condensed-matter physics: Glasses made from pure metals.

The free energy of mechanically unstable phases.

Metabolite concentrations, fluxes and free energies imply efficient enzyme usage.

Free energies of molecular clusters determined by guided mechanical disassembly.

Blind prediction of solvation free energies from the SAMPL4 challenge.

Bent molecules with a 60° central core angle that form B7 and B2 phases.

Problems of robustness in Poisson-Boltzmann binding free energies.

Sensitivity in Binding Free Energies Due to Protein Reorganization.

Real single ion solvation free energies with quantum mechanical simulation.

Unusual excess free energies of mixing in mixtures of partially fluorinated and hydrocarbon surfactants at the air-water interface: correlation with the structure of the layer.

Predicting hydration free energies with chemical accuracy: the SAMPL4 challenge.

Superconductivity in room-temperature stable electride and high-pressure phases of alkali metals.

Parity-Breaking Phases of Spin-Orbit-Coupled Metals with Gyrotropic, Ferroelectric, and Multipolar Orders.