The generalised infrared structure of the gluon propagator

TL;DR

This paper investigates the infrared singularities of the gluon propagator in gauge theories, showing they imply on-shell states with zero norm relevant for confinement, and tests lattice QCD data for these features.

Contribution

It provides a theoretical framework linking infrared singularities to zero-norm states and confirms their presence in lattice QCD data for the gluon propagator.

Findings

Infrared singularities imply the existence of zero-norm on-shell states.

Lattice data for the gluon propagator is consistent with these generalized poles.

The results have implications for understanding confinement in QCD.

Abstract

Gauge theory correlators are potentially more singular in the infrared than those in non-gauge theories. We determine the implications that these singularities have on the spectrum of the theory, proving that the appearance of generalised poles implies the existence of on-shell states with fixed mass, but zero norm. For quantum chromodynamics these poles have direct relevance for the confinement of coloured states. Using lattice data for the Landau gauge gluon propagator we subsequently test for the presence of these poles, establishing that the data is indeed consistent with such a component.

Figures

Tables

Table 1. Table 1: Chi-squared fit results for D i ​ ( p ) subscript 𝐷 𝑖 𝑝 D_{i}(p) and D i ​ j ​ ( p ) subscript 𝐷 𝑖 𝑗 𝑝 D_{ij}(p) under the different systematic error scenarios with the 64 4 superscript 64 4 64^{4} lattice data.

	Stat. only	Stat. + Shape	Stat. + Poly.
	${\chi }_1^2$/d.o.f. ($p_{\text{max}}$)	${\chi }_2^2$/d.o.f. ($p_{\text{max}}$)	${\chi }_3^2$/d.o.f. ($p_{\text{max}}$)
$D_0\left(p\right)$	$>150$  (0.7)	$>30$ (0.6)	3.5 (0.7)
$D_1\left(p\right)$	$>90$  (1.0)	3.5  (0.9)	2.7  (1.0)
$D_2\left(p\right)$	$>45$  (1.3)	$>45$  (1.3)	3.1  (1.3)
$D_{01}\left(p\right)$	1.9  (1.0)	2.2  (1.0)	2.1  (1.0)
$D_{02}\left(p\right)$	—	2.7  (1.2)	2.8  (1.3)
$D_{12}\left(p\right)$	$>15$  (1.0)	3.4  (1.0)	2.5  (1.0)

Table 2. Table 2: Chi-squared fit results for D i ​ ( p ) subscript 𝐷 𝑖 𝑝 D_{i}(p) and D i ​ j ​ ( p ) subscript 𝐷 𝑖 𝑗 𝑝 D_{ij}(p) under the different systematic error scenarios with the 80 4 superscript 80 4 80^{4} lattice data.

	Stat. only	Stat. + Shape	Stat. + Poly.
	${\chi }_1^2$/d.o.f. ($p_{\text{max}}$)	${\chi }_2^2$/d.o.f. ($p_{\text{max}}$)	${\chi }_3^2$/d.o.f. ($p_{\text{max}}$)
$D_0\left(p\right)$	$>80$  (0.7)	7.3  (0.6)	2.4  (0.7)
$D_1\left(p\right)$	$>25$  (1.0)	1.8  (0.9)	1.7  (1.0)
$D_2\left(p\right)$	$>15$  (1.3)	$>15$  (1.3)	1.8  (1.3)
$D_{01}\left(p\right)$	1.1  (1.0)	1.2  (1.0)	1.2  (1.0)
$D_{02}\left(p\right)$	—	1.6  (1.2)	1.6  (1.3)
$D_{12}\left(p\right)$	6.3  (1.2)	1.7  (1.2)	1.5  (1.2)

Equations

$${\hat{D}}_{\mu \nu }^{abF}(p)$$

$$D(p)=3i{\int }_0^{\infty }\frac{ds}{2\pi }\frac{{\rho }_1(s)}{p^2-s+iϵ}+\sum _{n=0}^{N+1}{g}_n({\partial }^2{)}^n{\delta }^4(p),$$

$${\delta }^{(n)}(s-m_n^2)={\left(\frac{d}{ds}\right)}^n\delta (s-m_n^2),n\ge 1.$$

$$(P^2-m_n^2{)}^n|{\Psi }_n⟩\ne 0,$$

