Nose-Hoover Chain Research Articles

Enhanced H2 adsorption on 2D B4N monolayer modified with Na atoms: A density functional theory exploration

In this study, H2 adsorption properties on two-dimensional (2D) B4N monolayer modified with Na atoms are comprehensively investigated via density functional theory (DFT). Our results expose the intrinsic 2D B4N monolayer barely adsorb H2 molecules owes to a weak H2 adsorption energy. However, Na atom can be strongly anchored on the surface of 2D B4N monolayer with a large binding energy, which results in that the metal clustering is prevented. Meanwhile, H2 adsorption energy is drastically enhanced by introduced Na atom on B4N monolayer. Moreover, Na/B4N complex can adsorb up to 10 H2 molecules with average adsorption energy of −0.109 eV/H2 while 2Na/B4N complex can adsorb up to 20 H2 molecules with average adsorption energy of −0.116 eV/H2. The obtained H2 adsorption energy falls into an optimal range between −0.1 eV/H2 and −0.8 eV/H2 for mobile hydrogen storage applications. Also, the hydrogen storage gravimetric density reaches 7.27 wt% and exceeds the latest objective of 6.5 wt% from U.S. Department of Energy (DOE). In addition, Mulliken atomic charges analysis, charges density difference and projected density of states systematically reveal H2 adsorption mechanisms mainly from polarization effect and orbital hybridization. Ab initio molecular dynamic (AIMD) simulations using Nose-Hoover (NH) chain thermostat under canonical (NVT) ensemble manifest H2 molecules can be easily released from the 2Na/B4N complex under indoor temperature at 300 K. Our theoretical study guarantee 2D B4N monolayer modified with Na atoms can be recognized as a prospective hydrogen storage material.

Unification of nonequilibrium molecular dynamics and the mode-resolved phonon Boltzmann equation for thermal transport simulations

Nano-size confinement induces many intriguing non-Fourier heat conduction phenomena, such as nonlinear temperature gradients, temperature jumps near the contacts, and size-dependent thermal conductivity. Over the past decades, these effects have been studied and interpreted by nonequilibrium molecular dynamics (NEMD) and phonon Boltzmann transport equation (BTE) simulations separately, but no theory that unifies these two methods has ever been established. In this work, we unify these methods using a quantitative mode-level comparison and demonstrate that they are equivalent for various thermostats. We show that different thermostats result in different non-Fourier thermal transport characteristics due to the different mode-level phonon excitations inside the thermostats, which explains the different size-dependent thermal conductivities calculated using different reservoirs, even though they give the same bulk thermal conductivity. Specifically, the Langevin thermostat behaves like a thermalizing boundary in phonon BTE and provides mode-level thermal-equilibrium phonon outlets, while the Nose-Hoover chain thermostat and velocity rescaling method behave like biased reservoirs, which provide a spatially uniform heat generation and mode-level nonequilibrium phonon outlets. These findings explain why different experimental measurement methods can yield different size-dependent thermal conductivity. They also indicate that the thermal conductivity of materials can be tuned for various applications by specifically designing thermostats. The unification of NEMD and phonon BTE will largely facilitate the study of thermal transport in complex systems in the future by, e.g., replacing computationally unaffordable first-principles NEMD simulations with computationally less expensive spectral BTE simulations.

The melting limit in sodium clusters

Thermodynamic properties of the small sodium clusters $$\hbox {Na}_6,\, \hbox {Na}_8$$ and $$\hbox {Na}_{10}$$ have been studied by Born–Oppenheimer molecular dynamics (BOMD) simulations. The simulations were performed with auxiliary density functional theory as implemented in the deMon2k code. This approach has already proved accurate for the calculations of thermodynamic properties of larger sodium clusters. The Nose–Hoover chain thermostat was applied to control the temperature. BOMD simulations were performed in the temperature range from 250 to 1000 K. The obtained trajectories were analyzed using the multiple-histogram method in order to obtain continuous functions for the energies and heat capacities. For the $$\hbox {Na}_6$$ and $$\hbox {Na}_8$$ clusters, besides the fragmentation of the clusters at higher temperature, no other characteristic features in the heat capacity curves are found. On the other hand, a small peak at low temperature was found in the $$\hbox {Na}_{10}$$ heat capacity curve which is characteristic for molecular melting. Our analysis of the $$\hbox {Na}_{10}$$ melting shows that electronic structure parameters are better suited than geometrical ones to describe the melting process due to the fluctional nature of the clusters. We find that energetical resorting of the occupied cluster orbitals is characteristic for the $$\hbox {Na}_{10}$$ cluster melting.

