Investigating the topological structure of quenched lattice QCD with overlap fermions by using multi-probing approximation

The topological charge density and topological susceptibility are determined by multi-probing approximation using overlap fermions in quenched SU(3) gauge theory. Then we investigate the topological structure of the quenched QCD vacuum, and compare it with results from the all-scale topological density, the results are consistent. Random permuted topological charge density is used to check whether these structures represent underlying ordered properties. Pseudoscalar glueball mass is extracted from the two-point correlation function of the topological charge density. We study $3$ ensembles of different lattice spacing $a$ with the same lattice volume $16^{3}\times32$, the results are compatible with the results of all-scale topological charge density, and the topological structures revealed by multi-probing are much closer to all-scale topological charge density than that by eigenmode expansion.


Introduction
The topological susceptibility χ of Yang-Mills theory is an important physical quantity in explaining the U (1) A anomaly and the mass of η , e.g. the study of this problem in N c → ∞ limit brings up the famous Witten-Veneziano formula for η mass [1][2][3]. In the continuum limit, the topological susceptibility as a function of fourmomentum k, denoted by χ(k 2 ) in Euclidean version, is defined as here the quantity C qq (x) = 0|T (q(x)q(0))|0 is the twopoint correlator of the topological charge density q(x) and |0 denotes the QCD vacuum.
The contribution of the U (1) A anomaly to the η mass is given by [1,2,4] in which Y M is for pure Yang-Mills theory, N f is the number of flavors, and F π is the corresponding π decay constant, and will leads to χ Y M (0) ≈ (180MeV) 4 [2,4].
To calculate q(x) and χ(0) of SU(3) pure Yang-Mills theory on lattice, we have to define the topological density q(x) first. There are several ways to define it. One usually way is to define q(x) with the field strength tensor F µν on the lattice as [5] in which subscript L denotes lattice. This definition has the classical continuum limit q L (x)  lattice, since the direct calculated topological charge Q L = x q L (x) usually deviates from the integral values, and Q L = a 4 Z(β)Q + O(a 6 ) [6][7][8], with the multipliciative renormalization coefficient Z(β) dependent on the inverse coupling constant β. This dependence can be computed perturbatively, or some kind of smoothing steps could be used to make the effective inverse coupling constant β eff larger and such that Z(β ef f ) → 1 [8] is achieved. As to χ(0), it was found that multiplicative and additive renormalization terms both exist as χ L (0) = a 4 Z 2 (β)χ(0)+a 4 A(β) < S(x) > n.p. +P (β) [6], S(x) is the action density, n.p. indicates non-perturbative, A(β) and P (β) are ultra-violet effects that can be approximated in perturbation theory. We can also smooth configurations to approach integer value of Q L for calculating χ(0). But the smooth procedure would influence q(x).
Another way to define q(x) is using massless Dirac operator D(0, x, y) that obeys the Ginsparg-Wilson relation Dγ 5 + γ 5 D = aDγ 5 D. Here the first argument in D(m, x, y) denotes the mass of the quark while the second and third index are the space-time index on the lattice. For simplicity, we will also use D(0) to denote the massless Dirac operator with the space-time indices suppressed. With this notation, the topological charge density can be defined as where c and d being the index of color and Dirac spinor respectively. This is a proper approach that guarantees the topological charge Q L an integral value since Ginsparg-Wilson fermion satisfies Atiyah-Singer index theorem [9] in which n L and n R are number of left-and right-handed zero modes of D(0). Massless overlap fermion is exactly Ginsparg-Wilson fermion [10,11], so in this paper we will use massless overlap Dirac operator D ov (0) to calculate the topological charge density. The definition of D ov (0) reads, D w (s) is the Wilson-Dirac operator, with parameter |s| < 1 for the optimization of the locality of overlap operator [12], D w (s) = D w − (1 + s), D w is the massless Wilson-Dirac operator with mass parameter m = −(1+s).
Nevertheless, it will cost too much computer resources to calculate topological charge density using the definition in Eq. (4) directly, i.e., using point sources to calculate the trace, the resulting topological charge density q(x) is named all-scale topological charge density in Ref. [13]. An approximation method called eigenmode expansion was proposed in the literature [13][14][15]: , ψ λ is the eigenmode of D(0) with eigenvalue λ, λ cut is the supremum of the absolute eigenvalue of the corresponding eigenmodes included in the eigenmode expansion.
Only if all zero modes are included, which obey the relationships γ 5 ψ λ = ±ψ λ , D ov (0)ψ λ = 0, q λ cut (x) would sum the individual local chiralities, and the topological charge of eigenmode expansion Q λ cut = Q L . However, some features such as the positive core of two-point correlator of topological charge density C qq (x) would increases, and the height of peaks of C qq (x) would decrease, even the negative region of C qq (x) would disappear when λ cut is too small [13], therefore more eigenmodes should be considered which would also cost too much computer time, and it's not sure what has been discarded in the eigenmodes of larger λ cut .
In this paper we will use multi-probing method [16] to approximate q(x) owing to the finite range of D(0) [13,17], δQ, the deviation of Q L from an integer would be regarded as a parameter that signals the validity of this approximation. The lattice version of Eq. (1) is defined as with the lattice version of two-point correlation function C qq (r) = 1 V x∈V q(x + r)q(x) . Noting that the true value of expectation of topological charge is 0, the lattice version of topological susceptibility at k 2 = 0 can be expressed as obviously it has the same accuracy of topological charge. It is believed that the U (1) A problem [18], topological susceptibility, η mass, θ dependence and spontaneous chiral symmetry breaking are all related to the vacuum fluctuation of topological charge density [1,19]. So in Sect. 3 we attempt to investigate the local topological structures of the quenched QCD vacuum, hoping that this will reflect some ordered features of the quenched QCD vacuum. Since they can be explained by instantons model, which was first introduced by Polyakov [20,21], it is reasonable to try to describe the structure of topological charge density in terms of instanton liquid model, the topological charge density shows sign-coherent unitquantized lumps locally. In the past, several papers [22][23][24][25][26][27][28][29][30] have attempted to justify this idea, but these papers used methods such as cooling, smoothing and smearing, which will lead to configurations towards classi-cal solutions, therefore the results should be taken with some reservations. On the other hand, it was predicted that QCD vacuum doesn't resemble the instanton liquid model [19], and this idea was further supported by lattice evidence in Ref. [14]. In fact, there is no four dimension global structure of topological charge density but only lower dimensional structures survive [13][14][15]. In Sect. 3 we will discuss more details about the structure of topological charge density, and the results of all-scale q(x) in Ref. [13] will be compared with, which were built up from three ensembles only with 53, 5 and 2 configurations, due to the large amount of computer time for calculating all-scale topological charge density.

