Searching for Entropically Stabilized Phases - The Case of Silver Iodide

restricted predictive power since in the definition of the collective variables (CVs) the nature of the final phase is ... their equilibrium phases is...
0 downloads 5 Views 1MB Size
Subscriber access provided by UCL Library Services

Article

Searching for Entropically Stabilized Phases - The Case of Silver Iodide Dan Mendels, James J. McCarty, Pablo Miguel Piaggi, and Michele Parrinello J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.7b11002 • Publication Date (Web): 26 Dec 2017 Downloaded from http://pubs.acs.org on December 30, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry C is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Searching for Entropically Stabilized Phases - the Case of Silver Iodide Dan Mendels,†,¶ James McCarty,†,¶ Pablo M. Piaggi,‡,¶ and Michele Parrinello∗,†,¶ †Department of Chemistry and Applied Biosciences, ETH Zurich, Zurich, Switzerland ´ ‡Theory and Simulation of Materials (THEOS), Ecole Polytechnique F´ed´erale de Lausanne ¶Facolt`a di Informatica, Instituto di Scienze Computationali, Universit della Svizzera italiana, Via Giuseppe Buffi 13, CH-6900 Lugano, Switzerland E-mail: [email protected]

Abstract Computational crystal phase prediction means have become an important tool in materials science for the acceleration of material discovery. Material phases stabilized by substantial entropic effects, however, cannot be accounted for with today’s widely used methods. Here, we propose a new search method for such phases in the framework of Metadynamics using surrogates for entropy and enthalpy as collective variables. Applying the method to AgI at low pressure, the method is found effective in discovering the material’s entropically stabilized superionic phase as well as its enthalpically stabilized β phase. Through calculating the material’s free energy surfaces at different temperatures, we additionally estimate the system’s melting temperature as well as the transition temperature between the β and superionic phases. The calculation of the free energy is accelerated by the use of a path collective variable defined in the space of our surrogates for entropy and enthalpy.

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Introduction One of the main tasks of materials science is the search for materials with new properties. Since experiments are costly and time consuming, computational methods have been called in to speed up this search. Of course, materials’ properties depend crucially on their crystalline structures. For this reason several structural prediction methods have been proposed and their remarkable degree of success has made them very popular. 1–6 These methods are based on the search on the potential energy surface for low lying minima. Entropy effects are either neglected or added at most within the quasi-harmonic approximation. 7,8 This leaves out important classes of materials whose structure is stabilized by strong entropic effects like superionic or plastic crystals. We propose here a method to address this problem and we use as an example AgI, a compound that over the years has attracted considerable attention. This material is known to crystallize at low pressures into a super-ionic phase (α phase), characterized by high ionic conductance. In this phase the iodine ions form a bcc lattice while the silver ions diffuse through the iodine lattice with a diffusion coefficient comparable to that of the liquid. 9–11 At low pressure and sufficiently low temperatures AgI subsequently transforms into a Wurtzite structure (β phase). Thus, AgI is almost an ideal testing ground for methods aimed at establishing the phase diagram of such complex materials. We should, thus, pretend that nothing is known about the system except a reasonable good approximation for its interaction potential and develop methods to discover its different polymorphs. A naive way of approaching this problem would be to simulate the system and, starting from the liquid state, slowly reduce the temperature and observe the system’s spontaneously transform first into the α phase and later into the β Wurtzite structure. However, this strategy is hampered by the very long time scale that the nucleation of a crystal from a liquid takes. This is a well recognized problem and much effort has been devoted to its solution. 12,13 Several enhanced sampling methods have been proposed to study the nucleation process. Most methods rely on the introduction of suitable collective co2

ACS Paragon Plus Environment

Page 2 of 17

Page 3 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

