Free energy barriers to melting in atomic clusters

100 Years of the Lennard-Jones Potential

It is now 100 years since Lennard-Jones published his first paper introducing the now famous potential that bears his name. It is therefore timely to reflect on the many achievements, as well as the limitations, of this potential in the theory of atomic and molecular interactions, where applications range from descriptions of intermolecular forces to molecules, clusters, and condensed matter.

A universal signature in the melting of metallic nanoparticles

Predicting when phase changes occur in nanoparticles is fundamental for designing the next generation of devices suitable for catalysis, biomedicine, optics, chemical sensing and electronic circuits. The estimate of the temperature at which metallic nanoparticles become liquid is, however, a challenge and a standard definition is still missing. We discover a universal feature in the distribution of the atomic-pair distances that distinguishes the melting transition of monometallic nanoparticles. We analyse the solid-liquid change of several late-transition metals nanoparticles, i.e. Ni, Cu, Pd, Ag, Au and Pt, through classical molecular dynamics. We consider various initial shapes from 146 to 976 atoms, corresponding to the 1.5-4.1 nm size range, placing the nanoparticles in either a vacuum or embedded in a homogeneous environment, simulated by an implicit force-field. Regardless of the material, its initial shape, size and environment, the second peak in the pair-distance distribution function, expected at the bulk lattice distance, disappears when the nanoparticle melts. As the pair-distance distribution is a measurable quantity, the proposed criterion holds for both numerical and experimental investigations. For a more straightforward calculus of the melting temperature, we demonstrate that the cross-entropy between a reference solid pair-distance distribution function and the one of nanoparticles at increasing temperatures present a quasi-first order transition at the phase-change temperature.

Free energy calculations from adaptive molecular dynamics simulations with adiabatic reweighting

We propose an adiabatic reweighting algorithm for computing the free energy along an external parameter from adaptive molecular dynamics simulations. The adaptive bias is estimated using Bayes identity and information from all the sampled configurations. We apply the algorithm to a structural transition in a cluster and to the migration of a crystalline defect along a reaction coordinate. Compared to standard adaptive molecular dynamics, we observe an acceleration of convergence. With the aid of the algorithm, it is also possible to iteratively construct the free energy along the reaction coordinate without having to differentiate the gradient of the reaction coordinate or any biasing potential.

Decoding the energy landscape: extracting structure, dynamics and thermodynamics

Describing a potential energy surface in terms of local minima and the transition states that connect them provides a conceptual and computational framework for understanding and predicting observable properties. Visualizing the potential energy landscape using disconnectivity graphs supplies a graphical connection between different structure-seeking systems, which can relax efficiently to a particular morphology. Landscapes involving competing morphologies support multiple potential energy funnels, which may exhibit characteristic heat capacity features and relaxation time scales. These connections between the organization of the potential energy landscape and structure, dynamics and thermodynamics are common to all the examples presented, ranging from atomic and molecular clusters to biomolecules and soft and condensed matter. Further connections between motifs in the energy landscape and the interparticle forces can be developed using symmetry considerations and results from catastrophe theory.

New generalization of supersymmetric quantum mechanics to arbitrary dimensionality or number of distinguishable particles

We present here a new approach to generalize supersymmetric quantum mechanics to treat multiparticle and multidimensional systems. We do this by introducing a vector superpotential in an orthogonal hyperspace. In the case of N distinguishable particles in three dimensions this results in a vector superpotential with 3N orthogonal components. The original scalar Schrödinger operator can be factored using a 3N-component gradient operator and introducing vector "charge" operators: Q(1) and Q(1)(dagger). Using these operators, we can write the original (scalar) Hamiltonian as H(1) = Q(1)(dagger) x Q(1) + E(0)((1)), where E(0)((1)) is the ground-state energy. The second sector Hamiltonian is a tensor given by H(2) = Q(1)Q(1)(dagger) + E(0)((1)) and is isospectral with H(1). The vector ground state of sector 2, psi(0)((2)), canbe used with the charge operator Q(1)(dagger) to obtain the excited-state wave function of the first sector. In addition, we show that H(2) can also be factored in terms of a sector 2 vector superpotential with components W(2j) = -(partial partial differential ln psi(0j)((2)))/partial partial differentialx(j). Here psi(0j)((2)) is the jth component of psi(0)((2)). Then one obtains charge operators Q(2) and Q(2)(dagger) so that the second sector Hamiltonian can be written as H(2) = Q(2)(dagger)Q(2) + E(0)((2)). This allows us to define a third sector Hamiltonian which is a scalar, H(3) = Q(2) x Q(2)(dagger) + E(0)((2)). This prescription continues with the sector Hamiltonians alternating between scalar and tensor forms, both of which can be treated by the variational method to obtain approximate solutions to both scalar and tensor sectors. We demonstrate the approach with examples of a pair of separable 1D harmonic oscillators and the example of a nonseparable 2D anharmonic oscillator (or equivalently a pair of coupled 1D oscillators). We consider both degenerate and nondegenerate cases. We also present a generalization to arbitrary curvilinear coordinate systems in the Appendix.

