Saddle points and dynamics of Lennard-Jones clusters, solids, and supercooled liquids

Assessment of Anharmonicities in Clusters: Developing and Validating a Minimum-Information Partition Function

Precise thermodynamic calculations are essential for understanding the dynamics of cluster systems and new particle formation. However, the widely employed harmonic statistical mechanical approach often falls short in terms of accuracy. In this study, we present an improved statistical model that incorporates vibrational anharmonicity via a novel partition function that requires only one additional system-specific input parameter. In addition to considering vibrational aspects, we also account for anharmonicity related to the configurational space. The role of anharmonicities is thoroughly examined in the case of general clusters, where the complete sets of conformers, mechanically stable spatial arrangements, are known up to clusters composed of 14 monomers. By performing consistent Monte Carlo simulations on these systems, we benchmark the statistical model's efficacy in reproducing key thermodynamic properties (formation free energy and potential energy) in the classical limit. The model exhibits exceptional alignment with simulations, accurately reproducing free energies within a precision of 2kBT and reliably capturing cluster melting temperatures. Furthermore, we demonstrate the significance and applicability of the model by reproducing thermodynamic barriers in homogeneous gas-phase nucleation of larger clusters. The transferability of our developed approach extends to more complex molecular systems and bears relevance for atmospheric multicomponent clusters, in particular.

Heterogeneous-elasticity theory of instantaneous normal modes in liquids

Since decades, the concept of vibrational density of states in glasses has been mirrored in liquids by the instantaneous-normal-mode spectrum. In glasses instantaneous configurations are believed to be situated close to minima of the potential-energy hypersurface and all eigenvalues of the associated Hessian matrix are positive. In liquids this is no longer true, and modes corresponding to both positive and negative eigenvalues exist. The instantaneous-normal-mode spectrum has been numerically investigated in the past, and it has been demonstrated to bring important information on the liquid dynamics and transport properties. A systematic deeper theoretical understanding is now needed. Heterogeneous-elasticity theory has proven to be particularly successful in explaining many details of the low-frequency excitations in glasses, ranging from the thoroughly studied boson peak, to other anomalies related to the crossover between wave-like and random-matrix-like excitations. Here we present an extension of heterogeneous-elasticity theory to the liquid state, and show that the outcome of the theory agrees well to the results of extensive molecular-dynamics simulations of a model liquid at different temperatures. We find that the spectrum of eigenvalues [Formula: see text] has a sharp maximum close to (but not at) [Formula: see text], and decreases monotonically with [Formula: see text] on both its stable and unstable side. We show that the spectral shape strongly depends on temperature, being symmetric at high temperatures and becoming rather asymmetric at low temperatures, close to the dynamical critical temperature. Most importantly, we demonstrate that the theory naturally reproduces a surprising phenomenon, a zero-energy spectral singularity with a cusp-like character developing in the vibrational spectra upon cooling. This feature, known from a few previous numerical studies, has been generally overlooked in the past due to a misleading representation of the data. We provide a thorough analysis of this issue, based on both very accurate predictions of our theory, and computational studies of model liquid systems with extended size.

Thomson problem in the disk

We investigate the classical ground state of a large number of charges confined inside a disk and interacting via the Coulomb potential. By realizing the important role that the peripheral charges play in determining the lowest energy solutions, we have successfully implemented an algorithm that allows us to work with configurations with a desired number of border charges. This feature brings a consistent reduction in the computational complexity of the problem, thus simplifying the search of global minima of the energy. Additionally, we have implemented a divide and conquer approach which has allowed us to study configurations of size never reached before (the largest one corresponding to N=40886 charges). These last configurations, in particular, are seen to display an increasingly rich structure of topological defects as N gets larger.

Neural network atomistic potentials for global energy minima search in carbon clusters

The global energy optimization problem is an acute and important problem in chemistry. It is crucial to know the geometry of the lowest energy isomer (global minimum, GM) of a given compound for the evaluation of its chemical and physical properties. This problem is especially relevant for atomic clusters. Due to the exponential growth of the number of local minima geometries with the increase of the number of atoms in the cluster, it is important to find a computationally efficient and reliable method to navigate the energy landscape and locate a true global minima structure. Newly developed neural network (NN) atomistic potentials offer a numerically efficient and relatively accurate approach for molecular structure optimization. An important question that needs to be answered is "Can NN potentials, trained on a given set, represent the potential energy surface (PES) of a neighboring domain?". In this work, we tested the applicability of ANI-1ccx and ANI-nr NN atomistic potentials for the global minima optimization of carbon clusters Cn (n = 3-10). We showed that with the introduction of the cluster connectivity restriction and consequent DFT or ab initio calculations, ANI-1ccx and ANI-nr can be considered as robust PES pre-samplers that can capture the GM structure even for large clusters such as C20.

