Relation between the continuum threshold and the Polyakov loop with the QCD deconfinement transition

Using vector and axial-vector correlators within finite energy sum rules with inputs from a chiral quark model, coupled to the Polyakov loop, with nonlocal vector interactions, we extend our previous work to confirm the equivalence between the continuum threshold $s_0$ and the trace of the Polyakov loop $\Phi$ as order parameters for the deconfinement transition at finite temperature $T$ and quark chemical potential $\mu$. The obtained results are in agreement with our initial conclusion, where we showed that $s_0(T,\mu)$ and $\Phi(T,\mu)$ provide the same information for the QCD deconfinement transition.


I. INTRODUCTION
Since quark color charge in a medium is screened due to density and temperature effects, Quantum Chromodynamics (QCD) predicts that at very high temperatures (T Λ QCD ) and low baryon densities, matter appears in the form of a plasma of quarks and gluons [1].
If one of those variables increases beyond a certain critical value, the interactions between quarks no longer confine them inside hadrons. This is usually referred to as the deconfinement phase transition. Simultaneously [2], at small densities another transition takes place, the chiral restoration. For high values of chemical potential, these two transitions can arise at different critical temperatures. The result will be a quarkyonic phase, where the chiral symmetry is restored but the quarks and gluons remain confined.
In the confined region, QCD is strongly coupled and coupling-constant expansions become inapplicable. At finite density, lattice QCD (lQCD) methods based on large-scale Monte Carlo simulations are also not applicable because lQCD has the sign problem [3,4].
Therefore, predictions of the transition features at regions that are not accessible through lattice techniques arise from effective theories, like for instance, the nonlocal Polyakov−Nambu−Jona-Lasinio (nlPNJL) models (see Ref. [5] and references therein), where quarks interact through covariant nonlocal chirally symmetric couplings in a background color field. These approaches, which can be considered as an improvement over the local model [6][7][8][9][10][11][12], offer a common framework to study both chiral restoration and deconfinement transitions for hadronic systems at finite temperature and/or chemical potential (see e.g. Refs. [5,[13][14][15][16][17][18][19]). In fact, the nonlocal character of the interactions arises naturally in the con-text of several successful approaches to low-energy quark dynamics [20][21][22], and leads to a momentum dependence in the quark propagator that can be made consistent with lQCD results.
In this work we use a nlPNJL model with vector and axial-vector interactions [23]. In addition to the standard scalar and pseudoscalar quark-antiquark currents, we consider couplings between vector and axial-vector nonlocal currents, satisfying proper QCD symmetry requirements.
In order to study the properties of chiral and deconfinement phase transitions it has been customary to study the behavior of corresponding order parameters as functions of the temperature and quark chemical potential, namely the quark anti-quark chiral condensate qq and the trace of the Polyakov loop (PL) Φ, respectively.
In addition to Φ, another phenomenological QCD deconfinement parameter that has been introduced in the literature [24] is the continuum threshold s 0 , for the onset of perturbative QCD (PQCD) in hadronic spectral functions. Around this energy, and at zero temperature, the resonance peaks in the spectrum are either no longer present or become very broad.
The natural framework to determine s 0 has been that of QCD sum rules [25]. This quantum field theory framework is based on the operator product expansion (OPE) of current correlators at short distances, extended beyond perturbation theory, and on Cauchy's theorem in the complex s-plane.
In this article we reconsider the light-quark axialvector channel with an improved hadronic spectral function involving the a 1 (1260) resonance [26], in addition to the already considered pion pole approximation [27]. We also include the vector channel with a ρ-meson saturated spectral function.
Within this theoretical framework, using finite energy sum rules with inputs for the spectral functions obtained from a nonlocal quark model (masses, decay constants and widths), we compare the thermal behavior of s 0 and Φ in both channels, at zero and finite chemical potential.
The paper is organized as follows. In the next section, we briefly review the finite energy sum rules (FESR) program for the light-quark vector and axial-vector channel. In Sect. III we present the general formalism for a finite temperature and density system within the nlPNJL effective model. The numerical and phenomenological analyses at zero and finite density are included in Sect. IV. Finally, in Sect. V we summarize our results and present the conclusions.

II. FINITE ENERGY SUM RULES
Within the formalism of finite energy sum rules [25] we extend our previous work [27], reconsidering the lightquark axial-vector channel with an improved hadronic spectral function involving the a 1 (1260) resonance, in addition to the already considered pion pole approximation. Moreover, we also include in this analysis the vector channel with a ρ-meson saturated spectral function We begin considering the current-current correlation function where J(x) is a local quark current. Let us consider, for our study, the correlator of light-quark vector and axialvector currents :] is the electric charge neutral conserved vector current in the chiral limit, and q µ = (ω, q) is the four-momentum carried by the current. The transverse parts Π A,V 0 are related to the vector meson resonances, whereas the longitudinal axial contribution Π A 1 correspond to the pion pole. Since we are working with non-strange current correlators, there is no longitudinal term in the vector channel [28].
In addition, through Operator Product Expansion of current correlators at short distances [29,30], one has where Ô 2N (µ 2 ) ≡ 0|Ô 2N (µ 2 )|0 , µ 2 is a renormalization scale, the Wilson coefficients C N depend on the Lorentz indexes and quantum numbers of the currents, and on the local gauge invariant operatorsÔ N built from the quark and gluon fields in the QCD Lagrangian. These operators are ordered by increasing dimensionality and the Wilson coefficients are calculable in PQCD. The unit operatorÎ has dimension d ≡ 2N = 0, and C 0 stands for the purely perturbative contribution. On the other hand, using Cauchy's theorem in the complex squared energy s-plane, we obtain the quark-hadron duality where the radius of the circle s 0 is large enough for QCD and the OPE to be used on the circle. Using the OPE, Eq. (3), the finite energy sum rules at finite temperature T and chemical potential µ become [24,31,32] ( For N = 1, the dimension d = 2 term in the OPE does not involve any condensate, as it is not possible to construct a gauge invariant operator of such dimension from the quark and gluon fields. Moreover, there is no evidence for such a term for low values of T [33][34][35].
The dimension d = 4 term, a renormalization group invariant quantity, is given in the chiral limit, for the vector and axial sector by As was mentioned previously, there is a difference in the currents involved in vector and axial processes, i.e. the former involves electrically neutral currents, while the latter involves electrically charged currents. Then, the normalization of these vector current correlators differs by a factor two.
The leading power correction of dimension d = 6 is the four-quark condensate, which in the vacuum saturation approximation [30] becomes which is channel dependent and has a very mild dependence on the renormalization scale. This approximation has no solid theoretical justification, other than its simplicity.
In the static limit ( q → 0), and for finite T and µ, the spectral function in PQCD, Π A,V (ω 2 , T, µ)| PQCD (to simplify the notation we shall omit the s, T and µ dependence), is given by where Li 2 (x) is the dilogarithm function, s = ω 2 , and n ± (x) = 1 e (x∓µ)/T + 1 (9) are the Fermi-Dirac thermal distributions for particles and antiparticles, respectively. In the chiral limit, the axial-vector hadronic sector in the spectral function can be approximated by the pion pole, followed by the a 1 (1260) resonance where f π is the pion decay constant, and a fit in the resonance region (s < 2 GeV 2 ) to the ALEPH data quoted in Ref. [36], gives with α = 0.1136 ± 0.0002, β = 0.302 ± 0.001 and δ = 0.8501 ± 0.0006. The mass, width and decay constants of the a 1 meson, m a1 , Γ a1 and f a1 , respectively, are obtained from a two flavor nonlocal PNJL model with vector interactions [23], that will be described in Sect. III. This fit together with the ALEPH data are shown in Fig. 1.
In our previous work [27] we only considered for the axial spectral function the pion pole. Therefore, incorporating the a 1 resonance improves the approximation and constitutes the best possible approach within the present theoretical advances.
For the vector sector, as usual [37][38][39], we will assume ρ-meson saturation of the spectral function in terms of a Breit-Wigner resonance. This parametrization has been normalized such that its area is equal to the area under a zero-width expression [37]. In addition, we will consider a contribution due to the coupling of the vector current to two pions in the thermal bath, the so called scattering term, and is given by where n B (z) = 1/(e z − 1) is the Bose thermal function. Therefore, we have where as before, m ρ , Γ ρ and f ρ are the ρ mass, width and decay constants, respectively, also calculated within the nlPNJL model.