Replica exchange statistical temperature Monte Carlo

The replica exchange statistical temperature Monte Carlo algorithm (RESTMC) is presented, extending the single-replica STMC algorithm [J. Kim, J. E. Straub, and T. Keyes, Phys. Rev. Lett. 97, 050601 (2006)] to alleviate the slow convergence of the conventional temperature replica exchange method (t-REM) with increasing system size. In contrast to the Gibbs-Boltzmann sampling at a specific temperature characteristic of the standard t-REM, RESTMC samples a range of temperatures in each replica and achieves a flat energy sampling employing the generalized sampling weight, which is automatically determined via the dynamic modification of the replica-dependent statistical temperature. Faster weight determination, through the dynamic update of the statistical temperature, and the flat energy sampling, maximizing energy overlaps between neighboring replicas, lead to a considerable acceleration in the convergence of simulations even while employing significantly fewer replicas. The performance of RESTMC is demonstrated and quantitatively compared with that of the conventional t-REM under varying simulation conditions for Lennard-Jones 19, 31, and 55 atomic clusters, exhibiting single- and double-funneled energy landscapes.

Thermodynamics of atomic clusters using variational quantum hydrodynamics

Small clusters of rare-gas atoms are ideal test cases for studying how quantum delocalization affects both the thermodynamics and the structure of molecular scale systems. In this paper, we use a variational quantum hydrodynamic approach to examine the structure and dynamics of (Ne)n clusters, with n up to 100 atoms, at both T = 0 K and for temperatures spanning the solid-to-liquid transition in bulk Ne. Finite temperature contributions are introduced to the grand potential in the form of an "entropy" potential. One surprising result is the prediction of a negative heat capacity for very small clusters that we attribute to the nonadditive nature of the total free-energy for very small systems.

Potential energy and free energy landscapes

Familiar concepts for small molecules may require reinterpretation for larger systems. For example, rearrangements between geometrical isomers are usually considered in terms of transitions between the corresponding local minima on the underlying potential energy surface, V. However, transitions between bulk phases such as solid and liquid, or between the denatured and native states of a protein, are normally addressed in terms of free energy minima. To reestablish a connection with the potential energy surface we must think in terms of representative samples of local minima of V, from which a free energy surface is projected by averaging over most of the coordinates. The present contribution outlines how this connection can be developed into a tool for quantitative calculations. In particular, stepping between the local minima of V provides powerful methods for locating the global potential energy minimum, and for calculating global thermodynamic properties. When the transition states that link local minima are also sampled we can exploit statistical rate theory to obtain insight into global dynamics and rare events. Visualizing the potential energy landscape helps to explain how the network of local minima and transition states determines properties such as heat capacity features, which signify transitions between free energy minima. The organization of the landscape also reveals how certain systems can reliably locate particular structures on the experimental time scale from among an exponentially large number of local minima. Such directed searches not only enable proteins to overcome Levinthal's paradox but may also underlie the formation of "magic numbers" in molecular beams, the self-assembly of macromolecular structures, and crystallization.

Equilibrium density of states and thermodynamic properties of a model glass former

This paper presents an analysis of the thermodynamics of a model glass former. We have performed equilibrium sampling of a popular binary Lennard-Jones model, employing parallel tempering Monte Carlo to cover the crystalline, amorphous, and liquid regions of configuration space. Disconnectivity graphs are used to visualize the potential energy landscape in the vicinity of a crystalline geometry and in an amorphous region of configuration space. The crystalline global minimum is separated from the bulk of the minima by a large potential energy gap, leading to broken ergodicity in conventional simulations. Our sampling reveals crystalline global minima that are lower in potential energy than some of the previous candidates. We present equilibrium thermodynamic properties based on parallel tempering simulations, including heat capacities and free energy profiles, which depend explicitly on the crystal structure. We also report equilibrium melting temperatures.