Conserving the linear momentum in stochastic dynamics: Dissipative particle dynamics as a general strategy to achieve local thermostatization in molecular dynamics simulations.

Stochastic dynamics is a widely employed strategy to achieve local thermostatization in molecular dynamics simulation studies; however, it suffers from an inherent violation of momentum conservation. Although this short-coming has little impact on structural and short-time dynamic properties, it can be shown that dynamics in the long-time limit such as diffusion is strongly dependent on the respective thermostat setting. Application of the methodically similar dissipative particle dynamics (DPD) provides a simple, effective strategy to ensure the advantages of local, stochastic thermostatization while at the same time the linear momentum of the system remains conserved. In this work, the key parameters to employ the DPD thermostats in the framework of periodic boundary conditions are investigated, in particular the dependence of the system properties on the size of the DPD-region as well as the treatment of forces near the cutoff. Structural and dynamical data for light and heavy water as well as a Lennard-Jones fluid have been compared to simulations executed via stochastic dynamics as well as via use of the widely employed Nose-Hoover chain and Berendsen thermostats. It is demonstrated that a small size of the DPD region is sufficient to achieve local thermalization, while at the same time artifacts in the self-diffusion characteristic for stochastic dynamics are eliminated. © 2016 Wiley Periodicals, Inc.

On the configurational temperature Nosè–Hoover thermostat

Abstract In this paper we reformulate the configurational temperature Nose–Hoover thermostat of Braga and Travis (2005) by means of a quasi-Hamiltonian theory in phase space Sergi and Ferrario (2001). The quasi-Hamiltonian structure is exploited to introduce a hybrid configurational-kinetic temperature Nose–Hoover chain thermostat that can achieve a uniform sampling of phase space (also for stiff harmonic systems), as illustrated by simulating the dynamics of one-dimensional harmonic and quartic oscillators. An integration algorithm, based on the symmetric Trotter decomposition of the propagator, is presented and tested against implicit geometric algorithms with a structure similar to the velocity and position Verlet. In order to obtain an explicit form for the symmetric Trotter propagator algorithm, in the case of non-harmonic and non-linear interaction potentials, a position-dependent harmonically approximated propagator is introduced. Such a propagator approximates the dynamics of the configurational degrees of freedom as if they were locally moving in a harmonic potential. The resulting approximated locally harmonic dynamics is tested with good results in the case of a one-dimensional quartic oscillator: The integration is stable and locally time-reversible. Instead, the implicit geometric integrator is stable and time-reversible globally (when convergence is achieved). We also verify the stability of the approximated explicit integrator for a three-dimensional N-particle system interacting through a soft Weeks–Chandler–Andersen potential.

Solvation Structure and Quasidynamics of Biomolecules Steered with Effective Solvation Forces Obtained From Molecular Theory of Solvation