Discovery of two-dimensional binary nanoparticle superlattices using global Monte Carlo optimization

Binary nanoparticle (NP) superlattices exhibit distinct collective plasmonic, magnetic, optical, and electronic properties. Here, we computationally demonstrate how fluid-fluid interfaces could be used to self-assemble binary systems of NPs into 2D superlattices when the NP species exhibit different miscibility with the fluids forming the interface. We develop a basin-hopping Monte Carlo (BHMC) algorithm tailored for interface-trapped structures to rapidly determine the ground-state configuration of NPs, allowing us to explore the repertoire of binary NP architectures formed at the interface. By varying the NP size ratio, interparticle interaction strength, and difference in NP miscibility with the two fluids, we demonstrate the assembly of an array of exquisite 2D periodic architectures, including AB-, AB₂-, and AB₃-type monolayer superlattices as well as AB-, AB₂-, A₃B₅-, and A₄B₆-type bilayer superlattices. Our results suggest that the interfacial assembly approach could be a versatile platform for fabricating 2D colloidal superlattices with tunable structure and properties.

The Loss Surface of Deep Linear Networks Viewed Through the Algebraic Geometry Lens

By using the viewpoint of modern computational algebraic geometry, we explore properties of the optimization landscapes of deep linear neural network models. After providing clarification on the various definitions of "flat" minima, we show that the geometrically flat minima, which are merely artifacts of residual continuous symmetries of the deep linear networks, can be straightforwardly removed by a generalized L₂-regularization. Then, we establish upper bounds on the number of isolated stationary points of these networks with the help of algebraic geometry. Combining these upper bounds with a method in numerical algebraic geometry, we find all stationary points for modest depth and matrix size. We demonstrate that, in the presence of the non-zero regularization, deep linear networks can indeed possess local minima which are not global minima. Finally, we show that even though the number of stationary points increases as the number of neurons (regularization parameters) increases (decreases), higher index saddles are surprisingly rare.

Instanton theory for Fermi's golden rule and beyond

Instanton theory provides a semiclassical approximation for computing quantum tunnelling effects in complex molecular systems. It is typically applied to proton-transfer reactions for which the Born-Oppenheimer approximation is valid. However, many processes in physics, chemistry and biology, such as electron transfers, are non-adiabatic and are correctly described instead using Fermi's golden rule. In this work, we discuss how instanton theory can be generalized to treat these reactions in the golden-rule limit. We then extend the theory to treat fourth-order processes such as bridge-mediated electron transfer and apply the method to simulate an electron moving through a model system of three coupled quantum dots. By comparison with benchmark quantum calculations, we demonstrate that the instanton results are much more reliable than alternative approximations based on superexchange-mediated effective coupling or a classical sequential mechanism. This article is part of the theme issue 'Chemistry without the Born-Oppenheimer approximation'.

Energy Landscape of State-Specific Electronic Structure Theory

State-specific approximations can provide a more accurate representation of challenging electronic excitations by enabling relaxation of the electron density. While state-specific wave functions are known to be local minima or saddle points of the approximate energy, the global structure of the exact electronic energy remains largely unexplored. In this contribution, a geometric perspective on the exact electronic energy landscape is introduced. On the exact energy landscape, ground and excited states form stationary points constrained to the surface of a hypersphere, and the corresponding Hessian index increases at each excitation level. The connectivity between exact stationary points is investigated, and the square-magnitude of the exact energy gradient is shown to be directly proportional to the Hamiltonian variance. The minimal basis Hartree-Fock and excited-state mean-field representations of singlet H2 (STO-3G) are then used to explore how the exact energy landscape controls the existence and properties of state-specific approximations. In particular, approximate excited states correspond to constrained stationary points on the exact energy landscape, and their Hessian index also increases for higher energies. Finally, the properties of the exact energy are used to derive the structure of the variance optimization landscape and elucidate the challenges faced by variance optimization algorithms, including the presence of unphysical saddle points or maxima of the variance.

Revisiting the single-saddle model for the β-relaxation of supercooled liquids

