Molecular dynamics integration and molecular vibrational theory. III. The infrared spectrum of water

The new symplectic molecular dynamics (MD) integrators presented in the first paper of this series were applied to perform MD simulations of water. The physical properties of a system of flexible TIP3P water molecules computed by the new integrators, such as diffusion coefficients, orientation correlation times, and infrared (IR) spectra, are in good agreement with results obtained by the standard method. The comparison between the new integrators' and the standard method's integration time step sizes indicates that the resulting algorithm allows a 3.0 fs long integration time step as opposed to the standard leap-frog Verlet method, a sixfold simulation speed-up. The accuracy of the method was confirmed, in particular, by computing the IR spectrum of water in which no blueshifting of the stretching normal mode frequencies is observed as occurs with the standard method.     

Explicit symplectic integrators of molecular dynamics algorithms for rigid-body molecules in the canonical, isobaric-isothermal, and related ensembles

The authors propose explicit symplectic integrators of molecular dynamics (MD) algorithms for rigid-body molecules in the canonical and isobaric-isothermal ensembles. They also present a symplectic algorithm in the constant normal pressure and lateral surface area ensemble and that combined with the Parrinello-Rahman algorithm. Employing the symplectic integrators for MD algorithms, there is a conserved quantity which is close to Hamiltonian. Therefore, they can perform a MD simulation more stably than by conventional nonsymplectic algorithms. They applied this algorithm to a TIP3P pure water system at 300 K and compared the time evolution of the Hamiltonian with those by the nonsymplectic algorithms. They found that the Hamiltonian was conserved well by the symplectic algorithm even for a time step of 4 fs. This time step is longer than typical values of 0.5-2 fs which are used by the conventional nonsymplectic algorithms.     

Molecular dynamics integration and molecular vibrational theory. I. New symplectic integrators

New symplectic integrators have been developed by combining molecular dynamics integration with the standard theory of molecular vibrations to solve the Hamiltonian equations of motion. The presented integrators analytically resolve the internal high-frequency molecular vibrations by introducing a translating and rotating internal coordinate system of a molecule and calculating normal modes of an isolated molecule only. The translation and rotation of a molecule are treated as vibrational motions with the vibrational frequency zero. All types of motion are thus described in terms of the normal coordinates. The method's time reversibility requirement was used to determine the equations of motion for internal coordinate system of a molecule. The calculation of long-range forces is performed numerically within the generalized second-order leap-frog scheme, in the same way as in standard second-order symplectic methods. The new methods for integrating classical equations of motion using normal mode analysis allow us to use a long integration step and are applicable to any system of molecules with one equilibrium configuration.     

Molecular dynamics integration meets standard theory of molecular vibrations

An iterative SISM (split integration symplectic method) for molecular dynamics (MD) integration is described. This work explores an alternative for the internal coordinate system prediction in the SISM introduced by JaneZic et al. (J. Chem. Phys. 2005, 122, 174101). The SISM, which employs a standard theory of molecular vibrations, analytically resolves the internal high-frequency molecular vibrations. This is accomplished by introducing a translating and rotating internal coordinate system of a molecule and calculating normal modes of an isolated molecule only. The Eckart frame, which is usually used in the standard theory of molecular vibrations as an internal coordinate system of a molecule, is adopted to be used within the framework of the second order generalized leapfrog scheme. In the presented MD integrator the internal coordinate frame at the end of the integration step is predicted halfway through the integration step using a predictor-corrector type iterative approach thus ensuring the method's time reversibility. The iterative SISM, which is applicable to any system of molecules with one equilibrium configuration, was applied here to perform all-atom MD simulations of liquid CO2 and SO2. The simulation results indicate that for the same level of accuracy, this algorithm allows significantly longer integration time steps than the standard second-order leapfrog Verlet (LFV) method.     

Symplectic molecular dynamics simulations on specially designed parallel computers