Multi-probing approximation
The massless overlap operator D ov (0) is not ultralocal but local. It is exponentially decaying with respect to lattice distance with a decay rate that has little independent on inverse coupling constant β [13,17] D ov (0, x, x + r) ∝ exp(−µr/a) , µ > 0, (9) in which r is the lattice Euclidean distance, a is the lattice spacing, so we have D ov (0, x, x + r) D ov (0, x, x) when r/a is large enough. Therefore we can approximate topological charge density q(x) with massless overlap fermions as in which {r i } are some distances satisfied exp(−µr i /a) 1, and c, d are indices for Dirac and color respectively. In actual computing, we choose sites on the lattice {x i } as a group of probes called multi-probe and these sites are distributed uniformly with interval δr, δr should be an integer number that is a common divisor of spacetime lengths of the lattice, or multiplied with √ 2 when the sites in the same group are separated by odd and even into two new groups [16]. The whole lattice sites are now divided into N mp groups,then we can produce 12 for each multi-probing group, which only has non-zero unity value on sites {x i } belong to the corresponding group. The element of (γ 5 D(0, x, x)) cd , x ∈ {x i } can now be probed by Thus, only 12N mp overlap operator multiplications with vector sources have to be performed with multi-probing method, which is to be compared with 12L x L y L z L t if calculating all-scale q(x) directly, saving a substantial amount of computations especially for large lattices.
In Tab. 1 we present the three ensembles of SU(3) pure gauge configurations in our study which are generated with periodic boundary using Iwasaki action. All have the same lattice size of N s = 16, N t = 32, and we use Wilson flow with w 0 = 0.1670(10)fm to set the scale of lattice spacing a [31][32][33]. In the calculations of topological charge density, we choose the interval δr = 4 √ 2a to approximate the all-scale topological density q(x), which ensures that almost all configurations have δQ < 0.1, in which δQ is the deviation of Q = x q mp (x) from an integer, since the massless overlap operator D ov (0) should satisfy Atiyah-Singer index theorem. Therefore only 4 4 × 2 × 3 × 4 = 6144 times of matrix multiplication to vector are needed, as for calculating all-scale q(x) directly, 16 3 × 32 × 3 × 4 = 1572864 times will be needed. In column 5 and 6 we list the percentage of the configurations that satisfy the condition δQ below some arbitrary cut values, 0.08 and 0.05 respectively, and the obtained values for χ 1 4 (0) using multi-probing method are tabulated in last column. It is noted that these values are compatible with phenomenological estimates In Fig. 1, we present topological charge density computed both by all-scale q(x) directly and multi-probing method with δr = 4 √ 2a on a 2D surface, Fig. 1(b) only differs a little bit from Fig. 1(a), the features of q(x) such as peaks, valleys and ridges are consistent, more detailed analyses of the structure of q(x) will be discussed in Sect. 3.2. The scheme of multi-probing with interval δr = 4 √ 2a is presented on 2D surface in Fig. 2(a), sites with same number belong to the same multi-probing group, and 3D sketch Fig. 2(b) shows that the diagonal points, which are on the surfaces of cube with side length equals 4a in the lattice, belong to the same multi-probing group. As to the scheme of multi-probing with interval δr = 4a, vertices of cube with side length equals 4a should belong to the same multi-probing group, in other words, sites with number 1 and 2 in Fig. 2(b) belong to the same group now, then the number of groups will be reduced by half. Therefore the schemes of multi-probing with interval δr = 4 √ 2a and δr = 4a are correlated strongly, same relationship exists between schemes with interval δr = n √ 2a and δr = na universally. 3 Analysing the topological charge density q(x) Noting that topological susceptibility can be computed from topological charge individually, and the topological charge Q should be an integer value, which can be used to check the degree of approximation. Although the multi-probing method seems to approximate Q rather well as indicated in Tab. 1, one is still concerned about the performance of topological charge density and especially the slope of the susceptibility. Thus, we need to find some other parameters for testing the validity of the approximation. From Eq. (??) we may compare the towpoint correlator C qq (r) from multi-probing and all-scale q(x). Other tools proposed in Ref. [13] will be utilized to analyze the vacuum structure of topological charge density, results will also be compared with that obtained from all-scale q(x).

