Extracting Physics from Topologically Frozen Markov Chains
Abstract
In Monte Carlo simulations with a local update algorithm, the autocorrelation with respect to the topological charge tends to become very long. In the extreme case one can only perform reliable measurements within fixed sectors. We investigate approaches to extract physical information from such topologically frozen simulations. Recent results in a set of models and gauge theories are encouraging. In a suitable regime, the correct value of some observable can be evaluated to a good accuracy. In addition there are ways to estimate the value of the topological susceptibility.
Extracting Physics from Topologically Frozen Markov Chains
Christoph P. Hofmann
Facultad de Ciencias, Universidad de Colima
Bernal Díaz del Castillo 340, Colima C.P. 28045, Mexico
\abstract@cs
1 Introduction
In many relevant models, the configurations are divided into topological sectors (for periodic boundary conditions). This includes the models in , all the 2d models, 2d and 4d Abelian gauge theory, and 4d YangMills theories. The topology persists if we include fermions, hence this class of models also includes the Schwinger model, QED and QCD.
In the continuum formulation, a continuous deformation of a configuration (at finite Euclidean action) cannot change the topological charge . On the lattice there are no topological sectors in this strict sense, but at fine lattice spacing the configurations of the above models occur in distinct sectors with local minima, separated by boundary zones of higher action. Thus it is possible — and often useful — to introduce topological sectors also in lattice field theory, although the definition of is somewhat ambiguous. For the models that we are going to consider, the geometric definition [2] has the virtue that it naturally provides integer values of .
Most simulations in lattice field theory are performed with local update algorithms, such as the Metropolis algorithm for spin models, the heatbath algorithm for pure gauge theories, and the Hybrid Monte Carlo algorithm for QCD with dynamical quarks. If there are wellseparated topological sectors, such simulations may face a severe problem: a Markov chain hardly ever changes . Thus the simulation tends to get stuck in one topological sector for an extremely long computation time. Such a tremendous topological autocorrelation time was observed e.g. by the JLQCD Collaboration in their QCD simulations with dynamical overlap quarks [3]. For QCD simulations with nonchiral quarks (e.g. given by Wilson fermions) the problem has been less severe so far, i.e. for lattice spacings that have typically been used. However, in the future even finer lattices will be employed, and then this problem will become manifest.
So how can we measure the expectation value of some observable, , or the topological susceptibility if only topologically frozen simulations can be performed? ^{1}^{1}1 is the volume, and we will deal with parity symmetric models, where .
Lüscher suggested open boundary conditions, so can change continuously [4]. This overcomes the problem, but giving up integer has disadvantages, like losing the link to aspects of field theory in the continuum, e.g. regarding the regime of QCD.
Here we investigate approaches where periodic boundaries, and therefore , are preserved. In the framework of nonlinear models, we test methods to extract physical results from Markov chains, which are permanently trapped in a single topological sector, hence numerical measurements are available only at fixed . We start with a procedure to determine from the correlation of the topological charge density, which was introduced by Aoki, Fukaya, Hashimoto and Onogi [5]. Then we probe a way to assemble an expectation value from topologically restricted results . That approach is based on the BrowerChandrasekharanNegeleWiese (BCNW) formula [6], which also yields an estimate for .
2 Correlation of the topological charge density
Ref. [5] derived an approximate formula for the correlation of the topological charge density , at topological charge and large separation (we now use lattice units),
(2.0) 
The derivation assumes to be large, and to be small.^{2}^{2}2Actually this formula also involves a kurtosis term (which vanishes for Gaussian distributions). However, its contribution is negligible in all examples that we considered, so here we skip that term. Therefore we will limit our considerations to the sectors with . We are going to consider the 1d model and the 2d model, and the explicit condition for will be tested.
In our simulations we use the Wolff cluster algorithm [8], which performs nonlocal cluster updates. Hence we can also measure directly, which is useful for testing this method in view of other models (in particular gauge theories), where no efficient cluster algorithm is available. Preliminary results were anticipated in Ref. [9], and Ref. [10] presented before a related study (with different densities) in 2flavour QCD.
The 1d model, or quantum rotor, describes a free quantum mechanical scalar particle on the circle . We use periodic boundary conditions in Euclidean time over the size . The continuum formulation deals with an angle , where . The lattice variables are the angles , , with . We define the nearest site difference as
(2.0) 
i.e. the modulus function acts such that it minimises the absolute value. This yields the (geometrically defined) topological charge density and the topological charge ,
(2.0) 
We now give the continuum action and the three lattice actions — standard action, Manton action [11] and constraint action [12] — that we studied,
(2.0) 
The parameter corresponds here to the moment of inertia, and is the constraint angle. The continuum limit is attained at and , respectively. In the limit , and the correlation length are known analytically for all four actions in eqs. (2.0) [7, 12].
Figures 1 and 2 show results for the standard action, the Manton action and the constraint action at different sizes and parameters and . They are all in excellent agreement with the prediction based on eq. (2) (horizontal lines), even down to .
Next we address the 2d model, or Heisenberg model, on square lattices of size , with classical spins . Here we simulated the standard action and the constraint action, which are analogous to the formulations (2.0),
(2.0) 
Also here we use the geometric definition of the topological charge, which is written down explicitly in Ref. [12].
Figure 3 shows results for the topological charge density correlation. Again the comparison to the prediction (2) works in all cases (in this context we don’t have to worry about the fact that diverges logarithmically in the continuum limit). However, we also observe that for increasing volume it becomes soon difficult to resolve a clear signal from the statistical noise (as required for the determination of ), even with the huge statistics provided by the cluster algorithm. The quantitative results will be given in Ref. [14].
Thus we confirm that formula (2) is a valid approximation over a broad set of parameters. Nevertheless, in view of 4d quantum field theory its application is not promising, since for large volume it becomes very statistics demanding. That limitation is in agreement with the conclusion of an earlier study in the 2flavour Schwinger model with dynamical overlap hypercube fermions (and with the plaquette gauge action) [13]: at , and fermion masses a statistics of configurations in one topological sector was insufficient to determine in this way; to achieve this to about 2 digits would take at least configurations.
3 Applications of the BCNW formula
We now turn to the more ambitious goal of evaluating an observable , when only some values — at various and volumes — are available. To this end, we use an approximate formula, which was derived in Ref. [6],
(3.0) 
Our input are measured values for the lefthandside in various and , and a fit determines , and , where the former two are of interest. We refer to a regime of moderate , where these three quantities practically take their infinitevolume values, but the are still well distinct.
This is the beginning of an expansion in , hence should be large, but what that means has to explored numerically. Moreover the assumption of a small value of is involved again (see also the rederivation in Ref. [13]), hence we only use sectors with . An extended expansion has been mentioned in Ref. [6] (first work), and explored in great detail in Refs. [15, 16, 17]. It involves further free parameters, and the crucial question if this improves the results for and was addressed in Refs. [15, 16, 17], and will be discussed further in Ref. [14].
As our observables we consider the action density and the magnetic susceptibility (where is the magnetisation, and ). Results for the 2d model in are shown in the plots of Figure 4, which reveal the aforementioned regimes of “moderate ”. The fitting results involving the sectors , and various ranges in in those regimes, are given in Table 1. In particular we see an impressive precision of the values for , and also the fitting results for and are quite good.
Standard action  directly measured  
fitting range for  in all sectors at  
1.24038(12)  1.24027(8)  1.24015(5)  1.24008(5)  
0.0173(6)  0.0169(5)  0.0164(5)  0.01721(4)  
Constraint action  directly measured  
fitting range for  in all sectors at  
36.56(4)  36.58(3)  36.57(2)  36.57(2)  
0.00262(17)  0.00256(16)  0.00259(14)  0.002790(5) 
4 Conclusions
In simulations with local update algorithms and fine lattices, the Monte Carlo history tends to be confined to a single topological sector for an extremely long (simulation) time. This rises questions about the ergodicity (even within one sector). Here we do not address this conceptual issue; we trust the topologically restricted measurements of , and try to interpret them physically.
In very large volumes , the restricted expectation values all coincide with the physical result, , cf. eq. (3.0), but in practical simulations such large volumes are often inaccessible. For smaller , where is well converged to its large limit, but the are still significantly distinct, the BCNW formula (3.0) often allows us to determine to a good accuracy, and it also provides useful results for . It is favourable to employ only the sectors with , and (roughly speaking) the method is successful if .
If we relax that requirement to , we can still measure from the topological charge density correlation . The (theoretical) condition of a large separation turns out to be harmless in practice, but for increasing the wanted signal decreases very rapidly. Therefore that method is hardly promising for 4d models, where — in reasonable volumes — the signal would most likely be overshadowed by statistical noise.
On the other hand, the BCNW formula is promising for applications
in QCD, where typical simulations take place at .
This observation is supported by studies in the Schwinger model
[18, 13, 16], the quantum rotor with a potential [15]
and in 4d gauge theory [17], which will be
reported in detail in Ref. [14].
Acknowledgements: We are indebted to Christopher Czaban, Arthur Dromard, Lilian Prado and Marc Wagner for valuable communication and collaboration. This work was supported by the Mexican Consejo Nacional de Ciencia y Tecnología (CONACyT) through project 155905/10 “Física de Partículas por medio de Simulaciones Numéricas”, as well as DGAPAUNAM. The simulations were performed on the cluster of the Instituto de Ciencias Nucleares, UNAM.
References
 [1]
 [2] B. Berg and M. Lüscher, Nucl. Phys. B 190 (1981) 412.
 [3] H. Fukaya et al., Phys. Rev. Lett. 98 (2007) 172001; Phys. Rev. D 76 (2007) 054503.
 [4] M. Lüscher, JHEP 1008 (2010) 071. M. Lüscher and S. Schaefer, JHEP 1107 (2011) 036.
 [5] S. Aoki, H. Fukaya, S. Hashimoto and T. Onogi, Phys. Rev. D 76 (2007) 054508.
 [6] R. Brower, S. Chandrasekharan, J.W. Negele and U.J. Wiese, Nucl. Phys. B Proc. Suppl. 106 (2002) 581; Phys. Lett. B 560 (2003) 64.
 [7] W. Bietenholz, R. Brower, S. Chandrasekharan and U.J. Wiese, Phys. Lett. B 407 (1997) 283.
 [8] U. Wolff, Phys. Rev. Lett. 62 (1989) 361.
 [9] I. Bautista, W. Bietenholz, U. Gerber, C.P. Hofmann, H. MejíaDíaz and L. Prado, arXiv:1402.2668 [heplat].
 [10] S. Aoki et al., Phys. Lett. B 665 (2008) 294.
 [11] N. Manton, Phys. Lett. B 96 (1980) 328.
 [12] W. Bietenholz, U. Gerber, M. Pepe and U.J. Wiese, JHEP 1012 (2010) 020. W. Bietenholz, M. Bögli, U. Gerber, F. Niedermayer, M. Pepe, F.G. RejónBarrera and U.J. Wiese, arXiv:1309.6278 [heplat].
 [13] W. Bietenholz, I. Hip, S. Shcheredin and J. Volkholz, Eur. Phys. J. C 72 (2012) 1938.
 [14] I. Bautista, W. Bietenholz, C. Czaban, A. Dromard, U. Gerber, C.P. Hofmann, H. MejíaDíaz and M. Wagner, in preparation.
 [15] A. Dromard and M. Wagner, arXiv:1309.2483 [heplat].
 [16] C. Czaban and M. Wagner, arXiv:1310.5258 [heplat].

[17]
A. Dromard and M. Wagner,
arXiv:1404.0247 [heplat].
C. Czaban, A. Dromard and M. Wagner, Acta Phys. Polon. Supp. 7 (2014) 3, 551.  [18] W. Bietenholz and I. Hip, PoS LATTICE2008 (2008) 079; J. Phys. Conf. Ser. 378 (2012) 012041.