We have developed a computer program for molecular dynamics (MD) simulation that implements the Split Integration Symplectic Method (SISM) and is designed to run on specialized parallel computers. The MD integration is performed by the SISM, which analytically treats high-frequency vibrational motion and thus enables the use of longer simulation time steps. The low-frequency motion is treated numerically on specially designed parallel computers, which decreases the computational time of each simulation time step. The combination of these approaches means that less time is required and fewer steps are needed and so enables fast MD simulations. We study the computational performance of MD simulation of molecular systems on specialized computers and provide a comparison to standard personal computers. The combination of the SISM with two specialized parallel computers is an effective way to increase the speed of MD simulations up to 16-fold over a single PC processor.     

Recent advances in polymer molecular dynamics simulation and data analysis

Significant advances in molecular simulation methodology over the past decade have greatly reduced the traditional size-timescale bottleneck in molecular dynamics calculations. The development of the geometric statement function method allows for systems up to several hundred thousand atoms to be simulated for up to several nanoseconds in reasonable times on standard workstations. For constant energy simulations, the use of symplectic integrators ensures accurate dynamics, even at long simulation times, without velocity or other artificial rescaling schemes. Finally, new methods of frequency estimation allow for accurate vibrational mode frequency calculations even in the presence of chaotic motion on time scales twenty times shorter than the standard fast Fourier transform, with an additional improvement in the sensitivity of the results when initial dynamics conditions are carefully chosen.     

Time reversible and symplectic integrators for molecular dynamics simulations of rigid molecules

Molecular dynamics integrators are presented for translational and rotational motion of rigid molecules in microcanonical, canonical, and isothermal-isobaric ensembles. The integrators are all time reversible and are also, in some approaches, symplectic for the microcanonical ensembles. They are developed utilizing the quaternion representation on the basis of the Trotter factorization scheme using a Hamiltonian formalism. The structure is similar to that of the velocity Verlet algorithm. Comparison is made with standard integrators in terms of stability and it is found that a larger time step is stable with the new integrators. The canonical and isothermal-isobaric molecular dynamics simulations are defined by using a chain thermostat approach according to generalized Nosé-Hoover and Andersen methods.     

New open modified trigonometrically-fitted Newton-Cotes type multilayer symplectic integrators for the numerical solution of the Schrödinger equation

In this paper we present new modified open Newton Cotes integrators and we develop a new modified trigonometrically-fitted open Newton-Cotes method. We study the connection between the new proposed schemes, the differential methods and the symplectic integrators. although The research on multistep symplectic integrators is very poor, although, much research has been done on one step symplectic integrators and several of then have obtained based on symplectic geometry. In this paper a new open modified numerical algorithm of Newton-Cotes type is produced. We present the new obtained method as symplectic multilayer integrator. The new obtained symplectic schemes are applied for the solution of the resonance problem of the radial Schrödinger Equation. The results show the efficiency of the new proposed algorithm.     

Symplectic splitting operator methods for the time-dependent Schrodinger equation

We present a family of symplectic splitting methods especially tailored to solve numerically the time-dependent Schrodinger equation. When discretized in time, this equation can be recast in the form of a classical Hamiltonian system with a Hamiltonian function corresponding to a generalized high-dimensional separable harmonic oscillator. The structure of the system allows us to build highly efficient symplectic integrators at any order. The new methods are accurate, easy to implement, and very stable in comparison with other standard symplectic integrators.     

Advanced gradient-like methods for rigid-body molecular dynamics

A novel approach is developed to integrate the equations of motion in many-body systems of interacting rigid polyatomic molecules. It is based on an advanced gradient-like decomposition technique in the presence of translational and orientational degrees of freedom. As a result, a new class of reversible phase-space volume preserving fourth-order algorithms for rotational motion is introduced. Contrary to standard nongradient decomposition integrators, the algorithms derived take into account additional analytically integrable terms in the exponential propagators, while the arising gradients are expressed in terms of forces and torques. This allows one to increase significantly the precision of the integration and, at the same time, reduce the increased computational costs. The optimized second-order integrator is also presented. The gradient-like and optimized algorithms are tested in molecular dynamics simulations of water versus well-established integrators known previously. It is shown that the new algorithms lead to the best efficiency in the rigid-body integration.     

Improvement of methods for molecular dynamics integration