III. THERMODYNAMICS AT FINITE DENSITY
In this section we present the formalism of a twoflavor quark model coupled to the Polyakov loop, that includes nonlocal vector and axial-vector quark-antiquark currents, in addition to the standard scalar and pseudoscalar nonlocal interactions [23]. The corresponding Euclidean effective action is given by [40] where ψ is the N f = 2 fermion doublet ψ ≡ (u, d) T , andm = diag(m u , m d ) is the current quark mass matrix. In what follows we consider isospin symmetry, , τ a are the Pauli matrices and the functions G(z) and F(z) in Eq. (15) are nonlocal covariant form factors characterizing the corresponding interactions. The scalar-isoscalar component of the nonlocal currents will generate a momentum dependent quark mass in the quark propagator, while the "momentum" current j M (x) leads to a momentum-dependent wave function renormalization of the quark propagator, in consistency with lQCD analyses.
To work with mesonic degrees of freedom, the fermionic theory is bosonized in a standard way [41] by considering the corresponding partition function Z = DψDψ exp[−S E ], introducing auxiliary bosonic fields σ 1 (x), σ 2 (x) (scalar, related to j S (x) and j M (x)), π a (x) (pseudoscalar), v a µ (x) (vector) and a a µ (x) (axial vector) and integrating out the quark fields. Details of this procedure can be found in Ref. [40].
Since we are interested in deconfinement and chiral restoration transitions, we extend the bosonized effective action to finite temperature T and chemical potential µ. This will be done using the standard imaginary time formalism. Concerning the gauge fields A µ , we assume that quarks move in a constant background field φ = A 4 = iA 0 = ig δ µ0 G µ a λ a /2. We will work in the so-called Polyakov gauge, in which the matrix φ is given a diagonal Although at finite chemical potential, the traced PL Φ and its conjugate Φ * could in principle have different values, the difference between the Φ and Φ * emerges beyond the mean field approximation [10,42,43]. Therefore we choose for our model to adopt the usual prescription Φ = Φ * (see Ref. [23]) implying that φ 8 = 0, leaving only φ 3 as an independent variable. Consequently, the traced For the light quark sector the trace of the Polyakov loop turns out to be an approximate order parameter, in the same way the chiral quark condensate is an approximate order parameter for the chiral symmetry restoration.
Thus, in the mean field approximation (MFA), and following the same prescriptions as in Refs. [5,16,19], the thermodynamic potential Ω MFA at finite temperature T and chemical potential µ is given by [44] where whith the shorthand notation The constantsσ 1,2 andω are the mean field values of the scalar fields and the isospin zero vector field. At nonzero quark densities, the flavor singlet term of the vector interaction develops a nonzero expectation valuē ω, while all other components of the vector and axial vector interactions have vanishing mean fields [45]. The mean field values can be calculated by minimizing Ω MFA , while the functions m(p) and z(p) -momentumdependent effective mass and wave function renormalization (WFR)-are related to the nonlocal form factors and the vacuum expectation values of the scalar fields by where g(p) and f (p) are the Fourier transforms of the form factors in Eq. (15). We have also defined the sums over color indices run over c = r, g, b, with the color background fields components being φ r = −φ g = φ, φ b = 0, q = q 2 + m 2 , and ω n = (2n + 1)πT are the fermionic Matsubara frequencies. The vector coupling generates a shifting in the chemical potential as [44] whereq The term Ω reg is the regularized expression with the thermodynamic potential of a free fermion gas Ω free , and finally the last term in Eq. (16) is just a constant fixed by the condition that Ω MFA vanishes at T = µ = 0.
The effective gauge field self-interactions are given by the Polyakov loop potential U(Φ, T ). At finite temperature, it is usual to take for this potential a functional form based on properties of pure gauge QCD. Among the most used effective potentials, the Ansatz that provides the best agreement with lQCD results [5,18] is the polynomial function based on a Ginzburg-Landau Ansatz [9,46]: Numerical values for the parameters can be found in Table I   In absence of dynamical quarks, from lattice calculations one expects a deconfinement temperature T 0 = 270 MeV. However, it has been argued that in the presence of light dynamical quarks this temperature scale, a further parameter of the model, should be adequately reduced to about 210 and 190 MeV for the case of two and three flavors, respectively, with an uncertainty of about 30 MeV [47].
In order to fully specify the model under consideration, we proceed to fix the model parameters as well as the nonlocal form factors g(q) and f (q). We consider here Gaussian functions which guarantee a fast ultraviolet convergence of the loop integrals and offer a momentum dependence in good agreement with lQCD, and other form factors [18,40].
Once the mean field values are obtained, the behavior of other relevant quantities as functions of the temperature and chemical potential can be determined. We concentrate, in particular, on the chiral quark condensate qq = ∂Ω MFA reg /∂m and the traced PL Φ, which will be taken as order parameters for the chiral restoration and deconfinement transitions, respectively. The associated susceptibilities will be defined as χ ch = ∂ qq /∂T and χ Φ = dΦ/dT .

