The determination of efficient collective variables is crucial to the success of many enhanced sampling methods. As inspired by previous discrimination approaches, we first collect a set of data from the different metastable basins. The data are then
projected with the help of a neural network into a low-dimensional manifold in which data from different basins are well discriminated. This is here guaranteed by imposing that the projected data follows a preassigned distribution. The collective variables thus obtained lead to an efficient sampling and often allow reducing the number of collective variables in a multi-basin scenario. We first check the validity of the method in two-state systems. We then move to multi-step chemical processes. In the latter case, at variance with previous approaches, one single collective variable suffices, leading not only to computational efficiency but to a very clear representation of the reaction free energy profile.
The development of enhanced sampling methods has greatly extended the scope of atomistic simulations, allowing long-time phenomena to be studied with accessible computational resources. Many such methods rely on the identification of an appropriate s
et of collective variables. These are meant to describe the systems modes that most slowly approach equilibrium. Once identified, the equilibration of these modes is accelerated by the enhanced sampling method of choice. An attractive way of determining the collective variables is to relate them to the eigenfunctions and eigenvalues of the transfer operator. Unfortunately, this requires knowing the long-term dynamics of the system beforehand, which is generally not available. However, we have recently shown that it is indeed possible to determine efficient collective variables starting from biased simulations. In this paper, we bring the power of machine learning and the efficiency of the recently developed on-the-fly probability enhanced sampling method to bear on this approach. The result is a powerful and robust algorithm that, given an initial enhanced sampling simulation performed with trial collective variables or generalized ensembles, extracts transfer operator eigenfunctions using a neural network ansatz and then accelerates them to promote sampling of rare events. To illustrate the generality of this approach we apply it to several systems, ranging from the conformational transition of a small molecule to the folding of a mini-protein and the study of materials crystallization.
We present an approach that extends the theory of targeted free energy perturbation (TFEP) to calculate free energy differences and free energy surfaces at an accurate quantum mechanical level of theory from a cheaper reference potential. The converg
ence is accelerated by a mapping function that increases the overlap between the target and the reference distributions. Building on recent work, we show that this map can be learned with a normalizing flow neural network, without requiring simulations with the expensive target potential but only a small number of single-point calculations, and, crucially, avoiding the systematic error that was found previously. We validate the method by numerically evaluating the free energy difference in a system with a double-well potential and by describing the free energy landscape of a simple chemical reaction in the gas phase.
We present a modified version of the nudged elastic band (NEB) algorithm to find minimum energy paths con-necting two known configurations. We show that replacing the harmonic band-energy term with a discretized version of the Onsager-Machlup action
leads to a NEB algorithm with adaptive spring lengths that automatically increase the resolution of the minimum energy path around the saddle point of the potential energy surface. The method has the same computational cost per optimization step of the standard NEB algorithm and does not introduce additional parameters. We present applications to the isomerization of alanine dipeptide, the elimination of hydrogen from ethane and the healing of a 5-77-5 defect in graphene.
The study of liquid-liquid phase transition has attracted considerable attention. One interesting example of such phenomenon is phosphorus for which the existence a first-order phase transition between a low density insulating molecular phase and a c
onducting polymeric phase has been experimentally established. In this paper, we model this transition by an ab-initio quality molecular dynamics simulation and explore a large portion of the liquid section of the phase diagram. We draw the liquid-liquid coexistence curve and determine that it terminates into a second-order critical point. Close to the critical point, large coupled structure and electronic structure fluctuations are observed.
Supersolid is a mysterious and puzzling state of matter whose possible existence has stirred a vigorous debate among physicists for over 60 years. Its elusive nature stems from the coexistence of two seemingly contradicting properties, long-range ord
er and superfluidity. We report computational evidence of a supersolid phase of deuterium under high pressure ($p >800$ GPa) and low temperature (T $<$ 1.0 K). In our simulations, that are based on bosonic path integral molecular dynamics, we observe a highly concerted exchange of atoms while the system preserves its crystalline order. The exchange processes are favoured by the soft core interactions between deuterium atoms that form a densely packed metallic solid. At the zero temperature limit, Bose-Einstein condensation is observed as the permutation probability of $N$ deuterium atoms approaches $1/N$ with a finite superfluid fraction. Our study provides concrete evidence for the existence of a supersolid phase in high-pressure deuterium and could provide insights on the future investigation of supersolid phases in real materials.
The sampling problem lies at the heart of atomistic simulations and over the years many different enhanced sampling methods have been suggested towards its solution. These methods are often grouped into two broad families. On the one hand methods suc
h as umbrella sampling and metadynamics that build a bias potential based on few order parameters or collective variables. On the other hand, tempering methods such as replica exchange that combine different thermodynamic ensembles in one single expanded ensemble. We instead adopt a unifying perspective, focusing on the target probability distribution sampled by the different methods. This allows us to introduce a new class of collective-variables-based bias potentials that can be used to sample any of the expanded ensembles normally sampled via replica exchange. We also provide a practical implementation, by properly adapting the iterative scheme of the recently developed on-the-fly probability enhanced sampling method [Invernizzi and Parrinello, J. Phys. Chem. Lett. 11.7 (2020)], which was originally introduced for metadynamics-like sampling. The resulting method is very general and can be used to achieve different types of enhanced sampling. It is also reliable and simple to use, since it presents only few and robust external parameters and has a straightforward reweighting scheme. Furthermore, it can be used with any number of parallel replicas. We show the versatility of our approach with applications to multicanonical and multithermal-multibaric simulations, thermodynamic integration, umbrella sampling, and combinations thereof.
Accurate prediction of a gas solubility in a liquid is crucial in many areas of chemistry, and a detailed understanding of the molecular mechanism of the gas solvation continues to be an active area of research. Here, we extend the idea of constant c
hemical potential molecular dynamics (C{mu}MD) approach to the calculation of the gas solubility in the liquid under constant gas chemical potential conditions. As a representative example, we utilize this method to calculate the isothermal solubility of carbon dioxide in water. Additionally, we provide microscopic insight into the mechanism of solvation that preferentially occurs in areas of the surface where the hydrogen network is broken.
We present an ab initio molecular dynamics (MD) investigation of the tautomeric equilibrium for aqueous solutions of glycine and acetone at realistic experimental conditions. Metadynamics is used to accelerate proton migration among tautomeric center
s. Due to the formation of complex water-ion structures involved the proton dynamics in the aqueous environment, standard enhanced sampling approaches may face severe limitations in providing a general description of the phenomenon. Recently, we developed a set of Collective Variables (CVs) designed to study protons transfer reactions in complex condensed systems [Grifoni et al. PNAS, 2019, 116(10), 4054-4057]. In this work we applied this approach to study proton dissociation dynamics leading to tautomeric interconversion of biologically and chemically relevant prototypical systems, namely glycine and acetone in water. Although relatively simple from a chemical point of view, the results show that even for these small systems complex reaction pathways and non-trivial conversion dynamics are observed. The generality of our method allows obtaining these results without providing any prior information on the dissociation dynamics but only the atomic species that can exchange protons in the process. Our results agree with literature estimates and demonstrate the general applicability of this method in the study of tautomeric reactions.
We present a method for performing path integral molecular dynamics (PIMD) simulations for fermions and address its sign problem. PIMD simulations are widely used for studying many-body quantum systems at thermal equilibrium. However, they assume tha
t the particles are distinguishable and neglect bosonic and fermionic exchange effects. Interacting fermions play a key role in many chemical and physical systems, such as electrons in quantum dots and ultracold trapped atoms. A direct sampling of the fermionic partition function is impossible using PIMD since its integrand is not positive definite. We show that PIMD simulations for fermions are feasible by employing our recently developed method for bosonic PIMD and reweighting the results to obtain fermionic expectation values. The approach is tested against path integral Monte Carlo (PIMC) simulations for up to 7 electrons in a two-dimensional quantum dot for a range of interaction strengths. However, like PIMC, the method suffers from the sign problem at low temperatures. We propose a simple approach for alleviating it by simulating an auxiliary system with a larger average sign and obtaining an upper bound to the energy of the original system using the Bogoliubov inequality. This allows fermions to be studied at temperatures lower than would otherwise have been feasible using PIMD, as demonstrated in the case of a three-electron quantum dot. Our results extend the boundaries of PIMD simulations of fermions and will hopefully stimulate the development of new approaches for tackling the sign problem.
