Ordering Behavior of the Two-Dimensional Ising Spin Glass with Long-Range Correlated Disorder
OOrdering Behavior of the Two-Dimensional Ising Spin Glass with Long-RangeCorrelated Disorder
L. M¨unster, C. Norrenbrock, A. K. Hartmann, ∗ and A. P. Young Institut f¨ur Physik, Universit¨at Oldenburg, 26111 Oldenburg, Germany University of California Santa Cruz, CA 95064, USA (Dated: February 3, 2021)The standard two-dimensional Ising spin glass does not exhibit an ordered phase at finite temper-ature. Here, we investigate whether long-range correlated bonds change this behavior. The bondsare drawn from a Gaussian distribution with a two-point correlation for bonds at distance r thatdecays as (1 + r ) − a/ , a ≥
0. We study numerically with exact algorithms the ground state anddomain wall excitations. Our results indicate that the inclusion of bond correlations does not leadto a spin-glass order at any finite temperature. A further analysis reveals that bond correlationshave a strong effect at local length scales, inducing ferro/antiferromagnetic domains into the system.The length scale of ferro/antiferromagnetic order diverges exponentially as the correlation exponentapproaches a critical value, a −→ a crit = 0. Thus, our results suggest that the system becomes aferro/antiferromagnet only in the limit a → PACS numbers: 75.40.Mg, 02.60.Pn, 68.35.Rh
I. INTRODUCTION
Spin glasses are disordered magnetic materials whichexhibit peculiar properties at very low temperatures [1].To understand these materials, the Edwards-Anderson(EA) model and the Sherrington-Kirkpatrick (SK) model[2, 3] have been developed. Spin glasses exhibit essen-tial aspects of complex behavior [4] and research on spinglasses [5–7] has stimulated progress in numerous otherfields, such as information processing [8], neuronal net-works [9], discrete optimization [10, 11] and Monte-Carlosimulation [12].In this work we study the two-dimensional EA spinglass model with Ising spins. This model has short-rangequenched random pair-wise interactions, described in de-tail in Sec. II. Its properties are well described in theframework of the scaling/droplet picture [13–15], as hasbeen confirmed by numerical calculations for large sys-tems using exact ground-state algorithms [16, 17]. Themodel exhibits no finite-temperature spin glass phas,e incontrast to the three or higher dimensional variants [18–20].At the zero-temperature phase transition, the distri-bution of the interaction disorder, in particular differ-ences between continuous Gaussian and discrete bimodal ± J disorder distributions have been the subject of in-tensive research [21–24]. Since the non-existence of afinite-temperature spin-glass phase for short range two-dimensional models is independent of the disorder distri-bution, we consider here only the Gaussian case. Previ-ous works have shown how the increase of the mean of theGaussian from zero to a sufficiently large value induces aferromagnetic phase [25–29].In this work we address the question how long-range ∗ [email protected] correlations in the interactions (bonds) affects the or-dering behavior, in particular whether it leads to a low-temperature spin-glass phase. Here, long-range meansthat the bond correlation decays with a power law andso does not have a characteristic length scale. For disor-dered ferromagnets the effects of long-range correlationshave already been taken studied [30].Our study was partially motivated by a correspondingnumerical study of the three-dimensional random-fieldIsing model with long-range correlation [31], where, forstrong correlation, an influence on the quantitative or-dering behavior has been observed for some critical ex-ponents. Note, however, that the case of the random-fieldmodel is a bit different, because the correlation acts onthe random fields which provide local randomness com-peting against long-range ferromagnetic order. Conse-quently, from an extended Imry-Ma argument it was pre-dicted [32] that the random-field Ising model with stronglong-range correlation will show an increase of the lowercritical dimension for ferromagnetic order. We note thatthere is no analogous prediction for the spin-glass case.The following content is structured into three parts.First, the model is introduced and it is outlined howground state computations under changing boundaryconditions are used to produce domain wall excitations.Second, the results of the simulations will be presented.Finally we give a discussion. II. MODEL AND METHODSA. The Ising Spin Glass with Correlated Bonds
The Ising spin glass consist of Ising spins s m ∈ {± } on the sites m ∈ Λ of a two-dimensional lattice, i.e. Λ ⊂ Z . In this study only square systems are considered,such that the spin glass has L spins in each direction and a r X i v : . [ c ond - m a t . d i s - nn ] F e b | Λ | = L . The Hamiltonian is given by H J ( s ) = − (cid:88) { m , n }∈M J m , n s m s n , (1)where the sum runs over all pairs M of nearest-neighborspin sites with periodic boundary conditions (BC) in oneand free boundary conditions in the other direction. Thebonds, J m , n , which represent the interaction betweentwo spins, are random in strength and sign but remainconstant over time. Hence, one speaks of a quenched dis-order where the system is investigated under a fixed real-ization of the bonds, J . Here, the bonds originate froma Gaussian random field, J ( x ), which is a function of a continuous position x , and has zero mean, (cid:104) J ( x ) (cid:105) = 0and a covariance given by (cid:104) J ( x ) J ( x + r ) (cid:105) = (1 + r ) − a/ , (2) x , r ∈ R , a ≥ r = (cid:107) r (cid:107) . The entries of J aregiven by, J m , n = J (( m + n ) / a , is the only parameter tocontrol the correlation. For a = 0 one obtains the Isingmodel of a ferromagnet or antiferromagnet, respectively,depending on the bond realization. When a −→ ∞ theuncorrelated Ising spin glass model with Gaussian disor-der is recovered.To generate the correlated bonds numerically we uti-lized the Fourier Filtering Method (FFM) [31, 33]. TheFFM is a procedure to create stationary correlated ran-dom numbers of previously independent random num-bers. Because it is based on the convolution theorem itis possible to benefit from the computationally efficientFast Fourier Transform Algorithm [34]. For its imple-mentation we relied upon the functions of the FFTWlibrary, version 3.3.5 [35]. Figure 1 shows the averagebond correlation along the main axes of a system with | Λ | = 46 spins calculated by the estimator, C ( r ) = 1 |M (cid:48) ( r ) | (cid:88) { n , m }∈M (cid:48) ( r ) (cid:104) J m , n J m + r , n + r (cid:105) J . (3)Here (cid:104) ... (cid:105) J denotes the average with respect to the disor-der. M (cid:48) ( r ) ⊂ M contains those bonds { m , n } for whichthe bond { m + r , n + r } is also on the lattice. The factthat M (cid:48) ( r ) does not contain all the bonds { n , m } in M is due to the free boundary conditions in one direction.The point is that for vectors r which are not exactly par-allel to the free boundary, the bond { m + r , n + r } doesnot exist for all bonds of M . For the correlations shownin Figure 1 we only investigated the directions paralleland perpendicular to the free boundary. B. Ground States and Domain Walls
The nature of the ground state (GS) of the two-dimensional Ising spin glass is an intriguing subject on r . . . . . C periodic BCfree BC FIG. 1. (color online) Correlation of the bonds of a spinglass, calculated with Eq. (3), with a = 1 and L = 46 spinsin each direction, x and y . The system has free boundaryconditions in one direction and periodic boundary conditionsin the other direction. The bond correlation was generatedaccording to the Fourier Filtering Method (FFM) [31, 33] withperiodic boundary conditions in both directions, but the sizein the directions of free boundaries was chosen much largerthan the corresponding system size L . The line follows a fit oftype B C (1 + r ) − A C / , yielding B C = 1 . A C =0 . its own [36, 37]. Furthermore, GS computations of finitesystems are a well established tool [10, 11] to investigatethe glassy behavior of the model in the zero-temperaturelimit [38]. The GS is the spin configuration which min-imizes Eq. (1) for a given realization of the bonds. Incase of two-dimensional planar lattices there exist exactprocedures to generate the GS with a polynomial worst-case running time. This is in contrast to the three orhigher-dimensional variants which belong to the class ofNP-hard problems [39].There is more than one approach to compute the GSof two-dimensional planar spin glasses, such as the algo-rithm of Bieche et al. [40] and that of Barahona et al. [41]. The key idea of these algorithms is to create a map-ping from the original problem defined on the underlyinglattice graph of the spin glass onto a related graph whichis constructed in such a manner that the GS can be ex-tracted from a minimum-weight perfect matching , whichis polynomially computable. In this work we applied anansatz which includes Kasteleyn-city subgraphs into themapping process and thus is more efficient [42, 43] interms of speed and memory usage than the above men-tioned algorithms. This allowed us to investigate sys-tems up to a linear system size of L = 724 spins ineach direction without needing excessive computationalresources. For the computation of the minimum-weightperfect matching the Blossom IV algorithm [44] imple-mented by A. Rohe [45] was utilized. a −→ ∞ a = 1 a = 0 . − − − − FIG. 2. (color online) The right hand side of the figure showsone realization of the Gaussian random field. This is a func-tion of continuous position but here we show it discretizedwith a spatial resolution of half a lattice spacing. From thiscontinuous Gaussian random field, the bonds are extracted,at “half lattice points”. We show the Gaussian field for threedifferent correlation exponents. The corresponding GSs areon the left. The correlations exponents are given by a −→ ∞ (top), a = 1 (center) and a = 0 . s m = − s m = 1. In ferromag-netic order two neighboring spins have same sign and in an-tiferromagnetic order the sign alternates. The system size is L = 100. To study the spin glass at nonzero temperature we cre-ate domain wall (DW) excitations in the system. Thiswas done by computing GSs under periodic and antiperi-odic BCs also referred to as P-AP [46]. It works as fol-lows. First, a spin glass with periodic BCs in one di- rection and free BCs in the other direction is generatedwith quenched disorder, J (p) . Then, its GS configura-tion, s (p)gs , is computed. Next, the periodic BCs are re-placed by antiperiodic BCs by reversing the sign of onecolumn of bonds parallel to the direction of periodicity,which leads to J (ap) . Afterwards, the new GS configura-tion, s (ap)gs , is calculated. The change of the BCs imposesa DW of minimal energy between the two spin configu-rations s (p)gs and s (ap)gs . The energy of the DW is givenby ∆ E = H J (ap) (cid:16) s (ap)gs (cid:17) − H J (p) (cid:16) s (p)gs (cid:17) . (4)The geometrical structure of a DW can be char-acterized by the number D L of bonds which are in-cluded in the surface. To avoid including those bondswhich are a direct result of the different BCs, the sur-face is defined to consist of those bonds which fulfill J (p) m , n s (p) m s (p) n J (ap) m , n s (ap) m s (ap) n <
0, where m , n runs overall unordered pairs of nearest neighbor lattice sites [21]. III. RESULTS
We have obtained exact ground states for systems withcorrelation exponents in the range a ∈ [10 − , ∞ ], where a = ∞ corresponds to independently sampled bonds.Since the GS calculation requires only polynomial timeas a function of the system size, we were able to studysizes up to a large value of N ≡ L = 724 for each valueof a . For each value of a and L , we performed an averageover many realizations of the disorder, ranging from 10 realizations for the smallest sizes, 10000 realizations for L = 512, to 2000 realizations for the largest system size.Figure 2 provides a first impression how the bond cor-relation impacts the ordering of the GS. It is apparentthat for strong correlations there are large areas wherethe spins are either in ferromagnetic or antiferromagneticorder. This is related to there being large areas where,due to the correlation, the bonds have identical sign. A. Domain Wall Energy
Next, we look at the influence of bond correlation onthe properties of the previously discussed DW excita-tions. The absolute value of the DW energy is propor-tional to the coupling strength between block spins in thezero-temperature limit [14]. A stable order is possible ifthe average absolute value of the DW energy increaseswith the system size. In the uncorrelated case, a −→ ∞ ,one obtains a power law behavior [14, 20, 21, 47], (cid:104)| ∆ E |(cid:105) J ∼ L θ , (5)where θ is the stiffness exponent with its current bestestimate θ = − . θ < L h | ∆ E | i J a − →∞ a = a = . a = . a = . a = . a = . a = . a = . a = . a − . − . − . θ FIG. 3. (color online) Scaling of the average of the absolutevalue of the DW energy, (cid:104)| ∆ E |(cid:105) J , as a function of system size L for different values of a . The full lines are guides to theeyes only. The broken lines are fits of type (cid:104)| ∆ E |(cid:105) J = A θ L θ .The inset shows the values of θ which were obtained by fitsfor values of a ≥ .
9. The red line marks the value of thestiffness exponent for a −→ ∞ according to [21]. Figure 3 shows the impact of bond correlation on thescaling of (cid:104)| ∆ E |(cid:105) J . For a ≥ . a ≤ . (cid:104)| ∆ E |(cid:105) J initially grows withsystem size, but starts to decrease for larger sizes, i.e. thecurves exhibit a peak. The system size at the peak, L ∗ ,shifts to larger system sizes on decreasing the correlationexponent a . This will be analyzed below.First we show in figure 4 the results for the averageDW energy (cid:104) ∆ E (cid:105) J , which behaves in a somewhat similarmanner as the average of the absolute value.If spin-glass order existed in this model for some valueof a one would observe an increase (cid:104)| ∆ E |(cid:105) J in the limit L → ∞ , while at the same time (cid:104) ∆ E (cid:105) J would remainat, or converge to, zero as a function of L . The lat-ter indicates the absence of ferromagnetic order, but notnecessarily the absence of spin glass order because for aspin glass the change of the boundary conditions fromperiodic to antiperiodic is symmetric, i.e., could eitherincrease or decrease the GS energy, so the average would L − − h ∆ E i J a − → ∞ a = a = . a = . a = . a = . a = . a = 0 . a = . a = . FIG. 4. (color online) A log-log plot of the average DW energy (cid:104) ∆ E (cid:105) J as a function of the system size L for some values ofcorrelation exponent a . The inset shows the same data withlinear energy scale for a = 0 .
4, 0.5, 0.6 and 0.9 to highlightthe peak structure. The lines are guides to the eyes only. be zero even in the case of spin glass order. An increaseof both (cid:104)| ∆ E |(cid:105) J and (cid:104) ∆ E (cid:105) J for small values of L and a corresponds to a ferro/antiferromagnetic ordering onlocal length scales, which is visible in Fig. 2 and will bediscussed more below. Whether there is a true orderedphase for very small values of a will be discussed next.To track how the length scale of local order changes asa function of the correlation strength we measure the sizeat the peak, L ∗ , for both the domain-wall energy and theabsolute value, as a function of the correlation exponent a . Numerically, L ∗ was computed by fitting a parabolain the vicinity of the peaks of (cid:104)| ∆ E |(cid:105) J and (cid:104) ∆ E (cid:105) J . Asthe fit in figure 5 demonstrates, the data for L ∗ ( a ) is welldescribed by an exponential function L ∗ ( a ) = A L exp (cid:110) b L ( a − a crit ) − c L (cid:111) . (6)A non-zero value of a crit would indicate that an orderedphase exists for a < a crit . A true spin-glass phase wouldbe possible if a c rit for (cid:104) ∆ E (cid:105) J is smaller than a c rit for (cid:104)| ∆ E |(cid:105) J . We obtained values of a crit = 0 . .
10) for (cid:104) ∆ E (cid:105) J and 0 . .
17) for (cid:104)| ∆ E |(cid:105) J with quality of the fit Q = 0 .
87 and Q = 0 .
58, respectively. Thus, a zero valuefor the critical correlation exponent parameter a c seemslikely.To consider this further, we set a crit to zero and obtain . . . . a L ∗ /a ) h| ∆ E |i J h ∆ E i J FIG. 5. (color online) The system size where the peak of (cid:104)| ∆ E |(cid:105) J and (cid:104) ∆ E (cid:105) J occurs, denoted by L ∗ , as a function of a . The lines are fits according to Eq. (6) in the range a ≤ . L ∗ on a logarithmic scale as a function of1 /a exhibiting a straight-line behavior, and thus confirmingthe behavior obtained from the fit. the values of the other fit parameters, which here are A L = 3 . b L = 0 . c L = 2 . Q = 0 . (cid:104) ∆ E (cid:105) J and A L = 3 . b L = 0 . c L = 1 . Q = 0 .
63) for (cid:104)| ∆ E |(cid:105) J . Since the qualities of the fitsremain almost identical in comparison to a crit (cid:54) = 0 thedata is considered to be consistent with a crit = 0, whichwould imply that there is no global order when a > L ∗ ( a ) is also compatible with c L = 2. Fitsof this type, with a crit = 0 and c L = 2 fixed, have qualityof the fit larger than 0 .
4, which is reasonable. Note thatan exponential dependence of the “breakup” length scaleof ferromagnetic order as a function of disorder strengthwas also found in the two-dimensional random field Isingmodel, both at low temperatures [48] and in the GS [49].
B. Domain Wall Surface
The behavior of DW surfaces is regarded as one ofthe essential parameters which describe the propertiesof random systems [50, 51]. DWs separate spins in GSand reversed GS. Their surface is defined as those bondswhich belong to the DW, and we denote the surface sizeby D .In the uncorrelated case, a −→ ∞ , the average DWsurface size exhibits a power law [52, 53], (cid:104)D(cid:105) J ∼ L d s , (7)where d s is the fractal surface dimension, for which thebest numerical estimate at present is d s = 1 . a = 0 the system is a ferro/antiferromagnet L h D i J a −→ ∞ a = 0 . a = 0 . a = 0 . L FIG. 6. (color online) Scaling of the DW surface area (cid:104)D(cid:105) J fordifferent values of a . The broken black line shows the scalingof the DW surface in case of a ferro/antiferromagnet. Thefull lines follow fits according to Eq. (8) with L (fit)min = 8. L . . . d s a −→ ∞ a = 0 . a = 0 . a = 0 . FIG. 7. (color online) The fractal surface dimension, d s , ofthe DW surface area as a function of the smallest system size L of a fit window. and (cid:104)D(cid:105) J = L , implying that d s = 1, i.e., the surface isnot fractal here. Figure 6 shows the scaling of the DWsurface size for different values of a . In general, it can beseen that the correlation decreases the number of bondsin the DW surface. Even for a = 0 . a strictly zero.To characterize this behavior we fitted pure power-lawsto the data. For this purpose, we did not do a full fit toall data points, but used the sliding-window approach in-stead. Here, four values of (cid:104)D(cid:105) J which are adjacent interms of system size were grouped together in one fit win-dow, respectively. The independent variables of such afit window were given by ( L , L , L , L ) with L i < L i +1 , i = 1 , ...,
4. The dependent variables corresponded to thedata, i.e. (cid:104)D L i (cid:105) J . The smallest independent variable ofeach fit window is denoted as L .For each window, we fitted the power law to the dataresulting in a value of d s . The dependence of d s as afunction of L for different values of a can be found infigure 7. In the uncorrelated case d s decreases as a func-tion of L , whereas for strong correlations d s increases.In any case, for system sizes L <
32 and small values a ≤ d s on thesystem size. This motivated us to include corrections toscaling the power law behavior in Eq. (7) by considering[16, 21] (cid:104)D(cid:105) J = A D L d s (1 + B D L − ω s ) . (8)By using Eq. (8) the quality of the fit is larger than 0 . a ≤ . L (fit)min = 8. For larger values of a , thepure linear fit was always fine, given our statistical accu-racy. Figure 8 shows the resulting fractal surface dimen-sion for all considered values of a . At a ≈ . d s = 1 . a the fractal structure of the cluster changes,although there is no phase transition. Note that, in or-der to extract d s ( a ), we also performed fits of the form d s ( L ; a ) = d s ( a ) + κ d L γ d (not shown; κ d and γ d alsodepend on a ) to the sliding window fractal dimensionsshown in Fig. 7. The behavior of this extrapolated frac-tal dimension also exhibits the same notable decrease of d s for a ≤ .
4. Of course, it can not be ruled out thaton sufficiently large length scales, i.e., for much largersizes than are currently accessible, the fractal dimensionof the uncorrelated model would be recovered again andthus d s = 1 . a > C. Spin Correlations in the Ground State
In this section the effects of the bond correlations onspin correlations in the GS will be discussed. Note thatat zero temperature there is no thermal disorder and,in our study which uses bonds with a continuous distri-bution, the ground state is non-degenerate apart fromoverall spin inversion.The two-point spin correlation is given by (cid:104) s (gs) m s (gs) m + r (cid:105) J , where s (gs) m denotes a spin in the GS con-figuration. In the uncorrelated model (cid:104) s (gs) m s (gs) m + r (cid:105) J = 0if r (cid:54) = 0. The correlated bonds induce a localferro/antiferromagnetic order into the GS. When a = 0the system is a ferro/antiferromagnet and the order isglobal. In the ferromagnetic case s (gs) m s (gs) m + r = 1 and inthe antiferromagnetic case s (gs) m s (gs) m + r alternates betweenplus and minus one. Therefore, the GS spin correlation − − − a . . . . d s .
25 0 .
50 0 . . . . FIG. 8. (color online) The fractal surface dimension as afunction of a . The values for a ≥ . (cid:104)D(cid:105) J = A D L d s with smallest system sizeused in the fits L (fit)min = 128. The values for a ≤ . L (fit)min = 8. The red line marks thevalue of d s for a −→ ∞ according to [21]. can be estimated by G gs ( r ) = 1 | Λ (cid:48) ( r ) | (cid:88) m ∈ Λ (cid:48) ( r ) (cid:18) ˆ σ ( r ) + 12 (cid:19) (cid:104) s (gs) m s (gs) m + r (cid:105) J , (9)ˆ σ ( r ) = σ ( r ) σ ( r ) with σ ( r i ) = (cid:40) r i is even − r i is unevenand r = ( r , r ) ∈ Z . Λ (cid:48) ( r ) ⊂ Λ, similar to Eq. (3),contains those sites m for which m + r is on the lattice,given the free BC in one direction. Note that this defini-tion means that the correlation is measured for each site m ∈ Λ (cid:48) on one of the two sublattices of a checkerboardpartition of the square lattice, such that the correlationis insensitive to whether the order is ferromagnetic orantiferromagnetic.In figure 9 one can see the GS correlation for differentvalues of a . The data is well described by a scaling formof type G gs ( r ) ∼ r υ exp (cid:26) − (cid:18) rξ gs (cid:19) ϕ (cid:27) . (10)From the perspective of ordering we are especially in-terested in the correlation length, ξ gs , that provides thedistance over which spins are notably correlated. Thestraightforward method to extract ξ gs is by fitting thefunction of Eq. (10) to the data. The problem with suchan approach is that neither υ nor ϕ are known. Alsofinite-size corrections reduce the match between scalingform and actual data. Thus, we used different approachesto obtain ξ g s .First, we perform a separate fit for each value of a downto small correlations where the error bars start to exceed r − − − G g s a = . a = . a = . a = . a = . a = . FIG. 9. (color online) Spin correlation of the GS for differentvalues of a . The correlation was computed by utilizing theestimator of Eq. (9) along the main axes. The black lines arefits according to Eq. (10) when using a multifit, i.e., fittingtwo exponents υ and ϕ and many values of ξ gs ( a ) simultane-ously to the correlation obtained for all values of a . one quarter of the correlation value. We observed thatfor a ≥
2, the values of the exponents υ and ϕ did notchange much, while for smaller values of a the exponentswere a bit smaller. Therefore, we fixed the exponents tothe (averaged) values seen for a ≥ ξ gs .Second, we also performed a multifit , i.e., we fitted thecorrelation function simultaneously [54] with one valuefor υ , one value for ϕ and values of the correlation lengthat many values of a . The results of this multifit are alsoshown in Fig. 9. The obtained values for ξ gs for thesetwo fitting approaches are shown in Fig. 10. As can beseen, the results from fixing the values of the exponentsand from using the multifit approach do not differ much.Before we discuss the behavior of ξ gs ( a ), we describethe third approach we have used to estimate the correla-tion length. Here, we used the integral estimator that wasintroduced in Ref. [55]. It presupposes that for r ≤ ξ gs the correlation function is dominated by a power law oftype r − υ , whereas for r > ξ gs the correlation is negligible.As a consequence, the integral I k = (cid:90) ∞ d r r k G gs ( r ) (11)is given by I k ∝ ( ξ gs ) k +1 − υ and thus ξ ( k,k +1)gs := I k +1 I k ∝ ξ gs . (12)This result would be exact and independent of k if thecorrelation actually was only a power law. For real cor-relations, the value of k dictates which part of the corre-lation function contributes most to the integral. Because such a scaling approach is only valid when r is muchlarger than the lattice constant a high value of k reducesthe systematic error of the method. On the other hand,large values of k increase the statistical error of I k . Fol-lowing the recommendation of [56] we used ξ (1 , as acompromise. Note that since the statistical error of themeasured correlation grows with distance r , one usuallydefines a cutoff distance up to which the data is directlyused for the integral. Similar to [56] we specified thiscutoff distance as the value of r where G gs is smallerthan three times its error. For values of r larger thanthis cutoff distance we computed the integral up to themaximal length of L from fits to Eq. (10). The startvalue of these fits were set to r (fit)min ≥
2. Because it wasobserved that the correlation decays slightly differentlyalong the directions with free and periodic BCs, the com-putations of the GS were also done with an independentset of simulations for full free BCs.Next, the correlation length ξ (1 , was extracted fromthe average of the GS correlation, G gs , along the mainaxes. The statistical error of ξ (1 , was estimated bybootstrapping [54, 57] and the integrals according to Eq.(11) were computed by utilizing the midpoint integrationrule. For small values of a , the contribution to I k fromthe integral beyond the cutoff gets increasingly large. Forinstance, when a = 1 .
15 this contribution made up ap-proximately 4% of the value of I . Hence, to estimatethe total error, we added an extra systematic error tothe statistical error. This was done by analyzing thevalue of ξ (1 , for two other choices of the cutoff distance,i.e., being the distance where the error of the correlationfunction is two or four times larger than its estimate, re-spectively. The maximal deviation of these two valuesfrom the standard definition, which uses a magnitude ofthree error bars to define cutoff, was set to be the system-atic error. Furthermore, because the statistical error of I k grows large for small values of the correlation exponent a and the system has to be sufficiently large to neglectboundary effects, values of a < . ξ (1 , as a function of the correla-tion exponent. In contrast to the data for the lengthscales L (cid:63) , the plot shows that the behavior of noneof the three correlation lengths that we have defined, ξ gs (fit) , ξ gs (multifit) and ξ (1 , , is close to exponen-tial. Nonetheless, for a → ξ gs = A ξ ( a − a crit ) − d ξ exp (cid:8) b ξ ( a − a crit ) − c ξ (cid:9) . (13)The exponential in Eq. (13) is consistent with the pre-vious results for the length scale L ∗ and will dominatefor a →
0. The power-law part is chosen to describe thebehavior for large values of a .We have fitted the data to this function. For example,for the data obtained from the multifit we again obtained /a ) L ∗ a nd ξ g s ξ (1 , ξ gs (multifit) ξ gs (fit) L ∗ of h| ∆ E L |i L ∗ of h ∆ E L i FIG. 10. (color online) The correlation length of the GS, ξ (1 , , as a function of 1 /a , obtained by three different ap-proaches. For comparison the length scales L (cid:63) of the maximaare shown again. a small value of a crit = 0 . a crit ≡ A ξ = 4 . d ξ = 0 . b ξ =1 . c ξ = 1 . c ξ , are found for the othertwo definitions of ξ gs that we used. This means that forlarge values of the correlation exponent a the behavior ofthe correlation length as function of a seems to be betterdescribed by a power law. Also, the behavior of the peaklengths and of the ferromagnetic correlation lengths differa lot, but it is still possible that for very small values of a , an exponential dependence of ξ gs on 1 /a would berecovered. Unfortunately, to investigate this issue muchlarger system sizes would have to be treated, well beyondcurrent numerical capabilities. IV. DISCUSSION
The standard two-dimensional Ising spins glass doesnot exhibit a finite-temperature spin glass phase in con-trast to the three or higher dimensional cases [18–20].This work deals with the question how long-range cor-related bonds influence this characteristic. Therefore,the ordering behaviour of the two-dimensional Ising spinglass with spatially long-range correlated bonds is studiedin the zero-temperature limit. The bonds are drawn froma standard normal distribution with a two-point correla-tion for bond distance r that decays as ( r +1) − a/ , a ≥ a = 0, the system is either aferromagnet or antiferromagnet, depending on the bondrealization. For a → ∞ the uncorrelated EA model isrecovered.For 0 < a < ∞ we observed that the correlationhas local effects on the zero-temperature ordering be-haviour. The correlation locally effects the average value of the bonds as well as their standard deviation for eachindividual realization of the disorder. These parame-ters are decisive to distinguish between a spin glass orferro/antiferromagnet in case of the uncorrelated model[25, 28]. In correspondence to that, the spin correla-tion of the GS reveals how the correlation induces a lo-cal ferro/antiferromagnetic order into the GS. This is re-flected by a growing correlation length ξ gs ( a ) when de-creasing a .Complementary results to the direct study of the GSspin configurations were obtained by investigating DWexcitations. The average of the absolute value of theDW energy can be interpreted as the coupling strengthbetween block spins at zero temperature [14]. We found,that for strong bond correlation, the average of absolutevalue of the DW energy initially increases as a function ofthe system size up to a peak, and then decreases. Sincewe made the same observation for the actual DW energyit shows that the increase of the absolute value of the DWenergy is a consequence of local ferro/antiferromagneticorder of the system in GS. The system size where thepeak occurs, L ∗ , is interpreted as the length scale of lo-cal order. For small values of the correlation exponent a ,both L ∗ ( a ) and the correlation length of the GS, ξ gs ( a ),can be described by an exponential divergence. Interest-ingly, a similar exponential length scale was also foundin the two-dimensional random field Ising model by GScomputations [49] and at low temperatures [48]. In thesestudies the length scale of ferromagnetic order was ex-amined as a function of the standard deviation of therandom magnetic field.For the two-dimensional Ising spin glass, the distri-bution of the absolute value of the domain wall energyis “universal” with respect to the initial bond distribu-tion. This means for any continuous, symmetric bonddistribution with sufficiently small mean and finite highermoments, the absolute value of the domain wall energyshould approach the same scaling function [13, 14, 58].Thus, we expect the same kind of universality for ourmodel.The stiffness exponent θ describes the scaling of thewidth of this distribution, and is related to the criticalexponent ν describing the divergence of the correlationlength as T → ν = − /θ [14, 24]. At this zero tem-perature transition ν is the only independent exponent.Therefore, any bond correlation which leaves the stiffnessexponent unchanged does not influence the universalityof the model. From the data of L ∗ it is expected thatthere is no global ordered phase for a >
0. This impliesthat the stiffness exponent is negative for a >
0. Further-more, for values of a in the range a ≥ . a less than 0 . (cid:104)D(cid:105) J ∼ L d s , where d s = 1 . a ≥ .
5. At a ≈ . a ≤ . (cid:104)D(cid:105) J = A D L d s (1 + B D L − ω s ). The de-crease in d s implies that for strong correlations DWs withshorter lengths are energetic favorable. In the extremecase when a = 0 the system is a ferro/antiferromagnetand thus d s = 1. Of course, it can not be ruled out thatthe decline in d s is local and on sufficiently long lengthscales the pure power law with d s = 1 . a > d s =1+3 / (4(3+ θ )) [59]. According to highly accurate numer-ical results [21] this equation is probably not exact. Ourresults for the fractal surface dimension d s = 1 . Q = 0 .
99) and the stiffness exponent θ = − . Q = 0 .
13) deviate by approximately 6.5 standard de-viations from the mentioned conjecture, since d s − − / (4(3 + θ )) = − . d s and θ which exist because both values were obtained from the same data set. Nonetheless, our results also support theconclusion that the proposed scaling relation is not exact.In conclusion, it is observed that correlation among thebonds has strong effect on the ordering on local lengthscales, inducing ferro/antiferromagnetic domains into theGS. The length scale of local ferro/antiferromagnetic or-der diverges exponentially when the correlation exponentapproaches zero. The fractal surface dimension decreasesfor strong correlations on the studied length scales. Nosignature of a spin-glass phase at finite temperature isobserved. ACKNOWLEDGMENTS
The simulations were performed at the HPC clusterCARL, located at the University of Oldenburg (Ger-many) and funded by the DFG through its Major Re-search Instrumentation Programme (INST 184/157-1FUGG) and the Ministry of Science and Culture (MWK)of the Lower Saxony State. [1] E. Vincent and V. Dupuis, in
Frustrated Materials andFerroic Glasses , edited by T. Lookman and X. Ren(Springer, Cham, 2018) pp. 31–56.[2] S. F. Edwards and P. W. Anderson, J. Phys. F , 965(1975).[3] D. Sherrington and S. Kirkpatrick, Phys. Rev. Lett. ,1792 (1975).[4] D. L. Stein and C. M. Newman, Spin Glasses and Com-plexity (Princeton University Press, Princeton, 2013).[5] M. M´ezard, G. Parisi, and M. Virasoro,
Spin Glass The-ory and Beyond (World Scientific, Singapore, 1987).[6] A. P. Young, ed.,
Spin Glasses and Random Fields (World Scientific, Singapore, 1998).[7] E. Bolthausen and A. Bovier, eds.,
Spin Glasses (Springer, Berlin, Heidelberg, 2007).[8] H. Nishimori,
Statistical Physics of Spin Glasses and In-formation Processing: An Introduction (Oxford Univer-sity Press, Oxford, 2001).[9] V. Dotsenko,
An introduction to the theory of spin glassesand neural networks (World Scientific, Singapore, 1994).[10] A. K. Hartmann and H. Rieger,
Optimization Algorithmsin Physics (Wiley-VCH, Weinheim, 2002).[11] A. K. Hartmann and H. Rieger, eds.,
New OptimizationAlgorithms in Physics (Wiley-VCH, Weinheim, 2004).[12] Z. Zhu, A. J. Ochoa, and H. G. Katzgraber, Phys. Rev.Lett. , 077201 (2015).[13] W. L. McMillan, J. Phys. C , 3179 (1984).[14] A. J. Bray and M. A. Moore, in Heidelberg Colloquium onGlassy Dynamics , edited by J. van Hemmen and I. Mor-genstern (Springer, Berlin, Heidelberg, 1987) pp. 121–153. [15] D. S. Fisher and D. A. Huse, Phys. Rev. B , 386 (1988).[16] A. K. Hartmann and M. A. Moore, Phys. Rev. Lett. ,127201 (2003).[17] A. K. Hartmann and M. A. Moore, Phys. Rev. B ,104409 (2004).[18] A. P. Young and R. B. Stinchcombe, J. Phys. C , 4419(1976).[19] B. W. Southern and A. P. Young, J. Phys. C , 2179(1977).[20] A. K. Hartmann and A. P. Young, Phys. Rev. B ,180404(R) (2001).[21] H. Khoshbakht and M. Weigel, Phys. Rev. B , 064410(2018).[22] F. Parisen Toldin, A. Pelissetto, and E. Vicari, Phys.Rev. E , 051116 (2011).[23] C. K. Thomas, D. A. Huse, and A. A. Middleton, Phys.Rev. Lett. , 047203 (2011).[24] L. A. Fernandez, E. Marinari, V. Martin-Mayor,G. Parisi, and J. J. Ruiz-Lorenzo, Phys. Rev. B ,024402 (2016).[25] D. Sherrington and W. Southern, J. Phys. F , L49(1975).[26] C. Amoruso and A. K. Hartmann, Phys. Rev. B ,134425 (2004).[27] R. Fisch and A. K. Hartmann, Phys. Rev. B , 174415(2007).[28] C. Monthus and T. Garel, Phys. Rev. B , 184408(2014).[29] H. Fajen, A. K. Hartmann, and A. P. Young, Phys. Rev.E , 012131 (2020). [30] A. Weinrib and B. I. Halperin, Phys. Rev. B , 413(1983).[31] B. Ahrens and A. K. Hartmann, Phys. Rev. B , 144202(2011).[32] T. Nattermann, in Spin Glasses and Random Fields ,edited by A. P. Young (World Scientific, Singapore, 1998)pp. 277–298.[33] H. A. Makse, S. Havlin, M. Schwartz, and H. E. Stanley,Phys. Rev. E , 5445 (1996).[34] J. W. Cooley and J. W. Tukey, Mathematics of Compu-tation , 297 (1965).[35] M. Frigo and S. G. Johnson, Proceedings of the IEEE , 216 (2005).[36] C. M. Newman and D. L. Stein, Phys. Rev. Lett. ,3966 (2000).[37] L.-P. Arguin and M. Damron, Annales de l’I.H.P. Prob-abilit´es et statistiques , 28 (2014).[38] A. K. Hartmann, in Rugged Free Energy Landscapes ,edited by W. Janke (Springer, Berlin, Heidelberg, 2008)pp. 67–106.[39] F. Barahona, J. Phys. A , 3241 (1982).[40] L. Bieche, J. P. Uhry, R. Maynard, and R. Rammal, J.Phys. A , 2553 (1980).[41] F. Barahona, R. Maynard, R. Rammal, and J. Uhry, J.Phys. A , 673 (1982).[42] G. Pardella and F. Liers, Phys. Rev. E , 056705 (2008).[43] C. K. Thomas and A. A. Middleton, Phys. Rev. B ,220406(R) (2007).[44] W. Cook and A. Rohe, INFORMS Journal on Computing , 077201 (2002).[47] A. K. Hartmann, A. J. Bray, A. C. Carter, M. A. Moore,and A. P. Young, Phys. Rev. B , 224401 (2002).[48] E. Pytte and J. F. Fernandez, J. Appl. Phys. , 3274(1985).[49] E. T. Sepp¨al¨a, V. Pet¨aj¨a, and M. J. Alava, Phys. Rev.E , R5217 (1998).[50] A. J. Bray and M. A. Moore, Phys. Rev. Lett. , 57(1987).[51] C. M. Newman and D. L. Stein, in Spin glasses , editedby E. Bolthausen and A. Bovier (Springer, Berlin, Hei-delberg, 2007) pp. 145–158.[52] N. Kawashima, J. Phys. Soc. Jpn. , 987 (2000).[53] O. Melchert and A. K. Hartmann, Phys. Rev. B ,174411 (2007).[54] A. K. Hartmann, Big Practical Guide to Computer Sim-ulations (World Scientific, Singapore, 2015).[55] F. Belletti, M. Cotallo, A. Cruz, L. A. Fernandez,A. Gordillo-Guerrero, M. Guidetti, A. Maiorano, F. Man-tovani, E. Marinari, V. Martin-Mayor, A. M. Sudupe,D. Navarro, G. Parisi, S. Perez-Gaviro, J. J. Ruiz-Lorenzo, S. F. Schifano, D. Sciretti, A. Tarancon,R. Tripiccione, J. L. Velasco, and D. Yllanes, Phys. Rev.Lett. , 157201 (2008).[56] F. Belletti, A. Cruz, L. A. Fernandez, A. Gordillo-Guerrero, M. Guidetti, A. Maiorano, F. Mantovani,E. Marinari, V. Martin-Mayor, J. Monforte, et al. , J.Stat. Phys. , 1121 (2009).[57] A. P. Young,
Everything You Wanted to Know AboutData Analysis and Fitting but Were Afraid to Ask (Springer, Cham, 2015).[58] C. Amoruso, E. Marinari, O. C. Martin, and A. Pagnani,Phys. Rev. Lett. , 087201 (2003).[59] C. Amoruso, A. K. Hartmann, M. B. Hastings, and M. A.Moore, Phys. Rev. Lett.97