Thermodynamics of water octamer in a uniform electric field

We study the water octamer in a uniform electric field using the all-exchanges parallel tempering Monte Carlo method in the canonical ensemble. The heat capacity, quenched energy configurations, and the order parameter Q(4) are employed to understand the phase changes observed as a function of temperature and the strength of the applied electric field. At a low field strength of 0.1 V A(-1) a solidlike to liquidlike "melting" transition is detected. The corresponding heat capacity peak appears around 206 K, where Q(4) shows a significant change of slope. For E> or =0.5 V A(-1) such features are absent. However, at E=0.5 V A(-1) we find a solidlike to solidlike transition between cubic and extended structures around T approximately 25 K.

Superheating and solid-liquid phase coexistence in nanoparticles with nonmelting surfaces

We present a phenomenological model of melting in nanoparticles with facets that are only partially wet by their liquid phase. We show that in this model, as the solid nanoparticle seeks to avoid coexistence with the liquid, the microcanonical melting temperature can exceed the bulk melting point and that the onset of coexistence is a first-order transition. We show that these results are consistent with molecular dynamics simulations of aluminum nanoparticles which remain solid above the bulk melting temperature.

Multiple structural transformations in Lennard-Jones clusters: Generic versus size-specific behavior

The size-temperature "phase diagram" for Lennard-Jones clusters LJn with sizes up to n=147 is constructed based on the analysis of the heat capacities and orientational bond order parameter distributions computed by the exchange Monte Carlo method. Two distinct types of "phase transitions" accompanied by peaks in the heat capacities are proven to be generic. Clusters with Mackay atom packing in the overlayer undergo a lower-temperature melting (or Mackay-anti-Mackay) transition that occurs within the overlayer. All clusters undergo a higher-temperature transition, which for the three-layer clusters is proven to be the 55-atom-core-melting transition. For the two-layer clusters, the core/overlayer subdivision is ambiguous, so the higher-temperature transition is better characterized as the breaking of the local icosahedral coordination symmetry. A pronounced size-specific behavior can typically be observed at low temperatures and often occurs in clusters with highly symmetric global minima. An example of such behavior is LJ135, which undergoes a low-temperature solid-solid transition, besides the two generic transitions, i.e., the overlayer reconstruction and the core melting.

A Self-Consistent Field Quantum Hydrodynamic Approach for Molecular Clusters

We present a novel self-consistent orbital-free method useful for quantum clusters. The method uses a hydrodynamical approach based on the de Broglie-Bohm description of quantum mechanics to satisfy an orbital-free density functional-like Euler-Lagrange equation for the ground state of the system. In addition, we use an information theoretical approach to obtain the optimal density function derived from a series of statistical sample points in terms of density approximates. These are then used to calculate an approximation to the quantum force in the hydrodynamic description. As a demonstration of the utility and flexibility of the approach, we compute the lowest-energy structures for small rare-glass clusters of argon and neon with 4, 5, 13, and 19 atoms. Extension to more complex systems is straightforward.

Transition from icosahedral to decahedral structure in a coexisting solid-liquid nickel cluster

We have used molecular dynamics simulations to construct a microcanonical caloric curve for a 1415 atom Ni icosahedron. Prior to melting, the Ni cluster exhibits static solid-liquid phase coexistence. Initially, a partial icosahedral structure coexists with a partially wetting melt. However, at energies very close to the melting point the icosahedral structure is replaced by a truncated decahedral structure that is almost fully wet by the melt. This structure remains until the cluster fully melts. The transition appears to be driven by a preference for the melt to wet the decahedral structure.

Static, transient, and dynamic phase coexistence in metal nanoclusters

Molecular dynamics simulations are used to examine static and dynamic coexistence between solid and liquid phases in nanoscale silver, copper, and nickel clusters. We find static coexistence in the 561-atom copper icosahedron, the 561-atom silver icosahedron, and the 923-atom nickel icosahedron, and in cluster sizes above these thresholds, but not in smaller clusters. Nonetheless, in smaller clusters we typically observe either dynamic coexistence between fully solid and liquid states or transient coexistence which is essentially dynamic coexistence between a fully solid state and a solid-liquid state.