Meson masses
Meson masses can be obtained from the terms in the Euclidean action that are quadratic in the bosonic fields. From the zero temperature action [40] one can obtain, using the imaginary time formalism, the finite temperature action, with p m ≡ (ν m , p). The functions G ab (ν m , p) are given by finite temperature one-loop integrals arising from the fermionic determinant in the bosonized action, and ν m = 2mπT are the bosonic Matsubara frequencies. Once cross terms have been eliminated, the resulting functions G M (ν m , p) stand for the inverses of the effective thermal meson propagators, in the imaginary time formalism. The functions G ρ,a1 (ν m , p) correspond to the transverse projections of the vector and axial-vector fields. Thus the masses of the physical ρ 0 , ρ ± (which are degenerate in the isospin limit) and a 1 (the transverse parts of the a µ fields do not mix with the pions) can be obtained by solving the equation G ρ,a1 (0, −im ρ,a1 ) = 0, where with D(q nc ) = q 2 nc + m 2 (q nc ) and q ± nc = q nc ± p m /2. In the case of the pseudoscalar sector, from Eq. (25) it is seen that there is a mixing between the pion fields and the longitudinal part of the axial vector fields [48,49]. Therefore, after the elimination of the mixing term, the pion mass can be then calculated from Gπ(0, −im π ) = 0, where with

