9512.net
甜梦文库
当前位置:首页 >> >>

A New Method to Calculate the Spin-Glass Order Parameter of the Two-Dimensional +-J Ising M


arXiv:cond-mat/9908370v3 [cond-mat.dis-nn] 24 Apr 2000

A new method to calculate the spin-glass order parameter of the two-dimensional ±J Ising model?
Hidetsugu Kitatani ? and Akira Sinada
Department of Electrical Engineering, Nagaoka University of Technology, Nagaoka, Niigata 940-2188, Japan Abstract. A new method to numerically calculate the nth moment of the spin overlap of the two-dimensional ±J Ising model is developed using the identity derived by one of the authors (HK) several years ago. By using the method, the nth moment of the spin overlap can be calculated as a simple average of the nth moment of the total spins with a modi?ed bond probability distribution. The values of the Binder parameter etc have been extensively calculated with the linear size, L, up to L = 23. The accuracy of the calculations in the present method is similar to that in the conventional transfer matrix method with about 105 bond samples. The simple scaling plots of the Binder parameter and the spin-glass susceptibility indicate the existence of a ?nite-temperature spin-glass phase transition. We ?nd, however, that the estimation of Tc is strongly a?ected by the corrections to scaling within the present data (L ≤ 23). Thus, there still remains the possibility that Tc = 0, contrary to the recent results which suggest the existence of a ?nite-temperature spin-glass phase transition.

PACS numbers: 75.50.Lk,02.70.Lq,64.60.Cn,05.50.+q 1. Introduction

Over the last two decades, investigations of spin-glass problems have been extensively performed[1-20]. It is now widely believed that the three-dimensional ±J Ising model shows a ?nite-temperature spin-glass phase transition [1-6], while the critical temperature of the two-dimensional ±J Ising model is zero [5-11]. Most of these studies have been done using Monte Carlo simulations, where the thermal relaxation time in the simulations becomes very large in the low-temperature region. This makes the investigations of the two-dimensional models rather di?cult, since the calculations of the physical quantities have to be performed at very low temperature. In previous studies, the data in the ?nite-size scaling analysis were in good agreement with a scaling function with the critical temperature, Tc = 0. The results, however, have not completely excluded the possibility of a ?nite-temperature spin-glass phase transition. Recently, Shirakura et al[12-15] have deduced the existence of a ?nite-temperature spin-glass phase transition of the two-dimensional models using extensive Monte Carlo simulations. To clarify the critical properties of the two-dimensional ±J Ising model, more precise results in the low-temperature region are necessary. The transfer matrix method is free from the problem of thermal equilibration, and has been very successfully used to determine the ferromagnetic-nonferromagnetic
? E-mail address:kitatani@vos.nagaokaut.ac.jp

2 phase boundary of the two-dimensional ±J Ising model in the p ? T plane [16, 17] (p is the concentration of the ferromagnetic bond). But, for the problem of the spin-glass phase transition, the transfer matrix method has only been used for the calculations of defect energies and correlation functions [9, 10], and has not been widely used for direct calculations of the nth moment of the spin ovelap (the spinglass susceptibility, the Binder parameter, etc) since, when we use real replicas in the calculations, the maximum linear size applicable is one-half of that in the calculations of the nth moment of the total spins. Thus, so far, no extensive result for the spinglass phase transition with direct calculation of the spin-glass susceptibility etc has been given by the transfer matrix method. In this paper, we present a new method to numerically calculate the nth moment of the spin overlap of the two-dimensional ±J Ising model, using the identity derived by one of the present authors several years ago [19]. In this identity, the nth moment of the spin overlap is transformed as a simple average of the nth moment of the total spins with a modi?ed bond probability distribution. Following a newly developed process, explained in section 2, we successively make the bond con?gurations according to the modi?ed bond probability distribution using the Monte Carlo technique. In each bond con?guration, we calculate the nth moment of the total spins by the transfer matrix method. We have performed extensive calculations of the nth moment of the spin overlap up to the linear size, L = 23. The accuracy of the calculations in the present method is similar to that in the conventional transfer matrix method with about 105 bond samples. Thus, the statistical errors in the present study are about an order of magnitude smaller than those in previous studies. Therefore, we can analyse the obtained data in detail using ?nite-size scaling analysis including the corrections to scaling. Our results show that the estimation of Tc is strongly a?ected by the corrections to scaling in the two-dimensional ±J Ising model. Thus, there still remains the possibility that Tc = 0, contrary to the recent results by Shirakura et al [12, 13, 15] which suggest the existence of a ?nite-temperature spin-glass phase transition. 2. New calculation method for the spin-glass order parameter

