Lattice gluodynamics computation of Landau-gauge Green's functions in the deep infrared

We present recent results for the Landau-gauge gluon and ghost propagators in SU(3) lattice gluodynamics obtained on a sequence of lattices with linear extension ranging from L=64 to L=96 at $\beta = 5.70$, thus reaching"deep infrared"momenta down to 75 MeV. Our gauge-fixing procedure essentially uses a simulated annealing technique which allows us to reach gauge-functional values closer to the global maxima than standard approaches do. Our results are consistent with the so-called decoupling solutions found for Dyson-Schwinger and functional renormalization group equations.


Introduction
The infrared behaviour of gauge-variant Green's functions of Yang-Mills theories has become an increasingly interesting topic during the last decade. The interest was originally stimulated by confinement scenarios proposed a long time ago by Kugo and Ojima [1,2], Gribov [3] and Zwanziger [4,5]. The recent progress is due to the discovery of consistent asymptotic solutions of the whole tower of Dyson-Schwinger (DS) equations and more recently of functional renormalization group (FRG) equations in the deep infrared (IR) limit [6,7,8,9,10,11,12]. These solutions called scaling or conformal solutions behave according to power laws with well-determined exponents [8,9,11]. The expectation is that they respect global BRST invariance. Since these solutions lead to a vanishing gluon propagator and correspondingly to an IRsingular ghost dressing function they fit nicely with the aforementioned scenarios. The running coupling related to the ghost-ghost-gluon vertex exhibits an IR fixed point [6] as also proposed by Shirkov and Solovtsov [13,14].
There is a different set of so-called decoupling solutions as proposed in [15,16,17,18,19]. These solutionsthoroughly discussed also in [12] -are characterised by a non-zero IR gluon propagator as well as by an IR-finite ghost dressing function, i.e., they do not agree with the Kugo-Ojima criterion. The name decoupling reflects the fact that the corresponding running coupling decreases towards zero in the limit of vanishing momenta. But this behaviour does not mean that the decoupling solutions contradict gluon and quark confinement (see also the discussion in [20]). Both sets of solutions demonstrate the expected positivity violation of the gluon propagator as well as providing the expected Polyakov loop behaviour at the deconfinement transition in pure SU (2) and SU (3) gauge theories [21]. Moreover, both types of solutions, when interpolated from the infrared asymptotics to the perturbative region by solving numerically the (properly truncated) system of DS or FRG equations, behave quite similarly in the momentum range relevant for hadron phenomenology.
Even if it might appear to be an academic question, it is important to ask which set of solutions is the correct one. Ab-initio lattice gauge theory computations are expected to solve this issue. SU (2) and SU (3) lattice computations of Landau-gauge gluon and ghost propagators have been carried out by several groups. For a recent review see [22] and papers cited therein. On the one hand, as long as in the four-dimensional case the linear lattice sizes did not reach far beyond O(5 fm) the ghost dressing function was observed to rise towards the infrared limit. When fitted with a power law the corresponding exponent was found to be much smaller than predicted by the scaling solutions. On the other hand, the gluon propagator has not yet been found to decrease towards the IR limit but rather approaching a finite plateau value at p = 0 (for SU (3) see [23,24,25]). However, investigations of the DS equations on a 4D torus [26] have demonstrated that linear box sizes of O(10 fm) or even more, might be necessary in order to see the correct asymptotic behaviour.
In the meantime simulations of SU (2) pure lattice gauge theory have reached huge lattice sizes [22,27] confirming nothing but an IR-plateau behaviour for the gluon propagator. Moreover, in [28] Cucchieri and Mendes derived infrared bounds for the gluon and ghost propagators indicating, e.g., a non-zero value for the zero-momentum gluon propagator in the infinite-volume extrapolation.
However, the aforementioned results of DS and FRG studies as well as most of the lattice results obtained with the standard minimal Landau gauge did not take into account the effect of Gribov copies. For SU (2) some of us have shown [29,30,31,32,33] that the influence of Gribov copies is more severe than many authors might have expected (see also the recent investigation in [34]). The effect is even stronger if non-periodic Z(2) gaugetransformations (flips of all link variables orthogonal to fixed 3D sheets with a factor -1) are taken into account. Already on modest lattice volumes the gluon propagator as well as the ghost dressing function were seen to run into IR plateaus indicating that lattice results seem to support the decoupling solution of DS and FRG equations [33].
In the SU (3) case the effect of Gribov copies on the Landau gauge gluon and ghost propagators has been studied without Z(3) flips so far [23,35,36] showing that only the ghost propagator seems to be systematically affected within 5 to 10 %.
In the meantime, we have extended the computation of Landau gauge gluon and ghost propagators in SU (3) gluodynamics to linear lattice sizes of around O(16 fm). Preliminary results for lattice sizes up to 80 4 have been presented already in [37] clearly showing that an infrared plateau of the gluon propagator evolves and demonstrating -to our knowledge for the first time -a flattening of the ghost dressing function as well. Here we go a step further in order to get more confidence by increasing the lattice size for the gluon propagator up to as much as 96 4 as well as by increasing considerably our statistics for both propagators. In our investigation we have put some emphasis on careful gauge fixing with the simulated annealing method which wins in efficiency in comparison with standard over-relaxation, the more degrees of freedom the system has (see below and discussions in [31]).

