Energy gap of the bimodal two-dimensional Ising spin glass

An exact algorithm is used to compute the degeneracies of the excited states of the bimodal Ising spin glass in two dimensions. It is found that the specific heat at arbitrary low temperature is not a self-averaging quantity and has a distribution that is neither normal or lognormal. Nevertheless, it is possible to estimate the most likely value and this is found to scale as L^3 T^(-2) exp(-4J/kT), for a L*L lattice. Our analysis also explains, for the first time, why a correlation length \xi ~ exp(2J/kT) is consistent with an energy gap of 2J. Our method allows us to obtain results for up to 10^5 disorder realizations with L<= 64. Distributions of second and third excitations are also shown.

In spite of its comparative simplicity, the twodimensional bimodal short-range Ising spin glass model remains an interesting source of controversy. Although it is now widely accepted that the spin glass only exists at zero temperature [1,2], the nature of excitations, in particular the energy gap with the first excited state, has commanded much interest in the literature.
The bimodal model has bond (nearest-neighbor) interactions of fixed magnitude J and random sign. If we think of an infinite square lattice, without open boundaries, it is easy to appreciate that any finite number of spin flips cannot result in an excitation energy of less than 4J. Nevertheless, some 20 years ago, Wang and Swendsen [3] gave credible evidence that the energy gap with the first excited state should be 2J in the thermodynamic limit. These excitations must involve an infinite number of spin flips. The issue here is the noncommutativity of the zero-temperature and thermodynamic limits. In such a situation it is imperative to perform the thermodynamic limit first. Support for the 2J energy gap has included studies that have involved exact computations of partition functions [4], a worm algorithm [5] and Monte Carlo simulation [6]. Challenges to 2J have also appeared. Saul and Kardar [7] maintained that the energy gap should be 4J as naive analysis suggests. More recently [8], it has been reported that the specific heat should follow a power law in temperature. In particular, it was proposed that the critical exponents must be the same as those for the model with a Gaussian distribution of bond interactions, indicating universality with respect to the type of disorder. Fisch [9], using a steepest descent approximation, has also suggested power law behaviour albeit with a different exponent.
An important quantity, intimately involved with studies of the thermodynamic limit at fixed temperature, is the correlation length. Essentially this measures the spatial extent of the influence of one spin on others. The correlation length is infinite at a critical temperature. For the spin glass model of interest here, the correlation length has been determined by reliable Monte Carlo tech-niques [2,10] to be ξ ∼ exp(2J/kT ). This is in agreement with Ref. [7] and consistent with a qualitative study [11]. If this is true then hyperscaling predicts that the energy gap should be 4J. A 2J gap would be consistent with ξ ∼ exp(J/kT ) as proposed in Ref. [6]. Another scenario is power law behaviour [8].
A good current review of the issues involved here has been given in Ref. [12]. The power law behaviour of Ref. [8] is discussed in the light of new Monte Carlo data. The conclusion is that the suggested universality cannot be reliably proven with the computational facilities currently available. Further, it is stated that extrapolation of the data of Ref. [6] to zero temperatures may not be plausible. The main message of this letter is to argue, for the first time, why a correlation length ξ ∼ exp(2J/kT ) is perfectly consistent with an energy gap 2J, in apparent violation of hyperscaling.
We have performed exact calculations of the degeneracies of excited states at a fixed arbitrarily low temperature. Each disorder realization consisted of a frustrated L × L patch with periodic boundary conditions in one dimension, embedded in an infinite unfrustrated environment in the second dimension. This choice of boundary conditions, with even L, definitely does not allow any first excitation with energy gap less than 4J. There are no open boundaries and no diluted bonds. A 2J energy gap can only arise in the thermodynamic limit. Any planar Ising model is isomorphic to a system of interacting fermions. The Pfaffian formalism [13] is particularly convenient for disordered systems [14,15]. Each bond is decorated with two fermions, one either side, so that a square plaquette has four; left, right, top and bottom. The partition function is given by Z ∼ (det D) 1/2 where D is a real skew-symmetric 4N × 4N matrix for a lattice with N sites. Perturbation theory is used to determine ground state properties. Basically we require the low-temperature behaviour of the defect (meaning zero at zero temperature) eigenvalues of D. Exactly, The temperature dependence appears in δ only. D 0 is singular in the ground state and has eigenvectors | d corresponding to zero eigenvalue, that is D 0 | d = 0, localized on each frustrated plaquette. These eigenvectors form the subbasis for the perturbation theory.
At first order, the matrix D 1 , which is 2 × 2 block diagonal across bonds, is diagonalized in the subbasis. To continue at second order requires the continuum Green's function [15,16] The first term g c1 is 4 × 4 block diagonal in the four fermions in each plaquette and allows us to connect frustrated plaquettes across two bonds. At second order we diagonalize D 2 = D 1 g c1 D 1 . The matrix g c2 is, just like D 1 , 2 × 2 block diagonal across bonds and only relevant for excited states. For higher orders we require Green's functions G r constructed at previous orders. At third order the matrix to be diagonalized is until all degeneracy is lifted. Here we define, for r ≥ 1, and there are N (r) such pairs at order r. Note that G r is real here, as are all matrices; in contrast to Ref. [15] where they are imaginary. The internal energy of the system is where σ i is an Ising spin, ζ ij represents the sign of the bond and the nearest-neighbor correlation function can be written [14] ζ where G −+ means the matrix element of G = −D −1 between the two nodes decorating the bond ij [17]. The full Green's function G can be expanded for a finite system in powers of δ where r max is the highest order of perturbation theory required. It is obvious that K m , for m > 1, has no direct physical meaning. We provide here a brief outline of how the perturbation theory for the Green's function allows us to expand the internal energy exactly. Equating powers of δ in GD = −1, we obtain (for m < r max ) Since P D 0 = 0, we can show that for m ≤ 0, It also can be proven that (1+K 1 D 1 ) is idempotent, that is Since both g c2 and D 1 are bond diagonal and g c2 D 1 = 1 2 [16], we get K 0 = (1+K 1 D 1 )g c1 (1+D 1 K 1 )+g c2 + 1 2 K 1 , and for m < 0 Substituting t = 1 − δ and the Green's function in Eq. (3), the correlation function can be expressed in terms of K 1 as We can now use the binomial theorem to expand in exp(−2J/kT ) instead of δ, which eliminates the explicit effect of g c2 , and the internal energy can be expressed as where Since for any skew-symmetric matrix A. We can then show under the trace that for m > 0, (12) This can be rewritten using the idempotence relation (Eq. (7)) to get Finally, with recursive expansion, this can be expressed as where This is sensible since it includes all the features of the ground-state calculation. It is also exact. We can imagine coloring plaquettes black and white; like a chess board. Matrices D r and G r are color diagonal for even r: otherwise off-diagonal. Since g c1 is color diagonal, R is off-diagonal and it follows that U m = 0 for odd m. This explicitly excludes any 2J excitations involving a finite number of spins.
The specific heat per spin can be derived in terms of the internal energy as We denote the degeneracy of the i th excited state as M i . Expanding ln Z, we obtain, for example, As a means of establishing bearings, we first report a study of the first excitations of the fully frustrated Villain model [18,19]. Fig. 1 shows clearly that the ratio M1 M0 ∼ 2 π L 2 ln L. The correlation length is known to be ξ ∼ δ −1 [20]. Replacing ln L with − ln δ ∼ T −1 , the correct form for the specific heat is obtained, that is c v ∼ T −3 exp(−4J/kT ). The degeneracy per spin of the first excited state is infinite in the thermodynamic limit, although only weakly (logarithmically) so. In Fig. 2 we show the distribution of M1 M0 for the spin glass. It is clear that the most likely value scales as L 3 . In this case we have an extra factor L, unlike ln L for the Villain model. In consequence, taking ξ ∼ exp(2J/kT ), the specific heat varies like c v ∼ T −2 exp(−2J/kT ) explaining why a 2J energy gap arises from the form of the correlation length obtained from Monte Carlo calculations [2,10]. Hyperscaling fails here. We also emphasize that the Villain model is not a spin glass and comparisons, such as in Ref. [7], are not meaningful.
The distributions are scaled with a factor of L, probably indicating that the amount of effort required for an experiment to find the mode of the distribution scales like L. This is consistent with the known difficulties that arise when trying to extrapolate Monte Carlo data to zero temperature. Also, the distributions are obviously neither normal or lognormal. Further, the data is not self averaging. In fact, the relative variance grows quickly with L which is unusually severe; convergence to a constant is normally expected [21]. We can discuss our choice of correlation length in the light of Ref. [12]. Two crossover temperatures are defined. Below T * (L) the ground state behaviour dominates. Above the finite-size crossover temperature T ξ (L) no size limitations are expected. It is clearly stated that, in all situations, T ξ > T * . This immediately rules out ξ ∼ exp(J/kT ). For our case, we can in fact make the definite prediction T ξ = 1.5T * . A power law behaviour for ξ must also be ruled out if we also expect similar behaviour for the specific heat.
Pfaffians are also computed, for the partition function, in Ref. [4] at definite finite temperatures. For comparison, we have computed the mean of ln( M1 M0 ) and find fairly close agreement. Nevertheless, we do not believe that this is physically meaningful. The distributions of M1 M0 are not lognormal and the most likely value scales as L 3 , not L 4 . We emphasize here that we are not computing the entire Pfaffians.
We have also checked distributions for second excitations. The most likely value of M2 M0 − N scales as L 6 . We have subtracted off here the (infinite) total number of spins N . Single spin flips can occur anywhere and it is sensible to measure the internal energy relative to the unfrustrated system. However, the appropriate contribu- tion to the internal energy is M2 M0 − 1 2 ( M1 M0 ) 2 −N . As shown in Fig. 3, the most likely value is close to zero. We note that the mode of 1 2 ( M1  Fig. 4. Here, the mode of M3 M0 − 2N scales as L 9 and the unfrustrated system has M 3 = 2N . It is unlikely that higher excitations will interfer with our arguments based on first excitations. In conclusion, we have reported exact results for the excitations of the bimodal two-dimensional Ising spin glass by expanding in arbitrary temperature from the ground state. All other treatments, except Ref. [7], have required extrapolation from definite finite temperatures. We have argued that an energy gap of 2J is consistent with ξ ∼ exp(2J/kT ) as found from Monte Carlo simulations. The manner in which our model is arranged excludes the possibility of any 2J excitation for a finite system. Nevertheless, it should also be interesting to study systems with obvious 2J excitations, for example the hexagonal lattice or square lattice with open boundaries or diluted bonds, and investigate the distributions. W. A. thanks the Commission on Higher Education Staff Development Project, Thailand for a scholarship. J. P. acknowledges fruitful conversations with J. A. Blackman. Some of the computations were performed on the Tera Cluster at the Thai National Grid Center.