On the contribution of different coupling constants in the infrared regime of Yang-Mills theory: a Curci-Ferrari approach

TL;DR

This paper explores how different interaction vertices affect two-point correlation functions in the infrared regime of Yang-Mills theory, using a phenomenological massive model and one-loop calculations, aligning well with lattice data.

Contribution

It introduces a phenomenological approach to analyze the role of various vertices in the infrared regime of Yang-Mills theory using a Curci-Ferrari model.

Findings

Good agreement with lattice data.

Different vertices have distinct impacts in the infrared.

One-loop calculations effectively capture infrared behavior.

Abstract

We investigate the influence of the different vertices of two-point correlation functions in the infrared regime of Yang-Mills theory using a phenomenological description. This regime is studied in Landau-gauge and using perturbation theory within a phenomenological massive model. We perform a one-loop calculation for two-point correlation functions taking into account the different role of the various interactions in the infrared. Our results show a good agreement with the lattice data.

Equations

$$\mathcal{L}=\frac{1}{4}{F}_{\mu \nu }^a{F}_{\mu \nu }^a+{\partial }_{\mu }{\bar{c}}^a(D_{\mu }c{)}^a+\frac{m^2}{2}{A}_{\mu }^a{A}_{\mu }^a+\text{i}{b}^a{\partial }_{\mu }{A}_{\mu }^a,$$

$$F_{\mu \nu }^a={\partial }_{\mu }{A}_{\nu }^a-{\partial }_{\nu }{A}_{\mu }^a+g_A{f}^{abc}{A}_{\mu }^b{A}_{\nu }^c,$$

$$(D_{\mu }c{)}^a={\partial }_{\mu }{c}^a+g_C{f}^{abc}{A}_{\mu }^b{c}^c.$$

$${\Gamma }_{\mu \nu }^{ab}(p)={\delta }^{ab}\left({\Gamma }^{\perp }(p)P_{\mu \nu }^{\perp }(p)+{\Gamma }^{\parallel }(p)P_{\mu \nu }^{\parallel }(p)\right),$$

$$P_{\mu \nu }^{\parallel }(p)=\frac{p_{\mu }{p}_{\nu }}{p^2},$$

$$P_{\mu \nu }^{\perp }(p)={\delta }_{\mu \nu }-P_{\mu \nu }^{\parallel }(p).$$

$$\begin{array}{rl}{D}_1=& \frac{N(g_C^B{)}^2{\delta }^{ab}}{(4\pi {)}^{d/2}}\frac{\Gamma (2-d/2)}{d-2}{p}^{d-4}\\ & \times ((d-2)p_{\mu }{p}_{\nu }+p^2{\delta }_{\mu \nu })\frac{{\Gamma }^2(d/2)}{\Gamma (d)},\end{array}$$

$$\begin{array}{r}{D}_2=\frac{N(g_A^B{)}^2{\delta }^{ab}}{(4\pi {)}^{d/2}}\frac{2(d-1{)}^2}{d(d-2)}{m}^{d-2}\Gamma (2-d/2){\delta }^{\mu \nu }.\end{array}$$

$$\begin{array}{rl}& A(m)=\int \frac{d^{d}q}{(2\pi {)}^d}\frac{1}{q^2+m^2},\\ & B(m_1,m_2)=\int \frac{d^{d}q}{(2\pi {)}^d}\frac{1}{q^2+m_1^2}\frac{1}{(q+p{)}^2+m_2^2},\end{array}$$

$$\begin{array}{rl}& D_3=(g_A^B{)}^{2}N{\delta }^{ab}[P_{\mu \nu }^{\parallel }(p)\left\{A(m)\left(d+\frac{1}{d}-\frac{m^2}{4p^2}-\frac{7}{4}\right)-\frac{{\left(m^2+p^2\right)}^{2}B(m,0)}{4p^2}\right\}\\ & +P_{\mu \nu }^{\perp }(p)\{-\frac{p^{6}B(0,0)}{8(d-1)m^4}-\frac{\left(4m^2+p^2\right)B(m,m)\left(4(d-1)m^4+4(3-2d)m^2{p}^2+p^4\right)}{8(d-1)m^4}\\ & +\frac{{\left(m^2+p^2\right)}^{2}B(m,0)\left(2(3-2d)m^2{p}^2+m^4+p^4\right)}{4(d-1)m^4{p}^2}+\frac{A(m)\left(-4d^2{p}^4+d\left(m^4+5m^2{p}^2+7p^4\right)-4m^2{p}^2\right)}{4(d-1)d{m}^2{p}^2}\}].\end{array}$$

