Empirical relations between static and dynamic exponents for Ising model cluster algorithms

We have measured the autocorrelations for the Swendsen-Wang and the Wolff cluster update algorithms for the Ising model in two, three, and four dimensions. The data for the Wolff algorithm suggest that the autocorrelations are linearly related to the specific heat, in which case the dynamic critical exponent is {ital z}{sub int,}{ital E}{sup W}={alpha}/{nu}. For the Swendsen-Wang algorithm, scaling the autocorrelations by the average maximum cluster size gives either a constant or a logarithm, which implies that {ital z}{sub int,}{ital E}{sup SW}={beta}/{nu} for the Ising model.


Introduction
The Monte Carlo cluster update algorithms of Swendsen and Wang (SW) [1] and Wolff [2] can dramatically reduce critical slowing down in computer simulations of spin models, and thus greatly increase the computational efficiency of the simulations (for reviews of cluster algorithms, see refs. [3] [4] ). There is little theoretical understanding of the dynamics of these algorithms. In particular, little is known as to why they seem to eliminate critical slowing down completely in some cases, and not others. There is no known theory which can predict the value of the dynamic critical exponent z for any spin model, although a rigorous bound on z for the SW algorithm for Potts models has been derived [5] . Another problem which is not well understood is why the SW and Wolff algorithms give similar values of z for the 2-d Potts model [6] , but have very different behavior for other models, such as the Ising model in more than two dimensions [7] [8] .
The measurement of dynamic critical exponents is notoriously difficult, and both very good statistics and very large lattices are required in order to obtain accurate results. This is certainly the case for the Ising model, where a number of different measurements have given conflicting results. For the two dimensional Ising model, initial results suggested z ≈ 1/3 for both the SW and Wolff algorithms [1] [8]. Further work [7] gave z ≈ 1/4, and it was later shown that the data were consistent with a logarithmic divergence, suggesting that z = 0 [9] . Recent results show that it is very difficult to distinguish between a logarithm and a small power [6].
Measurements on the three dimensional model have proven to be just as difficult, with values of z for the SW algorithm ranging from 0.339(4) to 0.75(1) [1][7] [10] . For the Wolff algorithm, Tamayo et al. [8] obtained 0.44(10), while Wolff found a value of 0.28(2) for the energy autocorrelations [7]. We have examined Wolff's data and found that it also fits well to a logarithm, so that z = 0 is also a possibility.
In four dimensions only one result is known, which is z = −0.05 (15) for the Wolff algorithm [8]. Simulations have also been done on the mean-field Ising model, which is expected to give the same exponents as the Ising model in four or more dimensions [11] . The mean-field data are consistent with z being 0 for the Wolff algorithm [8] and 1 for SW [12] , with the latter result being supported by theoretical arguments.

Simulations
Due to the discrepancies between the various measurements of the dynamic critical exponents, we have done numerical simulations of the Ising model in 2, 3 and 4 dimensions using the SW and Wolff algorithms, with the aim of obtaining good statistics on fairly large lattices, in order to get reliable values for the dynamic exponents. We measured the time correlation function ρ(t) for the energy, and extracted the integrated autocorrelation time [3] The dynamic critical exponent z is given by τ ∼ L z , where τ for the different lattice sizes is measured at the infinite volume critical point. We have used the Potts formulation of the Ising model, for which the critical point in two dimensions is known to be β c = log(1 + √ 2) ≈ 0.8813736 [13] . For the 3-d model we used the value 0.443308 [14] , while in the 4-d case we have used 0.29962 [15] . A detailed account of the methods we used to do the measurements, fits and error estimates, is given in ref. [6].
Autocorrelations are traditionally measured between each update of the entire lattice, so for the single cluster Wolff update, where only a fraction of the lattice sites are updated at each iteration, the measured autocorrelation time τ needs to be scaled by the ratio of the average Wolff cluster size < |c W | > and the number of lattice sites L d . The scaled autocorrelation time is what we present for the Wolff autocorrelations. Since this scaling ratio is an estimator for the susceptibility [2], the dynamic critical exponent z for the unscaled autocorrelations is given by z = z + (d − γ/ν), where ν is the critical exponent for the correlation length, and γ is the critical exponent for the susceptibility, which diverges as L γ/ν . For the SW algorithm on the larger lattice sizes in two and three dimensions, we used a parallel cluster labeling algorithm which we have developed [16] in order to run on large parallel machines. For the other lattice sizes, we ran multiple simulations in parallel using smaller shared memory machines and networks of workstations.

