Two-point correlation functions of QCD in the Landau gauge
aa r X i v : . [ h e p - t h ] S e p Two-point correlation functions of QCD in the Landau gauge
Marcela Peláez a,b , Matthieu Tissier a , and Nicolás Wschebor b a LPTMC, Laboratoire de Physique Théorique de la Matière Condensée,CNRS UMR 7600, Université Pierre et Marie Curie,boîte 121, 4 place Jussieu, 75252 Paris Cedex 05, France. b Instituto de Física, Facultad de Ingeniería, Universidad de la República,J. H. y Reissig 565, 11000 Montevideo, Uruguay. (Dated: September 9, 2014)We investigate the gluon, ghost and quark propagators in the Landau gauge with dynamic quarks.We perform a one-loop calculation in a model where the standard Faddeev-Popov Lagrangian iscomplemented by a mass term for the gluons which is seen as a minimal way of taking into accountthe effect of the Gribov copies. The analytic results are compared with lattice data obtained infour dimension and for two, three and four quark flavors. The gluon and ghost propagators arereproduced with a few percent accuracy in the whole range of accessible momenta. The scalarpart of the quark propagator is found to be in good agreement with the lattice data. However,the quark renormalization is poorly described. We attribute this discrepancy to the fact that theone-loop corrections to this quantity are unusually small so that the two loop contribution can notbe discarded. The results are expressed in terms of the coupling, the gluon mass and the light quarkmass at 1 GeV.
PACS numbers: 12.38.-t, 12.38.Aw, 12.38.Bx,11.10.KkKeywords:
I. INTRODUCTION
It is well-known that the Faddeev-Popov construc-tion, which is the standard analytic method for fixingthe gauge, is not sufficient in order to treat the infraredregime of QCD. The problem originates in the existenceof Gribov copies [1] which are ignored in the Faddeev-Popov construction. Gribov copies do not play any rolein the ultraviolet regime. For that reason, the Yang-Millsaction complemented by the Faddeev-Popov action is agood starting point to analyze both gauge-invariant andnon-invariant quantities at momentum scales much big-ger than GeV. However, at low momenta the question ofhow to include the effect of Gribov copies in a Lagrangianstill remains unsolved even though some ideas have beenalready developed. Among those, the most popular isprobably the Gribov–Zwanziger proposal [2–5]. It con-sists in restricting the functional integral to the first Gri-bov region by adding several new fields. Unfortunately itwas shown in [6] that the first Gribov region also includesmany Gribov copies so that the Gribov ambiguity is notcompletely removed. Moreover, the Gribov-Zwanzingerprocedure relies on some formal manipulations that arenot fully justified from first principles.In the past decades, lattice simulations improved con-siderably and became the most reliable technique to de-scribe the infrared region of QCD [7]. In particular,quenched lattice simulations in Landau gauge [8–14] andunquenched simulations [15, 16] have demonstrated thatthe gluon propagator acquires a finite value at low mo-mentum, contrary to the original belief. Numerical sim-ulations succeeded in convincing the community that, inthis gauge, the gluon propagator behaves in the IR as ifit was massive. Even before that, the infrared regime was studiedthrough several semi-analytical methods, among whichthe more popular were based on Dyson-Schwinger (DS)equations. Depending on the precise implementation,two solutions were observed, called scaling and decou-pling (or massive) solutions. The scaling solution is char-acterized by a gluon propagator that vanishes as a powerlaw at low momentum [4, 17–21]. On the other hand, thedecoupling solution corresponds to a finite gluon propa-gator [21–28], in agreement with lattice results.Acknowledging the fact that no analytic method candeal with the Gribov ambiguity in a fully consistent way,two of the authors have proposed an alternative strat-egy to study analytically the IR regime of the theory.The idea relies on working with the simplest Lagrangianwhich reproduces lattice results. To account for the de-coupling solution observed on the lattice in a minimalway, the Faddeev-Popov Lagrangian is complemented bya mass term for the gluons. In euclidean space, the La-grangian density reads: L = 14 F aµν F aµν + ∂ µ ¯ c a ( D µ c ) a + ih a ∂ µ A aµ + m A aµ A aµ + N f X i =1 ¯ ψ i ( γ µ D µ + M i ) ψ i , (1)where g is the coupling constant, γ µ are euclidean Diracmatrices satisfying { γ µ , γ ν } = 2 δ µ,ν , the flavor index i runs over the N f quark flavors and F aµν = ∂ µ A aν − ∂ ν A aµ + gf abc A bµ A cν , ( D µ c ) a = ∂ µ c a + gf abc A bµ c c ,D µ ψ = ∂ µ ψ − igA aµ t a ψ. The latin indices correspond to the SU ( N ) gauge group, t a are the generators of the algebra in the fundamentalrepresentation and f abc are the structure constants.The model described by Eq. (1) is a particular case ofthe Curci-Ferrari model [29]. It is well known that themass term violates the BRST symmetry of the Faddeev-Popov Lagrangian. However, it is still symmetric undera modified BRST symmetry responsible for its renormal-izability. This Lagrangian can also be motivated fromfirst principles by taking into account Gribov copies [30].A great advantage of the phenomenological model de-scribed by Eq. (1) is that it is very simple and allowsto perform perturbative calculations very easily. More-over, we found [31, 32] that there exist renormalizationschemes in which the coupling constant remains finite atall momentum scales. The absence of Landau pole allowsus to implement perturbation theory even in the infraredregime. The quenched one-loop calculations for the two-point [31, 32] and three-point [33] correlation functionscompare very well with lattice simulations. The previousmodel was also used for studying the two-point correla-tion functions at finite temperature [35] where it repro-duces at a qualitative level the properties of the gluonand ghost propagators.In this article, we pursue our systematic comparisonof the correlation functions obtained within the model(1) with those extracted from lattice simulations. We in-clude here dynamical quarks and compute, at one loop,the two-point correlation function for the gluon, ghostand quark, for arbitrary dimension, number of colors ( N )and flavors ( N f ). The rest of the article is organized asfollows. In Sect. II, we discuss the 1-loop calculationfor the gluon and ghost propagators which are expressedin terms of Passarino-Veltman integrals [36]. We presentour results in arbitrary dimension and give explicit ex-pressions in the physically relevant case d = 4 . Sect. IIIis devoted to the quark propagator. In Sect. IV we in-troduce the renormalization schemes and discuss the im-plementation of the renormalization group. In Sect. Vwe perform a comparison of lattice correlation functionswith our unquenched perturbative results in the gluonand ghost sector. We focus on N = 3 with two light de-generates quarks ( N f = 2 ) and with two light quarks andtwo heavier quarks ( N f = 2 + 1 + 1 ). In Sect. VI we workwith N f = 2 + 1 . The quark propagator was extractedfrom lattice simulations in this case, which enables us tomake a comparison with our analytical results. Finally, inSect. VII we estimate the two-loop contributions whichgives an indication of the error bars on our results. II. UNQUENCHED GLUON AND GHOSTPROPAGATOR
The gluon and ghost propagators have been extensivelystudied in lattice simulations, both in the quenched andunquenched case. It is found that the addition of thesea quarks does not change qualitatively the behavior of both propagators. It however tends to lower the plateauobserved at small momenta for the gluon propagator.We parametrize the gluon and ghost 2-point vertexfunctions in the following way Γ (2) A aµ A bν ( p ) = δ ab (cid:2) Γ ⊥ ( p ) P ⊥ µν ( p ) + Γ k ( p ) P k µν ( p ) (cid:3) , (2) Γ (2) c a c b ( p ) = δ ab p J ( p ) , (3)where P ⊥ µν and P k µν are the transverse and parallel pro-jector respectively, defined as: P k µν ( p ) = p µ p ν p and P ⊥ µν ( p ) = δ µν − P k µν ( p ) . In linear covariant gauges and in particular in the lin-ear version of the Landau gauge without the inclusion ofthe mass term, the gluon self-energy is transverse. Here,the longitudinal part is non-zero but it is controlled bya non-renormalization theorem (see [32]). It is impor-tant to stress, however, that this longitudinal part has noimpact on the gluon propagator, which is transverse inLandau gauge even in presence of the mass term. Similarremarks apply to higher vertex functions, that containslongitudinal parts, but which do not contribute to thecorresponding correlation functions. The function J ( p ) (usually called the dressing function of the ghost) andthe Γ ⊥ ( p ) (the transverse part of the two-point gluonvertex, related to the gluon dressing function p Γ ⊥ ( p ) ) ex-plicitly appear in the ghost and gluon propagators andare therefore of special interest in order to compare ourcalculations with lattice results. As the longitudinal partof the two-point gluon vertex is not directly accessible inlattice simulations, we do not describe it here.At one-loop, the gluon and ghost two-point vertexfunctions are given by the diagrams shown in Fig. 1.Most of these diagrams were already computed in [31]. In FIG. 1: First line: four diagrams contributing to the gluonself-energy. Second line: diagram contributing to the ghostself-energy what concerns the gluon propagator, we just need here tocompute the diagram with a quark loop (fourth diagramof Fig. 1), which can be expressed in terms of Passarino-Veltman integrals as Γ ⊥ , ( p ) = 2 g T f ( d − (cid:8) [4 M − ( d − p ] B ( p, M, M )+2( d − A ( M ) } , (4)where the index represents the fourth diagram in Fig. 1and T f is defined by Tr( t a t b ) = T f δ ab (in the fundamen-tal representation, T f = 1 / ). The A and B functionsare the analogue of Passarino-Veltman integrals [36] ineuclidean space: A ( m ) = Z d d q (2 π ) d q + m B ( p, m , m ) = Z d d q (2 π ) d q + m q + p ) + m Each quark flavor contributes to the gluon vertex and it istherefore necessary to sum this diagram over the flavors,with M replaced by the corresponding quark mass. Wehave checked that (4) coincides with the expression of [37]when the quark mass is set to zero.In d = 4 − ǫ , the diagram can be expressed in a com-pletely analytical form: Γ ⊥ , ( p ) = g T f p π n − ǫ + log (cid:18) M e γ/ √ π (cid:19) −
56 + 2 t + (1 − t ) √ t + 1 coth − (cid:0) √ t + 1 (cid:1) o + O ( ǫ ) where t = M p .At one loop, the ghost propagator is the same as inthe quenched situation due to the non-existence of aghost-quark vertex. However, the ghost dressing func-tion will be indirectly influenced by the quarks throughthe renormalization-group flow of the coupling constantand gluon mass. III. QUARK PROPAGATOR
The one-loop contribution to the quark two-point ver-tex involves only one diagram which is shown in Fig. 2.It has two independent structures in Dirac indices and istherefore parametrized by two scalar functions: Γ (2) ψ ¯ ψ ( p ) = Z − ( p ) (cid:0) i/p + M ( p ) (cid:1) Both M ( p ) and Z ( p ) have been determined in latticesimulations. FIG. 2: diagram contributing to the quark two point vertex
This diagram is expressed in terms of Passarino- Veltman integrals as follows: Γ (2) ψ ¯ ψ, ( p ) = g ( d − M B ( p, m, M ) − i/p g C f m p n(cid:2) (2 − d ) m + ( d − m ( M − p )+ ( M + p ) (cid:3) B ( p, m, M )+ ( M + p ) B ( p, , M ) + A ( m )[(2 − d ) m − M − p ] + ( d − m A ( M ) o , where C f is defined by t a t a = C f
1. (In the fundamentalrepresentation, C f = N − N .) The previous expressioncoincides with that of [38] when the gluon mass m is setto zero.In d = 4 − ǫ , the diagram takes the analytical form: Γ (2) ψ ¯ ψ, ( p ) = − i/pg C f π m p n k (cid:2) m + m ( p − M ) − ( M + p ) (cid:3) Q − m p ( − m + M + p ) − m + 3 m ( p − M ) + ( M + p ) ] log (cid:18) Mm (cid:19) − M + p ) log (cid:18) M + p M (cid:19) o + 3 g C f M π n − ǫ + log (cid:18) me γ/ √ π (cid:19) − − k p Q + 12 p ( m − M + p ) log (cid:18) Mm (cid:19) o + O ( ǫ ) with: k = p m + 2 m ( p − M ) + ( M + p ) ,Q = log (cid:20) ( k − p ) − ( M − m ) ( k + p ) − ( M − m ) (cid:21) . (5)It is worth to be mentioned that, at one loop, the fieldrenormalization (the part proportional to /p in the previ-ous expression) is finite. In fact, it even vanishes in thelimit of vanishing gluon mass. IV. RENORMALIZATION ANDRENORMALIZATION GROUPA. Renormalization scheme
In four dimensions most of the expressions presentedabove are divergent. In order to absorb these divergences,we redefine the coupling constant, masses and fields byintroducing renormalization factors: A a µ = p Z A A a µ , c a = p Z c c a , ψ a = p Z ψ ψ a , ¯ c a = p Z c ¯ c a ,g = Z g g m = Z m m M = Z M M The index represents the bare quantities and for nowon, when not specified, all quantities are the renormalizedones.The renormalization constants are redefined in theinfrared-safe (IS) scheme: Γ ⊥ ( p = µ ) = m + µ , J ( p = µ ) = 1 ,Z ( p = µ ) = 1 , M ( p = µ ) = M,Z m Z A Z c = 1 , Z g p Z A Z c = 1 . (6)The IS scheme is convenient because it does not presenta Landau pole [32]. It combines the Taylor scheme [40]and a non-renormalization theorem for the gluon mass,conjectured in [39], and proved in [41–43].The explicit expressions for the renormalization factorsare given in [32] for N f = 0 and for arbitrary N f theirdivergent parts match in Landau gauge with the resultspresented in [39]. B. Renormalization Group
The renormalization procedure leads to finite expres-sions for the vertex functions. However, as is well known,these expressions are hampered by large logarithms andwe have to implement the renormalization-group proce-dure to control perturbation theory. The β functions andanomalous dimensions of the fields are: β g ( g, m , { M i } ) = µ dgdµ (cid:12)(cid:12)(cid:12) g ,m ,M i, ,β m ( g, m , { M i } ) = µ dm dµ (cid:12)(cid:12)(cid:12) g ,m ,M i, ,γ A ( g, m , { M i } ) = µ d log Z A dµ (cid:12)(cid:12)(cid:12) g ,m ,M i, ,γ c ( g, m , { M i } ) = µ d log Z c dµ (cid:12)(cid:12)(cid:12) g ,m ,M i, ,β M i ( g, m , { M i } ) = µ dM i dµ (cid:12)(cid:12)(cid:12) g ,m ,M i, ,γ ψ i ( g, m , { M i } ) = µ d log Z ψ i dµ (cid:12)(cid:12)(cid:12) g ,m ,M i, . For completeness, we give here the contribution of thequarks to the various β and γ functions: γ quarks A = N f X i =1 g T f π n t − t − (cid:16) √ t +4 −√ t √ t + √ t +4 (cid:17)p t ( t + 4) o γ ψ = g C f π m µ n m µ (cid:0) m − M + µ (cid:1) + (cid:2) − (cid:0) m − m M + M (cid:1) − µ (cid:0) M + m (cid:1) + µ (cid:3) × log (cid:18) Mm (cid:19) + (cid:0) − M − M µ + µ (cid:1) log (cid:18) M + µ M (cid:19) + 12 k (cid:2) m − M ) (2 m + M )+ µ ( m − M )(7 m + 6 m M + 5 M )+ µ (cid:0) m − m M − M (cid:1) + µ (cid:0) m + M (cid:1) + µ (cid:3) Q o β M = M γ ψ + 3 g C f M π (cid:26) − u + ( m − M ) log (cid:18) Mm (cid:19) − k (cid:0) m + m (cid:0) u − M (cid:1) + M (cid:0) M + u (cid:1)(cid:1) Q (cid:27) where Q is the expression appearing in Eq. (5) with p replaced by the renormalization-group scale µ . As dis-cussed above, the ghost 2-point function is not affectedby the quarks at one loop so that γ c is not modified andwe recall that, in the IS scheme, β m = m ( γ A + γ c ) and β g = g ( γ A / γ c ) .We can then use the RG equation for the vertex func-tion with n A gluon legs, n c ghost legs and n ψ quark legs: (cid:16) µ∂ µ −
12 ( n A γ A + n c γ c + n ψ γ ψ )+ β g ∂ g + β m ∂ m + X i β M i ∂ M i (cid:17) Γ ( n A ,n c ,n ψ ) = 0 , to relate these functions at different scales, giving: Γ ( n A ,n c ,n ψ ) ( { p i } , µ, g ( µ ) , m ( µ ) , { M i ( µ ) } ) = z A ( µ ) n A / z c ( µ ) n c / z ψ ( µ ) n ψ / × Γ ( n A ,n c ,n ψ ) ( { p i } , µ , g ( µ ) , m ( µ ) , { M i ( µ ) } ) . Here g ( µ ) , m ( µ ) and M i ( µ ) are obtained by integrationof the beta functions with initial conditions given at somescale µ and: log z A ( µ ) = Z µµ dµ ′ µ ′ γ A (cid:0) g ( µ ′ ) , m ( µ ′ ) , M ( µ ′ ) (cid:1) , log z c ( µ ) = Z µµ dµ ′ µ ′ γ c (cid:0) g ( µ ′ ) , m ( µ ′ ) , M ( µ ′ ) (cid:1) , log z ψ ( µ ) = Z µµ dµ ′ µ ′ γ ψ (cid:0) g ( µ ′ ) , m ( µ ′ ) , M ( µ ′ ) (cid:1) . Note that each quark mass M i has its own β M i functionthat must be integrated. In our 1-loop calculation, theflow of M i does not depend on all the quark masses butonly on M i itself. As the infrared safe scheme does notpresent a Landau pole the RG scale µ will be chosen as µ ≃ p for a correlation function with typical momentum p in order to avoid large logarithms. V. RESULTS FOR THE GLUON AND GHOSTSECTORS
In this section we compare the gluon and ghost prop-agator obtained in our one-loop analytical expressionsand in lattice simulations, for SU (3) with N f = 2 and N f = 2 + 1 + 1 .When comparing our findings with the lattice data, wehave to fix the initial conditions of the renormalization-group flow. These quantities were considered here as fit-ting parameters and were chosen to minimize simultane-ously the relative error for the gluon and ghost propaga-tors, respectively defined as: χ AA = 1 N X i Γ ⊥ lt . ( p i ) (cid:18) ⊥ lt . ( p i ) − ⊥ th . ( p i ) (cid:19) χ cc = 1 N X i J − . ( p i ) ( J lt . ( p i ) − J th . ( p i )) There are only few parameters to choose in the fittingprocedure. They correspond to the initial conditions atsome scale µ of the coupling constant and masses in therenormalization-group flow. In the N f = 2 + 1 + 1 casethere are a priori three quark masses to fit. Howeverwe fixed, in the initial condition, the middle and heavyquarks to be twice and 20 times heavier than the lightestone (to obtain the same mass ratios as obtained from thelattice at two GeV). Therefore we have to fit only threeparameters: the coupling constant, the gluon mass andthe light quark mass, all at some renormalization scale(we choose µ = 1 GeV here).[47]For N f = 2 , the best fits were obtained for g = 4 . , m = 0 . GeV and M u,d = 0 . GeV. The correspondinggluon and ghost propagators are depicted in Fig. 3. For N f = 2 + 1 + 1 , the best fits were obtained for g = 5 . , m = 0 . GeV and M u,d = 0 . GeV. The correspondinggluon and ghost propagators are shown in Fig. 4. Thecomparison is very satisfactory, with an error of at mostthree percent. Similar precisions were already obtainedin the quenched situation, see [31–33].It is important to stress that a priori we have threeparameters to fix but, in fact, in QCD, once the scale isfixed, the coupling constant should not be arbitrary. Ata given scale, its value should be compatible to the RGevolution from the known value (for example, at the Z mass scale). We verified that the RG evolution from thecoupling g = 5 . to the Z mass scale gives a couplingat that scale of g = 1 . to be compared to g = 1 . (see, for example [34]). In order to arrive to this scale wehad to use the bottom mass scale 4.2 GeV above wherewe took into account a fifth quark in the running of thecoupling. This gives a 17% error on the coupling that ismore or less the overall 1-loop calculation. This is similarto the corresponding calculation done in the quenchedcase [32]. In this sense, the coupling constant is not afree adjustable parameter. However, in order to studythe infrared, it is more convenient to fit the coupling atscales of the order of 1 GeV. H GeV L p (cid:144) G ¦ H p L H GeV L J H p L FIG. 3: Gluon dressing function (top) and ghost dressingfunction (bottom) as a function of momentum in d = 4 for N f = 2 . The points are lattice data of [44]. The contribution of the new diagram with an internalloop of quark is not strong enough to modify consider-ably the infrared behavior of the propagators due to theIR-safe structure of the quark propagator. Moreover theoptimal values of the coupling constant and the gluonmass do not depend on slight changes in the quark mass.That means that the error for the gluon and ghost prop-agators do not change significantly if we change the valueof the light quark mass slightly. Fig. 5 shows the contourregions for the errors. The contour regions are vertical,showing that the quality of the fit is almost insensitiveto the quark mass (vertical axes). As expected, if thequark mass considerably increases at values of the orderof the GeV, the propagators tend to those obtained inthe quenched approximation.It is interesting to compare the results obtained for thegluon propagator and the ghost dressing function at oneloop for different number of quarks. In Fig. 6 we presentthe gluon propagator and the ghost dressing function for N f = 2 and N f = 2 + 1 + 1 . To compare their infraredbehaviors, we normalized the curves such that they co-incide at GeV. As was observed in lattice simulation[15, 45], the addition of heavy quarks leads to a suppres-sion of the IR saturation point [48]. We also see that theghost dressing function is enhanced. H GeV L p (cid:144) G ¦ H p L H GeV L J H p L FIG. 4: Gluon dressing function (top) and ghost dressingfunction (bottom) as a function of momentum in d = 4 for N f = 2 + 1 + 1 . The points are lattice data of [45]. VI. RESULTS IN THE QUARK SECTOR
We computed the quark propagator for arbitrary num-ber of flavors but as the only available data for the quarksector corresponds to two degenerates light quarks plus aheavy one ( N f = 2 + 1 ), we show the results only for thiscase. As in the previous section, we fixed the ratio of thelight and heavy quarks to the ratio of the bare massesused in lattice simulations. Accordingly, in the presentcase, we fix the mass of the heavy quark to be five timesheavier than the light one at µ = 1 GeV.As discussed in the previous section, the fits of thegluon and ghost propagators are rather insensitive to thechoice of the quark mass (see Fig. 5). However our resultsfor M u,d ( p ) depend strongly on this parameter. There-fore we fixed the fitting parameters by minimizing simul-taneously the error on the gluon propagator and on thefunction M u,d ( p ) , χ M . (Note that the ghost propagatorwas not extracted from lattice simulations for N f = 2 + 1 and could not be used for fixing the parameters, see be-low.) It is defined as χ M = 1 N X i M − . ( µ ) ( M lt . ( p i ) − M th . ( p i )) (7) H GeV L M H G e V L H GeV L M H G e V L H GeV L M H G e V L H GeV L M H G e V L FIG. 5: Contour levels for the quantities χ AA (left) and χ c ¯ c (right) for d = 4 , both for N f = 2 (above) and N f = 2 + 1 + 1 (below). The contourlines correspond to 8%, 10% and 12 %(top left), 1.5%, 2% and 2.2 % (top right), 3%, 6% and 10 %(bottom left), 1%, 2% and 4% (bottom right).The quality ofthe fits is almost insensitive to the value of the quark mass. H GeV L p (cid:144) G ¦ H p L N f = + + N f = H GeV L J H p L N f = + + N f = FIG. 6: Gluon dressing function (top) and ghost dressingfunction (bottom) for different number of flavors in d = 4 Note that we used a slightly different definition of theerror. The reason is that the mass function becomesrapidly close to zero when the momentum grows. Ac-cordingly a relative error point by point would be badlyestimated because the signal to noise ratio of lattice databecomes of order one. In order to avoid this problem wenormalized the data with a typical value of the mass inthe infrared (at the scale µ ).The best fit parameters are g = 4 . , m = 0 . GeV and M u,d = 0 . GeV at µ = 1 GeV. We insist on that thecoupling and gluon mass are fixed from the gluon propa-gator and once this is done, we fixed the mass parameterfrom the scalar part of the quark propagator. The resultsobtained are shown in Fig. 7. The gluon propagator isagain reproduced with an accuracy similar to that ob-tained for N f = 2 and N f = 2 + 1 + 1 .The one-loop results for M u,d ( p ) compare correctlywith lattice simulations. In particular we retrieve thedramatic increase of the mass at a momentum scale ofthe order of 2 GeV. We want to stress that we couldhave a perfect fit of M u,d ( p ) , but in doing so we wouldspoil the quality of the gluon propagator. Instead, wechoose to maintain the quality of the gluon propagatorat the prize of a not-so-good quark mass curve. In theFig. 8 it is shown that there is a region of parameterswhere the fit is excellent (see Fig. 9), but we insist inusing a value of the coupling compatible with a good fitfor the gluon dressing function. The ghost dressing func-tion was not extracted from lattice data and we can nottest our findings in this case. However, we expect thatthis function is rather insensitive to the inclusion of aheavy quarks. Under this hypothesis, it is interesting tocompare our findings with the ghost dressing function for N f = 2 + 1 + 1 . We find fits (not shown) of the samequality as in the previous section. For completeness, wealso compare the mass functions for the light and heavyquarks in Fig. 10.If the one-loop correlation functions compare very wellfor the ghosts, gluons propagators and correctly for thefunction M u,d ( p ) , the results for the quark renormaliza-tion function Z u,d ( p ) is not of the same quality [49], ascan be seen in Fig. 11. We attribute this mismatch tothe fact that the one-loop contribution to this functionis unusually small. This is a consequence of the factthat, when the mass of the gluon is zero, Z u,d ( p ) has nocontribution at one loop but it does have contributionsat two-loops. In this situation, the two-loop contribu-tions are not negligible and we do not expect a one-loopcalculation to be sufficient for describing the character-istics of the field renormalization of the quarks. In thenext section, we estimate the contribution of the two-loop diagrams and show that they have a small impacton our findings for the ghost and gluon propagators andon M u,d ( p ) but that they give a large contribution to Z u,d ( p ) . H GeV L p (cid:144) G ¦ H p L H GeV L M u , d H p LH G e V L FIG. 7: Gluon dressing function (top) and quark mass M u,d ( p ) (bottom) in d = 4 for N f = 2 + 1 . The parame-ters are g = 4 . , m = 0 . GeV and M = 0 . . The pointsare lattice data of [15, 46]. H GeV L g FIG. 8: Contour levels for the quantity χ M . The contourlinescorrespond to 10%, 20% and 30 %. VII. PRELIMINARY ESTIMATE OFTWO-LOOP CONTRIBUTION
The results presented in Fig. 11 show that first orderperturbation theory is not enough to reproduce the be-havior of the Z ( p ) in the Landau gauge. For this reasonwe want to estimate the importance of two loop correc-tions. Ideally, we should compute all the two-loop dia- H GeV L M H p LH G e V L FIG. 9: Quark mass M u,d ( p ) in d = 4 for N f = 2 + 1 . Theparameters are g = 9 . , m = 1 . GeV and M = 0 . GeV.The points are lattice data of [15, 46]. H GeV L M H p LH G e V L M s M u , d FIG. 10: Quark masses of the light (full line, red) and heavy(dashed, black) quark. grams and repeat the analysis performed above. Obvi-ously this work exceeds the purpose of the present article.Instead we use a cartoon expression of the two-loop β functions (let us call it hybrid expression) obtainedby complementing the β functions and anomalous di-mensions derived here at one loop with the ultraviolettwo-loop contribution computed in [39] with a standard(massless) Faddeev-Popov Lagrangian. To take into ac-count the suppression of massive contribution at low mo-menta, the ultraviolet two-loop contributions have to beappropriately regularized in that regime [50]. Here wechoose to modify the anomalous dimension presented in[39] in a minimal way by multiplying them by: µ µ + s This factor goes to one in the UV limit where the wholeexpression matches with [39]. On the other hand, in theinfrared regime it goes to zero as µ , as expected. In thefollowing, we present our results for s = 1 GeV but wechecked that we obtain the same conclusions by changing H GeV L Z u , d H p L FIG. 11: Z u,d ( p ) in d = 4 . The points are lattice data of [46]. the parameter in the range 0.5—2 GeV.We have performed a fit of the lattice data with thishybrid flow equations. The best fits were obtained for g = 3 . , m = 0 . GeV and M = 0 . GeV, again de-fined at µ = 1 GeV. These fits are depicted in Fig. 13,in the center of the shaded areas. We estimate the er-ror coming from neglecting the 2-loop contribution as thedifference between the best fit in the hybrid model and inthe purely one-loop calculation. In the Fig. 13 the hybridcalculation is taken as the central value (that includes arough estimate of the two-loop contributions) and theerror band estimated as explained just before is repre-sented by the shaded area. We clearly see that highercorrections for the gluon propagator are small and com-patible with the discrepancy between our one-loop resultsand the lattice data. The two-loop contribution to thefunction M ( p ) is small but does not seem sufficient tomake the one-loop results compatible with lattice data.This could be related with the spontaneous breaking ofchiral symmetry which was not treated here. Finally, thecorrections to Z ( p ) are large and can explain the discrep-ancy between the lattice data and the 1-loop results. Asmentioned in the previous section, this last function istherefore more difficult to describe and it is necessary togo to next to leading order to find a satisfactory match-ing. VIII. CONCLUSIONS
In this article, we presented a one-loop calculation ofthe two-point correlation functions of QCD in the Lan-dau gauge. The effect of the Gribov copies is minimallyencoded in a mass term for the gluons. We find that thegluon and ghost propagators are reproduced with highprecision. We estimate the higher order corrections tobe small which seems to indicate that the Curci-Ferrarimodel reproduces well the lattice data.On the other hand, one loop calculations are notenough to describe properly all the properties of the H GeV L p (cid:144) G ¦ H p L H GeV L M u , d H p LH G e V L H GeV L Z u , d H p L FIG. 12: Gluon dressing function (top), quark mass M ( p ) (middle) and Z ( p ) (bottom) in d = 4 using the hybridapproximation. The points are lattice data of [15, 46]. quark sector. The one loop contribution to Z u,d ( p ) func-tion is very small and it vanishes if the gluon mass goesto zero in Landau gauge. That is why the two loopscontributions are important and we have to take theminto account if we want to quantitatively reproduce the Z u,d ( p ) function behavior.The mass of the quark is correctly reproduced eventhough the accuracy is not as good as it is for the gluonsector. It is even unclear if a perfect matching for thisquantity by such a simple one-loop calculation is possible H GeV L p (cid:144) G ¦ H p L H GeV L M u , d H p LH G e V L H GeV L Z u , d H p L FIG. 13: Gluon dressing function (top), quark mass M u,d ( p ) (middle) and Z u,d ( p ) (bottom) as a function of momentumin d = 4 . In red the hybrid calculation results (see text) andin black the one loop calculation and the symmetrized oneloop results with respect to hybrid calculation. The pointsare lattice data of [44] and [15, 46]. given that we have not studied the influence of the chiralsymmetry breaking. This analysis remains for a futurework.As we discuss in the article, computing the two loopcontribution can help to estimate the validity of pertur-bation theory in this model. That is why we are consid-ering in doing the two loop calculation completely.0 Acknowledgments
The authors want to thank A. Sternbeck, J. Rodriguez-Quintero and B. El-Bennich for kindly making availablethe lattice data. We also want to thank J. A. Gracey for providing the two loop massless results for each dia-gram independently. The authors want to acknowledgepartial support from PEDECIBA and ECOS programs.M. Peláez wants to thank the ANII for financial support. [1] V. N. Gribov, Nucl. Phys. B (1978) 1.[2] D. Zwanziger, Nucl. Phys. B , 513 (1989).[3] D. Zwanziger, Nucl. Phys. B , 477 (1993).[4] D. Zwanziger, Phys. Rev. D (2002) 094039.[5] D. Dudal, J. A. Gracey, S. P. Sorella, N. Vandersickeland H. Verschelde, Phys. Rev. D (2008) 065047.[6] P. van Baal, Nucl. Phys. B , 259-275 (1992).[7] I. Montvay and G. Munster, “Quantum fields on a lat-tice,” Cambridge, UK: Univ. Pr. (1994).[8] A. Sternbeck, L. von Smekal, D. B. Leinweberand A. G. Williams, PoS LAT , 340 (2007)[arXiv:0710.1982 [hep-lat]].[9] A. Cucchieri and T. Mendes, Phys. Rev. Lett. (2008)241601 and also arXiv:1001.2584 [hep-lat].[10] A. Cucchieri and T. Mendes, Phys. Rev. D (2008)094503.[11] A. Sternbeck and L. von Smekal, Eur. Phys. J. C , 487(2010)[12] A. Cucchieri and T. Mendes, Phys. Rev. D (2010)016005.[13] I. L. Bogolubsky et al. , Phys. Lett. B , 69 (2009).[14] D. Dudal, O. Oliveira and N. Vandersickel, Phys. Rev. D (2010) 074505.[15] P. O. Bowman, U. Heller, D. Leinweber, M. Parappilly,A. Williams, and others, Phys. Rev. D , 034509(2004).[16] M. B. Parappilly, P. O. Bowman, U. M. Heller,D. B. Leinweber, A. G. Williams and J. B. Zhang, AIPConf. Proc. , 237 (2006) [hep-lat/0601010].[17] L. von Smekal, R. Alkofer and A. Hauck, Phys. Rev. Lett. (1997) 3591.[18] R. Alkofer and L. von Smekal, Phys.Rept. (2001) 281.[19] C. S. Fischer and R. Alkofer, Phys. Rev. D (2003)094020.[20] J. C. R. Bloch, Few Body Syst. (2003) 111.[21] C. S. Fischer, A. Maas and J. M. Pawlowski, AnnalsPhys. (2009) 2408.[22] A. C. Aguilar and A. A. Natale, JHEP (2004) 057.[23] Ph. Boucaud et al. , JHEP (2006) 001.[24] A. C. Aguilar and J. Papavassiliou, Eur. Phys. J. A (2008) 189.[25] A. C. Aguilar, D. Binosi and J. Papavassiliou, Phys. Rev.D (2008) 025010.[26] P. Boucaud, J. P. Leroy, A. Le Yaouanc, J. Micheli,O. Pene and J. Rodriguez-Quintero, JHEP (2008) 099.[27] J. Rodriguez-Quintero, JHEP (2011) 105.[28] M. Q. Huber and L. von Smekal, JHEP , 149 (2013)arXiv:1211.6092 [hep-th].[29] G. Curci and R. Ferrari, Nuovo Cim. A , 151 (1976).[30] J. Serreau and M. Tissier, Phys. Lett. B (2012) 97 [arXiv:1202.3432 [hep-th]].[31] M. Tissier and N. Wschebor, Phys. Rev. D (2010)101701 [arXiv:1004.1607 [hep-ph]].[32] M. Tissier and N. Wschebor, Phys. Rev. D (2011)045018 [arXiv:1105.2475 [hep-th]].[33] M. Pelaez, M. Tissier and N. Wschebor, Phys. Rev. D , 125003 (2013).[34] B. Blossier, P. Boucaud, M. Brinet, F. De Soto, X. Du,M. Gravina, V. Morenas and O. Pene et al. , Phys. Rev.D , 034503 (2012)[35] U. Reinosa, J. Serreau, M. Tissier, and N. Wschebor,Phys. Rev. D , 105016 (2014).[36] G. Passarino and M. J. G. Veltman, Model,” Nucl. Phys.B (1979) 151.[37] A. I. Davydychev, P. Osland and O. V. Tarasov, Phys.Rev. D , 101 (2003)[hep-th/0211144].[40] J. C. Taylor, Nucl. Phys. B (1971) 436.[41] D. Dudal, H. Verschelde and S. P. Sorella, Phys. Lett. B (2003) 126.[42] N. Wschebor, Int. J. Mod. Phys. A (2008) 2961.[43] M. Tissier and N. Wschebor, Phys. Rev. D , 065008(2009).[44] A. Sternbeck, K. Maltman, and M. Muller-Preussker andL. von Smekal , PoS LATTICE2012, 243 (2012).[45] A. Ayala, A. Bashir, D. Binosi, M. Cristoforetti and J.Rodriguez-Quintero, Phys. Rev. D , 74512 (2012).[46] P. O. Bowman, U. Heller, D. Leinweber, M. Parappilly,A. Williams, and others, Phys. Rev. D , 54507 (2005).[47] The propagators are, generally speaking, not normalizedas in (6). In consequence, by using multiplicative renor-malizability, we multiplied each lattice propagator by aconstant in order to impose the renormalization condi-tion (6). In the case of unrenormalized propagators, thisalso implements the required field renormalization.[48] This effect is a consequence of both varying N f and changing the renormalization-group trajectory. If wechange only N f but use the same initial conditions (at,say, µ = 4 GeV to obtain similar correlations in theUV), the quarks tend to suppress the gluon propagatoraround 1 GeV but enhance it at very low momenta.[49] For this reason, we did not use the function Z ( p ) forchoosing the best fitting parameters, but concentratedon the function M u,d ( p ))