Recursive-search method for ferromagnetic ising systems: combination with a finite-size scaling approach

Recursive-Search Method for Ferromagnetic Ising Systems: Combination with a Finite-Size Scaling Approach

P. C. da Silva, U. L. Fulco^*, F. D. Nobre, L. R. da Silva, and L. S. Lucena

Departamento de Física Teórica e Experimental

Universidade Federal do Rio Grande do Norte

Campus Universitário, Caixa Postal 1641, 59072-970, Natal, RN, Brazil

Received on 19 October, 2001; Revised version received on 2 February, 2002

A method for obtaining critical properties of physical systems is presented. Based on a recursive relation involving a physical parameter of the system, it drives the system spontaneously to the critical point, providing an efficient way to estimate critical properties. The method is illustrated for several ferromagnetic Ising systems on well-known Bravais lattices. A finite-size scaling approach is performed, by applying the method on lattices of different sizes. The efficiency of the method is confirmed by evaluating critical temperatures, as well as critical exponents, that turn up to be in good agreement with those available in the literature, with a relatively small computational effort.

I Introduction

Although statistical mechanics represents one of the most successful physical theories nowadays, only a few simple theoretical models have been solved exactly within such a framework. The evaluation of thermal averages, that are, in principle, to be performed over the whole phase space, becomes a hard task in most of the cases. As a consequence, many approximation methods have been proposed in order to deal with complicated systems, characterized by many-interacting constituents. Due to the recent improvements in computer technology, the computer simulations [1, 2, 3, 4] became nowadays one of the most important tools for studying physical systems. Among many different types of computer simulations, one may single out the Monte Carlo (MC) simulations, that are probably the most commonly employed of all numerical simulations. The MC method consists in performing the usual averages of statistical mechanics over a restricted part of phase space, i.e., only those configurations which contribute significantly to the thermal averages are considered; such a procedure reduces a lot of computing time, in such a way that one may study large - but finite - physical systems. In a standard MC simulation, each dynamical variable (which may be defined on sites of regular lattices) is visited either at random or in well-defined sequences, to be afterwards updated according to certain dynamical rules; depending on the size of the sistem, a MC simulation may require a large computational effort. The main drawback of any numerical simulation is that one is restricted to work with finite-size systems, and sometimes, the finite-size effects may disguise important physical results. A typical case is when one is working with systems at criticality, where in order to obtain reliable results, one needs to extrapolate the data obtained for finite sizes to the infinite-size limit (the so-called thermodynamic limit). The most commonly used extrapolation technique for systems exhibiting critical behavior is the finite-size scaling (FSS) approach [1, 2]. Through the FSS method one is able to extract critical properties of a physical system, e.g., critical exponents, from finite-size data; however, a good estimate of critical properties requires a precise determination of the critical point. Although some FSS methods are able to produce critical exponents, as well as the location of the critical point, the application of such a procedure becomes much simpler if one knows a priori the location of the critical point.

An important step in the theory of critical phenomena occurred through the concept of self-organized criticality [5], according to which certain dynamical systems evolve spontaneously towards the critical state, i.e., the critical state is an attractor of the dynamics. Recently, the concept of self-organized criticality has been applied to the determination of critical properties in polymers [6], percolation [7], and magnetic systems [8, 9, 10]. The method uses an algorithm based on a recursive relation,

involving two dimensionless variables (X_n,Y_n), associated with parameters of a given physical system. The variables (X_n, Y_n) change at each iteration step n, in such a way that after a sufficient number of steps, X_n will converge to a stationary value X^*, compatible with the stationary value Y^*º Y(X^*), assumed by Y_n. The desired stationary state (X^*,Y^*) may be previously selected by an appropriate choice of the quantity Y^*; the rate of convergence to the stationary state is controlled by the parameter a. As an example, for a ferrromagnetic system, such quantities may be related to the temperature and magnetization [8, 9] (or to the temperature and inverse of magnetic susceptibility [10]), respectively; instead of controlling the temperature T, one may reach criticality by taking the magnetization m ® 0 (or the inverse of susceptibility 1/ c® 0), which is equivalent to approaching the critical temperature, T ® T_c. It is important to mention that the idea of keeping a physical parameter of the system, associated with the critical state, close to a small positive value (pushing the system automatically to the vicinity of the critical point), has been guessed by Sornette et al. [11], although it was not operationally worked out in terms of a recursive relation.

