# Calculation of ground states of four-dimensional ±𝐽 Ising spin glasses

## I Introduction

Optimization methods have found widespread application in computational physics. Among these the investigation of the low-temperature behavior of spin glasses  attracted most of the attention within the statistical physics community. The reason is that despite its simple definition (see below) its behavior is far from being understood. From the computational point of view the calculation of spin-glass ground states is very demanding, because it belongs to the class of the NP-hard problems . This means that only algorithms are available, for which the running time on a computer increases exponentially with the system size. In this work a method recently proposed, the cluster-exact approximation (CEA)  is applied to four-dimensional Ising spin glasses.

The model under investigation here consists of $`N`$ spins $`\sigma _i=\pm 1`$, described by the Hamiltonian

$$H\underset{i,j}{}J_{ij}\sigma _i\sigma _j$$
(1)   
where $`\mathrm{}`$ denotes a sum over pair of nearest neighbors. In this report simple 4d lattices are considered, i.e. $`N=L^4`$. The nearest neighbor interactions (bonds) take independently $`J_{ij}=\pm 1`$ with equal probability. Periodic boundary conditions are applied to the systems. No kind of external magnetic field is present here.

Four-dimensional Ising spin glasses have been investigated rather rarely. Most of the results were obtained via Monte-Carlo (MC) simulations at finite temperature, see e.g. . Here the $`T=0`$ behavior is investigated, i.e. ground states are calculated. This has the advantage, that one does not encounter ergodicity problems or critical slowing down like in algorithms which base on MC methods. Only one attempt  to address the 4d spin-glass ground-state problem is known to the author. But, as we will see later, the former results suffer from the problem, that not the true global minima of the energy were obtained. Furthermore, no analytic predictions of the ground-state energy have been noted by the author.

The question whether finite-dimensional Ising spin glasses show an ordered phase below a non-zero transition temperature $`T_c`$ is of crucial interest. By MC simulations around the (expected) transition temperature this question is hard to solve. Another way to address this question is to calculate the stiffness or domain wall energy $`\mathrm{\Delta }=E^aE^p`$ which is the difference between the ground-state energies $`E^a,E^p`$ for antiperiodic and periodic boundary conditions in one direction. Here the antiperiodic boundary conditions for calculating $`E^a`$ are realized by inverting one plane of bonds. For the other directions periodic boundary conditions are applied always. This treatment introduces a domain wall into the system. If a model exhibits an ordered low-temperature phase, the domain wall increases with growing system size, which becomes visible through the behavior of $`\mathrm{\Delta }`$: the disorder-averaged stiffness energy shows a finite size dependence

$$|\mathrm{\Delta }|L^{\mathrm{\Theta }_S}$$
(2)   
A positive value of the stiffness exponent $`\mathrm{\Theta }_S`$ indicates the existence of an ordered phase for non-zero temperature. For example a simple $`d=2`$ Ising ferromagnet has $`\mathrm{\Theta }_S=1`$. For spin glasses, the stiffness exponent plays additionally an important role within the droplet-scaling theory , where it describes the finite-size behavior of the basic excitations (the droplets).

Using this kind of analysis is was proven that the 2d spin glass exhibits no ordering for $`T>0`$ . For the three-dimensional problem in a recent calculation  by applying genetic CEA a value of $`\mathrm{\Theta }_S=0.19(2)`$ was found, which shows, that indeed the $`d=3`$ model has a spin-glass phase for nonzero temperature. For $`d=4`$ the existence of a finite $`T_c2.1`$ was proven rather early even by MC simulations , but the value for the stiffness-exponent $`\mathrm{\Theta }_S`$ is of interest on its own. In  recently a value of $`\mathrm{\Theta }_S=0.82(6)`$ was found by performing a MC simulation near $`T_c`$. In the work presented here the value is obtained via ground-state calculations.

The paper is organized as follows: In the next section the algorithm applied here is briefly presented. The main section contains the results for the ground-state energy and the stiffness exponent. Finally a summary is given.

## II Algorithm

The technique for the calculation bases on a special genetic algorithm  and on cluster-exact approximation  which is an optimization method designed especially for spin glasses. Now a brief description of the method is given.

