Critical properties of the Hintermann-Merlini model
Many critical properties of the Hintermann-Merlini model are known exactly through the mapping to the eight-vertex model. Wu [J. Phys. C 8, 2262 (1975)] calculated the spontaneous magnetizations of the model on two sublattices by relating them to the conjectured spontaneous magnetization and polarization of the eight-vertex model, respectively. The latter conjecture remains unproved. In this paper, we numerically study the critical properties of the model by means of a finite-size scaling analysis based on transfer matrix calculations and Monte Carlo simulations. All analytic predictions for the model are confirmed by our numerical results. The central charge is found for critical manifold investigated. In addition, some unpredicted geometry properties of the model are studied. Fractal dimensions of the largest Ising clusters on two sublattices are determined. The fractal dimension of the largest Ising cluster on the sublattice A takes fixed value , while that for sublattice B varies continuously with the parameters of the model.
pacs:05.50.+q, 64.60.Cn, 64.60.Fr, 75.10.Hk
The exact solutionsIsing (); YangIsing () of the two-dimensional (2D) Ising model significantly promote the research of phase transitions and critical phenomena. After that, the Ising model becomes one of the most famous lattice models in statistical physics. Another famous Ising system, known as the Baxter-Wu modelbw (), is defined on the triangular lattice with pure three-spin interactions. The model was firstly proposed by Wood and Griffithswood () and exactly solved by Baxter and Wubw () by relating the model with the coloring problem on the honeycomb lattice. The solution gives the critical exponents () and (), which are exactly the same as those of the 4-state Potts modelpotts (); wfypotts (), which can be derived via the Coulomb gas theorydenNijs1 (); denNijs2 (); cg (). This means that the Baxter-Wu model belongs to the universality class of the 4-state Potts model. The critical properties of the 4-state Potts model are modified by logarithmic correctionslog2 (); log () due to the second temperature field, which is marginally irrelevantNBRS (); log2 (). With the two leading temperature fields simultaneously vanishing, the leading critical singularities of the Baxter-Wu modeldengbw () do not have logarithmic factors. Deng et al. generalize the Baxter-Wu model in dengbw () where the spins are allowed to be states ( can be larger than 2) and the up- and down-triangles can have different coupling constants. Both generalizations lead to discontinuous phase transitions.
where the sum takes over all the square unit cells and () are the coupling constants. The model is similar to the Baxter-Wu model in the sense that both are Ising systems with pure three-spin interactions. However the critical properties of the HM model are much more complicated. For the ferromagnetic case (), the model has four-fold degenerate ground states, as shown in Fig. 2. The nature of the phase transition breaking the symmetry of the order parameter cannot be determined by the dimensionality of the system and the symmetry properties of the ground statesSuzuki (). The exponents may vary with some tuning parameter without changing the symmetry of the order parameter. Such behavior has been found, e.g., in the Ashkin-Teller model (AT) AT (); ATMC (), the eight-vertex model 8vertex (); exactbook (), the 2D model in a four-fold anisotropic fieldXYh (), and the ferromagnetic Ising model with antiferromagnetic next-nearest neighbor couplings on square latticej1j2-henk (); j1j2-sandvik (); j1j2-guo ().
Through mapping to the eight-vertex model8vertex (); exactbook (), the free energy of the HM model has been found exactlyhintermann (). The critical manifold and critical exponent varying with the ratio of couplings were obtained. Based on the conjectured spontaneous magnetization 8v-mag () and the spontaneous polarization 8v-polar () of the eight-vertex model, Wu calculated the spontaneous magnetizations of the model on two sublattices. The results show that the two magnetizations possess different critical exponentswubw (). The spontaneous magnetization of the eight-vertex model was derived later 8v-mag-proof (), however, the spontaneous polarization remains a conjecture. It is thus very necessary to verify the results numerically.
In current paper, we study the critical behavior of the ferromagnetic HM model with and numerically. The numerical procedure includes transfer matrix calculations and Monte Carlo simulations. The critical properties we studied includes not only the verification of the analytic predictions, but also the fractal structure of spin clusters, which is studied for the first time.
The paper is organized in the following way: In Sec. II, we summarize the theoretical results of the model. In Sec. III, we introduce the transfer matrix method and present the associated numerical results. In Sec. IV, we describe the algorithm used in Monte Carlo simulations and give the associated numerical results. In Sec. V, we define Ising clusters on two sublattices and numerically determine the fractal dimensions of the corresponding largest clusters, respectively. We summarize in Sec. VI.
Ii Exact solutions
Assign arrows between the nearest-neighboring Ising spins of sublattice A (Fig. 3): if the two spins beside the arrow are the same, the arrow is rightward (or upward), otherwise it is leftward (or downward). There is a two-to-one correspondence between the Ising spin configurations and the arrow configurations.
After taking sum over the spin (in the center of the unit cell), the Boltzmann weights of the vertices are
Thus the model is mapped to the symmetric eight-vertex model8vertex (); exactbook (), which has been exactly solved by Baxter. For the ferromagnetic case considered in current paper, the critical manifold is given by
which leads to for the uniform case ().
The singularity of the free energy density is governed by
when is not an integer (The case that = an integer is not considered in current paper). is the critical temperature and is given by wubw ()
For the ferromagnetic case that we consider, is determined by (6). Since and , the critical exponent is thus found
For the uniform case, (7) gives .
A very interesting critical property of the HM model is that the spontaneous magnetizations and of the sublattice A and B possess different critical exponents:
The critical exponents and were obtained by Wuwubw ()
(It should be noted that the author reversely wrote the two critical exponents in Eq. (16) of wubw ().) According to the scaling law, this leads to two magnetic exponents and
Here has fixed value , while varies continuously with the parameters of the model.
In following sections, we will numerically study the critical properties of the model. The critical points and critical exponents will be numerically verified. Especially, Wu’s result of () is based on the unproved conjecture of the spontaneous polarization8v-polar () of the eight-vertex model, thus numerical verification is very necessary. Our numerical studies are focused on the subspace that and . By defining , critical manifold and exponents are expressed as functions of . It is worthy to notice that there is a symmetry for the transformation . At a special point and its dual , , the analytic results give . This is the 4-state Potts point.
Iii Transfer matrix calculations
As shown in Fig. 4, we define the row-to-row transfer matrix
where and are the states of two neighboring rows, respectively. Here periodic boundary conditions are applied. For a system with rows, the partition sum is found to be
with the periodic boundary conditions applied. In the limit , the free energy density is determined by the leading eigenvalue of
For the HM model, the dimension of the matrix is . To numerically calculate the eigenvalues of , we used the sparse matrix technique, which sharply reduces the requirement of computer memory for storing the matrix elements. (For details of this technique, see henk1982 (); henk1989 (); henk1993 (); Qian2004 ().) We are able to calculate the eigenvalues of with up to 22. In our calculations, we restrict the system size to even values, because the ordered configurations, as shown in Fig. 2, do not fit well in odd systems.
The critical properties can be revealed by calculating three scaled gaps :
where are three subleading eigenvalues of the matrix , respectively. According to the conformal invariance theoryconformal (), the scaled gap is related with a correlation length
where governs the decay of a correlation function . According to the finite-size scalingfss (), the gap in the vicinity of a critical point scales as
where is the corresponding scaling dimension, which is related to an exponent according to the conformal invarianceconformal (). is an irrelevant field, and is the corresponding irrelevant exponent. ,, and are unknown constants.
We focus on the the energy-energy correlation and two magnetic correlations and . Generally speaking, the magnetic correlation function is defined as . For the HM model, since the spontaneous magnetizations on the two sublattices behave differently, we look two types of magnetic correlation functions: with and on the A sublattice and on the B sublattice.
Let be the largest eigenvalue in the subspace that breaks the spin up-down symmetry, which means that the associated eigenvector satisfies
where is the operator flipping spins. Thus, the scaled gap is identified as and the corresponding correlation is .
Let and be the second and the third largest eigenvalue in the subspace keeping the spin up-down symmetry. It is not a priori clear which of the two corresponding gaps is the thermal one or the magnetic one . Based on the magnitudes of and , we identify and as and , respectively.
We then numerically solve the finite-size scaling equation
for , respectively. The solution satisfies
Here is an unknown constant. Since and , converges to the critical point with increasing system sizes. (19) is used to determine the critical point in our numerical procedure. For the uniform case, we obtain , which is in good agreement with the exact solution . We have also estimated the critical points of the cases with , and 5, where and . Both and have been used to estimate the critical points, we list the best estimations in Table 1. Our numerical estimations of are consistent with the theoretical results in a high accuracy.
Exactly at the critical point , (16) reduces to
which is used to determine the scaling dimensions and . For the uniform case , it gives , and , which are consistent with the analytic results (Table 1). We applied the same procedure to . The numerical estimations and the associated analytic results are consistent, see Table 1.
We also calculate the scaled gap in the vicinity of critical point. The scaling behavior of seems different from that of or , as shown in Fig. 7 for the uniform case. This is the result of in Eq. (16). The location of the minimum of for a given can be denoted as , which converges to the critical point when . This property can also be used to estimate , e.g., as was done in henkcubic (). But we don’t do this in current paper. By calculating at the estimated , we obtain the scaling dimension according to (20). The results are listed in Table 1.
Fitting the data of the free energy density according to (21), we obtain the central charge for the uniform case. For the other cases, the numerical estimations of the central charge are listed in Table 1. The central charge for all takes the fixed value as expected for transitions breaking four-fold symmetry.
Furthermore, we consider the 4-state Potts point . A finite-size scaling analysis based on transfer matrix calculations is performed at this point. The numerical estimations obtained are listed in Table 1. No logarithmic corrections are found, as in the Baxter-Wu model.
Iv Monte Carlo simulations
The Baxter-Wu model has been simulated by the Metropolis algorithmbwMC (), the Wang Landau algorithmwanglandau (), and a cluster algorithmhenkbw (); dengbw (). The cluster algorithm is similar to the Swendsen-Wang algorithmsw () for the Potts model. In this algorithm, the triangular lattice is divided as three triangular sublattices. By randomly freezing the spins on one of the sublattices, the other two sublattices compose a honeycomb lattice with pair interactions, then a Swendsen-Wang type algorithm can be formed to update the spins on this honeycomb lattice. In current paper, we suitably modify this cluster algorithm to simulate the HM model.
The algorithm is divided into three steps:
Step 1, update the spins on sublattice A.
Let the spins on sublattice B be unchanged (frozen), the interactions of the spins reduce to pair interactions on sublattice A.
A. bonds. An vertical (horizontal) edge on sublattice A belongs to two triangles, the interactions of the left (upper) and right (down) triangles can be denoted and respectively, which may take the values , or . Place a bond on this edge with probability if both the products of spins in the left (upper) triangle and in the right (down) triangle are 1, if only the product in the left (upper) triangle is 1, if only the product in the right (down) triangle is 1, otherwise.
B. clusters. A cluster is defined as a group of sites on sublattice A connected through the bonds.
C. update spins. Independently flip the spins of each cluster with probability 1/2.
Step 2, update the spins on sublattice B and half of the spins on sublattice A.
This step is very similar to step 1, but we freeze only half of the spins on sublattice A, which are labeled (Fig. 8). The other spins (on sublattice B+) are updated, whose interactions also reduce to pair interactions when spins are frozen.
Step 3, update the spins on sublattice B and the other half of the spins on sublattice A.
In this step, the spins labeled (in Fig. 8) are frozen, other spins (on sublattice B+) are updated.
In a complete sweep, all the spins are updated twice.
In the simulations, the sampled variables include the specific heat , the magnetizations on sublattice A, and on sublattice B.
The specific heat is calculated from the fluctuation of energy density
where is the linear size of the system, and means the ensemble average. The magnetizations are defined as
where is the number of total sites on sublattice A or B.
In order to demonstrate our numerical procedure based on the Monte Carlo simulations, we take the uniform case () as an example. The cluster algorithm is very efficient, which easily allows us to do meaningful simulations for system with linear size up to . All the simulations are performed at the critical point . After equilibrating the system, samples were taken for each system size. Using the Levenberg-Marquardt algorithm, we fit the data of according to the finite-size scaling formulafss ()
where are the leading correction-to-scaling term with the leading irrelevant exponent. is the dimensionality of the lattice. , , and are unknown parameters. The fitting yields , which is in good agreement with the exact result .
The finite-size scaling behaviors of the magnetizations and are
respectively. The log-log plot of and versus are shown in Fig. 9. The fits yield , and , where and are in good agreement with the analytic predictions and , respectively.
We also simulate the , and cases. The numerical estimations of and the corresponding analytic predictions are listed in Table 1. They are in good agreement. In these fits, no logarithmic corrections to scaling are found.
V Fractal structure of the model
In addition to calculate the critical exponents , , and , we also investigate the geometry properties of the model.
The configurations of the HM model are represented by Ising spins. Thus we can define “Ising clusters” for the model as in the Ising modelIsingCluster (). To do so, each sublattice is viewed as a square lattice. For two nearest-neighboring Ising spins on a sublattice, they are considered to be in the same cluster if they have the same sign. Here, in saying “nearest-neighboring”, the neighborhood of the site is the same to that of a site in a square lattice. This definition is applied to the two sublattices respectively, thus there are two types of Ising clusters. For clarity, we call the Ising cluster on sublattice A as “A cluster”, while the Ising cluster on sublattice B as “B cluster”.
The simulations are also performed at the critical points. The size of the largest A cluster and the largest B cluster are denoted and , respectively. They satisfy the finite-size scaling
which means that the largest A cluster and the largest B cluster are fractals, with and the fractal dimensions, respectively. Figure 10 is an illustration of and versus system size at the uniform case. Fitting the data according to (28) and (29), we find and .
The same procedure were applied to and . The results obtained are listed in Table. 2. The fractal dimension of the A cluster has a fixed value , while that of the B cluster varies with the ratio .
Vi Conclusion and discussions
In conclusion, we have numerically studied the critical properties of the HM model by means of finite-size scaling analysis based on the transfer matrix calculations and the Monte Carlo simulations. For the critical points and the critical exponents (), (), and (), our numerical estimations are in good agreement with the corresponding analytic predictionshintermann (); wubw (). The analytic prediction for is based on the conjectured spontaneous polarization of the eight-vertex model, which remains unproved. Our numerical results verified the correctness of the prediction.
In addition the central charge of the model is found to be . This is consistent with the fact that , and vary continuously with the parameters, which is related to the four-fold degeneracy of the ground states of the model. Usually, logarithmic correctionslog (); lvAT (); j1j2-guo () show up in the 4-state Potts point for the model with 4-fold degenerate ground states due to the second temperature field. However, we do not see such logarithmic corrections in the 4-state Potts point of the HM model. This is the same as the BaxterâWu model, in which the amplitude of the marginally irrelevant operator is zero.
Furthermore, two unpredicted critical exponents are found to describing geometry properties of the model. We define clusters based on the Ising-spin configurations on the two sublattices of the model. The fractal dimension of the largest cluster on sublattice A takes fixed value , while the fractal dimension of the largest cluster on sublattice B varies continuously with the parameter of the model.
Ding thanks F. Y. Wu for valuable discussion in terms of understanding the analytical results of HM model. This work is supported by the National Science Foundation of China (NSFC) under Grant No. 11205005 (Ding) and 11175018 (Guo).
- (1) L. Onsager, Phys. Rev. 65, 117 (1944).
- (2) C. N. Yang, Phys. Rev. 85, 808 (1952).
- (3) R. J. Baxter and F. Y. Wu, Phys. Rev. Lett. 31, 1294 (1973).
- (4) D. W. Wood and H. P. Griffiths, J. Phys. C 5, L253 (1972).
- (5) R. B. Potts, Proc. Camb. Phys. Soc. 48, 106 (1952).
- (6) F. Y. Wu, Rev. Mod. Phys. 54, 235 (1982).
- (7) M. P. M. den Nijs, J. Phys. A 12, 1857 (1979).
- (8) M. P. M. den Nijs, Phys. Rev. B 27, 1674 (1983).
- (9) B. Nienhuis, J. Stat. Phys. 34, 731 (1984).
- (10) B. Nienhuis, A. N. Berker, E. K. Riedel, M. Schick, Phys. Rev. Lett. 43, 737 (1979).
- (11) J. Salas and A. D. Sokal, J. Stat. Phys. 88, 567 (1997).
- (12) M. Nauenberg, D. J. Scalapino, Phys. Rev. Lett. 44, 837 (1980).
- (13) Y. Deng, W.-A. Guo, J. R. Heringa, H. W. J. Blöte, and B. Neinhuis, Nucl. Phys. B 827[FS], 406 (2010).
- (14) A. Hintermann and D. Merlini, Phys. Lett. A, 41, 208 (1972).
- (15) M. Suzuki, Prog. Theor. Phys. 51, 1992 (1974).
- (16) J. Ashkin and E. Teller, Phys. Rev. 64, 178 (1943).
- (17) S. Wiseman and E. Domany, Phys. Rev. E 48, 4080 (1993).
- (18) J.-P. Lv, Y. Deng, and Q.-H. Chen, Phys. Rev. E 84, 021125(2011).
- (19) M. P. Nightingale and H. W. J. Blöte, Physica A 251, 211 (1998).
- (20) S. Jin, A. Sen, and A. W. Sandvik, Phys. Rev. Lett. 108, 045702 (2012).
- (21) S. Jin, A. Sen, W.-A. Guo, and A. W. Sandvik, Phys. Rev. B 87, 144406(2013).
- (22) R. J. Baxter, Phys. Rev. Lett. 26, 832 (1971).
- (23) R. J. Baxter, Exactly Solved Models in Statistical Mechanics (Academic, London, 1982).
- (24) J. V. José, L. P. Kadanoff, S. Kirkpatrick, and D. R. Nelson, Phys. Rev. B 16, 1217 (1977).
- (25) M. N. Barber and R. J. Baxter, J. Phys. C 6, 2913 (1973).
- (26) R. J. Baxter and S. B. Kelland, J. Phys. C 7, L403 (1974).
- (27) F. Y. Wu, J. Phys. C 8, 2262 (1975).
- (28) R. J. Baxter, J. Stat. Phys. 15, 485 (1976).
- (29) H. W. J. Blöte and M. P. Nightingale, Physica A (Amsterdam) 112, 405 (1982).
- (30) H. W. J. Blöte and B. Nienhuis, J. Phys. A 22, 1415 (1989).
- (31) H. W. J. Blöte and M. P. Nightingale, Phys. Rev. B 47, 15046 (1993).
- (32) X. F. Qian, M. Wegewijs, and H. W. J. Blöte, Phys. Rev. E 69, 036127 (2004).
- (33) J. L. Cardy, in Phase Transitions and Critical Phenomena, edited by C. Domb and J. L. Lebowitz (Academic Press, London, 1987), Vol. 11, p. 55, and references therein.
- (34) For reviews, see e.g. M. P. Nightingale in Finite-Size Scaling and Numerical Simulation of Statistical Systems, ed. V. Privman (World Scientific, Singapore 1990), and M. N. Barber in Phase Transitions and Critical Phenomena, Vol. 8, eds. C. Domb and J. L. Lebowitz (Academic, New York, 1983).
- (35) H. W. J. Blöte, M. P. Nightingale, Physica A (Amsterdam) 129, 1 (1984).
- (36) H. W. J. Blöte, J. L. Cardy, and M. P. Nightingale, Phys. Rev. Lett. 56, 742 (1986).
- (37) I. Affleck, Phys. Rev. Lett. 56, 746 (1986).
- (38) L. N. Shchur and W. Janke, arXiv:1007.1838.
- (39) N. Schreiber and J. Adler, J .Phys. A 38, 7253 (2005).
- (40) H. W. J. Blöte, J. R. Heringa, and E. Luijten, Comp. Phys. Comm. 147, 58 (2002).
- (41) R. H. Swendsen and J. S. Wang, Phys. Rev. Lett. 58, 86 (1987).
- (42) W. Janke and A. M. J. Schakel, Phys. Rev. E 71, 036703 (2005).