Systematic enumeration of configuration classes for entropic sampling of Ising models

Systematic enumeration of configuration classes for entropic sampling of Ising models


Departamento de Física, Instituto de Ciências Exatas, and National Institute of Science and Technology for Complex Systems,

Universidade Federal de Minas Gerais

C.P. 702, 30123-970, Belo Horizonte, MG, Brazil


We describe a systematic method for complete enumeration of configuration classes (CCs) of the spin-1/2 Ising model in the energy-magnetization plane. This technique is applied to the antiferromagnetic Ising model in an external magnetic field on the square lattice, which is simulated using the tomographic entropic sampling algorithm. We estimate the number of configurations, , and related microcanonical averages, for all allowed energies and magnetizations for to 30, with . With prior knowledge of the CCs, we can be sure that all allowed classes are sampled during the simulation. Complete enumeration of CCs also enables us to use the final estimate of to obtain good initial estimates, , for successive system sizes () through a two-dimensional interpolation. Using these results we calculate canonical averages of the thermodynamic quantities of interest as continuous functions of temperature and external field . In addition, we determine the critical line in the - plane using finite-size scaling analysis, and compare these results with several approximate theoretical expressions.

Keywords: Ising model; antiferromagnet; configuration classes; Monte Carlo simulation; phase diagram.

1 Introduction

The most important task of entropic sampling algorithms [1][6] is to visit the full configuration space (CS) to obtain good estimates of the number of configurations, , as functions of the energy E and other quantities of interest. In studies of the Ising model in an external field, for example, we require with the magnetization; each allowed pair defines a class of configurations (CC). Only if we know beforehand the possible values of for a given system size, can we be sure that all CCs are sampled during the simulation.

Figure 1 shows the CCs in the plane for the spin-1/2 Ising model with nearest-neighbor (NN) interactions, on a square lattice of sites with periodic boundaries. [We use to denote the number of NN pairs of spins with the same orientation; the interaction energy of the antiferromagnetic (AF) Ising model is .] Although the CCs tend to fill in a triangular region, some gaps are evident near the lower apex and along the upper edge. Knowing just which values are allowed for a given lattice size is important if we are to implement entropic sampling with confidence. In this paper we present a method for systematically enumerating all CCs of Ising models in the plane.

Figure 1: Allowed configuration classes for system size ; and .

We study the spin-1/2 antiferromagnetic Ising model in an external magnetic field, whose energy is given by


where , indicates a sum over NN pairs of spins, is the external field, and is the number of spins; the model is defined on a square lattice of sites, with periodic boundary conditions. Unlike the ferromagnetic Ising model (), which exhibits a unique critical point in the plane and has an exact solution [7], the AF model () possesses a critical line, which is not completely understood.

Various approximate methods have been applied to determine the critical line of the AF Ising model on the square lattice [8][18]; these results, however, do not agree altogether. Binder and Landau [19] estimated the critical line via Monte Carlo simulation, obtaining very good agreement with the approximate closed-form expression of Müller-Hartmann and Zittartz [8], raising the possibility that the latter expression was in fact exact. (It was later shown that this is not so [12].) Hwang et al. [20] studied the AF Ising model on the square and triangular lattices using a microcanonical transfer matrix method and the Wang-Landau algorithm [3]. They performed an exact enumeration of the number of configurations, , and found CSs in the plane similar to that shown in Fig. 1.

Using the tomographic entropic sampling (TES) algorithm [6] we estimate , and associated microcanonical averages, for lattice sizes to 30, with an increment . We then calculate the canonical averages of the thermodynamic quantities of interest. Using these data we map out the critical line in the plane, and compare our results with several theoretical expressions.

Prior determination of the set of allowed CCs is an important tool to verify the quality of the sampling: we want to be sure that all CCs are visited during the simulation. Since this algorithm uses an initial guess, , to begin the study, it is convenient to use the final estimate (after the -th iteration) to obtain the initial guess of the next system size to be studied . As we will show, good initial estimates, , can be obtained using a two-dimensional interpolation because we know a priori the set of allowed CCs.

This paper is organized as follows. In Sec. 2 we define the basic CCs and the respective allowed values of for the spin-1/2 Ising model on the square lattice; the main goal of this section is to find all gaps in the plane. In Sec. 3, this information is used in simulations of the AF Ising model via the TES algorithm. There we describe the method used to determine via two-dimensional interpolation of the final estimate, , of the previous system size studied. Simulation results are reported in Sec. 4, for the order parameter and the staggered susceptibility as functions of and . Points along the critical line in the plane are obtained using finite size scaling analysis, and the results compared with several theoretical expressions. We summarize our findings in Sec. 5.

2 Configuration Classes

As pointed out above, one of the main problems in entropic sampling methods is the prior determination of the complete set of configuration classes for a given system size. Let us denote by and the number of up and down spins, respectively, on a square lattice with spins. The number of pairs of opposite spins, , and are related to and , respectively, via