$$(P^2-m_n^2{)}^{n+1}|{\Psi }_n⟩=0,$$

$$|{\tilde{\Psi }}_n⟩\equiv (P^2-m_n^2{)}^n|{\Psi }_n⟩.$$

$$(P^2-m_n^2)|{\tilde{\Psi }}_n⟩=(P^2-m_n^2{)}^{n+1}|{\Psi }_n⟩=0,$$

$$⟨{\tilde{\Psi }}_n|{\tilde{\Psi }}_n⟩=\left(⟨{\Psi }_n|(P^2-m_n^2{)}^n,(P^2-m_n^2{)}^n|{\Psi }_n⟩\right)$$

$$=⟨{\Psi }_n|(P^2-m_n^2{)}^{n-1}(P^2-m_n^2{)}^{n+1}|{\Psi }_n⟩=0,$$

$$D(p)=3{\int }_0^{\infty }\frac{ds}{2\pi }\frac{{\rho }_1(s)}{p^2+s}+\sum _{n=0}^{N+1}{g}_n(-{\nabla }^2{)}^n{\delta }^4(p).$$

$$D_i(p)=\frac{Z_i}{(p^2+m_i^2{)}^{i+1}},i=0,1,2.$$

$$Z_1=48\begin{array}{c}+12\\ -8\end{array}{\text{GeV}}^2,m_1=0.93\begin{array}{c}+0.07\\ -0.06\end{array}\text{GeV},$$

$$D_{ij}(p)=\frac{Z_i}{(p^2+m_i^2{)}^{i+1}}+\frac{Z_j}{(p^2+m_j^2{)}^{j+1}},$$

$$Z_1=51\begin{array}{c}+8\\ -6\end{array}{\text{GeV}}^2,m_1=0.88\begin{array}{c}+0.09\\ -0.06\end{array}\text{GeV},$$

$${\chi }_1^2=(\mathbf{G}-\mathbf{D}{)}^{\text{T}}{\mathbf{C}}^{-1}(\mathbf{G}-\mathbf{D}),$$

$${\chi }_2^2=(\mathbf{G}\cdot \mathbf{f}-\mathbf{D}{)}^{\text{T}}{\mathbf{C}}^{-1}(\mathbf{G}\cdot \mathbf{f}-\mathbf{D})+\frac{(a-\alpha {)}^2}{{\sigma }_a^2}+\frac{(b-\beta {)}^2}{{\sigma }_b^2},$$

$$f_i(p_i,a,b)=\frac{a+p_i^2}{ab+p_i^2},$$

$${\chi }_3^2=(\mathbf{G}\cdot \mathbf{g}-\mathbf{D}{)}^{\text{T}}{\mathbf{C}}^{-1}(\mathbf{G}\cdot \mathbf{g}-\mathbf{D})+\frac{|\mathbf{g}-\mathbf{\gamma }{|}^2}{{\sigma }_g^2},$$

Full text

1 Introduction

As with any quantum field theory (QFT), the fundamental field correlators in quantum chromodynamics (QCD) completely encode both the mathematical and physical properties of the theory. Analysing the structure of these objects is therefore essential for understanding these characteristics. Although the confinement of coloured states is a well-established feature of QCD, the precise mechanism which governs this phenomenon remains a deep open problem. Since the non-observability of colour occurs irrespective of whether the theory contains quarks or not, the structure of the gluon field correlators must therefore play a significant role in determining why confinement occurs in QCD, as opposed to other QFTs. Besides their purely theoretical relevance, gluon correlators are also of importance phenomenologically, entering into the calculation of many non-perturbative observables, including bound-state properties of hadrons such as decay widths and mass spectra [1, 2, 3, 4, 5].

