We present a simple and powerful method for extrapolating finite-volume Monte Carlo data to infinite volume, based on finite-size-scaling theory. We discuss carefully its systematic and statistical errors, and we illustrate it using three examples: the two-dimensional three-state Potts antiferromagnet on the square lattice, and the two-dimensional $O(3)$ and $O(\infty)$ $\sigma$-models. In favorable cases it is possible to obtain reliable extrapolations (errors of a few percent) even when the correlation length is 1000 times larger than the lattice

Caracciolo, S

Edwards, RG

Ferreira, SJ

Pelissetto, A

Sokal, AD

UCL Discovery

VOLUME 74, NUMBER 15 PHYSICAL REVIEW LETTERS 10 APRIL 1995Extrapolating Monte Carlo Simulations to Infinite Volume: Finite-Size Scaling at g/L )) 1Sergio Caracciolo*Dipartimento di Fisica, Universita di Lecce and INFRA —Sezione di Lecce, I-73100 Lecce, ItalyRobert G. Edwards~Supercomputer Computations Research Institute, Florida State University, Tallahassee, Florida 32306Sabino Jose Ferreira~Minas Gerais 30161, BrasilAndrea Pelissetto~ and Alan D. SokalDepartment of Physics, ¹wYork University, 4 Washington Place, New York, New York l0003(Received 9 September 1994)We present a simple and powerful method for extrapolating finite-volume Monte Carlo data to infinitevolume, based on finite-size-scaling theory. We discuss carefully its systematic and statistical errors,and we illustrate it using three examples: the two-dimensional three-state Potts antiferromagnet on thesquare lattice, and the two-dimensional O(3) and O(~) cr models. In favorable cases it is possible toobtain reliable extrapolations (errors of a few percent) even when the correlation length is 1000 timeslarger than the lattice.PACS numbers: 02.70.Lq, 05.70.Jk, 11.15.Ha, 64.60.FrNo matter how powerful computers become, physicistswill always want to study problems that are too difficultfor the computers at hand. For example, in statisticalmechanics and quantum field theory, physicists want topush to ever larger correlation lengths g. But MonteCarlo simulations must perforce be carried out on latticesof finite linear size L (limited by computer memory andspeed); the data are then extrapolated to the infinite-volume limit L = ~. Obviously this extrapolation—which is based on the theory of finite-size scaling (FSS)[1]—is feasible in practice only if g/L is not too large.But how largerIn this Letter we present a simple and powerful methodfor performing the extrapolation to L = ~, and discusscarefully its systematic and statistical errors. We illustratethe method using three examples: the two-dimensionalthree-state Potts antiferromagnet on the square lattice [2],and the two-dimensional O(3) and O(~) o. models [3,4].We found —much to our surprise —that in favorablecases it is possible to obtain reliable extrapolations (errorsof a few percent) at g/L as large as 10—1000.Consider, for a start, a model controlled by arenormalization-group (RG) fixed point having one rele-vant operator. Let us work on a periodic lattice of linearsize L Let se(P, L) be a suit. ably defined finite-volumecorrelation length [5], and let 6 be any long-distance ob-servable (e.g. , the correlation length or the susceptibility).Then finite-size-scaling theory [1]predicts that [7]6,I.6(P, ) = f (F(P, )/L) + O(6, L "),where fthm is a universal function and co is a correction-to-scaling exponent. It follows that if s is any fixed scalefactor (usually we take s = 2) then6(P, L) = Fr (g(P, L)/L) + O(g ",I. "),where F@ can be expressed in terms of f@ and f~Our method proceeds as follows [8]. Make MonteCarlo runs at numerous pairs (P, L) and (P, sL). Plot6(p, sL)/6(P, L) versus $(P, L)/L, using those pointssatisfying both s (P, L) ~ some value g;„and L ~ somevalue L;„. If all these points fall with good accuracy on asingle curve —thus verifying the ansatz (2) for P ~ $;„,L ~ L;„—choose a smooth fitting function Fg. Then,using the functions Ft and Fto, extrapolate the pair (g, 6)successively from L ~ sL ~ s2L, ~ . .We have chosen to use functions Fg of the form(2)Fg(x) = I +a~e 't +a2e '+ . . +a e "'. (3)This form is partially motivated by theory, which tells usthat F(x) I exponentially fast as x ~ 0 [11]. Typicallya fit of order 3 ~ n ~ 12 is sufficient; we increase nuntil the y of the fit becomes essentially constant. Theresulting y value provides a check on the systematicerrors arising from corrections to scaling and/or from theinadequacies of the form (3).The statistical error on the extrapolated value of6 (P) —= 6(P, ~) comes from three sources: (i) erroron 6(P, L), which gets multiplicatively propagated to6; (ii) error on se(P, L), which affects the argumentx = $(P, L)/L of the scaling functions Fe and F@., and(iii) statistical error in our estimate of the coefficientsa~, . . . , a„ in Fe and F@. The errors of types (i) and(ii) depend on the statistics available at the single point(P, L), while the error of type (iii) depends on the0031-9007/95/74(15)/2969(4)$06. 00 Oc 1995 The American Physical Society 2969VOLUME 74, NUMBER 15 PHYSICAL REVIEW LETTERS 10 APRIL 1995statistics in the whole set of runs. Errors (i) + (ii) [(i)+ (ii) + (iii)] can be quantified by performing a MonteCarlo experiment in which the input data at (P, L) [thewhole set of input data] are varied randomly within theirerror bars and then extrapolated [12—14].The discrepancies between the extrapolated values fromdifferent lattice sizes at the same P —to the extent thatthese exceed the estimated statistical errors —indicate thepresence of systematic errors and thus the necessity ofincreasing L;„and/or g;„and/or nA figure of (de)merit of the method is the relativevariance on the extrapolated value 6 (P), multiplied bythe computer time needed to obtain it [15]. We expectthis relative variance time -product (RVTP) [for errors (i)+ (ii) only] to scale asRVTP(P, L) = 6-(/3)"""' GG(r-(p)/L), (4)2.00, 0 0.2I I I I I I I0.8where d is the spatial dimension and z;„,@ is the dy-namic critical exponent of the Monte Carlo algorithm be-ing used; here Gg is a combination of several static anddynamic finite-size-scaling functions, and depends bothon the observable 6 and on the algorithm but not onthe scale factor s. As g /L tends to zero, we expectGo to diverge as (g /L) " (it is wasteful to use a lat-tice L &) g ). As $ /L tends to infinity, we expectGg —(g /L)" [16], but the power p can be either positive or negative. If p ) 0, there is an optimum value ofg /L; this determines the best lattice size at which to per-form runs for a given p. If p ( 0, it is most efficientto use the smallest lattice size for which the correctionsto scaling are negligible compared to the statistical errors.(This neglects errors of type (iii) [17]; the optimizationbecomes much more complicated if they are included. )Our first example [2] is the two-dimensional three-state Potts antiferromagnet on the square lattice, whichis believed to have a critical point at P = ~ [18].We used the Wang-Swendsen-Kotecky cluster algorithm[19], which appears to have no critical slowing down( 5 uniformly in P and L) [2]. We ran on lat-tices L = 32, 64, 128, 256, 512, 1024,1536 at 153 differentpairs (P, L) in the range 5 ( g ( 20000. Each run wasbetween 2 X 10 and 2.2 ~ 10 iterations, and the totalCPU time was modest by our standards (about 2 yr on anIBM RS-6000/370). We took g;„= 10 and L;„=128and used a quintic fit in (3); the result for Ft is shown inFig. 1 [y = 75.41, 66 DF (No. of degrees of freedom),level = 20%]. The extrapolated values from different lat-tice sizes at the same p agree within the estimated statisti-cal errors (g2 = 43.03, 75 DF, level ) 99%); see Table Ifor an example. The result for G~ is shown in Fig. 2; theerrors are roughly constant for $o/L ) 0.4 but rise sharplyfor smaller g /L. The theoretical exponent p computed[16] from the fitted Ft is equal to p = —0.10 ~ 0.06; thecurve suggests that p = 0. In practice we were able toobtain $ to an accuracy of about 1% (2%,3%,5%) at= 1000 (2000,5000, 10000).FIG. 1. $(p, 2L)/$(p, L) vs $(p, L)/L for the two-dimensional three-state Potts antiferromagnet. Symbolsindicate L = 32 (+), 64 (X), 128 ( ), 256 (C), and 512 (Q).Error bars are 1 standard deviation. Curve is a quintic fit in(3), with $;„=10 and L;„=128.TABLE I. Raw and extrapolated correlation lengths for thetwo-dimensional three-state Potts antiferromagnet at p = 3.5.Extrapolation based on g;„= 10 and L;„=128 and a quinticfit. For each extrapolated value we reported the standarddeviation of the estimate, including errors of all three types.The mean value and the g2 have been computed, taking intoaccount the full covariance matrix [14].32641282565121024mean(L. ~ 128)Iterations2.2 X 10'1061061065 x 10'2 X 10'ExtrapolatedF(p, )90.68 (1.03)92.47 (0.81)93.01 (0.44)93.01 (0.31)92.85 (0.39)93.78 (1.17)93.01 (0.29)~' = 0.68 (3 DF, level = 88%)Raw Dataa(p, L)19.02 (0.01)35.52 (0.04)60.61 (0.09)84.69 (0.15)92.48 (0.33)93.78 (1.17)Next let us consider [3,4] the two-dimensional O(3) o.model. We used the Wolff embedding algorithm withstandard Swendsen-Wang updates [6,20,21]; again criticalslowing down appears to be completely eliminated. Weran on lattices L = 32, 48, 64, 96, 128, 192, 256, 384,512 at180 different pairs (p, L) in the range 20 ( $ ( 10~.Each run was between 10 and 5 && 10 iterations, and thetotal CPU time was 7 yr on an IBM RS-6000/370. Wetook $;„=20 and used a tenth-order fit. There appearto be weak corrections to scaling (of order (1.5%) inthe region 0.3 ~ (I /L ( 0.7 for lattices with L ( 64—96. We therefore chose L;„=128 for $L/L ( 0.7, andL;„=64 for $1 /L ) 0.7. The result for Ft is shown inFig. 3 (~ = 72.91, 73 DF, level = 48%). The result for2970VOLUME 74, NUMBER 15 PH YS ICAL REVIEW LETTERS 10 ApRIL 1995I I I I I I I I I I I I I I I I010 I I I I I III I I I I I III0 20.04—100.02—100.000. 1I I I I I I I I101010 100( /L1000I I I I I III I I I I I III I I I I I IIIFIG. 2. Relative variance-time product [for errors (i) +(ii) only] divided by g„(P)~, plotted vs g„(P)/L, for thetwo-dimensional three-state Potts antiferromagnet. Symbolsindicate L = 128 ( ), 256 (0), 512 (Q), 1024 (+), and1536 (+).FIG. 4. Relative variance-time product [for errors (i) + (ii)only] divided by $„(P)2, plotted vs g (P)/L, for the two-dimensional O(3) o. model. Symbols indicate L = 64 (X), 96(+), 128 ( ), 192 (&), 256 (0), 384 (+), and 512 (Q). Forcomparison, the line shows the theoretical limiting slope —2.G~ is shown in Fig. 4; at large g /L it decreases sharply,with a power p = —2, in agreement with theory [16]. Inpractice we obtained g to an accuracy of about 0.2%(0.7%,1.1%,1.6%) at g = 10 (10,10,10').We also carried out a "simulated Monte Carlo" experi-ment for the O(N) rr model at N = ~, by generating dataI I I I I I I I I I I I I I I I I I2.01.00.0I I I I I I I I I I I I I I I I I I I0.2 0,4 0.6 0.8(-(L)/L1.0FIG. 3. $(P, 2L)/s(P, I ) vs g(P, L)/L for the two-dimensional 0 (3) o. model. Symbols indicate L = 32(+), 48 (+), 64 (x), 96 (X), 128 ( ), 192 (&), and 256 (0).Error bars are 1 standard deviation. Curve is a tenth-order fitin (3), with g;„= 20 and L;„=128 (64) for g(L)/L ~ 0.7(~ 0.7).from the exact finite-volume solution plus random noiseof 0.1% for I. = 64,96,128, 0.2% for I = 192,256, and0.5% for L = 384,512 [which is the order of magnitudewe attain in practice for O(3)]. We considered 35 valuesof P in the range 20 ~ s ~ 10 . We used s,„=20 andL;„=64 (in fact much smaller values could have beenused, as corrections to scaling are very small here) and atenth-order fit; for two different data sets we obtain g~ =159.18 and 179.39 with 165 DF. In practice we obtainwith an accuracy of 0.6% (1.2%,2%,3%) at g = 10'(10,10,10 ). Here we can also compare the extrapolatedvalues $'""(P) with the exact values g'""'(P). Defin-ing R = gp[g"'"(P) —g'""'(P)]2/o. 2(P), we find forthe two data sets R = 28.05 (32.53) with 35 DF. Only13 (10) points differ from the exact value by more than 1standard deviation, and only 1 (1) by more than 2. Detailson all of these models will be reported separately [2,4].The method is easily generalized to a model controlledby a RG fixed point having k relevant operators. Itsuffices to choose k —1 dimensionless ratios of long-distance observables, call them R = (RI, . . . , Rk I), thenthe function F@ will depend parametrically on R(P, L). Inpractice one can divide R space into "slices" within whichF@ is empirically constant within error bars, and performthe fit (3) within each slice. We used this approach tostudy the mixed isovector-isotensor o- model, taking R tobe the ratio of isovector to isotensor correlation length[3,4].The method can also be applied to extrapolate theexponential correlation length (inverse mass gap). For2971VOLUME 74, NUMBER 15 PHYSICAL REVIEW LETTERS 10 ApRrr 1995this purpose one must work in a system of size I" ' X Twith T )) $,„~(P,L) (compare [9]).We wish to thank Martin Hasenbusch and especiallyJae-Kwon Kim for sharing their data with us, andfor challenging us to push ever larger values of g/L.The authors' research was supported by CNR, INFN,CNPq, FAPEMIG, and DOE Contracts No. DE-FG05-85ER250000 and No. DE-FG05-92ER40742, NSF GrantNo. DMS-9200719, and NATO CRG 910251.*Electronic address: CARACCIO LE.INFN. ITElectronic address: EDWARDS SCRI.FSU.EDU~Electronic address: SABINO@FISICA.UFMG. BR~Permanent address: Dipartirnento di Fisica and INFN-Sezione di Pisa, Universita degli Studi di Pisa,I-56100 Pisa, Italy. Electronic address:PELISSET@SUNTHPI1.DIFI.UNIPI. ITSOKAL@ACF4.NYU. EDU[I] Finite Size Scaling and Numerical Simulation of Statistical Systems edited by V. Privman (World Scientific,Singapore, 1990).[2] S.J. Ferreira and A. D. Sokal, Phys. Rev. B (to bepublished); (to be published).[3] S. Caracciolo, R. G. Edwards, A. Pelissetto, and A. D.Sokal, Phys. Rev. Lett. 71, 3906 (1993).[4] S. Caracciolo, R. G. Edwards, A. Pelissetto, and A. D.Sokal, Reports No. hep-lat/9411009 and No. 9411064; (tobe published).[5] It is necessary to choose a definition of g that makessense in a fully finite system. We use the second-momentcorrelation length defined by Eqs. (4.11)—(4.13) ofRef. [6).[6] S. Caracciolo, R. G. Edwards, A. Pelissetto, and A. D.Sokal, Nucl. Phys. 8403, 475 (1993).[7] This form of finite-size scaling assumes hyperscaling, andthus is expected to hold only below the upper criticaldimension of the model. See, e.g. , Ref. [1], Chap. I,Sec. 2.7.[8] Our method has many features in common with thoseof Liischer, Weisz, and Wolff [9] and Kim [10]. Inparticular, all these methods share the property of workingonly with observable quantities (g, 6, and L) and not withbare quantities (P). Therefore, they rely on "scaling" andnot on "asymptotic scaling;" they differ from other FSS-based methods such as phenomenological renormalization.[9)[10)[11][12][13][14][15)[16][17][18][19][201[21]M. Luscher, P. Weisz, and U. Wolff, Nucl. Phys. B359,221 (1991).J.-K. Kim, Phys. Rev. Lett. 70, 1735 (1993); Nucl. Phys.B (Proc. Suppl. ) 34, 702 (1994); Phys. Rev. D 50, 4663(1994); Europhys. Lett. 28, 211 (1994).H. Neuberger, Phys. Lett. B 233, 183 (1989);S. Caracciolo and A. Pelissetto (to be published).In principle, g and 6 should be generated from a jointGaussian with the correct covariance. We ignored thissubtlety and simply generated independent fluctuations on$ and 6.Errors of types (i) and (ii) can be exactly computedin terms of ft and fo. For $ we have b, g„/gKt(g /L)Ag&/g&, where 5$L (Ag„) is the standard de-viation on the raw (extrapolated) value, and K&(z) =ft(z)/[ft(z) + zft(z)]. For a deviation, see Ref. [4].The extrapolated estimates from different (P, L) arecorrelated. In setting the error bars we have kept accountof the full covariance matrix between extrapolations at thesame P but different L, but we have ignored correlationsbetween extrapolations at different values of P.This variance-time product tends to a constant as theCPU time tends to infinity. However, if the CPU timeused is too small, then the variance-time product canbe significantly larger than its asymptotic value, due tononlinear cross terms between error sources (i) and (ii).The exponent p, which does not depend on s, can be com-puted in terms of F~. Define an exponent h by assumingthat b gt. /$L = $„'H($ /L)/N;„, —g„'+"/N;„„where Apt.is the error on g& for a run of length N;„,. If jt/Lx. ~ ~, when g„/L ~ ~ (this is the case for models withcritical exponent 7/ ) 0) and R —= 1 + x„F&(x„)/s ~ 1,then p = —d + 21ogR/logs + 2h. If R = 1, logarithmicterms appear and G&(z) —z "+2"(logz)~, with q ) 0. Ifinstead g&/L ~, when s /L ~ (this is the case forasymptotically free theories), then Ge(z) —z "+ "(Iogz) .In our three examples, the total error (i) + (ii) + (iii) isnever more than twice the error (i) + (ii), and is usuallymuch less.R. J. Baxter, Proc. R. Soc. London A 383, 43 (1982).J.-S. Wang, R. H. Swendsen, and R. Kotecky, Phys. Rev.Lett. 63, 109 (1989); Phys. Rev. B 42, 2465 (1990);M. Lubin and A. D. Sokal, Phys. Rev. Lett. 71, 1778(1993).U. Wolff, Phys. Rev. Lett. 62, 361 (1989).R. H. Swendsen and J.-S. Wang, Phys. Rev. Lett. 58, 86(1987).2972

Extrapolating Monte Carlo Simulations to Infinite Volume: Finite-Size Scaling at ξ/L ≫1

http://discovery.ucl.ac.uk/1323422/1/1323422.pdf

Extrapolating Monte Carlo Simulations to Infinite Volume: Finite-Size Scaling at ξ/L ≫1

Abstract

Similar works

Full text

Available Versions

UCL Discovery