$$\begin{split}\Gamma_{c^{a}\bar{c}^{b}}^{(2)}(p)=&\delta^{ab}p^{2}-\frac{\delta^{ab}N(g_{C}^{B})^{2}}{4m^{2}}\Big{(}(m^{2}-p^{2})A(m)\\ &-p^{4}B(0,0)+(m^{2}+p^{2})^{2}B(0,m)\Big{)}.\end{split}$$

$$\begin{array}{rl}& A_B^{a\mu }=\sqrt{Z_A}{A}^{a\mu },c_B^a=\sqrt{Z_c}{c}^a,\\ g_A^B=& Z_{g_A}{g}_A,g_C^B=Z_{g_C}{g}_C,m_B^2=Z_{m^2}{m}^2.\end{array}$$

$$\begin{array}{rl}& {\Gamma }_{AA}^{\perp }(\mu )=m^2+{\mu }^2,\\ & {\Gamma }_{c\bar{c}}(\mu )={\mu }^2,\\ & {\Gamma }_{AA}^{\parallel }(\mu )=m^2.\end{array}$$

$$\begin{array}{r}{\Gamma }_{A_{\alpha }^a{C}^b{\bar{C}}^c}^{(3)}(\mu ,-\mu ,0)=-i{f}^{abc}{g}_C{\mu }_{\alpha }.\end{array}$$

$$Z_C\sqrt{Z_A}{Z}_{g_C}=1.$$

$$\begin{array}{rl}{\Gamma }_{A_{\alpha }^a{A}_{\beta }^b{A}_{\gamma }^c}^{(3)}(\mu ,-\mu ,0)& =-i{f}^{abc}{g}_A{\Gamma }_{A_{\alpha }{A}_{\beta }{A}_{\gamma }}^{\text{tl}}(\mu ,-\mu ,0)\end{array}$$

$${\Gamma }_{A_{\alpha }{A}_{\beta }{A}_{\gamma }}^{\text{tl}}(\mu ,-\mu ,0)=({\mu }_{\alpha }{\delta }_{\beta \gamma }-2{\mu }_{\beta }{\delta }_{\alpha \gamma }+{\mu }_{\gamma }{\delta }_{\alpha \beta }).$$

$$\begin{array}{rl}{\Gamma }_{A_{\alpha }^a{A}_{\beta }^b{A}_{\gamma }^c}^{(3)B}& (p,-p,0)=f_1(p^2)(p_{\gamma }{\delta }_{\alpha \beta }+p_{\alpha }{\delta }_{\beta \gamma })\\ & +f_2(p^2)p_{\beta }{\delta }_{\gamma \alpha }+f_3(p^2)p_{\alpha }{p}_{\beta }{p}_{\gamma }.\end{array}$$

$$\begin{array}{rl}{g}_A& =\frac{-Z_A^{\frac{3}{2}}}{2}{f}_2({\mu }^2,g_A^B,g_C^B).\end{array}$$

$$\begin{array}{rl}{\Gamma }_{\text{div}}^{\perp }(p)& =\frac{N}{96{\pi }^2ϵ}\left((g_C^B{)}^2{p}^2+(g_A^B{)}^2(25p^2-9m^2)\right),\\ {\Gamma }_{\text{div}}^{\parallel }(p)& =\frac{N}{32{\pi }^2ϵ}\left(\left((g_C^B{)}^2-(g_A^B{)}^2\right)p^2-3(g_A^B{)}^2{m}^2\right).\end{array}$$

$$\frac{N{p}^2\left((g_C^B{)}^2-(g_A^B{)}^2\right)}{348{\pi }^2}\left[\frac{12}{ϵ}+6(\gamma +2)+6\text{Log}\left(\frac{p^2}{4\pi }\right)\right]$$