ordinates or order parameters whose sampling is accelerated to observe nucleation in an affordable computation time. This approach has been extremely successful and much has been learned about nucleation. 14 While illuminating, these simulations have a somewhat restricted predictive power since in the definition of the collective variables (CVs) the nature of the final phase is taken into account. However, very recently Piaggi et al. have circumvented this problem and proposed the use of two CVs that are surrogates for enthalpy and entropy. 15 The rational for this choice is that in a liquid-solid phase transition there is a trade off between enthalpy and entropy. The enhanced sampling method used, metadynamics, was able to enhance the fluctuations of these variables and thus lead spontaneously to a transition to the thermodynamic stable state. In the work of Piaggi et al. it was shown how this scheme lead to predict, in agreement with experiments, that Na and Al crystallize into the bcc and fcc forms respectively. Clearly, Na and Al are very simple cases and predicting their equilibrium phases is not a very remarkable achievement and could have been done more simply using enthalpy based methods. However, the case of AgI is much more complicated being a multi-component system and the α phase being favored over the enthalpic minimum of the β phase by the strong entropic contributions coming from the diffusing Ag ions. We find that by using the entropy and enthalpy inspired collective variables of ref., 15 the existence of the α phase can be predicted. Furthermore, one can estimate the melting temperature and compute the entropy of fusion. In addition, one can go from the α to the β Wurtzite phase and also in this case calculate the transition temperature, thus completing the low pressure region of the phase diagram of this important substance and providing a model for further investigation. In the latter case, however, being generic the two CVs were found not to be particularly efficient and in order to obtain an accurate evaluation of the free energy difference between the α and β phases we supplemented our calculation with a path variable metadynamics calculation using a path from α to β in the space of entropy and enthalpy.

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 17

Choice of Collective Variables The use of metadynamics as an enhancement sampling method requires the introduction of CVs. As suggested earlier and following the footsteps of Piaggi et al. we shall use surrogates for enthalpy and entropy. The former can be naturally defined while defining the latter is not obvious. Piaggi et al. have suggested to use the expression

S2k

Z∞ = −2πkB ρk [gk (r) ln gk (r) − gk (r) + 1]r2 dr

(1)

0

where ρk is the system’s density, kB the Boltzmann constant, and gk (r) the radial distribution function, that has been shown to provide an excellent approximation to the entropy of a single component uniform liquid. 16 Despite AgI being obviously a two component system, we found instrumental to use the entropy of either one of the two species for the construction of our CVs, and hence in the remainder of the paper S2k , ρk and gk (r) will refer to the entropy, density and radial distribution of the system’s kth species, respectively, where k = I, Ag. In practice we shall limit the upper bound of the integral in Eq. 1 and use a mollified version of the radial distribution function as

gm (r) =

X 1 (r−rij )2 1 2σ 2 √ e 4πN ρr2 i6=j 2πσ 2

(2)

where rij is the distance between particles i and j, N the number of particles in the system and σ a broadening parameter. The mollification is required to ensure that the derivatives of gm (r) relative to the atomic positions are continues. Finally, while both H, i.e. enthalpy, and S2k have thermodynamic meaning only when averaged over time, for the purpose of defining the CVs in the present study, we utilized their instantaneous values. In the reminder of the paper we, henceforth, refer to the CVs corresponding to the above definitions as skS and sH where

sH = 4

H . N

ACS Paragon Plus Environment

(3)