An application of symplectic implicit Runge–Kutta (RK ) integration schemes, the s-stage Gauss–Legendre Runge–Kutta (GLRK ) methods of order 2s, for the numerical solution of molecular dynamics (MD ) equation is described. The two-stage fourth-order GLRK method, the implicit midpoint rule, and the three-stage diagonally implicit RK method of order four are studied. The fixed-point iteraction was used for solving the resulting nonlinear system of equations. The algorithms were applied to a complex system of N particles interacting through a Lennard-Jones potential. The proposed symplectic methods for MD integration permit a wide range of time steps, are highly accurate and stable, and are thus suitable for the MD integration. © 1994 John Wiley & Sons, Inc.     

ORAC: A Molecular dynamics program to simulate complex molecular systems with realistic electrostatic interactions

In this study, we present a new molecular dynamics program for simulation of complex molecular systems. The program, named ORAC, combines state-of-the-art molecular dynamics (MD) algorithms with flexibility in handling different types and sizes of molecules. ORAC is intended for simulations of molecular systems and is specifically designed to treat biomolecules efficiently and effectively in solution or in a crystalline environment. Among its unique features are: (i) implementation of reversible and symplectic multiple time step algorithms (or r-RESPA, reversible reference system propagation algorithm) specifically designed and tuned for biological systems with periodic boundary conditions; (ii) availability for simulations with multiple or single time steps of standard Ewald or smooth particle mesh Ewald (SPME) for computation of electrostatic interactions; and (iii) possibility of simulating molecular systems in a variety of thermodynamic ensembles. We believe that the combination of these algorithms makes ORAC more advanced than other MD programs using standard simulation algorithms. © 1997 John Wiley & Sons, Inc. J Comput Chem 18 : 1848–1862, 1997     

High order closed Newton-Cotes exponentially and trigonometrically fitted formulae as multilayer symplectic integrators and their application to the radial Schr?dinger equation

In this paper we study the connection between: (i) closed Newton-Cotes formulae of high order, (ii) trigonometrically-fitted and exponentially-fitted differential methods, (iii) symplectic integrators. Several one step symplectic integrators have been produced based on symplectic geometry during the last decades (see the relevant literature and the references here). However, the study of multistep symplectic integrators is very poor. In this paper we investigate the High Order Closed Newton-Cotes Formulae and we write them as symplectic multilayer structures. We develop trigonometrically-fitted and exponentially-fitted symplectic methods which are based on the closed Newton-Cotes formulae. We apply the symplectic schemes in order to solve the resonance problem of the radial Schr?dinger equation. Based on the theoretical and numerical results, conclusions on the efficiency of the new obtained methods are given.     

Application of symplectic algorithms to QCT calculation:H+H_2 system

The basic theory of symplectic algorithm was introduced. A comparison between Runge-Kutta method and symplectic integration method was preformed in the simulation of the long time behavior of H + H2 system on BKMP potential energy surface. Our results reveal a dis-sipative behavior in the integral of ordinary differential equation by the fourth order Runge-Kutta method, which causes incorrect simulation results in QCT calculations. However, when the symplectic integration method is applied, the dissipative behavior is not found in the same system. When the initial state is the same, the energy deviation of fourth order symplectic integral method is almost one percent of that of fourth order Runge-Kutta method in a 60000-step simulation, and that of sixth order symplectic integral method is much less. These results show that the symplectic integral methods are always the better choice in the integral calculation of the long time behavior in maintaining energy conservation.     

Closed Newton–Cotes trigonometrically-fitted formulae of high order for the numerical integration of the Schrödinger equation

In this paper we investigate the connection between (i) closed Newton–Cotes formulae, (ii) trigonometrically-fitted differential methods, (iii) symplectic integrators and (iv) efficient solution of the Schrödinger equation. In the last decades several one step symplectic integrators have been produced based on symplectic geometry, (see the relevant literature and the references here). However, the study of multistep symplectic integrators is very poor. In this paper we investigate the closed Newton–Cotes formulae and we write them as symplectic multilayer structures. We also develop trigonometrically-fitted symplectic methods which are based on the closed Newton–Cotes formulae. We apply the symplectic schemes to the well known radial Schrödinger equation in order to investigate the efficiency of the proposed method to these type of problems.     