Three-dimensional reference interaction site model with Kovalenko-Hirata closure (3D-RISM-KH molecular theory of solvation) is Ornstein-Zernike type integral equation theory of liquids. Based on first principles of statistical mechanics, 3D-RISM-KH consistently accounts for effect of chemical specificities of solvent, co-solvent, ions, and ligands on biomolecules solvation structure and thermodynamics, including steric forces, hydrophobicity, hydrogen bonding, and other effective interactions. 3D-RISM-KH supplemented with the partial molar volume correction, a.k.a. “Universal Correction” (UC), provides an excellent agreement with experimental data on a large set of small compounds for solvation free energies in octanol and water, and octanol-water partition coefficients. 3D-RISM-KH structural water detection and placement are implemented in Amber Tools and Molecular Operating Environment (MOE) packages. Structural water plays critical role in protein interactions and functions; case studies include HET-s prion, Aβ oligomer and fibril formation, folding promotion in GroEL chaperonin, and ligand binding to maltose-binding protein. Multi-time-step MD of biomolecules steered with mean solvation forces obtained from 3D-RISM-KH at outer steps and treated with generalized solvation force extrapolation (GSFE) at inner steps is efficiently stabilized with the optimized isokinetic Nose-Hoover chain (OIN) thermostat. The accuracy of GSFE and efficiency of OIN allow picoseconds outer steps while accurately reproducing conformational properties, as validated on hydrated alanine dipeptide, miniprotein 1L2Y, and protein G. Due to 3D-RISM-KH statistical averaging over rare events of solvent exchange and localization in biomolecular spaces, this quasidynamics results in time scale compression of protein conformational changes coupled with solvent and so in huge acceleration of conformational sampling. This provides up to 1000-fold effective speedup of sampling, compared to conventional MD with explicit solvent, and enables to fold the miniprotein from a fully denatured state in 60ns quasidynamics, cf. 4-9μs folding time in experiment.

An efficient and stable hybrid extended Lagrangian/self-consistent field scheme for solving classical mutual induction.

We have adapted a hybrid extended Lagrangian self-consistent field (EL/SCF) approach, developed for time reversible Born Oppenheimer molecular dynamics for quantum electronic degrees of freedom, to the problem of classical polarization. In this context, the initial guess for the mutual induction calculation is treated by auxiliary induced dipole variables evolved via a time-reversible velocity Verlet scheme. However, we find numerical instability, which is manifested as an accumulation in the auxiliary velocity variables, that in turn results in an unacceptable increase in the number of SCF cycles to meet even loose convergence tolerances for the real induced dipoles over the course of a 1 ns trajectory of the AMOEBA14 water model. By diagnosing the numerical instability as a problem of resonances that corrupt the dynamics, we introduce a simple thermostating scheme, illustrated using Berendsen weak coupling and Nose-Hoover chain thermostats, applied to the auxiliary dipole velocities. We find that the inertial EL/SCF (iEL/SCF) method provides superior energy conservation with less stringent convergence thresholds and a correspondingly small number of SCF cycles, to reproduce all properties of the polarization model in the NVT and NVE ensembles accurately. Our iEL/SCF approach is a clear improvement over standard SCF approaches to classical mutual induction calculations and would be worth investigating for application to ab initio molecular dynamics as well.

Nonergodicity of the Nose-Hoover chain thermostat in computationally achievable time.

The widely used Nose-Hoover chain (NHC) thermostat in molecular dynamics simulations is generally believed to impart the canonical distribution as well as quasi- (i.e., space-filling) ergodicity on the thermostatted physical system (PS). Working with the standard single harmonic oscillator, we prove analytically that the two-chain Nose-Hoover thermostat with unequal thermostat masses approaches the standard Nose-Hoover dynamics, and hence the PS loses its canonical and quasiergodic nature. We also show through numerical simulations over substantially long times that for certain Poincaré sections, for both the equal and unequal thermostat mass cases, the bivariate distribution function of position and momentum (x,p) and of reservoir degrees of freedom (ξ,η) lose their Gaussian nature. Further, the four-dimensional x-p-ξ-η extended phase space exhibits two holes of nonzero measure. The NHC thermostat therefore does not generate the canonical distribution or preserve quasiergodicity for the PS.

Efficient stochastic thermostatting of path integral molecular dynamics

