Inferring free-energy barriers and kinetic rates from molecular dynamics via underdamped Langevin models.

Rare events include many of the most interesting transformation processes in condensed matter, from phase transitions to biomolecular conformational changes to chemical reactions. Access to the corresponding mechanisms, free-energy landscapes and kinetic rates can in principle be obtained by different techniques after projecting the high-dimensional atomic dynamics on one (or a few) collective variable. Even though it is well-known that the projected dynamics approximately follows - in a statistical sense - the generalized, underdamped or overdamped Langevin equations (depending on the time resolution), to date it is nontrivial to parameterize such equations starting from a limited, practically accessible amount of non-ergodic trajectories. In this work we focus on Markovian, underdamped Langevin equations, that arise naturally when considering, e.g., numerous water-solution processes at sub-picosecond resolution. After contrasting the advantages and pitfalls of different numerical approaches, we present an efficient parametrization strategy based on a limited set of molecular dynamics data, including equilibrium trajectories confined to minima and few hundreds transition path sampling-like trajectories. Employing velocity autocorrelation or memory kernel information for learning the friction and likelihood maximization for learning the free-energy landscape, we demonstrate the possibility to reconstruct accurate barriers and rates both for a benchmark system and for the interaction of carbon nanoparticles in water.

MetalWalls: Simulating electrochemical interfaces between polarizable electrolytes and metallic electrodes.

Electrochemistry is central to many applications, ranging from biology to energy science. Studies now involve a wide range of techniques, both experimental and theoretical. Modeling and simulations methods, such as density functional theory or molecular dynamics, provide key information on the structural and dynamic properties of the systems. Of particular importance are polarization effects of the electrode/electrolyte interface, which are difficult to simulate accurately. Here, we show how these electrostatic interactions are taken into account in the framework of the Ewald summation method. We discuss, in particular, the formal setup for calculations that enforce periodic boundary conditions in two directions, a geometry that more closely reflects the characteristics of typical electrolyte/electrode systems and presents some differences with respect to the more common case of periodic boundary conditions in three dimensions. These formal developments are implemented and tested in MetalWalls, a molecular dynamics software that captures the polarization of the electrolyte and allows the simulation of electrodes maintained at a constant potential. We also discuss the technical aspects involved in the calculation of two sets of coupled degrees of freedom, namely the induced dipoles and the electrode charges. We validate the implementation, first on simple systems, then on the well-known interface between graphite electrodes and a room-temperature ionic liquid. We finally illustrate the capabilities of MetalWalls by studying the adsorption of a complex functionalized electrolyte on a graphite electrode.

Anharmonic spectral features via trajectory-based quantum dynamics: A perturbative analysis of the interplay between dynamics and sampling.

The performance of different approximate algorithms for computing anharmonic features in vibrational spectra is analyzed and compared on model and more realistic systems that present relevant nuclear quantum effects. The methods considered combine approximate sampling of the quantum thermal distribution with classical time propagation and include Matsubara dynamics, path integral dynamics approaches, linearized initial value representation, and the recently introduced adaptive quantum thermal bath. A perturbative analysis of these different methods enables us to account for the observed numerical performance on prototypes for overtones and combination bands and to draw qualitatively correct trends for the numerical results obtained for Fermi resonances. Our results prove that the unequal performances of these approaches often derive from the method employed to sample initial conditions and not, as usually assumed, from the lack of coherence in the time propagation. Furthermore, as confirmed by the analysis reported in Benson and Althorpe, J. Chem. Phys. 130, 194510 (2021), we demonstrate, both via the perturbative approach and numerically, that path integral dynamics methods fail to reproduce the intensities of these anharmonic features and follow purely classical trends with respect to their temperature behavior. Finally, the remarkably accurate performance of the adaptive quantum thermal bath approach is documented and motivated.

Nuclear Quantum Effects in Liquid Water at Near Classical Computational Cost Using the Adaptive Quantum Thermal Bath.

We demonstrate the accuracy and efficiency of a recently introduced approach to account for nuclear quantum effects (NQEs) in molecular simulations: the adaptive quantum thermal bath (adQTB). In this method, zero-point energy is introduced through a generalized Langevin thermostat designed to precisely enforce the quantum fluctuation-dissipation theorem. We propose a refined adQTB algorithm with improved accuracy and report adQTB simulations of liquid water. Through extensive comparison with reference path integral calculations, we demonstrate that it provides excellent accuracy for a broad range of structural and thermodynamic observables as well as infrared vibrational spectra. The adQTB has a computational cost comparable to that of classical molecular dynamics, enabling simulations of up to millions of degrees of freedom.