Meson decay constants
The pion weak decay constant f π is given by the matrix elements of axial currents between the vacuum and the physical one-pion states at the pion pole, On the other hand, the matrix elements of the electromagnetic current J em between the neutral vector meson state and the vacuum determine the vector decay constant f ρ [40], with p 2 = −m 2 ρ , where e is the electron charge. And finally, the axial-vector decay constant f a1 is defined by the matrix elements of the electroweak charged currents J ew between the axial-vector meson state and the vacuum, at p 2 = −m 2 a1 , as [23] 0|J ew µ (0)|a 1ν (p) = f a1 (p 2 )(g µν p 2 − p µ p ν ) .
In order to obtain these matrix elements within our model, we have to gauge the effective action through the introduction of gauge fields, and then we have to calculate the functional derivatives of the bosonized action with respect to the currents and the renormalized meson fields. In addition, due to the nonlocality of the interaction, the gauging procedure requires the introduction of gauge fields not only through the usual covariant derivative in the Euclidean action, but also through a transport function that comes with the fermion fields in the nonlocal currents (see e.g. Refs. [41,50,51]).
After a lengthy calculation, the decay constants at finite T and µ, are given by where the mesons wave function renormalization can be obtained from The analytical expressions for F 0 , F 1 , J (I,II) V,A and the thermal behavior of f π , f v and f a can be found in Ref. [23].
In absence of vector meson fields, the mixing term F 1 in f π vanishes, and the expression reduces to that quoted in Ref. [52].
The resulting one-loop contributions are diagrammatically schematized in Fig. 2. Tadpole-like diagrams, which are not present in the local PNJL model, arise from the occurrence of gauge fields.