Structural Behavior and Self-Assembly of Lennard-Jones Clusters on Rigid Surfaces

The phase behavior and surface pattern formation for intermediate size Lennard-Jones clusters on rigid surfaces are examined. We use a parallel tempering Monte Carlo algorithm, in the canonical ensemble. Tempering is done over the temperature domain in most of the calculations. A two-dimensional temperature and Hamiltonian tempering algorithm is also implemented, to examine its usefulness in investigating this type of problem. In general, we observe gas phase systems as they undergo a condensation transition on the surface, followed by a freezing transition. The final solid state pattern formed by the cluster on the surface is the result of a number of competing effects. First, there is a competition between attraction within the cluster and that between cluster and surface atoms. Second, a monolayer of Lennard-Jones atoms tends to pack in a hexadic geometry. This geometry is frustrated on a surface with a different symmetry. The molecular organization of the substrate has a serious impact on the cluster packing. The surface morphology and the size mismatch between cluster and surface atoms, along with the relative interaction strengths, determine which of the effects prevail. When the surface atoms are small enough, the interactions within the cluster determine the symmetry of the pattern. In such a case, the substrate behaves similarly to a continuous surface, and the low-temperature pattern is a hexadic monolayer. When the sizes of the surface and cluster atoms are comparable, the low-temperature adsorbed geometry mimics the substrate symmetry. On a face-centered cubic surface, face-centered cubic monolayers or droplets are obtained.

New results for phase transitions from catastrophe theory

Catastrophe theory predicts that in certain limits universal relations should exist between barrier heights, curvatures and the positions of local maxima and minima on a potential or free energy surface. In the present work we investigate these relations for both first- and second-order phase transitions, revealing that the ideal ratios often hold quite well over a wide range of conditions. This elementary catastrophe theory is illustrated for the melting transition of an atomic cluster, the isotropic-to-nematic transition in a liquid crystal, and the ferromagnetic-to-paramagnetic phase transition in the two-dimensional Ising model.

Liquid-solid transition in nuclei of protein crystals

It is generally assumed that crystallization begins with a small, crystalline nucleus. For proteins this paradigm may not be valid. Our numerical simulations show that under conditions typically used to produce protein crystals, small clusters of model proteins (particles with short-range, attractive interactions) cannot maintain a crystalline structure. Protein crystal nucleation is therefore an indirect, two-step process. A nucleus first forms and grows as a disordered, liquid-like aggregate. Once the aggregate grows beyond a critical size (about a few hundred particles) crystal nucleation becomes possible.

Transient backbending behavior in the Ising model with fixed magnetization

The physical origin of the backbendings in the equations of state of finite but not necessarily small systems is studied in the Ising model with fixed magnetization (IMFM) by means of the topological properties of the observable distributions and the analysis of the largest cluster with increasing lattice size. Looking at the convexity anomalies of the IMFM thermodynamic potential, it is shown that the order of the transition at the thermodynamic limit can be recognized in finite systems independently of the lattice size. General statistical mechanics arguments and analytical calculations suggest that the backbending in the caloric curve is a transient behavior which should not converge to a plateau in the thermodynamic limit, while the first-order transition (in the Ehrenfest sense) is still signaled by a discontinuity in the magnetization equation of state.

Approach to ergodicity in monte carlo simulations

The approach to the ergodic limit in Monte Carlo simulations is studied using both analytic and numerical methods. With the help of a stochastic model, a metric is defined that enables the examination of a simulation in both the ergodic and nonergodic regimes. In the nonergodic regime, the model implies how the simulation is expected to approach ergodic behavior analytically, and the analytically inferred decay law of the metric allows the monitoring of the onset of ergodic behavior. The metric is related to previously defined measures developed for molecular dynamics simulations, and the metric enables the comparison of the relative efficiencies of different Monte Carlo schemes. Applications to Lennard-Jones 13-particle clusters are shown to match the model for Metropolis, J-walking, and parallel tempering based approaches. The relative efficiencies of these three Monte Carlo approaches are compared, and the decay law is shown to be useful in determining needed high temperature parameters in parallel tempering and J-walking studies of atomic clusters.