C qq (r) from multi-probing method
The reflection positivity [34] requires that two-point correlator of topological charge density should be nega-tive [4,35,36]: Together with massless overlap operator D ov (0) has a finite range, C qq (r) with overlap fermions will have a positive core, which has been studied extensively already. In Figs. 3(a), 3(b) the two-point correlator from our work displays positive core and negative tails, just as all-scale q(x) does in Figs. 3(c), 3(d) ( results of Fig.14 in Ref. [13]). The details showed in Figs. 3(b), 3(d) demonstrate that the ranges of positive cores are consistent with β = 9.45 (a ≈ 0.1f m) in Fig. 3(b) and of β = 8.45, 8.60(a ≈ 0.1f m) in Fig. 3(d) with the same lattice size 16 3 ×32. From Fig. 3(b), ensemble with different lattice spacing a has different positive core in physical unit, but almost the same in lattice unit, which means it is independent of β in lattice unit, this effect comes from the fact that the finite range of overlap operator D ov (0) has little dependence on β. In other words, as lattice spacing a decreases, the positive core of C qq (x) will disappear in physical unit, when a → 0, only the contact term C qq (0) > 0, and the reflection positivity condition Eq. (12) will be restored in the continuum limit.
Inspired by the instanton model, a functional form of the two-point correlator C qq (x) in the negative region can be used to extract pseudoscalar particle mass [37]: where K 1 (z) is the modified Bessel function which has asymptotic form as We can use Eqs. (13) (14) to fit the data of C qq (r). In the fitting procedure the amplitude and mass are treated as free parameters and then we can extract mass of pseudoscalar glueball in quenched QCD, since only glueballs exist in pure gauge theory. The fitting range is deter-mined by the effective mass m ef f (r) which is defined by C qq (r + ∆r) C qq (r) = r r + ∆r where we set ∆r = 0.5a and C qq (r) will be averaged over [r − ∆r 2 , r + ∆r 2 ]. In Fig. 4 we present the effective mass plateaus of 2 ensembles in Tab. 1 with errors computed by Jackknife method, except ensemble a = 0.1161fm, whose lattice spacing may be too coarse. The starting and ending points of fitting range would be varied in [3.25a, 5.25a] for a = 0.0769fm ensemble and [3.5a, 4.5a] for a = 0.0963fm ensemble respectively, they are the plateaus shown in Fig. 4, and fitting range is set to be greater than or equal to ∆r = 0.5a, then the final fitting range is decided by Tab. 2 presents the fitted results for the masses for ensembles a = 0.0963, 0.0769fm, both of which are compatible with the value 2560(35)Mev from Ref. [39], though the fitting range is determined somewhat subjectively. It is also seen that the mass has little dependence on the lattice spacing a.