The path integral molecular dynamics (PIMD) method provides a convenient way to compute the quantum mechanical structural and thermodynamic properties of condensed phase systems at the expense of introducing an additional set of high frequency normal modes on top of the physical vibrations of the system. Efficiently sampling such a wide range of frequencies provides a considerable thermostatting challenge. Here we introduce a simple stochastic path integral Langevin equation (PILE) thermostat which exploits an analytic knowledge of the free path integral normal mode frequencies. We also apply a recently developed colored noise thermostat based on a generalized Langevin equation (GLE), which automatically achieves a similar, frequency-optimized sampling. The sampling efficiencies of these thermostats are compared with that of the more conventional Nosé-Hoover chain (NHC) thermostat for a number of physically relevant properties of the liquid water and hydrogen-in-palladium systems. In nearly every case, the new PILE thermostat is found to perform just as well as the NHC thermostat while allowing for a computationally more efficient implementation. The GLE thermostat also proves to be very robust delivering a near-optimum sampling efficiency in all of the cases considered. We suspect that these simple stochastic thermostats will therefore find useful application in many future PIMD simulations.

The Jarzynski identity derived from general Hamiltonian or non-Hamiltonian dynamics reproducing NVT or NPT ensembles

The Jarzynski identity (JI) relates nonequilibrium work averages to thermodynamic free energy differences. It was shown in a recent contribution [M. A. Cuendet, Phys. Rev. Lett. 96, 120602 (2006)] that the JI can, in particular, be derived directly from the Nose-Hoover thermostated dynamics. This statistical mechanical derivation is particularly relevant in the framework of molecular dynamics simulation, because it is based solely on the equations of motion considered and is free of any additional assumptions on system size or bath coupling. Here, this result is generalized to a variety of dynamics, along two directions. On the one hand, specific improved thermostating schemes used in practical applications are treated. These include Nose-Hoover chains, higher moment thermostats, as well as an isothermal-isobaric scheme yielding the JI in the NPT ensemble. On the other hand, the theoretical generality of the new derivation is explored. Generic dynamics with arbitrary coupling terms and an arbitrary number of thermostating variables, both non-Hamiltonian and Hamiltonian, are shown to imply the JI. In particular, a nonautonomous formulation of the generalized Nose-Poincare thermostat is proposed. Finally, general conditions required for the JI derivation are briefly discussed.

Heat conduction in one-dimensional oscillator lattices using Nosé–Hoover chain thermostats

In this work, we numerically study the dynamical evolution and heat transport properties of a system that consists of two time-reversible thermostats connected either by a one-dimensional Fermi–Pasta–Ulam or a Frenkel–Kontorova oscillator lattice, which are representative models of momentum-conserving and nonconserving systems, respectively. The thermostats were described by a chain of variables, Nose–Hoover chains, which enhances the ergodicity of the thermostats in comparison to the Nose–Hoover method. The time evolution of both lattices is not significantly altered by the dynamics of the thermostats. The temperature profile and heat flux of the Fermi–Pasta–Ulam model are more sensitive to the dynamics of the extended variables than those corresponding to the Frenkel–Kontorova model. Nevertheless we reproduce the scaling properties of the thermal conductivity with system size obtained by other authors.

A proof of Jarzynski’s nonequilibrium work theorem for dynamical systems that conserve the canonical distribution

We present a derivation of the Jarzynski [Phys. Rev. Lett. 78, 2690 (1997)] identity and the Crooks [J. Stat. Phys. 90, 1481 (1998)] fluctuation theorem for systems governed by deterministic dynamics that conserves the canonical distribution such as Hamiltonian dynamics, Nose-Hoover dynamics, Nose-Hoover chains, and Gaussian isokinetic dynamics. The proof is based on a relation between the heat absorbed by the system during the nonequilibrium process and the Jacobian of the phase flow generated by the dynamics.

An algorithm for the uniform sampling of iso‐energy surfaces and for the calculation of microcanonical averages