The dynamics of glass-forming liquids display several outstanding features, such as two-step relaxation and dynamic heterogeneities, which are difficult to predict quantitatively from first principles. In this work, we revisit a simple theoretical model of the β-relaxation, i.e., the first step of the relaxation dynamics. The model, first introduced by Cavagna et al. [J. Phys. A: Math. Gen. 36, 10721 (2003)], describes the dynamics of the system in the neighborhood of a saddle point of the potential energy surface. We extend the model to account for density–density correlation functions and for the four-point dynamic susceptibility. We obtain analytical results for a simple schematic model, making contact with related results for p-spin models and with the predictions of inhomogeneous mode-coupling theory. Building on recent computational advances, we also explicitly compare the model predictions against overdamped Langevin dynamics simulations of a glass-forming liquid close to the mode-coupling crossover. The agreement is quantitative at the level of single-particle dynamic properties only up to the early β-regime. Due to its inherent harmonic approximation, however, the model is unable to predict the dynamics on the time scale relevant for structural relaxation. Nonetheless, our analysis suggests that the agreement with the simulations may be largely improved if the modes’ spatial localization is properly taken into account.

Modeling the instantaneous normal mode spectra of liquids as that of unstable elastic media

We present a theory, based on an unstable elasticity model with thermally fluctuating elastic constants, which explains the salient features of the instantaneous normal mode (INM) spectrum of a numerically simulated liquid over a large range of temperatures. The INM spectrum of a liquid is obtained in a molecular dynamics (MD) simulation by considering the harmonic eigenvalue spectrum of the potential energy, taken at a given instant (snapshot). Because the INM spectrum records the curvatures of the fluctuating potential landscape of the liquid, an understanding of this spectrum may pave the way toward understanding the liquid to glass transformation process.

We study the instantaneous normal mode (INM) spectrum of a simulated soft-sphere liquid at different equilibrium temperatures T. We find that the spectrum of eigenvalues ρ(λ) has a sharp maximum near (but not at) λ=0 and decreases monotonically with |λ| on both the stable and unstable sides of the spectrum. The spectral shape strongly depends on temperature. It is rather asymmetric at low temperatures (close to the dynamical critical temperature) and becomes symmetric at high temperatures. To explain these findings we present a mean-field theory for ρ(λ), which is based on a heterogeneous elasticity model, in which the local shear moduli exhibit spatial fluctuations, including negative values. We find good agreement between the simulation data and the model calculations, done with the help of the self-consistent Born approximation (SCBA), when we take the variance of the fluctuations to be proportional to the temperature T. More importantly, we find an empirical correlation of the positions of the maxima of ρ(λ) with the low-frequency exponent of the density of the vibrational modes of the glasses obtained by quenching to T=0 from the temperature T. We discuss the present findings in connection to the liquid to glass transformation and its precursor phenomena.

Spin Crossover of Thiophosgene via Multidimensional Heavy-Atom Quantum Tunneling

The spin-crossover reaction of thiophosgene has drawn broad attention from both experimenters and theoreticians as a prime example of radiationless intramolecular decay via intersystem crossing. Despite multiple attempts over 20 years, theoretical predictions have typically been orders of magnitude in error relative to the experimentally measured triplet lifetime. We address the T₁ → S₀ transition by the first application of semiclassical golden-rule instanton theory in conjunction with on-the-fly electronic-structure calculations based on multireference perturbation theory. Our first-principles approach provides excellent agreement with the experimental rates. This was only possible because instanton theory goes beyond previous methods by locating the optimal tunneling pathway in full dimensionality and thus captures "corner cutting" effects. Since the reaction is situated in the Marcus inverted regime, the tunneling mechanism can be interpreted in terms of two classical trajectories, one traveling forward and one backward in imaginary time, which are connected by particle-antiparticle creation and annihilation events. The calculated mechanism indicates that the spin crossover is sped up by many orders of magnitude due to multidimensional quantum tunneling of the carbon atom even at room temperature.

Solution landscapes of the diblock copolymer-homopolymer model under two-dimensional confinement

We investigate the solution landscapes of the confined diblock copolymer and homopolymer in two-dimensional domain by using the extended Ohta-Kawasaki model. The projection saddle dynamics method is developed to compute the saddle points with mass conservation and construct the solution landscape by coupling with downward and upward search algorithms. A variety of stationary solutions are identified and classified in the solution landscape, including Flower class, Mosaic class, Core-shell class, and Tai-chi class. The relationships between different stable states are shown by either transition pathways connected by index-1 saddle points or dynamical pathways connected by a high-index saddle point. The solution landscapes also demonstrate the symmetry-breaking phenomena, in which more solutions with high symmetry are found when the domain size increases.

Orbital Optimized Density Functional Theory for Electronic Excited States

