Universal dependence on disorder of 2D randomly diluted and random-bond \pm J Ising models

Universal dependence on disorder of 2D randomly diluted and random-bond Ising models

Martin Hasenbusch Institut für Theoretische Physik, Universität Leipzig, Postfach 100 920, D-04009 Leipzig, Germany.    Francesco Parisen Toldin Max-Planck-Institut für Metallforschung, Heisenbergstrasse 3, D-70569 Stuttgart, Germany
and Institut für Theoretische und Angewandte Physik, Universität Stuttgart, Pfaffenwaldring 57, D-70569 Stuttgart, Germany
   Andrea Pelissetto Dipartimento di Fisica dell’Università di Roma “La Sapienza” and INFN, Roma, Italy.    Ettore Vicari Dipartimento di Fisica dell’Università di Pisa and INFN, Pisa, Italy.
July 8, 2019

We consider the two-dimensional randomly site diluted Ising model and the random-bond Ising model (also called Edwards-Anderson model), and study their critical behavior at the paramagnetic-ferromagnetic transition. The critical behavior of thermodynamic quantities can be derived from a set of renormalization-group equations, in which disorder is a marginally irrelevant perturbation at the two-dimensional Ising fixed point. We discuss their solutions, focusing in particular on the universality of the logarithmic corrections arising from the presence of disorder. Then, we present a finite-size scaling analysis of high-statistics Monte Carlo simulations. The numerical results confirm the renormalization-group predictions, and in particular the universality of the logarithmic corrections to the Ising behavior due to quenched dilution.

75.10.Nr, 64.60.Fr, 75.40.Cx, 75.40.Mg

I Introduction and summary

Random Ising systems represent an interesting theoretical laboratory in which one can study general features of disordered systems. Among them, the two-dimensional (2D) random-site and random-bond Ising models have attracted much interest. In particular, the effects of quenched disorder on the critical behavior at the paramagnetic-ferromagnetic transitions, which are observed for sufficiently small disorder, have been much investigated and debated, see Refs. DD-81, ; Shalaev-84, ; Cardy-86, ; Shankar-87, ; Ludwig-87, ; LC-87, ; Mayer-89, ; Ludwig-90, ; Ziegler-90, ; WSDA-90, ; Heuer-92, ; Shalaev-94, ; KP-94, ; Kuhn-94, ; QS-94, ; TS-94, ; TS-94-2, ; MS-95, ; DPP-95, ; AQS-96, ; JS-96, ; CHMP-97, ; BFMMPR-97, ; AQS-97, ; SSLI-97, ; RAJ-98, ; SSV-98, ; AQS-99, ; MK-99, ; LSZ-01, ; Nobre-01, ; SV-01, ; MC-02, ; COPS-04, ; Queiroz-06, ; LQ-06, ; PHP-06, ; MP-07, ; KJJ-07, .

Renormalization-group (RG) and conformal field theory Shalaev-84 (); Shankar-87 (); Ludwig-90 (); Shalaev-94 () predict the marginal irrelevance of random dilution at the paramagnetic-ferromagnetic transition. Therefore, the asymptotic behavior is controlled by the standard Ising fixed point, characterized by the critical exponents and ; disorder gives only rise to logarithmic corrections. The marginality of quenched disorder coupled to the energy density, as it is the case for random dilution, is already suggested by the Harris criterion,Harris-74 () which states that the relevance or irrelevance of quenched dilution depends on the sign of the specific-heat exponent of the pure system; in the case of the 2D Ising model, the specific heat diverges only logarithmically at the transition, i.e. . The marginal irrelevance of disorder has also been supported by numerical studies of lattice models; see, e.g., Refs. WSDA-90, ; BFMMPR-97, ; AQS-97, ; RAJ-98, ; MK-99, ; SV-01, ; MC-02, ; COPS-04, ; LQ-06, ; PHP-06, ; MP-07, (see, however, Refs. KP-94, ; Kuhn-94, ; LSZ-01, for different scenarios). We recall that in three dimensions random dilution is a relevant perturbation of the pure Ising fixed point, leading to a new three-dimensional randomly diluted Ising (RDI) universality class, which is characterized by different critical exponents; see, e.g., Refs. PV-02, ; HPPV-07, .