We consider the two-dimensional ±J Ising model on an L × L square lattice with only nearest neighbour interactions. The Hamiltonian is written as follows: H=?
(ij )

τij Si Sj ,

(1)

where Si = ±1, and the summation of (ij ) runs over all the nearest neighbours. We take the skew boundary condition in one direction, and the free boundary condition in the other direction. Each τij is determined according to the following probability distribution: P (τij ) = pδ (τij ? 1) + (1 ? p)δ (τij + 1). (2)

In this paper, we make J = 1. We de?ne the spin overlap, Q, between two replicas in each bond con?guration as
N

Q=
i=1

α β Si Si ,

(3)

3 where α and β denote the two replicas, and N is the total number of spins. When we de?ne Kp as p exp(2Kp ) = , (4) 1?p the 2nth moment of the spin overlap is written as [< Q2n >T,{S α ,S β } ]p = 1 (2 cosh(Kp ))NB exp(Kp
τij =±1 (ij )

τij ) < Q2n >T,{S α ,S β } ,

(5)

where < · · · >T,{S α ,S β } denotes the thermal average both for the ”S α ” and ”S β ” spins in a bond con?guration, {τ } at temperature, T . [· · ·]p denotes the con?gurational average at the ferromagnetic bond concentration, p, and NB is the number of bonds [18]. By the use of a local gauge transformation, an identity has been derived [19]: [< Q2n >T,{S α ,S β } ]p = 1 (2 cosh(Kp ))NB exp(K
τij =±1 (ij )

τij )

Z (Kp , {τ }) < M 2n >T,{S } ,(6) Z (K, {τ })

where M denotes the total spins,
N

M=
i=1

Si ,

(7)

< · · · >T,{S } denotes the thermal average for the ”S ” spins at temperature, T , and Z (K, {τ }) is the partition function at the inverse temperature, K (= 1/T ), with the bond con?guration, {τ } . We show the summary of the derivation of equation (6) in the appendix. (For the details of the derivation, see [19].) When we de?ne the modi?ed probability distribution, P2 (K, Kp , {τ }), for the bond con?guration, {τ }, as P2 (K, Kp , {τ }) = we can then write [< Q2n >T,{S α ,S β } ]p = {< M 2n >T,{S } }K,Kp , (9) 1 exp(K (2 cosh(Kp ))NB τij )
(ij )

Z (Kp , {τ }) , Z (K, {τ })

(8)

where {· · ·}K,Kp denotes the con?gurational average by the modi?ed bond probability distribution. That is, [< Q2n >T,{S α ,S β } ]p at temperature, T , with the ferromagnetic bond concentration, p, is transformed into the con?gurational average of < M 2n >T,{S } by the modi?ed bond probability disribution, P2 (K, Kp , {τ }). Similarly, we can get the following identity [19]: [< M 2n >T,{S } ]p = {< M 2n >Tp ,{S } }K,Kp , (10)

where Tp = 1/Kp . (Note that the above argument applies to any dimension.) Hereafter, we explain a new approach to numerically calculate the values of [< Q2n >T,{S α ,S β } ]p , using equation (9). To realize the bond con?guration with the modi?ed bond probability distribution, P2 (K, Kp , {τ }), we use the conventional Monte Carlo technique. We de?ne W (τij → ?τij , {τ }′ ) as the transition probability by which the value of the bond, τij , changes. To guarantee that the stationary probability distribution becomes P2 (K, Kp , {τij }), the following detailed balance must be satis?ed: P2 (K, Kp , τij , {τ }′ )W (τij → ?τij , {τ }′ ) = P2 (K, Kp , ?τij , {τ }′ )W (?τij → τij , {τ }′ ) (11)

4 Using equation (8), we obtain cosh(2Kp ) ? sinh(2Kp ) < τij Si Sj >Tp ,{S } W (τij → ?τij , {τ }′ ) = exp(?2Kτij ) W (?τij → τij , {τ }′ ) cosh(2K ) ? sinh(2K ) < τij Si Sj >T,{S } (12)