Density functional theory (DFT) based modeling of electronic excited states is of importance for investigation of the photophysical/photochemical properties and spectroscopic characterization of large systems. The widely used linear response time-dependent DFT (TDDFT) approach is, however, not effective at modeling many types of excited states, including (but not limited to) charge-transfer states, doubly excited states, and core-level excitations. In this perspective, we discuss state-specific orbital optimized (OO) DFT approaches as an alterative to TDDFT for electronic excited states. We motivate the use of OO-DFT methods and discuss reasons behind their relatively restricted historical usage (vs TDDFT). We subsequently highlight modern developments that address these factors and allow efficient and reliable OO-DFT computations. Several successful applications of OO-DFT for challenging electronic excitations are also presented, indicating their practical efficacy. OO-DFT approaches are thus increasingly becoming a useful route for computing excited states of large chemical systems. We conclude by discussing the limitations and challenges still facing OO-DFT methods, as well as some potential avenues for addressing them.

Critical Point-Finding Methods Reveal Gradient-Flat Regions of Deep Network Losses

Despite the fact that the loss functions of deep neural networks are highly nonconvex, gradient-based optimization algorithms converge to approximately the same performance from many random initial points. One thread of work has focused on explaining this phenomenon by numerically characterizing the local curvature near critical points of the loss function, where the gradients are near zero. Such studies have reported that neural network losses enjoy a no-bad-local-minima property, in disagreement with more recent theoretical results. We report here that the methods used to find these putative critical points suffer from a bad local minima problem of their own: they often converge to or pass through regions where the gradient norm has a stationary point. We call these gradient-flat regions, since they arise when the gradient is approximately in the kernel of the Hessian, such that the loss is locally approximately linear, or flat, in the direction of the gradient. We describe how the presence of these regions necessitates care in both interpreting past results that claimed to find critical points of neural network losses and in designing second-order methods for optimizing neural networks.

Rare events and first passage time statistics from the energy landscape

We analyze the probability distribution of rare first passage times corresponding to transitions between product and reactant states in a kinetic transition network. The mean first passage times and the corresponding rate constants are analyzed in detail for two model landscapes and the double funnel landscape corresponding to an atomic cluster. Evaluation schemes based on eigendecomposition and kinetic path sampling, which both allow access to the first passage time distribution, are benchmarked against mean first passage times calculated using graph transformation. Numerical precision issues severely limit the useful temperature range for eigendecomposition, but kinetic path sampling is capable of extending the first passage time analysis to lower temperatures, where the kinetics of interest constitute rare events. We then investigate the influence of free energy based state regrouping schemes for the underlying network. Alternative formulations of the effective transition rates for a given regrouping are compared in detail to determine their numerical stability and capability to reproduce the true kinetics, including recent coarse-graining approaches that preserve occupancy cross correlation functions. We find that appropriate regrouping of states under the simplest local equilibrium approximation can provide reduced transition networks with useful accuracy at somewhat lower temperatures. Finally, a method is provided to systematically interpolate between the local equilibrium approximation and exact intergroup dynamics. Spectral analysis is applied to each grouping of states, employing a moment-based mode selection criterion to produce a reduced state space, which does not require any spectral gap to exist, but reduces to gap-based coarse graining as a special case. Implementations of the developed methods are freely available online.

Perspective: new insights from loss function landscapes of neural networks

We investigate the structure of the loss function landscape for neural networks subject to dataset mislabelling, increased training set diversity, and reduced node connectivity, using various techniques developed for energy landscape exploration. The benchmarking models are classification problems for atomic geometry optimisation and hand-written digit prediction. We consider the effect of varying the size of the atomic configuration space used to generate initial geometries and find that the number of stationary points increases rapidly with the size of the training configuration space. We introduce a measure of node locality to limit network connectivity and perturb permutational weight symmetry, and examine how this parameter affects the resulting landscapes. We find that highly-reduced systems have low capacity and exhibit landscapes with very few minima. On the other hand, small amounts of reduced connectivity can enhance network expressibility and can yield more complex landscapes. Investigating the effect of deliberate classification errors in the training data, we find that the variance in testing AUC, computed over a sample of minima, grows significantly with the training error, providing new insight into the role of the variance-bias trade-off when training under noise. Finally, we illustrate how the number of local minima for networks with two and three hidden layers, but a comparable number of variable edge weights, increases significantly with the number of layers, and as the number of training data decreases. This work helps shed further light on neural network loss landscapes and provides guidance for future work on neural network training and optimisation.