Erratum: Fluctuation relations for systems in a constant magnetic field [Phys. Rev. E 102, 030101(R) (2020)].

This corrects the article DOI: 10.1103/PhysRevE.102.030101.

Fluctuation Relations for Dissipative Systems in Constant External Magnetic Field: Theory and Molecular Dynamics Simulations.

We illustrate how, contrary to common belief, transient Fluctuation Relations (FRs) for systems in constant external magnetic field hold without the inversion of the field. Building on previous work providing generalized time-reversal symmetries for systems in parallel external magnetic and electric fields, we observe that the standard proof of these important nonequilibrium properties can be fully reinstated in the presence of net dissipation. This generalizes recent results for the FRs in orthogonal fields-an interesting but less commonly investigated geometry-and enables direct comparison with existing literature. We also present for the first time a numerical demonstration of the validity of the transient FRs with nonzero magnetic field via nonequilibrium molecular dynamics simulations of a realistic model of liquid NaCl.

Fluctuation relations for systems in a constant magnetic field.

The validity of the fluctuation relations (FRs) for systems in a constant magnetic field is investigated. Recently introduced time-reversal symmetries that hold in the presence of static electric and magnetic fields and of deterministic thermostats are used to prove the transient FRs without invoking, as commonly done, inversion of the magnetic field. Steady-state FRs are also derived, under the t-mixing condition. These results extend the predictive power of important statistical mechanics relations. We illustrate this via the nonlinear response for the cumulants of the dissipation, showing how the alternative FRs enable one to determine analytically null cumulants also for systems in a single magnetic field.

Adiabatic motion and statistical mechanics via mass-zero constrained dynamics.

In recent work [Coretti et al., J. Chem. Phys., 2018, 149, 191102], a new algorithm to solve numerically the dynamics of the shell model for polarization was presented. The approach, broadly applicable to systems involving adiabatically separated dynamical variables, employs constrained molecular dynamics to strictly enforce the condition that the force on the fast degrees of freedom, modeled as having zero mass, is null at each time step. The algorithm is symplectic and fully time reversible, and results in stable and efficient propagation. In this paper we complete the discussion of the mechanics of mass-zero constrained dynamics by showing how to adapt it to problems where the fast degrees of freedom must satisfy additional conditions. This extension includes, in particular, the important case of first principles molecular dynamics. We then consider the statistical mechanics of the mass-zero constrained dynamical system demonstrating that the marginal probability sampled by the dynamics in the physical phase space recovers the form of the Born-Oppenheimer probability density. The effectiveness of the approach and the favorable scaling of the algorithm with system size are illustrated in test calculations of solid Na via orbital-free density functional dynamics.

Charge fluctuations from molecular simulations in the constant-potential ensemble.

We revisit the statistical mechanics of charge fluctuations in capacitors. In constant-potential classical molecular simulations, the atomic charges of electrode atoms are treated as additional degrees of freedom which evolve in time so as to satisfy the constraint of fixed electrostatic potential for each configuration of the electrolyte. The present work clarifies the role of the overall electroneutrality constraint, as well as the link between the averages computed within the Born-Oppenheimer approximation and that of the full constant-potential ensemble. This allows us in particular to derive a complete fluctuation-dissipation relation for the differential capacitance, that includes a contribution from the charge fluctuations (around the charges satisfying the constant-potential and electroneutrality constraints) also present in the absence of an electrolyte. We provide a simple expression for this contribution from the elements of the inverse of the matrix defining the quadratic form of the fluctuating charges in the energy. We then illustrate numerically the validity of our results, and recover the expected continuum result for an empty capacitor with structureless electrodes at large inter-electrode distances. By considering a variety of liquids between graphite electrodes, we confirm that this contribution to the total differential capacitance is small compared to that induced by the thermal fluctuations of the electrolyte.

Sampling the thermal Wigner density via a generalized Langevin dynamics.

The Wigner thermal density is a function of considerable interest in the area of approximate (linearized or semiclassical) quantum dynamics where it is employed to generate initial conditions for the propagation of appropriate sets of classical trajectories. In this paper, we propose an original approach to compute the Wigner density based on a generalized Langevin equation. The stochastic dynamics is nontrivial in that it contains a coordinate-dependent friction coefficient and a generalized force that couples momenta and coordinates. These quantities are, in general, not known analytically and have to be estimated via auxiliary calculations. The performance of the new sampling scheme is tested on standard model systems with highly nonclassical features such as relevant zero point energy effects, correlation between momenta and coordinates, and negative parts of the Wigner density. In its current brute force implementation, the algorithm, whose convergence can be systematically checked, is accurate and has only limited overhead compared to schemes with similar characteristics. We briefly discuss potential ways to further improve its numerical efficiency.