Page 5 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Methods Simulation Implementation All simulations were carried out using LAMMPS 17 patched with PLUMED 2 18 in a simulation box consisting of 250 silver ions and 250 iodine ions using the interaction potential of ref. 10 A time-step of 1fs was used in all simulations and constant temperatures were maintained using the velocity rescaling thermostat of Bussi et al. 19 with a relaxation time of 0.1ps. A constant pressure of 1000 bar was additionally kept in all simulations using the Parrinello-Rahman method 20,21 with a relaxation time of 1ps to which a pressure correction term was added. This correction term was also used in ref. 10 to set the densities close to experimental values. In the study of the liquid to α phase we constrained the diagonal terms of the h matrix of the Parrinello-Rahman Lagrangian to be coupled and the off-diagonal ones to zero. This avoided the occurrence of inconveniently oblong simulation boxes in the liquid phase. 15 Similarly, in the study of the α → β solid-solid transition we found convenient to impose that only the h matrix components h11 , h32 and h33 could change. This choice was made after having observed the α → β transition taking place without any restraint. In these calculations due to modular invariance the h matrix could take different forms and the crystal align in different orientations with respect to the Cartesian axis. Using the restriction described above the transformation took place without hindrance and the analysis of the system became also more simple. Long range electrostatic interactions were calculated with a standard Ewald summation with a real space cutoff of 10˚ A and a root mean squared force accuracy of 10−4 and skS were calculated using a radial cutoff of rmax =12.5˚ A and rmax =10˚ A in simulations initialized in the liquid phase and simulations initialized in the α phase, respectively. All simulations were run with bias factors spanning between γ = 80 and γ = 160 with a hill deployment time interval of τG = 500f s and initial hill heights spanning between W = 0.05eV and 0.1eV . The employed hill widths were σskS = 0.15 and σsH = 0.005, for skS and sH , respectively. 5

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 17

Paths Construction and Path CV Simulation Implementation For the implementation of path based Well Tempered Metadynamics (WTMD) simulations, paths between the β and α states were constructed in the sH and sAg S plane using Eqs. 4 and 5, 22 R1

te−λ|g·(R−R(t))| dt s(R) = lim R0 1 λ→∞ e−λ|g·R−R(t)| dt 0 1 z(R) = lim − ln λ→∞ λ

Z

(4)

1

e−λ|g·(R−R(t))| dt

(5)

0

where s(R) is the path coordinate corresponding to the position along the parameterized path R(t) closest to the simulation’s position R in CV space (where R(0) = Rα and R(1) = Rβ ) and z the geometrical distance between the two. To account for sH ’s and sAg S ’s different units, a scaling vector g = (g1 , g2 ) was additionally included into the Eqs. 4 and 5, where g1 and g2 corresponded to the inverse values of the typical scales of each of the CVs (e.g. Ag Ag −1 −1 with δsAg g1 = [sAg or g1 = [ 12 (hδsAg S (α)i + hδsS (β)i)] S being the typical S (α) − sS (β)]

(unbiased) fluctuations of sAg S in each of the states). This, thus, insured that the path’s components exerted on average a similar influence in determining the system’s position on the path. As an initial approach, path contours were constructed to follow the paths taken by the system in (unbiased) temperature induced transitions and for simplicity were made of two straight line segments, each connecting α or β with a carefully selected point positioned in the middle section of the path. All paths were constructed using 18 to 34 frames (i.e. reference points), with path parameter values corresponding to the relation λ ≈

2.3 , ∆s

∆s

being the distance between adjacent frames. Simulations were run with the parametrization used in the exploration simulations, albeit with bias factors spanning between γ = 70 and γ = 120 and a hill widths of σs = 0.2. To limit the system’s entry into undesired regions of the CV plane, lower and upper parabolic constraints were also added near the path CV’s 6

ACS Paragon Plus Environment

Page 7 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

edges and in similarity an upper constraint was set on z at z ≈ 0.8. Additionally, upper and lower parabolic constraints were set on sIS at -4 and -12.5, respectively, and the fcc parameter CV introduced in ref. 23 was employed to enforce a constraint on the iodine ions structure at sf cc = −0.128. The calculation of sf cc was implemented with respect to the first sphere of I neighbors of each I, where the CV’s parameters were kept to their default values, namely φ = 0, θ = 0, ψ = 0 and α = 3.

Results From liquid to α superionic We attempt here to test our ability to discover crystalline phases straight from the liquid phase. To this effect, WTMD simulations employing enthalpy and either sIS or sAg S were initiated from the system’s liquid phase at different temperatures in the vicinity of the system’s putative melting temperature (technical details can be found in the methods’ section). I While utilizing sAg S did not induce crystallization, using sS clear structural transitions were