$$\begin{array}{rl}{Z}_C& =1+\frac{3g_C^{2}N}{32ϵ{\pi }^2}+\mathcal{O}(g^4),\\ Z_{m^2}& =1-\frac{34g_A^{2}N+g_C^{2}N}{96ϵ{\pi }^2}+\mathcal{O}(g^4),\\ Z_A& =1+\frac{25g_A^{2}N+g_C^{2}N}{96{\pi }^2ϵ}+\mathcal{O}(g^4).\end{array}$$

$$\begin{gathered}\gamma_{A}=\mu\dfrac{\partial\log(Z_{A})}{\partial\mu},\\ \gamma_{C}=\mu\dfrac{\partial\log(Z_{C})}{\partial\mu},\\ \beta_{m^{2}}=\mu\dfrac{\partial m^{2}}{\partial\mu}\Big{|}_{g_{C}^{B},g_{A}^{B},m^{2}_{B}}=-m^{2}\gamma_{m^{2}},\\ \beta_{g_{A}}=\mu\dfrac{\partial g_{A}}{\partial\mu}\Big{|}_{g_{C}^{B},g_{A}^{B},m^{2}_{B}}=\left(\dfrac{3}{2}\gamma_{A}+\mu\dfrac{\mathrm{d}\log(B(\mu^{2}))}{\mathrm{d}\mu}\right)g_{A},\\ \beta_{g_{C}}=\mu\dfrac{\partial g_{C}}{\partial\mu}\Big{|}_{g_{C}^{B},g_{A}^{B},m^{2}_{B}}=\left(\gamma_{C}+\dfrac{1}{2}\gamma_{A}\right)g_{C}.\end{gathered}$$

$${\beta }_{g_A}^{UV}(g_C,g_A)$$

$${\beta }_{g_C}^{UV}(g_C,g_A)$$

$$\begin{split}\Big{(}\mu\dfrac{\partial}{\partial\mu}+&\dfrac{1}{2}(n_{A}\gamma_{A}+n_{C}\gamma_{C})+\beta_{g_{A}}\dfrac{\partial}{\partial g_{A}}\\ &+\beta_{g_{C}}\dfrac{\partial}{\partial g_{C}}+\beta_{{m^{2}}}\dfrac{\partial}{\partial{m^{2}}}\Big{)}\Gamma^{(n_{A},n_{C})}_{R}=0,\end{split}$$

$$\begin{array}{rl}{\Gamma }_R^{(n)}(\{p\},\mu ,g_A,& g_C)=z_A(\mu {)}^{\frac{n_A}{2}}{z}_C(\mu {)}^{\frac{n_C}{2}}\\ & \times {\Gamma }_R^{(n)}\left(\{p\},{\mu }_0,g_A,g_C,\right).\end{array}$$

$$\begin{array}{r}\log \left(z_A(p)\right)={\int }_{\mu }^p\frac{\mathrm{d}{\mu }^{\prime }}{{\mu }^{\prime }}{\gamma }_A\left(g_A({\mu }^{\prime }),g_C({\mu }^{\prime }),m^2({\mu }^{\prime })\right),\\ \log \left(z_C(p)\right)={\int }_{\mu }^p\frac{\mathrm{d}{\mu }^{\prime }}{{\mu }^{\prime }}{\gamma }_C\left(g_A({\mu }^{\prime }),g_C({\mu }^{\prime }),m^2({\mu }^{\prime })\right).\end{array}$$

$$\begin{array}{r}{\Gamma }_{AA}^{(2)}(p,{\mu }_0,g_A,g_C,m^2)=\frac{p^2+m^2(p^2)}{z_A(p)}=\frac{p^2}{D(p)},\\ {\Gamma }_{c\bar{c}}^{(2)}(p,{\mu }_0,g_A,g_C,m^2)=\frac{p^2}{z_C(p)}=\frac{p^2}{J(p)}\end{array}$$

$$g_A=g_C+\frac{37N}{384{\pi }^2}{g}_C^3+\mathcal{O}\left(g_C^5\right).$$