Construction of a Pathway Map on a Complicated Energy Landscape

How do we search for the entire family tree of possible intermediate states, without unwanted random guesses, starting from a stationary state on the energy landscape all the way down to energy minima? Here we introduce a general numerical method that constructs the pathway map, which guides our understanding of how a physical system moves on the energy landscape. The method identifies the transition state between energy minima and the energy barrier associated with such a state. As an example, we solve the Landau-de Gennes energy incorporating the Dirichlet boundary conditions to model a liquid crystal confined in a square box; we illustrate the basic concepts by examining the multiple stationary solutions and the connected pathway maps of the model.

Structural, electronic and magnetic properties of stoichiometric cobalt oxide clusters (CoO)nq (n=3−10,q=0,+1): A modified basin-hopping Monte Carlo algorithm with spin-polarized DFT

The structural, electronic and magnetic properties of the neutral and cationic cobalt oxide clusters (CoO)[Formula: see text] ([Formula: see text], [Formula: see text]) have been studied using a modified basin-hopping Monte Carlo (BHMC) algorithm refined by spin-polarized DFT. A systematic search of global minimum structures predicts new global minima of (CoO)[Formula: see text] and reproduced other minima that are in excellent agreement with previous works. For most low-spin and high-spin states, the structural transition from planar-like to compact structure occurs at (CoO)[Formula: see text], which is in contrast with the general notion that the structural changes at (CoO)[Formula: see text]. Supported by the results of the binding energy, second-order total energy difference, chemical hardness, chemical potential and HOMO-LUMO gap confirms the stability of (CoO)4. Results of the spin magnetic moments for the global minima show that (CoO)4 and (CoO)8 spin configurations exhibit a fully antiferromagnetic (AFM) ordering, while (CoO)9 spin displays the highest ferromagnetic (FM) ordering. Interestingly, elongation of Co–Co bond in (CoO)4 causes O being polarized by the neighboring Co atoms that accordingly follows the Goodenough-Kanamori-Anderson rule of FM super-exchange coupling for the Co-O-Co structural rearrangement to 90∘ ([Formula: see text] structure) in order to accommodate the spin magnetic ordering changes. This rearrangement is a result of the valence band being shifted away from the Fermi level to lower energy causing high population of the spin-up density of state and leading to the asymmetrical polarization of the whole (CoO)4 structure. As far as the dissociation energy surfaces are concerned, the first ever such surfaces are constructed corresponding to [Formula: see text], which identify a complete dissociation pathway linking the cationic and neutral clusters and finally confirm (CoO)[Formula: see text] as the most stable cluster compared to the rest.

Exploring Energy Landscapes

Recent advances in the potential energy landscapes approach are highlighted, including both theoretical and computational contributions. Treating the high dimensionality of molecular and condensed matter systems of contemporary interest is important for understanding how emergent properties are encoded in the landscape and for calculating these properties while faithfully representing barriers between different morphologies. The pathways characterized in full dimensionality, which are used to construct kinetic transition networks, may prove useful in guiding such calculations. The energy landscape perspective has also produced new procedures for structure prediction and analysis of thermodynamic properties. Basin-hopping global optimization, with alternative acceptance criteria and generalizations to multiple metric spaces, has been used to treat systems ranging from biomolecules to nanoalloy clusters and condensed matter. This review also illustrates how all this methodology, developed in the context of chemical physics, can be transferred to landscapes defined by cost functions associated with machine learning.

Loss surface of XOR artificial neural networks

Training an artificial neural network involves an optimization process over the landscape defined by the cost (loss) as a function of the network parameters. We explore these landscapes using optimization tools developed for potential energy landscapes in molecular science. The number of local minima and transition states (saddle points of index one), as well as the ratio of transition states to minima, grow rapidly with the number of nodes in the network. There is also a strong dependence on the regularization parameter, with the landscape becoming more convex (fewer minima) as the regularization term increases. We demonstrate that in our formulation, stationary points for networks with N_{h} hidden nodes, including the minimal network required to fit the XOR data, are also stationary points for networks with N_{h}+1 hidden nodes when all the weights involving the additional node are zero. Hence, smaller networks trained on XOR data are embedded in the landscapes of larger networks. Our results clarify certain aspects of the classification and sensitivity (to perturbations in the input data) of minima and saddle points for this system, and may provide insight into dropout and network compression.

From sticky-hard-sphere to Lennard-Jones-type clusters