Decay widths
In general, various transition amplitudes can be calculated by expanding the bosonized action to higher orders in meson fluctuations. The decay amplitudes A ρ (ρ → ππ) and A a1 (a 1 → ρπ) are obtained by calculating the corresponding functional derivatives of the effective action.
For the vector sector, only the transverse piece contributes to ρ → ππ decay [40], while for the axial-vector sector both transverse and longitudinal parts contribute to the a 1 → ρπ decay [23].
In order to study the thermal dependence of these decay widths, it is necessary to modify the two-body phase space to include finite temperature effects.
Following Refs. [53,54], the decay of a particle at rest of mass M , into particles of masses m 1 and m 2 in equilibrium with the heat bath, is given by where A M is evaluated within the effective model. For the decays ρ → ππ (M = m ρ , m 1 = m 2 = m π ) and a 1 → ρπ (M = m a1 , m 1 = m ρ , m 2 = m π ) we have The factors g ρππ , g aρπ and f aρπ are one-loop functions that arise from the expansion up to order three of the effective action. And, due to the π − a 1 mixing, receive contributions from the diagrams sketched in Fig. 3. For the explicit forms of these functions and the temperature dependence of the widths, we refer the reader to Refs. [23,40]. ρ π π ρ a 1 π ρ a 1 a 1 Figure 3: Diagrams contributing to ρ and a 1 decays due to the π − a 1 mixing.

IV. RESULTS
To determine the relation between the perturbative QCD threshold s 0 and the trace of the Polyakov loop Φ, we study the finite energy sum rules first at zero density, where chiral restoration and deconfinement occurs simultaneously. Then we move into finite values of the chemical potential, where chiral symmetry is restored through a first order phase transition at a certain critical temperature whereas deconfinement has not yet been achieved, occurring the deconfinement transition at a higher temperature.
Therefore, for large chemical potentials, the critical temperatures for the restoration of the chiral symmetry and deconfinement transition begin to separate. The region between them denotes a phase where the chiral symmetry is restored but quarks remains confined. This splitting is strongly dependent on the functional form of the Polyakov loop effective potential and also on the parameter T 0 entering in the PL potential [5]. If we consider for this parameter an explicit dependence with µ [47,55,56], both transitions are always simultaneous, and therefore there is no such mixed phase.
It should be mentioned that the FESR program was performed at one-loop order. In the thermal perturbative QCD sector, only the leading one-loop contributions can be taken into account, since the problem of the appearance of two scales, i.e. the short-distance QCD scale and the critical temperature, remains unsolved. However, at zero temperature it is possible to extend the calculations up to five-loop order. The estimations at one-loop and five-loop order for s 0 , C 4 Ô 4 and C 6 Ô 6 differ considerably [26,37,57] (even with changes in the sign of the coefficients). Therefore, only the thermal behavior of the condensates obtained through the higher order FESR should be taken into account. The nlPNJL effective model includes six free parameters, namely the current quark mass m and the coupling constants G S , G V , G 0 , G 5 and κ p (see Eq. (14)). In addition, one has to determine the cutoffs Λ 0 and Λ 1 introduced in the form factors, Eq. (24).
Through a fit to lQCD results quoted in Ref. [58] for the functions m(p) and z(p), we obtain Λ 0 = 1092 ± 22 MeV and Λ 1 = 1173 ± 60 MeV. Furthermore, by requiring that the model reproduces the value of z(p = 0) and the empirical values of three physical quantities, namely the masses of mesons π and ρ and the pion weak decay constant f π , one can determine the model parameters quoted in Table II  Regarding the vector coupling constants G 0 and G 5 , in the MFA only the former should be fixed (see Eq. (17)). Therefore, we will follow the prescription used in Ref [44], parameterizing the isoscalar vector coupling as G 0 = η G V . Hence, the strength of the vector coupling can be evaluated by considering different values for η.
The influence of the vector coupling increases with the chemical potential. At zero density,ω vanishes for all temperatures and therefore the vector interactions do not contribute to the mean field thermodynamic potential.
Once the model parametrization is defined, one can calculate several meson properties at finite temperature and/or chemical potential. For the numerical results of meson masses, decay constants, decay widths and other observables we refer the reader to Ref. [23].
The inputs used in the FESR calculated within the nlPNJL model are the masses, decay constants and decay widths of the π, ρ and a 1 mesons.