Genetic algorithms are biologically motivated. An optimal solution is found by treating many instances of the problem in parallel, keeping only better instances and replacing bad ones by new ones (survival of the fittest). The genetic algorithm starts with an initial population of $`M_i`$ randomly initialized spin configurations (= individuals), which are linearly arranged in a ring. Then $`\nu \times M_i`$ times two neighbors from the population are taken (called parents) and two offspring are created using the so called triadic crossover . Then a mutation with a rate of $`p_m`$ is applied to each offspring, i.e. a fraction $`p_m`$ of the spins is reversed.

Next for both offspring the energy is reduced by applying CEA. The algorithm bases on the concept of frustration . The method constructs iteratively and randomly a non-frustrated cluster of spins, whereas spins with many unsatisfied bonds are more likely to be added to the cluster. The non-cluster spins act like local magnetic fields on the cluster spins. For the spins of the cluster an energetic minimum state can be calculated in polynomial time by using graph-theoretical methods : an equivalent network is constructed , the maximum flow is calculated  and the spins of the cluster are set to the orientations leading to a minimum in energy. This minimization step is performed $`n_{\mathrm{min}}`$ times for each offspring.

Afterwards each offspring is compared with one of its parents. The pairs are chosen in the way that the sum of the phenotypic differences between them is minimal. The phenotypic difference is defined here as the number of spin where the two configurations differ. Each parent is replaced if its energy is not lower (i.e. better) than the corresponding offspring.

After this creation of offspring is performed $`\nu \times M_i`$ times the population is halved: From each pair of neighbors the configuration which has the higher energy is eliminated. If not more than 4 individuals remain the process is stopped and the best individual is taken as result of the calculation.

The whole algorithm is performed $`n_R`$ times and all configurations which exhibit the lowest energy are stored, resulting in $`n_g`$ statistical independent ground state configurations. The method was already applied for the investigation of the ground-state landscape of 3d Ising spin glasses .

## III Results

In this section at first the values for the simulation parameters, which are defined above, are presented. Then the finite-size behavior of the ground-state energy is investigated. Finally results for the stiffness energy are discussed.

The simulation parameters were determined in the following way: For the system sizes $`L=2,4,6,7`$ several different combinations of the parameters $`M_i,\nu ,n_{min},p_m`$ were tested. For the final parameter sets it is not possible to obtain lower energies even by using parameters where the calculation consumes four times the computational effort. For $`L=3,5`$ the parameter sets for $`L+1`$ were used. Using parameter sets chosen this way genetic CEA calculates true ground states, as shown in . It should be pointed out that it is relatively easy to obtain states, which exhibit an energy slightly above the true ground state energy. The hard task is to obtain really the global minimum of the energy.

Here $`p_m=0.1`$ and $`n_R=5`$ were used for all system sizes. Table I summarizes the parameters. Also the typical computer time $`\tau `$ per ground state computation on a 80 MHz PPC601 is given.

Ground states were calculated for system sizes up to $`7\times 7\times 7\times 7`$ for $`N_L`$ independent realizations (see table I) of the random variables. For each realization the ground states with periodic and antiperiodic boundary condition in one direction were calculated. The remaining three directions are always subjected to periodic boundary conditions. One can extract from the table that for small system sizes $`L4`$ ground states are rather easily to obtain, while the $`L=7`$ systems alone required 6560 CPU-days. Using these parameters on average $`n_g>2.7`$ ground states were obtained for every system size $`L`$ using $`n_R=5`$ runs per realization.

The average ground-state energy $`e_0`$ per spin is shown in Fig.1 as a function of the system size $`L`$. Using a fit to $`e_0(L)=e_0^{\mathrm{}}+aL^b`$ the value for the infinite system is extrapolated, resulting in $`e_0^{\mathrm{}}=2.095(1)`$ ($`a=7.1(7),b=4.2(1)`$). This value is compatible with the lower bound of $`e_0=\sqrt{2d\mathrm{ln}2}2.35`$ given by the random energy model . The value calculated here is substantially smaller than the result $`e_0^{\mathrm{}}=2.054(3)`$, which was obtained in  using a pure genetic algorithm. This shows that in  not the true global minima were found, which can be concluded also from the fact, that there $`e_0(L)`$ increases with growing system size. Because the periodic boundary conditions impose additional constraints on the systems, the opposite behavior is expected, as found for the results presented here. For further comparison additionally some calculations were performed by the author by simply rapidly quenching from random chosen spin configurations. By executing an analogous fit, a value of $`e_0^{\mathrm{}}=2.04(2)`$ is obtained. This shows, that the result from  seems to be only slightly better than the data obtained by applying a very simple minimization method.