General setup
We compute the SU (3) gluon and ghost propagators with Monte Carlo (MC) techniques on a lattice with periodic boundary conditions. The standard Wilson singleplaquette action and the lattice definition for the gauge potentials are adopted. In order to fix the Landau gauge for each lattice gauge field {U } generated by means of the MC procedure, the gauge functional is iteratively maximised with respect to a gauge transformation {g x } which is taken as a periodic field as well. In order to approach the global maximum (related to the fundamental modular region) as closely as possible, we use the simulated annealing (SA) algorithm [38], in combination with subsequent standard over-relaxation (OR). The latter is applied in the final stage of the gauge-fixing procedure in order to finalise the transformation to any required precision of the transversality condition ∇ µ A µ = 0. More than a decade ago, SA was shown to be very efficient when dealing with the maximally Abelian gauge [39,40].
In case of the present approach the SA algorithm generates a field of gauge transformations g x by MC iterations with a statistical weight proportional to exp(F U [g]/T ). The temperature T is an auxiliary parameter which is gradually decreased in order to maximise the gauge functional F U [g]. In the beginning, the T chosen must be large enough to allow the traversing of the configuration space of g x fields in large steps. An initial value T init = 0.45 was found to be sufficient for that [41]. After each sweep, including one heatbath and four micro-canonical update steps at each lattice site, T is decreased until g x is captured in a particular basin of attraction. We choose the lowest temperature value to be T final = 0.01 and use a finetuned simulated annealing schedule before applying final over-relaxation steps to reach Landau gauge with good ). An infinitely slow and long simulated annealing process would definitely lead to the global maximum. In practice, however, we chose schedules according to which the temperature is reduced in steps of finite but varying size. These schedules are motivated by the following observation. When the thermal average F U [g] is considered as a function of T , one easily notices its monotonous but non-linear dependence. In fact, F U [g] increases when T decreases, but its gradient reaches a strong maximum at a "critical" value T ≃ 0.405 ± 0.01, which resembles very much a phase transition [41]. This means that after starting with small ∆T steps, as required for proper ther-malisation close to T init = 0.45, it is important to keep the step size especially small within the narrow critical region T ≃ 0.40 . . . 0.41. Indeed, the distribution of the final gauge-functional values (those of the local maxima after simulated annealing and over-relaxation) is shifted to noticeably higher values if the temperature step size ∆T is reduced according to the rise of F U [g] , i.e., if ∆T is taken relatively small within the critical interval. This is in particular an efficient approach, as almost no improvement was observed, when further reducing ∆T outside this interval [41].
As an example, Fig. 1 shows the step sizes we use for the different temperature intervals on a 96 4 lattice, where, especially around T ≃ 0.40 . . . 0.41, ∆T is taken to be very small. Such extremely long SA "chains", with an order O(10 4 ) of iterations, allow us with high probability to reach local maxima of the gauge functional F U [g] close to the global maximum (as far as possible for the computer resources available), i.e., to the fundamental modular region, even with only one gauge-fixing attempt ("first-copy approach").
Note that we do not apply here Z(3) flips, the so-called FSA method [30,32] which in principle is capable of providing even larger F U [g] values than our SA algorithm does. But we have seen the difference between SA and FSA results decreases when L increases [32].
The computations presented are carried out at rather strong coupling, β = 6/g 2 0 = 5.70. The reason for this choice is to get access to comparatively large physical volumes at the price of a rather coarse lattice (a ≈ 0.17 fm). In order to study the volume dependence we calculate the gluon propagator for linear lattice sizes L ranging from 64 to 96. Thus, our largest lattice size corresponds to (16 fm) 4 . The ghost propagator is studied for linear lattice sizes L = 64 and L = 80. At larger lattice sizes and the lowest momenta the inversion of the Faddeev-Popov matrix turns out to converge only in rare cases. Obviously, this is due to the occurrence of very small eigenvalues which generate some algorithmic problems. There, a modification of the used matrix inversion method (see Sec. 4) or an even better gauge fixing, driving configurations further away from the Gribov horizon [36], could be valuable means of reducing such problems in future, but it is beyond the scope of the current work.