Symplectic integrators and the conservation of angular momentum

In this article we observe that generally symplectic integrators conserve angular momentum exactly, whereas nonsymplectic integrators do not. We show that this observation extends to multiple timesteps and to constrained dynamics. Both of these devices are important for efficient molecular dynamics simulations. © 1995 by John Wiley & Sons, Inc.     

Application of symplectic algorithms to QCT calculation: H + H2 system

The basic theory of symplectic algorithm was introduced. A comparison between Runge-Kutta method and symplectic integration method was preformed in the simulation of the long time behavior of H + H₂ system on BKMP potential energy surface. Our results reveal a dissipative behavior in the integral of ordinary differential equation by the fourth order Runge-Kutta method, which causes incorrect simulation results in QCT calculations. However, when the symplectic integration method is applied, the dissipative behavior is not found in the same system. When the initial state is the same, the energy deviation of fourth order symplectic integral method is almost one percent of that of fourth order Runge-Kutta method in a 60000-step simulation, and that of sixth order symplectic integral method is much less. These results show that the symplectic integral methods are always the better choice in the integral calculation of the long time behavior in maintaining energy conservation.     

Application of symplectic algorithms to QCT calculation: H + H2 system

The basic theory of symplectic algorithm was introduced. A comparison between Runge-Kutta method and symplectic integration method was preformed in the simulation of the long time behavior of H + H₂ system on BKMP potential energy surface. Our results reveal a dissipative behavior in the integral of ordinary differential equation by the fourth order Runge-Kutta method, which causes incorrect simulation results in QCT calculations. However, when the symplectic integration method is applied, the dissipative behavior is not found in the same system. When the initial state is the same, the energy deviation of fourth order symplectic integral method is almost one percent of that of fourth order Runge-Kutta method in a 60000-step simulation, and that of sixth order symplectic integral method is much less. These results show that the symplectic integral methods are always the better choice in the integral calculation of the long time behavior in maintaining energy conservation.     

Molecular latent space simulators

Small integration time steps limit molecular dynamics (MD) simulations to millisecond time scales. Markov state models (MSMs) and equation-free approaches learn low-dimensional kinetic models from MD simulation data by performing configurational or dynamical coarse-graining of the state space. The learned kinetic models enable the efficient generation of dynamical trajectories over vastly longer time scales than are accessible by MD, but the discretization of configurational space and/or absence of a means to reconstruct molecular configurations precludes the generation of continuous atomistic molecular trajectories. We propose latent space simulators (LSS) to learn kinetic models for continuous atomistic simulation trajectories by training three deep learning networks to (i) learn the slow collective variables of the molecular system, (ii) propagate the system dynamics within this slow latent space, and (iii) generatively reconstruct molecular configurations. We demonstrate the approach in an application to Trp-cage miniprotein to produce novel ultra-long synthetic folding trajectories that accurately reproduce atomistic molecular structure, thermodynamics, and kinetics at six orders of magnitude lower cost than MD. The dramatically lower cost of trajectory generation enables greatly improved sampling and greatly reduced statistical uncertainties in estimated thermodynamic averages and kinetic rates.

Latent space simulators learn kinetic models for atomistic simulations and generate novel trajectories at six orders of magnitude lower cost.     

An efficient method for computing excess free energy of liquid

We present a new theoretical method for efficient calculation of free energy of liquid. This interaction entropy method allows one to compute entropy and free energy of liquid from standard single step MD (molecular dynamics) simulation directly in liquid state without the need to perform MD simulations at many intermediate states as required in thermodynamic integration or free energy perturbation methods. In this new approach, one only needs to evaluate the interaction energy of a single (fixed) liquid molecule with the rest of liquid molecules as a function of time from a standard MD simulation of liquid and the fluctuation of distribution of this interaction energy is then used to calculate the interaction entropy of the liquid. Explicit theoretical derivation of this interaction entropy approach is provided and numerical calculations for the benchmark liquid water system were carried out using three different water models. Numerical analysis of the result was performed and comparison of the computational result with experimental data and other theoretical results were provided. Excellent agreement of calculated free energies with the experimental data using TIP4P model is obtained for liquid water.