The Fluctuation-Dissipation Theorem as a Diagnosis and Cure for Zero-Point Energy Leakage in Quantum Thermal Bath Simulations.

Quantum thermal bath (QTB) simulations reproduce statistical nuclear quantum effects via a Langevin equation with a colored random force. Although this approach has proven efficient for a variety of chemical and condensed-matter problems, the QTB, as many other semiclassical methods, suffers from zero-point energy leakage (ZPEL). The absence of a reliable criterion to quantify the ZPEL without resorting to demanding comparisons with path integral-based calculations has so far hindered the use of the QTB for the simulation of real systems. In this work, we establish a quantitative connection between ZPEL in the QTB framework and deviations from the quantum fluctuation-dissipation theorem (FDT) that can be monitored along the simulation. This provides a rigorous general criterion to detect and quantify the ZPEL without any a priori knowledge of the system under study. We then use this criterion to build an adaptive QTB method that strictly enforces the quantum FDT at all frequencies via an on-the-fly, spectrally resolved fine-tuning of the system-bath coupling coefficients. The validity of the adaptive approach is first demonstrated on a simple two-oscillator model. It is then applied to two more realistic problems: the description of the vibrational properties of a model aluminum crystal at low temperature and the simulation of the liquid-solid phase transition in a 13-atom neon cluster. In both systems, the standard QTB results are strongly altered by the ZPEL, which can be essentially eliminated using the adaptive approach.

Communication: Constrained molecular dynamics for polarizable models.

A new algorithm to solve numerically the evolution of empirical shell models of polarizable systems is presented. It employs constrained molecular dynamics to satisfy exactly, at each time step, the crucial condition that the gradient of the potential with respect to the shell degrees of freedom is null. The algorithm is efficient, stable, and, contrary to the available alternatives, it is symplectic and time reversible. A proof-of-principle calculation on a polarizable model for NaCl is presented to illustrate its properties in comparison with the current method, which employs a conjugate-gradient procedure to enforce the null gradient condition. The proposed algorithm is applicable to other cases where a minimum condition on a function of an auxiliary set of driven dynamical variables must be satisfied.

Thermal Diffusion in Binary Mixtures: Transient Behavior and Transport Coefficients from Equilibrium and Nonequilibrium Molecular Dynamics.

Equilibrium and nonequilibrium molecular dynamics simulations are combined to compute the full set of coefficients that appear in the phenomenological equations describing thermal transport in a binary mixture subject to a constant thermal gradient. The Dynamical Non-Equilibrium Molecular Dynamics approach (D-NEMD) is employed to obtain the microscopic time evolution of the density and temperature fields, together with that of the mass and energy fluxes. D-NEMD enables one to study not only the steady state, but also the evolution of the fields during the transient that follows the onset of the thermal gradient, up to the establishment of the steady state. This makes it possible to ensure that the system has indeed reached a stationary condition, and to analyze the transient mechanisms and time scales of the mass and energy transport. A local time averaging procedure is applied to each trajectory contributing to the calculation to improve the signal-to-noise ratio in the estimation of the fluxes and to obtain a clear signal with the, relatively limited, statistics available.

Fermi resonance in CO2: Mode assignment and quantum nuclear effects from first principles molecular dynamics.

Vibrational spectroscopy is a fundamental tool to investigate local atomic arrangements and the effect of the environment, provided that the spectral features can be correctly assigned. This can be challenging in experiments and simulations when double peaks are present because they can have different origins. Fermi dyads are a common class of such doublets, stemming from the resonance of the fundamental excitation of a mode with the overtone of another. We present a new, efficient approach to unambiguously characterize Fermi resonances in density functional theory (DFT) based simulations of condensed phase systems. With it, the spectral features can be assigned and the two resonating modes identified. We also show how data from DFT simulations employing classical nuclear dynamics can be post-processed and combined with a perturbative quantum treatment at a finite temperature to include analytically thermal quantum nuclear effects. The inclusion of these effects is crucial to correct some of the qualitative failures of the Newtonian dynamics simulations at a low temperature such as, in particular, the behavior of the frequency splitting of the Fermi dyad. We show, by comparing with experimental data for the paradigmatic case of supercritical CO2, that these thermal quantum effects can be substantial even at ambient conditions and that our scheme provides an accurate and computationally convenient approach to account for them.