For many years it has been understood that the gluon propagator, and in particular its low-momentum infrared behaviour, is important for unravelling the non-perturbative structure of QCD [6]. Motivated by the idea that gluons could directly give rise to a linearly rising confinement potential, in the late 1970’s Mandlestam [7] hypothesised that the non-perturbative gluon propagator may have a double massless pole as $p^{2}\rightarrow 0$. This idea was supported by phenomenological models at the time [8, 9], as well as initial calculations of the Dyson-Schwinger equations [10]. Nevertheless, the earliest lattice data suggested that the infrared propagator was in fact describable by an ordinary massive pole, corresponding to some sort of effective excitation [11, 12]. From a different perspective, motivated by the issues surrounding gauge copies, Gribov suggested that the gluon propagator could instead have a pair of complex conjugate poles [13]. Based on these ideas it was further argued by Zwanziger [14, 15] that the propagator must also strictly vanish in the limit $p^{2}\rightarrow 0$, and up until a decade ago lattice evidence seemed to support this hypothesis. However, once it became feasible to perform calculations using large lattice volumes, results started to suggest that the propagator was actually non-vanishing at zero momentum [16, 17]. Despite further ideas being put forward to account for this behaviour [18, 19], including several recent studies that propose the existence of complex mass poles [20, 21, 22], it still remains an open question as to how the propagator should behave in the infrared regime, and ultimately what theoretical implications this has for QCD.

Although the precise non-perturbative structure of correlators in QCD remains largely unknown, progress can be made by first establishing the general properties of these objects. In local formulations of QFT quantum fields are defined to be operator-valued distributions, which among other things implies that any field correlator must be a distribution [23, 24, 25]. It turns out that by combining this characteristic with the Poincaré transformation properties of the fields, and the physical assumption that the states in the theory have positive energy, this imposes significant analytic constraints on the structure of correlators, including the existence of a spectral representation [25]. Despite the complications that can arise in gauge theories [26], analogous results also hold for the correlators in QCD, and in particular the gluon propagator [27, 28]. The major goal of this work will be to fully elaborate the connection between the infrared structure of this representation, and the spectrum of the theory. Focussing on the specific case of Landau gauge, we will then use the Euclidean generalisation of this representation in order to test the feasibility of different infrared ansätze using the high precision pure Yang-Mills lattice data from [29], and ultimately help to shed new light on the non-perturbative gluon spectrum.

2 The spectral structure of the gluon propagator

Irrespective of the matter content of the theory, it turns out that the gluon propagator has the following spectral representation [27, 28]:

[TABLE]

where $c_{n}^{ab}$ and $d_{n}^{ab}$ are complex coefficients, and $\partial_{\mu}=\tfrac{\partial}{\partial p^{\mu}}$. Defining: $D(p)=g^{\mu\nu}\widehat{D}_{\mu\nu}^{aa\,F}(p)$ and $\rho_{i}=\rho_{i}^{aa}$, it follows from Eq. (2.1) that the trace component of the Landau gauge gluon propagator $D(p)$ can be written

[TABLE]

with $g_{n}=4c_{n}^{aa}+d_{n}^{aa}$. As detailed in [27] and references within, the representation in Eq. (2.2) follows from very broad assumptions, including the Poincaré covariance of the fields and the positivity of the energy-momentum spectrum. In particular, this implies that Eq. (2.2) holds independently of the non-perturbative subtleties associated with local quantisations such as Landau gauge. One can immediately recognise the first term in Eq. (2.2) as having the familiar-looking Källén-Lehmann spectral form, whereas the second term is purely singular. The potential appearance of non-vanishing singular terms involving derivatives of $\delta^{4}(p)$ arises because in gauge theories such as QCD the gauge symmetry provides an obstacle to the locality of the theory [26]. In order to consistently quantise the theory one therefore has two options: either allow non-local fields, or explicitly preserve locality. A general feature of local quantisations such as Landau gauge is that additional degrees of freedom are introduced into the theory, resulting in a space of states which includes unphysical negative-norm states such as ghosts [30]. It is precisely this loss of state-space positivity which prevents one from ruling out the existence of $\delta^{4}(p)$-derivative components [25].

2.1 Generalised pole terms

Just as the loss of state-space positivity allows for the possibility of purely singular terms, this feature also implies that the gluon spectral density $\rho_{1}(s)$ can potentially have more singular types of components. An important example of such components are generalised pole terms [31, 32]:

[TABLE]

Due to the appearance of derivatives in this expression these terms give rise to contributions in the propagator of the form $(p^{2}-m_{n}^{2}+i\epsilon)^{-n-1}$, which are more singular than ordinary massive poles. One can also demonstrate that each such term is associated with the existence of a finite-norm state $|\Psi_{n}\rangle$, satisfying the following conditions [31]:

[TABLE]