FESR program at zero density
At zero chemical potential, for the above set of parameters and for the polynomial Polyakov loop potential, Eq. (22), with T 0 = 210 MeV, we obtain through the corresponding susceptibilities, almost the same chiral and deconfinement critical temperatures T c = 202 MeV (less than 3% of difference), as expected, since at µ = 0 chiral restoration and deconfinement takes place simultaneously as crossover phase transitions. This behavior was verified by lQCD calculations [2], in nlPNJL models [5,17,18] and also obtained by finite energy sum rules [32].
Moreover, in Ref [2], the deconfinement temperature defined at the peak of the entropy of a static quark (which is related to the Polyakov loop) is located at the same temperature, within errors, as the chiral susceptibility even at finite lattice spacing. In several works (see [2] and references therein), the deconfinement transition in lQCD with light dynamical quarks has been studied in terms of the inflection point of the renormalized Polyakov loop and fluctuations of conserved charges. Usually, these critical deconfinement temperatures are equal or larger than the restoring chiral transition critical temperature. In addition, these approaches have the disadvantage of being lattice scheme dependent and therefore the obtained values may differ considerably between them.
Furthermore, as it is discussed in Refs. [59][60][61], the strict comparison between our results and lattice data for the traced Polyakov loop has to be taken with some care, owing to the difference between the definitions of Φ in the continuum and on the lattice. From Eq. (5), for N = 1 and µ = 0, we obtain the first FESR as function of the temperature T , where n F (x) = 1/(1 + e x ) is the Fermi thermal function, and the spectral functions Π A,V | HAD are given by Eqs. (10) and (13). The continuum threshold s A,V 0 can be calculated, as function of the temperature, by solving this equation with the corresponding spectral function.
In Fig. 4 we plot, as function of the reduced temperature T /T c , the continuum threshold for the vector (axial) channel in solid (dashed) line, together with the trace of the PL and the quark condensate normalized by its value at T = 0 in dotted and dash-dotted line, respectively.
In addition, in dot-dashed line, we quote the continuum threshold for the axial-vector channel in the pion pole approximation s A,π 0 , its thermal behavior is equivalent to that found in our previous work, Ref. [27].
The FESR have solutions up to T ∼ 0.9 T c , a temperature at which s A,V 0 reaches its minimum. A short extrapolation, denoted in the figure by a thin dotted line, should be understood for all results in the sequel.
As we expected for both channels, the PQCD threshold vanishes at critical temperatures T V c = 202 MeV and T A c = 208 MeV, located almost at the chiral critical temperature T ch c = 202 MeV and the PL deconfinement temperature T Φ c = 196 MeV. From the figure one can see that the thermal behavior of s A 0 (T ) and s V 0 (T ) close to T c is similar, even when the hadronic spectral functions are very different in these two cases. This result is pointing to an approximate universality of the deconfinement transition in light-quark systems.
The higher order FESR, from where it is possible to estimate the thermal dependence of the gluon and four-quark condensate, can be analogously obtained from Eq. (5) with N = 2 and N = 3. Both condensates show the expected behavior with a finite value at zero temperature, and decreasing monotonically as function of the temperature.
In order to avoid the mentioned discrepancies at T = 0 due to different loop-order calculations, we plot in Fig. 5 for the axial and vector channel in dashed and solid line, respectively, a normalized C 4 Ô 4 as function of the reduced temperature, where we have defined such quantity as It can be seen from the figure that the thermal evolution of the gluon condensate is quite similar in both channels (as it is expected), even when the spectral functions, Eqs. (10) and (13), are completely different.
The thermal behavior of the meson masses directly affects the decay widths, since the kinematic condition in Eq. (34) tends to zero as T increases and therefore, even when the phase space is increased due to the Bose enhancement, the width decreases [23].
For the process ρ → ππ, the decay width starts to drop above the chiral critical temperature, since beyond this temperature the π mass grows faster than the ρ mass.
For the a 1 decay, the width begins to diminish before the chiral critical temperature, vanishing close to T c . This is caused by the chiral partner mass degeneration. Near above T c , vector mesons have approximately the same mass.
The physical decay process is a 1 → πππ, however in our approach we are only considering the main channel of the partial decay a 1 → ρπ. Other partial widths contribute approximately with 40% of the total width [54]. These processes, not considered here, will contribute to the total width and could modify the decreasing behavior of the a 1 width.
The consequence of this temperature dependence in Γ a1 is an small increment in the value of s A 0 (see Fig. 4). This indicates that the non considered decay processes for the a 1 could be relevant for temperatures close to the critical temperature T c .
Nevertheless, it should be noticed that although we are just considering the a 1 main decay channel, the general thermal dependence of the continuum threshold s A 0 (T ) not only is in agreement with other QCD sum rules results at one-loop order [57], but also with the behavior found in this work for the vector channel.