In this article an algorithm is proposed to efficiently perform the uniform sampling of an iso-energy surface corresponding to a fixed potential energy U of a molecular system, and for calculating averages of certain quantities over microstates having this energy (microcanonical averages). The developed sampling technique is based upon the combination of a recently proposed method for performing constant potential energy molecular dynamics simulations [Rapallo, A. J Chem Phys 2004, 121, 4033] with well-established thermostatting techniques used in the framework of standard molecular dynamics simulations, such as the Andersen thermostat, and the Nose-Hoover chain thermostat. The proposed strategy leads to very accurate and drift-free potential energy conservation during the whole sampling process, and, very important, specially when dealing with high-dimensional or complicated potential functions, it does not require the calculation of the potential energy function hessian. The technique proved to be very reliable for sampling both low- and high-dimensional surfaces.

Canonical Sampling Method for Initial Conditions for Reactive Flux Calculations Using Nose-Hoover Chains

Canonical sampling method has been presented to generate the initial conditions for reactive flux studies of organic reactions in water. Velocity Verlet version of Nose-Hoover chain dynamics algorithm has been employed to sample the initial conditions according to canonical distribution. The unstable normal mode of a transition state has been introduced to define a dividing plane separating reactant and product regions in reaction processes. This method has been implemented and tested for the case iels-Alder reaction of methyl vinyl ketone (MVK) and cyclopentadiene (CPD) in water, providing a reliable tool for further reactive flux molecular dynamics studies in condensed media.

Nose-Hoover chain method for nonequilibrium molecular dynamics simulation

It is demonstrated that the Nose-Hoover chain (NHC) thermostat method influences the properties of a system out of equilibrium. A simple modification of the NHC equations is proposed which enables both steady-state averages and time correlation functions to be obtained correctly for nonequilibrium states, and consistent with those computed using Gaussian isokinetic and Nose-Hoover dynamics.

Generalized Gaussian moment thermostatting: A new continuous dynamical approach to the canonical ensemble

A new method for generating the canonical ensemble via continuous dynamics is presented. The new method is based on controlling the fluctuations of an arbitrary number of moments of the multidimensional Gaussian momentum distribution function. The equations of motion are non-Hamiltonian, and hence have a nonvanishing phase space compressibility. By applying the statistical mechanical theory of non-Hamiltonian systems recently introduced by the authors [M. E. Tuckerman, C. J. Mundy, and G. J. Martyna, Europhys. Lett. 45, 149 (1999)], the equations are shown to produce the correct canonical phase space distribution function. Reversible integrators for the new equations of motion are derived based on a Trotter-type factorization of the classical Liouville propagator. The new method is applied to a variety of simple one-dimensional example problems and is shown to generate ergodic trajectories and correct canonical distribution functions of both position and momentum. The new method is further shown to lead to rapid convergence in molecular dynamics based calculations of path integrals. The performance of the new method in these examples is compared to that of another canonical dynamics method, the Nosé–Hoover chain method [G. J. Martyna, M. L. Klein, and M. E. Tuckerman, J. Chem. Phys. 97, 2635 (1992)]. The comparison demonstrates the improvements afforded by the new method as a molecular dynamics tool. Finally, when employed in molecular dynamics simulations of biological macromolecules, the new method is shown to provide better energy equipartitioning and temperature control and to lead to improved spatial sampling over the Nosé–Hoover chain method in a realistic application.