where $P^{2}$ is the squared energy-momentum operator. Eq. (2.4) clearly emphasises that these states are off shell, and therefore represent unphysical degrees of freedom. The appearance of these types of states was first mentioned in [31], although it was never established whether their presence is significant for the on-shell spectrum of the theory. In the remainder of this section we will address this question.

Given $|\Psi_{n}\rangle$ satisfying Eqs. (2.4) and (2.5), consider the state $|\widetilde{\Psi}_{n}\rangle$ defined by the repeated action ($n$ times) of the mass-shell operator $(P^{2}-m_{n}^{2})$:

[TABLE]

Due to Eq. (2.5) and the definition of $|\widetilde{\Psi}_{n}\rangle$ it follows that:

[TABLE]

and therefore $|\widetilde{\Psi}_{n}\rangle$ defines an on-shell state with mass $m_{n}$. Moreover, taking the inner product of $|\widetilde{\Psi}_{n}\rangle$ with itself one finds

[TABLE]

where the last equality is implied by Eq. (2.5). From this one can conclude that the appearance of an unphysical off-shell state $|\Psi_{n}\rangle$ necessarily implies the existence of an on-shell zero-norm state $|\widetilde{\Psi}_{n}\rangle$, with mass $m_{n}$. Since $|\Psi_{n}\rangle$ are intrinsically connected to the appearance of generalised pole terms as in Eq. (2.3), this proves that despite the unconventional structure of these components, they nevertheless have a definite impact on the spectrum of the theory. For many years it has been understood that $\delta(s-m_{0}^{2})$ components in the spectral density with weights of differing signs correspond to on-shell states with either positive or negative norm, depending on the normalisation. Since $\delta^{(n)}(s-m_{n}^{2})$ with $n\geq 0$ are the only discrete mass terms that can appear in any spectral density, this analysis completes the classification of all possible on-shell states. Although the main focus of this work is the gluon propagator, these results apply generally to any QFT that has a space of states with an indefinite inner product.

2.2 Asymptotic behaviour

The relevance of the generalised pole terms in Eq. (2.3) to the asymptotic behaviour of correlators was first recognised in [33, 34]. In particular, it was proven that these terms can potentially cause the correlation strength between clusters of states to grow with distance, a feature which does not occur in ordinary non-gauge theories [23]. This characteristic is often referred to as a violation of the cluster decomposition property (CDP). For clusters of coloured states this characteristic provides a mechanism which can guarantee their absence from the asymptotic spectrum, since a growth in correlation strength between coloured states would prevent the independent measurement of either of these states at large distances [30]. Establishing the spectral structure of the gluon propagator, or in fact the correlators involving any coloured fields, is therefore of particular importance for understanding confinement. In order to violate the CDP the spectral density must not only contain a component as in Eq. (2.3), but this component also needs to influence the behaviour of the correlator at $p^{2}=0$ [32, 33, 34], which requires that $m_{n}=0$. The purely singular terms in Eq. (2.1) involving derivatives of $\delta^{4}(p)$ can also cause the CDP to be violated. However, in this case this occurs not because of the spectral characteristics, but due to the singular properties of the correlator itself.

3 Lattice data fits

Although many different parametrisations of the gluon propagator have been considered over the years, the potential appearance of generalised poles has not yet been fully investigated. In the remainder of this section we will outline infrared fits to the lattice data of [29] testing for the existence of these terms. Since lattice calculations are inherently Euclidean, one must first consider the Euclidean generalisation of Eq. (2.2). Under the assumption of the validity of the Euclidean analytic continuation, the Landau gauge gluon propagator has the form

[TABLE]

Although the pole structure of the Euclidean and Minkowski propagators differ to one another, the spectral density itself remains unchanged. Nevertheless, since lattice simulations inherently possess a finite momentum resolution, the singular terms at $p=0$, if they exist, are presumably not detectable. We therefore disregarded these terms for the purpose of this analysis.

In this analysis we used the $\beta=6.0$ data of [29] with both $64^{4}$ and $80^{4}$ lattices, corresponding to a lattice spacing of $a=0.1016(25)\,\text{fm}$ and physical volumes of $(6.57\,\text{fm})^{4}$ and $(8.21\,\text{fm})^{4}$, respectively. Specific details of the sampling, gauge fixing and other definitions regarding the data can be found in [29] and references within. The first goal of our analysis was to test whether the following single-pole infrared propagator ansätze could fit the lattice data up to some scale $p_{\text{max}}$

