θ -dependence in the small- N limit of 2d C P N−1 models
aa r X i v : . [ h e p - l a t ] S e p θ -dependence in the small- N limit of d CP N − models Mario Berni, ∗ Claudio Bonanno, † and Massimo D’Elia ‡ Universit`a di Pisa and INFN Sezione di Pisa,Largo Pontecorvo 3, I-56127 Pisa, Italy (Dated: September 30, 2020)We present a systematic numerical study of θ -dependence in the small- N limit of 2 d CP N − models, aimed at clarifying the possible presence of a divergent topological susceptibility in thecontinuum limit. We follow a twofold strategy, based on one side on direct simulations for N = 2and N = 3 on lattices with correlation lengths up to O (10 ), and on the other side on the small- N extrapolation of results obtained for N up to 9. Based on that, we provide conclusive evidence for afinite topological susceptibility at N = 3, with a continuum estimate ξ χ = 0 . N = 2 are still inconclusive: they are consistent with a logarithmicallydivergent continuum extrapolation, but do not yet exclude a finite continuum value, ξ χ ∼ .
4, withthe divergence taking place for N slightly below 2 in this case. Finally, results obtained for thenon-quadratic part of θ -dependence, in particular for the so-called b coefficient, are consistent witha θ -dependence matching that of the Dilute Instanton Gas Approximation at the point where ξ χ diverges. PACS numbers: 12.38.Aw, 11.15.Ha,12.38.Gc,12.38.Mh
1. INTRODUCTION
The CP N − models in two space-time dimensions havebeen extensively studied in the literature because theyrepresent an interesting theoretical laboratory for thestudy of gauge theories [1–3]. As a matter of fact, manyintriguing non-perturbative properties, such as confine-ment, the existence of field configurations with non-trivial topology and the related θ -dependence are featuresthat these models share with 4 d Yang-Mills theories.The Euclidean action of these models, includingthe topological term, can be written through a non-propagating Abelian field A µ as S ( θ ) = Z d x (cid:20) Ng ¯ D µ ¯ z ( x ) D µ z ( x ) − iθq ( x ) (cid:21) , (1)where N is the number of components of the complexscalar field z , which is subject to the constraint ¯ zz = 1, D µ = ∂ µ + iA µ and Q = Z d x q ( x ) = 12 π ǫ µν Z d x ∂ µ A ν ( x ) (2)is the topological charge. The θ -dependent vacuum en-ergy (density) is defined through the path integral as E ( θ ) = − log Z ( θ ) V = − V log Z [ d ¯ z ][ dz ][ dA ] e − S ( θ ) (3)and it can be parametrized in terms of the cumulants h Q n i c of the topological charge distribution at θ = 0 ∗ Electronic address: [email protected] † Electronic address: [email protected] ‡ Electronic address: [email protected] P ( Q ) as: f ( θ ) ≡ E ( θ ) − E (0) = 12 χθ ∞ X n =1 b n θ n ! , (4)where the topological susceptibility χ ≡ h Q i c V = h Q i V (5)parametrizes the leading θ term while the coefficients b n ≡ ( − n n + 2)! h Q n +2 i c h Q i (6)parametrize the non-quadratic part.One of the most interesting features of the CP N − models is the possibility of performing a systematic ex-pansion of any observable, including those related to θ -dependence, in the inverse of the number of field compo-nents 1 /N when N → ∞ while g is kept fixed [1]; thatclosely resembles the 1 /N expansion of 4 d SU ( N ) gaugetheories for large number of colors. The similarities arenon-trivial, since the relevant scaling quantity turns outto be θ/N in both cases [4–7], leading to the predictionthat the quantities ¯ b n ≡ N n b n are finite in the large- N limit: this scaling behavior has been verified by numeri-cal lattice simulations both for the CP N − [7, 8] and for SU ( N ) gauge theories [7, 9] up to the quartic coefficient b . From a quantitative point of view, CP N − modelsare more predictive, since the leading term in the 1 /N expansion is known [1, 6, 7, 9] for all coefficients in the θ -expansion of the free energy in Eq. (4), and also the firstsubleading term in the case of the topological susceptibil-ity [10], while in the case of 4 d SU ( N ) gauge theories onehas just phenomenological predictions, based on the spec-trum of pseudoscalar mesons, for the leading term in the1 /N expansion of the topological susceptibility [11, 12].Numerical simulations underlined also some differencesbetween the two theories. Indeed, while in the Yang-Mills case the large- N expected scaling practically holdsalready for N ≥ CP N − models. Indeed, the large- N limit of b showssignificant deviations from large- N predictions even for N ∼
50, and an agreement with lattice data can be re-covered only by including large higher-order correctionsin 1 /N ; a similar behavior is observed for the topologicalsusceptibility [8, 13].Another important difference emerges when looking atthe small- N limit. Indeed, while this is predicted andobserved to be regular in the case of 4 d SU ( N ) gaugetheories, a singular behavior is expected for 2 d CP N − approaching N = 2. For instance, one expects a diver-gence of the topological susceptibility of the CP theory,which can be justified in perturbation theory on the basisof the ultraviolet (UV) divergence of the instanton sizedistribution [14–18] P N ( ρ ) ∝ ρ N − (7)occurring for N = 2. This result has been tested in manylattice studies and there seems to be a general consensusabout the singular behavior of χ for N = 2 and aboutits origin due to the presence of small instantons, see,e.g., Refs. [19–24]. However, there are still many aspectsdeserving a more careful investigation.First of all, as we discuss later on, the actual verifica-tion of the divergent behavior occurring for N = 2 re-quires to perform a continuum limit extrapolation whichis highly non-trivial, since the divergence is expected toappear as a logarithm of the UV cutoff, i.e., of the latticespacing a , which could be difficult to disentangle from aregular power-law behavior in a and requires numericalresults spanning over orders of magnitude in terms of thedimensionless correlation length of the system.The second issue is that a divergent behavior for thetopological susceptibility has been claimed to be observedalso for N = 3 in Ref. [22], where the authors employeda standard lattice action along with an overlap definitionof the topological charge, however this is in contradictionwith previous results obtained using a different actiondiscretization and a geometric lattice definition of Q [25]:the origin of this discrepancy may be due to the fact that,according to Eq. (7), also for N = 3 small instantonsare expected to dominate, so that also in this case thecontinuum limit has to be handled with care.Finally, one would like to understand whether thesmall- N divergent behavior regards just the topologicalsusceptibility or also other coefficients in the Taylor ex-pansion of the free energy in Eq. (4). In asymptoticallyfree theories, small size instantons are expected to beweakly interacting, so that a possible conjecture [13] isthat θ -dependence in the N = 2 limit be described bythe Dilute Instanton Gas Approximation (DIGA): f DIGA ( θ ) = χ (1 − cos θ ) (8) with the divergence appearing just in χ , while the b n coefficients are finite and have alternate sign: b DIGA n = ( − n n + 2)! . (9)A DIGA-like small- N behavior of f ( θ ) could explainwhy, unlike the case of SU ( N ) gauge theories, largecorrections are observed when studying the large- N limit of CP N − models, being the 1 /N series at large N not able the capture the small- N behavior of thetheory. A first evidence of near-DIGA behavior for thequartic coefficient b for N = 2 has been reported inRef. [26], where however no continuum extrapolation forthis quantity is reported. As for higher values of N , noresult is known for N <
N θ -dependence of CP N − models bylattice simulations, in order to go beyond the presentstate of the art. To do so, we have attacked the prob-lem from two different sides. On one hand, we haveperformed extensive numerical simulations for N = 2and N = 3 considering dimensionless correlation lengthsspanning over two orders of magnitude, namely reach-ing values of ξ going above 10 , and considering variousdifferent ansatzes for the continuum extrapolation in or-der to fairly assess our systematic uncertainties on thefinal results. On the other hand, we have also consid-ered numerical simulations for larger values of N , namely N ∈ [4 , N = 2 , N extrapolation of these results.We have applied this double-front strategy to the de-termination of the topological susceptibility and of thefirst coefficients of the θ expansion of f ( θ ), up to O ( θ ):consistency of results obtained in the two different waysprovides a solid way to assess the reliability of our finalstatements, among which, for instance, the fact that the θ -dependence is finite for N = 3.The paper is organized as follows. In Sec. 2 we describethe lattice setup adopted for discretizing the theory andfor determining the cumulants of the topological chargedistribution, as well as our strategy for taking the contin-uum limit of our results. In Sec. 3 we present and discussour numerical results and finally, in Sec. 4, we draw ourconclusions.
2. NUMERICAL SETUP
In this section we briefly discuss various issues relatedto the discretization of the model and of its observables,in particular those related to topology, and to the con-tinuum extrapolation of the numerical results.
A. Lattice discretization
We discretized the action in Eq. (1) at θ = 0 on a peri-odic square lattice of size L using the tree-level Symanzik-improved lattice discretization [27]: S L = − N β L X x,µ (cid:8) c ℜ (cid:2) ¯ U µ ( x )¯ z ( x + ˆ µ ) z ( x ) (cid:3) + c ℜ (cid:2) ¯ U µ ( x + ˆ µ ) ¯ U µ ( x )¯ z ( x + 2ˆ µ ) z ( x ) (cid:3)(cid:9) , (10)where c = 4 / c = − /
12 are improvement coeffi-cients, β L ≡ /g L is the inverse bare coupling and U µ ( x )are the U (1) gauge link variables. The adoption of theimproved action cancels out logarithmic corrections tothe leading O ( a ) behavior of the discretization errors,where a is the lattice spacing.Being CP N − models asymptotically free for all valuesof N , the continuum limit is approached as β L → ∞ .The a → ξ L . Our choice is for the secondmoment correlation length ξ , defined in the continuumtheory in terms of the two-point correlation function of P ij ( x ) ≡ z i ( x )¯ z j ( x ) G ( x ) ≡ h P ij ( x ) P ij (0) i − N (11)as ξ ≡ R G ( x ) d x Z G ( x ) | x | d x. (12)To define the lattice length ξ L , we adopted the followingdefinition [28], expressed through the Fourier transform˜ G L ( p ) of G L ( x ) (i.e., the discretization of Eq. (11)): ξ L = 14 sin ( π/L ) " ˜ G L (0)˜ G L (2 π/L ) − . (13) B. Discretization of the topological charge
Regarding the topological charge Q , several equivalentlattice discretizations Q L can be adopted, all having thesame continuum limit. However, at finite lattice spac-ing, these definitions and their correlations are related tothe continuum by a finite multiplicative renormalization Z [29, 30]: q L ∼ Za q + O ( a ) (14)where q is the topological charge density. We adopted a geometric definition of the lattice charge [27, 31], i.e., adefinition with Z = 1, meaning that it yields always aninteger number for every lattice configuration. Amongthe several geometric definitions, we chose one which in-volves only the link variables [27]: Q U = X x q U ( x ) = 12 π X x ℑ { log [Π ( x )] } , (15) where Π µν ( x ) ≡ U µ ( x ) U ν ( x + ˆ µ ) ¯ U µ ( x + ˆ ν ) ¯ U ν ( x ) is theplaquette operator. Despite the fact that Z = 1, renor-malization effects are still present, since in general Q L isrelated to the physical charge Q , configuration by config-uration, by a relation like Q = Z ( β L ) Q L + η (16)where η is a noise with zero average stemming fromfluctuations at the UV scale. For a geometric charge,such noise appears in the form of the so-called disloca-tions [32, 33], i.e., integer valued fluctuations at the scaleof the UV cutoff and proliferating over physical contri-butions as the continuum limit is approached. The pres-ence of a non-zero η gives rise to further additive renor-malizations as one considers cumulants of the topolog-ical charge, including the topological susceptibility (seeRef. [34] for a review).To suppress such noise, a smoothing method, such ascooling [32, 35–40] or the gradient flow [41, 42], can beadopted. Since it has been shown that they are all nu-merically equivalent [43, 44], we adopted cooling for itssimplicity and numerical cheapness. Indeed, this methodconsists in a sequence of local steps of action minimiza-tion that damp UV fluctuations: for simplicity, the actionthat is minimized is the non-improved one, i.e., that inEq. (10) with c = 1 and c = 0.Contrary to Refs. [8, 13], in this study we did not adoptneither simulations at imaginary values of θ in order toimprove the signal-to-noise ratio of cumulants [7, 45, 46],nor improved algorithms in order to defeat the criticalslowing down of topological modes [47–50]. This is dueto the relative ease in obtaining precise determinations ofthe cumulants of the topological charge for small values of N , even using standard algorithms such as heat-bath orover-relaxation. For these reasons, the coefficients b n aredetermined in this study by simply using the definitiongiven in Eq. (6), with the average taken over the pathintegral distribution at θ = 0, where the choice for Q isthe geometrical topological charge in Eq. (15) measuredafter a certain number of cooling steps, as discussed lateron. C. Continuum limit at small N Since ξ L = ξ/a diverges as 1 /a in the continuumlimit, finite lattice spacing corrections can be expressedas a function of 1 /ξ L . Thus, with the adoption of theSymanzik-improved action in Eq. (10), one expects ul-traviolet corrections to the lattice expectation value of ageneric observable O to have the form hOi latt ( ξ L ) = hOi cont + c ξ − L + O (cid:0) ξ − L (cid:1) . (17)However, when approaching N → a , henceon ξ L .An a priori estimate of such effects can be done onlywith some assumptions, nevertheless it can be a usefulguide. For instance, taking the perturbative estimatefor the instanton size distribution reported in Eq. (7)and assuming that topological fluctuations are domi-nated by a non-interacting gas of small instantons andanti-instantons, one has that the number of instantons n I and anti-instantons n A are distributed as two inde-pendent Poissonians with equal mean: h n I i = h n A i ∝ Z ρ a P N ( ρ ) dρ = Z ρ a ρ N − dρ, (18)where the integral is carried over sizes ranging from theUV scale, set by the lattice spacing a , up to a certain IRlength scale ρ , which is proportional to La . Since, withthese hypotheses, χ ∝ h ( n I − n A ) i = 2( h n I i − h n I i ) = 2 h n I i (19)we have the following predictions: χ ∝ ρ N − − a N − N − , if N > (cid:16) ρ a (cid:17) , if N = 2 . (20)Ultraviolet corrections predicted by Eqs. (20) becomenegligible as N grows, in particular one expects them todisappear for N ≥
4, where the contribution from smallinstantons becomes negligible. On the other hand, thecontribution of small instantons becomes dominant for N = 2 and 3, where it leads either to a logarithmicallydivergent continuum limit for N = 2, or to linear, insteadof quadratic, corrections in the lattice spacing for N = 3.Notice that, under the assumptions of independent Pois-son distributions for n I and n A , the b n coefficients arefinite with values as predicted by DIGA in Eq. (9), sothat no further corrections are expected, with respect toEq. (17), in their approach to the continuum limit.Such considerations will be used with caution in thefollowing. In particular, in order to correctly assess oursystematics on the continuum limit, we will consider thepossible presence of generic power law corrections in thelattice spacing for N ≤
4, both for χ and for the otherterms in the θ -expansion. D. Continuum limit and smoothing
Since we adopt a smoothing method in order to re-move field fluctuations at the UV cutoff scale, which areresponsible for unphysical contributions to lattice topo-logical observables, we need to fix how the amount ofsmoothing is changed as one approaches the continuumlimit. Smoothing algorithms work in general as diffusive pro-cesses, affecting field correlations up to a given distance,i.e., up to a given smoothing radius , which is propor-tional, in dimensionless lattice units, to the square rootof the amount of smoothing, i.e., to √ n cool for cooling,where n cool is number of cooling steps, or to √ t for thegradient flow, where t is the flow time. It is a standardprocedure to change the amount of smoothing so thatthe smoothing radius is kept fixed in physics units: thatwould correspond to change n cool proportionally to ξ L .However, for a model like the one we are exploring, wheresmall distance physical contributions are expected to bequite significant, one should be careful and take care, ad-ditionally, of the dependence of continuum results on thephysical smoothing radius, eventually sending it to zero.An alternative to this difficult double limit procedureis to take the continuum limit at fixed number of coolingsteps n cool , so that the smoothing radius goes to zero pro-portionally to a and there is no possibility that physicalcontributions at small scales are smoothed away. Thereare good reasons to believe that such procedure workscorrectly in the present case.As we have discussed above, for a geometric chargelike the one used in this study renormalization effects areessentially due to dislocations leading to a wrong and/orambiguous counting of the topological winding number.Such dislocations consist of exceptional field fluctuationsliving at the lattice spacing scale, hence it is reasonableto expect that they will be suppressed by a given andfixed number of cooling steps n cool , independently of thevalue of the lattice spacing a .Based on such considerations, in the following we willconsider results obtained by performing the continuumlimit at fixed n cool , then carefully checking the possi-ble systematics related to this procedure. In particular,we will show that while results obtained at finite latticespacing, but at different number of cooling steps, usuallydiffer from each other well beyond their statistical errors,the so obtained continuum-extrapolated results are notsignificantly dependent on n cool , and usually well withinstatistical errors.
3. NUMERICAL RESULTS
In Tables I-VII we summarize the parameters of theperformed simulations, along with the total accumulatedstatistics. Configurations were generated using standardlocal algorithms, in particular our elementary MonteCarlo step consisted in 4 lattice sweeps of over-relaxationfollowed by a sweep of over-heath-bath; measures weretaken every 10 Monte Carlo steps.We simulated CP N − models with N ranging from 2 to8. For each value of N we simulated several runs at dif-ferent values of the correlation length (i.e., at differentvalues of β L ); for each ξ L we measured ξ χ , b and b inorder to be able to extrapolate these quantities towardsthe continuum limit.Concerning the choice of the lattice size, we performedsimulation at fixed L/ξ L , ensuring that L/ξ L ∼ ξ L runs for N = 2), several latticeswith different sizes were simulated and then the infinitevolume limit was performed by fitting the L -dependenceof every observable O with the law O ( L ) = O ∞ (cid:16) − a e − bL/ξ L (cid:17) , (21)where O ∞ is the desired quantity and a and b are addi-tional fit parameters. An example of extrapolation to-wards the thermodynamic limit is shown in Fig. 2. . . . × − ξ χ L/ξ L − − × − b FIG. 1: Example of finite size scaling of ξ χ and b as afunction of L/ξ L for N = 6 and β L = 0 .
80 with measurestaken after n cool = 50 cooling steps. A. Results for ξ χ , N > First, we consider the case
N >
3, for which a finitecontinuum limit is surely expected for the topologicalsusceptibility, with no qualitative differences in the con-tinuum scaling compared to the large- N case. For thisreason, in order to extrapolate the quantity ξ χ towardsthe continuum limit we have fitted its dependence on ξ L according to the ansatz f ( x ) = a + a x + a x , x = 1 /ξ L . (22)Only for N = 4, due to its proximity to N = 2 and 3 (seelater discussion), we have also considered the possiblepresence of further power law corrections. An example . . ξ χ − . − . b L/ξ L ξ L FIG. 2: Examples of extrapolation towards the infinite vol-ume limit of ξ χ , b and ξ L for N = 2 and β L = 1 .
65 withmeasures taken after n cool = 50 cooling steps. Best fits ofEq. (21) performed with 4 degrees of freedom (dof) yield, re-spectively, χ / dof = 1 . , . , . of continuum extrapolation is depicted, for N = 4, inFig. 3.Several sources of systematic errors have been checked:first, the extrapolation was performed fitting data in sev-eral ranges of ξ L to check that the obtained extrapola-tions were all consistent with each other and that, asthe fit range is restricted, the O ( x ) corrections becamenegligible. Second, as anticipated in Sec. 2 D, we extrap-olated the continuum limit at fixed number of coolingstep n cool for several values of n cool , checking that thisprocedure does not introduce significant systematics incontinuum-extrapolated values.When changing the number of cooling steps, we ob-serve that, while measures at coarse lattice spacing dif-fer, the dependence on n cool is less and less visible asthe continuum limit is approached, making the contin-uum extrapolation stable as n cool is varied. In Fig. 4 weshow an example of the continuum extrapolation at twodifferent values of n cool for N = 5 while, in Fig. 5, weshow, for the same case, that any variation in the contin-uum extrapolation observed when changing n cool is wellcontained inside our statistical errors. The final, con-tinuum determinations obtained for ξ χ are reported inTab. VIII. β L L ξ L L/ξ L Stat. (M)0.80 22 1.704(13) 12.9 3.20.85 26 2.037(16) 12.8 3.20.90 30 2.4764(72) 12.1 20.80.95 36 3.0309(28) 12.0 2011.00 46 3.7726(48) 12.3 1301.05 58 4.7345(59) 12.3 1501.10 74 6.0142(91) 12.2 1141.15 94 7.735(16) 12.3 66.71.20 122 9.958(51) 12.4 12.21.25 160 13.10(14) 12.2 3.41.30 210 17.16(19) 12.2 3.91.35 280 22.77(34) 12.3 3.61.40 98 29.031(40) 3.4 3.5198 31.115(55) 6.4 13.5298 31.32(19) 9.5 6.0396 31.41(50) 12.6 3.41.45 122 38.167(55) 3.2 3.5240 41.352(89) 5.8 8.7360 41.48(28) 8.7 4.0480 42.30(76) 11.3 2.21.50 168 51.43(11) 3.3 1.5336 55.75(21) 6.0 4.5504 57.26(78) 8.8 1.9672 56.4(1.5) 11.9 2.01.55 214 67.95(21) 3.1 1.5400 74.97(31) 5.3 3.8500 74.60(57) 6.7 2.4600 78.0(1.0) 7.7 1.6700 76.9(1.7) 9.1 1.21.60 260 88.13(34) 3.0 1.1400 97.92(33) 4.1 3.8500 98.95(59) 5.1 2.4600 102.23(94) 5.9 1.6700 98.6(1.5) 7.1 1.21.65 200 89.44(16) 2.2 2.4306 112.57(63) 2.7 0.77400 123.54(60) 3.2 1.5500 129.82(52) 3.9 3.7600 131.98(80) 4.5 2.5700 132.1(1.2) 5.3 1.81024 131.4(1.5) 7.8 6.3TABLE I: Simulations summary for N = 2. Statistics is ex-pressed in millions (M) and every measure is taken after 10elementary Monte Carlo steps (see the text for further de-tails). B. Results for ξ χ , N = 3 In the N = 3 case we fit the dependence of ξ χ on ξ L according to the following function: g ( x ) = a + a x + a x + a x c , x = 1 /ξ L , (23)where we consider both the extrapolations obtained withfixed c = 1, as suggested by the ansatz in Eq. (20), andwith c left as a free parameter. In both cases, the ansatzin Eq. (23) provides a very good description of our nu-merical results, and the x corrections turn out to be β L L ξ L L/ξ L Stat. (M)0.70 18 1.4797(59) 12.2 3.50.75 22 1.8227(76) 12.1 3.50.80 28 2.298(10) 12.2 3.50.85 36 2.906(14) 12.4 3.50.90 46 3.750(31) 12.3 1.50.95 58 4.930(23) 11.8 3.50.975 68 5.710(29) 11.9 3.21.00 80 6.588(25) 12.1 7.81.025 92 7.589(28) 12.1 8.31.05 106 8.828(32) 12.1 8.31.075 122 10.231(37) 11.9 8.81.10 146 11.834(39) 12.3 13.51.15 200 16.019(55) 12.5 15.51.20 264 21.717(71) 12.2 19.41.25 374 29.342(82) 12.7 37.81.366 720 58.69(19) 12.3 44.9TABLE II: Simulations summary for N = 3. β L L ξ L L/ξ L Stat. (M)0.70 22 1.7842(59) 12.3 2.80.75 30 2.3137(97) 13.0 2.80.80 38 3.039(12) 12.5 3.10.85 50 4.009(17) 12.5 3.10.90 66 5.398(22) 12.2 3.10.95 90 7.314(33) 12.3 3.11.00 120 9.997(45) 12.0 3.11.05 170 13.700(71) 12.4 3.11.10 226 18.55(10) 12.2 3.51.15 320 25.00(20) 12.8 3.5TABLE III: Simulations summary for N = 4. β L L ξ L L/ξ L Stat. (M)0.70 32 2.1263(98) 15.0 3.50.75 44 2.830(13) 15.5 4.60.80 58 3.882(20) 14.9 3.50.85 84 5.270(35) 15.9 3.50.90 106 7.149(41) 14.8 3.50.95 146 9.849(59) 14.8 3.51.00 198 13.40(11) 14.8 2.21.05 278 18.117(93) 15.3 7.41.10 368 24.77(19) 14.9 3.6TABLE IV: Simulations summary for N = 5. necessary only when considering in the fit range correla-tion lengths as small as ξ L .
5. Moreover, even when c is treated as a free parameter, its values turn out tobe compatible with 1 within errors, thus giving furthernumerical support to the ansatz in Eq. (20).Examples of continuum extrapolations are shown inFig. 6, where we also report our final continuum estimatefor N = 3, ξ χ = 0 . a when changing the fitting ansatz, the fitting β L L ξ L L/ξ L Stat. (M)0.60 18 1.4035(36) 12.8 2.60.65 24 1.8608(53) 12.9 2.60.70 32 2.5032(76) 12.8 2.50.75 42 3.394(11) 12.4 2.20.80 58 4.671(15) 12.4 2.40.85 78 6.389(21) 12.2 2.60.90 108 8.760(31) 12.3 2.60.95 150 11.955(48) 12.5 2.81.00 200 16.112(41) 12.4 6.5TABLE V: Simulations summary for N = 6. β L L ξ L L/ξ L Stat. (M)0.70 44 2.876(12) 15.3 2.20.75 58 3.940(19) 14.7 1.70.80 80 5.403(26) 14.8 2.40.85 112 7.436(35) 15.1 3.10.90 154 10.096(48) 15.3 3.40.95 210 13.717(68) 15.3 3.61.00 304 18.76(11) 16.2 3.9TABLE VI: Simulations summary for N = 7. β L L ξ L L/ξ L Stat. (M)0.65 36 2.3175(76) 15.5 3.00.70 48 3.232(12) 14.8 2.20.75 68 4.452(17) 15.3 3.00.80 92 6.056(27) 15.2 2.40.85 124 8.242(42) 15.0 2.10.90 170 11.306(52) 15.0 2.70.95 240 15.164(73) 15.8 3.4TABLE VII: Simulations summary for N = 8. range (with ξ min varied between 2 and 7) and the numberof cooling steps (with n cool varied between 10 and 50).Let us discuss more in detail the systematics related tocooling. Fig. 6 shows that results obtained for n cool = 10and 50 are significantly different from each other, nev-ertheless, as also evident for one particular fit ansatz inFig. 6, their continuum limit shows very little variations,when compared to statistical errors on the fit parameters.There is a simple way to understand why results at dif-ferent n cool differ so much at finite ξ L , while coincidingin the ξ L → ∞ limit. As already discussed in Sec. 2 D,cooling acts as a diffusive process which smooths awayfield fluctuations (both physical and unphysical) belowan effective radius r = a ˆ r ( n cool ), where the radius inlattice spacing units, ˆ r , is a function of n cool only, i.e. in-dependent of the lattice spacing a . At fixed lattice spac-ing, different values of n cool lead to different values of r hence to different values of the topological susceptibil-ity, because a different amount of physical signal below r is removed: the effect can be particularly significant forsmall N , where topological fluctuations at small scales .
000 0 .
025 0 .
050 0 .
075 0 .
100 0 . /ξ L . . . . . . . ξ χ × − FIG. 3: Extrapolation towards the continuum limit of ξ χ for N = 4 using data in the range ξ L > n cool = 50cooling steps. The solid line represents the best fit obtainedusing fit function f ( x ) = a + a x (where x = 1 /ξ L ) while thedashed line represents the one obtained with f ( x ) = a + a x c .In the latter case the exponent is c = 1 . χ / dof = 2 . / χ / dof = 1 . / are more abundant. On the other hand, the effect mustfade away as a → a and different n cool , but such that r = a ˆ r ( n cool ) isthe same, should coincide: for instance, results shown inFig. 6 for n cool = 10 should go onto those at n cool = 50 ifthey are shifted along the horizontal axis by a constantfactor equal to ˆ r (10) / ˆ r (50). Such experiment has beenperformed for n cool = 10 , ,
50 in Fig. 7, showing thatit works perfectly: plotting all data as a function of aneffective variable x eff proportional to r/ξ = ˆ r ( n cool ) /ξ L ,where ˆ r ( n cool ) has been found empirically, data at dif-ferent values of n cool collapse perfectly onto each otherover the whole range of explored correlation lengths; asimilar collapse can be obtained also for other values of N . Therefore, the dependence of ξ χ on n cool observedat finite lattice spacing can indeed be simply interpretedin terms of a global effective rescaling of 1 /ξ L , so thatsuch dependence naturally fades away when ξ L → ∞ .To summarize, our results for N = 3 provide solid ev-idence that ξ χ is indeed finite in the continuum limitof this theory. On the other hand, that will be furtherchecked and supported by an independent small- N ex-trapolation based on N > .
000 0 .
005 0 .
010 0 .
015 0 . /ξ L . . . . . . . ξ χ × − n cool = 10 n cool = 40 FIG. 4: Extrapolation towards the continuum limit of ξ χ for N = 5 using data in the range ξ L >
7. The solid line repre-sents the best fit obtained using fit function f ( x ) = a + a x (where x = 1 /ξ L ) with measures taken after n cool = 40cooling steps while the dashed line represents the one ob-tained with the same fit function with n cool = 10 coolingsteps. The best fits yield, respectively, χ / dof = 1 . / χ / dof = 2 . /
3. Data obtained for different numbers of cool-ing steps have been plotted slightly shifted to improve read-ability. The diamond point represents our final estimation ofthe continuum limit. ξ min . . . . ξ χ × − n cool = 10 n cool = 20 n cool = 30 n cool = 40 FIG. 5: Example of study of systematic errors on the con-tinuum extrapolation of ξ χ for N = 5 and for 4 differentnumber of cooling steps n cool . Extrapolations obtained fordifferent values of n cool are plotted slightly shifted to improvereadability. This extrapolations are all obtained using the fitfunction f ( x ) = a + a x , where x = 1 /ξ L . The diamondpoint represents our final estimation of the continuum limit. .
00 0 .
05 0 .
10 0 . /ξ L . . . . . . . ξ χ × − n cool = 10 n cool = 50 FIG. 6: Extrapolation towards the continuum limit of ξ χ for N = 3 using data in the range ξ L >
6. The dashed and dottedlines represent, respectively, the best fits obtained for n cool =50 using fit function g ( x ) = a + a x + a x c (where x = 1 /ξ L )either leaving c as a free parameter ( χ / dof = 2 . / c = 1 χ / dof = 1 . / c = 0 . c as a free parameter for n cool = 10 cooling steps: in this case c = 0 . χ / dof =1 . / C. Results for ξ χ , N = 2 For the N = 2 case two possibilities are open. Theansatz based on Eq. (20) could be correct, leading to asusceptibility which diverges logaritmically in the contin-uum limit, i.e. has the following dependence on x = 1 /ξ L : ξ χ ( ξ L ) = a log( x/a ) + a x + O ( x ) . (24)On the other hand, corrections to such prediction leadingto a finite ξ χ cannot be excluded a priori, in this caseone should consider a dependence on ξ L like that used inEq. (23) for the N = 3 case, i.e.: ξ χ ( ξ L ) = a + a x c + a x + O ( x ) . (25)where c is a positive exponent. We have to say that, un-fortunately, despite the fact that our present results ex-tend over almost two orders of magnitudes in terms of thecorrelation length ξ L , we are still not able to clearly dis-tinguish between the two possibilities, in the sense thatboth ansatzes, Eqs. (24) and (25), return acceptable val-ues of the χ / dof test, even if marginally better for theconvergent ansatz.Let us consider, for instance, the data reported inFig. 8, where we consider only results obtained for ξ L > x eff ξ χ n cool = 50n cool = 25n cool = 10 FIG. 7: Results obtained for ξ χ in the N = 3 model and fordifferent number of cooling steps, plotted as a function of aneffective variable x eff proportional to ˆ r ( n cool ) /ξ L = r/ξ . Theperfect collapse of data proves that the dependence of ξ χ on n cool can be interpreted in terms of varying effective radius r = a ˆ r ( n cool ) below which the topological signal is lost, andas such disappears after continuum extrapolation a →
0. Theexact value of x eff has been fixed conventionally to x eff =(ˆ r ( n cool ) /ξ L ) / ˆ r ( n cool = 50) and we have found empiricallyˆ r (25) / ˆ r (50) ≃ . r (10) / ˆ r (50) ≃ .
10. In this range O ( x ) corrections turn out to be unnec-essary. Considering data obtained for n cool = 50 we ob-tain χ / dof = 7 . / χ / dof = 3 . / c = 0 . ξ χ = 0 . n cool = 20, the latter pre-diction changes to ξ χ = 0 . c = 0 . χ / dof = 3 . / c is at two standard deviationsfrom zero, which is the boundary where also this fit be-comes divergent: this is consistent with the fact that thelogarithmic fit is marginally acceptable.The situation does not change appreciably whenchanging the fit range, with both ansatzes remaining ac-ceptable even if with slightly lower values of the χ / doftest in the case of a finite continuum susceptibility: ifthe latter scenario could be assumed a priori we wouldconclude ξ χ ∼ . N = 2. However, the factthat no conclusive answer can be obtained based on ourpresent data for N = 2, is confirmed by looking at thebest fit curves reported in Fig. 8: the two curve profiles(convergent and divergent) are hardly distinguishable inthe explored range of ξ L and deviate from each otheronly for much larger values of ξ L . Another independentpossibility to discriminate between these two different be-haviors is to use data obtained for N > N = 2. This topic is covered in Sec. 3 D. .
01 0 .
02 0 .
03 0 .
04 0 .
05 0 .
06 0 .
07 0 . /ξ L ξ χ × − n cool = 20 n cool = 50 FIG. 8: Extrapolation towards the continuum limit of ξ χ for N = 2 using data in the range ξ L >
10. The solid and dashedlines represent, respectively, the best fits obtained using fitfunctions h ( x ) = a log( x/a ) + a x (where x = 1 /ξ L ) and g ( x ) = a + a x + a x + a x c with measure taken after n cool =50 cooling steps. The former best fit yields χ / dof = 7 . / a =0 . c = 0 . χ / dof = 3 . / g ( x ) = a + a x + a x + a x c but withmeasures taken after n cool = 20, which gives a = 0 . c = 0 . χ / dof = 3 . / N ξ χ · ξ χ as a function of N . The value marked with * is taken fromRef. [13]. D. N → limit of ξ χ In this section we aim at tackling the question aboutthe small N behavior of the topological susceptibilityfrom another independent front, by extrapolating con-tinuum results obtained for ξ χ at N > N = 2. A summary of such results, including all N < ξ χ when approaching the N → ξ χ ∼ N − , N > . (26)0Therefore, we fitted the N -dependence of ξ χ for N > F ( N ) = a ( N − N ∗ ) γ . (27) N . . . . . ξ χ fit with N ∗ = 2full fit FIG. 9: Best fit of the small- N behavior of ξ χ for N rangingfrom 3 to 9. The solid and dashed lines represent, respectively,the best fits obtained using fit function (27) with N ∗ = 2 orleft as a free parameter. The best fits yield, respectively, χ / dof = 5 . / χ / dof = 4 . / The best fit is quite good and is shown in Fig. 9, thecorresponding parameters are the following: a = 0 . ,N ∗ = 1 . ,γ = 0 . ,χ / dof = 4 . / . This best fit, yielding a central value for N ∗ <
2, tech-nically supports a finite extrapolation towards N = 2;however, N = 2 is well within one standard deviationfrom N ∗ , so that even this approach reveals to be incon-clusive, at least for this issue. As shown in Fig. 10, theerror on the best fit blows up as N → N ∗ very close to it, so that this extrapolation turnsout to be compatible with the hypothetical finite valuefound from the convergent continuum limit of Sec. 3 C, ξ χ = 0 . N ∗ = 2 yields a perfectly compatible results: a = 0 . γ = 0 . χ / dof = 5 . / N data also in narrower ranges, in particularconsidering only data for N >
N >
4. Best fitsare shown in Fig. 11. In both cases, N ∗ and γ turned N . . . . . ξ χ N > ξ χ ( N = 2) FIG. 10: Best fit of the small- N behavior of ξ χ for N rang-ing from 3 to 9. The dashed line represents the best fit ob-tained using fit function (27) leaving N ∗ as a free parameter,the shadowed band represents the fit error and the diamondpoints represent the hypothetical finite continuum limit ob-tained from the convergent continuum extrapolation of ξ χ for N = 2. out, again, to be compatible with 2 and 1 respectively: N ∗ = 1 . , γ = 0 . , N > ,N ∗ = 2 . , γ = 0 . , N > . Moreover, extrapolating these fits towards N = 3, weobtain values for ξ χ ( N = 3) which are, in both cases,in perfect agreement with the continuum limit obtainedfrom direct N = 3 simulations in Sec. 3 B. Based on theseresults, we conclude that the finiteness of the topologicalsusceptibility for N = 3 can be definitely assessed. E. Small- N behavior of the b n coefficients In this section we study the small- N behavior of the b n coefficients, in particular, we aim at checking the hy-pothesis that these coefficients are finite in the N → b , while already b always turnedout to be compatible with zero, after continuum extrap-olation, in all explored cases. Thus, we can presentlydiscuss only the small- N behavior of b .As we have discussed in Sec. 2 C, contrary to whathappens for the topological susceptibility, in this case wedo not expect a priori modifications to the standard formof UV corrections reported in Eq. (17). As a confirmationof this expectation, a continuum extrapolation performedaccording to the fit function in Eq. (22) turns out to workwell for all values of N , including N = 2 and N = 3, with O ( x ) corrections becoming irrelevant and not needed1 . . ξ χ N > N . . . . N >
FIG. 11: Best fits of the small- N behavior of ξ χ consideringjust N >
N > N ∗ as a freeparameter, while the shadowed bands represent the fit errors.Best fits give, respectively, N ∗ = 1 . χ / dof = 4 . / N ∗ = 2 . χ / dof = 3 / N = 3, both qualitatively and quantitatively. when restricting the fit range to large enough correlationlengths. Examples of such extrapolations are reported inFigs. 12, 13 and 14.For N = 2 and 3 we have also considered the possi-ble addition of generic power law corrections to Eq. (22)proportional to x c , as we have done for ξ χ , howeverthis turned out to be irrelevant in this case, with mod-ifications to the continuum extrapolation staying withinerrors. The only modification which can be noticed inthe N = 2 , ξ L available from our simulations in thesecases.All continuum extrapolations are reported in Tab. IX.Also in this case the quoted errors include systematics,which have been assessed by observing the variation ofcentral fit values when changing the fit function, the fitrange and the number of cooling steps n cool . Our datafor N = 2 confirm that b is finite, with a continuumestimate b ( N = 2) = − . ∼ σ off fromthe DIGA prediction b DIGA = − / ≃ − . N ∗ of N for which b reaches the DIGA prediction, to thataim we try to fit our data according to a critical behaviorlike: G ( N ) = b DIGA + a ( N − N ∗ ) γ . (28) .
000 0 .
005 0 .
010 0 .
015 0 . /ξ L − . − . − . − . − . b × − n cool = 10 n cool = 40 FIG. 12: Extrapolation towards the continuum limit of b for N = 5 using data in the range ξ L >
7. The solid anddashed lines represent the best fits obtained using fit function f ( x ) = a + a x (where x = 1 /ξ L ) with measures taken, re-spectively, after n cool = 10 and n cool = 40 cooling steps. Bestfits yield, respectively, χ / dof = 2 . / χ / dof = 3 . / .
00 0 .
25 0 .
50 0 .
75 1 .
00 1 . /ξ L × − − . − . − . − . − . b × − n cool = 25 n cool = 50 FIG. 13: Extrapolation towards the continuum limit of b for N = 3. The solid and dashed lines represent the best fits ob-tained using the fit function g ( x ) = a + a x (where x = 1 /ξ L )using measures taken after n cool = 50 cooling steps in theranges ξ L > ξ L >
15. The best fits yield, respectively, a = − . χ / dof = 4 . /
5) and a = − . χ / dof = 1 . / n cool = 20 coolingsteps in the range ξ L >
8. The best fit yields a = − . χ / dof = 2 . /
5. Data obtained for different numbersof cooling steps have been plotted slightly shifted to improvereadability. The diamond point represents our final estima-tion of the continuum limit. . . . . . /ξ L × − − − − − b × − DIGA prediction n cool = 20 n cool = 50 FIG. 14: Extrapolation towards the continuum limit of b for N = 2. The solid and dashed lines represent the bestfits obtained using the fit function g ( x ) = a + a x (where x = 1 /ξ L ) using measures taken after n cool = 50 coolingsteps in the ranges ξ L >
20 and ξ L >
30. The best fitsyield, respectively, a = − . χ / dof = 6 . /
5) and a = − . χ / dof = 2 . / n cool = 20 cooling steps in the range ξ L >
20. The bestfit yields a = − . χ / dof = 4 . /
5. Data ob-tained for different numbers of cooling steps have been plottedslightly shifted to improve readability. The diamond pointrepresents our final estimation of the continuum limit.
N b · b asa function of N . The value marked with * is taken fromRef. [13]. While we do not have any argument to support such anansatz, it turns out to work quite well: a best fit in therange [2 ,
9] is reported in Fig. 15 and returns the followingparameters: a = 0 . ,N ∗ = 1 . ,γ = 0 . ,χ / dof = 1 . / . It is interesting to observe that the N ∗ found in this wayfrom the analysis of b is perfectly compatible with the N − − − − b × − DIGA prediction
FIG. 15: Small- N behavior of the quartic coefficient b .The solid lines represent the best fit obtained using fit func-tion (28). The best fit yields N ∗ = 1 . χ / dof =1 . /
5. The black point represents the found interval for N ∗ ,for which b ( N ∗ ) = b DIGA . one found for the critical fit of ξ χ , and still compatiblewith N ∗ = 2.
4. CONCLUSIONS
In this work, we have presented a systematic numericalstudy of the peculiar features of the θ -dependence of thevacuum energy f ( θ ) of 2 d CP N − models in the small- N limit. To that aim, we have performed numerical simu-lations for N ∈ [2 , N , if any, at which the topological susceptibility χ di-verges: this is predicted to be N = 2 by perturbativecomputations of the instanton size distribution. A sec-ond interesting question regards the fate of the b n coef-ficients, which parametrize the non-quadratic part of θ -dependence and which could possibly approach the valuespredicted by the Dilute Instanton Gas Approximation atthe point where χ diverges, if the theory can be approx-imated in this case by a gas of small and non-interactinginstantons and anti-instantons.Our strategy has been twofold. On one hand, we havededicated particular efforts to correctly assess the con-tinuum limit of ξ χ for N = 2 and N = 3 from directsimulations of these theories: to that aim, we have per-formed simulations on lattices with correlations lengthsranging from a few units up to O (10 ). On the otherhand, we have exploited results obtained for larger val-ues of N , where the continuum extrapolation is easier, inorder to perform a small- N extrapolation. Based on thisdouble strategy, we have obtained consistent and conclu-sive evidence that the topological susceptibility is finite3for N = 3, providing the estimate ξ χ = 0 . N = 2 are still incon-clusive: results obtained directly at N = 2 are con-sistent with a logarithmically-divergent continuum ex-trapolation, but do not yet exclude a finite continuumvalue with ξ χ ∼ .
4, which is even marginally favoredfrom the point of the χ / dof test. A similar pictureemerges from the extrapolation from results obtained for N >
2, which provides evidence for a critical behavior ξ χ ∝ / ( N − N ∗ ) γ , with N ∗ = 1 . b coefficient, we have provided continuumextrapolations down to N = 2. While there is no com-pelling reason to expect that the DIGA prediction bevalid at the point where ξ χ diverges, it is interestingto observe that our numerical results are consistent with that. For N = 2 we obtain b ( N = 2) = − . σ off the DIGA value b DIGA = − / ≃− . N (see Eq. (28)) provides evidence for b reach-ing b DIGA for N = 1 . N ∗ at which ξ χ diverges reported above. Acknowledgments
The authors thank C. Bonati for useful discussions.Numerical simulations have been performed at the Scien-tific Computing Center at INFN-PISA and on the MAR-CONI machine at CINECA, based on the agreement be-tween INFN and CINECA (under projects INF19 npqcdand INF20 npqcd). [1] A. D’Adda, M. L¨uscher, and P. Di Vecchia, Nucl. Phys.B , 63 (1978).[2] M. Shifman,
Advanced topics in Quantum Field The-ory (Cambridge University Press, Cambridge, 2012), pp.171–268, 361–367.[3] E. Vicari and H. Panagopoulos, Phys. Rept. , 93(2009), 0803.1593.[4] E. Witten, Annals Phys. , 363 (1980).[5] E. Witten, Phys. Rev. Lett. , 2862 (1998), hep-th/9807109.[6] P. Rossi, Phys. Rev. D94 , 045013 (2016), 1606.07252.[7] C. Bonati, M. D’Elia, P. Rossi, and E. Vicari, Phys. Rev.
D94 , 085017 (2016), 1607.06360.[8] M. Berni, C. Bonanno, and M. D’Elia, Phys. Rev.
D100 ,114509 (2019), 1911.03384.[9] L. Del Debbio, G. M. Manca, H. Panagopoulos, A. Sk-ouroupathis, and E. Vicari, JHEP , 005 (2006), hep-th/0603041.[10] M. Campostrini and P. Rossi, Phys. Lett. B , 305(1991).[11] E. Witten, Nucl. Phys. B , 285 (1979).[12] G. Veneziano, Nucl. Phys. B , 213 (1979).[13] C. Bonanno, C. Bonati, and M. D’Elia, JHEP , 003(2019), 1807.11357.[14] A. Jevicki, Nucl. Phys. B127 , 125 (1977).[15] D. Forster, Nucl. Phys.
B130 , 38 (1977).[16] B. Berg and M. L¨uscher, Commun. Math. Phys. , 57(1979).[17] V. A. Fateev, I. V. Frolov, and A. S. Shvarts, Nucl. Phys. B154 , 1 (1979).[18] J.-L. Richard and A. Rouet, Nucl. Phys.
B211 , 447(1983).[19] M. D’Elia, F. Farchioni, and A. Papa, Nucl. Phys. B ,313 (1995), hep-lat/9505004.[20] M. D’Elia, F. Farchioni, and A. Papa, Phys. Rev.
D55 ,2274 (1997), hep-lat/9511021.[21] M. Blatter, R. Burkhalter, P. Hasenfratz, and F. Nieder-mayer, Phys. Rev.
D53 , 923 (1996), hep-lat/9508028.[22] Y. Lian and H. B. Thacker, Phys. Rev.
D75 , 065031(2007), hep-lat/0607026.[23] W. Bietenholz, U. Gerber, M. Pepe, and U.-J. Wiese, JHEP , 020 (2010), 1009.2146.[24] W. Bietenholz, P. de Forcrand, U. Gerber, H. Meja-Daz,and I. O. Sandoval, Phys. Rev. D , 114501 (2018),1808.08129.[25] D. Petcher and M. L¨uscher, Nuclear Physics B , 53(1983), ISSN 0550-3213.[26] W. Bietenholz, K. Cichy, P. de Forcrand, A. Dro-mard, and U. Gerber, PoS LATTICE2016 , 321 (2016),1610.00685.[27] M. Campostrini, P. Rossi, and E. Vicari, Phys. Rev. D , 2647 (1992).[28] S. Caracciolo and A. Pelissetto, Phys. Rev. D , 105007(1998), hep-lat/9804001.[29] M. Campostrini, A. Di Giacomo, and H. Panagopoulos,Phys. Lett. B212 , 206 (1988).[30] F. Farchioni and A. Papa, Phys. Lett. B , 108 (1993).[31] B. Berg and M. L¨uscher, Nucl. Phys. B , 412 (1981).[32] B. Berg, Phys. Lett. B , 475 (1981).[33] M. Campostrini, P. Rossi, and E. Vicari, Phys. Rev. D , 4643 (1992), hep-lat/9207032.[34] M. D’Elia, Nucl. Phys. B661 , 139 (2003), hep-lat/0302007.[35] Y. Iwasaki and T. Yoshie, Phys. Lett. B , 159 (1983).[36] S. Itoh, Y. Iwasaki, and T. Yoshie, Phys. Lett. B ,141 (1984).[37] M. Teper, Phys. Lett. B , 357 (1985).[38] E.-M. Ilgenfritz, M. Laursen, G. Schierholz, M. M¨uller-Preussker, and H. Schiller, Nucl. Phys. B , 693(1986).[39] M. Campostrini, A. Di Giacomo, H. Panagopoulos, andE. Vicari, Nucl. Phys. B , 683 (1990).[40] B. Alles, L. Cosmai, M. D’Elia, and A. Papa, Phys. Rev.D , 094507 (2000), hep-lat/0001027.[41] M. L¨uscher, Commun. Math. Phys. , 899 (2010),0907.5491.[42] M. L¨uscher, JHEP , 071 (2010), [Erratum:JHEP03,092(2014)], 1006.4518.[43] C. Bonati and M. D’Elia, Phys. Rev. D89 , 105005 (2014),1401.2441.[44] C. Alexandrou, A. Athenodorou, and K. Jansen, Phys.Rev.
D92 , 125014 (2015), 1509.04259. [45] H. Panagopoulos and E. Vicari, JHEP , 119 (2011),1109.6815.[46] C. Bonati, M. D’Elia, and A. Scapellato, Phys. Rev. D93 , 025028 (2016), 1512.01544.[47] L. Del Debbio, G. M. Manca, and E. Vicari, Phys. Lett.B , 315 (2004), hep-lat/0403001.[48] A. Laio, G. Martinelli, and F. Sanfilippo, JHEP , 089(2016), 1508.07270. [49] J. Flynn, A. Juttner, A. Lawson, and F. Sanfilippo(2015), 1504.06292.[50] M. Hasenbusch, Phys. Rev. D , 054504 (2017),1706.04443.[51] P. Rossi and E. Vicari, Phys. Rev. D48