Namely, when we can calculate < Si Sj >T,{S } in a particular bond con?guration, we can estimate the transition probability. The processes to calculate [< Q2n >T,{S α ,S β } ]p are as follows: 1)We start from an arbitrary bond con?guration. 2)Using the transfer matrix method, we exactly calculate the value of < Si Sj >T,{S } , and successively ?ip the bond, τij , according to the transition probability, W (τij → ?τij , {τ }′ ), using the conventional Monte Carlo technique. 3)We continue process 2) until the modi?ed bond probability distribution, P2 (K, Kp , {τ }) is realized. 4)We calculate the value of < M 2n >T,{S } for each bond con?guration using the transfer matrix method. 5)We repeat processes, 2) and 4). 6)Finally, the simple average of < M 2n >T,{S } gives the value of [< Q2n >T,{S α ,S β } ]p with the bond probability distribution, P (τij ). We now show the e?ciency of this method. We de?ne na , nb and nc as the number of initial Monte Carlo skip steps, the number of Monte Carlo steps where we calculate < M 2n >T,{S } , and the number of independent Monte Carlo runs, respectively. In all the calculations, we have evaluated {τij (0)τij (t)}K,Kp using the statistical dependence time method [22], and ?nd that the relaxation time of {τij (0)τij (t)}K,Kp is always very small even when compared with one Monte Carlo step time. That is, only several tens of initial skip steps are enough to realize the stationary bond probability distribution. For example, we have compared the values of the spin-glass susceptibility χSG (= [< Q2 >T,{S α ,S β } ]p /N ) calculated by the present method and that by the conventional transfer matrix method using real replicas. Table 1 shows the results at T = 0.1 for L = 7. The calculations by the conventional transfer matrix method have been done with 105 independent bond con?gurations. The error bars of the present methods have been estimated in the same way as those of conventional Monte Carlo simulations. From table 1, we can see that all the data are consistent within the error bars, and the size of the error bars of all the calculations are of the same magnitude. Consequently, we ?nd that only 20 steps are enough for the initial Monte Carlo skip steps. Furthermore, we have examined whether equation (9) holds or not at p = 0.5 ? 0.95, T = 0.1 ? 0.5 for L = 7. We have also examined whether equation (10) holds or not at p = 0.8 ? 0.9, T = 0.1 ? 0.4 for L = 15. All the results are consistent in a statistical sense, from which we conclude that the present method is usable and not a?ected by systematic errors. 3. The spin-glass phase transition of the two-dimensional ±J Ising model

We have extensively investigated the two-dimensional ±J Ising model for p = 0.5. The results of the asymmetric case (p > 0.5) will be given in a subsequent paper [20]. For p = 0.5, we have calculated [< Q2n >T,{S α ,S β } ]p at T = 0.1 ? 0.5 with the linear size L = 7 ? 23. The calculations have been performed with (na , nb , nc ) = (200, 1000, 100) for L ≤ 21 and (200,200,480) for L = 23.

5
Table 1. The values of χSG with various (na , nb , nc ) at T = 0.1 for L = 7. The value of ? is calculated by the conventional transfer matrix method using real replicas with 105 bond samples. (na , nb , nc ) (2000,10000,10) (200,1000,100) (20,100,1000) (20,20,5000) ? χSG 29.076(24) 29.084(38) 29.094(35) 29.045(34) 29.037(35)

The energy gap between the ground state and the ?rst excitation state is two in the unit of the interaction strength. Thus, in ?nite systems, any physical quantity at a very low temperature must saturate to its value at T = 0. We show the temperature dependence of the spin-glass susceptibility, χSG , in ?gure 1. We ?nd, indeed, that the values of χSG for each L show the strong saturation near T = 0, and the tendency becomes clearer as the system size becomes smaller, as has already been pointed out by several authors [8, 12, 14]. The Binder parameter [21] gL = [< Q4 >T,{S α ,S β } ]p 1 (3 ? ) 2 [< Q2 >T,{S α ,S β } ]2 p (13)