[TABLE]

The motivation for performing these specific fits is that it provides a way of testing whether the data is consistent with the appearance of different types of simple isolated poles in the gluon spectral density $\rho_{1}(s)$. If there does indeed exist such a pole at $s=m_{i}^{2}$, together with a gap $(m_{i}^{2},s_{c})$ in the spectrum up to some continuum threshold $s_{c}$, the infrared structure of the propagator will be dominated by these poles111A simple argument for this characteristic can be made by performing a small momentum expansion of the spectral representation, and observing that the continuum contributions are increasingly suppressed as the gap $(m_{i}^{2},s_{c})$ becomes larger.. In particular, $D_{i}(p,m_{i})$ would provide a good fit to the data up to some $p_{\text{max}}<\sqrt{s_{c}}$. For these fits to make physical sense one requires that $p_{\text{max}}>m_{i}$. It is important to emphasise here that the ansätze in Eq. (3.2) correspond to real mass poles of differing orders, in contrast to the complex mass pole ansätze often considered in the literature [20, 21, 22].

Due to the difficulty in precisely accessing the systematic uncertainties of the lattice data we decided to perform the fits using three different choices of uncertainties, with increasing levels of conservatism. For the first case we considered the statistical errors only, which are $\mathcal{O}(10^{-3})$. For the second case, we took into account both the statistical errors as well as an energy-dependent systematic error. This shape uncertainty was derived using gluon propagator data with different lattice spacings in [35], and used an empirical function that can describe the differences between these results and hence capture the general shape changes which occur due to systematic effects. This approach is particularly well-suited to cases where the data is heavily influenced by spectrum shape effects, such as in the assignment of nuclear reactor flux errors for neutrino experiments [36]. Lastly, we adopted a data-independent choice for the systematic errors: we assumed that the propagator data could be modified by a fourth-order polynomial in $p^{2}$, with deviations at the $2.5$% level. This was the most conservative choice for the systematic error since it allows for significant shape changes. In all of these fitting scenarios the goodness of fit was assessed using a Chi-squared minimisation procedure. The Chi-squared statistic definitions together with the fit results for the $64^{4}$ and $80^{4}$ lattice data are provided in the appendix. Since the $p=0$ data points in both the $64^{4}$ and $80^{4}$ data sets have a considerable systematic uncertainty relative to the other points, due to their calculation, we decided to omit $p=0$ from each of the final fits. Nevertheless, the omission of this point did not result in a significant modification to the parameter values in the fits.

Comparing all of the fits in the different systematic error scenarios we found that $D_{1}(p)$ was the only ansatz which could provide a consistent fit to the data. In order to assess the robustness of these fits we performed a two-dimensional Chi-squared minimisation of the ansätze parameters, and also tested the sensitivity of these fits to the momentum cutoff $p_{\text{max}}$. We found that neither $D_{0}(p)$ nor $D_{2}(p)$ provided a reliably stable fit for both data sets as $p_{\text{max}}$ varied, and that the fitted parameter values were highly sensitive to the choice of systematic errors. The $D_{1}(p)$ fits on the other hand resulted in a good fit to the data for both choices of (non-vanishing) systematic errors. Although the overall quality of the fits improved when using the larger volume data, the parameter values remained stable across the $64^{4}$ and $80^{4}$ fits. The Chi-squared values in $(Z_{1},m_{1})$ parameter space for the $80^{4}$ data with the polynomial shape choice of systematics are plotted in Fig. 1. The parameter values at the global minimum of this fit are

[TABLE]

where the uncertainties indicate a 1$\sigma$ variation. We found that this fit remained convergent and physically consistent ($p_{\text{max}}>m_{1}$) up to $p_{\text{max}}=1.5\,\text{GeV}$.

In principle it is also possible that the gluon spectral density could contain multiple types of poles. To investigate this possibility we therefore considered the simplest case of the two-term linear combinations of single pole components

[TABLE]

where $i,j=0,1,2$ and $i\neq j$. Overall we found that $D_{01}(p)$ provided a very good description of the lattice data, even when only statistical errors were included. The best fit was obtained with the $80^{4}$ data, yielding the corresponding generalised pole parameter values

[TABLE]