FESR program at finite density
In general, one can find regions in the QCD phase diagrams where chiral symmetry is either broken or restored through a crossover or a first order phase transition and regions where the system remains either in confined or deconfined states.
For relatively low densities, chiral restoration takes place as a smooth crossover, whereas for high values of chemical potential the order parameter has a discontinuity at a given critical temperature T µ c signaling a first order phase transition. This gap in the quark condensate induces also a jump in the trace of the PL, and the PL susceptibility present a divergent behavior at the chiral critical temperature. Therefore, as in Ref [17], when the chiral phase transition is first order, we define the deconfinement critical temperature T µ Φ requiring that Φ = 0.4, which could be taken as large enough to denote deconfinement.
The value of Φ at both sides of the discontinuity indicates, for a chiral symmetric state, whether the system remains confined or not. The region where the chiral symmetry is restored but the quarks and gluons remain confined, is usually referred as the quarkyonic phase [62][63][64].
If we move, in the T − µ plane, along the first order phase transition curve, the critical temperature rises from zero up to a critical endpoint (CEP) temperature T CEP , while the critical chemical potential decreases from its value at zero temperature µ c to a critical endpoint chemical potential µ CEP . Beyond this point, the chiral restoration phase transition proceeds as a crossover.
To determine these temperatures and densities we need to fix the value of the coefficient η in the definition of G 0 . Here, as in Ref. [23], we choose η = 0 and η = 0.3 as representative cases. Since the former leads to a mean field theory without vector interactions, and the later provides the best phenomenological agreement with other effective models (see Ref. [44] and references therein).
In Fig. 6 we plot the phase diagram for a mean field theory with and without vector interactions. Specifically, in the upper and lower panel we quote the reduced critical temperatures as function of the reduced chemical potential for η = 0 and η = 0.3, respectively.
In the crossover region the deconfinement temperatures are determined with the PL susceptibility, whereas for µ > µ CEP are obtained by requiring that 0.4 < Φ < 0.6, quoted in dotted line and with the color shaded area, respectively. In addition, the chiral critical temperatures are represented by the solid and dashed line, for the first order and crossover phase transition curves, respectively. Finally, the dot indicates the position of the critical endpoint. These phase diagrams were taken from Ref. [23].
At zero density the FESR have solutions up to 0.9 T c where the continuum threshold reaches its minimum. Therefore, to define the value of T where s A,V 0 = 0, it is usual to perform a short extrapolation.
However, in the crossover region, when µ increases the FESR equations stop having solutions at lower temperatures and consequently the extrapolation is no longer well defined. This led us to set a range of critical temperatures T µ s A,V 0 , for each value of µ, that enlarges as the density grows.
We quote this, in both panels of Fig. 6, with pattern filled areas. The critical temperatures within those bands, correspond to a family of extrapolations of the obtained results, for both channels, with similar χ 2 .
Despite this broadening, it can be seen that the density dependence of T µ s A,V 0 is equivalent to that found for T µ Φ and T µ c , showing that the continuum threshold and the traced Polyakov loop provide analogous information about the QCD deconfinement transition in a wide range scenario.
On the other hand, for chemical potentials larger than the µ CEP , the thermal dependence of the continuum threshold in the available temperature range, does not allow to extrapolate the s 0 curve beyond the last obtained value, usually located at T T µ c (see Fig. 7). In a finite density scenario, when µ > µ CEP , the two transitions takes places separately at different critical temperatures, and therefore provides unique conditions to identify the phenomenological equivalence between s 0 and Φ as deconfinement order parameters.
In Table III we summarize, for these two scenarios, the critical temperatures, critical chemical potentials and the CEP coordinates.   In the upper and lower panel of Fig. 7 we plot for η = 0 and η = 0.3, respectively, the continuum threshold for the vector (axial) channel in solid (dashed) line, together with the trace of the PL and the normalized quark condensate in dotted and dash-dotted line, respectively.
As before, we also quote in dot-dashed line, the continuum threshold for the axial-vector channel in the pion pole approximation, s A,π 0 . In both situations, η = 0 and η = 0.3, we see that for bigger densities than the critical end point chemical potential, the thermal equation has not solution beyond the critical temperature T µ c . The continuum threshold, for both channels, stops with a finite value at this temperature, signaling that the system continues in a confined state.
As in the zero density case s A 0 has a slight increase before decreasing, consequence of the approximation used for the a 1 width decay.
Regarding the Polyakov loop, the value of Φ at both sides of T /T µ c = 1 indicates that the system, at this temperature, remains in a confined state, even when the chiral symmetry has been restored.
In this way, we see that the Polyakov loop and the continuum threshold provide the same information. When the chiral symmetry is restored, s 0 and Φ show that we are still in a confined phase. This characterizes the occurrence of a quarkyonic phase.