Gluon propagator
The gluon propagator is defined by where whereas at q = 0 the zero-momentum propagator D(0) is defined as Our data for the gluon propagator obtained for various lattice sizes is presented 1 in Fig. 2. One can clearly see a flattening of the gluon propagator as a function of q 2 for small momenta. Note also the weak volume dependence of the results. To illustrate the latter, in Fig. 3 we present also the dependence of the zero-momentum propagator D(0) (acc. to Eq. (5)) on the inverse lattice size 1/L. From this point of view it is natural to conclude that in the infinite volume limit the gluon propagator will approach some constant value as q 2 → 0, i.e., the gluon dressing function Z(q 2 ) = q 2 D(q 2 ) seems to decrease linearly with q 2 . The observed IR behaviour related to the zero-momentum modes of the lattice gauge field A µ (x) can be associated with a massive gluon. Our SU (3) IR gluon plateau results are in close agreement with analogous SU (2) results found recently on huge 4D symmetric lattices [22,27].

Ghost propagator
The Landau-gauge ghost propagator is defined by where M denotes the lattice Faddeev-Popov operator, being the Hessian of the gauge functional (2) with respect to g x , in the background of the gauge-fixed links U xµ with To invert M we use the conjugate gradient (CG) algorithm with plane-wave sources ψ c with colour and position components ψ a c (x) = δ ac exp(2πi k·x/L). In fact, we apply a pre-conditioned CG algorithm (PCG) to solve M ab xy φ b (y) = ψ a c (x) where as pre-conditioning matrix we use the inverse Laplacian operator ∆ −1 with diagonal colour substructure (for details see [23,42]).
For the large lattice sizes as considered here, we are confident that finite-volume distortions for all lattice momenta besides the two minimal ones do not change considerably with increasing L (see Fig. 4 and [37] for details). In this figure the ghost dressing function J(q 2 ) = q 2 G(q 2 ) is presented in a log-log scale. We do not see any powerlike singular behaviour in the limit q 2 → 0. Instead, we have a good indication that J(q 2 ) reaches a plateau just as the decoupling solution of DS and FRG equations does (see also [28,33]).

Running coupling
Finally, let us present the running coupling defined as the renormalization group (RG) invariant product of the Landau-gauge gluon and ghost dressing functions. This definition is based on the ghost-gluon vertex in a momentum-subtraction scheme with the vertex renormalisation constant (in Landau gauge) set to one. This is possible [43], since the vertex is known to be regular [44] (see also the lattice studies [45,46]). Note that the relation of α s in this scheme to the running coupling in the MS scheme is known to four loops and it can provide a valuable alternative to the MS coupling in phenomenological applications [43]. Beyond perturbation theory, the behaviour of α s differs at low scales for the scaling and decoupling solutions. Based on our propagator data we can calculate α s for intermediate and lower scales, and it clearly shows a decrease towards q 2 → 0 (see Fig. 5). This is again consistent with the decoupling DS and FRG solutions.

Conclusions
The progress achieved on the lattice during the last two years in studying the IR limit of gluodynamics and checking the well-known scenarios of confinement in terms of Landau-gauge Green's functions leads us to the following conclusions. Within the standard lattice approach as described above only the decoupling-type solution of DS and FRG equations seems to survive. Since for this solution the gluon propagator tends to a non-zero IR value, it corresponds to a massive gluon. It has been argued that this behaviour contradicts global BRST invariance [12]. But the lattice approach as discussed here has a few weak points. The choice of the gauge potentials A µ (x) and correspondingly of the gauge functional F U (g) is in no way unique. As long as we are reaching the infrared limit by employing rather large lattice couplings the continuum limit is not under control. Moreover, we have used standard periodic boundary conditions which certainly have an impact on the IR limit. The fact that under these conditions the gluon propagator does not tend to zero is related to the behaviour of the zero-momentum modes, which do not become sufficiently suppressed as the lattice size increases. Changing the definition of lattice Landau gauge, and correspondingly the lattice definitions of A µ (x) and M , modifying the boundary conditions and further improving the gauge-fixing procedure, e.g., by taking Z(3) flips as mentioned in Sec. 1 into account, may essentially suppress the zero-momentum modes and correspondingly change the behaviour of both the gluon and ghost propagators. Therefore, a final conclusion still cannot be drawn.
Note that one of us (A.S.) has recently carried out a lattice computation in the strong-coupling limit. For the gluon and the ghost propagator at larger a 2 q 2 it was possible to extract the right exponents as expected for the scaling solution [47]. At asymptotically small momenta, however, results were shown to depend strongly on the lattice definitions of A µ (x) and M . This corresponds to observations in studies of DS and FRG equations, namely that for decoupling-like solutions any IR-asymptotic values of the gluon propagator and the ghost dressing function can be considered as boundary conditions [12].
It has been argued in [48] that a BRST-invariant gaugefixing prescription is possible on the lattice. It remains to be seen, whether the preferred scaling behaviour of Landau-gauge gluon and ghost propagators can be achieved consistently also in lattice Yang-Mills theories.