and: $Z_{0}=-3.4\begin{subarray}{c}+0.5\\ -0.5\end{subarray}\,\text{GeV}^{2}$, $m_{0}=0.36\begin{subarray}{c}+0.04\\ -0.04\end{subarray}\,\text{GeV}$. We found that this fit remained convergent and physically consistent ($p_{\text{max}}>\{m_{0},m_{1}\}$) up to $p_{\text{max}}=1.4\,\text{GeV}$. Although the $D_{02}(p)$ and $D_{12}(p)$ ansätze appeared to give a reasonable fit to the data based on the $\chi^{2}$/d.o.f. values alone, the fitted mass parameters were found to be degenerate, possessing a sensitivity to the momentum cutoff $p_{\text{max}}$. An important conclusion from this analysis is that the best fit $Z_{1}$ and $m_{1}$ parameter values from the $D_{01}(p)$ fit coincide (within errors) with those obtained in the $D_{1}(p)$ only fit. Moreover, adding a $D_{1}(p)$ component to either $D_{0}(p)$ or $D_{2}(p)$ ended up leading to a significant improvement in the quality of these fits. Taken together, these results support the hypothesis that the gluon propagator data is consistent with the existence of an infrared $D_{1}(p)$ contribution. The $D_{1}(p)$ and $D_{01}(p)$ ansätze with the best fit parameter values in Eqs. (3.3) and (3.5), respectively, are plotted together with the $80^{4}$ lattice data points in Fig. 2.

Due to the spectral structure in Eq. (3.1), a $D_{1}(p)$ contribution means that $\rho_{1}(s)$ contains a generalised component $\mathcal{Z}_{1}\delta^{\prime}(s-m_{1}^{2})$, where $\mathcal{Z}_{1}=\tfrac{2}{3}\pi Z_{1}$. With the parameters given in Eq. (3.5) it follows from the previous theoretical analysis that the appearance of such a component implies the existence of an on-shell zero-norm state in the spectrum, with mass $m_{1}=0.88\,\text{GeV}$ and weight $\mathcal{Z}_{1}$. An important characteristic of this component is that it directly satisfies the Landau gauge superconvergence constraint [28, 37, 38]: $\int ds\,\rho_{1}(s)=0$, not due to any continuous oscillatory-type behaviour [39], but because of the generalised singular properties of the $\delta^{\prime}$ distribution. Since $m_{1}>0$ this also suggests that the gluon spectral density does not induce a singular enough behaviour near $p^{2}=0$ to result in a violation of the CDP. It could well be the case that a violation is indeed caused by the non-vanishing of one of the coefficients $g_{n}$ in Eq. (3.1), but this remains to be seen. A non-violation of the CDP would imply that clusters of states created from single gluonic fields should decorrelate as they are separated further apart. Whilst this might initially appear surprising, one should note that the gluon correlator itself is not gauge-invariant, and so a non-violation of the CDP in this specific case does not contradict the physical expectation that asymptotic coloured states are prohibited.

4 Conclusions

In locally quantised gauge theories the field correlators can potentially contain more singular generalised-pole contributions in their spectral densities. In this work we establish for the first time the connection between these types of components and the spectrum of the theory. In particular, we prove that the presence of generalised poles implies the existence of on-shell states in the spectrum with fixed mass, but zero norm. Despite their non-physicality, it turns out that these states can fundamentally alter the large-distance behaviour of correlators, which is important in the context of confinement. With these theoretical motivations in mind we analysed the high-precision Landau gauge gluon propagator lattice data from [29] in order to assess whether the data can accommodate generalised-pole contributions. Overall, we found that the data was consistent with the appearance of a massive single-derivative of delta component in the spectral density, and hence the existence of a corresponding zero-norm state in the spectrum. The non-vanishing mass of this state suggests that the spectral structure of the propagator alone is not sufficiently singular to guarantee the confinement of states involving single gluonic fields. Since generalised spectral components are also potentially present in other correlators, this work opens up a new direction for understanding the infrared structure of QCD.

Acknowledgements