clusters of topological charge density from multi-probing method
We use cluster analysis proposed in Ref. [13] to investigate the sign-coherent structure of topological charge density (the word "sign-coherent" will be omitted below for simplicity): 1. The cluster on the lattice is composed of sites that link-connected of same sign of the topologi-  Fig.14(b) in Ref. [13]    cal charge density q(x), q cut means that only sites with |q| > q cut are taken into account, and q max is the maximum value of |q(x)|.

The connectivity is defined by
r max = L x 2 + L y 2 + L z 2 + L t 2 /2 is the largest distance available on the periodic lattice, Θ c (x) is a character function of the cluster c, namely Θ c (x) = 1 , if x ∈ c , else Θ c (x) = 0. The value q cut = q perc for which f (r max ) = 0 is called the onset of percolation. When q cut < q perc , f (r max ) > 0 indicating that global clusters exist.

The Euclidean distance between the two largest
clusters c and c is defined by : 4. Random walker method is utilized to define the fractal dimension d * of a cluster. Starting from the site that has the maximum q max inside a cluster, and after τ steps of independent random walkers, the probability of returning to the starting point P (τ ) = (2πτ ) −d * /2 . The walker jumps to linked site in the same cluster with equal probability every step.
Since the physics of vacuum is thought to be mostly affected by intense fields, and the super-long-distance might be relevant to the physics of the QCD vacuum such as Goldstone boson propagation [15], therefore we are concerned about the large topological charge density (q cut is large), at the same time global structures still survive, i.e. f (r max ) > 0. Left side from multi-probing method, right side from all-scale q(x) in Ref. [13]. 1st row: the total number of clusters. 2nd row: the size of the largest cluster relative to all clusters.
In Fig. 5, the number of clusters, which we call cluster multiplicity, increases as q cut increases till it reaches the maximum, then decreases; and the fractional size of the largest cluster, which denote the size of the largest cluster relative to total volume occupied by all clusters, reaches the minimum near the maximum point of cluster multiplicity , then rises. The trend of fractional size of the largest cluster is just opposite to cluster multiplicity. At q cut = 0 the two largest clusters with opposite sign both have fractional size of the largest cluster ≈ 0.5. Except the two largest clusters, the dozens of clusters left are fragments only with volume of 1 or 2, and their total volume compared with the lattice volume 16 3 ×32 are so small that can be ignored. Left side from multi-probing method, right side from all-scale q(x) in Ref. [13]. 1st row: the distance between the two largest clusters. 2nd row: the connectivity f (rmax) of the largest cluster.
In Fig. 6, the distance d(c, c ) ≈ √ 3a between two largest clusters at q cut = 0 is independent of β, like results shown in 6(b), and just as the positive core of twopoint correlator of topological charge density does. This means that at q cut = 0 the two largest clusters of opposite sign tangle with each other closely, when a → 0, they are just like two tangled 3D papers dominating the whole 4D space. After the onset of percolation, the distance between two largest clusters reaches a plateau, which is also independent of β, this had been proved in Fig. 6(b) from Ref. [13]. The turning point of cluster multiplicity, fractional size of the largest cluster and the plateau of the distance between two largest clusters all happen near the onset of percolation at q cut ≈ 0.2q max .
In Fig. 7 we present effective fractal dimension of the largest cluster using random walker algorithm, lower than 4 dimensional structure of topological charge density can be found. As the q cut increases, the effec-tive fractal dimension drops down, but even q cut at 0.20q max − 0.25q max , when percolation is just starting, global clusters with lower effective dimension d * = 1−1.4 have already appeared, which is consistent with result  from all-scale q(x) showed in Fig. 7(d), and this structure had been first found as one dimensional skeleton in Ref. [15].
The trends of the properties of clusters with respect to q cut /q max from our work are shown on left side of Figs. 5, 6, 7 are consistent with the trends shown on the right sides in the same corresponding figures which we quote from Ref. [13]. Most of the quantitative values are also compatible, except for the cluster multiplicity. Especially when q cut = 0 there are still dozens of clusters survive, meanwhile in all-scale q(x) there will be only two global clusters with opposite sign, and noting that the magnitude of cluster multiplicity in Fig. 5(a) is about twice in Fig. 5(b), which shows a scene much more chaotic than all-scale q(x) does, this extra chaos or clusters might come from the multi-probing approximation method itself. After doing more study on the data and considering Fig. 5(c) is consistent with Fig. 5(d), we found that these extra small clusters add to the total number of clusters, but they do not affect the fraction of the largest cluster to the total occupied volume substantially. Especially at q cut = 0, there are few small fragments composed of only one or two sites, the total volume of these fragments can be ignored with respect to the whole lattice space or the two largest clusters with opposite sign, which are the clusters that all-scale q(x) only have at q cut = 0.
Noting that in Figs. 5(b), 5(d), 6(b), 6(d), lattice ensembles at β = 8.10 with 53 configurations, at β = 8.45 with 5 configurations, and 2 configurations at β = 8.60, the numbers of configurations are so small owing to that the direct all-scale q(x) calculating consumes too much computer resource. All the 3 ensembles used in this work consist of 200 configurations.
In Fig. 8 we present the fractional size, which is separated by different sections of size of clusters, relative to the total occupied volume by clusters, and the clusters are ranked according to the volume of clusters. From Fig. 8 we can find that most volume are occupied by the first half largest clusters in three different ensembles below the onset of the percolation q cut = 0.2q max , then the fractional size of extra clusters should not be concerned. But according to Figs. 5, 6, 7, above the onset of percolation the extra clusters from multi-probing approximation still do not influence the role of the largest clusters substantially.  The average fractional size of largest clusters relative to occupied volume by total clusters, separated by sections of size of clusters: the first 2 largest clusters, the 3rd to 10% largest clusters, 10% to 30% largest clusters, 30% to 50% largest clusters, and the left 50%.
It is easy to understand the reason why multi-probing method will bring extra clusters. Because we use the clusters of topological charge density as tools to investigate the structures of the quenched QCD vacuum, with the hope that the structures would reflect ordered properties of q(x), and the multi-probing method will bring superposition from different probing sites to topological charge density, which can be regarded as disordered disturbance, then noise will be introduced to the ordered structure of clusters.
Noting that the properties of clusters of topological charge density from multi-probing method showed above are much closer to the all-scale topological charge density than that from eigenmode expansion showed in Figs. 20, 21 from Ref. [13].

