Thermodynamics of QCD at vanishing density

We study the phase structure of QCD at finite temperature within a Polyakov-loop enhanced quark-meson model. Such a model describes the chiral as well as the confinement-deconfinement dynamics. In the present investigation, based on the approach and results put forward in [1-4], both, matter as well as glue fluctuations are included. We present results for the order parameters as well as some thermodynamic observables and find very good agreement with recent results from lattice QCD.


I. INTRODUCTION
The study of strongly interacting matter under extreme conditions is a very active field of research. Experiments conducted at CERN, RHIC and the future FAIR and NICA facilities aim at probing the phase structure of Quantum Chromodynamics (QCD).
From the theoretical side, calculating the phase structure from first principles is a hard task which requires the use of non-perturbative methods. Over the recent years a lot of progress has been made in this direction. In particular it has been shown that, apart from lattice QCD, also continuum methods, such as the Functional Renormalisation Group (FRG) [5][6][7][8][9][10][11][12] are well suited to study the QCD phase diagram. This has been demonstrated in e.g. [1,3,[13][14][15][16] at vanishing as well as finite temperature and chemical potential. Complementary to first-principles studies, low-energy QCD has been studied successfully within effective models. Especially the use of Polyakov-loop extended chiral models makes it possible to study the interrelation of the chiral and deconfinement phase transitions, e.g. [2,3,. However, the confinement sector in these models is not fully constrained, resulting in various parametrisations of the corresponding order-parameter potential, the glue or Polyakov-loop potential. Furthermore, the important unquenching effects on the glue potential are usually not included. Ideally, this potential is derived from QCD directly, leaving no ambiguity. This has recently been accomplished with the FRG for two-flavour QCD in the chiral limit [1] and for 2 + 1 flavours in [16] and puts us in the position to make use of these results to improve the effective description of the gauge sector. In summary, these effective models can be systematically improved towards full QCD, using input from the lattice and other first-principles studies, see e.g. [1][2][3]10]. In [3,40] this approach has already been tested in a mean-field approximation.
In the present work we aim at quantitative results for the thermodynamics of QCD. To achieve this goal, we combine the previous efforts of [2,4] and include quantum and thermal fluctuations with the FRG in an effective Polyakov-quark-meson (PQM) model with 2+1 flavours. Furthermore, we apply the augmentation of the gauge sector by QCD results as in [3]. In combination, this gives us a good handle on the chiral and confinementdeconfinement transitions and thermodynamics of QCD. This work is structured as follows: In Sec. II we briefly review the FRG approach to QCD and its connection to low-energy effective models. In particular, we discuss how to augment low-energy effective models with firstprinciples results from QCD. In Sec. III we provide the details of our truncation and present the resulting flow equation in Sec. III B. Results for the order parameters and thermodynamic observables for 2+1 flavours are presented in Sec. IV A and Sec. IV C, respectively. Sec. IV B contains our prediction for the thermodynamics in the two-flavour case. Concluding remarks and a summary are presented in Sec. V. We discuss the dependence of our results on various parameters in the appendix.

II. FUNCTIONAL RENORMALISATION GROUP APPROACH TO LOW-ENERGY QCD
The mapping of QCD degrees of freedom to low-energy effective models is discussed in depth in, e.g. [1][2][3]10]. Here, we only briefly recapitulate the main points. Fig. 1 shows the pictorial representation of the FRG flow of QCD, where the first two loops represent the gluon and ghost contributions, respectively, whereas the third loop denotes the quark degrees of freedom. The fourth loop corresponds to mesonic degrees of freedom which have been introduced via the dynamical hadronisation technique [7,[41][42][43].
It is well-established that the ghost-gluon sector decouples from the matter dynamics below the chiral and Partially hadronised version of the FRG flow for QCD. The loops denote the gluon, ghost, quark and meson contributions, respectively. The crosses mark the FRG regulator term. deconfinement temperatures, see e.g. [44]. In terms of the flow equation, Fig. 1, this means that in this regime we are only left with the dynamical matter sector given by the last two loops, explicitly The full field content is collected in φ = (a, c,c, q,q, H) . In the non-perturbative domain of QCD the spectrum is gapped and only light constituent quarks (q,q) and the corresponding hadrons (H) do not decouple, whereas the ghost (c,c) and gluon (A =Ā + a) fields, as well as the heavy matter sector act as spectators at low densities. Here we have decomposed the gauge fields into a constant backgroundĀ and a fluctuating part a.
The effective action of full QCD can then be written as where V is the spatial volume and β = 1/T the inverse temperature. In Eq. (2), the first term denotes the QCD glue potential, encoding the ghost-gluon dynamics in the presence of matter fields. The second term contains the matter contribution coupled to a background gluon fieldĀ 0 . This contribution is well-described in terms of low-energy chiral models, such as the Nambu-Jona-Lasinio (NJL) and quark-meson (QM) models, coupled to Polyakov loops. In this work we make use of a Polyakovquark-meson (PQM) truncation [25,29,31] for the matter sector at low energies. It is important to notice that the glue potential V [A 0 ] in full QCD is different from its Yang-Mills counterpart due to unquenching effects, see e.g. the discussion in [2,3]. The glue potential used in effective models, on the other hand, is usually fitted to pure Yang-Mills lattice results [18,21,23,27,45,46]. To approximate unquenching effects we formulate the glue potential in terms of the reduced temperature and To be more precise, there are also two reduced temperatures, defined in terms of the critical temperatures T cr = T YM cr and T cr = T glue cr . One important effect of dynamical matter fields is to lower the scale T glue cr as compared to T YM cr , which can be used to model the unquenching effects [25,29]. In the present work we remedy the scale mismatch with the help of first-principles QCD results, see [3] for a detailed discussion. There, the FRG results for the glue potential in YM theory [13,15,47] and QCD with two massless quark flavours [1] have been compared, see also [16] for results with 2 + 1 flavours. It was found that, apart from a rescaling, the shape of the glue potential in both theories is very similar close to T cr , see Figs. 5 and 7 in [3]. The simple linear relation is already capable of connecting the scales of both theories. In this manner, a potential V glue [Φ, is defined, where Φ,Φ denote the Polyakov loop and its conjugate. Note that the relation Eq. (4) holds only for small and moderate temperatures, as the slope saturates at high scales, where the perturbative limit is reached. In the following, Eq. (4) is used to account for the scale mismatch introduced by the fit of the PQM glue potential to YM lattice data. The only quantity left to fix is then the critical temperature of the glue sector, T glue cr . This value can in principle also be deduced from the QCD glue potential, see [1], and yields T glue cr (N f = 2) = 203 MeV. Since the absolute scale in [1] was not computed in a chiral extrapolation of the theory with physical quark masses, we consider T glue constrained by the estimates in [25,29].

III. POLYAKOV-QUARK-MESON MODEL
In the following we provide some details of the Polyakov-quark-meson (PQM) model [25,29,31] and discuss the corresponding FRG flow equation at leading order in an expansion in derivatives.
The chiral sector of this model is given by the wellknown quark-meson model [48][49][50][51][52]. The integration of the gluonic degrees of freedom results in a potential for the Polyakov-loops (Φ(A 0 ),Φ(A 0 )). They are coupled to the matter sector via the quark fields.

A. Setup
The Euclidean Lagrangian for the PQM model reads with a flavour-blind Yukawa coupling h and the covariant derivative D / (Φ) = γ µ ∂ µ − i gγ 0 A 0 (Φ) coupling the quark fields to the Polyakov loop. In this work we assume isospin symmetry in the light sector and use a flavourblind chemical potential µ . The mesonic Lagrangian is given by [4,53,54] Here, the field Σ is a complex (3 × 3)-matrix where σ a denotes the scalar and π a the pseudo-scalar meson nonets and the Hermitian generators of the flavour U (3) symmetry are defined via the Gell-Mann matrices as T a = λ a /2 . It is advantageous to rotate into the nonstrange-strange basis via with σ x the non-strange and σ y the strange condensate. Then, the explicit symmetry breaking term consistent with isospin symmetry takes the simple form Tr C(Σ + Σ † ) → c x σ x + c y σ y , with c x , c y governing the bare light and strange quark masses, respectively. The meson potential U can be expressed via the chiral invariants ρ i = Tr (ΣΣ † ) i , i = 1, . . . , N f [55]. In the (2 + 1)-flavour approximation, where σ 3 = 0, ρ 3 can be expressed in terms of the other invariants and we use the set {ρ 1 ,ρ 2 , ξ} with Here, ξ represents the 't Hooft determinant [56,57], rewritten in the mesonic language [58,59], and as such implements the chiral U A (1) anomaly. The strength of its coupling, c, determines the mass splitting between the η, η ′ and pions, see e.g. [4,54,60] for a detailed discussion. Furthermore, the quasi-particle energies of the quarks and mesons are given by for the light and strange quarks, respectively, and for the mesons. Here, H Σ (.) denotes the Hessian w.r.t. Σ and eig{.} denotes the set of eigenvalues of the given operator. For further details on this model we refer the reader to [4,53,54]. The two-flavour case considered in Sec. IV B is obtained by omitting the strange quark sector as well as all mesons except the sigma and pions.
What is now left is to specify the glue potential V glue . We have argued above that we can use the YM-based parametrisations U of the glue potential and modify the scale according to QCD FRG results, Eq. (4). Several parametrisations of the Polyakov-loop potential have been put forward in the recent years [18, 21,23,27,36,45,46]. In the main text we only show results for a polynomial version, introduced in [18, 21] The temperature-dependent coefficient, expressed in terms of the reduced temperature, is given by The parameters a i , b i of Eqs. (13) and (14) have been determined in [21] by a fit to pure Yang-Mills lattice results to be a 0 = 6.75 , a 1 = −1.95 , a 2 = 2.625 , a 3 = −7.44 (15) and We use the lattice result for the pressure to fix the open parameter T glue cr = 210 MeV in U(Φ,Φ; t). A discussion of the dependence of our results on this choice and on the parametrisation of U can be found in the appendix.

B. Fluctuations in the PQM model
In the present work we go beyond the mean-field approximation used in [3] and apply the FRG to include quantum and thermal fluctuations of the PQM model. This provides us with a more realistic description of the chiral and deconfinement phase transitions. In fact it has been shown previously, see e.g. [52], that fluctuations smear out the phase transition, yielding smoother transitions that are in better agreement with lattice results.
The flow equation for the two-flavour PQM model has been derived previously in [29,31] while the flow of the (2 + 1)-flavour quark-meson model is discussed in depth in [4]. It is then straight-forward to deduce the flow equation of the full (2 + 1)-flavour PQM model The Polyakov-loop extended quark/anti-quark occupation numbers are given by Note that we restrict ourselves to leading order in a derivative expansion and neglect the running of any couplings involving quark interactions. The RG running of the mesonic couplings, on the other hand, is encoded in the scale-dependent effective potential Ω k .
In order to solve the flow Eq. (17), we have to specify an initial potential at the cutoff scale Λ. In this work we have chosen Λ = 1 GeV, in accordance with our interpretation of the quark-meson model as a low-energy effective description. We keep only renormalisable terms in the mesonic potential at the cutoff scale and restrict ourselves to only two chiral invariants ρ 1 ,ρ 2 , cf. the discussion in [4] The parameters are fixed to a 10 = (972.63 MeV) 2 , a 01 = 50, a 20 = 2.5 which, together with the choices h = 6.5 for the Yukawa coupling, c = −4807.84 MeV for the 't Hooft-determinant coupling and explicit breaking strengths c x = (120.73 MeV) 3 and c y = (336.41 MeV) 3 , reproduces the physical spectrum as well as the pion and kaon decay constants in the vacuum [61]. In particular, we have chosen a sigma-meson mass of m σ = 400 MeV. In Appendix 3 we discuss the dependence of our results on this choice. At temperatures 2πT Λ , thermal fluctuations become important also at scales above the cutoff k > Λ. These thermal fluctuations are, however, not taken into account in the solution to the flow Eq. (17) with finite cutoff Λ. Therefore, the initial potential Ω Λ is not fully independent of temperature, which is quantitatively important in the region 2πT Λ. On the other hand, this temperature dependence of the initial potential Ω Λ is also governed by the flow Eq. (17) and can be obtained by integrating the vacuum flow from the cutoff Λ up to a scaleΛ ≫ 2πT and subsequently integrating the finite temperature flow down to Λ again This procedure is equivalent to a change of the initial scale from Λ toΛ, while keeping the infrared physics fixed, i.e. a change in the renormalisation scale. However, as we expect mesonic fluctuations to be less important at scales k > Λ = 1 GeV, we approximate the difference by the purely fermionic contribution to Eq. (17).
Since the fermionic contribution to the flow is independent of Ω k , the approximate temperature dependence of Ω Λ is given by the simple integral [29,62] Here, we have chosenΛ = ∞, since the fermionic difference flow is finite. Finally we obtain for the initial potential at the cutoff scale Λ, including fermionic temperature corrections.

A. QCD Crossover
From the solution to the flow Eq. (17) we can determine the phase structure and thermodynamics of the PQM model. For the time being, we restrict ourselves to vanishing chemical potential. This has the advantage that in this limit the Polyakov loop and its conjugate coincide,Φ = Φ. Hence, the numerical effort to solve the equations of motion (EoM) which determine the order parameters χ 0 = (σ x , σ y , Φ) for given temperature and chemical potential, is drastically reduced. A discussion of the numerical method used to solve this multi-dimensional system of partial differential equations can be found in [4,64]. In Fig. 2 our result for the subtracted chiral condensate, is shown in comparison with the lattice result by the Wuppertal-Budapest collaboration [63]. Due to the finite quark masses, we find a smooth crossover and there is no exact definition of the transition temperature. Nevertheless, it is customary to associate a transition temperature with the peak position of the temperature derivative of the order parameter, d∆ l,s /dT . Using this definition, we obtain T χ = 172 MeV for the chiral crossover temperature, and similarly T d = 163 MeV for the Polyakov-loop related transition via dΦ/dT . Both values agree roughly with the transition region on the lattice (147 − 165) MeV [63] and the pseudocritical temperatures T χ = 157 ± 3 MeV [63] and T χ = 154 ± 9 MeV [67]. Note that, apart from a shift along the temperature axis,the slope of our FRG result for the subtracted condensate coincides with the lattice one, cf. Figs. 2 and 3. This indicates that the relative strength of the relevant dynamics is included properly. However, there is a difference in the absolute scale, T χ , in our calculation and the lattice. This is to some extent related to our choice of the sigma meson mass. From experiment it is known that the σ (f 0 (500)) is a broad resonance, (400 − 550) − i(200 − 350) MeV [68]. It has been shown previously, see e.g. [54], that a lower sigma mass results in a lower chiral transition temperature, while the slope of the condensate is only changed marginally. For the curve shown in Fig. 2 the value m σ = 400 MeV has been used. We have checked that for lower values, entailing stronger mesonic fluctuations, the result is shifted even closer to the lattice one. However, the value of m σ = 400 MeV is already at the lower experimental boundary, hence we refrain from using a lower mass in the following. The impact of this mass parameter on thermodynamics is also discussed in Appendix 3.
The axial anomaly similarly influences the transition. It has been demonstrated in [4] that the transition temperature is reduced for vanishing anomaly coupling, c = 0. In fact, with our choice of m σ = 400 MeV and c = 0, the resulting condensate lies almost exactly on top of the lattice points. Note that in this case the η ′ meson would be an additional pseudo-Goldstone boson, again leading to enhanced fluctuations. However, this choice is unphysical and results in, e.g., a too high pressure at low temperatures. The use of a temperature-dependent anomaly coupling, c(T ), is expected to resolve this issue. In summary, we have found that for a correct description of the absolute scale, further mesonic fluctuations need to be included. Within our FRG treatment this would correspond to higher mesonic operators in our cutoff potential, Ω Λ . In full QCD such contributions are dynamically generated at higher scales, but we have omitted them in the present work since we restrict our cutoff action to contain renormalisable operators only.
We conclude that while the absolute scale, T χ , differs from the lattice one by about 10% in our FRG calculation, the relative strength of the relevant dynamics of the transition is captured well. This is due to the inclusion of unquenching effects as well as matter fluctuations. We postpone the improvement of our scale-setting procedure to future work and concentrate on the discussion of the dynamics of the transition in the following. To this end, all results are expressed in terms of the reduced temperature t = (T − T χ )/T χ . This choice allows us to compare the overall shape -and thereby the proper inclusion of the relevant dynamics -of the observables, while a mismatch of the critical temperatures is scaled out. Fig. 3 shows the subtracted chiral condensate (left) as well as the Polyakov loop (right) in terms of the reduced temperature. As argued above, we observe excellent agreement between the FRG (solid line) and lattice (symbols) result for the chiral condensate, especially at temperatures below T χ . In turn, for temperatures above the transition the present model overestimates the importance of mesonic fluctuations, and the FRG result for the order parameter is above the lattice result. The use of dynamical hadronisation, [7,[41][42][43], should compensate this effect.
For comparison, in Fig. 3 we also show results from a PQM mean-field calculation without ("MF", dotted line) and with ("eMF", dashed line) the fermionic vacuum loop contribution [32,[69][70][71]. Note that we fix the remaining parameters of the model, m σ and T glue cr , by comparing the pressure to the lattice pressure, cf. discussion in Appendix 2. In this manner, the effect of fluctuations is partially included in the model parameters. This results in different parameter values for the meanfield and FRG calculations. We use m σ = 500 MeV and T glue cr = 210 MeV for the standard mean-field calculation, as discussed in [3,40] and m σ = 400 MeV, T glue cr = 260 MeV for the extended mean-field calculation. However, it is clear from, e.g., Fig. 3, that a modification of the parameters is not sufficient to describe the full dynamics of the transition. The inclusion of fluctuations, as done in our FRG setup, is crucial to reproduce the slope of the order parameter as well as thermodynamic observables correctly.
In fact, the pseudocritical temperature of the standard mean-field approximation is closer to the lattice one than our FRG result, T χ = T d = 158 MeV. However, this approach neglects mesonic fluctuations, and the transition comes out too steep, see Fig. 3. Including the fermionic vacuum fluctuations yields too high pseudocritical temperatures, T χ = 181 MeV and T d = 173 MeV. Compared to the standard MF result, on the other hand, the slope of the condensate is reduced.
A word of caution needs to be added concerning the Polyakov loop, Fig. 3 (right). It is well-known that the definitions of this quantity used on the lattice, Φ , and the present continuum formulation, Φ[ A 0 ], differ and a direct comparison is not possible. In view of this, we do not expect agreement of these two observables. However, it can be shown that the continuum definition serves as an upper bound for the lattice one, Φ[ A 0 ] ≥ Φ , up to renormalisation issues, cf. [13,72]. Hence, an approximate coincidence of the respective crossover regions is still anticipated. Indeed, we find that our transition temperature, defined by the inflection point of the Polyakov loop, roughly agrees with the transition region found on the lattice. Within the FRG framework, the full quantum effective potential is defined by the effective average potential Ω k in the infrared, evaluated on the solution of the EoM, The pressure of the system is then given by the negative of the effective potential, normalised in the vacuum and serves as a thermodynamic potential, from which we can deduce other bulk thermodynamic quantities in the standard way. In particular, we are interested in the free energy density with the entropy density s = ∂P/∂T and the quark number densities n f = ∂P/∂µ f for f = u, d, s . Moreover, we consider the interaction measure which quantifies the deviation from the equation of state of an ideal gas, ǫ = 3P . We compare these quantities to results of the HotQCD collaboration, [65], using the HISQ action and temporal lattice extents of N τ = 8, 12 as well as to the continuum extrapolated results of the Wuppertal-Budapest collaboration [66].
We start our discussion of thermodynamic quantities by studying the two-flavour case. The flow equation of the two-flavour PQM model has previously been discussed in [29,31]. Here, we use the parameter set given in [2], with m σ = 540 MeV. To our knowledge, there are no recent two-flavour lattice results for thermodynamics available. This entails that we cannot fix the remaining parameter, T glue cr , to, e.g., the lattice pressure. Instead we have chosen the same value as for 2 + 1 flavours, T glue cr = 210 MeV, see also our discussion in Sec. IV C below. This value is close to the phenomenological HTL estimate, T glue cr (N f = 2) = 208 MeV put forward in [25,29] and the FRG estimate, T glue cr = 203 MeV of [1]. Furthermore, this choice results in almost degenerate chiral and Polyakov-loop critical temperatures. Having fixed all parameters, the results of the present section serve as a prediction for the thermodynamics of two-flavour QCD.
In Fig. 4 we show the pressure (left) and interaction measure (right), both normalised by T 4 . For comparison we also show the lattice results for 2 + 1 flavours. Despite the fact that we only consider two quark flavours here, the overall agreement is rather good. At low temperatures, the lightest mesonic degrees of freedom, the pions, are expected to dominate the pressure. These are already included in the two-flavour model. At high temperatures, on the other hand, the two-flavour FRG result underes-timates the (2 + 1)-flavour lattice value. This is expected due to the additional third quark species contributing to the lattice pressure at high T .
The right panel of Fig. 4 shows the interaction measure. Deviations from the lattice are more pronounced in this quantity due to the presence of derivatives in its definition. Of course we do not expect perfect agreement between our two-flavour computation and the N f = 2 + 1 lattice result. However, the strongest modifications are expected around the phase transition, where there are more light degrees of freedom contributing to the thermodynamics for 2 + 1 flavours. In the low and high temperature regimes the quarks and mesons are heavy, respectively, and hence contribute less to the thermodynamic observables. This explains the surprisingly good agreement between the two and three flavour results. In fact, we find reasonable agreement with the lattice data below the phase transition, t ≤ 0. While the peak height is underestimated, the increase in ∆/T 4 around T c is similar to the lattice. Above T c , the two-flavour curve lies below the lattice result. C. Thermodynamics: N f = 2 + 1 Next, we turn to the (2+1)-flavour model. Here, we can directly compare to the available (2+1)-flavour lattice results and fix our open parameter, T glue cr by comparison of the pressure. Fig. 5 (left, solid line) then shows our result for the pressure, normalised by T 4 , which agrees very well with the lattice result in the continuum limit. Near t = 0 this is a consequence of our choice of T glue cr . The nice agreement with lattice data away from T χ , on the other hand, indicates that we have included all relevant degrees of freedom, especially below the transition temperature. The grey band gives an error estimate of our FRG result, which is obtained from the change of the threshold functions with respect to the temperature, at vanishing mass at the ultraviolet cutoff Λ. This results in The propagation of uncertainty in the interaction measure as a derived quantity has been taken into account via Also shown in Fig. 5 are the mean-field results. To achieve a better description of the thermodynamics at low T , we have augmented the MF and eMF results by the contribution of a thermal pion gas, where the pion in-medium mass is determined by the mean-field potential. Our two flavour FRG calculation, see Sec. IV B, confirms that these are the relevant degrees of freedom below the phase transition. To highlight the impact of this pion contribution, we also show results for the pure eMF calculation in Fig. 5 (yellow, solid line). Strictly speaking, this contribution picks up a field-dependence via the in-medium pion mass, which would modify the equations of motion. Here, however, we consider it as a correction to the thermodynamic potential only, and hence neglect its backcoupling on the equations of motion. For consistency, we also neglect all terms containing field derivatives, ∂P π /∂φ i , in higher thermodynamic observables.
While the pressure of the mean-field approximation including pions (dotted line) lies above the FRG and lattice ones, the inclusion of the vacuum term (dashed line) results in an additional increase in P at low t and a decrease at high t. As expected, the omission of the pion contribution in the eMF calculation (solid, yellow line) yields a pressure that is too low, especially below the phase transition, where pions are expected to dominate the thermodynamics.
The interaction measure is displayed in the right panel of Fig. 5. Similarly to the other observables, also the interaction measure is too steep within the standard meanfield approximation. Including the vacuum term, the transition is smoothened out and already agrees quite well with the Wuppertal-Budapest results. Turning to our FRG result (red, solid curve), we find remarkably good agreement with the continuum extrapolated lattice result from the Wuppertal-Budapest collaboration. There is a stronger deviation from the HotQCD data, but we attribute this to the lacking continuum limit of their data. In fact it is observed that the peak height of ∆/T 4 goes down as the continuum is approached, cf. [65,66]. Although not shown explicitly, we want to stress that the drastic reduction of the peak height in the interaction measure towards the lattice results is due the inclusion of t YM (t glue ) in both the FRG and mean-field approaches, see also [3,40].
In comparison to our two-flavour result, we see that the inclusion of the heavier strange quark and especially the full scalar and pseudoscalar meson nonets increases the peak height of the interaction measure and puts our curve right on top of the lattice result. At high temperatures we find that ∆/T 4 decreases too slowly in our calculation. However, this is the region where our scale matching procedure, Eq. (4), ceases to be valid and corrections are expected.

V. CONCLUSION AND OUTLOOK
We have investigated order parameters and thermodynamic observables of two and 2+1 flavour QCD within effective Polyakov-quark-meson models. This type of models can be systematically related to full QCD, as e.g. discussed in [2,3]. Thus far, the glue sector of these models is badly constrained. One often resorts to a Ginzburg-Landau-like ansatz for the glue potential obtained from fits to lattice Yang-Mills theory.
Recently, first-principles continuum results for the unquenched glue potentials have become available. These have been used to augment the glue sector considerably in a mean-field approach to the PQM model [3]. It was shown that by a simple rescaling of the temperature in the standard Yang-Mills based Polyakov-loop potentials one can already capture the essential glue dynamics of the unquenched system In the present work, we have extended the previous investigation [3] by additionally including thermal and quantum fluctuations via the functional renormalisation group. A comparison to lattice QCD simulations with 2 + 1 flavours shows excellent agreement up to temperatures of approximately 1.3 times the critical temperature. Therefore, we conclude that most of the relevant dynamics for the QCD crossover can already be captured within the PQM model. Additionally, we have put forward a prediction for the pressure and interaction measure for two quark flavours, where no recent lattice results with physical masses are available.
The present work serves as a benchmark of our system at vanishing chemical potential, which allows us to conclude that we have all relevant fluctuations included. Since our approach is not restricted to the zero chemical potential region we can now aim at the full phase diagram, µ ≥ 0 .

APPENDIX: PARAMETER DEPENDENCE
In this appendix we estimate the parameter dependence of our results. In particular, we discuss the influence of the Polyakov-loop potential chosen, our choice of the glue critical temperature T glue cr and the sigma-meson mass.

Polyakov-loop Potential
In the Polyakov-loop extended chiral models one is free to choose a parametrisation of the Polyakov-loop potential. It is customary to employ a Landau-Ginzburg-like ansatz and fit the parameters to available lattice data. However, in this manner only the region close to the minimum is constrained, not the overall shape. This is the reason why several different functional forms have been chosen in the past. In practise, when the Polyakov-loop is coupled to the matter sector, regions away from the minimum are probed and one should not expect to find exactly the same results with different versions of the potential.
In the main text we have presented results for a polynomial parametrisation [21]. Here, we want to compare these results to those using a logarithmic version of the potential [23] U log (Φ,Φ; t) In this variant, the logarithmic form arises from the integration of the Haar measure and constrains the Polyakovloop variables Φ,Φ to values smaller than one. Furthermore, we show results with the parametrisation recently proposed in [46] U log−II (Φ,Φ; t) The parameters of Eqs. (31) and (32) have been fixed in [23] and [46], respectively. Tab. I summarises our FRG results for the transition temperatures using the initial values given in Sec. III B and these different Polyakov-loop potentials. Using the polynomial Polyakov-loop potential, the chiral and deconfinement transitions lie closer to each other than with the logarithmic ones.
In Fig. 6, the pressure (left) and interaction measure (right) are shown for N f = 2+1 and the three parametrisations. The pressure for the two logarithmic versions lies below the lattice result at low temperature, while the overall agreement is quite good also in this case. In the interaction measure, however, differences are seen more clearly. The trace anomaly resulting from the standard logarithmic potential rises much steeper than the one obtained with the polynomial version. The peak height, on the other hand, agrees well with the results of the Wuppertal-Budapest collaboration, independent of the parametrisation of the potential.  as obtained from fixing the pressure has been used.

Glue Critical Temperature
Next, we discuss the impact of the glue critical temperature, T glue cr . In the main text we have chosen T glue cr = 210 MeV for 2+1 flavours with physical anomaly strength. This value has been obtained by a comparison of the resulting pressure to the lattice results of the Wuppertal-Budapest collaboration. Fig. 7 (left) summarises the dependence of the pressure on T glue cr . In terms of the reduced temperature t, the effect of the glue criti-cal temperature is to stretch the transition region, i.e. the transition becomes less steep for larger T glue cr .

Sigma-Meson Mass
The sigma meson, f 0 (500) is a rather broad resonance, (400 − 550) − i(200 − 350) MeV, which leaves us some freedom to fix this mass in our setup. The choice of this mass influences, e.g. the position of the phase transition and the location of a possible critical endpoint in the phase diagram [54]. For 2+1 flavours we have fixed m σ = 400 MeV, see also our discussion in Sec. IV A. In Tab. II we show the critical temperatures for different choices of the sigma-meson mass. Although the absolute value T χ is quite susceptible to m σ , we find that the thermodynamic observables in terms of the reduced temperature are not. This is demonstrated for the interaction measure in the right panel of Fig. 7. We have to stress, however, that this m σ independence is partially due to the fact that a change in the sigma-meson mass can be compensated to some degree by a change in T glue cr , see also [3,40].