S. Li was supported by the U.S. Department of Energy under contract DE–AC02–76SF00515, and P. Lowdon by the Agence Nationale de la Recherche under project No. ANR–18–ERC1–0002. P. J. Silva acknowledges support from FCT under the contracts SFRH/BPD/40998/2007 and SFRH/BPD/109971/2015, and both O. Oliveira and P. J. Silva acknowledge support from FCT under the contracts UID/FIS/04564/2016 and CERN/FIS-COM/0029/2017. This work was granted access to the HPC resources of the PDC Center for High Performance Computing at the KTH Royal Institute of Technology, Sweden, made available within the Distributed European Computing Initiative by the PRACE-2IP, receiving funding from the European Community’s Seventh Framework Programme (FP7/2007-2013) under grand agreement no. RI-283493. The use of Lindgren has been provided under DECI-9 project COIMBRALATT. We acknowledge that the results of this research have been achieved using the PRACE-3IP project (FP7 RI312763) resource Sisu based in Finland at CSC. The use of Sisu has been provided under DECI-12 project COIMBRALATT2. We also acknowledge the Laboratory for Advanced Computing at the University of Coimbra (http://www.uc.pt/lca) for providing access to the HPC resource Navigator.

Appendix A Goodness of fit results

The goodness of fit results for the $D_{i}(p)$ and $D_{ij}(p)$ ansätze with the $64^{4}$ and $80^{4}$ lattice data of [29] are found in Tables 1 and 2, respectively. In each fit a value of $p_{\text{max}}=1\,\text{GeV}$ was initially chosen to determine whether a convergent and physically consistent ($p_{\text{max}}>m_{i}$) fit could be achieved. If not, the value of $p_{\text{max}}$ was lowered or raised until these conditions were satisfied. The values of $p_{\text{max}}$ in the tables reflect this final choice of cutoff. The $D_{02}(p)$ fits to both data sets were found to not be consistent when only statistical errors were included, regardless of the value of $p_{\text{max}}$.

Appendix B Chi-squared definitions

For each propagator ansatz tested we used three different choices for the systematic error, each of which requires a separate Chi-squared statistic definition in order to determine the goodness of fit. In the least conservative scenario for which we only consider the statistical errors, we used the following Chi-squared definition:

[TABLE]

where $\bm{G}$ is the $n$-dimensional vector of lattice data at the discrete momenta $\{p_{1},...,p_{n}\}$, and $\bm{D}$ is the vector of the fitted propagator functional form evaluated at these same momenta. The covariance matrix $\bm{C}$ quantifies the statistical correlations between different field configurations in the data, and is highly diagonal. In the limit of zero cross-correlations $\bm{C}$ becomes the matrix of statistical variances at each lattice point, and Eq. (B.1) reduces to the standard Chi-squared form. The minimum $\chi_{1}^{2}$ follows a $\chi^{2}$-distribution with $n-m$ degrees of freedom, where $m$ is the number of free parameters in the propagator ansatz.

In the second scenario we took into account both the statistical error and an energy-dependent systematic error. The corresponding Chi-squared statistic is defined by

[TABLE]

where now the lattice data is modified by an effective shape function $\bm{f}=\bm{f}(\{p_{i}\},a,b)$. This fitting approach is particulary well-suited to cases where the data is heavily influenced by energy-dependent effects, for example in assigning a shape error to a nuclear reactor flux for neutrino experiments [36]. Since the purpose of $\bm{f}$ is to simulate the sort of shape changes that can occur due to systematic effects, such as the appearance of Gribov copies [29], we used the propagator data from [35] in order to establish a reasonable functional form. In particular, $\bm{f}$ was established such that the shape of the data sets for the different lattice spacings in [35] could be transformed into one another after multiplication by this function. With this procedure the components of $\bm{f}$ were determined to have the functional form

[TABLE]

and the priors: $\left(\alpha=0.21,\sigma_{a}=0.016\right)$ and: $\left(\beta=1.0,\sigma_{b}=0.1\right)$ were found to cover the full range of shape differences observed in [35]. Due to the definition in Eq. (B.2), the minimum $\chi_{2}^{2}$ has $n-(m+2)$ degrees of freedom. In the most conservative systematic error scenario we used the following Chi-squared statistic

[TABLE]

where the systematic shape function $\bm{g}$ was chosen to be a forth-order polynomial in $p^{2}$ (with 4 free coefficients), with priors $\bm{\gamma}=1$ and $\sigma_{g}=0.025$ simulating a $2.5$% variation from the data. The corresponding minimum $\chi_{3}^{2}$ has $n-(m+4)$ degrees of freedom.