The distribution of the stiffness energy, which is obtained from performing ground-state calculations for systems with either periodic or antiperiodic boundary conditions in one direction, are shown in Fig. 2 for $`L=5`$ and $`L=7`$. With increasing system size the distribution broadens. This means that larger domain walls become more and more likely. To study this effect more quantitatively, in Fig. 3 the disorder-averaged absolute value $`|\mathrm{\Delta }|`$ of the stiffness energy is plotted as a function of the system size $`L`$. Also shown is a fit $`|\mathrm{\Delta }(L)|L^{\mathrm{\Theta }_S}`$ which results in $`\mathrm{\Theta }_S=0.64(5)`$. Here, the system sizes $`L=2,3`$ were left out of the analysis, since they are below the scaling regime. Because of the large sample sizes the error bars are small enough, so we can be pretty sure that $`\mathrm{\Theta }_S>0`$. It confirms earlier results from MC simulations  that the 4d EA spin glass exhibits a non-zero transition temperature $`T_c`$. The value $`\mathrm{\Theta }_S=0.64(5)`$ is comparable to a recent result from MC simulations $`\mathrm{\Theta }_S=0.82(6)`$ , given the facts that the system sizes are rather small and the other result was obtained at finite temperature near the transition point $`T_c2.1`$. Additionally, the prediction from droplet-scaling theory $`\mathrm{\Theta }_S<(d1)/2=1.5`$  is fulfilled.

It should be pointed out, that the method described above does not guaranty to find exact ground states, although the method for choosing the parameters makes it very likely. If states with a slightly higher energy are obtained, the result for $`e_0^{\mathrm{}}`$ is not affected very much. For the stiffness energy, it was shown in  that the result is very reliable as well, as long as the energies of the states are not too far away from the true ground-state energies.

## IV Conclusion

Results have been presented from calculations of a large number of ground states of 4d Ising spin glasses. They were obtained using a combination of cluster-exact approximation and a genetic algorithm. Using a huge computational effort it was ensured that true ground states have been obtained with a high probability.

The finite size behavior of the ground-state energy and the stiffness energy have been investigated. By performing a $`L\mathrm{}`$ extrapolation, the ground-state energy per spin for the infinite system is estimated to be $`e_0^{\mathrm{}}=2.095(1)`$. The absolute value of the stiffness energy increases with system size and shows a $`|\mathrm{\Delta }(L)|L^{\mathrm{\Theta }_S}`$ behavior with $`\mathrm{\Theta }_S=0.64(5)`$. For systems with a Gaussian distribution of the bonds qualitatively similar results are expected, since the ordering behavior depends only on the sign of the interactions and not on their magnitudes.

A more detailed study of the ground-state landscape of 4d systems, similar to , requires more than $`n_G3`$ ground states per realization to be calculated. Since this requires a substantial higher computational effort, it remains to be done for the future.

## V Acknowledgements

The author thanks A.P. Young for interesting discussions, critical reading of the manuscript and for the allocation of computer time on his workstation cluster at the University of California in Santa Cruz. This work was suggested by him during the “Monbusho Meeting” held at the Fondation Royaumont near Paris. The author was supported by the Graduiertenkolleg “Modellierung und Wissenschaftliches Rechnen in Mathematik und Naturwissenschaften” at the Interdisziplinäres Zentrum für Wissenschaftliches Rechnen in Heidelberg and the Paderborn Center for Parallel Computing by the allocation of computer time. Financial support was provided by the DFG (Deutsche Forschungsgemeinschaft) and the organizers of the “Monbusho Meeting”.