$${\chi }_{\mathrm{tot}}^2=\frac{{\chi }_{\mathrm{abs}}^2+{\chi }_{\mathrm{rel}}^2}{2}.$$

$${\chi }_{\mathrm{rel}}^2$$

Full text

Quantum chromodynamics, infrared correlation functions

I Introduction

Several efforts have been made in order to understand the infrared (IR) behavior of correlation functions in Yang-Mills theories. The IR is usually known as the non-perturbative regime since standard perturbation theory based on Faddeev-Popov Lagrangian is no longer valid. For this reason other semi-analytical techniques have been implemented (see e.g. Pawlowski et al. (2004); Alkofer and von Smekal (2001); Zwanziger and Daniel (2002); Fischer and Alkofer (2003); Bloch (2003); Fischer and Gies (2004); Natale (2007); Fischer et al. (2009); Rodriguez-Quintero (2011); Huber and von Smekal (2013); Aguilar et al. (2008); Athenodorou, Andreas et al. (2018); Aguilar and Papavassiliou (2008); Alkofer and von Smekal (2001); Aguilar and Natale (2004); Wetterich (1993); Berges et al. (2002)). Correlation functions are important since they can be related to scattering amplitudes and used as a testing ground for different approaches in the study of IR Yang-Mills or QCD. However, these quantities are gauge dependent and it is therefore necessary to compare analytic results with gauge fixed simulations.

At present, one of the main approaches to nonperturbative problems is the lattice simulation. In the last decades Landau gauge-fixed lattice calculations have improved considerably and a consensus has been reached on several unexpected features. First, they found that gluon two-point correlation function behaves as a massive propagator in the IR Cucchieri and Mendes (2009, 2010, 2008); Bogolubsky et al. (2009); Dudal et al. (2010). Second, that the coupling constants remain moderate even at low momenta Bogolubsky et al. (2009); Blossier et al. (2012); Boucaud et al. (2014, 2017, 2018). Both facts motivate the use of a phenomenologically modified gauge-fixed-Lagrangian with the addition of a gluon mass term Tissier and Wschebor (2010, 2011). This model is a phenomenological modification of Faddeev-Popov Lagrangian in Landau gauge which matches with Curci-Ferrari Lagrangian in the same gauge Curci and Ferrari (1976).

The advantage of this phenomenological model is that it is renormalizable de Boer et al. (1996) and IR safe, in the sense that it is possible to find a renormalization scheme with no Landau pole Reinosa et al. (2017). Within this model, one-loop calculations of the two and three point vertex functions have shown good agreement with lattice simulations Tissier and Wschebor (2010, 2011); Peláez et al. (2013). More recently a two-loop calculation have been performed for two-point correlation functions in quenched Landau-gauge Curci-Ferrari model Gracey et al. (2019), and again an overall satisfactory description of the lattice data is obtained. Moreover, this model was also extended to finite temperature and chemical potential Reinosa et al. (2016); Reinosa et al. (2015a, b); Reinosa et al. (2014). In this work we present an improvement to the one-loop calculation of two-point correlation functions within the massive model presented in Tissier and Wschebor (2011). In previous works, as is common in the field, only one running coupling constant was considered. The logic behind this choice is that this is enough to renormalize the theory (in order to preserve gauge symmetry the interactions in different channels must be related to the same bare coupling constant). However, lattice simulations show that the strength of the interaction in different channels can be very different for momenta below 1 GeV (see e.g. Bogolubsky et al. (2009); Blossier et al. (2012); Boucaud et al. (2014); Athenodorou et al. (2016); Boucaud et al. (2017, 2018); Zafeiropoulos et al. (2019)). In order to take into account this effect, we propose to consider different renormalization factors for couplings extracted from different interactions.

The one-loop effect on the infrared behaviour of the different coupling constants in two-point correlation functions can also be obtained with a two-loop calculation in the equal coupling case. In this way, we can interpret our results as if some two-loops effects were taken into account. We stress that, even though these are two-loop corrections in the equal coupling case, the present computation has still the simplicity of a one-loop calculation. We concentrate on the Landau gauge which has been largely studied in the past.