In the present work we illustrate the algorithm based on recursive relation (1) by applying it to the ferromagnetic Ising model, defined on some well-known Bravais lattices, namely, the square, triangular, honeycomb and cubic lattices. In each case, the recursive approach is applied for several finite sizes, and a FSS approach is performed. In spite of a small computational effort, the critical temperatures and critical exponents obtained are in good agreement with those available in the literature. In the next section we describe the numerical procedure. In section 3 we present and discuss our results.

II The Numerical Procedure

Let us consider the nearest-neighbor interaction ferromagnetic Ising model, defined through the Hamiltonian,

with J > 0 , and S_i = ±1. The model will be studied on some well-known Bravais lattices, namely, the square, triangular, honeycomb, and cubic ones. In each case, several linear sizes L will be considered; the total number of spins is N = L², for the two-dimensional lattices, and N = L³, for the cubic lattice.

For the present problem, the parameter X of Eq. (1) will be related to the temperature, X º K = J/(k_BT). The recursive method is carried out through several steps, as described below [9, 10].

(a) One previously defines a (which controls the rate of convergence) and Y^* (which defines the desired stationary state to be accessed).

(b) At the first iteration, one chooses the initial value K₀. This will define the range of parameters to be investigated [(K_n, Y_n) will vary from (K₀, Y₀) up to (K^*, Y^*)].

(c) A particular initial configuration is assigned to the spin variables and the system is let to evolve dynamically according to a given MC prescription. Herein, we have used a single spin-flip updating, following the Glauber dynamics [12], according to which,

and

In the equations above, z_i(t) is a uniform random number in the interval [0, 1] and p_i(t) is the probability

where

is the local field acting on site i, at time t , and for the first iteration, K º K₀ = J/(k_BT₀).

(d) After equilibration is attained (t₀ MC steps) one may calculate thermodynamic averages (associated with the particular choice of K₀) over t₁ MC steps. To improve the statistics, this procedure is repeated over N_s samples, i.e., N_s distinct sequences of random numbers. The average value Y₀ is computed, and from Eq. (1) one obtains K₁.

(e) Steps (c) and (d) are performed for parameter K₁, and so on, in such a way that one gets, iteratively, (K₀, Y₀) ® (K₁, Y₁) ® (K₂, Y₂) ¼ .

(f) The process converges when Y_n and K_n present small oscillations around the values Y^* [defined in step (a)] and K^* (the desired stationary value of the parameter K), respectively.

(g) After the stationary regime is attained, one may consider a number Dn of oscillations around (K^*,Y^*), in order to get a statistics for the stationary temperature T^*.

It is important to mention that Eq. (1) has two parameters to be adjusted a priori, namely, a and Y^*, in order to get a proper convergence of the recursive approach. The parameter a must be set to a small value; for an appropriate choice of Y^*, we found an optimal value for a, below which K^* does not change within the error bars: a = 10^-2. The choice of Y^* is somewhat more subtle; it must be chosen according to the desired stationary state. If one chooses Y as the order parameter of the system, for achieving a convergence towards criticality, Y^* should be set to a small value (typically, Y^* = 10^-2), whereas for a convergence to a low-temperature state, Y^* must be set close to its maximum value. In each problem a few attempts are required in order to find the appropriate Y^*, in such a way as to provide convergence to the desired stationary state.

In the models investigated herein, one may choose in Eq. (1), Y_nº m_n [8, 9] as the dimensionless magnetization per spin (m=N^-1å_iáS_iñ, where á ñ stands for a thermodynamic average). However, for physical systems exhibiting strong finite-size effects, e.g., disordered magnets, the magnetization may present pronounced tails close to the critical temperature, and such a choice may lead to a large error in the location of the critical temperature. For cases like that, one may choose Y_n as the inverse of a quantity which diverges at the critical temperature [10]; an appropriate quantity may be the magnetic susceptibility,

In the present work we estimated the stationary temperatures T^* by using Eq. (1) with X_nº K_n = J/(k_BT_n) and Y_nº 1/(J c_n) as the dimensionless parameter associated with the magnetic susceptibility. In the next section we present and discuss our results.

III Results and Discussion

For the results that follow, we have simulated the model defined through Eq. (2) on two-dimensional lattices (square, triangular, and honeycomb lattices) of linear dimensions L = 20, 40, 50, and 60, and on cubic lattices of linear sizes L = 10, 12, 14, and 16. At each iteration n, our simulations always started with a completely ordered configuration (all spins up) and a number of runs, corresponding to a time t₀ = N MC steps, were discarded before calculating averages. After that, we have computed thermodynamic averages over t₁ = N MC steps. All simulations were repeated over N_s = 200 samples, i.e., different sequences of random numbers.