Thus, once the possible values of are determined, so are those of . Note that , and can only take even values. The gaps in CS fall near the maximum and minimum values of ( and , respectively) for a given . Therefore, we will identify the possible values of near its maximum, , and minimum, , for a given ; note that the number of configurations is symmetric under interchange of and .

2.1 Determining and nearby classes

2.1.1 Compact configurations

Compact configurations consist of a square or rectangular cluster with up spins. To begin, consider the case of a square cluster of size with , as is illustrated in Fig. 2. It is evident that this configuration corresponds to the minimum value of , . A configuration with the same number of up spins, and is obtained by transferring an up spin from one of the corners of the square to an edge, as shown in Fig. 3. From this configuration, further rearrangements leading to , etc., are possible. When is not a square number, the most compact configuration (i.e., with the smallest perimeter) is a rectangle, or a square or rectangle with an incomplete layer of sites along one edge. For we have , while for , . In all cases, moving a corner site to an edge, one constructs a configuration with , and further arrangements yield additional increases in .

Figure 2: Basic compact configuration. Up and down spins are represented by “” and “”, respectively; wavy lines represent pairs of opposite spins.
Figure 3: Modified compact configuration. The new up and down spins are represented by “” (previously “”) and “” (previously “”), respectively; double straight and double wavy lines represent new pairs of identical and opposite NNs, respectively.

2.1.2 Extended configurations

Thus there are no gaps in the large- region due to compact configurations. Such compact configurations, however, are not necessarily the ones that minimize for a given value of . Consider for example the case with an integer. The cluster of up spins can be arranged to wrap around the lattice, yielding , independent of . We call such a configuration extended. For greater than a certain value, on the order of , the configurations that minimize are extended rather than compact. Suppose now that for a given , is obtained with an extended configuration, and that all compact configurations have strictly greater than . If lies between and , then the configuration that minimizes has at least two corners, and by the previous construction, configurations with exist. But if , the minimizing configuration has no corners and any modification yields a configuration with . This is how the gaps near the maximum values of arise.

Summarizing, if and are such that is obtained with a compact configuration, then there are configurations with , , …, etc., and no gap exists. If, on the other hand, is realized only for an extended configuration, and is an integer multiple of , then there are no configurations with .

2.2 Determining and nearby classes

The largest possible value of , , occurs in a configuration such that with spins arranged in a chess board (CB) configuration, such that all up spins have down spins as NNs and vice versa. One readily verifies that for , the maximum number of unlike NN pairs is . Starting from the CB configuration, we can reduce by exchanging an up and a down spin. If the exchanged spins are NNs, the resulting configuration has ; otherwise one has . Thus, for , there are no configurations with or . One readily verifies that configurations with , , etc., can be obtained via further exchanges of spins.

For , we have . Such configurations can be constructed by flipping one up spin in the CB configuration. Starting from this configuration, exchanging a pair of spins, one can reduce by 4, 6 or 8, but there is no rearrangement which reduces by just two. Thus for , there is no configuration having ; configurations with , , etc. do exist.

Finally, we note that for , there are no gaps in the neighborhood of . To verify this, consider a configuration obtained by flipping up spins in the CB configuration, so that , where . By hypothesis there are now at least four more sites with down spins than with up spins. Starting from a configuration with , in which each up spin is surrounded by down spins, we can create a single NN pair of up spins, with all six neighbors down, and with all remaining up spins completely surrounded by down spins. In this manner, is reduced by two. Configurations with , , etc. can be obtained by further exchanges of up and down spins.

Using the facts summarized above, it is straightforward to construct an algorithm that determines which values of are possible, for a given and , and thereby which values are accessible for a given system size. In the simulations reported below, we have verified that our entropic sampling scheme converges to visit all allowed classes.

3 Implementation

Using tomographic sampling, we study the antiferromagnetic Ising model in an external field on the square lattice; we consider periodic boundary conditions and NN interactions. The CCs of the systems are defined in the energy-magnetization space .

The TES method is applied in order to generate estimates of . For the smallest system size we begin with a guess, , obtained using a mean-field approximation. For subsequent system sizes, however, we use a two-dimensional interpolation of (the final estimate of for the smaller system size, ) to construct , where . For most studies we use five iterations, each one with initial configurations, which are simulated for lattice updates or Monte Carlo steps.

Let us denote by and the CCs that contain configurations and , respectively. The simulation uses a single spin-flip dynamics, so that the possible variations of and are and . At iteration , the acceptance probability for the transition is


These probabilities are stored in a table. For each configuration generated, be it a new one (if it is accepted) or the same (if it is rejected), we update the sums used to calculate the microcanonical and canonical averages of , , , and , where