observed. Given the large number of transitions observed, the difference in free energy between liquid and solid can be accurately determined. The resulting free energy surface (FES) shows two clear metastable states that correspond to the liquid and to the α phase that is thus predicted to be the dominant crystal phase in this temperature. This is not so surprising since the liquid to α transition can be regarded in rather loose terms as a crystallization of the I ions’ subsystem. In fact, in the α phase the I ions formed a bcc lattice while the Ag ions exhibited a diffusive behaviour with a diffusion coefficient comparable to that of the liquid. Here, in addition, the liquid-like form of the Ag ions in the α phase was found to be reflected in the Ag ions’ subsystem entropy, which was found to retain its value from the liquid phase, being sAg S = −0.95 ± 0.2[kB ]. In the calculation of the liquid-solid free energy difference ∆G we used the WTMD relation 24 between the bias V (s) and the free energy: 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

G(s) = −

 γ  V (s). γ−1

Page 8 of 17

(6)

where γ is the WTMD bias factor and ∆G is as usual given by: R ∆G = −kB T log R

dse A

G(s) BT

−k

dse B

G(s) BT

−k

(7)

where in this case A and B correspond to the α phase and liquid basins, respectively. We repeated the calculation at different temperatures. We find our data compatible with a linear dependence as shown in Fig. 1. The intercept with the ∆G = 0 line gives for the freezing temperature the value T = 905◦ ± 10◦ K. This value is bracketed between the values T ≈ 1050◦ K and T ≈ 750◦ K obtained by overcooling or overheating (with a cooling/heating rate of 0.1◦ K/ps) the liquid and the solid, respectively. We note that the utilized model’s transition temperature is somewhat lower than the experimental one which is T = 831◦ K.

Figure 1: (a)Gibbs free energy as function of sH and sIS at T = 900◦ K. Also shown, snapshots of the iodine ions in the two phases colored using the common neighbor analysis tool in OVITO. 25 As expected, the iodine ions are distributed in a liquid like form within the liquid phase (white) and form a bcc lattice (blue) in the α phase. (b)Temperature dependence of ∆G. Shown also, a linear fit to the data points adjacent to the ∆G = 0 line, indicating the estimated phase transition temperature.

8

ACS Paragon Plus Environment

Page 9 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Simulations Initialized in the System’s α Phase To assess the ability of our approach to find more stabilized crystalline phases in the system at lower temperatures, WTMD simulations employing enthalpy and either sIS or sAg S were initiated at the system’s α phase at different temperatures below the freezing point. The use of sIS proved ineffective, instead sAg S promoted the appearance of a new crystalline phase (see Fig. 2 for one example of such a transition). In this new phase which we identify as the β phase, the Ag ions cease to diffuse and the system crystallizes into a Wurtzite structure. Here, the nature of the transition is somewhat specular relative to the previously discussed liquid → α transformation. There, in first approximation, we could regard it as a freezing of the I subsystem, here instead in a rather superficial way of speaking the transition can be regarded as a freezing of the Ag subsystem. Hence, the different efficacy of sIS and sAg S in driving the two transformations. While sAg S was able to induce the α → β transition, it was difficult to obtain an accurate FES given that during the simulations several spurious non diffusive states were visited, thus hampering efficient sampling of the relevant regions in the system’s phase space. We examined the stability of these states and found them unstable. In order to get around this difficulty, we decided to resort to the path CV. We first drew a reference path in the 22 {sAg we introduced two variables s and z. The first one measures S , sH } space and as in ref.

the progression along the reference path while the other gives the distance from the guessed path (see methods for details). We found the results to be insensitive to small changes in the definition of the path. We also checked whether the addition of an extra CV that describes a local close packed arrangement of the I would have an effect on the rate of convergence, since it combines information on the crystal structure. Quite satisfyingly it did, however, in the spirit of this work we confined ourselves only to the use of sS and sH . As can be seen in Fig. 3 these two variables alone were able to focalize the configurational space exploration to the relevant regions. In the same picture we also monitor the Ag ions’ diffusion that is very low 9

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

