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 $\left(x,p\right)$ and of reservoir degrees of freedom $\left(\xi ,\eta \right)$ lose their Gaussian nature. Further, the four-dimensional $x\text{−}p\text{−}\xi \text{−}\eta$ 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.

Figure 1

JPDF plot of $x\text{−}p$ corresponding to different Poincaré sections for four different pairs of thermostat masses with initial conditions as $x=1.1,p=1.1,\eta =0,\xi =0$. The different figures correspond to (a) $Q_{\eta }=1,Q_{\xi }=1$ at the Poincaré section $\eta =-1,\xi =1$, (b) $Q_{\eta }=10,Q_{\xi }=0.1$ at the Poincaré section $\eta =-3.16,\xi =0.03$, (c) $Q_{\eta }=50,Q_{\xi }=0.02$ at the Poincaré section $\eta =-7.07,\xi =0.003$, (d) $Q_{\eta }=100,Q_{\xi }=0.01$ at the Poincaré section $\eta =-10,\xi =-0.003$. For each case we can observe that there is a deviation from normal distribution.

Figure 2

The $\eta \text{−}\xi$ joint probability density plots for a Poincaré section at $x=0,p=0$. The different figures correspond to (a) $Q_{\eta }=1,Q_{\xi }=1$, (b) $Q_{\eta }=10,Q_{\xi }=0.1$, (c) $Q_{\eta }=50,Q_{\xi }=0.02$, (d) $Q_{\eta }=100,Q_{\xi }=0.01$. The initial conditions remain the same as before. The figures indicate the presence of two holes in the system dynamics.

Figure 3

Convergence of Hellinger (a) and modified KL distances (b) for various cases. Bins of size $\Delta x=0.1$ and $\Delta p=0.1$ are used. None of the distributions converge to the standard uncorrelated bivariate normal.

Figure 4

Convergence of the (a) joint fourth moment and (b) joint sixth moment of position and velocity for two values of thermostat parameters. It is evident that the moments have converged to values that deviate from the desired standard normal values (shown as dashed pink line). A similar behavior was observed for other moments as well.

Figure 5

Inability of projected variables $n_1\text{−}{n}_2$ to capture the four-dimensional hole of radius 0.25 forcefully embedded within a four-dimensional joint standard normal. (a) Projected values of $n_1\text{−}{n}_2$ onto the $n_3=0,n_4=0$ plane. (b) Marginal distribution of $n_1$ and (c) marginal distribution of $n_2$. No difference in marginal distributions from that of standard normal can be observed.

Figure 6

Phase space plot of $x\text{−}p$ corresponding to different Poincaré sections for four different pairs of thermostat masses. The different figures correspond to (a) $Q_{\eta }=1,Q_{\xi }=1$ at the Poincaré section $\eta =-1,\xi =1$, (b) $Q_{\eta }=10,Q_{\xi }=0.1$ at the Poincaré section $\eta =-3.16,\xi =0.03$, (c) $Q_{\eta }=50,Q_{\xi }=0.02$ at the Poincaré section $\eta =-7.07,\xi =0.003$, (d) $Q_{\eta }=100,Q_{\xi }=0.01$ at the Poincaré section $\eta =-10,\xi =-0.003$. Existence of unoccupied regions near origin for all the cases suggests that holes are present in phase space.

Figure 7

Phase space plot of $\eta \text{−}\xi$ corresponding to different Poincaré sections for four different pairs of thermostat masses. The different figures correspond to (a) $Q_{\eta }=1,Q_{\xi }=1$, (b) $Q_{\eta }=10,Q_{\xi }=0.1$, (c) $Q_{\eta }=50,Q_{\xi }=0.02$, (d) $Q_{\eta }=100,Q_{\xi }=0.01$. Two holes can be easily seen.

Figure 8

Nonzero measure holes in each of the four cases using 20 billion time steps each of 0.01. In each of the figures, $-0.1\le x\le 0.1,-0.1\le p\le 0.1$. The different figures correspond to (a) $Q_{\eta }=1,Q_{\xi }=1$, (b) $Q_{\eta }=10,Q_{\xi }=0.1$, (c) $Q_{\eta }=50,Q_{\xi }=0.02$, (d) $Q_{\eta }=100,Q_{\xi }=0.01$.