Random permutation of topological charge density
Though we have multi-probing observed structures of topological charge density, it is still a question where these stuctures come from. For example, are these structures arise from the ordered properties of the quenched QCD vacuum or they can be established by disorders. To get some information about this, we permute the topological charge density randomly through space-time to generate stochastically distributed q(x), which should no longer contain information of the ordered properties of the quenched QCD vacuum, and then investigate the consequences of this permutation. Fig. 9 shows results of the permuted topological charge density. It is found that the cluster multiplicity becomes much larger than in Fig. 5(a). Even at q cut = 0 there are still hundreds of clusters. The onset of percolation is found near 0.15q max now, the fractional size of the largest cluster and the connectivity drop quickly, the distance are √ 2a, which is less than results above in Fig. 6(a). These results suggest that the clusters tangle with each other more closely. Above the onset of percolation, the distance between the two largest clusters drops down, when it should reach a plateau before. The effective dimension shown in Figs. 9(e), 9(f), 9(g) drops faster than those shown in Figs. 7(a), 7(b), 7(c) and become negligible at q cut = 0.4q max already. On the other hand, the changes shown in Fig. 9 also demonstrate that, the differences between results of multi-probing and allscale q(x) seen above really come from disordered disturbance on the lattice. Therefore we can be sure to say that the properties of clusters from multi-probing method do reflect the underlying ordered structure of QCD vacuum, especially the structure recognized between q cut = 0.15q max − 0.20q max where the structure of clusters is similar to all-scale q(x) but different from the permuted case.

Summary
In this paper, using the multi-probing method with overlap fermion to approach the topological charge density, the topological aspects of quenched QCD are investigated. To be specific, we study the topological charge, topological susceptibility and the structure of topological charge density in space-time. It is found that the cluster properties are consistent with the results obtained from all-scale topological charge density rather well. After comparing them with clusters of randomly permuted topological charge, we conclude that, the global structure that appear near the onset of percolation with effective dimension near 1 may well reflect the underlying ordered properties of the quenched QCD vacuum. Since there are still global clusters survive when the sites that have topological charge density of absolute value below q cut have been discarded.