in the β phase but that becomes very significant every time the system visits the α phase. Given the high number of recrossings, the free energy can be accurately determined and from the dependence of ∆Gαβ on temperature, the transition temperature estimated to be at T = 560◦ ± 8◦ K (see Fig. 4), where the obtained value is bracketed between the values T ≈ 630◦ K and T ≈ 480◦ K obtained by overcooling or overheating (with a cooling/heating rate of 0.1◦ K/ps) the α and the β phases, respectively. We note that the utilized model’s transition temperature is higher than the experimental value which is T = 420◦ K.

Figure 2: A WTMD exploration simulation run at T = 550◦ K. Top: The Ag ions’ entropy CV. Middle: The enthalpy CV. Bottom: the mean squared displacement of the Ag ions. Errors of all displayed quantities are contained within their corresponding displayed symbols.

10

ACS Paragon Plus Environment

Page 10 of 17

Page 11 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 3: Excerpts from the time dependence behaviour of a WTMD simulation run at T = 530◦ K. a)The path CV, b)The Ag ions entropy CV, c)the system’s enthalpy CV, d)mean squared displacement of the Ag ions. Frames (17-22) of the path CV correspond to the β phase, whereas frames (3-7) to the α one. Errors of all displayed quantities are contained within their corresponding displayed symbols.

11

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4: a)Gibbs free energy as function of the path CV at T = 530◦ K. In solid line, attained average FES, and in dashed lines, the bounds of its statistical uncertainty. Also shown, snapshots of the I ions in the two phases colored using the common neighbor analysis tool in OVITO. As expected, in the α phase the iodine ions form a bcc lattice (blue). In contrast, the β phase the iodine ions form an hcp (red) lattice. (The observed stacking fault (green) was found to appear in both biased and unbiased simulations.) (b)Temperature dependence of ∆G. Shown also, a linear fit to the data points adjacent to the ∆G = 0 line, indicating the estimated phase transition temperature.

12

ACS Paragon Plus Environment

Page 12 of 17

Page 13 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Conclusions In conclusion, an approach for searching for uncharted entropically stabilized phases was introduced. The approach, applied within the Metadynamics framework, relies on the employment of surrogates for an investigated system’s enthalpy and entropy as CVs, enabling the scanning of the system’s FES in the space spanned by these CVs. In applying it to AgI, the approach was found to be highly effective in inducing the exploration of the CV space, in enabling the ’discovery’ of the material’s entropically stabilized super-ionic α phase within simulations initialized in the system’s liquid phase and in turn the calculation of the system’s melting temperature. In utilizing this approach to search for additional stable phases starting from the discovered α phase, the approach was found effective in finding also the system’s enthalpically stabilized β phase. Here, attaining a converged estimate of the system’s corresponding FES required combining the introduced approach with the path CV scheme, which in addition in turn enabled the calculation of the transition temperature between the α and β phases. In light of the efficacy exhibited by the approach in the present study, we envision it to constitute a new and effective tool for the discovery and characterization of material phases which stabilization relies on entropic contributions as well as ones which are solely enthalpically stabilized.

Acknowledgments The research was supported by the European Union Grant No. ERC-2014-AdG-670227 / VARMET. The authors also acknowledge funding from NCCR MARVEL, funded by the Swiss National Science Foundation. Calculations were carried out on the M¨onch cluster at the Swiss National Super-computing Center (CSCS).

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

