Oct 20, 2016 - that do not assure the model is in the BKT universality class. In the Fisher zeros approach one starts with a discrete canonical partit...

0 downloads 9 Views 636KB Size

arXiv:1507.02231v3 [cond-mat.stat-mech] 20 Oct 2016

J.C.S. Rochaa,b,∗, L.A.S. M´ola , B.V. Costaa a

Laborat´orio de Simula¸ca˜ o, Departamento de F´ısica, ICEx Universidade Federal de Minas Gerais, 31720-901 Belo Horizonte, Minas Gerais, Brazil b Departamento de F´ısica, ICEB, Universidade Federal de Ouro Preto, 35400-000 Ouro Preto, Minas Gerais, Brazil

Abstract Using the two dimensional XY − (S (O(3)) model as a test case, we show that analysis of the Fisher zeros of the canonical partition function can provide signatures of a transition in the Berezinskii-Kosterlitz-Thouless (BKT ) universality class. Studying the internal border of zeros in the complex temperature plane, we found a scenario in complete agreement with theoretical expectations which allow one to uniquely classify a phase transition as in the BKT class of universality. We obtain T BKT in excellent accordance with previous results. A careful analysis of the behavior of the zeros for both regions Re(T ) ≤ T BKT and Re(T ) > T BKT in the thermodynamic limit show that Im(T ) goes to zero in the former case and is finite in the last one. Keywords: Phase transitions: general studies, Classical spin models, Monte Carlo methods PACS: 05.70.Fh, 75.10.Hk, 05.10.Ln 1. Introduction The Berezinskii-Kosterlitz-Thouless (BKT ) transition already has more than 40 years of history [1] and is still intriguing the scientific community. The nature of this transition is completely different from the common discontinuous (first order) or continuous (second order) phase transitions. Long range order does not exist and the two point correlation function has an algebraic decay at low temperature (T ≤ T BKT ) and an exponential decay for T > T BKT [2]. Here T BKT is known as the BKT temperature, and a model displaying a BKT transition has an entire line of critical points for T ≤ T BKT . In addition to that, the corresponding free energy, which is a C ∞ function, is not analytical in this region. Besides these striking features, the correlation function has a characteristic universal exponent decay η(T BKT ) = 1/4 at the transition temperature. Its phenomenology relies on the belief that it is driven by a vortex-antivortex unbinding mechanism [3, 4]. Another proposition that also describes the transition is based on a polymerization of domain walls [5]. Many systems, e.g., superfluid films, ∗ Corresponding

author Email addresses: [email protected] (J.C.S. Rocha), [email protected] (L.A.S. M´ol), [email protected] (B.V. Costa) Preprint submitted to Computer Physics Communications

Coulomb gases and crystal surface roughening, undergo transitions that can be classified as belonging to this universality class [6]. More recently, Berker et al [7] found that a BKT transition can occur in a scale free network. Although the BKT transition is well known, the characterization of an unknown phase transition as being in the BKT universality class is not an easy task since there is no standard method to do so. The lack of a criterion capable of determining the nature of a transition beyond any reasonable doubt is a problem discussed by Bramwell and Holdsworth [8]. They pointed out that to be able to see the transition the system under investigation should be very large. They estimate that for “a system with atomic spacing of 3Å the area should correspond to the size of a postage stamp”. From the analytical point of view, the renormalization group approach is able to describe the main features associated with the transition. However, the approximations involved in the study of a given model may hinder the discovery of its real nature. In this paper we present an algorithm for the systematic analysis of the Fisher zeros of the canonical partition function [9, 10, 11] for the 2D XY model (with spin symmetry O(3)) looking for possible signatures of the BKT transition. It was Fisher [11], in 1964, who proposed considering the partition function zeros October 21, 2016

1

in the complex temperature plane to study phase transitions. As known, the thermodynamic behavior of a given physical system is encoded by its partition function, Z, and all thermodynamic quantities can be obtained as derivatives of the free energy, F = −kBT ln Z. The basic assumption of the Fisher zeros approach to study phase transitions is that a system undergoes a phase transition at a given (real) temperature, T c , if Z(T c ) = 0, reflecting the non-analyticity of F at T c [12]. Since in a BKT transition there is an entire line of critical points for T ≤ T BKT , one should expect that a map of the Fisher zeros in the complex temperature plane exhibits an entire line of zeros in the real temperature axis for T ≤ T BKT , signaling the BKT behavior of the transition, while for T > T BKT the zeros should not touch the real positive axis, emphasizing the analytical behavior of thermodynamic functions at high temperatures. Of course, these considerations apply only to the thermodynamic limit. We warn the reader that the Fisher zeros studied here should not be confused with the Yang-Lee zeros [9, 10] defined on the complex fugacity plane instead of the temperature plane [13, 14, 15]. In what follows we will analyze the Fisher zeros map for the XY-model.

0.5

Im(x)

0

-0.5

0.04

0.02

rz 0

0

-1 -1

0.2

-0.5

0.4

0

0.5

1

Re(x) Figure 1: Fisher zeros map on the x = e−βε plane for the 2D classical XY-model in a 50 × 50 lattice. The inset shows a zoom on the inner region for five distinct simulations, represented by different symbols. The cusp, rz , is indicated.

or the more reliable helicity modulus [6, 18], quantities that do not assure the model is in the BKT universality class. In the Fisher zeros approach one starts with a discrete canonical partition function that can be written as Z = P E g(E) exp (−βE), where g(E) stands for the density of states (DOS). For a continuous system one may perform a discretization of the DOS [22]. This can be done by dividing the energy range into M bins of size ε. We can organize the energies inside the interval [E0 , E M−1 ] by counting the energy of the nth bin as En = E0 + nε. By defining a variable x ≡ e−βε , we get a discretized version,

2. Model and Methods Here we study a “fruit fly” model of the BKT transition [6], the classical two-dimensional XY-model on a square lattice, defined by X y y H = −J (S ix S xj + S i S j ). (1) hi, ji

The sum runs over the nearest neighbors, J stands for the exchange coupling constant and S iα stands for the component α = (x, y, z) of the ith spin. The same Hamiltonian also defines the Planar-Rotator (O(2)) model [16, 17], whose spins have only two components (one degree of freedom), and can also be viewed as an example belonging to the BKT universality class. In spite of the lack of long-range order for the model [2], there is a non-zero magnetization for any finite volume [18], resulting in a thermodynamic behavior very similar to that observed at continuous phase transitions. In fact, the behavior can be easily confused with a second order phase transition or something else [19, 20, 21]. Thus, distinguishing between continuous and BKT behavior in finite systems is difficult and may require system sizes beyond those feasible. Usually, the BKT temperature is estimated by using the Binder cumulant, the divergence of the magnetic susceptibility, the correlation functions

ZD = e−βE0

M−1 X n=0

gn xn = e−βE0

M−1 Y

(x − x j ),

(2)

j=1

where gn = g(En ), and x j is the jth zero of the polynomial, usually called Fisher zero. Once gn is obtained, any reliable zero finder can be used 1 to find the x j ’s. As long as gn ∈ R+ , Z is an analytical function and has no real positive zeros for finite systems. The zeros show up as complex conjugate pairs. The analysis of Fisher zeros in finite systems is done by considering the special set of zeros {x⋆ (L) = a(L) + ib(L)} ∈ x j (L), where L is the linear size of the system, called first or leading zeros. They have the following properties: b(L) → 0 as L → ∞ 1 In this work we decided to use: E.W. Weisstein, “Polynomial roots”, in MathWorld – A Wolfram Web Resource, http://mathworld.wolfram.com/PolynomialRoots.html.

2

while limL→∞ a(L) = a(∞), a constant value. The leading zeros are very stable against statistical fluctuations, in contrast to non-leading zeros, and are, in general, featured in the map (see, for example, Ref. [23]). They are related to the transition temperature of the system in the thermodynamic limit as 1/kB T c = − ln(a(∞))/ε [23] and their impact angle on the real positive axis is directly related to the order of the phase transition [24]. Recent results also suggest that the pattern of zeros as a whole can be used to characterize the order of the transition [25]. In a paper of 1983 Ytzykson et al. [26], analysing the zeros distribution for small lattices, found that some properties of the Ising and gauge models seem to exhibit a universal behavior close to complex singularities. However, it is not clear how to extend their arguments to the BKT transition. In order to obtain the DOS of the XY-model we used the Replica Exchange Wang-Landau (REWL) method [27, 28, 29, 30], a parallel version of WangLandau (WL) sampling [31, 32, 33], capable of sampling the entire configuration space efficiently in a single simulation. In this scheme the energy range is split into smaller overlapping windows. We considered e0 = −1.9J, e M−1 = 0, and an overlap of 75%, where e = E/L2 stands for the energy per spin. Several different random walkers are allowed to run in each of these windows following the original WL scheme2 . We use regular square lattices with sizes ranging from L = 10 up to 200. In this work we chose J = 1, S = 1, kB = 1, and the lattice parameter a = 1.

0.2 0.1

Im(x)

0 L = 10 L = 20 L = 30 L = 40 L = 50 L = 60 L = 80 L = 100

-0.1 -0.2 0

0.2

0.4

0.6

0.8

1

Re(x) Figure 2: Zoom on the real positive semi-axis of the zeros maps in the x plane for L = 10 → 100. The zeros on the internal border are highlighted.

low temperatures. The line for Re(x) > rz on the other hand should not touch the real axis. The cusp rz should, then, give T BKT . In order to systematically investigate the finite size effects and confirm the above expectations, we opted to work in the complex temperature (T ) plane instead of the complex x = e−ε/kB T plane. This choice is justified by the fact that finite size scaling in terms of the temperature is known [15] and that by doing so our results are almost independent on the chosen bin size, ε. This last point is especially important when dealing with bigger lattices. In this work we used ε = 1 for L ≤ 100 and ε = 3 for L = 200; we noted that small modifications in the bin size did not change our results. To determine the limits of the internal border we divided the real temperature axis into small bins of size ∆T x centered at T x = Re(T ). The border line was identified by looking inside a given bin for the smallest value of T y = Im(T ) that appears. At least five different simulations were used for each lattice size. The resulting curves are presented in Fig. 3. Fig. 4 shows T y as a function of L−1 for some typical T x ≤ T BKT values. The solid lines are linear regressions showing that T y → 0 for L → ∞, i.e., the internal border coalesces with the real positive axis in the thermodynamic limit. Other scaling functions and exponents were also tested (not shown here), but they do not describe our data as well as this ansatz.

3. Results and discussion In Fig. 1 we show a typical zeros map of the imaginary and real parts of all x j ’s, for a 50 × 50 lattice. In a continuous phase transition a single leading zero is expected. However, we cannot identify a single point featured in comparison to others. The inset shows in different symbols the zeros for five different simulations for L = 50, in order to analyze statistical fluctuations. A leading zero cannot be identified in this picture. Instead, a cusp at rz is evident and the border line is quite stable against fluctuations. The size dependence of the internal border of the map pattern as a function of the lattice size is shown in Fig. 2. One should expect the internal border to coalesce into the real axis for Re(x) ≤ rz and L → ∞, in accordance with the existence of an entire line of critical points at

To estimate the critical temperature we used the location of the cusp position, T z (L). Since the scaling function for the “pseudocritical” temperature is given by [15] T BKT (L) ∼ [ln(L)]−2 we plotted T z (L) as a function of [ln(L)]−2 in the Fig. 6. A linear regression, discarding the point corresponding to L = 10, gives T BKT = 0.709(2), and discarding the points for L < 40

2 The considered flatness criteria is p = 0.7, the final ln f value is 10−9 , and the acceptance ratio is 60%.

3

0.25

2.5

0.2

2 L = 10 L = 20 L = 30 L = 40 L = 50 L = 60 L = 70 L = 80 L = 90 L = 100 L = 200

0.15

Ty 0.1 0.05 0 0

1

0.5

1.5

Ty 1

0 0

2

1.5

ω = 2.0 2.3 2.5

0.5

Tx

0.1

0.2

0.3 -ω

L

Figure 3: Internal border of zeros obtained by using the binning process for L = 10 → 200. Here ∆T x = 0.1 and the error bars represent statistical fluctuations.

0.4

0.5

0.6

0.7

-3

(10 )

Figure 5: Typical T y × L−ω behavior. Here T x = 3.0 > T BKT . Linear fit gives support to a power law behavior with exponent ω ≈ 2.4.

0.3 0.25 0.2

Ty

0.15

tion at high temperatures.

Re(T) = 0.2 Re(T) = 0.3 Re(T) = 0.4 Re(T) = 0.5 Re(T) = 0.6 Re(T) = 0.7

0.85

0.8

0.1

Tz

0.05 0 0

0.75 0.02

0.04

0.06

0.08

0.1

0.12

1/L 0.7

Figure 4: Finite size scaling analysis of the imaginary part of the internal border according to the ansatz T y ∼ L−1 . The lines represent a linear regression of the data, showing good agreement with the ansatz and the convergence to zero or small negative values for L → ∞. We use ∆T x = 0.1.

0

0.05

0.1

[ln(L)]

0.15

0.2

-2

Figure 6: Finite size scaling of the cusp position, T z (L), according to the prediction for the BKT “pseudocritical” temperature T BKT (L) ∼ [ln(L)]−2 . Discarding the point corresponding to L = 10 the BKT temperature is estimated by a linear regression as 0.709(2) (dashed blue line) and discarding L < 40, T BKT = 0.704(3) (solid red line).

gives T BKT = 0.704(3), which agrees very well with previous results [34, 18], 0.700(5) and 0.700(1) respectively. More precise results could be obtained by increasing the number of zeros in the map by reducing ε, increasing the precision in the location of the cusp. However, the zeros finder task may become a problem for such a high degree polynomial (> 3 × 104 ). In any case, conventional methods to locate the BKT may be used together with this one to get even more precise results for T BKT .

To analyse the behavior of T y in the non-critical region we have to be careful. First we have to consider the following. An analysis based in finite size scaling, as we did in the critical region, is meaningless since the basic assumption of any FS S is that the free energy behaves as a homogeneous function. Following Ytzykson we identify three regions. (1) The region far from T BKT where L >> ξ ≈ 1, where finite size lattice effects are not relevant. Here ξ is the correlation length which for the XY model diverges exponentially when T BKT is approached from above, remaining infinite in the critical region. (2) The region T & T BKT where L >> ξ >> 1. The finite lattice exhibits the same scaling behavior as the infinite system. And, (3) the region

From Figs. 2 and 3 one can see that the curves diverge from the positive real axis, in accordance with the expectation that for T > T BKT the imaginary part of the zeros should remain finite. Moreover, since the free energy has to be an analytic function at high temperatures, there can be no real positive zero of the partition func4

where L ≈ ξ >> 1 for which are observed severe finite size effects. Regions (2) and (3) are scaling regions. With that in mind we content ourselves with an analysis in a region far enough from T BKT where we may expect that the border of the zeros map should converge fast to its asymptotic limit. For T > T BKT the free energy is an analytical function of T so that it can be approximated by a polynomial in T with the coefficients depending on the lattice size in powers of 1/L. Farthest we are from T BKT faster will decay the higher coefficients of the polynomial, so that, it is reasonable to suppose that in this region T y has a power law T y (∞) + AL−ω behavior. Asymptotic convergence is assured by the Brouwer fixed point theorem [35]. Reasoning in this way we can plot T y for the largest lattices of our simulations as a function of L−ω . By varying ω we can find the best value that adjusts a linear function to our data. This procedure is shown in Fig. 5. It is noteworthy the different “scaling” behavior of Im(T ) in both regions, T < T BKT and T > T BKT . As should be expected the exponent ω is not universal, depending on T .

tional Physics at UGA where part of this work was realized. This work was partially supported by CNPq and Fapemig, Brazilian Agencies. JCSR thanks CNPq for the support under grant PDJ No. 150503/2014-8 and Grant CNPq 402091/2012-4. References References [1] J. V. Jos´e, ed., 40 Years History of Berezinskii-KosterlitzThouless Theory (World Scientific, 2013). [2] N. D. Mermin and H. Wagner, Phys. Rev. Lett. 17, 1133 (1966). [3] V. L. Berezinskii, Journal of Experimental and Theoretical Physics - JETP 32, 493 (1971). [4] J. M. Kosterlitz and D. J. Thouless, Journal of Physics C: Solid State Physics 6, 1181 (1973). [5] A. Patrascioiu and E. Seiler, Phys. Rev. Lett. 60, 875 (1988). [6] P. Minnhagen and P. Olsson, Phys. Rev. B 44, 4503 (1991). [7] A. N. Berker, M. Hinczewski, and R. R. Netz, Phys. Rev. E 80, 041118 (2009). [8] S. T. Bramwell and P. C. W. Holdsworth, Journal of Applied Physics 75, 5955 (1994). [9] C. N. Yang and T. D. Lee, Phys. Rev. 87, 404 (1952). [10] T. D. Lee and C. N. Yang, Phys. Rev. 87, 410 (1952). [11] M. E. Fisher, in Lectures in Theoretical Physics: Volume VII C - Statistical Physics, Weak Interactions, Field Theory : Lectures Delivered at the Summer Institute for Theoretical Physics, University of Colorado, Boulder, 1964, edited by W. Brittin (University of Colorado Press, Boulder, 1965), v. 7. [12] B.-B. Wei, S.-W. Chen, H.-C. Po, and R.-B. Liu, Sci. Rep. 4, 5202 (2014). [13] R. Kenna and A. Irving, Nuclear Physics B 485, 583 (1997). [14] A. C. Irving and R. Kenna, Phys. Rev. B 53, 11568 (1996). [15] R. Kenna and A. Irving, Physics Letters B 351, 273 (1995). [16] P. Olsson, Phys. Rev. B 52, 4526 (1995). [17] M. Hasenbusch, Journal of Physics A: Mathematical and General 38, 5869 (2005). [18] B. V. Costa, P. Z. Coura, and S. A. Leonel, Physics Letters A 377, 1239 (2013), ISSN 0375-9601. [19] L. A. S. M´ol and B. V. Costa, Phys. Rev. B 79, 054404 (2009). [20] L. A. S. M´ol and B. V. Costa, Journal of Physics: Condensed Matter 22, 046005 (2010). [21] L. A. S. M´ol and B. V. Costa, Journal of Magnetism and Magnetic Materials 353, 11 (2014). [22] J. C. S. Rocha, S. Schnabel, D. P. Landau, and M. Bachmann, Phys. Rev. E 90, 022601 (2014). [23] J. C. S. Rocha, S. Schnabel, D. P. Landau, and M. Bachmann, Physics Procedia 57, 94 (2014). [24] W. Janke, D. Johnston, and R. Kenna, Nuclear Physics B 682, 618 (2004), ISSN 0550-3213. [25] M. P. Taylor, P. P. Aung, and W. Paul, Phys. Rev. E 88, 012604 (2013). [26] R. P. C. Itzykson and J. Zuber, Nuclear Physics B 220, 415 (1983). [27] T. Vogel, Y. W. Li, T. W¨ust, and D. P. Landau, Phys. Rev. Lett. 110, 210603 (2013). [28] T. Vogel, Y. W. Li, T. W¨ust, and D. P. Landau, Journal of Physics: Conference Series 487, 012001 (2014). [29] Y. W. Li, T. Vogel, T. W¨ust, and D. P. Landau, Journal of Physics: Conference Series 510, 012012 (2014). [30] T. Vogel, Y. W. Li, T. W¨ust, and D. P. Landau, Phys. Rev. E 90, 023302 (2014).

4. Conclusion In summary, we show that a method of investigating the Fisher zeros of the partition function can identify whether or not the model exhibits a BerezinskiiKosterlitz-Thouless (BKT) phase transition. By studying the 2D XY-model we found a qualitative picture that is completely consistent with expectations for the BKT transition, i.e., the zeros map is consistent with the existence of an entire line of zeros in the real positive axis in the thermodynamic limit, with different behaviors for points below and above the transition temperature that could be used to signalize the BKT transition. Moreover, the BKT transition temperature was successfully obtained by considering the location of the cusp that splits regions with different behaviors giving T BKT = 0.704(3), in excellent accordance with previous results [18, 34]. Our quantitative analysis shows that for temperatures above T BKT the imaginary part of the zeros converge to finite values. This behavior is in accordance with the fact that the free energy must be an analytical function over the real axis at high temperatures, which prevents the existence of any real positive zero at high T . Acknowledgments We would like to thank Prof. David P. Landau for very fruitful discussions and the Center for Simula5

[31] F. Wang and D. P. Landau, Phys. Rev. Lett. 86, 2050 (2001). [32] F. Wang and D. P. Landau, Phys. Rev. E 64, 056101 (2001). [33] D. P. Landau, S.-H. Tsai, and M. Exler, American Journal of Physics 72, 1294 (2004). [34] H. G. Evertz and D. P. Landau, Phys. Rev. B 54, 12302 (1996). [35] M. H. Stone, Mathematics Magazine 21, 237 (1948).

6