---

Critical Behaviour of the 3D XY -Model:

arXiv:cond-mat/9305020v1 18 May 1993KL-TH-93/10CERN-TH.6885/93Critical Behaviour of the 3D XY -Model:A Monte Carlo StudyAloysius P. GottlobUniversit¨at Kaiserslautern, D-6750 Kaiserslautern, GermanyandMartin HasenbuschCERN, Theory Division, CH-1211 Gen`eve 23, SwitzerlandAbstractWe present the results of a study of the three-dimensional XY -model ona simple cubic lattice using the single cluster updating algorithm com-bined with improved estimators. We have measured the susceptibilityand the correlation length for various couplings in the high temperaturephase on lattices of size up to L = 112.

At the transition temperature westudied the fourth-order cumulant and other cumulant-like quantities onlattices of size up to L = 64. From our numerical data we obtain for thecritical coupling Kc = 0.45420(2), and for the static critical exponentsγ/ν = 1.976(6) and ν = 0.662(7).KL-TH-93/10CERN-TH.6885/93May 1993

I. INTRODUCTIONThree-dimensional classical O(N) vector models are of great interest, both as simplest sta-tistical models with a continuous symmetry and as a lattice version of the scalar quantumfield theory with ( ⃗φ2)2-interaction. In particular, the 3D O(2) model, also called the XY -model, is relevant to the critical behaviour of a number of physical systems, such as thephase transition of superfluid 4He and magnetic systems with planar spin Hamiltonians.Quantitative knowledge of the critical behaviour of the O(N) vector models is mostly basedon the field theoretic renormalization group techniques at dimension D = 3 [1] and theǫ-expansion [2,3].

Very accurate values of the critical exponents are among the most suc-cessful predictions of these approaches. In addition, the analysis of high temperature seriesexpansions [4] provides estimates for the critical temperature of particular lattice models.Monte Carlo simulations have also succeeded in providing detailed information aboutthe critical behaviour of the 3D O(N) vector models, but only in the case of the threedimensional Ising model (N = 1) [5,6] is an accuracy close to that of analytic calculationsreached.The difficulty encountered in Monte Carlo simulations with local updates is the criticalslowing down near a phase transition.Considerable progress has been achieved duringthe last 6 years with the development of efficient non-local Monte Carlo algorithms whichovercome critical slowing down to a large extent [7].In the present paper we extend previous Monte Carlo studies of the 3D XY model [8,9]where cluster algorithms [10,11] were first applied to simulate this model.

These studieswere performed on vector computers, using a moderate amount of CPU-time. Since optimalvectorization cannot be reached for cluster algorithms, we used modern RISC stations for thepresent study.

Using about two months of CPU-time we were able to simulate larger latticesand reached a considerably better statistical accuracy than in the previous studies. Thisaccuracy allowed us to control systematic errors in the estimates for the critical couplingand the critical exponents.This paper is organized as follows.In section 2 we give the definition of the modeland describe the cluster updating algorithm, in section 3 we give our results for the hightemperature phase, while section 4 contains our data obtained in the critical region.

Insection 5 we compare our results with those of previous studies.II. CLUSTER UPDATE MONTE CARLO OF THE 3D XY -MODELWe study the XY -model in three dimensions defined by the partition functionZ =Yi∈ΛZS1dsi exp(KX⟨i,j⟩⃗si · ⃗sj) ,(1)where ⃗si is a two dimensional unit vector, the summation is taken over all nearest neighbourpairs of sites i and j on a simple cubic lattice Λ and K =JkbT is the coupling, or moreprecisely, the reduced inverse temperature.1

For ferromagnetic interactions J > 0, the XY -model has a second-order phase transitionseparating a low temperature phase with non-zero magnetization from a massive disorderedphase at high temperatures. This phase transition can be viewed alternatively as due toBose condensation of spin waves [12] or the unbinding of vortex strings [13,14].A major difficulty encountered in Monte Carlo simulation at second-order phase transi-tions is critical slowing down.

The autocorrelation time τ, which is roughly the time neededto generate statistical independent configurations, grows as τ ∝Lz at criticality, where Lis the linear size of the system and z is the dynamical critical expontent. Random walkarguments indicate that local updates like the Metropolis algorithm result in z = 2, whichis consistent with the numerical finding for the 3D XY -model [8].In the case of O(N) vector models, critical slowing down can be drastically reduced,using cluster algorithms [8,11,9].In the present work we employ the single cluster algorithm which was introduced byWolff[11].

Let us shortly recall the steps of the update. First choose randomly a reflectionaxis in the IR2 plane.

Denote the component of the spin ⃗si which is parallel to this reflectionaxis by s∥i and that which is orthogonal by s⊥i . Then choose randomly a site i of the latticeas a starting point for the cluster C. Visit all neighbour sites j of i.

These sites join thecluster with the probabilityp(i, j) = 1 −exp(−K(s⊥i s⊥j + |s⊥i s⊥j |)). (2)After this is done, visit the neighbours of the new sites in the cluster and add them to thecluster with probability p(i, j) which is given above.

Iterate this step until no new sites enterthe cluster. Now flip the sign of all s⊥contained in the cluster.III.

NUMERICAL RESULTS IN THE HIGH TEMPERATURE PHASE.A. Observables to be measuredLet us first summarize the definitions of the observables that we studied.

The energydensity is given by the two-point correlation function G(xi, xj) = ⟨⃗si⃗sj⟩at distance oneE =13L3X⟨i,j⟩⟨⃗si⃗sj⟩. (3)The specific heat of the system at constant external field is defined by the derivative ofthe energy density with respect to the inverse temperature.

It can be obtained from thefluctuations of the energy H = −P⟨i,j⟩⃗si⃗sjCh = 1L3⟨H2⟩−⟨H⟩2. (4)The magnetic susceptibility χ gives the reaction of the magnetization m = Pi∈Λ ⃗si to anexternal field.

In the high temperature phase one gets2

χ = 1L3⟨m2⟩,(5)since ⟨m⟩= 0.Cluster algorithms enable one to reduce the variance of the expectation values in thehigh temperature phase by using improved estimators [15,16]. The improved estimator ofthe magnetic susceptibility is given byχimp = ⟨2|C| (Xi∈Cs⊥i )2 ⟩,(6)where |C| denotes the number of spins in the cluster C.There are two common definitions of a correlation length ξ.

The exponential correlationlength ξexp is defined via the decay of the two-point correlation function at large distancesξexp =lim|xi−xj|→∞−|xi −xj|log G(xi, xj),(7)which is equal to the inverse mass gap. For the measurement of the exponential correlationlength we consider the correlation functionG(t) ≡⟨O0Ot⟩∝ exp( −tξexp) + exp(−(L −t)ξexp)!,(8)of the translational invariant time slice magnetization Ot = Pi ⃗s(xi, t).The second-moment correlation length is defined byξ2nd = (χ/F) −14 sin2(π/L)!

12,(9)with F = bG(k)||k|=2π/L, where bG(k) = Pj∈Λ⟨exp(ikxj)⃗s0⃗sj⟩is the Fourier transform of thetwo-point correlation function and χ the magnetic susceptibility. For more details see forexample ref.

[17]. The two definitions of the correlation length do not coincide, since in ξexponly the first excited state enters, while in the case of ξ2nd a mixture of the full spectrum istaken into account.

However, near the critical point the two quantities should scale in thesame way. As for the magnetic susceptibility there exist improved estimators for the twodefinitions of the correlation length.

The improved estimator of the two-point correlationfunction is given by⟨⃗si⃗sj⟩imp =* 2|C| δij(C) s⊥i s⊥j+,(10)where δij(C) = 1 if i and j belong to the same cluster C, otherwise δij(C) = 0 [15,16]. For ξ2ndone has to provide a Fimp.

This is given by the Fourier transform of the improved two-pointcorrelation function [16,18]Fimp = bG(k)imp|k|= 2πL,(11)3

withbG(k)imp =* 2|C| Xi∈Cs⊥i cos(kxi)!2+ Xi∈Cs⊥i sin(kxi)!2+. (12)The helicity modulus describes the reaction of the system to a suitable phase twisting field[19].

The lattice definition of the helicity modulus is given byΥµ = 1L3*X⟨i,j⟩sisj(ǫ⟨i,j⟩µ)2+−KL3*(X⟨i,j⟩(s1i s2j + s2i s1j)ǫ⟨i,j⟩µ)2+,(13)where µ is a unit vector in x, y or z direction and ǫ⟨i,j⟩the unit vector connecting the sitesi and j [20].B. Monte Carlo SimulationsIn order to obtain an estimate of the critical coupling Kc and determine static criticalexponents, we have done 15 simulations at couplings from K = 0.4 up to K = 0.452 onlattices of linear size L = 24 up to L = 112.

The simulation parameter and the results of theruns are given in Tables I and II. The statistics is given in terms of N measurements takenevery N0 update steps.

N0 is chosen such that approximately N0 × ⟨C⟩= L3 and hence thewhole lattice is updated once for a measurement. We estimated the statistical errors σA ofexpectation values ⟨A⟩fromσA2 = ⟨A2⟩−⟨A⟩2N/(2τ)(14)and from a binning analysis.

These error estimates were consistent throughout. The statisti-cal error of quantities which contain several expectation values we calculated from Jackknife-blocking [21].1.

Finite-size effectsWe tried to avoid a finite-size scaling analysis. Hence we had to choose our lattices largeenough to ensure that deviations of the values of the observables from the thermodynamiclimit values are negligible.We therefore have measured the energy density E, the specific heat Ch, the helicitymodulus Υ , the exponential correlation length ξexp and the second-moment correlationlength ξ2nd for fixed coupling K = 0.435 and increasing system size L = 4 up to L = 32.The results are summarized in Table III.

The values of the observables obtained for L = 24and L = 32 are consistent within error bars. Furthermore, the values of the helicity modulusΥ for L = 24 and L = 32 are consistent with 0, which is the thermodynamic limit value ofthe helicity modulus in the high temperature phase.

The correlation length at K = 0.435is approximately 4. Hence we conclude, assuming scaling, that the systematical deviationsfrom the thermodynamic limit are smaller than our statistical errors for L/ξ ≥6.

Thiscondition is fulfilled by all the simulation parameters of our runs given in Table I.4

2. Energy density and specific heatThe energy density E shows, as expected, no singular behaviour close to the criticaltemperature.

In the scaling region the specific heat Ch should followCh = Creg + C0Kc −KKc−α,(15)where Creg denotes the regular part of the specific heat and α is the critical exponent of thespecific heat. In order to estimate α we did a four-parameter least-square fit.

However, it wasnot possible to extract meaningful estimates. The best fit to the data leads to Kc = 0.456(2),α = 0.23(13) with χ2/d.o.f.

≈0.93 and the relative errors of the constants are about 100%.If we fix the critical coupling to Kc = 0.45420 (this is our estimate obtained at criticality)the quality of the fit gets worse. Therefore we assumed α = 0 and fitted the data followingCh = Creg + C0 logKc −KKc.

(16)The best three-parameter fit to our data leads to Kc = 0.4543(4), Creg = −0.49(20) andC0 = −1.61(7) with χ2/d.o.f. ≈0.61, where data with ξ > 2.5 are taken into account.This result shows that our data for the specific heat, combined with extended ans¨atze, arecompatible with an α = −0.007(6) obtained from the hyperscaling relation α = 2 −Dν andthe estimate ν = 0.669(2) from resummed perturbation series [1], but have no predictivepower for the exponent α.3.

Magnetic susceptibilityFor comparison we give in Table I the results for the standard and the improved sus-ceptibilities. The statistical error of χimp is about 3.5 to 8 times smaller than the error ofthe standard susceptibility.

But one should remark that the statistical error of the standardestimator depends very much on how often one measures. In the following we only discussthe results obtained with the improved estimator.

In order to estimate the critical couplingand the susceptibility exponent γ we performed a three-parameter least-square fit followingthe scaling lawχ = χ0Kc −KKc−γ. (17)We obtained γ = 1.324(1), Kc = 0.454170(7) and χ0 = 1.009(2) with χ2/d.o.f.

= 0.65, whenall data are taken into account. In order to test the stability of the results we successivelydiscarded data points with small K. The results of these fits are summarized in Table IV.χ2/d.o.f.

remains small and the results for γ, Kc and χ0 are consistent within the errorbars for all data-sets that we used. But the small χ2/d.o.f.

of the fits discussed above, is ofcourse, no proof for the absence of corrections to the scaling. From renormalization groupconsiderations [22] one expects confluent and analytical corrections of the type5

χ(K) = χ0Kc −KKc−γ+ χconf.Kc −KKc−γ+∆1+ χanal.Kc −KKc−γ+1,(18)with ∆1 = ων, where ν is the critical exponent of the correlation length and ω denotesthe correction-to-scaling exponent. We fitted our data according to the scaling law withcorrections.Since a fit with 6 free parameters is hard to stabilize, we fixed the criticalexponents to the values γ = 1.3160(25), ω = 0.780(25) and ν = 0.669(2) which are obtainedfrom resummed perturbation series [1].

Including all the data points in the fit we get Kc =0.454162(9), χ0 = 1.058(7), χconf. = −0.16 = (6) and χanal.

= 0.18(10) with χ2/d.o.f. ≈0.73.

The χ0 which is obtained from the simple scaling fit (17), and that obtained from thefit allowing corrections to the scaling, differ by a larger amount than their statistical errors.This shows that one cannot interpret a small χ2/d.o.f. as the absence of systematic errorsdue to an incomplete fit ansatz.One can also write the scaling relations in terms of the temperature T = 1K.

This leadstoeχ(T) = eχ0T −TcTc−γ(19)andeχ(T) = eχ0T −TcTc−γ+ eχconf.T −TcTc−γ+∆1+ eχanal.T −TcTc−γ+1(20)with corrections. We repeated the analysis as done above.

Taking all 15 data points intoaccount we get for the simple scaling fit a χ2/d.o.f ≈61.2. We again subsequently discardeddata points with small K. A summary of the results is given in Table V. Starting from 5discarded data points the χ2/d.o.f.

becomes approximately 1. But the results obtained forKc, γ and χ0 are not consistent with those obtained from the fit according to eq.

(17).Finally, we performed a four-parameter fit to the scaling relation with corrections andfixed values for the exponents.Taking all data points into account we obtain Kc =0.454163(9), eχ0 = 1.059(7), eχconf. = −0.17(6) and eχanal.

= 1.59(10) with χ2/d.o.f. ≈0.71.The results for Kc, χ0 and χconf of the fits according to the ans¨atze (18) and (20) are consis-tent within the error bars.

The ambiguity between the ansatz with K as variable and thatwith T as variable is covered by the analytic corrections.We conclude that the scaling ansatz (17) fits well if one chooses the coupling K as thevariable. But we also learned that a small χ2/d.o.f.

does not exclude systematic errors, dueto corrections to the scaling, that are larger than the statistical ones. Hence it is hard to givefinal estimates obtained from the simple scaling ansatz that also take systematic errors intoaccount.

From the ansatz with corrections to the scaling we obtain, assuming γ = 1.3160(25)and ∆1 = 0.52182, the results Kc = 0.454162(13), χ0 = 1.058(22) and χconf. = −0.16(11),where the uncertainty of γ is taken into account.We also like to emphasize that the Wegner amplitude χconf.

is negative for the fits thattake corrections into account. This is in agreement with a field-theoretical renormalizationgroup calculation of Esser and Dohm, which predicts the confluent correction-to-scalingamplitude to be negative for a finite cut-off[23].6

4. Correlation lengthWe extracted ξexp from the large distance behaviour of the improved time slice correlationfunction eq.(8).

We therefore considered the effective correlation length, defined byξeff(t) = −ln G(t −1)G(t),(21)where for brevity we have suppressed the contribution due to periodic boundary conditions.As an example, we show in Fig. 1 the results for ξeff(t) obtained on a 1123 lattice atK = 0.452.

A single state dominates the correlation function and a plateau sets in aroundt = ξexp/2 and extends to t = 3ξexp, with no visible degradation due to increasing statisticalerrors at large t. As our final estimate for ξexp we took self-consistently ξeff(t) at the distancet = 2ξexp.In order to calculate ξ2nd we used the improved version of eq.(9). The advantage of thisdefinition is that no fit is needed to obtain the correlation length.

The data of ξexp and ξ2ndare given in Table II. The deviation of ξ2nd from ξexp is about 1% for K = 0.40 and becomessmaller than 0.1% for K ≥0.448.Since the difference of ξexp and ξ2nd is so small, we will discuss only the results of ξexp inthe following.

First we did a three-parameter fit for ξexp following the simple scaling ansatzξ(K) = ξ0Kc −KKc−ν. (22)The results are given in Table VI.

Taking all data into account, the fit has a large χ2/d.o.fof about 9. Starting with three data points with small K being discarded, the χ2/d.o.f.

isclose to 1. But still the critical exponent ν and the critical coupling systematically tend tosmaller values.If we fit the data to the simple scaling ansatz (22), where the coupling is replaced asvariable by the temperature, a similar behaviour is observable.

The results are shown inTable VII. Here one also has to discard three data points to obtain a χ2/d.o.f.

close to 1.But now the estimates of the critical exponent and the critical coupling start at lower valuesand tend to larger ones.This indicates that corrections to the simple scaling ansatz have to be taken into account.Therefore we have fitted all data to the scaling relation with corrections given byξ(K) = ξ0Kc −KKc−ν+ ξconf.Kc −KKc−ν+∆1+ ξanal.Kc −KKc−ν+1,(23)whereas, in the case of the magnetic susceptibility, ∆1 = ων. The four-parameter fit to alldata points with the critical exponents fixed to the resummed perturbation series estimatesleads to Kc = 0.454167(10), ξ0 = 0.498(2), ξconf.

= −0.10(2) and ξanal. = −0.07(4) withχ2/d.o.f.

≈0.63.We also made a four-parameter fit to the scaling relation with corrections where thecoupling is replaced by the temperature. This leads to Kc = 0.454165(10), eξ0 = 0.498(2),eξconf.

= −0.09(2), eξanal. = 0.24(3) with χ2/d.o.f.

≈0.63.7

Taking the uncertainty of ν into account, we leave at Kc = 0.454166(15), ξ0 = 0.498(8)and ξconf. = −0.10(6).In summary, we conclude that systematic deviations from the simple scaling ansatz (22)due to corrections to scaling are important for the analysis of the correlation length datain the coupling range that is accessible to Monte Carlo simulations.Thus it is hard toobtain accurate estimates for the critical exponents and the critical coupling from such anapproach.IV.

NUMERICAL RESULTS AT CRITICALITYOn lattices of the size L = 4, 8, 16, 32 and 64 we performed simulations at K0 = 0.45417which is the estimate for the critical coupling obtained in the previous section. As in thehigh temperature simulation the statistics are given in terms of N measurements takenevery N0 update steps.

We have chosen N0 such that on the average the lattice is updatedapproximately twice for one measurement. The results of the runs are summarized in TableVIII.A.

Phenomenological Renormalization GroupFirst we determined the critical coupling Kc and the critical exponent ν employingBinder’s phenomenological renormalization group method [24]. In addition to the fourth-order cumulant defined on the whole lattice we also studied cumulants defined on subblocksof the lattice.

Therefore let us first introduce blockspinsSB = L1/2(D−2)B1LBDXi∈B⃗si ,(24)where LB is the size of the block and 1/2(D −2) is the canonical dimension of the field. Inparticular we studied the fourth-order cumulantULB = 1 −< (S2B)2 >3 < S2B >2(25)for LB = L, L/2 and a nearest neighbour interaction on subblocksNN = < SB1SB2 >< S2B >(26)for LB = L/2.For the extrapolation to couplings K other than the simulation coupling K0, we usedthe reweighting formula⟨A⟩(K) =Pi Ai exp((−K + K0)Hi)Pi exp((−K + K0)Hi) ,(27)8

where i labels the configurations generated according to the Boltzmann-weight at K0. Wecomputed the statistical errors from Jackknife binning on the final result of the extrapolatedcumulants.

The extrapolation only gives good results within a small neighbourhood of thesimulation coupling K0. This range shrinks with increasing volume of the lattice.For sufficiently large LB the cumulants have a non-trival fixed point at the critical cou-pling.

When one considers the cumulants as a function of the coupling, the crossings of thecurves for different L provide an estimate for the critical coupling Kc. As an example weshow in Fig.

2 the fourth-order cumulant in a neighbourhood of Kc. The figure shows thatthe crossings of the cumulant are well covered by the extrapolation (27).

The error bars ofUL with L = 64 blow up for |K −K0| > 0.001 while |Kcross −K0| = 0.00003 for L = 32and L = 64. The results for the crossings are summarized in Table IX.

The given errorsare taken from the size of the crossings of the error bars. The convergence of the crossingcoupling Kcross towards Kc should followKcross(L) = Kc (1 + const.L−(ω+ 1ν )),(28)where ω is the correction to scaling exponent [24].

Our data for the crossings of the cumulantsdid not allow us to perform a two-parameter fit, keeping the exponents fixed, following theabove formula. Within the statistical errors the results of the crossings of the fourth-ordercumulants on L = 8 and L = 16 up to L = 32 and L = 64 are consistent.

The convergence ofthe crossings of NN towards Kc seems to be slower than that of the fourth-order cumulants,but it is interesting to note that the Kcross for the fourth-order cumulant and NN comefrom different sides with increasing L.This is shown in Fig. 3, where the estimates ofKcross versus the lattice size L are plotted.

Our final estimate for the critical coupling isKc = 0.45420(2) obtained from the L = 32 and L = 64 crossing of the fourth-order cumulanton the full lattice. Taking into account the fast convergence of the crossings towards Kc,that is predicted by (28), we conclude that the systematic error of our estimate for Kc issmaller than the given statistical error.At the critical coupling Kc the cumulants converge with increasing lattice size L to auniversal fixed point.

The convergence rate is given by [24]UL(Kc) = U∞(1 + const.L−ω) . (29)The results for the cumulants at K = 0.45420 , which is our estimate of critical coupling,are given in Table X.

The data did not allow us to perform a two parameter fit with ω beingfixed. Hence we take the value UL(Kc) = 0.589(2) form L = 64 as our final estimate for thefixed point of the fourth-order cumulant on the full lattice, where we now have taken intoaccount the uncertainty of the estimated critical coupling.We extracted the critical exponent ν of the correlation length from the L dependence ofthe slope of the fourth-order cumulant at criticality [24].

According to Binder, the scalingrelation for the slope of the fourth-order cumulant is given by∂U(L, K)∂KKc∝L1/ν. (30)9

We evaluated the slopes of the observables A entering the cumulant U according to∂⟨A⟩∂K = ⟨AH⟩−⟨A⟩⟨H⟩,(31)where A is an observable and H is the energy. The statistical errors are calculated froma Jackknife analysis for the value of the slope.

First we estimated the exponent ν fromdifferent lattices viaν =ln (L2) −ln (L1)ln ∂A(L2, K)∂KKc!−ln ∂A(L1, K)∂KKc!. (32)The results are given in Table XI.

The estimates for ν stemming from UL and UL/2 are stablewith increasing L and consistent with each other for L ≥16. Therefore we performed a fitaccording to eq.

(30) with UL from lattices of the size L = 16 up to L = 64. We consider theresult ν = 0.662(7) as our final estimate for the critical exponent of the correlation length.B.

Magnetic SusceptibilityIn order to estimate the ratio γ/ν of the critical exponents we studied the scaling be-haviour of the magnetic susceptibility defined on the full lattice and on subblocks. Thedependence of the susceptibility on the lattice size is given byχ ∝Lγ/ν(33)at the critical coupling.

We have estimated γ/ν from pairs of lattices with size L1, L2. Theratio is then given byγν = ln(χ(L1, Kc)) −ln(χ(L2, Kc))ln(L1) −ln(L2).

(34)The second column of Table XII shows the estimates of the ratio obtained from the sus-ceptibility defined on the full lattice, while the third column shows the estimates obtainedfrom the blockspin-susceptibility with subblocks of the size L/2. The estimates for γ/ν ob-tained from the subblocks monotonically increase with increasing lattice size L, while thoseobtained from the full lattice decrease.

The results obtained from the full lattice for L ≥16and the result from the subblocks of the largest lattices are consistent within error bars.Hence we take γ/ν = 1.976(6) as our final result, where statistical as well as systematicerrors should be covered. Using the scaling relation η = 2 −γν we obtain for the anomalousdimension η = 0.024(6).C.

Hyperscaling and specific heatIn ref. [24] a dangerous irrelevant scaling field u is proposed as explanation for a possibleviolation of hyperscaling.Dangerous means that the scaling function of the correlation10

length vanishes with some power q of the vanishing irrelevant scaling field.Hence thecorrelation length should scale asξ ∝L1+qyu(35)at the critical point.Remember that yu is negative for an irrelevant scaling field.AtKc = 0.45420, which we obtained from the analysis of the fourth-order cumulant, we havefitted ξ2nd to this relation. The reweighted estimates of ξ2nd are shown in Table XIII.

Takinglattices of size L = 16 up to L = 64 into account we estimate qyu = 0.007(2) with χ2/d.o.f. ≈0.3 and only statistical errors considered.

This indicates that there is no or only very smallhyperscaling violation due to a dangerous irrelevant field.At criticality the specific heat should scale asCh(L) = Creg + const. Lαν ,(36)where Creg denotes the regular part of the specific heat and α is the critical exponent ofthe specific heat.

Using the critical exponent ν = 0.662(7) obtained from the analysis abovethe hyperscaling relation α = 2 −Dν gives α = 0.014(21). We also tried to estimate α viaa three-parameter fit, following the finite-size scaling relation.

However, we are not able togive a stable estimate for α.D. Helicity modulusThe 3D XY model is assumed to share the same universality class as an interacting Bosefluid, and the helicity modulus Υ should be proportional to the superfluid density ̺s of theBose fluid [19].

Near the critical coupling the superfluid density, resp. the helicity modulusshould scale as̺s ∝Υ ∝|K −Kc|v ,(37)with v the critical exponent of the superfluid density.

Assuming hyperscaling the Josephsonrelation reads v = (D −2)ν [19]. Hence the productΥ · L = const(38)should stay constant at the critical point in 3D.

To check this prediction we have measuredthe helicity modulus Υ on lattices of size L = 4 to L = 32. The estimator of ΥL becomesnoisy with increasing lattice size.

We tried to overcome this problem by measuring moreoften, which did not remove the problem completely. Hence we skipped the measurement ofΥ for L = 64.

The results, shown in Table XIII, indicate that the above relation holds.E. Performance of the AlgorithmThe efficiency of a stochastic algorithm is characterized by the autocorrelation time11

τ = 12∞Xt=−∞ρ(t) ,(39)where the normalized autocorrelation function ρ(t) of an observable A is given byρ(t) = ⟨Ai · Ai+t⟩−⟨A⟩2⟨A2⟩−⟨A⟩2. (40)We calculated the integrated autocorrelation times τ with a self-consistent truncation win-dow of width 6τ for the energy density E and the magnetic susceptibility χ for lattices withL = 4 up to L = 64 at the coupling K = 0.45417.

In Fig. 4 we show a log-log plot of theintegrated autocorrelation times τ of the energy density E and the magnetic susceptibilityχ versus the lattice size L given in units of the average number of clusters that is neededto cover the volume of the lattice.

Our estimates for the critical dynamical exponents arezE = 0.21(1) and zχ = 0.07(1) taking only statistical errors into account. These results areconsistent with those of Janke [9].Finally let us briefly comment on the CPU time: 160 single cluster updates of the 643lattice at the coupling K = 0.45417 plus one measurement of the observables took on average26 sec CPU time on a IBM RISC 6000-550 workstation.

All our MC simulations of the 3DXY model together took about two months of CPU-time.V. COMPARISON OF OUR RESULTS WITH PREVIOUS STUDIESIn Table XIV we display estimates of critical properties of the 3D XY -model obtained byvarious methods.

Our estimates of Kc from the scaling fit to the high temperature data andfrom the penomenological RG approach are consistent within 2 standard deviations. Butonly for the result from the phenomenological RG approach are the systematical errors fullyunder control.

Our error of Kc is about 4 times smaller than that of previous MC studies[8,9], and also about 4 times smaller than that obtained recently [26] from the analysis ofa 14th order high temperature series expansion [28]. Recently Butera et al.

[27] extendedthe high temperature series expansion for the sc lattice to the order 17. Their value for thecritical coupling is by three times their error estimate smaller than our value.The error of γ obtained from a fit to the simple scaling ansatz is about 5 times smallerthan those of previous MC studies [8,9]; however, the systematical errors are not undercontrol.

The value of γ is, within two standard deviations, consistent with the estimate ofRef. [9].

Our estimate of γ is consistent with the value obtained from the high temperatureseries expansion [25,26,27] and, within two standard deviations, consistent with the value ofthe ǫ-expansion [3] while the very accurate estimate from the resummed perturbation series[1] is smaller than our estimate by three times the quoted error.Our estimate for ν is consistent within error bars with all other estimates we quote inTable XIV. Our quoted error bars are 3.5 times larger than that of ref.

[9]. Janke used finitedifferences to determine the slope of the cumulant [29], while we used fluctuations at a singletemperature (31).

Furthermore the smallest lattice size L = 16 included in our fit is chosen12

to be rather conservative. The most accurate number for ν stems from the measurement ofthe superfluid fraction of 4He [30].In this work we give for the first time an accurate direct MC estimate for the exponentη, the anomalous dimension of the field.

The uncertainty of the estimate is comparable withthose obtained with field theoretical methods. Our value of the exponent η is consistentwith the estimates from the high temperature series expansion [25,26] and with that of theresummed perturbation series [1], but is smaller than the ǫ-expansion [3] result by more thantwice our error estimate.Our result for the critical fourth-order cumulant, is consistent with previous MC results[8,9].

But the value obtained from ǫ-expansion [31] is offby about 20 times our error estimatethat also takes into account systematic errors. Furthermore we provide estimates for thecritical fourth order cumulant on subblocks and a nearest neighbour blockspin product NN.These numbers might be useful in testing other models sharing the XY universality class.VI.

CONCLUSIONSThe application of the single cluster algorithm [11], which is almost free of critical slowingdown for the 3D XY model, and the extensive use of modern RISC workstations allowedus to increase the statistics as well as the studied lattices sizes considerably compared withprevious MC simulations [8,9]. In the high temperature phase of the model we measuredcorrelation length up to 17.58 with an accuracy of about 0.1%.

But the analysis of our datafor the correlation length and the magnetic susceptibility showed that it is hard to controlsystematic errors due to confluent and analytic corrections.It seems to be much easierto fight the systematic errors in the phenomenological RG approach. Analytic correctionsare absent at the critical point and corrections to the scaling are less harmful, since therelevant length scale at criticality is the lattice size, which can be chosen much larger thanthe correlation length in the thermodynamic limit of the high temperature phase.

From thecrossings of the fourth-order cumulant we obtain Kc = 0.45420(2), which reduces the errorby a factor of about 4 compared with previous MC studies [8,9]. Further improvements of theaccuracy of the estimates of the critical coupling and critical exponents seem to be reachableby just increasing the statistics, while keeping the present lattice sizes.

The accurate valuesobtained for critical cumulants could be very useful for testing whether other models sharethe XY universality class. Here of course a proper block-spin definition is essential.ACKNOWLEDGMENTSThe numerical simulations were performed on an IBM RISC 6000 cluster of the Re-gionales Hochschulrechenzentrum Kaiserslautern (RHRK).

The work was supported in partby Deutsche Forschungsgemeinschaft (DFG) under grant Me 567/5-3. It is a pleasure tothank S. Meyer for discussions.

We would like to thank W. Janke for communicating to ussome unpublished details of his study.13

REFERENCES[1] J.C. Le Guillou and J. Zinn-Justin, Phys.

Rev. B 21, 3976 (1980).

[2] K.G. Wilson and M.E.

Fisher, Phys. Rev.

Lett. 28, 240 (1972).

[3] J.C. Le Guillou and J. Zinn-Justin, J. Phys.

Lett. (Paris) 46, L137 (1985).

[4] Phase Transitions and Critical Phenomena vol. 3 , C. Domb and M.S.

Green eds. (Aca-demic Press, London, 1974).

[5] A.M. Ferrenberg and D.P. Landau, Phys.

Rev. B 44, 5081 (1991).

[6] C.F. Baillie, R. Gupta, K.A.

Hawick and G.S. Pawley, Phys.

Rev. B 45, 10438 (1992).

[7] For reviews see, e.g.:A.D. Sokal, Nucl.

Phys. B (Proc.

Suppl.) 20, 55 (1990);R.H. Swendsen, J.S.

Wang and A.M. Ferrenberg, in The Monte Carlo Method in Con-densed Matter Physics, ed. K. Binder (Spinger-Verlag, Berlin-Heidelberg, 1992).

[8] M. Hasenbusch and S. Meyer, Phys. Lett.

B 241, 238 (1990). [9] W. Janke, Phys.

Lett. A 148, 306 (1990).

[10] R. H. Swendsen and J.-S. Wang, Phys. Rev.

Letters 58, 86 (1987). [11] U. Wolff, Phys.

Rev. Lett.

62, 361 (1989);U. Wolff, Nucl. Phys.

B322, 759 (1989). [12] J. Fr¨ohlich, B. Simon and T. Spencer, Commun.

Math. Phys.

50, 79 (1976). [13] J. Fr¨ohlich and T. Spencer, Commun.

Math. Phys.

83, 411 (1982). [14] G. Kohring, R.E.

Shrock and P. Wills, Phys. Rev.

Lett. 57, 1358 (1986).

[15] U. Wolff, Nucl. Phys.

B 334, 581 (1990). [16] M. Hasenbusch, Nucl.

Phys. B333, 581 (1990).

[17] F. Cooper, B. Freedman and D. Preston, Nucl. Phys B 210, 210 (1982);R.G.

Edwards, S.J. Ferreira, J. Goodman and A.D. Sokal, Nucl.

Phys. B 380, 621 (1992).

[18] C. Holm and W. Janke, HLRZ preprint 77/92 J¨ulich (1992). [19] M.E.

Fisher, M.N. Barber and D. Jasnow, Phys.

Rev. A 8, 1111 (1973).

[20] Y.-H. Li and S. Teitel, Phys. Rev.

B 40, 9122 (1989). [21] R.G.

Miller, Biometrica 61, 1 (1974);B. Efron, The Jackknife, the Bootstrap and other Resampling Plans (SIAM, Philadelphia,PA, 1982). [22] F.J. Wegner, Phys.

Rev. B5, 4529 (1972)[23] A. Esser and V. Dohm, in The 18th IUPAP International Conference on Statistical Physicsp.159, ed.

W. Loose (TU Berlin, Germany, 1992). [24] K. Binder, Phys.Rev.Lett 47, 693 (1981) ;K. Binder, Z. Phys.

B - Condensed Matter 43, 119 (1981);K. Binder, M. Nauenberg, V. Privman and A.P. Young, Phys.

Rev. B 31, 1498 (1985).

[25] M. Ferer, M.A. Moore and W. Wortis, Phys.

Rev. B 8, 5205 (1973).

[26] J. Adler, C. Holm and W. Janke, preprint FUB-HEP 18/92 , HLRZ 76/92 Berlin/J¨ulich(1992). [27] P. Butera, M. Comi and A.J.

Guttmann, preprint IFUM 444/FT Milano (1992). [28] M. L¨uscher and P. Weisz, Nucl.

Phys. B 300, 325 (1988).

[29] W. Janke, private communication.14

[30] L.S. Goldner and G. Ahlers, Phys.

Rev. B 45, 13129 (1992).

[31] E. Br´ezin and J. Zinn-Justin, Nucl. Phys.

B 257, 867 (1985).15

FIGURESFIG. 1.The effective correlation length ξeff(t) as a function of separation for couplingK = 0.452 on a lattice of size L = 112.FIG.

2. Reweighting plot of the Binder cumulant UL of the full lattice from the simulation atK = 0.45417.

The dashed lines give the statistical errors obtained by a binning analysis.FIG. 3.

Plot of the convergence of the critical coupling obtained by the cumulant crossingmethod. Because of the small statistical errors one is able to see systematic convergence of thecritical coupling.FIG.

4. Intergated autocorrelation times τ of the energy density E and the magnetic suscep-tibility χ versus the lattice size L. The dynamical critical exponent is given by the slopes of thefits.16

TABLESTABLE I. Results of the energy density E, the specific heat at constant external field Ch, theimpoved susceptibility χimp and the standard susceptibility χ obtained from the simulations in thehigh temperature phase.

The parameters of the runs are given in terms of the simulation couplingK, the linear size of the system L, and the statistics, with N the number of measurements takenevery N0 update steps.KLNN0EChχimpχ0.4002420k0.9k0.24697(6)3.17(3)16.848(14)16.70(12)0.4103210k2k0.25779(5)3.39(6)22.108(19)22.24(32)0.4203210k1.6k0.26970(5)3.65(6)30.994(29)31.18(32)0.4253210k1.2k0.27605(6)3.96(7)38.21(5)38.47(38)0.4303220k0.8k0.28286(4)4.22(4)49.12(6)48.87(36)0.4353220k0.6k0.29023(7)4.57(8)66.57(15)66.6(7)0.4373220k0.4k0.29329(5)4.72(6)77.11(16)76.9(5)0.4403220k0.4k0.29818(5)5.11(6)99.50(20)99.2(7)0.4434821k1k0.30336(3)5.39(7)136.41(20)135.8(9)0.4454840k0.5k0.30705(3)5.84(6)177.04(30)176.2(9)0.4486425k1k0.31303(2)6.43(8)299.14(53)301.0(2.0)0.4496425k1k0.315214(18)6.61(7)377.75(64)378.9(2.5)0.4509612k3k0.317497(15)7.02(11)503.85(85)506.5(4.6)0.4519612k2k0.319854(13)7.33(10)722.1(1.3)724.6(6.5)0.45211212k2k0.322446(12)8.10(14)1193.0(3.0)1201.0(10.0)17

TABLE II. Results for the correlation lengths ξexp and ξ2nd obtained from the simulations inthe high temperature phase.

ξ2nd/ξexp gives the ratio between the values of the two correlationlengths.KLξ2ndξexpξ2nd/ξexp0.400241.876(1)1.898(2)0.98840.410322.182(1)2.202(2)0.99090.420322.624(2)2.639(2)0.99430.425322.938(2)2.953(2)0.99490.430323.361(2)3.375(3)0.99590.435323.947(5)3.959(5)0.99690.437324.262(5)4.273(6)0.99740.440324.875(5)4.885(6)0.99790.443485.746(5)5.756(6)0.99830.445486.582(6)6.593(7)0.99830.448648.638(9)8.645(10)0.99910.449649.738(9)9.747(10)0.99910.4509611.288(10)11.295(10)0.99940.4519613.587(16)13.594(16)0.99950.45211217.570(19)17.580(20)0.9994TABLE III. Results of the energy density E, the specific heat Ch, the helicity modulus Υ, theimproved susceptibility χimp, the correlation lengths ξexp and ξ2nd and the ratio L/ξ2nd obtainedfrom simulations at K = 0.435 with linear system size L. The statistics is given in terms of Nmeasurements taken every N0 update steps.LNN0EChΥχimpξ2ndξexpL/ξ2nd425k0.1k0.3621(26)6.43(13)0.2180(44)15.83(21)2.015(19)2.08(2)1.98820k0.2k0.3054(26)6.56(13)0.0608(16)43.00(07)3.228(3)3.272(4)2.481620k0.2k0.29114(9)4.92(6)0.0056(17)64.31(16)3.886(6)3.901(6)4.122020k0.3k0.29049(7)4.71(5)0.0031(16)66.21(14)3.938(4)3.949(5)5.082420k0.3k0.29024(7)4.65(5)−0.0001(17)66.55(15)3.946(6)3.959(6)6.083210k0.6k0.29013(7)4.52(9)0.0009(25)66.64(15)3.949(5)3.960(5)8.1018

TABLE IV. Estimates for the critical coupling Kc, the static critical exponent γ and theamplitude χ0 obtained from a fit of the improved susceptibility χimp to eq.(17).

χ2/d.o.f gives thequality of the fit. # denotes the number of discarded data points at small couplings.#Kcγχ0χ2/d.o.f.00.454170(7)1.3241(10)1.0090(24)0.6510.454168(8)1.3238(12)1.0099(30)0.6920.454175(9)1.3252(14)1.0057(40)0.4930.454173(10)1.3248(18)1.0069(53)0.5340.454170(11)1.3239(22)1.0101(66)0.5250.454179(14)1.3264(31)1.0016(100)0.4160.454176(15)1.3256(35)1.0043(117)0.4470.454174(17)1.3251(43)1.0060(145)0.5280.454174(19)1.3250(52)1.0063(180)0.6590.454180(24)1.3272(74)0.9979(267)0.81100.454197(42)1.3337(156)0.9736(576)1.11110.454208(58)1.3384(240)0.9557(885)2.14TABLE V. Estimates for the critical coupling Kc, the static critical exponent γ and the am-plitude χ0 obtained from a fit of the improved susceptibility to eq.(19).

χ2/d.o.f gives the qualityof the fit. # denotes the number of discarded data points at small couplings.#Kcγbχ0χ2/d.o.f.00.453871(6)1.2351(8)1.3972(27)77.110.453932(7)1.2471(10)1.3500(34)40.020.453995(8)1.2604(13)1.2970(45)12.830.454028(9)1.2683(16)1.2648(59)6.9840.454040(10)1.2733(20)1.2440(74)5.2750.454091(13)1.2849(29)1.1959(111)1.3960.454098(14)1.2871(33)1.1870(128)1.3170.454110(16)1.2907(40)1.1720(158)1.0880.454120(18)1.2938(49)1.1587(196)1.0290.454135(23)1.2991(70)1.1359(289)0.99100.454166(42)1.3112(150)1.0845(618)1.06110.454182(59)1.3177(231)1.0570(945)1.9819

TABLE VI. Estimates for the critical coupling Kc, the static critical exponent ν and theamplitude ξ0 obtained from a fit of the exponential correlation length ξexp to eq.(22).

χ2/d.o.fgives the quality of the fit. # denotes the number of discarded data points at small couplings.#Kcνξ0χ2/d.o.f.00.454325(9)0.7029(7)0.4294(8)9.210.454301(9)0.7003(8)0.4327(9)4.920.454286(10)0.6985(9)0.4351(11)3.630.454269(11)0.6964(10)0.4381(13)2.040.454247(12)0.6933(13)0.4426(18)0.6650.454243(14)0.6927(18)0.4436(25)0.7160.454235(16)0.6914(21)0.4457(31)0.5970.454223(18)0.6895(24)0.4487(37)0.3180.454216(20)0.6882(30)0.4509(48)0.2590.454210(24)0.6870(39)0.4529(64)0.26100.454208(46)0.6866(90)0.4537(157)0.39110.454218(66)0.6890(140)0.4493(246)0.74TABLE VII.

Estimates for the critical coupling Kc, the static critical exponent ν and theamplitude eξ0 obtained from a fit of the exponential correlation length ξexp to eq. (22), where thecoupling is replaced by the inverse temperature.

χ2/d.o.f gives the quality of the fit. # denotesthe number of discarded data points at small couplings.#Kcνeξ0χ2/d.o.f.00.454079(8)0.6612(6)0.5024(8)8.510.454098(8)0.6632(7)0.4994(9)5.320.454118(9)0.6656(8)0.4957(11)1.730.454125(10)0.6665(9)0.4944(13)1.540.454134(11)0.6677(12)0.4923(19)1.450.454157(14)0.6711(16)0.4866(26)0.1860.454160(15)0.6717(19)0.4857(31)0.1670.454161(17)0.6717(23)0.4856(38)0.1980.454164(20)0.6723(29)0.4846(49)0.2290.454165(24)0.6724(38)0.4844(66)0.29100.454177(45)0.6750(90)0.4796(160)0.38110.454192(65)0.6783(136)0.4734(250)0.6520

TABLE VIII. Results of the energy density E, the specific heat Ch, the susceptibility χ and thesecond moment correlation length ξ2nd obtained from simulations at the fixed coupling K = 0.45417near the final estimate of the critical coupling.

τ denotes the integrated autocorrelation time of thespecified observable, given in units of the average number of clusters needed to cover the volumeof the lattice. The statistics is given in terms of N measurements taken every N0 update steps.LNN0EτEChχτχξ2nd4100k100.40440(44)2.0(1)6.561(27)19.095(34)1.84(5)2.3104(37)895k200.35585(20)2.4(1)8.890(39)77.80(15)1.97(5)4.6852(65)16100k400.338945(7)2.6(1)10.757(66)309.95(60)1.96(7)9.447(15)3283k800.332815(3)3.1(1)12.520(73)1216.0(2.7)2.11(5)18.922(38)6472k1600.330628(2)3.7(1)14.35(11)4732(12)2.32(7)37.793(77)TABLE IX.

Estimates for Kc(L) obtained via Binder’s cumulant crossing technique of thereweighted fourth-order cumulants UL and UL/2 and nearest neighbour observable NN. L1 −L2gives the pair of linear lattice sizes which determine the intersection point.Kc(L)L1 −L2ULUL/2NN4 −80.4565(4)0.4617(3)0.4378(4)8 −160.4544(2)0.45457(8)0.4517(1)16 −320.45423(5)0.45424(4)0.45393(4)32 −640.45420(2)0.45421(2)0.45415(2)TABLE X.

Results for the fourth-order cumulants UL, UL/2 and the nearest neighbour ob-servable NN at K = 0.45420 obtained with the reweighting technique from the simulations atK = 0.45417. The errors are obtained by a Jackknife-blocking procedure.LULUL/2NN40.59640(42)0.56860(25)0.70557(61)80.59134(42)0.55270(31)0.77439(45)160.58966(43)0.55040(32)0.79925(44)320.58907(50)0.54974(37)0.80640(48)640.58909(44)0.54925(33)0.80963(49)21

TABLE XI. Estimates for the static critical exponent ν obtained using eq.

(32), where A isreplaced by the fourth-order cumulants UL and UL/2 and the nearest neighbour observable NNwith the critical coupling is set to Kc = 0.45420, the final estimate of the critical coupling.latticeνL1 −L2ULUL/2NN4 −80.6496(93)0.5807(50)0.8443(86)8 −160.6799(111)0.6576(74)0.7519(76)16 −320.6649(126)0.6694(84)0.6977(68)32 −640.6584(154)0.6565(103)0.6779(83)TABLE XII. Estimates for the ratio of the static critical exponents γ/ν obtained using eq.

(34).The first column gives the results of the ratio obtained from the susceptibility of the full latticewhile the second column is obtained from the susceptibility of the subblocks. L1 −L2 gives thepair of lattices, which is used to calculate the ratio of the exponents.latticeγ/νL1 −L2full latticesubblocks4 −82.027(4)1.898(2)8 −161.996(4)1.954(3)16 −321.978(4)1.966(3)32 −641.979(5)1.974(4)TABLE XIII.

Expectation values of the specific heat Ch, the second moment correlation lengthξ2nd, and the product of the helicity modulus Υ times the linear size of the system that arereweighted to the final estimate of the critical coupling Kc = 0.45420. The errors are calulatedfrom a Jackknife analysis.LChξ2ndΥ · L46.561(27)2.3112(34)1.090(2)88.877(33)4.6856(68)1.091(4)1610.704(72)9.4639(142)1.12(1)3212.564(63)19.003(34)1.13(2)6414.406(102)38.25(69)-22

TABLE XIV. Comparison of critical properties determined from various methods.

The resultsgiven for the simulations of the model in the high temperature phase are obtained from fits ac-cording to the simple scaling ansatz with the coupling K as parameter. For γ and ν from Ref.

[25]we took the estimates of the fcc lattice, since the errors are smaller than those obtained from thesc lattice.MethodRef.KcγνηULPhenomenological RGthis work0.45420(2)-0.662(7)0.024(6) 0.589(2)High temperature MCthis work 0.454170(7)1.324(1)---Phenomenological RG[9]0.4542(1)-0.670(2)≈0.020.586(1)High temperature MC[9]0.45408(8)1.316(5)---Phenomenological RG[8]--≈0.67-0.590(5)High temperature MC[8]0.45421(8)1.327(8)---ǫ-expansion[3]-1.315(7)0.671(5)0.040(3)-ǫ-expansion[31]----0.552Resummed perturbation series[1]-1.3160(25)0.669(2)0.033(4)-High temperature series[25]0.4539(12)1.323(15)0.670(7)0.028(5)-High temperature series[26]0.45414(7)1.3250.6730.030-High temperature series[27]0.45406(5)1.315(9)0.68(1)--Experiment 4He[30]--0.6705(6)--23

---

출처: arXiv:9305.020 • 원문 보기