The outline of the article is the following. In Section II we describe the model and present our one-loop calculation. We then describe in Section III the IR safe renormalization scheme and introduce the renormalization group in Section IV. In Section V we show the comparison between our results and lattice data, and discuss the relevance of treating differently the interactions.

II Model

In order to study the IR regime of Yang-Mills theory we consider a model based on phenomenological observations. We take as a starting point the Curci-Ferrari model in the Landau gauge Curci and Ferrari (1976). It consists in adding a gluon mass term in Faddeev-Popov Lagrangian. The addition of this term is motivated by lattice simulations that show a gluon propagator that saturates in the IR behaving as a massive propagator. The advantage of this model is that it admits renormalization schemes without Landau pole and therefore perturbation theory is potentially under control even in the IR. It is important to mention that the inclusion of a mass scale does not modify ultraviolet results since for momentum $p\gg m$ the mass scale can be neglected and, in this limit, the Faddeev-Popov Lagrangian is recovered.

In order to trace out the influence of different couplings in the Lagrangian, we name the coupling constants for the three-gluons vertex, $g_{A}$, and for the ghost-gluon vertex, $g_{C}$. For simplicity, we identify the coupling extracted from the four-gluon vertex with the one associated with the three-gluon vertex.

In the Euclidean space the Landau-gauge Curci-Ferrari Lagrangian density reads:

[TABLE]

where

[TABLE]

The different labels in the couplings are in order to trace them out when renormalizing even though they are the same in the bare Lagrangian.

Accordingly, we depict the different vertices as follow:

We will show that, even though coupling constants behave very differently in the infrared, this difference barely modifies the results of two-point functions.

II.1 Propagators

In this section we compute one-loop correction for the gluon and ghost two-point vertex tracing out the influence of the different interactions.

The gluon two-point vertex function can be described by two scalar functions as

[TABLE]

where the projectors are defined as

[TABLE]

The diagrams contributing to $\Gamma_{\mu\nu}^{ab}(p)$ at one-loop are

The first two diagrams can be computed analytically for arbitrary dimension in terms of elementary functions Tissier and Wschebor (2011).

[TABLE]

and

[TABLE]

For the last one, it is convenient to describe it in terms of one-loop master integrals

[TABLE]

that can be computed analytically for integer dimension in terms of elementary functions,

[TABLE]

The diagram contributing to the ghost two-point vertex function, $\Gamma_{c\bar{c}}(p)$, is

[TABLE]

III Renormalization scheme

In order to absorb the divergent part of one-loop diagrams, we introduce renormalization factors:

[TABLE]

The index $B$ represents the bare quantities and for now on, when not specified, all quantities are the renormalized ones. In contrast with previous one-loop calculation for the equal coupling case, in this work we renormalize differently the ghost-gluon and the three-gluon couplings, even though they are related to the same bare value, $g^{B}_{A}=g^{B}_{C}$.

The field renormalization factors are defined in the scheme presented in Weber and Dall’Olio (2014), which is a reformulation of the IS-scheme of Tissier and Wschebor (2011). In this renormalization scheme, renormalization factors for the fields and the gluon mass respect the conditions:

[TABLE]

The ghost-gluon coupling constant, $g_{C}$, is defined through the ghost-gluon vertex when the ghost momentum vanishes,

[TABLE]

As this vertex is finite, we have a non-renormalization theorem Taylor (1971) that allows to define the renormalization factor $Z{g_{C}}$ in terms of the fields renormalization factor as

[TABLE]

On the other hand, the three-gluon coupling is defined through the three-gluon vertex, $\Gamma_{A_{\alpha}^{a}\hskip 0.85358ptA_{\beta}^{b}A_{\gamma}^{c}}^{(3)}(p,r,k)=-if^{abc}g_{A}\Gamma_{A_{\alpha}\hskip 0.85358ptA_{\beta}A_{\gamma}}(p,r,k)$. In particular, we define $g_{A}$ at the scale $\mu$ so the vertex has a similar form as tree-level when one momentum is zero. The renormalization scheme can be summarized as:

[TABLE]