In this paper we return to the issue of the critical behavior of 2D randomly diluted Ising systems. By using the RG results reported in Refs. Cardy-86, ; Ludwig-87, ; LC-87, ; MK-99, , we show that their critical behavior can be derived from the RG equations


where is the flow parameter (the logarithm of a length scale), , , and are the scaling fields associated with the leading operators of the three different conformal families of the 2D Ising model, i.e. the identity, energy, and spin families, and is the marginally irrelevant scaling field associated with disorder. Higher-order terms in Eqs. (1) are not necessary, because they can be reabsorbed by appropriate analytic redefinitions of the scaling fields. The appearance of the term in the first equation, where is a nonuniversal constant, is due to the resonance of the identity operator with the thermal operator, which already occurs in the pure Ising model.Wegner-76 () It is interesting to note that randomness couples only to the thermal scaling field . It would be interesting to understand if these conclusions also apply to the irrelevant operators, i.e., if the only operators that couple disorder are those that belong to the conformal family of the energy.

The analysis of the RG equations shows that random dilution gives rise to logarithmic corrections which are universal after an appropriate normalization of the scaling field associated with disorder. Additional scaling corrections due to the irrelevant operators are suppressed by power laws as in standard continuous transitions. For these reasons, we prefer to distinguish the randomly dilute Ising (RDI) critical behavior characterized by the RG equations (1) from the standard 2D Ising universality class of pure systems, although they share the same 2D Ising fixed point.

The RG equations (1) allow us to determine the scaling behavior of any thermodynamic quantity. Denoting with , , , and the reduced temperature, the magnetic external field, the disorder parameter, and the lattice size, respectively, the free energy satisfies the scaling equation


(we consider models defined on square lattices with periodic boundary conditions), where , , , and are the solutions of the RG equations. From Eq. (2) one can derive the scaling behavior of the relevant thermodynamic quantities and determine the logarithmic corrections due to the quenched disorder. At the critical point , we obtain the asymptotic behaviorsMK-99 ()


for the specific heat, and


for the magnetic susceptibility , any RG invariant quantity , such as the quartic Binder cumulant and the ratio , and its derivative with respect to the temperature. Here indicates the Ising fixed-point value and is the solution of Eq. (1) with . For , behaves as


where is a length scale. The functions , , and are normalized such that and are universal once the scaling field is appropriately normalized. In the above-reported equations we have neglected scaling corrections which are suppressed by power laws. They are due to the analytic dependence of the scaling fields on the Hamiltonian parameters, to the background (i.e., the contribution of the identity operator in the RG language), and to the irrelevant operators. CCCPV-00 (); CHPV-02 () In particular, we expect scaling corrections associated with the leading irrelevant operator appearing in the pure Ising model (the corresponding exponent is ).

Moreover, in this paper we compare the theoretical predictions with a finite-size scaling (FSS) analysis of numerical Monte Carlo (MC) results for the randomly site-diluted Ising model and for the random-bond Ising model, also known as Edwards-Anderson model. Our main results can be summarized as follows. Our FSS analyses provide a robust evidence that the paramagnetic-ferromagnetic transitions in these models present the same RDI critical behavior. Note that this implies that frustration in the random-bond Ising model is irrelevant at the paramagnetic-ferromagnetic transition line. The FSS behaviors are in agreement with the predictions of the RG equations (1). The asymptotic critical behavior appears to be controlled by the Ising fixed point. The logarithmic corrections and their universal behavior are consistent with the theoretical results obtained from Eqs. (1), cf. Eqs. (2), (3), (4), (6), (7).

The paper is organized as follows. In Sec. II we define the randomly site-diluted model and the Ising model, we briefly discuss their phase diagrams, and define the quantities that we consider in the paper. In Sec. III we discuss the RG flow at a 2D Ising fixed point in the presence of a marginally irrelevant operator, and the implications for the infinite-volume and finite-size critical behavior. In particular, we focus on the universal features of the logarithmic corrections due to disorder. Finally, in Sec. IV we present our FSS analysis of high-statistics MC results for the randomly site-diluted and random-bond Ising models.

Ii Randomly site-diluted and random-bond Ising models

ii.1 The models and their phase diagrams

The randomly site-diluted Ising model (RSIM) is defined by the Hamiltonian