First of all, in order to find the stationary temperature , associated with each linear size L of a given lattice, it is important to determine with a good accuracy the stationary parameter Y^*_L=1/(). In Fig. 1 we present the evolution of the temperature with the iteration step n, for different choices of , for a square lattice of linear size L=20; the same is done for the magnetic susceptibility in Fig. 2. From Figs. 1 and 2, one sees clearly that there is an optimal value of . For choices above the optimal value, Eq. (1) will not converge to the stationary values, in such a way that the temperature will diverge [e.g., see the curve for the choice 1/()=0.050 in Fig. 1], whereas the magnetic susceptibility will go down after a finite number of iteration steps [e.g., see the curves for the choices 1/()=0.050 (Fig. 2(a)) and 1/()=0.052 (Fig. 2(b))]. For choices below the optimal value, Eq. (1) will converge to stationary values that are not those associated with criticality; indeed, one finds a convergence to temperature values below the critical temperature [e.g., see the curves for the choices 1/()=0.056 and 1/()=0.058 in Fig. 1], whereas the susceptiblity will converge to values that are lower than the one at criticality. The optimal value of should be the highest choice for which Eq. (1) converges to the corresponding stationary values. For the cases considered in Figs. 1 and 2, our closest estimate to the optimal value is 1/()=0.054, which leads to a stable convergence up to iteration step n=10⁴, as exhibited in Fig. 2(c). Such a choice is associated with the stationary temperature k_B/J=2.380 ±0.005, representing a discrepancy of about 5 % with respect to the well-known square-lattice exact critical temperature (k_BT_c/J = 2.269185...).

By carrying such a procedure for different linear sizes of the lattices considered herein, one may find stationary temperatures for each lattice size. In Fig. 3 we exhibit the evolution of the temperature T_L with the iteration step n for several sizes of a honeycomb lattice; each curve of Fig. 3 is produced with its corresponding optimal choice . One sees that the stationary temperature approaches the exact critical temperature for increasing lattice sizes. Therefore, one may extrapolate the finite-size stationary temperatures to the limit L ® ¥, in order to find the stationary temperatures in the thermodynamic limit, T^*. The results obtained in Fig. 3, for the honeycomb lattice, are extrapolated towards the thermodynamic limit in Fig. 4. Our extrapolated stationary temperatures are presented in Table 1 for the lattices considered herein. One sees that two of our extrapolations (honeycomb and cubic lattices) agree, within the error bars, with the values available in the literature, whereas for the other two cases (square and triangular lattices), we find a small discrepancy (less than 1 % ) with respect to the well-known exact values. Considering the modest lattice sizes investigated, the accuracy of the temperatures estimated in Table 1 show the potential of the recursive approach presented herein.

Let us now consider the calculation of critical exponents through the present method; there are two straightfoward ways to carry out such a computation, as we mention below.

(i) For a given system size, one may evaluate the optimal value for and consequently the corresponding stationary temperature , as described above. Then, one may compute the critical exponent associated with a given thermodynamic property by investigating how such a thermodynamic property behaves when one approaches the stationary temperature . By repeating such a procedure for different system sizes, one may compute critical exponents for different sizes and carry out an extrapolation to the limit L ® ¥. This approach has been applied successfully for pure [9] and site-diluted [10] ferromagnetic Ising models, showing that the recursive method defined by Eq. (1) approaches criticality through the correct thermodynamic path.

(ii) By knowing the optimal value for , one deliberatily chooses a above such an optimal value. Obviously, there will be no convergence of the recursion relation in Eq. (1), and if one starts the iteration from a sufficiently low temperature, the desired temperature range - including the critical region - may be completely explored. In this case, the recursive relation of Eq. (1) is used in marginal way, determining the sequence of temperatures at which the susceptibility will be computed. By repeating such a procedure for different lattice sizes, one may obtain the critical exponents, as well as the stationary temperatures, in the thermodynamic limit, by implementing a standard FSS approach [1,2].

In the present work we shall consider procedure (i) for the location of the critical point, i.e., for computing stationary temperatures, and procedure (ii) for computing both stationary temperatures and critical exponents. In Figs. 5 and 6 we exhibit the magnetization and magnetic susceptibility curves, respectively, obtained by scheme (ii), for several linear sizes of the square-lattice model (similar plots hold for the other lattices). One may now use the standard scaling functions,