is the order parameter (staggered magnetization); are the magnetizations of the two sublattices. At the end of each iteration the estimate of is refined according to


where is the histogram containing the number of times the CC is visited during the sampling, and is the average of over all accessible CCs; the acceptance probability [Eq. (4)] is updated using the new estimate of and the histogram, , is set to zero.

A single iteration of our method consists of ten independent simulations, each involving lattice updates, and each beginning from a different initial configuration. (By a lattice update we mean one attempted flip per spin; the initial configurations include, high, low, and intermediate interaction energies and both signs of the magnetization.)

3.1 Determining : mean field approach

As noted above, for the smallest system size we use an estimate of obtained via a mean-field approximation, specifically


where , and


To obtain this expression, we first note that , and use Stirling’s approximation. The dependence on is then obtained by estimating and for a given and , using a random-mixing approximation, and supposing that follows a Gaussian distribution.

In Fig. 4 we plot the estimate of given by Eq. (7). We present in Fig. 5 the final estimated value of after the 5 iteration of the simulation; this result is quite similar to that presented by Hwang et al. [20] for the square lattice, which was obtained using Wang-Landau sampling [3]. We note that the differences between the initial estimate obtained via mean field approximation and the final simulation result of are more evident near the edges of the CS, specially near to the maximum values of . To have a better idea of how close this initial estimate is to the final simulation result, we plot in Fig. 6 the difference between that estimate and the final result of the simulation for . Analyzing this figure, it is again clear that the differences are larger near the edges of the CS.

Figure 4: Estimate of via mean-field approximation.
Figure 5: Final estimate of after the 5 iteration of the simulation.
Figure 6: Difference between , estimated via mean-field approximation, and the final simulation result after the 5th iteration, .

3.2 Determining : interpolation

We expect that the closer the estimate is to the (unknown) exact value, , the faster the simulation will converge, and the more accurate our final estimate will be. Since the mean-field estimate worsens as the system size grows, we only use it for the smallest system size studied; initial estimates of subsequent system sizes are obtained by interpolating the final estimate, , of the previous system size studied. Our procedure is based on the existence of the limiting microcanonical entropy density as a function of the intensive parameters and (the bond and magnetization densities, respectively),


where and . (Since and are restricted to even integers we have and similarly for .) The idea is then to write


where , , and the r.h.s. is evaluated by extending to noninteger and via extrapolation and interpolation. Using this approach, we obtain better estimates, as is shown in Fig. 7. (Note that the largest differences continue to fall along the edges of the CS.)

Figure 7: Difference between , estimated by interpolating the final result of , and the final simulation result after the 5th iteration, .

3.3 Extrapolation and interpolation

Suppose we wish to construct the initial estimate on the basis of the simulation results for a smaller system, . We could do this via interpolation if every CC of the larger system were surrounded by four points of the smaller one in the - plane. Fig. 8 shows, however, that along the edges of the CS, the points corresponding to classes of the larger system are not surrounded by four points of the smaller one. For those points, one could in principle use extrapolation rather than interpolation. We found, however, that direct extrapolation yields poor estimates for .

Figure 8: Region of the plane with every possible configuration classes for system sizes and .

We obtain better estimates by first extrapolating the points along the edges of the CS for the smaller system, such that every point of the larger system is surrounded by four points of the smaller. Figure 9 shows the CS with accessible and extrapolated CCs for . Following this extrapolation we perform a linear two-dimensional interpolation as per Eq. (10).

Figure 9: Accessible and extrapolated configuration classes of a system of size .

4 Results

In this section we present results of the AF Ising model on the square lattice in an external field; periodic boundary conditions are employed. We use TES to simulate systems of sizes to 30, with . To calculate the uncertainties we perform five independent studies for each system size. We plot in Fig. 10 the staggered magnetization (order parameter) per site, , as a function of at . We can see that decreases considerably between and ; this behavior suggests a critical point, , marking a phase transition from the AF to the paramagnetic state. Figure 11 shows the staggered susceptibility per site,


as a function of at . As grows the peaks tend to the critical point, . The specific heat per site, , (not shown) has a similar behavior.

Figure 10: Staggered magnetization per site as a function of at , for to 30. The absolute uncertainties are plotted in the inset.
Figure 11: Staggered susceptibility per site as a function of at , for to 30. The highest peaks correspond to the largest system sizes. The absolute uncertainties are plotted in the inset.

4.1 Phase diagram

Using finite size scaling analysis [21] we estimate the critical line, , or, equivalently, , via the relations:




where is the field at which (the specific heat or the staggered susceptibility) takes its maximum for a given temperature and system size; is defined in an analogous manner. The estimates obtained using the maximum of and are averaged to yield and ; Figure 12 illustrates the procedure for estimating , for . Our estimated points for the phase boundary are shown in Fig. 13 along with the approximate expression derived by Müller-Hartmann and Zittartz [8]:


Our simulation results are in good agreement with their expression. For the critical field, the greatest relative difference between theory and simulation is about 0.9%, which occurs at .

Figure 12: Critical field determination using finite size scaling analysis: .
Figure 13: Phase diagram of the Ising AF on the square lattice. Comparison between simulation and the theoretical expression of Müller-Hartmann and Zittartz. Asterisks denote points obtained varying with fixed; circles denote points obtained varying , with fixed. The error bars of our results are smaller than the symbols.

In Table 1 we compare our simulation estimates for the critical magnetic field with some theoretical approximations. For temperature we find good agreement with the estimates of Monroe [17] (whose analysis involves a free parameter ), Wu and Wu [12], and Blöte and Wu [13], whereas there are significant discrepancies in relation to the other approximations. At higher temperatures, differences appear between simulation and the predictions of Monroe, Wu and Wu, and Blöte and Wu. These may reflect a small systematic error or an underestimate of simulation uncertainties. We intend to examine this issue in greater detail in future work.

TES Monroe Monroe MHZ WW BW WK
0.1 3.93304(16) 3.93307 3.93318 3.93069 3.93329 3.93330 3.93372
0.5 3.6648(8) 3.66506 3.66561 3.65309 3.66611 3.66614 3.67589
1.0 3.2906(14) 3.29303 3.29391 3.26843 3.29200 3.29261 3.31764
1.5 2.7258(14) 2.73243 2.73396 2.70401 2.73094 2.73176 2.75099
2.0 1.696(2) 1.71629 1.71872 1.69490 1.71492 1.71499 1.71512
TES: Tomographic entropic sampling [6]
MHZ: Müller-Hartmann and Zittartz [8]
WW: Wu and Wu [12]
BW: Blöte and Wu [13]
WK: Wang and Kim [16]
Table 1: Comparison between our simulation estimates of with some theoretical approximations.

5 Conclusions

The complete enumeration of CCs is of fundamental importance to study the antiferromagnetic Ising model using tomographic entropic sampling. The determination of CCs for entropic sampling of Ising models also enables us to obtain good initial estimates for the configuration numbers , using two-dimensional linear interpolation; the initial estimate is reasonably close to the final estimate . Despite the relatively small system sizes used in this study, we obtain a good estimate for the critical line in the temperature-magnetic field plane. Further details on critical behavior will be published elsewhere [22].


We are grateful to CNPq and CAPES, Brazil, for financial support.


  • [1] B. A. Berg and T. Neuhaus, Phys. Rev. Lett. 68, 9 (1992).
  • [2] J. Lee, Phys. Rev. Lett. 71, 211 (1993).
  • [3] F. Wang and D.P. Landau, Phys. Rev. Lett., 84, 10 (2001).
  • [4] F. Wang and D.P. Landau, Phys. Rev. E, 63, 056101 (2001).
  • [5] M.E.J Newman and G.T. Barkema. Monte Carlo Methods in Statistical Physics (Oxford University Press, New York, 2001), p. 161 - 169.
  • [6] R. Dickman and A.G. Cunha-Netto, Phis. Rev. E, 84, 026701 (2011).
  • [7] L. Onsager, Phys. Rev., 65, 117 (1944).
  • [8] E. Müller-Hartmann and J. Zittartz, Z. Phys., 73, 261 (1977).
  • [9] L. Sneddon, J. Phys. C, 12, 3051 (1979).
  • [10] R.R. Santos, J. Phys. C, 18, L1067 (1985).
  • [11] M. Kaufman, Phys. Rev. B, 36, 3697 (1987).
  • [12] X.N. Wu and F.Y. Wu, Phys. Lett. A, 144, 123 (1990).
  • [13] H.W.J. Blöte and X-N. Wu, J. Phys. A, 23, L627 (1990).
  • [14] M. Badehdah, A. Benyoussef and M. Touzani, J. Magn. Magn. Mater., 172, 254 (1997).
  • [15] X-Z. Wang and J.S. Kim, Phys. Rev. Lett., 78, 413 (1997).
  • [16] X-Z. Wang and J.S. Kim, Phys. Rev. E, 56, 2793 (1997).
  • [17] J.L. Monroe, Phys. Rev. E, 64, 016126 (2003).
  • [18] S.J. Penney, V.K. Cumyn and D.D. Betts, Physica A, 330, 507 (2003).
  • [19] K. Binder and D.P. Landau, Phys. Rev. B, 21, 5 (1980).
  • [20] C-O. Hwang, S-Y. Kim, D. Kang and J.M Kim, J. Stat. Mech: Theory Exp., 2007, L05001 (2007).
  • [21] V. Privman. Finite Size Scaling Analysis and Numerical Simulations of Statistical Systems (World Scientific, London, 1990).
  • [22] B.J. Lourenço and R. Dickman, in preparation.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description