where the sum is extended over all pairs of nearest-neighbor sites of a square lattice, are Ising spin variables, and are uncorrelated quenched random variables, which are equal to one with probability (the spin concentration) and zero with probability (the impurity concentration). The RSIM is expected to undergo a continuous transition for any , whereNZ-00 () corresponds to the site-percolation point of the impurities; moreover, for , see, e.g., Ref. Cardy-book, . For the ferromagnetic phase is absent. Thus, the paramagnetic-ferromagnetic transition line starts from the pure Ising point , where is the critical temperature of the 2D Ising model, and ends at . Along this line the critical behavior is expected to be universal, i.e. independent of dilution, and to be characterized by the RG equations (1). As we shall see, this is supported by the analysis of our MC results.

The random-bond Ising model, also known as Edwards-Anderson model,EA-75 () is defined by the lattice Hamiltonian


where , the sum is over all pairs of nearest-neighbor sites of a square lattice, and the exchange interactions are uncorrelated quenched random variables, taking values with probability distribution


In the following we set without loss of generality. For we recover the standard Ising model, while for we obtain the bimodal Ising spin-glass model. The Ising model is a simplified model EA-75 () for disordered spin systems showing glassy behavior in some region of their phase diagram. Its phase diagram in two dimensions is sketched in Fig. 1 (it is symmetric for ). For sufficiently small values of the probability of the antiferromagnetic bonds, the model presents a paramagnetic phase and a ferromagnetic phase. The paramagnetic-ferromagnetic transition line starts from the Ising point and extends up to the multicritical Nishimori point (MNP) , located along the so-called Nishimori line ( line) defined by Nishimori-81 (); GHDB-85 (); LH-88 (); LH-89 () with HPPV-08 () and . The critical behavior is expected to be in the same universality class as that of the transition in the RSIM. As we shall see, our FSS analysis strongly supports this scenario. A detailed discussion of the phase diagram can be found in Ref. HPPV-08, and references therein.

Figure 1: Phase diagram of the square-lattice random-bond Ising (Edwards-Anderson) model in the - plane.

ii.2 Observables

In our FSS analyses we consider models defined on a square lattice with periodic boundary conditions. The two-point correlation function is defined as


where the angular and square brackets indicate the thermal average and the quenched average over disorder, i.e. over in the case of RSIM and over in the case of the Ising model. We define the magnetic susceptibility and the correlation length ,


where , , and is the Fourier transform of . We also consider quantities that are invariant under RG transformations in the critical limit. Beside the ratio


we consider the quartic cumulants and defined by




The above RG invariant quantities , , and are also called phenomenological couplings. In the critical () 2D pure Ising model, they converge for large to the universal valuesSS-00 ()


Finally, we consider the derivatives


which can be computed by measuring appropriate expectation values at fixed and .

Iii Renormalization-group flow and finite-size scaling

iii.1 Ising RG flow in the presence of a marginally irrelevant scaling field associated with disorder

In this section we discuss the RG flow close to the 2D Ising fixed point in the presence of a marginally irrelevant scaling field associated with disorder.

Let us consider a system with a marginal scaling field and with a set of scaling fields , , with RG dimensions . Close to the fixed point , the RG equations have the form


If there are no degeneracies ( for all ) and no resonancies (i.e., there is no combination with integer coefficients of the RG dimensions that vanishes), one can redefine the scaling fields in such a way to simplify the RG equations. We define


With a proper choice of the coefficients , , , we can simplify the RG equations, obtaining the simple canonical form:


By normalizing appropriately the scaling field , we can also set . In the case we are considering, is marginally irrelevant so that (we assume that is defined such that ). We can thus set . Once this choice has been made, and all coefficients are universal.

The simple form we have derived above does not strictly apply to the RSIM. Indeed, in the Ising model the RG operators belong to three different conformal families and within each family the RG dimensions differ by integers (see Ref. CHPV-02, for a discussion of the irrelevant operators in the pure Ising model). Thus, in the present case there are both degeneracies and resonancies. If we limit our considerations to the relevant scaling fields, we must only consider the resonance between the identity operator and the thermal operator, which is responsible for the logarithmic divergence of the specific heat in the pure Ising model.Wegner-76 () By a proper redefinition of the nonlinear scaling fields, one can show that in this case the RG equations (20) for the relevant scaling fields can be written as:


where the couplings to the irrelevant scaling fields due to the additional resonancies have been neglected. The scaling field is associated with the identity operator. The additional term which appears in Eq. (25) is due to the resonance with the thermal operator, as in the pure Ising model.Wegner-76 () The scaling fields and are the relevant scaling fields associated with the reduced temperature and the external field , respectively; and are the corresponding RG dimensions. Finally, is the marginally irrelevant operator associated with randomness. The coefficients , , , and are universal, being independent of the normalization of the scaling fields. By using conformal field theory, , , and have been computed: Cardy-86 (); Ludwig-87 (); LC-87 (); DPP-95 (); MK-99 ()


Let us now integrate the RG equations. Since , Eq. (28) has two fixed points with : one is and is stable; the second one is and is unstable. Thus, the basin of attraction of the Ising FP corresponds to ; for , flows to infinity. It is important to note that Eq. (28) is only valid within the basin of attraction of the stable fixed point . The redefinitions of the scaling fields that we have used to obtain the simple canonical form (28) cannot be extended outside the basin of attraction since they are expected to become singular at the unstable fixed point. The presence of an unstable fixed point indicates that the behavior for large values of the disorder is not controlled by the Ising fixed point. The RG flow could be attracted by a new fixed point—thus, for large values of the disorder the transition would be continuous and in a new universality class—or could go to infinity, indicating the absence of a continuous transition. A similar phenomenon was conjectured in three dimensionsCPPV-04 () on the basis of a perturbative field-theoretical analysis of the RG flow.

If , the function is implicitly given by (we do not replace with its theoretical value , in order to obtain general expressions that can be tested numerically)


The solution can be simplified if we introduce


which satisfies the implicit equation


This equation can be inverted to give


The function cannot be computed analytically. However, it is easy to determine it in the large- limit. We obtain


where . Since the equation for gives


where . In order to determine , we rewrite the corresponding equation as


which gives (, )


where . For large the function behaves as


Let us finally consider the identity operator. If the solution can be obtained as in the case of . Thus, we write


where is an unknown function of , such that . Substituting in the equation for and using the result for , we obtain


and therefore


The behavior of for depends on the value of . Since the integral appearing in diverges as for , as for , and is finite for , we obtain


The RG equations do not fix completely the normalization of the scaling fields. First, one can redefine , , and by a multiplicative constant;footnote-norm () such a redefinition is not possible for , since a multiplicative constant would break the condition . Beside these trivial redefinitions there is also a nonlinear set of transformations that leave the equations invariant. Given a constant , we define as the solution of the equation


Then, for any we have


Note that the transformation is analytic in a neighborhood of . If is defined as in Eq. (31), we obtain


Analogously, if we define


then satisfies the same equation of with replacing , as it can be seen from Eq. (37). A similar redefinition can be made for . This invariance implies that, beside fixing the normalizations of , , and , we must also appropriately fix . In practical terms, is completely arbitrary and must be fixed in order to define unambigously. Finally, note that there are no analytic redefinitions of that map Eq. (28) in an identical equation with , proving the universality of .

Neglecting scaling corrections that are suppressed by power laws, we write the free energy in the scaling formWegner-76 ()


for any . Note that the whole dependence on , , and is encoded in the constants , , , , and . Of course, and vanish at the critical point, while vanishes for . The independence of Eq. (47) on allows us to simplify the general expression for the free energy. We choose such that and thus


where . Substituting these expressions in Eq. (47), we obtain the general dependence of the free energy on and .

In order to determine , we consider the specific heat . The leading singular behavior is due to the temperature dependence of the scaling field . Using Eq. (42) we obtain


The asymptotic behavior of the specific heat of 2D randomly diluted Ising systems has been determined in Refs. DD-81, ; Shalaev-84, ; Shankar-87, , obtaining


Comparing with Eq. (49), we obtain . In this case we have


It is interesting to note that, since , randomness couples only to the thermal scaling field . This result appears quite natural from the point of view of the Landau-Ginzburg-Wilson approach to critical phenomena. Indeed, in field theory randomly diluted models are obtained by coupling disorder to the energy operator: Lubensky-75 (); GL-76 ()