where $\Gamma_{A_{\alpha}\hskip 0.85358ptA_{\beta}A_{\gamma}}^{\text{tl}}(\mu,-\mu,0)$ represents the tensor structure of tree-level vertex with one vanishing momentum. In particular,

[TABLE]

One loop diagrams contributing to $\Gamma_{A_{\alpha}^{a}\hskip 0.85358ptA_{\beta}^{b}A_{\gamma}^{c}}^{(3)}$ are depicted in Fig.2 and where computed in Peláez et al. (2013).

For one momentum vanishing the vertex has only three scalar functions, $f_{1}(p^{2})$, $f_{2}(p^{2})$ and $f_{3}(p^{2})$, defined as:

[TABLE]

Using the relation $\Gamma_{A_{\alpha}^{a}\hskip 0.85358ptA_{\beta}^{b}\hskip 0.85358ptA_{\gamma}^{c}}^{(3)\,R}=Z_{A}^{\frac{3}{2}}\Gamma_{A_{\alpha}^{a}\hskip 0.85358ptA_{\beta}^{b}\hskip 0.85358ptA_{\gamma}^{c}}^{(3)\,B}$ we obtain

[TABLE]

III.1 The renormalization process with multiple couplings

We observe that once the couplings related to different interactions in the Lagrangian, namely $g_{A}^{B}$ and $g_{C}^{B}$, are treated differently, the two-point vertex functions seems to have quadratic divergences. This is coherent with the fact that we need gauge symmetry, and therefore it is necessary to respect $g_{A}^{B}=g_{C}^{B}$ in order to avoid this kind of divergence. This problem does not happen for the ghost vertex function at one-loop since it has only one kind of interaction.

However, in the gluon vertex function we find that the divergent parts are:

[TABLE]

Note that divergences in the parallel part proportional to $(g_{C}^{B})^{2}-(g_{A}^{B})^{2}$ spoil renormalization. In particular, a term of the form:

[TABLE]

where $\gamma$ is the Euler-gamma, appears in the bare parallel part of gluon vertex function.

Obviously this kind of terms are zero since at bare level the coupling constants are the same. However, if we want the procedure to respect renormalizability even in perturbation theory, this would imply that the renormalized coupling constants at one-loop are the same (at one-loop all bare couplings can be replace by the renormalized ones). Thus, since we are mainly interested in describing phenomenologically the different behaviour of the coupling constants in the infrared with a simple calculation, we propose an ad-hoc procedure for this. Namely, we impose the equality of bare couplings in the ultraviolet divergences, but we still take into account the difference of the renormalized couplings in finite parts.

For practical purposes, in order to ensure the renormalization, we set term (6) to zero but we still keep the various couplings different in the rest of the expression. This ad-hoc procedure allows to obtain a first study of the different behaviours of the coupling constants with a simple one-loop calculation.

In light of this, the divergent part of the renormalization factors in four dimensions, using $d=4-\epsilon$, are given by

[TABLE]

IV Renormalization group

Once we find the renormalization factors we can compute the $\beta$ and $\gamma$ functions, defined as:

[TABLE]

Once we solved this system of differential equations we can find the flow for the different coupling constants and masses.

The $\beta$ functions in the UV for $g_{A}$ and $g_{C}$, keeping track of the couplings, are

[TABLE]

These expressions match with standard one-loop $\beta$-function, $\beta_{g}=-\frac{11Ng^{3}}{3(16\pi^{2})}$, when $g_{A}=g_{C}$. In order to obtain the propagators renormalized at a different scale we use the renormalization group equation:

[TABLE]

which is a generalization of the standard renormalization group equation for this model.

The vertex function with different energies scales are related as

[TABLE]

In equation (10) we took the energy $\mu$ as the reference point in order to resolve the RG equation, with $z_{A}$ and $z_{C}$ defined as

[TABLE]

Therefore two-point functions are computed as:

[TABLE]

where $D(p)$ and $J(p)$ denote the gluon and ghost dressing functions respectively.

V Results