Time-reversal symmetry for systems in a constant external magnetic field.

The time-reversal properties of charged systems in a constant external magnetic field are reconsidered in this paper. We show that the evolution equations of the system are invariant under a new symmetry operation that implies a new signature property for time-correlation functions under time reversal. We then show how these findings can be combined with a previously identified symmetry to determine, for example, null components of the correlation functions of velocities and currents and of the associated transport coefficients. These theoretical predictions are illustrated by molecular dynamics simulations of superionic AgI.

Probabilistic Derivation of Spatiotemporal Correlation Functions in the Hydrodynamic Limit.

In this paper, we use probability theory to prove in suitable conditions the equivalence of equilibrium time correlation functions of microscopic density fields with the time correlation functions of local macroscopic density fields evolved by hydrodynamics in (approximate) phenomenological continuum theories of matter. We further discuss a useful and rigorous numerical algorithm, derived from this framework, to compute macroscopic space- and time-dependent behaviors (such as the hydrodynamical one) via molecular dynamics simulations.

Photophysics and photochemistry of a DNA-protein cross-linking model: a synergistic approach combining experiments and theory.

The photophysical and photochemical properties of 5-benzyluracil and 5,6-benzyluracil, the latter produced by photocyclization of the former through irradiation with femtosecond UV laser pulses, are investigated both experimentally and theoretically. The absorption spectra of the two molecules are compared, and the principal electronic transitions involved are discussed, with particular emphasis on the perturbation induced on the two chromophore species (uracil and benzene) by their proximity. The photoproduct formation for different irradiation times was verified with high-performance liquid chromatography and nuclear magnetic resonance measurements. The steady-state fluorescence demonstrates that the fluorescence is a distinctive physical observable for detection and selective identification of 5- and 5,6-benzyluracil. The principal electronic decay paths of the two molecules, obtained through TDDFT calculations, explain the features observed in the emission spectra and the photoreactivity of 5-benzyluracil. The order of magnitude of the lifetime of the excited states is derived with steady-state fluorescence anisotropy measurements and rationalized with the help of the computational findings. Finally, the spectroscopic data collected are used to derive the photocyclization and fluorescence quantum yields. In obtaining a global picture of the photophysical and photochemical properties of the two molecules, our findings demonstrates that the use of 5-benzyluracil as a model system to study the proximity relations of a DNA base with a close-lying aromatic amino acid is valid at a local level since the main characteristics of the decay processes from the excited states of the uracil/thymine molecules remain almost unchanged in 5-benzyluracil, the main perturbation arising from the presence of the close-lying aromatic group.

Quantum dynamical structure factor of liquid neon via a quasiclassical symmetrized method.

We apply the phase integration method for quasiclassical quantum time correlation functions [M. Monteferrante, S. Bonella, and G. Ciccotti, Mol. Phys. 109, 3015 (2011)] to compute the dynamic structure factor of liquid neon. So far the method had been tested only on model systems. By comparing our results for neon with experiments and previous calculations, we demonstrate that the scheme is accurate and efficient also for a realistic model of a condensed phase system showing quantum behavior.

A statistical mechanics handbook for protein-ligand binding simulation.

In this work, the fundamental elements of statistical mechanics underlying the simulation of the protein-ligand binding process, such as statistical ensembles and the concept of microscopic estimators of macroscopic observables and free energy, are summarized in a self consistent fashion. Particular attention is then devoted to the introduction of some mathematical tools that are used in atomistic simulations aimed at estimating binding affinities and free energy profiles, and to the illustration of the origins of the difficulties encountered in this endeavor.

The quantum free energy barrier for hydrogen vacancy diffusion in Na3AlH6.

The path integral single sweep method is used to assess quantum effects on the free energy barrier for hydrogen vacancy diffusion in a defective Na(3)AlH(6) crystal. This process has been investigated via experiments and simulations due to its potential relevance in the H release mechanism in sodium alanates, prototypical materials for solid state hydrogen storage. Previous computational studies, which used density functional methods for the electronic structure, were restricted to a classical treatment of the nuclear degrees of freedom. We show that, although they do not change the qualitative picture of the process, nuclear quantum effects reduce the free energy barrier height by about 18% with respect to the classical calculation improving agreement with available neutron scattering data.