in order to obtain, from our finite-size data, the stationary temperature T^*, as well as the critical exponents b, g, and n in the thermodynamic limit. In Figs. 7 and 8 we exhibit the data colapse of the magnetization and susceptibility curves shown in Figs. 5 and 6, respectively, for the square-lattice model (similar data colapses apply for the other lattices). The results of our FSSs are exhibited in Table 2, and compared with the well-known values available in the literature. All the estimated stationary temperatures and critical exponents agree, within the error bars, with the values available in the literature. In two cases (square and triangular lattices), the stationary temperatures estimated through a simple extrapolation to the limit L ® ¥ (cf. Table 1), and those obtained by means of the FSS approach (cf. Table 2), present a small discrepancy (including the respective error-bar ranges). We believe that such discrepancies are consequences of the small lattice sizes investigated.

To conclude, we have investigated the ferromagnetic Ising model defined on several well-known Bravais lattices, by combining a recursive method that drives the system spontaneously towards criticality, with a standard FSS approach. The recursive method allows for the detection of the critical temperatures spontaneously, through a convergence towards a fixed point of a recursion relation involving a pair of dimensionless variables (X,Y), constructed in such a way that X is associated with the temperature, whereas Y may be related, in principle, to any physical parameter displaying a nontrivial behavior at criticality. In spite of the small computational effort involved, the effectiveness of the method has been confirmed by the evaluation of critical temperatures and critical exponents that are in good agreement with those available in the literature. The agreement of our critical-exponent estimates with well-known values indicates that the recursive method yields a convergence towards criticality following the correct thermodynamic path. Certainly, the accuracy of present results may be significantly improved by using larger lattice sizes. However, we have illustrated that the power of the method - which has already proven its efficiency in finding critical properties for branched polymers [6], percolation [7], pure [8, 9] and site-diluted [10] Ising ferromagnets - may be substantially enhanced if it is combined with a FSS approach.

Acknowledgments

The partial financial supports from CNPq, CAPES, Pronex/MCT, and ANP/CTPETRO/CENPES/Petrobras (Brazilian agencies) are acknowledged.

• [1] D. P. Landau and K. Binder, A Guide to Monte Carlo Simulations in Statistical Physics (Cambridge University Press, Cambridge, 2000).
• [2] M. E. J. Newman and G. T. Barkema, Monte Carlo Methods in Statistical Physics (Oxford University Press, Oxford, 1999).
• [3] H. Gould and J. Tobochnik, An Introduction to Computer Simulation Methods, Second Edition, (Addison-Wesley Publishing Company, Reading, Massachusetts, 1996).
• [4] K. Binder and D. W. Heermann, Monte Carlo Simulation in Statistical Physics (Springer-Verlag, Berlin, Heidelberg, 1988).
• [5] P. Bak, C. Tang, and K. Wiesenfeld, Phys. Rev. A 38, 364 (1988).
• [6] J. S. Andrade Jr., L. S. Lucena, A. M. Alencar, and J. E. Freitas, Physica A 238, 163 (1997).
• [7] A. M. Alencar, J. S. Andrade Jr., and L. S. Lucena, Phys. Rev. E 56, R2379 (1997).
• [8] U. L. Fulco, L. S. Lucena, and G. M. Viswanathan, Physica A 264, 171 (1999).
• [9] U. L. Fulco, F. D. Nobre, L. R. da Silva, L. S. Lucena, and G. M. Viswanathan, Physica A 284, 223 (2000).
• [10] U. L. Fulco, F. D. Nobre, L. R. da Silva, and L. S. Lucena, Physica A 297, 131 (2001).
• [11] D. Sornette, A. Johansen, and I. Dornic, J. Phys. I France 5, 325 (1995).
• [12] R. J. Glauber, J. Math. Phys. 4, 294 (1963).
• [13] R. J. Baxter, Exactly Solved Models in Statistical Mechanics (Academic Press, London, 1982).
• [14] H. G. Ballesteros, L. A. Fernández, V. Martín-Mayor, A. Muñoz Sudupe, G. Parisi, and J. J. Ruiz-Lorenzo, J. Phys. A 32,1 (1999).

^*

Present address: Departamento de Física, Universidade Federal de Alagoas, 57072-970, Maceió, AL, Brazil

Publication Dates

• Publication in this collection
  20 Nov 2002
• Date of issue
  June 2002

History

• Reviewed
  02 Feb 2002
• Received
  19 Oct 2001