In order to obtain the behavior of the coupling constants and consequently the propagators, we have to give the initial condition of the renormalization group flow. In this case we choose the initial value at $\mu_{0}=10$ GeV. A priori we can think that we have three free parameters to fit, $g_{A}(\mu_{0})$, $g_{C}(\mu_{0})$ and $m(\mu_{0})$. However, as both coupling constants have the same bare value, we can relate them at large energy ($\mu_{0}\approx 10$ GeV) using standard perturbation theory. We specify the calculation of this relation in appendix A. Perturbativelly $g_{A}$ and $g_{C}$ can be related as,

[TABLE]

For large energy $g_{A}$ is barely larger that $g_{C}$. Using this relation, at the end, we have to fix only $m_{0}=m(10\,\textnormal{GeV})$ and $g_{C,0}=g_{C}(10\,\textnormal{GeV})$ in order to fit both propagators at the same time. The parameters are chosen in order to minimize the error function $\chi_{\mathrm{tot}}$ defined as

[TABLE]

Relative and absolute error are defined as follow:

[TABLE]

where $M$ is the number of lattice data and A is the multiplicative factor that appears due to the difference between the IS-renormalization scheme and lattice renormalization. It takes the form

[TABLE]

in order to minimize the quadratic error. The parameters that minimize $\chi$ are $g_{C,0}=1.67$, $m_{0}=0.24$ GeV. In Fig. 3 we show the comparison of the gluon propagator of this work, computed with $g_{C,0}=1.67$, $m_{0}=0.24$ GeV, to the lattice data in $SU(3)$ from Bogolubsky et al. (2009).

From Fig. 3 we infer that a good fitting is obtained from this model. The fitting is only a bit better than in the equal coupling case Tissier and Wschebor (2011). The fact that there is no great difference with the equal coupling case is consistent with the fact that the perturbative expansion is under control.

Apart from the study of two-point functions we are able to study in detail the infrared behaviour of the coupling constants. With this purpose, we study the low momenta limit $\mu\ll m$ of Eqs. (8). Therefore, the one-loop $\beta$ functions in the infrared regime read,

[TABLE]

The analytic solution of the system (15) shows that both coupling constants go to zero in the very deep IR, see appendix B. While the ghost coupling goes to zero by positives values, the gluon coupling arrives from negatives values. The zero crossing is also observed by lattice simulations Duarte et al. (2016); Athenodorou et al. (2016); Boucaud et al. (2017). This is due to the fact that the three-gluon vertex in the infrared is dominated by massless ghost while gluons are frozen Peláez et al. (2013); Aguilar et al. (2019).

Still we can wonder why even though coupling constants are extremely different in the infrared this has almost no effects in the results for two-point functions. We think that it is because they start differing below 1 GeV as it is shown in Fig. 4. Moreover, at one-loop, the expansion parameter related to the ghost-gluon coupling is in fact $\alpha_{C}^{\tiny\text{1-loop}}=\frac{Ng_{C}^{2}}{16\pi^{2}}$ while the expansion parameter of the three-gluon coupling is $\alpha_{A}=\frac{Ng_{A}^{2}}{16\pi^{2}}\frac{\mu^{2}}{\mu^{2}+4m^{2}}$. In the last case, the extra factor appears due to massive gluons and to the fact that, since we are coupling 1PI functions, the internal gluons come at least in pairs. The three-gluon expansion parameter is suppressed in the infrared and therefore is almost equivalent to consider $g_{A}$ or $g_{C}$ in that quantity. Let us note that for higher loops the expansion parameters is for both couplings $\alpha_{A,C}=\frac{Ng_{A,C}^{2}}{16\pi^{2}}\frac{\mu^{2}}{\mu^{2}+4m^{2}}$, since higher loop diagrams have internal gluons. The latter is less than $0.3$ as it can be appreciated in Fig. 5.

VI Conclusion

In this work we present a one-loop calculation of two-point correlation functions using a modified Lagrangian in Landau gauge. Based on the fact that coupling constants related to different vertices behave differently in the infrared, we do not identify vertices with only one coupling as presented in Tissier and Wschebor (2011).