where , and is a spatially uncorrelated random field with Gaussian distribution. The 2D RDI critical behavior has been also investigated by using this field-theoretical approach and the so-called replica trick. The analysis of the corresponding five-loop perturbative expansionsCOPS-04 (); Mayer-89 () gives results which are substantially consistent with the marginal irrelevance of disorder.

iii.2 Finite-size scaling

Let us now discuss the implications of the above RG analysis for the FSS of thermodynamic quantities at the critical point. We start from the scaling behavior of the free energy


where the contributions of the irrelevant scaling fields have been neglected. By choosing , we can write


If we set in Eq. (32), for we have


The free energy can then be written as


The constants , , , and depend on , , and . Moreover, and close to the critical point. The terms proportional to , , and are due to the identity operator, whose dependence on is given in Eq. (51). Eq. (56) is valid up to contributions of the irrelevant operators, which are expected to scale as inverse powers of .

From Eq. (56) we can compute zero-momentum quantities that involve disorder averages of a single thermal average. For instance, for the specific heat at and we obtain


For the susceptibility at we obtain


where, as before, we neglect power-law scaling corrections. A similar scaling equation holds for :


The determination of the scaling behavior of and requires an extension of the scaling Ansatz (56). A detailed discussion is presented in Sec. 3.1 of Ref. HPPV-07, . It shows that both quantities behave as , apart from scaling corrections. Thus, if or , we also have


Derivatives of the phenomenological couplings have a simple behavior as well, the leading term being of the form


At the critical point we can set , so that we can write the scaling behaviors


The functions and are universal once an appropriate normalization is chosen for , which is independent of the model. For this purpose, let us consider a phenomenological coupling . For we can expand


Now we normalize by requiring . It is easy to prove that this is a correct normalization condition. Indeed, imagine that has been normalized arbitrarily so that . Then, redefine by using Eq. (45). By properly choosing , it is easy to see that one can set . This condition fixes uniquely the scale .

Note that, in the pure Ising model, we have , so that we expect at the critical point


for . It is natural to invert this relation to express in terms of . Then, we obtain the scaling forms


where , , are universal scaling functions that are normalized such that , and have a regular expansion in powers of . Note that these scaling equations are much simpler than those in terms of , since they are independent of the scale and of the normalization of .

iii.3 Finite-size scaling at a fixed phenomenological coupling

Instead of computing the various quantities at fixed Hamiltonian parameters, one may study FSS keeping a phenomenological coupling fixed at a given value , as proposed in Ref. Hasenbusch-99, and also discussed in Refs. HPV-05, ; HPPV-07, . This means that, for each , one considers such that


All interesting thermodynamic quantities are then computed at . The pseudocritical temperature converges to as .

In the next section we report a FSS analysis of MC simulations keeping the phenomenological coupling fixed. The value can be specified at will as long as it is between the corresponding high-temperature and low-temperature values. Since we wish to check the hypothesis that the asymptotic critical behavior is governed by the Ising fixed point, we choose , where is the universal value of at the critical point in the 2D Ising universality class SS-00 () for square lattices with periodic boundary conditions. Note, however, that this choice does not bias our analysis in favor of the Ising nature of the transition. By fixing to the critical Ising value, we can perform the following consistency check: if the transition belongs to the Ising universality class, then any other RG-invariant quantity must converge to its critical-point value in the Ising model.

In the plane, the line is obtained by solving the equation


It gives a relation


where has a regular expansion in powers of . Moreover, since we have chosen , we have . Substituting relation (70) in the above-reported scaling equations for the susceptibility , the phenomenological couplings , and their derivative, we obtain at fixed


The scaling functions are universal, have a regular expansion in powers of , and are normalized such that , . The additional corrections due to the irrelevant operators decay as powers of .

The large- behavior of follows from Eq. (70). Since , we obtain


where is computed at the critical point .

We finally mention that Eqs. (65), (66), (67) hold at fixed as well. The corresponding universal scaling functions depend on the values of at fixed, i.e., , (we denote them by , , and , respectively) and have a regular expansion in powers of .

Iv Finite-size scaling analysis of Monte Carlo simulations

iv.1 Monte Carlo simulations