V. SUMMARY AND CONCLUSIONS
Along this article, as in our previous work [27], we compare the behavior of two vastly used phenomenological order parameters for the deconfinement transition: the continuum threshold s 0 and the trace of the Polyakov loop Φ.
In Ref. [27] we study the finite energy sum rules for the axial-vector current correlator saturating the spectral function with the pion pole approximation. Here, we have extended that analysis in two complementary directions: we improve the approximation for the axial spectral function including the a 1 resonance, and we consider the vector current correlator assuming ρ-meson saturation for the spectral function.
In this way, both on the side of the FESR formalism and the nlPNJL model we have considered the best possible phenomenological approach. Since there are no further possible corrections to the spectral function or to the effective model in the light quark sector, our results seem to be strong and conclusive. The input parameters used in the FESR, namely the masses, decay constants and decay widths for the π, ρ and a 1 mesons, were obtained from a SU(2) PNJL model with nonlocal vector and axial-vector interactions [23].
At zero density, we determine that the continuum threshold vanishes, for both channels, at approximate the same temperature where the Polyakov susceptibility has its maximum value.
At finite density, beyond the critical end point chemical potential, we found that for both deconfinement parameters, the system remains in its confined phase even when the chiral symmetry is restored.
Therefore, based in our previous study and the results obtained here, we can conclude that both quantities, s 0 (T, µ) and Φ(T, µ), provide the same kind of physical information about the QCD deconfinement transition.