As it was mentioned, this one-loop difference on the couplings can be also taken into account in a two-loop calculation using the model with only one coupling constant. Therefore, the fact that this results are of the same quality as the equal-coupling case show that perturbation theory seems to be under control. A complete two loop calculation of these quantities has been recently reported Gracey et al. (2019).

We can appreciate in Fig. 4 the change of sign in the gluon coupling constant, $g_{A}$, associated to the three gluon vertex. The change of sign is due to the contribution of the diagram with one-loop of ghosts (left of Fig. 2). We find that, even though both coupling constants differ considerably in the infrared, the gluon and ghost propagators slightly improve with respect to the one coupling case of Tissier and Wschebor (2011). This can be understood since the expansion parameters differ below the gluon mass scale, and in that regime gluon fluctuations are almost frozen.

Even though we do not know how to justify this massive Lagrangian from first principles, we think that the simplicity of this model makes it useful to study phenomenologically different IR properties of Yang-Mills theory.

Acknowledgements.

We are very thankful to Nicolás Wschebor, Matthieu Tissier and Emerson Luna for useful discussions and reading previous version of this manuscript. We also thank the support of ECOS program, ANII “$FCE\_126412$” research project and the support of a CNRS-PICS project “irQCD”. M. Fernández wants to thank the Comisión Académica de Posgrado de la Universidad de la República for financial support.

Appendix A Relation between $g_{A}$ and $g_{C}$ in the ultraviolet

For $\mu\gg m$, the calculation is identical to the one in standard Yang-Mills. We can use perturbation theory, in this way we can relate the renormalize couplings with the bare ones as follows:

[TABLE]

were $\Lambda$ is the UV scale. At one-loop order we can write:

[TABLE]

From this we conclude:

[TABLE]

Because $g_{A}$ and $g_{C}$ are renormalized couplings both have a limit for $\Lambda\rightarrow\infty$. This implies that $f_{A}\left(\frac{\mu}{\Lambda}\right)-f_{C}\left(\frac{\mu}{\Lambda}\right)$ has a limit for $\Lambda\rightarrow\infty$. For dimensional reasons, this subtraction must be a pure number. The corresponding result for $f_{A}$ and $f_{C}$ is extracted from the definition of $g_{A}$ y $g_{C}$ using the one-loop vertices computed in Peláez et al. (2013).

Appendix B Coupling constants in the deep IR

The IR regime is characterized by energies that satisfy $\mu\ll m$, being $\mu$ the energy and $m$ the gluon mass. The functions $\beta$ in this regime are,

[TABLE]

We rewrite these equations, at order $g^{3}$, using the changes of variables: $t=\mathrm{log}\left(\frac{\mu}{\bar{\mu}}\right)$ where $\bar{\mu}$ is an arbitrary constant and $y=\frac{g_{C}^{3}}{g_{A}}$.

[TABLE]

Using the last equation of system (20) we get

[TABLE]

Substituting the relation (21) in the second equation of the system (20) we will have a non autonomous equation:

[TABLE]

As we are interested in scales $\mu\ll m$, $t$ will be negative, $t<0$. We transform Eq.(22) into an autonomous equation using $\lambda=g_{C}(-t)^{1/2}$ and $x=\mathrm{log}|t|$,

[TABLE]

Because the differential equation (23) is homogeneous in $x$, when $x\rightarrow\infty$, that is, for small energies, $\lambda\rightarrow$ constant. We could find the fixed points of $\lambda$ imposing $\partial_{x}\lambda=0$. This is equivalent to solve

[TABLE]

with $\sigma=\frac{\lambda^{2}}{(8\pi)^{2}}$. We can find numerically that the only real and positive root is $\sigma_{0}\sim 21.709$,

[TABLE]

The analytic result seen in the system (25) shows us that both coupling constants goes slowly to zero when the energy decrease. The ghost coupling constant does this by the positives values and the gluon one by the negatives. We know that the latter sign comes from ghost loop of the three-gluon vertex.

As a verification we can observe when we impose the approximation of $g_{C}\approx g_{A}$ at first order we get the beta function $\beta_{g_{c}}$:

[TABLE]

This match with the infrared limit of the beta function, $\beta_{g}$, calculated for the infrared scheme in Peláez (2015).