Singly-Thermostated Ergodicity in Gibbs' Canonical Ensemble and the 2016 Ian Snook Prize Award

TL;DR

This paper confirms the effectiveness of a new hyperbolic tangent thermostat in achieving ergodicity for various one-dimensional systems, advancing the development of single-thermostat algorithms for Gibbs' canonical ensemble.

Contribution

It introduces a stiff hyperbolic tangent thermostat force and demonstrates its ability to produce Gibbs' distribution in several systems, addressing the challenge of ergodic algorithms with a single thermostat.

Findings

Successfully reproduces Gibbs' distribution for multiple systems

Confirms the ergodic properties of the new thermostat

Highlights an interesting feature with adaptive integrators in the Mexican Hat problem

Abstract

The 2016 Snook Prize has been awarded to Diego Tapias, Alessandro Bravetti, and David Sanders for their paper -- Ergodicity of One-Dimensional Systems Coupled to the Logistic Thermostat. They introduced a relatively stiff hyperbolic tangent thermostat force and successfully tested its ability to reproduce Gibbs' canonical distribution for the harmonic oscillator, the quartic oscillator, and the Mexican Hat potentials. Their work constitutes an effective response to the 2016 Ian Snook Prize Award goal -- Finding ergodic algorithms for Gibbs' canonical ensemble using a single thermostat variable. We confirm their work here and highlight an interesting feature of the Mexican Hat problem when it is solved with an adaptive integrator.

Equations

$$\{(q^2/2);(q^4/4);(q^4/4)-(q^2/2)\}.$$

$$\{\dot{q}=p;\dot{p}=-q-\zeta p;\dot{\zeta }=[(p^2/T)-1]/{\tau }^2\}⟶$$

$$f(q,p,\zeta )\propto e^{-(q^2/2T)}{e}^{-(p^2/2T)}{e}^{-({\zeta }^2{\tau }^2/2)}[\mathrm{Nos}\acute{\mathrm{e}}-\mathrm{Hoover}].$$

$$(\partial f/\partial t)=-\nabla \cdot (fv)=0,$$

$$\{\dot{q}=p;\dot{p}=-q-\zeta p-\xi p^3;\dot{\zeta }=p^2-1;\dot{\xi }=p^4-3p^2\}⟶$$

$$f(q,p,\zeta )\propto e^{-(q^2/2)}{e}^{-(p^2/2)}{e}^{-({\zeta }^2/2)}{e}^{-({\xi }^2/2)}[\mathrm{Hoover}-\mathrm{Holian}].$$

$$\{\dot{q}=p;\dot{p}=-q-\zeta (0.05p+0.32p^3);\dot{\zeta }=0.05(p^2-1)+0.32(p^4-3p^2)\}.$$

$$q_{n+1}=c{q}_n(1-q_n)\mathrm{and}\dot{q}=q(1-q).$$

$$\dot{q}=\frac{1}{[e^{+t/2}+e^{-t/2}{]}^2}⟷q=\frac{e^{+t/2}}{[e^{+t/2}+e^{-t/2}]}⟷2q=1+\tanh (+t/2).$$

$$\{\dot{q}=p;\dot{p}=q-q^3-50p\tanh (25\zeta );\dot{\zeta }=p^2-1\}.$$

$$\{\dot{q}=p;\dot{p}=q-q^3-\alpha p\tanh (\alpha \zeta );\dot{\zeta }=p^2-1\},$$

$$10^{-12}>\sqrt{\delta q^2+\delta p^2+\delta {\zeta }^2}>10^{-14}\mathrm{to}10^{-15}>\sqrt{\delta q^2+\delta p^2+\delta {\zeta }^2}>10^{-17}.$$

Full text

[TABLE]

Their work constitutes an effective response to the 2016 Ian Snook Prize Award goal, “finding ergodic algorithms for Gibbs’ canonical ensemble using a single thermostat”. We confirm their work here and highlight an interesting feature of the Mexican Hat problem when it is solved with an adaptive integrator.

Ergodicity, Chaos, Algorithms, Dynamical Systems

I Nosé and Nosé-Hoover Canonical Dynamics Lack Ergodicity

In 1984 Shuichi Nosé used “time scaling”b1 ; b2 to relate his novel Hamiltonian ${\cal H}$ to an extended version of Gibbs’ canonical phase-space distribution $f$ , proportional to $e^{-{\cal H}/kT}$ . Hoover’s simpler “Nosé-Hoover” motion equationsb3 dispensed with Hamiltonian mechanics and time scaling, reducing the dimensionality of the extended phase space by one. For the special case of a harmonic oscillator the Nosé-Hoover motion equations and the corresponding modified Gibbs’ distribution are :