Results
Results for τ int,E , the integrated autocorrelation time for the energy, are shown in Figures 1(a), (b) and (c) for d = 2, 3 and 4 respectively. For d = 3 we have used a log-log plot, with the straight lines representing χ 2 fits to a power law, while for d = 2 and d = 4 we have used a log-linear plot, with the straight lines representing χ 2 fits to a logarithm. Note however that for d = 4 we plot log τ int,E rather than τ int,E for the SW algorithm, since the SW autocorrelations increase as a power of L. The measured values of the exponents from the fits to the data are shown in Table 1. For the Wolff algorithm in all dimensions, and the SW algorithm in two dimensions, it is very difficult to distinguish between a small exponent and a logarithmic increase in the autocorrelations (which would imply that z = 0).
In Fig. 1 we also include the measured value of the specific heat C H , scaled by an appropriate factor, in order to show that the bound of Li and Sokal [5] is indeed satisfied by the SW algorithm. Here α is the critical exponent for the specific heat (C H ∼ L α/ν ). No such bound has been proven for the Wolff algorithm [17] , although it appears from the figures that not only does the bound hold, but that there may actually be equality in the exponents. If we compare the results of fits to C H and τ W int,E (the autocorrelations in the energy for the Wolff algorithm), which correspond to the measured values of α/ν and z W int,E respectively, then for d = 3 we find 0.32(1) and 0.33(1). In two and four dimensions α = 0, and z W int,E is also consistent with zero. Hence the Wolff algorithm for the Ising model seems to satisfy the surprisingly simple relations where a and b are constants. In Figure 2 we plot the difference τ W int,E − (a + b × C H ) for the various dimensions, with a and b chosen to minimize χ 2 over a certain range of lattice sizes (smaller values of L are excluded from the fit). We can see that in all cases, values of a and b can be found such that the difference is zero within errors. Note that all the errors shown here are purely statistical (one standard deviation). In two dimensions the best fit is obtained with a ≈ −0.474 and b ≈ 0.957 (the data does not exclude the possibility that b = 1, which would imply that τ W int,E is just a constant plus C H ). For the 3-d model the additive constant a is consistent with zero, so that the autocorrelation time may be just a multiple of the specific heat, with b ≈ 0.148. In four dimensions we find a ≈ 0.167 and b ≈ 0.050.
The surprising simplicity of the result (3) led us to look for a similar relation for the SW algorithm. The power of cluster update algorithms comes from the fact that they flip large clusters of spins at a time. The relative average size of the largest SW cluster, m =< |c max SW | > /L d , is an estimator of the magnetization [18] , and the exponent β/ν characterizing the divergence of the magnetization has values which are similar to our measured values for the dynamic exponents of the SW algorithm. To investigate this further, we have scaled the SW autocorrelations by m, in a similar manner to the scaling of the Wolff autocorrelations in equation (1). If this gives a constant or a logarithm, then τ SW int,E diverges like the magnetization, and so we have z SW int,E = β/ν. The SW autocorrelations scaled by m (and by an additional arbitrary constant, so that these points are not entangled with others in the plots) are also shown in Fig. 1. For d = 4 the results are very close to a constant, while for d = 3 they seem to approach a constant as L increases. In two dimensions the scaled autocorrelations are not constant, but they fit very much better to a logarithm than does the unscaled data, as can be seen in Fig. 1(a), and fit very poorly to a power law. The data therefore support the assertion that Our measurements of z W int,E in 3-d and z SW int,E in 4-d give results which are very different from the accepted values of α/ν (≈ 0.10 − 0.20) [19] [20] and β/ν (= 1) [15] [21] , since corrections to scaling are known to be important for these quantities. If we do a simple power law fit to the specific heat in 3-d and the magnetization in 4-d, we get results which are also very different from the actual exponents, but which are very close to the measured values of the corresponding dynamic exponents, as expected from relations (3) and (4). Thus although we may not be able to measure the asymptotic behavior of the autocorrelation times, finding simple relations between the autocorrelations and static quantities whose asymptotic behavior is known enables us to infer the true values of the dynamic critical exponents. This is especially useful for the 2-d model, for which the autocorrelations grow so slowly that any corrections to scaling could have a big effect. It is therefore very difficult to say with any confidence that z = 0, even with data on very large lattices. The apparent relation (3) seems to be the most compelling evidence so far that z W int,E is in fact zero for the 2-d Ising model, while the relation (4) would imply that z SW int,E is actually 1/8, which is not apparent from the usual fits to either a logarithm or a power law.

Conclusions
We have measured the autocorrelations and dynamic critical exponents of the SW and Wolff cluster algorithms for the Ising model in 2, 3 and 4 dimensions. We have found what appear to be surprisingly simple empirical relations between the autocorrelation times of these algorithms and simple static quantities (the magnetization and specific heat). These relations could perhaps stem from the fact that the dynamics of cluster algorithms are closely linked to the physical properties of the system, since the Swendsen-Wang clusters are just the Coniglio-Klein-Fisher droplets [22] , or "physical clusters" [18], from which the critical behavior of the system may be described.
The relations (3) and (4) are certainly not general results, since for the 2-d q = 3 Potts model we find that z W > α/ν and z SW > β/ν [5] [6]. Also, it is quite surprising that these empirical relations imply that z SW is not equal to z W for the 2-d Ising model, whereas the two appear to be equal for the 2-d q = 3 Potts model. It is of course possible that these relations are not exact, but merely good approximations. We are currently collecting more data in order to check whether these results hold up with larger lattices and better statistics, and we will present more detailed results in a future publication [23] . Measured dynamic critical exponents for Ising model cluster algorithms. Asterisks indicate that the data is also consistent with a logarithmic divergence (z int,E = 0).