Comment on “Simple reversible molecular dynamics algorithms for Nosé–Hoover chain dynamics” [J. Chem. Phys. 107, 9514 (1997)

The leap-frog based integration algorithms for Nosé–Hoover chain dynamics recently presented by Jang and Voth [J. Chem. Phys. 107, 9514 (1997)] are shown to give rise to nonsecular growth in the conserved quantity at long times. Therefore, their use in molecular dynamics simulations should be approached with caution.

Modified nonequilibrium molecular dynamics for fluid flows with energy conservation

The nonequilibrium molecular dynamics generated by the SLLOD algorithm [so called due to its association with the DOLLS tensor algorithm (D. J. Evans and G. P. Morriss, Statistical Mechanics of Nonequilibrium Liquids (Academic, New York, 1990)] for fluid flow is considered. It is shown that, in the absence of time-dependent boundary conditions (e.g., shearing boundary conditions via explicit cell dynamics or Lees–Edwards boundary conditions), a conserved energy, H exists for the equations of motion. The phase space distribution generated by SLLOD dynamics can be explicitly derived from H. In the case of a fluid confined between two immobile boundaries undergoing planar Couette flow, the phase space distribution predicts a linear velocity profile, a fact which suggests the flow is field driven rather than boundary driven. For a general flow in the absence of time-dependent boundaries, it is shown that the SLLOD equations are no longer canonical in the laboratory momenta, and a modified form of the SLLOD dynamics is presented which is valid arbitrarily far from equilibrium for boundary conditions appropriate to the flow. From an analysis of the conserved energy for the new SLLOD equations in the absence of time-dependent boundary conditions, it is shown that the correct local thermodynamics is obtained. In addition, the idea of coupling each degree of freedom in the system to a Nosé–Hoover chain thermostat is presented as a means of efficiently generating the phase space distribution.

Centroid path integral molecular dynamics simulation of lithium para-hydrogen clusters

The real-time quantum dynamics of a series of lithium para-hydrogen clusters, Li(p-H2)n (n=13, 55, and 180), has been investigated at 2.5 and 4.0 K by means of normal mode centroid path integral molecular dynamics (NMCMD) simulation, following the methodology originally proposed by Cao and Voth [J. Chem. Phys. 101, 6168 (1994)]. The Li(p-H2)34 and neat (p-H2)34 clusters have also been simulated at 2.5 K to see the effect of doping of a Li atom on the cluster dynamics. We explicitly display both the microcanonical and the constant-temperature equations of motion for NMCMD simulations using the Nosé–Hoover chain thermostats and the reference system propagator algorithm (RESPA). In addition to the energetic and structural properties, the real-time semi-classical dynamics of the centroids of the Li atom and p-H2 molecules in the clusters has been explored to investigate the diffusive and vibrational properties. In general, quantization of the nuclear motion enhances the ease of melting and diffusion, and also causes the power spectra of the velocity autocorrelation functions of hydrogen molecules to be dramatically shifted to lower frequency than in the classical regime. The self-diffusion coefficient of p-H2 molecules in the Li(p-H2)13 cluster, which is melted completely at both 2.5 and 4.0 K, is comparable to the experimental value of bulk liquid p-H2 just above the melting temperature (14 K). No significant differences have been found for the structural properties, molecular diffusion, and the quantum vibrational spectra of p-H2 molecules between the neat (p-H2)34 and Li(p-H2)34 clusters at 2.5 K. The n=34, 55, and 180 clusters show a solidlike core and a tendency towards surface melting at 2.5 K, with diffusion restricted to the outer region. The core of the Li(p-H2)55 cluster exhibits slow liquidlike diffusion at 4.0 K, whereas that of the Li(p-H2)180 cluster remains solidlike at this temperature. The quantum power spectrum of the velocity autocorrelation function of the core molecules in the Li(p-H2)180 cluster reveals vibrational frequencies in good accord with experimental neutron and Raman data on the bulk solid. The power spectrum of the Li atom on the cluster surface is peaked at low frequency (<20 cm−1).

Generalized valence bond molecular dynamics at constant temperature

The Nose–Hoover chain thermostat of Martyna, Klein and Tuckerman has been adapted to work with our generalized valence bond ab initio molecular dynamics (AIMD) method, allowing us to perform AIMD in the canonical ensemble (constant number of particles temperature and volume). The canonical ensemble is ideal for simulated annealing which can be used to find global potential energy minima. The efficacy of constant-temperature-based AIMD simulated annealing is demonstrated for a Si6 cluster and is compared with that of a non-constant-temperature-based method. It is found that in comparison with the velocity-scaling approach to annealing, constant-temperature based annealing gives better control over the annealing process and, at least in the case under study, results in attainment of the preferred geometry of the system more efficiently.