A relation M_{SHS→LJ} between the set of nonisomorphic sticky-hard-sphere clusters M_{SHS} and the sets of local energy minima M_{LJ} of the (m,n)-Lennard-Jones potential V_{mn}^{LJ}(r)=ɛ/n-m[mr^{-n}-nr^{-m}] is established. The number of nonisomorphic stable clusters depends strongly and nontrivially on both m and n and increases exponentially with increasing cluster size N for N≳10. While the map from M_{SHS}→M_{SHS→LJ} is noninjective and nonsurjective, the number of Lennard-Jones structures missing from the map is relatively small for cluster sizes up to N=13, and most of the missing structures correspond to energetically unfavorable minima even for fairly low (m,n). Furthermore, even the softest Lennard-Jones potential predicts that the coordination of 13 spheres around a central sphere is problematic (the Gregory-Newton problem). A more realistic extended Lennard-Jones potential chosen from coupled-cluster calculations for a rare gas dimer leads to a substantial increase in the number of nonisomorphic clusters, even though the potential curve is very similar to a (6,12)-Lennard-Jones potential.

Structure, thermodynamics, and rearrangement mechanisms in gold clusters-insights from the energy landscapes framework

We consider finite-size and temperature effects on the structure of model AuN clusters (30 ≤ N ≤ 147) bound by the Gupta potential. Equilibrium behaviour is examined in the harmonic superposition approximation, and the size-dependent melting temperature is also bracketed using molecular dynamics simulations. We identify structural transitions between distinctly different morphologies, characterised by various defect features. Reentrant behaviour and trends with respect to cluster size and temperature are discussed in detail. For N = 55, 85, and 147 we visualise the topography of the underlying potential energy landscape using disconnectivity graphs, colour-coded by the cluster morphology; and we use discrete path sampling to characterise the rearrangement mechanisms between competing structures separated by high energy barriers (up to 1 eV). The fastest transition pathways generally involve metastable states with multiple fivefold disclinations and/or a high degree of amorphisation, indicative of melting. For N = 55 we find that reoptimising low-lying minima using density functional theory (DFT) alters their energetic ordering and produces a new putative global minimum at the DFT level; however, the equilibrium structure predicted by the Gupta potential at room temperature is consistent with previous experiments.

Pathways for diffusion in the potential energy landscape of the network glass former SiO2

We study the dynamical behaviour of a computer model for viscous silica, the archetypal strong glass former, and compare its diffusion mechanism with earlier studies of a fragile binary Lennard-Jones liquid. Three different methods of analysis are employed. First, the temperature and time scale dependence of the diffusion constant is analysed. Negative correlation of particle displacements influences transport properties in silica as well as in fragile liquids. We suggest that the difference between Arrhenius and super-Arrhenius diffusive behaviour results from competition between the correlation time scale and the caging time scale. Second, we analyse the dynamics using a geometrical definition of cage-breaking transitions that was proposed previously for fragile glass formers. We find that this definition accurately captures the bond rearrangement mechanisms that control transport in open network liquids, and reproduces the diffusion constants accurately at low temperatures. As the same method is applicable to both strong and fragile glass formers, we can compare correlation time scales in these two types of systems. We compare the time spent in chains of correlated cage breaks with the characteristic caging time and find that correlations in the fragile binary Lennard-Jones system persist for an order of magnitude longer than those in the strong silica system. We investigate the origin of the correlation behaviour by sampling the potential energy landscape for silica and comparing it with the binary Lennard-Jones model. We find no qualitative difference between the landscapes, but several metrics suggest that the landscape of the fragile liquid is rougher and more frustrated. Metabasins in silica are smaller than those in binary Lennard-Jones and contain fewer high-barrier processes. This difference probably leads to the observed separation of correlation and caging time scales.

Mean-field dynamic criticality and geometric transition in the Gaussian core model

We use molecular dynamics simulations to investigate dynamic heterogeneities and the potential energy landscape of the Gaussian core model (GCM). Despite the nearly Gaussian statistics of particles' displacements, the GCM exhibits giant dynamic heterogeneities close to the dynamic transition temperature. The divergence of the four-point susceptibility is quantitatively well described by the inhomogeneous version of the mode-coupling theory. Furthermore, the potential energy landscape of the GCM is characterized by large energy barriers, as expected from the lack of activated, hopping dynamics, and display features compatible with a geometric transition. These observations demonstrate that all major features of mean-field dynamic criticality can be observed in a physically sound, three-dimensional model.

Comparison of minimum-action and steepest-descent paths in gradient systems