is widely used for the determination of the critical temperature. The value of the Binder parameter becomes asymptotically size independent for large systems at the critical temperature. Therefore, the point where this quantity becomes asymptotically size independent gives an estimation of the critical point. The simple plot of the Binder parameter versus temperature is shown in ?gure 2. We can clearly see that the data for di?erent sizes intersect at almost the same temperature, T = 0.3, and the size dependence of the intersection points is very small. We cannot, however, immediately conclude that the spin-glass phase transition occurs at T ? 0.3, since the intersection might be due to the strong saturation of the data near T = 0 [8, 12, 14]. Therefore, we perform scaling analyses for gL and χSG . There is no general rule to avoid the disturbance from the saturation mentioned above in the scaling analyses. Here, we adopt a criterion that every data point is not used all through the scaling analyses, when the value of χSG increases less than 3% in the temperature interval, ?T = 0.05. Although the criterion, which we have determined from the observation in ?gure 1, seems rather arti?cial, we believe that this criterion systematically removes the saturation to T = 0 in a certain sense. Consequently, we use, for example, the data with T ≥ 0.35 for L = 7, and with T ≥ 0.2 for L = 23. First, we perform the scaling analyses without including the corrections to scaling. In this case, the Binder parameter, gL , has the scaling form gL = g ?(L1/ν (T ? Tc )), and that of the spin-glass susceptibility, χSG , is χSG = L2?η χ ?(L1/ν (T ? Tc )), (15) (14)

where ν is the critical exponent of the spin-glass correlation length, and η is the critical exponent which describes the decay of the correlation at the critical temperature. Figure 3 shows the best-?t scaling plot of the Binder parameter, which indicates that

6 Tc ? 0.3 and the critical exponent ν ? 1.3. We can see that the data at T < Tc = 0.3 and T ≥ Tc ?t rather well on one scaling function, which indicates that the spinglass phase transition of the two-dimensional ±J Ising model is a conventional phase transition, and there exists a ?nite long range order at T < Tc . The best-?t scaling plot of the spin-glass susceptibility is also shown in ?gure 4, which indicates that the critical exponent η ? 0.225. The estimated values of Tc and the critical exponents are similar to those determined by Shirakura et al[12]. Figures 5 and 6 show the scaling plots of the Binder parameter and the spin-glass susceptibility, assuming Tc = 0, ν = 2.6 and η = 0.2 [8], where we clearly see the systematic deviations. Namely, the scaling plots without including the corrections to scaling strongly indicate the existence of a ?nite-temperature spin-glass phase transition. Next, we perform the scaling analyses including the corrections to scaling. We take the scaling forms of the Binder paremeter and the spin-glass susceptibility as a (16) gL = g ?(L1/ν (T ? Tc ))(1 + ω ), L and b χSG = L2?η χ ?(L1/ν (T ? Tc ))(1 + ω′ ). (17) L There is little quantitative change in ?gures 3 and 4, even though we ?t the data using equations (16) and (17). Thus, when we assume Tc = 0.3, the e?ect of the corrections to scaling is rather small. On the other hand, assuming that Tc = 0, ν = 2.6 and η = 0.17 [8], the data of the Binder parameter and the spin-glass susceptibility ?t very well on one scaling function, respectively, as shown in ?gures 7 and 8 with ω = ω ′ = 0.5 and a = b = ?0.3, although the data with a small linear size, L = 7, deviate from the scaling function. (To ?t the data, we use η = 0.17, which is, for example, consistent with the result in [8], η = 0.2 ± 0.05.) Thus, including the corrections to scaling, both Tc = 0 and Tc ? 0.3 are consistent with the scaling analyses. Furthermore, we ?nd that every temperature for 0 ≤ T ≤ 0.3 might become the critical temperature, Tc , in this scaling form. Consequently, we ?nd that the estimation of the value of Tc is strongly a?ected by the corrections to scaling in the two-dimensional ±J Ising model within the present data (L ≤ 23). 4. Conclusions

