Monte Carlo simulation of classical spin models with chaotic billiards

It has recently been shown that the computing abilities of Boltzmann machines, or Ising spin-glass models, can be implemented by chaotic billiard dynamics without any use of random numbers. In this paper, we further numerically investigate the capabilities of the chaotic billiard dynamics as a deterministic alternative to random Monte Carlo methods by applying it to classical spin models in statistical physics. First, we verify that the billiard dynamics can yield samples that converge to the true distribution of the Ising model on a small lattice, and we show that it appears to have the same convergence rate as random Monte Carlo sampling. Second, we apply the billiard dynamics to finite-size scaling analysis of the critical behavior of the Ising model and show that the phase-transition point and the critical exponents are correctly obtained. Third, we extend the billiard dynamics to spins that take more than two states and show that it can be applied successfully to the Potts model. We also discuss the possibility of extensions to continuous-valued models such as the $XY$ model.

Figure 1

Absolute errors of the empirical distribution for the Ising model on a two-dimensional lattice of size $L=4$ with the periodic boundary condition at temperature $T=2.4$. The absolute error from the true distribution at time $t$ is calculated for the empirical distribution obtained from a trajectory until time $t$. The solid line indicates the average absolute error for 96 different trajectories, and the dashed lines indicate the minimum and maximum absolute errors among the trajectories. All the lines decrease at a gradient of nearly $-1/2$.

Figure 2

Absolute errors of the empirical distribution obtained by the billiard dynamics for the Ising model on a two-dimensional lattice with the periodic boundary condition. (a) The average absolute errors for different sampling algorithms. The lines indicated by “billiard1,” “billiard2,” “billiard3,” and “random” correspond to Eqs. (5), (7), (8), and random Monte Carlo sampling, respectively. (b) The average absolute errors for different lattice sizes $L=2,3,4$ at $T=2.4$. (c) The average absolute errors for different temperature values $T=2.0,2.4,2.8$ with $L=4$. All the lines show the average absolute error for 96 different trajectories. Overall, the errors of billiard sampling dynamics decay very similarly to those of random sampling.

Figure 3

Binder cumulant as a function of the temperature $T$ around $T_{\mathrm{c}}$ for $L=32$, 40, 48, 56, 64, 72, and 80. (a) The points show the values of the Binder cumulant calculated numerically by the billiard dynamics. The statistics are calculated during ${10}^7$ unit time, after initial ${10}^6$ unit time is skipped, for 96 different initial values. The solid lines indicate least-square fittings with cubic functions for each lattice size. (b) Magnification of the graphs in (a) around the critical temperature $T_{\mathrm{c}}$. The intersections of the cubic fittings provide an estimate of $T_{\mathrm{c}}=2.2690\pm 0.0002$, which is consistent with the theoretical critical temperature of $T_{\mathrm{c}}=2/\text{log}(1+\sqrt{2})=2.2691...$ of the Ising model on a two-dimensional lattice indicated by the vertical line. Note that the graphs for larger lattices tend to give better estimates.

Figure 4

Dependencies on lattice size $L$ of (a) the gradient of the cubic fittings to the Binder cumulant $\partial U/\partial T$, (b) the mean absolute magnetization $\langle \left|m\right|\rangle$, and (c) the magnetic susceptibility $\chi$, all calculated at the theoretical critical temperature $T=T_c$. The crossings show the values obtained by numerical calculations. These log-log plots give the estimates (a) $1/\nu =0.993\pm 0.016$, (b) $\beta /\nu =0.1248\pm 0.0004$, and (c) $\gamma /\nu =1.751\pm 0.004$. These results are in agreement with the theoretical values of the critical exponents $\nu =1$, $\beta =1/8$, and $\gamma =7/4$ of the Ising model on a two-dimensional lattice. The lines show least-square fits by linear functions with the theoretical gradient values.

Figure 5

Statistics of the Potts model on a two-dimensional lattice of size $L=24$ with the periodic boundary condition for $q=4$ (red lines) and $q=6$ (blue lines) calculated by the heat-bath algorithm (thick lines) and the oscillator sampling dynamics (thin lines). (a) Average order parameter $\langle m\rangle$. (b) Average energy per site $\langle \epsilon \rangle$. (c) Magnetic susceptibility $\chi$. (d) Specific heat $c$. The statistics are calculated for ${10}^5$ unit time, after initial ${10}^4$ unit time is skipped, for 100 different initial values.

Figure 6

Sampling dynamics for the $XY$ model on a two-dimensional lattice of $L=8$ at $T=1.2$ (top row) and $T=0.6$ (bottom row). The absolute errors of empirical energy distributions from the distribution obtained by the Metropolis algorithm for $t={10}^8$ are shown. The empirical energy distributions are constructed from samples until time $t$ generated by the (a) Metropolis algorithm, (b) irrational-rotation sampling dynamics, and (c) coupled-oscillator sampling dynamics. The solid lines indicate the average absolute error for 96 different sample sequences, and the dashed lines indicate the minimum and maximum absolute errors among the sequences.

Figure 7

Observation of a spin echo in the Ising model on a two-dimensional lattice with size $L=64$ at $T=2.6$. A time series of the magnetization $m$ sampled uniformly at every $0.1$ unit time is shown. At time $t=0$, all the spins are aligned. At time $t=50$, all the internal states are flipped. Subsequently, at time $t=100$, the spins become aligned again.