{'fraction_non_alphanumeric': 0.03955168358653131, 'fraction_numerical': 0.011159910720714235, 'mean_word_length': 4.421285454072339, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Covariant variational approach to Yang-Mills Theory: effective potential of the Polyakov loop

We compute the effective action of the Polyakov loop in SU(2) and SU(3) Yang-Mills theory using a previously developed covariant variational approach. The formalism is extended to background gauge and it is shown how to relate the low order Green's functions to the ones in Landau gauge studied earlier. The renormalization procedure is discussed. The self-consistent effective action is derived and evaluated using the numerical solution of the gap equation. We find a clear signal for a deconfinement phase transition at finite temperatures, which is second order for SU(2) and first order for SU(3). The critical temperatures obtained are in reasonable agreement with high precision lattice data.


I. INTRODUCTION
The low energy sector of quantum chromodynamics (QCD) and, in particular, its phase diagram continue to be one of the most actively researched topics in elementary particle physics.
Recently, large experimental facilities such as the large hadron collider (LHC) offer the possibility to study strongly interacting matter under extreme conditions, such as large temperatures and baryon densities. On the theoretical side, lattice simulations can be used to obtain numerical ab-initio solutions of QCD in a variety of settings, but restricted to zero or small chemical potential due to the sign problem. Alternative continuum or functional methods, which are not plagued by this issue, are thus of particular intereset. These methods comprise the functional renormalization group (FRG) flow equations [1,2] Dyson-Schwinger equations (DSE) [3][4][5] and variational methods [6][7][8]. In all of these approaches, the simplest objects to be studied are the low-order Green's functions of gluons, ghosts and quarks. Since these quantities are, however, gauge-dependent, one has to fix a gauge in order to get meaningful results. The DSE and FRG approach have a covariant setup and are thus predominantly studied in Landau gauge, while the variational approach has been formulated and worked out in the Hamiltonian approach using Coulomb gauge.
In the covariant case, the Becchi-Rouet-Stora-Tyutin (BRST) symmetry can be used to guide the analsysis, and the Kugo-Ojima criterion [9,10] claims a direction connection (based on an unbroken BRST symmetry) between the propagators in Landau gauge and color confinement. In reality, however, the functional methods in Landau gauge allow for a whole family of solutions for the propagator which are parametrized by different infrared boundary conditions [11]. In particular, the so-called decoupling solution shows excellent agreement with the lattice data, even though it (softly) breaks BRST symmetry. The Kugo-Ojima criterion is no longer applicable in this case, and the question of color confinement must be resolved in other ways. One particularly simple approach is to study the effective action for the Polyakov loop, a maximally extended temporal Wilson line that closes via the periodic boundary conditions in the euclidean time direction. The Polyakov loop is a direct order parameter for confinement [12,13] and can therefore be used to study the deconfinement phase transition that occurs in the QCD phase diagram at increasing temperatures. Naturally, the FRG [14][15][16], the DSE [17] and the Hamiltonian approach in Coulomb gauge [18][19][20] have been used to study the properties of this phase transition and compared it with findings from high precision lattice studies.
In the pure Yang-Mills case, both approaches predict a transition which is second order for the color group G = SU(2) and first order for G = SU (3). Moreover, the phase transition temperatures are in reasonable quantitative agreement with the lattice.
In Refs. [6] we have proposed an alternative covariant variational approach which combines the theoretical simplicity of the covariant setup with the efficiency of a variational method.
It yields a closed set of integral equations that can be renormalized by conventional counterterm techniques. In addition, the variational character automatically optimizes any ansatz for the low-order Green's functions, which therefore are predicted quite accurately even when based on simple Gaussian ansätze for the trial measure. In particular, the quantitative agreement with the lattice propagators at T = 0 is excellent [6], and also the corrections due to finite temperatures are reproduced in accord with the lattice data [21]. Due to the inherent breaking of BRST symmetry, the Kugo-Ojima ciriterion does not apply and the question of color confinement is inconclusive. As mentioned above, it is therefore important to study the effective potential of the Polyakov loop to see if our method, combined with the simple Gaussian ansatz, is strong enough to predict the properties of the QCD phase diagram in this region. In the present paper, we study this question at vanishing chemical potential and for pure Yang-Mills theory, i.e. in the absence of dynamical quarks. Both these restrictions are intended to be lifted in future studies.
Our paper is organized as follows: In section II, we recall the physical interpretation of the Polyakov loop and its relation to center symmetry and to color confinement. Section III extends our variational approach to the background gauge which is most convenient for the study of the Polyakov loop in our setup. In particular, the relevant variational kernels need not be re-computed, but instead are inherited from the known kernels in Landau gauge, as is demonstrated in section IV. The following section V computes the self-consistent effective action for the Polyakov loop in our approach and shows how to turn these results in an efficient numerical calculation. Section VI describes our numerical setup and presents our results, including a comparision with high-precision lattice simulations. Finally, we conclude in section VII with a short summary and an outlook on the future development of our approach.

II. CENTER SYMMETRY AND THE POLYAKOV LOOP
Let us recall the role of the center symmetry and the Polyakov loop as an order parameter for the deconfinement phase transition in pure Yang-Mills theory. At finite temperatures, the coupling to the heat bath excites higher energy states and the static quark potential V β (r) at inverse temperatuer β must be defined as the change in free energy due to a heavy qq-pair immersed in the thermal background [13]. Within the imaginary time formalism, V β (r) can extracted from non-contractible Polyakov loops winding around the compactified Euclidean time direction: where P is path ordering and the trace in the respective quark representation is normalised such that tr1 = 1. In a pure Yang-Mills theory, the insertion of a static heavy quark at position x in the thermal background state is, in fact, equivalent to the insertion of a Polyakov loop in the finite-temperature partition function [12]. More precisely, the free energy of a static qq-pair located at positions x and y is obtained from the correlator of two Polyakov lines, where r = |x − y| is the distance of the quarks and the star on L denotes complex conjugation.
Furthermore, H phys is the gauge invariant Yang-Mills Hamiltonian and the partition function Z(β) serves to normalize the expectation value. Since the pure phases on both sides of the deconfinement transition are expected to show asymptotic clustering, we find From eqs. (2) and (3), we conclude that a vanishing average for the Polyakov line immediately leads to L(x) · L * (y) → 0 for large qq separations, i.e. a static potential increasing with the distance between the charges. Thus, L = 0 implies color confinement, while L = 0 implies that the static qq-potential remains finite, and a finite energy is sufficient to separate the qq-pair (deconfinement). 1 The confined phase is characterized as an unordered phase with a higher degree of symmetry. In fact, the vanishing of the Polyakov loop can be seen as a 1 The formal association of a free energy of a fictious single quark, exp(−βF q ) = L , is, however, doubtful: since a single color charge in the fundamental representation cannot be screened by gluons, Gauss' law cannot be satisfied within the vacuum sector of gauge-invariant states, and the physical traces are void.
consequence of center symmetry which maps each Polyakov loop L → z L by a center element z of the color group, but leaves the action invariant. If unbroken, L = z k L for each center element z k , and thus L = 0 because ∑ k z k = 0. An unbroken center symmetry hence implies confinement. Conversely, a broken center symmetry implies L = z k L for each center element, and thus L = 0, which entails deconfinement.
In a lattice discretisation, the center symmetry can be realized by multiplying all temporal links emenating from a fixed time slice by the same center element z. A continuum version of this construction employs the so-called Polyakov gauge, where {T a } are a maximal set of commuting generators of the Lie algebra of the color group G, which span the Cartan subgroup (maximal torus) H ⊂ G, while the remaining generators span the coset G/H. The Polyakov loop becomes and requires no path ordering. For the color group G = SU(2), the maximal torus is spanned by T 3 = σ 3 /(2i) (where σ 3 is the diagonal Pauli matrix) and the Polyakov loop can be expressed as The fundamental modular region is φ(x) ∈ [0, 1] and the effective potential for the observable φ is even in φ (due to Weyl reflection symmetry) and periodic outside the range [−1, 1]. The In Refs. [15] and [14], it was argued that not only the Polyakov loop, but also the expectation value of its logarithm (i.e. φ ∼ A 3 0 ) can be an order parameter for confinement. Let us briefly repeat this reasoning: In the confined phase, the traced Polyakov loop vanishes and center symmetry is unbroken. This implies that configurations φ(x) and 1 − φ(x) are equally likely, i.e. φ(x) = 1 − φ(x) or φ(x) = 1/2 at any space position x. As a consequence confined (center symmetric): In the deconfined phase, L = cos(πφ) > 0 within the fundamental domain, and hence | cos(πφ)| > 0. (The latter conclusion follows because | cos(πφ)| = 0 would imply | cos(πφ)| = 0 almost everywhere, hence cos(πφ) = 0 almost everywhere and thus cos(πφ) = 0 in contradiction to the assumption.) Since the function φ → | cos(πφ)| is concave on the fundamental domain, Jensen's inequality can be applied [14,15] and we have deconfined (center broken): Thus, for the gauge group G = SU(2) at least, we can also use φ ∼ A 0 itself as an order parameter for confinment in Polyakov gauge. 2 The considerations above were valid in the (non-covariant) Polyakov gauge. Since the effective action for a gauge-variant operator such as A 0 is also gauge-variant, Γ[ φ ] will differ when changing gauges. For instance, φ ∼ A 0 = 0 in Landau gauge due to Lorentz and global color invariance, and nothing can be learned from A 0 . In order to have a covariant theory with A 0 = 0, we have to adopt the background gauge formalism described below. Even then, the effective action for the background field a 0 = A 0 does not necessarily agree with the one in Polyakov gauge discussed above. However, it has been argued in Refs. [14,15] that the background gauge formalism can be used to detect gauge invariant features such as the location and the order of the phase transition, if the background field a µ itself is taken to have only an a 0 component which obeys the Polyakov gauge condition. Moreover, the corresponding potentials in background and Polyakov gauge are very similar numerically, both on the lattice and in the functional renormalization group approach [15]. In the latter case, it was argued that this coincidence might be expected, as the relevant quantum fluctuations in background gauge include the ones in Polyakov gauge, if the backgrond field a 0 itself is taken in Polyakov gauge. 3 We will therefore adopt this formulation and employ the covariant variational approach [6,21] to compute the effective action for a time-independent, Abelian background field a 0 (x), and check if the formalism is strong enough to predict the known properties of the deconfinement phase transition for G = SU(2) and G = SU(3).

III. THE COVARIANT VARIATIONAL APPROACH IN BACKGROUND GAUGE
In the background field formalism, we choose an arbitrary background field a µ and decompose the full gauge connection A µ = a µ + Q µ , where Q µ is the fluctuating field. We fix the quantum gauge symmetry (θ = infinitesimal gauge angle) by imposing the background g.f. condition where d µ ≡ D µ (a) is the covariant background deriative. It is then easy to see that the gaugefixed action, the corresponding path integral and hence the quantum effective action Γ a [Q] for the classical fluctuation field Q µ = Q µ a (with the fixed background field a µ as a parameter) are all invariant under the background transformation Also Γ a [Q] = Γ[a + Q], where Γ is now the usual effective action, i.e. the generating functional of proper functions for the original gauge field A µ , in the (slightly unusual) gauge If we set the classical field Q µ = 0, we find that Γ[a] = Γ a [0] is a gauge-invariant functional of the background field a µ .
Next, we have to adapt this formalism to the variational approach of Refs. [6,21]. We keep a µ as a parameter and define the variational approach to the path integral for the fluctuation field, Here, the variation is over all normalised path integral measures µ for the fluctuation field, S fix is the Yang-Mills action with the gauge fixing term corresponding to the background gauge condition (10), and is the entropy of the measure dµ ≡ dQ ρ(Q) relative to the natural measure on the gauge orbit [6], i.e. the Faddeev-Popov determinant J . Before investigating this further, it should be stressed once again that the effective action eq. (13) is not the same as the effective action for the original gauge field A µ (in the unusual gauge eq. (12)), because the vev of the fluctuation Q µ , rather than the full A µ , is held fixed. As mentioned earlier, and we set the classical fluctuation field Q µ = 0 from this point on. In particular, this means that all trial measures for the variation principle (13) must be centered in the origin, Q µ = 0.
To further work out the variational approach, we write the measure µ as a modification of the flat measure, dµ(Q) = dQ ρ(Q) and denote the Faddeev-Popov determinant in background gauge by where the hat denotes the adjoint representation of the color group, i.e.D ab µ = ∂ µ δ ab +Â ab µ = ∂ µ δ ab − f abc A c µ with the antisymmetric structure coefficients f abc , and likewise for the covariant background derivatived µ =D µ (a). Then the free action F(µ) for a measure µ is given explicitly by where ξ is a gauge-fixing parameter. We cannot do the variation over the full set of measures dµ(Q). Instead, we will restrict our investigations to Gaussian trial measures, since these have already proven to be successful in the description of the low-order Green's functions, both at zero [6] and finite temperature [21]. More specifically, our ansatz is where ω is a variational kernel to be determined by minimizing the free action eq. (16). Furthermore, the Gaussian is centered at A µ = a µ (or Q µ = 0) as discussed above. Note that the FP determinant in the first line is evaluated at the original field A as in eq. (15), since this is the natural measure on the original gauge orbit. After the shift to the fluctuation field in the second step, the measure remains centered, Q µ = 0, even when the FP prefactor is taken into account. This is clear because fluctuations about the background field with opposite signs should be equally likely on the gauge orbit, and it will be verified explicitly when J is treated in curvature approximation [23], cf. eq. (19) below. A side effect of the curvature approximation is that the variational parameter α, which controls the localization of the measure near the Gribov horizon, becomes immaterial and the trial measure is exactly Gaussian for all α. This is benificial, because it allows to use Wick's theorem to evaluate expectation values. Since the free action (16) depends implicitly on the background field a, we must allow for the variational kernel ω in the ansatz (17) to be non-transversal and non-diagonal in color.
Next we insert the ansatz (17) in the free action eq. (16). After a short calculation, we find where an irrelevant field-independent constant has been dropped. To proceed we employ the curvature approximation where both sides of eq. (19) vanish at Q = 0 and the curvature kernel χ is determined to give the optimal approximation of the lhs by a Gaussian (up to two-loop order). The explicit calculation for χ leads to an integral equation that relates it to the ghost form factor η(k), which in turn is coupled to the gluon propagator through its own DSE, cf. eq. (A4). The detailed treatment of the ghost sector in our approach has been described in Refs. [6,23] and the resulting integral equations can also be found in appendix A. It should be noted at this point that the expression on the lhs of eq. (19) arises from the corresponding Landau gauge expression by the replacement 1 ∂ µ →d µ . This will allow us to relate the curvature in eq. (19) to the one in Landau gauge further below.
The normalization constant N can now be computed explicitly and we end up with the final form of our trial measure, whereω ≡ ω + (1 − 2α)χ. This is now suitable for variation. The free action becomes This is exactly the free action as found in Landau gauge [6], except for some shifts by the background field in strategic places. We will discuss the effect of these shifts in the next section.

IV. RELATING BACKGROUND GAUGE AND LANDAU GAUGE
Up to this point, the background field has been arbitrary. Since we are only interested in the effective potential for the Polyakov loop, we can take the background field to have only a temporal component, which can be chosen constant and Abelian, a µ = δ µ0 a with a = a c T c .
We also observe that the entire free action (21) including the curvature approximation (19) depends on the background field only through its (adjoint) covariant derivatived = ∂ +â.
This was already observed for the FP operator above, and it also holds for the gauge fixing term and the YM action because It is therefore convenient to go to an adjoint color base in which the matrixd orâ ab = (2), this is simply the spherical basis introduced in Ref. [20]: The corresponding eigenvectors are We use greek letters σ, τ . . . ∈ {−1, 0, 1} to denote color components in the new spherical basis, and latin letters a, b ∈ 1, 2, 3 for the usual Cartesian components of adjoint SU(2) color.
Any matrix M ab in the adjoint representation can then be transformed to the new basis using 4 In particular, this unitary transformation diagonalisesT 3 . The advantage of the new color basis is that the covariant background derivative is diagonal, In momentum space, i.e. when acting on e ipx , this operator becomes and we find the d'Alembertian A similar root decomposition exists for G = SU (3), where the single index σ ∈ {−1, 0, 1} must be replaced by a two-component vector (the so-called root vector) σ = (σ 3 , σ 8 ), because SU (3) has rank 2 and hence two Cartan generators. This is worked out in appendix B. One finds that the SU(3) algebra can be decomposed into three SU(2) sub-algebras, and the SU(3) potential can thus be written as a suitable "skew" superposition of SU(2) potentials. For simplicity, we will therefore continue to use SU(2) notation until we present explicit SU(3) results.
Returning to our variational ansatz eq. (20), we note that the kernelω cannot be taken color diagonal in the original (Cartesian) basis due to the background field, but it can be taken diagonal in the spherical basis because this diagonalisesd µ , which is the only way in which the background field enters. Furthermore, the constant background field does not break translational invariance and we can Fourier transform the kernel as usual, so that We can relate the kernelω σ µν (p) to the corresponding kernel in Landau gauge. To see this, recall that (i) the kernels reduce to their Landau gauge counterpart when the background field vanishes (a = 0), and (ii) the background field a enters the Yang-Mills action and the Faddeev-Popov determinant only through its (adjoint) covariant derivatived µ . Furthermore, there are no derivatives in the free action which are not part of a covariant derivatived µ . This observation make it intuitively clear that we only need to replace ∂ µ →d µ to obtain a solution of the gap equation for a = 0, and this amounts to a shift in the momentum of the Fourier transformed kernel. We will confirm this expectation in the next section.
More precisely, the replacement ∂ µ δ ab →d ab µ is a matrix equation and we need to go to the spherical color basis in whichd µ is diagonal. For every root vector σ, the replacement ∂ µ →d µ amounts to the shift eq. (24) of the momentum, and we haveω σ µν (p) =ω µν (p σ ). Note that this shift includes the Lorentz structure as well:ω µν is (ordinary) 4-transversal, whileω σ µν is background transversal, i.e. proportional to the projector The same argument holds for the curvature, since the lhs of the defining equation (19) contains the background field only throughd µ . More details on a similar argument in Coulomb gauge can be found in Ref. [20].
It should finally be noted that the relation between Landau and background gauge derived here only holds up to two-loop order in the free action. In higher orders, it is fairly easy to see that the free action must contain terms in which ordinary derivatives enter in other combinations than justd µ . Because the effective action is gauge invariant (in background gauge), such derivatives must act on, or appear in, gauge invariant operators such as L. This was already observed in Ref. [14] where it was also argued that such corrections are small or negligable.

V. THE EFFECTIVE ACTION FOR THE POLYAKOV LOOP
Let us go back to the expression eq. (21) for the free action. To keep the formulars readable, we employ an obvious shorthand notation where a roman digit stands for the combination of spactime, Lorentz and adjoint color index, 1 ≡ (x, µ, a), 2 ≡ (y, ν, b) etc. The expectation values of the Yang-Mills and gauge fixing action in the Gaussian trial measure eq. (20) can be worked out using Wick's theorem. After some algebra, we obtain Here, the kernelsω and χ are matrices in color, spacetime and Lorentz space; as for the latter, they multiply the transversal Lorentz projector (27). The kernels γ are the bare proper vertices from the Yang-Mills and gauge fixing action, Furthermore, the second to last term in eq. (28) is the contribution from the longitudinal degree of freedom (which receives no radiative correction), and the last term comes from the two perturbative ghost degrees of freedom, J [a] = det −d µdµ , because the curvature gives only the ratio of the full and free background FP determinant, cf. eq. (19). The last two terms in eq. (28) decouple from the dynamics, but they are required to get the counting of the degrees of freedom right.
The optimal gluon propagator can be found by minimising eq. (28) with respect to the variation kernel, We define the tadpole contraction of the 4-gluon vertex by The permutations in the free energy eq. (30) can also be worked out using the Jacobi identity for the structure constants, and they are found to give just an overall factor of 2 as compared to the tadpole contraction. The result can therefore be written in the form This is formally identical to the expression in Landau gauge. The difference is that the kernel γ(1, 2) and the FP determinant (and hence the curvature χ) now depend on the background gauge field a. More precisely, they are obtained from their Landau gauge counterpart by replacing each partial derivative by a background-covariant derivative, 1 ∂ µ →d µ . For a constant Abelian background field a, this amounts to the shift eq. (24) in momentum spaces.
As a consequence, we now conclude that the solutionω of the gap equation in background gauge is also obtained from its Landau gauge counterpart by the same momentum shift. This is very beneficial as we do not have to re-solve the gap equation when varying the background field.
One remark concerning the last term in the gap equation (32) is in order. In Refs. [6] and [21], we have simply discarded this term as a higher loop effect. A more systematic scheme would introduce a formal loop counting parameter λ as a prefactor in the exponent of the trial measure (20), which amounts to counting (fully dressed) internal gluon lines. In this scheme, our results (including the curvature approximation) are correct to order O(λ) in the gap equation, and O(λ 2 ) in the free energy. However, it is then intuitively clear that the last term in the gap equation (32) is of the same order as the previous one and cannot be discarded.
More precisely, a direct calculation (cf. appendix A) reveals that where χ 0 is a field-and momentum-independent quadratic divergence that is eventually removed by renormalization. 5 Since the tadpole term M 2 0 is a similar quadratic divergence, it is prudent to combine the two in the form This results in the free action and gap equation, respectively, Inserting the solution of the gap equation into F[a] yields the (unrenormalized) self-consistent effective action of the background field a, The subtraction (χ − χ 0 ) in the first term already removes the quadratic divergence in the curvature, and the sole effect of the renormalization is to further cancel the subleading logarithmic divergence in the curvature.
To see how this comes about, we need to introduce counter terms. The gap equation only requires a gluon mass and field counter term, while a third (ghost wave function) counter term is necessary for the ghost form factor which enters the integral equation for the curvature. 6 The renormalized gap equation thus reads The mass counter term δM 2 serves to cancel the quadratic divergence χ 0 in the curvature (δM 2 + χ 0 = 0), while the subleading logarithmic divergence is compensated by the wave function renormalization δZ A . (The finite pieces in the counter terms are fixed by precise renormalization conditions on the gluon propagator as described in Ref. [21].) With these arrangements, we findω were Z is a finite part of δZ A determined by the renormalization conditionω(µ) = Z µ 2 at a large scale µ 1, cf. next section. Note that the renormalized curvature in the gap equation is unambiguously determined by the counter terms. This would not have been possible if we had discarded the last term in eq. (32), as we could not combine the curvature and tadpole contribution in this case. Since only that combination has a unique (mass) counter term, we would then be faced with the problem of how to "distribute" the finite parts of the mass counter term onto M 2 0 and χ to define these quantities individually. The renormalized free action using the definition eq. (39) is If we finally employ the renormalized gap equation (38), we find the final expression for the self-consistent effective action of background field, Here, the spacetime volume V 4 has factorized because of translational invariance, and we are really computing the effective potential Γ/V 4 . Furthermore, the variational kernelω(k) (and likewise the curvature) is color diagonal and transversal in Landau gauge,ω ab µν (k) = δ ab t µν (k)ω(k) with a scalar variation kernelω(k). The color trace yields the prefactor (N 2 − 1), while the Lorentz trace gives the factors of 3 for the transversal degrees of freedom, one for the longitudinal degree of freedom, and (−2) for the ghost degrees of freedom. Since the curvature describes the deviation from the free-field FP determinant only, we must also include the contribution from the free ghost degrees of freedom J [0], which is the last term in eq. (42). Furthermore, the longitudinal gluon receives no radiative corrections beyond one loop and we have ω (k) = k 2 .
In the next step, we must put the system at finite temperatures since we want to study the deconfinement phase transition. This amounts to imposing periodic boundary conditions along the compactified euclidean time direction of length β, which is the inverse temperature.
As a consequence, the momentum integrals are i.e. the integral over k 0 is always understood as a discrete sum over the Masubara frequencies k 0 = ν n = 2πn/β. Moreover, the heat bath singles out a rest frame and the overall O(4) invariance of the theory is broken. As a consequence, the gluon propagator (and also the curvature) comes in two distinct Lorentz structures, which are both 4-dimensionally transversal, but also 3-dimensionally longitudinal or transversal, respectively. The corresponding kernels ω ⊥ andω (and χ ⊥ , χ ) were computed in Ref. [21]: as the temperature increases, there is a moderate suppression of the gluon propagators, and a slight enhancement of the ghost form factor. The temperature sensitivity is larger in the components longitudinal to the heat bath, but this affects only the infrared region. Generally, the ghost and gluon propagators are only moderatley changed by temperatures well up to twice the critical temeprature 2T * , and we can even discard their implicit temperature dependence and use the zero-temperature kernels throughout. On the one hand, this reduces the numerical effort dramatically, since the finite temperature propagators are easily three orders of magnitude harder to calculate. 7 It is also justified a posteriori by our numerical results, since the dominant contribution to the (Poisson-resummed) Matsubara series comes from frequencies in the mid-momentum regime around 1 GeV, where the Green's functions are only very slightly affected by finite temperature. As long as we are only interested in the demonstration of the basic phenomenon without excessive accuracy goals, the use of the T = 0 propagators is thus justified.
We have seen above that the (zero-temperature) propagators in Landau gauge can be re-used in the presence of a constant Abelian background field. All we have to do is to replace the color trace factor (N 2 − 1) in eq. (41) by a sum over all root vectors, and for each root σ shift the momentum argument in the kernels by (This is for G = SU(2); we will discuss G = SU ( Here, q ≡ |k|β/(2π) is the rescaled dimensionless momentum and 7 The reason is that the T = 0 model can use O(4) invariance to reduce the momentum integral to a double integral over the momentum norm and one angle. The same can be done for the spatial momentum integral at finite temperature, which is thus of the same complexity. On top of this, however, the finite temperature solution must sum up to 40 − 50 Matsubara frequencies which represent a coupled channel problem that must be solved by iteration. This easily adds three orders of magnitude to the numerical effort.
is the norm of the shifted 4-momentum. In the first term in eq. (44), we have subtracted the perturbative contribution (ω(k) = k 2 ) for the three transversal gluons, and added it back in the last term, combined with similar contributions from the longitudinal gluon and the perturbative ghost. These contributions in the last term of eq. (44) are known as the Weiss potential [24], and it can be calculated explicitly (see appendix C), Even with the Weiss potential subtracted, eq. (44) is not well suited for numerical evaluation.
Instead, it is prudent to first do some arithmetic manipulations such as Poisson resummation, introduction of spherical 4-coordinates etc. This is presented in detail in appendix C. The result is the equivalent formula where J 1 is a regular Bessel function and ξ = λk = mβk is dimensionless. We will use this formula as the basis for our numerical investigation, with the curvature from eq. (39) and the inverse gluon propagatorω(k) from the renormalized T = 0 gap equation (38), see also Ref. [21] for details.

VI. NUMERICAL SETUP AND RESULTS
We use eq. (47) directly to compute the effective action of the Polyakov loop for G = SU(2).
The integral over the rescaled momentum ξ is tricky because of the oscillating nature of the Bessel function. Since the argument of J 1 (ξ) is temperature independent (which was the purpose of the rescaling q → ξ), we can precompute the roots of J 1 and break up the integral at these roots into many small contributions from each half-oscillation of J 1 . The remainder of the integrand is positive, and so the contributions come with alternating sign and a series accelerator can be used to accurately estimate the remainder. We have checked spot values against Mathematica to ensure that this procedure is reliable. The remaining series over the frequencies m converges quite well due to the factor 1/m 4 . It can also be estimated by an accelerator, or simply summed up to convergence.
For the inverse gluon propagatorω(k), we take the numerical T = 0 solution of the gap equation (38), with the renormalization conditions as discussed in Ref. [21]. Since we no longer discard the derivative term in eq. (32), the gap equation itself differs slightly from the one considered in Ref. [21]. In a first step, we must therefore recompute the T = 0 solution in Landau gauge and determine the renormalization constants which best fit the lattice data. We employ the same high-precision lattice data [25] used for Refs. [6,21] and leave the renormalization conditions unchanged. The solutions for the improved gap equation eq. (38) with the renormalized curvature taken from eq. (39) (and the integral equations for χ(k) and the ghost form factor η(k) as in Ref. [21], see also appendix A) is presented in Fig. 1.
As can be seen, the improved gap equation describes the lattice propagators with the same (or even slightly better) accuracy than in the previous studies. The optimal value for the renormalizaltion scale µ is slightly larger than the µ = 5 GeV used previously; this also affects This group has rank 2 and the effective potential V eff (x, y) is thus a function of two parameters which are the rescaled Cartan components of the background field The potential must be computed as a sum over the SU(2) potentials corresponding to the three non-trivial positive root vectors. The appropriate shift in the momentum p σ is determined in detail in appendix B.
In the left panel of Fig. 3 we have plotted the slice y = a 8 0 = 0 of the effective potential as a function of x ∈ [0, 1] for several temperatures. The potential is no longer symmetric around Weyl alcove center symmetric (x, y) : We now observe a phase transition that is clearly first order, in accordance with lattice findings.
Our best estimate for the phase transition temperature is This should be compared to the lattice estimate [26] of T * ≈ 284 MeV. Again, we observe that the qualitative features of the deconfinement phase transition (such as its order) are correctly predicted, and the numerical estimate of the transition temperature is in reasonable agreement with the lattice data. The accuracy of the agreement is actually better than could be expected, since we have not bothered to re-optimise the renormalization constants of our approach for the case of SU(3), i.e. the entire calculation is still based on SU (2) propagators. This is also one obvious way for numerical improvement, in addition to what was already suggested for To better visualize the Polyakov loop potential V eff (x, y), we have also plotted it as a function of both SU(3) Cartan parameters. Figure 4 shows the result for a single temperature T = 141 MeV below the phase transition, i.e. in the confined phase. Both the surface and density plot clearly indicate that the minima occur at the center symmetric points (51) with L = 0.
By contrast, Fig. 5 shows the same plot for a temperature T = 400 MeV above the phase transition, i.e. in the confined phase. Now the minima are clearly at the symmetry breaking points eq. (52) with |L| = 1.

VII. CONCLUSIONS
In this paper, we have extended the covariant variational approach to Yang-Mills theory to background gauge and studied the effective action for the Polyakov loop. Our findings demonstrate that the formalism is not only well suited to describe low-order Green's functions, but it can also accurately capture the physics of the deconfinement phase transtition.
All known qualitative ascpects of this transition (in particular its order for different color groups) are correctly reproduced, and there is even reasonable numerical agreement with the transition temperatures measured on the lattice. Given the simplicity of our ansatz without adjustable parameters except for the T = 0 renormalization constants, this agreement is quite remarkable. Based on these encouraging findings, it is planned to extend the present formulation to include fermions and study the QCD phase diagram also in the region of non-vanishing chemical potentials. Preliminary studies in this direction are currently underway.
If we emplyo the loop counting scheme explained before eq. (33), this can be solved in the form and we find Next we consider the integral equation for the curvature, Taking the derivative w.r.t.ω −1 (p) under the integral and using eq. (A3) yields, after some algebra, This is precisely the O(λ) contribution if the expansion eq. (A2) is used in the curvature eq. (A4). The term of order O(λ 0 ) is where the numerical prefactor was computed in d = 4 using an O(4)-invariant momentum cutoff. (Other schemes will give different numerical prefactors.) Thus, we have shown This is the momentum space representation of eq. (33) in the main text. Although the derivation was carried out explicitly for Landau gauge, it generalizes to the case of a constant Abelian background field, since the necessary shift in the momentum arguments has no consequence at this point.

Appendix B: Root decomposition of SU(N)
The semi-simple Lie algebra SU(N) has rank r = (N − 1) and there are hence r mutually commuting generators H k which span the Cartan subalgebra of SU(N). As explained in the main text, the background field a µ = a δ µ0 must be chosen in the Cartan subalgebra, Since the H k are anti-hermitean and mutually commuting, they can be simultaneously diagonalized with purely imaginary eigenvalues(−iµ k ). The real numbers µ k are called the weights of H k , and the collection of one eigenvalue from each H k forms a weight vector µ = (µ 1 , . . . , µ r ). The number of such vectors, i.e. the number of the eigenvalues of H k depends on the representation.
In the present paper, we are mainly concerned with the background field in the adjoint representation,â ab = − f abc a c . The weights σ k in the adjoint representation are called the roots, This is the analogon to the cyclic SU(2) basis used in the main text. The complete argument about the diagonalization ford µ and the shift p → p − p σ in the momentum to go from Landau to background gauge remains valid, if only p σ is generalized to For G = SU(2), the root and weight vectors are pure numbers since the Cartan subalgebra is r = 1-dimensional. There are two weights ± 1 2 and three roots {−1, 0, 1}, of which only two are non-vanishing.
This structure easily generalizes to G = SU(3), which has rank r = 2. The two Cartan generators are usually taken as H 1 = T 3 = λ 3 /(2i) and H 2 = T 8 = λ 8 /(2i) in terms of Gell-Mann matrices. The root and weight vectors are both r = 2-dimensional. In fact, the weights read µ : 0 , To evaluate eq. (C1), we start by Poisson resumming the Matsubara series Next we shift the z-integral by introducing s ≡ (z + x) and combine the terms in the series with index m and (−m). After interchanging the order of summation and integration, ds e 2πims f 2πβ −1 s 2 + q 2 .
The q-dependent terms in the integrand are even in q so that we can extend the q-integral to all of R. For the resulting double integral in (q, s) over R 2 , we use polar coordinates (r, ϕ), dr r 2π 0 dϕ r 2 sin 2 ϕ · e 2πimr cos ϕ · f (2πr/β) .
The ϕ-integral leads to 2π 0 dϕ sin 2 ϕ · e 2πimr cos ϕ = J 1 (2πmr) mr , where J 1 is a regular Bessel function. Finally, we rescale the integration variable in the remaining integral r → ξ ≡ 2πmr and collect all pieces. The result is For later convenience, we write this formula as Note that the temperature only enters through the factor h(λ = βm), and the series now converges much more quickly due to the 1/m 4 term. (For implementation details, see section
We use the proper-time representation of the logarithm, where Λ → ∞ is a cutoff to be lifted at the end of the calculation, and γ E is Euler's constant.
The second term is λ-independent and only leads to a temperature-independent constant in the potential, which can be dropped. For the first term, we can do the ξ-integral to obtain h(λ) = λ 4 16 ∞ Λ −2 ds s 3 e −λ 2 /(4s)−sµ 2 .
After changing variables s → t = λ 2 4s , the cutoff can be lifted, where K 2 is an irregular modified Bessel function. We can put this in the main formula (C3) and find The rhs depends on the temperature only in the combination (µβ), i.e. increasing the mass is equivalent to decreasing the temperature, and vice versa. The Weiss potential follows by including a prefactor 9 2/3 and sending the mass µ → 0. From (βµ) 2 K 2 (βµ m) → 2/m 2 in this limit, we obtain eventually the SU(2) Weiss potential [24] β 4 W(x) = 4