We have developed a new method to numerically calculate [< Q2n >T,{S α ,S β } ]p of the two-dimensional ±J Ising model, where, using a local gauge transformation, [< Q2n >T,{S α ,S β } ]p can be calculated as a simple average of < M 2n >T,{S } with a modi?ed bond probability distribution. By using the present method, we have extensively calculated the values of [< Q2n >T,{S α ,S β } ]p , where the statistical errors become about an order of magnitude smaller than in previous studies, and we have investigated the scaling analyses including the corrections to scaling. By using the scaling analyses without including the corrections to scaling, our data strongly indicate a ?nite-temperature spin-glass phase transition. We ?nd, however, that the estimation of Tc is strongly a?ected by the corrections to scaling within the data with L ≤ 23, and that every temperature for 0 ≤ T ≤ 0.3 might be able to become the critical temperature. Consequently, our results indicate that there still remains the possibility that Tc = 0, contrary to the recent results of Shirakura et al [12, 13, 15] which suggest the existence of a ?nite-temperature spin-glass phase transition.

7 Acknowledgments The authors thank Dr. Shirakura for his useful comments and discussions. The calculations were made on the FACOM VPP500 at the Institue for Solid State Physics, University of Tokyo and on the HITAC SR8000 at the University of Tokyo. This work was partially supported by a Grant-in-Aid for Exploratory Research from the Ministry of Education, Science, Sports and Culture. Appendix A. In this appendix, we brie?y explain the derivation of equation (6). We show that both sides of equation (6) coincide with each other. Using equation (5), [< Q2n >T,{S α ,S β } ]p is written as [< Q2n >T,{S α ,S β } ]p = 1 C
N

exp(Kp
τij =±1 (ij )

τij ) < (
i=1

α β 2n Si Si ) >T,{S α ,S β } ,

(A1)

and we abbreviate (2 cosh(Kp ))NB as C from now on. Here, we perform the following local gauge transformation:
β β α α τij → τij σi σj , Si → Si σi , Si → Si σi , (σi = ±1)

(A2)

where each σi arbitrarily takes +1 or ?1. Since < ( invariant under this transformation, we obtain [< Q2n >T,{S α ,S β } ]p = 1 C
N

N i=1

α β 2n Si Si ) >T,{S α ,S β } is

exp(Kp
τij =±1 (ij )

τij σi σj ) < (
i=1

α β 2n Si Si ) >T,{S α ,S β } , (A3)

where we note that the summation over τij σi σj = ±1 is eqivalent to that over τij = ±1. As each σi arbitarily takes +1 or ?1, therefore, we take all the summations of σi and divide by 2N , namely [< Q2n >T,{S α ,S β } ]p = 1 C 1 1 C 2N
N

exp(Kp
σi =±1 τij =±1 N (ij )

τij σi σj ) < (
i=1

α β 2n Si Si ) >T,{S α ,S β }

=

τij

Z (Kp , {τ }) α β 2n Si Si ) >T,{S α ,S β } <( N 2 =± 1 i=1

(A4)

Next, we consider the r.h.s. (we denote it as A) of equation (6). The r.h.s. of equation (6) is written as A= 1 = C 1 C exp(K
τij =±1 (ij )

τij )

Z (Kp , {τ }) Si )2n >T,{S } <( Z (K, {τ }) i=1
S i =± 1

N

exp(K
τij =±1 (ij )

Z (Kp , {τ }) τij ) Z (K, {τ })

exp(K

(ij ) τij Si Sj )(

N i=1

Si )2n

Z (K, {τ })

(A5)

Here, we perform the same local gauge transformation. Then, we obtain 1 A= C exp(K
τij =±1 (ij )

Z (Kp , {τ }) τij σi σj ) Z (K, {τ })

S i =± 1

exp(K

(ij ) τij Si Sj )(

N i=1

Si σi )2n

Z (K, {τ })

8 1 = C 1 = C = 1 C Z (Kp , {τ })
τij =±1

exp(K

(ij ) τij σi σj )

S i =± 1

exp(K

(ij ) τij Si Sj )(

N i=1

Si σi )2n

Z (K, {τ })2
Si ,σi =±1

τij

Z (Kp , {τ }) 2N =± 1

exp(K

(ij ) τij σi σj ) exp(K Z (K, {τ })2

(ij ) τij Si Sj )(

N i=1

Si σi )2n

τij

Z (Kp , {τ }) Si σi )2n >T,{S,σ} . <( N 2 =± 1 i=1

N

(A6)

Thus, from equations (A4) and (A6), we conclude that both sides of equation (6) coincide with each other. References
[1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] Bhatt R N and Young A P 1985 Phys. Rev. Lett. 54 924 Ogielski A T and Morgenstern I 1985 Phys. Rev. Lett. 54 928 Ogielski A T 1985 Phys. Rev. B 32 7384 Kawashima N and Young A P 1996 Phys. Rev. B 53 R484 Singh R R P and Chakravarty S 1986 Phys. Rev. Lett. 57 245 Bray A J and Moore M A 1984 J. Phys. C 17 L463 Swendsen R H and Wang J -S 1986 Phys. Rev. Lett. 57 2607 Bhatt R N and Young A P 1988 Phys. Rev. B 37 5606 Morgenstern I and Binder K 1980 Phys. Rev. B 22 288 Ozeki Y 1990 J. Phys. Soc. Jpn. 59 3531 Kawashima N and Rieger H 1997 Europhys. Lett. 39 85 Shirakura T and Matsubara F 1996 J. Phys. Soc. Jpn. 65 3138 Shiomi M, Matsubara F and Shirakura T private communication Shirakura T and Matsubara F 1997 Phys. Rev. Lett. 79 2887 Matsubara F, Shirakura T and Shiomi M 1998 Phys. Rev. B 58 R11821 Ozeki Y and Nishimori H 1987 J. Phys. Soc. Jpn. 56 3265 Kitatani H and Oguchi T 1992 J. Phys. Soc. Jpn. 61 1598 Nishimori H 1981 Prog. Theor. Phys. 66 1169 Kitatani H 1992 J. Phys. Soc. Jpn. 61 4049 Kitatani H and Seike T in preparation Binder K 1991 Phys. Rev. Lett. 47 119 Kikuchi M and Ito N 1993 J. Phys. Soc.Jpn. 62 3052

9
300 250 200
1.1

?

3 2
1

= )