[TABLE]

[TABLE]

Here $q$ and $p$ are the oscillator coordinate and momentum. $\zeta$ is a “friction coefficient”, or “control variable”. In all that follows we choose the equilibrium temperature $T$ equal to unity to simplify notation. The timescale of the thermal response to the imposed equilibrium temperature is governed by the relaxation time $\tau$ . For simplicity we choose force constants, masses, and Boltzmann’s constant all equal to unity.

Hoover used the steady-state phase-space continuity equation :

[TABLE]

to show that Gibbs’ canonical distribution is consistent with the Nosé-Hoover motion equations. Here the phase-space flow velocity is $v\equiv(\dot{q},\dot{p},\dot{\zeta})$ . Hoover’s numerical work showed that only a portion of the three-dimensional Gaussian distribution (typically just a two-dimensional torus) is generated. That is, solutions of these three-dimensional motion equations are not ergodic. Particular solutions fail to cover the entire $(q,p,\zeta)$ phase space.

Considerable numerical work, following the comprehensive analyses of Kusnezov, Bulgac, and Bauerb4 ; b5 , suggested that using two thermostat variables rather than one was the simplest route to oscillator ergodicity. Including another thermostat variable requires a four-dimensional $(q,p,\zeta,\xi)$ phase space. A successful exampleb6 , ergodic in $(q,p,\zeta,\xi)$ space, controlled two velocity moments, $\langle\ p^{2}\ \rangle$ and $\langle\ p^{4}\ \rangle$ rather than just one :

[TABLE]

[TABLE]

In 2015 a single-thermostat approachb7 with simultaneous weak control of $\langle\ p^{2}\ \rangle$ and $\langle\ p^{4}\ \rangle$ was found to generate Gibbs’ entire distribution for the harmonic oscillator :

[TABLE]

Straightforward generalizations of this single-thermostat approach failed to thermostat the quartic and Mexican Hat potentials, leading to the posing of the 2016 Snook Prize problem solved by Tapias, Bravetti, and Sandersb8 .

II Tapias, Bravetti, and Sanders’ “Logistic” Thermostat

The Logistic Map and the Logistic Flow are two simple models for chaotic behavior :

[TABLE]

A solution of the logistic flow equation is

[TABLE]

With these logistic equations in mind Tapias, Bravetti, and Sandersb8 suggested a hyperbolic tangent form for the thermostat variable, and showed, with a variety of numerical techniques, convincing evidence for the ergodicity of their “Logistic Thermostat” motion equations for the quartic and Mexican Hat potentials as well as the simpler harmonic oscillator problem.

In the most challenging case, the Mexican Hat potential, the ergodic set of motion equations found by Tapias, Bravetti, and Sanders was feasible to solve, but relatively stiff :

[TABLE]

In replicating their work we also characterized solutions of a slight variant :

[TABLE]

where values of the parameter $\alpha$ in the neighborhood of seven lead to apparent ergodic behavior in $(q,p,\zeta)$ space.

One of the simplest and most useful tests for ergodicity in three dimensions is the lack of holes in the two-dimensional cross-sections (as opposed to projections) of the three-dimensional flow. For stiff equations it is convenient to use “adaptive” integrations of the motion equations where the timestep varies to maintain the accuracy of the integratorb9 .

In our own numerical work we integrated for a time of 10,000,000 using timesteps which maintained the rms difference between a fourth-order or fifth-order Runge-Kutta step of $dt$ and two such steps with $(dt/2)$ to lie within a band varying from

[TABLE]

We generated about 3,000,000 $\{\ q,0,\zeta\ \}$ double-precision cross-section points in laptop runs taking about an hour each. Typical timesteps were in the range from 0.0001 to 0.001 .

Figure 1 shows portions of the $p=0$ cross section with $\alpha=6.5$ which has evident holes at $(q=\pm 1.92,\zeta=\pm 0.39)$ . The holes disappear if $\alpha$ is increased to 6.9. But a look at the $(q,0,\zeta)$ section with an error band of $10^{-13\pm 1}$ reveals not only “normal” (irregularly-dotted) regions but also a few striped regions. In Figure 2 we see that the stripes using RK4 differ from those using RK5 showing that the stripes are artefacts. Tightening the error band to $10^{-15\pm 1}$ confirms this diagnosis, as shown in Figure 3. The interesting structure of these striped regions is a thoroughly unexpected fringe benefit of the new logistic thermostat.

We thank Drs Tapias, Bravetti, and Sanders for their stimulating prize-winning work.