RSIM, 8 5361 1.16476(3) 0.05083(3) 36.1853(9) 6.5911(11) 2.7285(5)
16 2560 1.16463(4) 0.04170(4) 122.367(4) 12.608(5) 3.4282(12)
32 1280 1.16507(4) 0.03618(4) 412.573(14) 24.029(9) 4.0283(14)
64 640 1.16563(6) 0.03237(6) 1389.57(8) 45.91(3) 4.550(3)
128 640 1.16597(6) 0.02947(6) 4677.0(3) 87.84(7) 5.014(3)
256 653 1.16619(5) 0.02704(5) 15741.3(7) 168.50(9) 5.431(3)
512 633 1.16656(4) 0.02522(5) 52962(2) 324.18(19) 5.815(3)
RSIM, 8 640 1.14557(10) 0.09561(10) 25.841(3) 1.2536(13) 0.2155(3)
16 2176 1.15206(6) 0.07526(6) 86.941(5) 2.6583(11) 0.30474(14)
32 1280 1.15682(6) 0.06297(7) 293.888(15) 4.967(2) 0.35203(12)
64 658 1.15996(7) 0.05491(8) 993.26(6) 9.140(4) 0.38283(10)
128 843 1.16185(6) 0.04871(7) 3351.92(14) 16.903(6) 0.40516(7)
256 1288 1.16313(4) 0.04368(5) 11299.1(3) 31.501(9) 0.42262(5)
Is, 8 3200 1.16399(2) 0.04026(3) 40.9962(8) 6.0887(15) 2.6027(7)
16 3200 1.16405(3) 0.04023(3) 139.318(5) 11.163(2) 3.1527(8)
32 3200 1.16439(2) 0.03847(3) 470.511(11) 20.924(3) 3.6299(5)
64 812 1.16482(4) 0.03592(5) 1585.27(7) 39.570(10) 4.0371(7)
128 658 1.16527(4) 0.03331(6) 5337.0(3) 75.12(2) 4.3865(5)
Table 1: MC data at fixed . For each model and lattice size , we report the number of samples , the quartic cumulants and , the magnetic susceptibility , the derivative , and the specific heat . If the asymptotic behavior is controlled by the Ising fixed point, for we should have and .

We perform high-statistics MC simulations of the RSIM at , and of the Ising model at . We consider square lattices of linear size with periodic boundary conditions. In the MC simulations of the RSIM we use a mixture of Metropolis and Wolff clusterWolff-89 () updates as we did in the three-dimensional case reported in Ref. HPPV-07, . In the case of the Ising model, the Wolff cluster update is expected to be slowHPPV-07-pmj () so that we only use Metropolis updates with multispin coding.

Instead of computing the different quantities at fixed Hamiltonian parameters, we compute them at fixed . This means that, given a MC sample generated at , we determine the value such that . All interesting observables are then computed at . The pseudocritical temperature converges to as . This method has the advantage that it does not require a precise knowledge of the critical value (an estimate is only needed to fix that should be close to ). Moreover, for some observables the statistical errors at fixed are smaller than those at fixed (close to ).HPV-05 (); HPPV-07 () In order to compute any quantity at , we determine its Taylor expansion around , as we did in our previous work.HPPV-07-pmj () Particular care has been taken to avoid any bias due to the finite number of iterations for each sample: we use the method described in Ref. HPPV-07, and extended to correlated data in Ref. HPPV-07-pmj, .

The results at fixed are reported in Table 1. For each model and lattice size , we report the number of samples, the MC estimates of the quartic cumulants and at fixed (we denote them with and , respectively), the magnetic susceptibility , the derivative , and the specific heat .

iv.2 Results

iv.2.1 Approach to the 2D Ising fixed-point values

Figure 2: The phenomenological coupling vs . The lines show the results of fits to Eq. (87). For the RSIM at and the Ising model we fit all data, while for the RSIM at , we use data satisfying . Note that the asymptotic slope as of the resulting curves is identical in the three cases, confirming the universality of , defined by , see Sec. IV.2.3.

Since we perform our FSS analysis keeping fixed, if the critical behavior is controlled by the Ising fixed point, in the large- limit we should have


where SS-00 () is the universal large- limit of the quartic (Binder) cumulant at the critical point in the 2D Ising model. Since disorder is expected to be marginally irrelevant, see Sec. III.1, the approach of and to their large- Ising limit is expected to be logarithmic.

The MC data of and , reported in Table 1, clearly approach the Ising values (75). In the case of