1(

4
7 9 11 13 15 17

χ ??
150 100 50 0 0 0.1 0.2 0.3
T

7 9 11 13 15 17 19 21 23

ν= η =

0( ) 8 556 7

0.9

χ
0.8

19 21 23

0.7

0.4

0.5

0.6

0.6 -1.5

-1

-0.5

0

$&

?

&' #%"! ν

0.5

1

1.5

2

2.5

Figure 1. A plot of χSG versus T.

Figure 4. χSG .

A scaling plot for

0.88 0.86


?
7 11 15 19 23

0.84 0.82 0.8 0.78 0 0.1 0.2 0.3 0.4

0.5

T
Figure 2. A plot of gL versus T.
0.9 0.9

0.85

   =   ν = 

?
7 9 11 13 15 17 19 21 23 0.85

G H = RP ν= Q I

9
7 9 11 13 15 17 19 21 23

S
0.8

0.8

0.75

0.75

0.7 0.7 -1.5 -1 -0.5 0 0.5 0 0.5 1  ? ¨§ ν 1.5 2 2.5



CE

?

?

1 EF B DA@ ν

1.5

2

Figure 3. A scaling plot for gL .

Figure 5. A scaling plot for gL , assuming Tc = 0.

10
1
 ?

1

0.9

ν= η=

=

¨  §

B A

   



7 9 11 13 15 17

0.9

η=
0.8

ν=

=

@ C 9 D F

E I H G

7 9 11 13 15 17

0.8

χ
0.7

χ
0.7

19 21 23

19 21 23

0.6

0.6

0.5 0 0.5
??

0.5 1 1.5 2 0 0.5
57

1

1.5

2

?

?? ? ? ?¤

ν

?

32 7 8 46

ν

Figure 6. A scaling plot for χSG , assuming Tc = 0.

Figure 8. A scaling plot for χSG including the corrections to scaling with ω ′ = 0.5 and b = ?0.3, assuming Tc = 0.

1
' &

0.95

ν=

=

% (

$

0 )

0.9
1

0.85

0.8

7 9 11 13 15 17 19 21 23

0.75 0 0.5
"

1

?

 " # ! ν

1.5

2

Figure 7. A scaling plot for gL including the corrections to scaling with ω = 0.5 and a = ?0.3, assuming Tc = 0.



更多相关文章:
更多相关标签:

All rights reserved Powered by 甜梦文库 9512.net

copyright ©right 2010-2021。
甜梦文库内容来自网络,如有侵犯请联系客服。zhit325@126.com|网站地图