On high-dimensional and complex potential energy surfaces, the identification of the most likely mechanism for the transition between local minima is a challenging task. Usually the steepest-descent path is used interchangeably with the minimum-energy path and is associated with the most likely path. Here we compare the meaning of the steepest-descent path in complex energy landscapes to the path integral formulation of a trajectory that minimizes the action functional for Brownian dynamics. In particular, for energy landscapes with bifurcation points and multiple minima and saddle points, there can be several steepest-descent paths associated with specific saddles that connect two predetermined states but largely differ from the path of maximum likelihood. The minimum-action path, however, additionally takes into account the scalar work along the trajectory. Minimizing the scalar work can be less ambiguous in the identification of the most likely path in different gradient systems. It can also be used to distinguish between multiple steepest-descent paths that connect reactant and product states. We illustrate that in systems with complex energy landscapes a careful assessment of the steepest-descent path is thus advisable. Here the evaluation of the action can provide valuable information on the analysis and description of the most likely path.

Computer simulations of glasses: the potential energy landscape

We review the current state of research on glasses, discussing the theoretical background and computational models employed to describe them. This article focuses on the use of the potential energy landscape (PEL) paradigm to account for the phenomenology of glassy systems, and the way in which it can be applied in simulations and the interpretation of their results. This article provides a broad overview of the rich phenomenology of glasses, followed by a summary of the theoretical frameworks developed to describe this phenomonology. We discuss the background of the PEL in detail, the onerous task of how to generate computer models of glasses, various methods of analysing numerical simulations, and the literature on the most commonly used model systems. Finally, we tackle the problem of how to distinguish a good glass former from a good crystal former from an analysis of the PEL. In summarising the state of the potential energy landscape picture, we develop the foundations for new theoretical methods that allow the ab initio prediction of the glass-forming ability of new materials by analysis of the PEL.

Exploring the potential energy landscape of the Thomson problem via Newton homotopies

Locating the stationary points of a real-valued multivariate potential energy function is an important problem in many areas of science. This task generally amounts to solving simultaneous nonlinear systems of equations. While there are several numerical methods that can find many or all stationary points, they each exhibit characteristic problems. Moreover, traditional methods tend to perform poorly near degenerate stationary points with additional zero Hessian eigenvalues. We propose an efficient and robust implementation of the Newton homotopy method, which is capable of quickly sampling a large number of stationary points of a wide range of indices, as well as degenerate stationary points. We demonstrate our approach by applying it to the Thomson problem. We also briefly discuss a possible connection between the present work and Smale's 7th problem.

Potential energy landscape and inherent dynamics of a hard-sphere fluid

Hard-sphere models exhibit many of the same kinds of supercooled-liquid behavior as more realistic models of liquids, but the highly nonanalytic character of their potentials makes it a challenge to think of that behavior in potential energy landscape terms. We show here that it is possible to calculate an important topological property of hard-sphere landscapes, the geodesic pathways through those landscapes, and to do so without artificially coarse-graining or softening the potential. We show, moreover, that the rapid growth of the lengths of those pathways with increasing packing fraction quantitatively predicts the precipitous decline in diffusion constants in a glass-forming hard-sphere mixture model. The geodesic paths themselves can be considered as defining the intrinsic dynamics of hard spheres, so it is also revealing to find that they (and therefore the features of the underlying potential energy landscape) correctly predict the occurrence of dynamic heterogeneity and nonzero values of the non-Gaussian parameter. The success of these landscape predictions for the dynamics of such a singular model emphasizes that there is more to potential energy landscapes than is revealed by looking at the minima and saddle points.

Certification and the potential energy landscape

Typically, there is no guarantee that a numerical approximation obtained using standard nonlinear equation solvers is indeed an actual solution, meaning that it lies in the quadratic convergence basin. Instead, it may lie only in the linear convergence basin, or even in a chaotic region, and hence not converge to the corresponding stationary point when further optimization is attempted. In some cases, these non-solutions could be misleading. Proving that a numerical approximation will quadratically converge to a stationary point is termed certification. In this report, we provide details of how Smale's α-theory can be used to certify numerically obtained stationary points of a potential energy landscape, providing a mathematical proof that the numerical approximation does indeed correspond to an actual stationary point, independent of the precision employed.

Potential energy landscape of the two-dimensional XY model: higher-index stationary points