References (1) Lyakhov, A. O.; Oganov, A. R.; Stokes, H. T.; Zhu, Q. New developments in evolutionary structure prediction algorithm USPEX. Comput. Phys. Commun. 2013, 184, 1172–1182. (2) Woodley, S. M.; Catlow, R. Crystal structure prediction from first principles. Nat. Mater. 2008, 7, 937–46. (3) Wales, D. J.; Scheraga, H. A. Global Optimization of Clusters, Crystals, and Biomolecules. Science 1999, 285, 1368–1372. (4) Kirkpatrick, S.; Gelatt, C. D.; Vecchi, M. P. Optimization by simulated annealing. Science 1983, 220, 671–680. (5) Pannetier, J.; Bassas-Alsina, J.; Rodriguez-Carvajal, J.; Caignaert, V. Prediction of crystal structures from crystal chemistry rules by simulated annealing. Nature 1990, 346, 343–345. (6) Martoˇ na´k, R.; Laio, A.; Parrinello, M. Predicting Crystal Structures: The ParrinelloRahman Method Revisited. Phys. Rev. Lett. 2003, 90, 075503. (7) Togo, A.; Tanaka, I. First principles phonon calculations in materials science. Scr. Mater. 2015, 108, 1–5. (8) Oganov, A. R.; Glass, C. W. Crystal structure prediction using ab initio evolutionary techniques: Principles and applications. J. Chem. Phys. 2006, 124, 244704. (9) Vashishta, P.; Rahman, A. Ionic motion in α-AgI. Phys. Rev. Lett. 1978, 40, 1337. (10) Parrinello, M.; Rahman, A.; Vashishta, P. Structural Transitions in Superionic Conductors. Phys. Rev. Lett. 1983, 50, 1073–1076.

14

ACS Paragon Plus Environment

Page 14 of 17

Page 15 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(11) Wood, B. C.; Marzari, N. Dynamical Structure, Bonding, and Thermodynamics of the Superionic Sublattice in α- AgI. Phys. Rev. Lett. 2006, 97, 166401. (12) Mandell, M.; McTague, J.; Rahman, A. Crystal nucleation in a three-dimensional Lennard-Jones system: A molecular dynamics study. J. Chem. Phys. 1976, 64, 3699– 3702. (13) Stillinger, F. H.; Weber, T. A. Study of melting and freezing in the Gaussian core model by molecular dynamics simulation. J. Chem. Phys. 1978, 68, 3837–3844. (14) Giberti, F.; Salvalaglio, M.; Parrinello, M. Metadynamics studies of crystal nucleation. IUCrJ 2015, 2, 256–266. (15) Piaggi, P. M.; Valsson, O.; Parrinello, M. Enhancing entropy and enthalpy fluctuations to drive crystallization in atomistic simulations. Phys. Rev. Lett. 2017, 119, 015701. (16) Nettleton, R.; Green, M. Expression in terms of molecular distribution functions for the entropy density in an infinite system. J. Chem. Phys. 1958, 29, 1365–1370. (17) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1 – 19. (18) Tribello, G. A.; Bonomi, M.; Branduardi, D.; Camilloni, C.; Bussi, G. {PLUMED} 2: New feathers for an old bird. Comput. Phys. Commun. 2014, 185, 604 – 613. (19) Bussi, G.; Donadio, D.; Parrinello, M. Canonical sampling through velocity rescaling. J. Chem. Phys. 2007, 126, 014101. (20) Parrinello, M.; Rahman, A. Crystal structure and pair potentials: A moleculardynamics study. Phys. Rev. Lett. 1980, 45, 1196. (21) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant pressure molecular dynamics algorithms. J. Chem. Phys. 1994, 101, 4177–4189. 15

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(22) Branduardi, D.; Gervasio, F. L.; Parrinello, M. From A to B in free energy space. J. Chem. Phys. 2007, 126, 054103. (23) Angioletti-Uberti, S.; Ceriotti, M.; Lee, P. D.; Finnis, M. W. Solid-liquid interface free energy through metadynamics simulations. Phys. Rev. B 2010, 81, 125416. (24) Barducci, A.; Bussi, G.; Parrinello, M. Well-tempered metadynamics: a smoothly converging and tunable free-energy method. Phys. Rev. Lett. 2008, 100, 020603. (25) Stukowski, A. Visualization and analysis of atomistic simulation data with OVITO - the Open Visualization Tool. Modelling and Simulation in Materials Science and Engineering 2009, 18, 1.

16

ACS Paragon Plus Environment

Page 16 of 17

Page 17 of 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

TOC Graphic

17

ACS Paragon Plus Environment