The application of numerical techniques to the study of energy landscapes of large systems relies on sufficient sampling of the stationary points. Since the number of stationary points is believed to grow exponentially with system size, we can only sample a small fraction. We investigate the interplay between this restricted sample size and the physical features of the potential energy landscape for the two-dimensional XY model in the absence of disorder with up to N = 100 spins. Using an eigenvector-following technique, we numerically compute stationary points with a given Hessian index I for all possible values of I. We investigate the number of stationary points, their energy and index distributions, and other related quantities, with particular focus on the scaling with N. The results are used to test a number of conjectures and approximate analytic results for the general properties of energy landscapes.

An inversion-relaxation approach for sampling stationary points of spin model Hamiltonians

Sampling the stationary points of a complicated potential energy landscape is a challenging problem. Here, we introduce a sampling method based on relaxation from stationary points of the highest index of the Hessian matrix. We illustrate how this approach can find all the stationary points for potentials or Hamiltonians bounded from above, which includes a large class of important spin models, and we show that it is far more efficient than previous methods. For potentials unbounded from above, the relaxation part of the method is still efficient in finding minima and transition states, which are usually the primary focus of attention for atomistic systems.

Clusters of coarse-grained water molecules

Global optimization for molecular clusters can be significantly more difficult than for atomic clusters because of the coupling between orientational and translational degrees of freedom. A coarse-grained representation of the potential can reduce the complexity of this problem, while retaining the essential features of the intermolecular interactions. In this study, we use a basin-hopping algorithm to locate putative global minima for clusters of coarse-grained water molecules modeled using a monatomic water potential for cluster sizes 3 ≤ N ≤ 55. We characterize these structures and identify structural trends using ideas from graph theory. The agreement with atomistic results and experiment is rather patchy, which we attribute to the tetrahedral bias in the three-body potential that results in too few nearest neighbor contacts and premature emergence of bulk-like structure. In spite of this issue, the results offer further useful insight into the relationship between the structure of clusters and bulk phases, and the mathematical form of a widely used model potential.

Potential energy landscapes for the 2D XY model: minima, transition states, and pathways

We describe a numerical study of the potential energy landscape for the two-dimensional XY model (with no disorder), considering up to 100 spins and central processing unit and graphics processing unit implementations of local optimization, focusing on minima and saddles of index one (transition states). We examine both periodic and anti-periodic boundary conditions, and show that the number of stationary points located increases exponentially with increasing lattice size. The corresponding disconnectivity graphs exhibit funneled landscapes; the global minima are readily located because they exhibit relatively large basins of attraction compared to the higher energy minima as the lattice size increases.

A strategy to find minimal energy nanocluster structures

An unbiased strategy to search for the global and local minimal energy structures of free standing nanoclusters is presented. Our objectives are twofold: to find a diverse set of low lying local minima, as well as the global minimum. To do so, we use massively the fast inertial relaxation engine algorithm as an efficient local minimizer. This procedure turns out to be quite efficient to reach the global minimum, and also most of the local minima. We test the method with the Lennard-Jones (LJ) potential, for which an abundant literature does exist, and obtain novel results, which include a new local minimum for LJ13 , 10 new local minima for LJ14 , and thousands of new local minima for 15≤N≤65. Insights on how to choose the initial configurations, analyzing the effectiveness of the method in reaching low-energy structures, including the global minimum, are developed as a function of the number of atoms of the cluster. Also, a novel characterization of the potential energy surface, analyzing properties of the local minima basins, is provided. The procedure constitutes a promising tool to generate a diverse set of cluster conformations, both two- and three-dimensional, that can be used as an input for refinement by means of ab initio methods.

Communication: Certifying the potential energy landscape

It is highly desirable for numerical approximations to stationary points for a potential energy landscape to lie in the corresponding quadratic convergence basin. However, it is possible that an approximation may lie only in the linear convergence basin, or even in a chaotic region, and hence not converge to the actual stationary point when further optimization is attempted. Proving that a numerical approximation will quadratically converge to the associated stationary point is termed certification. Here, we apply Smale's α-theory to stationary points, providing a certification serving as a mathematical proof that the numerical approximation does indeed correspond to an actual stationary point, independent of the precision employed. As a practical example, employing recently developed certification algorithms, we show how the α-theory can be used to certify all the known minima and transition states of Lennard-Jones LJ(N) atomic clusters for N = 7, ..., 14.

Energy landscape of the finite-size spherical three-spin glass model

We study the three-spin spherical model with mean-field interactions and Gaussian random couplings. For moderate system sizes of up to 20 spins, we obtain all stationary points of the energy landscape by means of the numerical polynomial homotopy continuation method. On the basis of these stationary points, we analyze the complexity and other quantities related to the glass transition of the model and compare these finite-system quantities to their exact counterparts in the thermodynamic limit.