QCD phase structure at finite temperature and density

We discuss the phase structure of QCD for $N_f=2$ and $N_f=2+1$ dynamical quark flavors at finite temperature and baryon chemical potential. It emerges dynamically from the underlying fundamental interactions between quarks and gluons in our work. To this end, starting from the perturbative high-energy regime, we systematically integrate out quantum fluctuations toward low energies by using the functional renormalization group. By dynamically hadronizing the dominant interaction channels responsible for the formation of light mesons and quark condensates, we are able to extract the phase diagram for ${\mu }_B/T\lesssim 6$. We find a critical endpoint at $(T_{\mathrm{CEP}},{\mu }_{B_{\mathrm{CEP}}})=(107,635)\mathrm{MeV}$. The curvature of the phase boundary at small chemical potential is $\kappa =0.0142(2)$, computed from the renormalized light chiral condensate ${\mathrm{\Delta }}_{l,R}$. Furthermore, we find indications for an inhomogeneous regime in the vicinity and above the chiral transition for ${\mu }_B\gtrsim 417\mathrm{MeV}$. Where applicable, our results are in very good agreement with the most recent lattice results. We also compare to results from other functional methods and phenomenological freeze-out data. This indicates that a consistent picture of the phase structure at finite baryon chemical potential is beginning to emerge. The systematic uncertainty of our results grows large in the density regime around the critical endpoint and we discuss necessary improvements of our current approximation toward a quantitatively precise determination of QCD phase diagram.

Figure 1

Phase diagram for $N_f=2+1$ flavor QCD in comparison to other theoretical approaches and phenomenological freeze-out data. Our result for the chiral crossover is depicted by the black dashed line. The crossover temperature has been determined through the peak position of the thermal susceptibility of the renormalized light chiral condensate, ${\partial }_T{\mathrm{\Delta }}_{l,R}$, at fixed baryon chemical potential ${\mu }_B$. For more details see Sec. 5a, and in particular Fig. 10. We show dotted black lines for ${\mu }_B/T=2$, 3 to indicate the reliability bounds for the lattice and functional methods. The phase boundary globally agrees well with recent lattice results. In particular the curvature of the phase boundary for small chemical potential, $\kappa =0.0142(2)$, is consistent with recent lattice results, $\kappa =0.0149(21)$ in [46], $\kappa =0.0144(26)$ in [49], and $\kappa =0.015(4)$ in [51], for an overview see [64]. We find a critical endpoint at $(T_{\mathrm{CEP}},{\mu }_{B_{\mathrm{CEP}}})=(107,635)\mathrm{MeV}$. Indications for an inhomogeneous regime (minimum of pion dispersion at nonvanishing spatial momentum and sizable chiral condensate) close to the chiral phase transition for ${\mu }_B\gtrsim 420\mathrm{MeV}$ are depicted by the hatched red area. For quantitative statements in this area the current approximation has to be upgraded systematically. Accordingly the hatched red area also serves as a reliability bound for the current approximation. For more details see Sec. 5b and Fig. 21. Other theoretical results: lattice QCD with analytic continuation from imaginary chemical potential [46] (WB), lattice QCD with Taylor expansion about vanishing chemical potential [51] (HotQCD), DSE approach with backcoupled quarks and a dressed vertex [37] (Fischer et al.), and DSE calculations with a gluon model [65] (Gao et al.). Freeze-out data: [2] (STAR), [66] (Alba et al.), [3] (Andronic et al.), [67] (Becattini et al.), [68] (Vovchenko et al.), and [69] (Sagun et al.). Note that freeze-out data from Becattini et al. with (light blue) and without (dark green) afterburner-corrections are shown in two different colors.

Figure 2

Flow of the effective action of QCD. The first three diagrams arise from the gluon, ghost, and the quark degrees of freedom respectively. The last diagram is that of the mesonic contribution. The double line with the up-down arrows indicates the nature of the mesons as quark-antiquark composites. The crossed circles indicate the regulator insertion in the flow equation.

Figure 3

$N_f=2+1$ gluon dressing function $1/Z_A(p^2)$ and gluon propagator $G_A=1/[Z_A(p^2)p^2]$ (inset) in the vacuum as functions of momenta in comparison to lattice results based on the $N_f=2+1$ configurations (domain wall fermions) of the RBC/UKQCD collaboration, e.g., [133]. The lattice results are continuum extrapolated from data with $\beta =2.25,a^{-1}=2.359(7)\mathrm{GeV},m_{\pi }=139.2(5)\mathrm{MeV}$ and $\beta =2.13,a^{-1}=1.730(4)\mathrm{GeV},m_{\pi }=139.2(4)\mathrm{MeV}$, see [134, 135, 136]. The momentum dependence of the fRG results is given by $p^2=k^2$, in accordance with the same identification of the $N_f=2$ flavor input data; for more details see Sec. 4a2.

Figure 4

Flow equation for the four-quark coupling in the present approximation. Further diagrams with the nonresonant part of the scalar-pseudoscalar channel and other four-quark tensor structures, as well as that with higher order couplings are dropped. The black dots denote the full vertices and the various lines represent for the full propagators. The derivative ${\tilde{\partial }}_t={{\partial }_t|}_{\{{\mathrm{\Gamma }}_k^{(n)}\}}$ only hits the $k$-dependence of the regulator, e.g., ${\tilde{\partial }}_t{G}_k=-G_k{\dot{R}}_k{G}_k$. The first line with the quark-gluon diagrams generates the four-quark interaction through gluon-exchange at large scales. The second line takes into account the self-interaction of the resonant four-quark interaction in terms of meson exchange diagrams. The last line comprises the mixed diagrams with meson and gluon exchange. Gluon-exchange dominates the flow for large and intermediate energy scales, while the meson-exchange diagrams dominate at low energy, see Appendix pp12. Accordingly, in the present work we keep the quark-gluon and quark-meson diagrams and drop the mixed ones.

Figure 5

Flow equations for the inverse propagators of quarks, gluons, mesons respectively, in the present approximation; in particular the tadpole contribution from two-quark–two-meson as well as the remnant four-quark scatterings are missing, as is the ghost-gluon tadpole. The derivative ${\tilde{\partial }}_t$ is defined in the caption of Fig. 4. The left-hand sides stand for, in a slight abuse of notation, ${\partial }_t{\mathrm{\Gamma }}_k^{(2)}={\partial }_t(G_k^{-1}-R_k)$, where we have dropped the $-R_k$ in the depiction.

Figure 6

Flow equations for the three-point functions computed here: quark-gluon, three-gluon and Yukawa couplings. In the present approximation, the diagrams with higher order vertices and the remnant four quark vertices are dropped. ${\tilde{\partial }}_t$ is defined in the caption of Fig. 4.

Figure 7

Left panel: $N_f=2+1$ constituent mass ${\overline{m}}_l$ of light quarks as a function of temperature $T$ at various baryon chemical potentials ${\mu }_B$. Right panel: Polyakov loops $L$, $\overline{L}$ as functions of temperature $T$ at different baryon chemical potentials ${\mu }_B$.

Figure 8

$N_f=2+1$ effective four-quark coupling ${\overline{h}}^2/(2{\overline{m}}_{\pi }^2)$ in the pseudoscalar channel as a function of the RG scale $k$ at different temperatures and ${\mu }_B=0$.

Figure 9

$N_f=2+1$ effective four-quark coupling ${\overline{h}}^2/(2{\overline{m}}_{\pi }^2)$ at vanishing cutoff scale, $k=0$, in the pseudoscalar channel as a function of temperature $T$ at various baryon chemical potentials ${\mu }_B$.

Figure 10

Left panel: $N_f=2+1$ renormalized light chiral condensate, ${\mathrm{\Delta }}_{l,R}$, see Eqs. (a3), (a7), as a function of temperature $T$ at various baryon chemical potentials ${\mu }_B$. The normalization constant in (a3), (a7) is chosen to match the scale in the lattice calculation [187]. Right panel: Thermal susceptibility of $N_f=2+1$ renormalized light chiral condensate, ${\partial }_t{\mathrm{\Delta }}_{l,R}$, as a function of temperature $T$ at various baryon chemical potentials ${\mu }_B$.

Figure 11

$N_f=2+1$ meson masses ${\overline{m}}_{\pi }$, ${\overline{m}}_{\sigma }$ as function of the temperature at different baryon chemical potentials ${\mu }_B$. They are rescaled by the multiplicative constant $z_{\varphi }\equiv ({\overline{Z}}_{\varphi }/Z_{\varphi }(0){)}_{T,{\mu }_B=0}^{1/2}$. As a result, $z_{\varphi }\overline{m}$ coincides with the pole mass defined in (88) and (98) in the vacuum. The in-medium behavior is that of the curvature masses defined in (43). For more details see Sec. 4a1.

Figure 12

Left panel: $Z_{\varphi }(0)$ at vanishing cutoff scale as a function of the temperature $T$ for various baryon chemical potentials ${\mu }_B$. Negative $Z_{\varphi }(0)$, see (90), entails a negative spatial momentum slope of the dispersion ${\mathrm{\Gamma }}_{\varphi \varphi }^{(2)}(p)$ at vanishing momentum. This indicates a minimum of the dispersion at a finite ${\mathit{p}}^2\ne 0$. Right panel: $|1/Z_{\varphi }(0)|$ for the data in the left panel. Between the spikes the wave function renormalization is negative at ${\mathit{p}}^2=0$.

Figure 13

$N_f=2+1$ gluon dressing function $1/Z_A$ as a function of spatial momenta $|\mathit{p}|$ at different temperatures in comparison to the $N_f=2+1+1$ lattice data from [197], based on the tmfT configurations [198], for the tmfT simulation setup see also [199, 200]. The lattice results depicted here are obtained for $\beta =2.1$, $a=0.0646\mathrm{fm}$ and an approximate pion mass of $m_{\pi }=369(15)\mathrm{MeV}$ leading to a pseudocritical temperature of $T_{\text{deconf}}=193(13)(2)\mathrm{MeV}$, for more details see [197].

Figure 14

$N_f=2+1$ gluon dressing function $1/Z_A$ at vanishing Matsubara frequency, $p_0=0$, as a function of spatial momenta $|\mathit{p}|$ at different temperatures (left panel) and baryon chemical potentials (right panel). As for the $N_f=2$ flavor input data and the $N_f=2+1$ vacuum results, see Fig. 3, we use the identification ${\mathit{p}}^2=k^2$.

Figure 15

$N_f=2+1$ quark (left panel, $Z_q$) and meson (right panel, ${\overline{Z}}_{\varphi }$) wave function renormalizations at vanishing cutoff scale, $k=0$, as functions of temperature $T$ at different baryon chemical potentials ${\mu }_B$.

Figure 16

Left panel: $N_f=2+1$ Yukawa coupling at vanishing cutoff scale, $k=0$, as a function of the temperature $T$ at different baryon chemical potentials ${\mu }_B$. Right panel: $N_f=2+1$ Yukawa coupling as a function of the RG scale $k$ at different temperatures $T$ and vanishing baryon chemical potential ${\mu }_B=0$.

Figure 17

$N_f=2+1$ quark-gluon coupling for light quarks (${\alpha }_{\overline{l}Al}$) and strange quarks (${\alpha }_{\overline{s}As}$) quarks, and three-gluon coupling (${\alpha }_{A^3}$) as functions of the RG scale $k$ for several values of the temperature $T$ and for vanishing baryon chemical potential, ${\mu }_B=0$.

Figure 18

Dimensionless four-quark single gluon exchange coupling ($g_{\overline{l}Al}^2$ and $g_{\overline{s}As}^2$) and four-quark single meson exchange couplings (${\overline{h}}^2/(1+{\tilde{m}}_{\pi }^2)$ and ${\overline{h}}^2/(1+{\tilde{m}}_{\sigma }^2)$) as functions of the cutoff scale $k$ in the vacuum. Gluons, quarks, and mesons decouple sequentially from the matter dynamics.

Figure 19

Dimensionless propagator gapping $1/(1+{\tilde{m}}_{{\mathrm{\Phi }}_i}^2)$ for ${\mathrm{\Phi }}_i=l,s,\sigma ,\mathit{\pi }$. For comparison, we also show the gluon dressing function $Z_A^{\mathrm{peak}}/Z_A$ with the solid blue line. Here, $1/Z_A$ is normalized by its peak position at $k_{\mathrm{peak}}$, that is $Z_A^{\mathrm{peak}}=Z_A(k_{\mathrm{peak}})$: for cutoff scales $k\gtrsim k_{\mathrm{peak}}$ the $k$-dependence of $Z_A$ is dominated by the running of the gluon wave function renormalization, for $k\lesssim k_{\mathrm{peak}}$ it shows the gluon gapping. In summary, the sequential decoupling of gluon, quark and meson dynamics toward the IR is evident here.

Figure 20

Phase diagram in the plane of the temperature and the baryon chemical potential. The gray and blue bands denote the crossover transitions for the $N_f=2$ and $2+1$ flavor QCD, respectively; and the red star and black circle are their relevant CEP. The bands are determined through the 80% peak height of $\partial {\mathrm{\Delta }}_{l,R}/\partial T$ at fixed ${\mu }_B$. The black and red dashed lines depict the peak positions for the $N_f=2+1$ and $N_f=2$, respectively. We also provide the two lines of ${\mu }_B/T=2$, 3 related to reliability bounds for both lattice and functional methods.

Figure 21

Phase diagram for $N_f=2+1$ flavor QCD in comparison to freeze-out data. The crossover temperature has been determined through the peak position of the thermal susceptibility of the renormalized light chiral condensate ${\mathrm{\Delta }}_{l,R}$, see (a3), (a7): ${\partial }_T{\mathrm{\Delta }}_{l,R}$ at fixed baryon chemical potential ${\mu }_B$. For more details see Sec. 5. The phase boundary globally agrees well with recent lattice results, and in particular the curvature of the phase boundary for small chemical potential, $\kappa =0.0142(2)$, is consistent with recent lattice results. We find a critical endpoint at $(T_{\mathrm{CEP}},{\mu }_{B_{\mathrm{CEP}}})=(107,635)\mathrm{MeV}$. fRG: inhom: The blue area depicts the regime with a negative slope of the mesonic dispersion at vanishing spatial momentum, that is $Z_{\varphi }(0)<0$, see Fig. 12. There, the meson dispersion has a minimum at a nonvanishing spatial momentum. This is a strong indication for an inhomogeneous regime. We also show the minimum of $Z_{\varphi }(0)$ in $T$ at a fixed value of ${\mu }_B$ with the blue dashed line. The red hatched area shows where the inhomogeneous regime overlaps with a sizable homogeneous chiral condensate. The latter is defined as the size of the (reduced) chiral condensate where the thermal susceptibility $\partial {\mathrm{\Delta }}_{l,R}/\partial T$ at fixed ${\mu }_B$ reaches 80% of its peak height (from above), and larger. In this region, a competition between homogeneous and inhomogeneous condensation is expected and the phase structure could be altered significantly. For a detailed discussion see Sec. 5b. Freeze-out data: [2] (STAR), [66] (Alba et al.), [3] (Andronic et al.), [67] (Becattini et al.), [68] (Vovchenko et al.), and [69] (Sagun et al.). Note that freeze-out data from Becattini et al. with (blue) and without (dark green) afterburning corrections are shown in two different colors.

Figure 22

$N_f=2+1$ renormalized light chiral condensate ${\mathrm{\Delta }}_{l,R}$, in Eqs. (a3) and (a7) as a function of the temperature with ${\mu }_B=0$, in comparison to the lattice result in [187]. The normalization constant in (a3) and (a7) is chosen to match the scale in the lattice calculation.

Figure 23

The reduced condensate in the approximation (a17) for the constituent quark mass difference $\mathrm{\Delta }{\overline{m}}_{sl}={\overline{m}}_s-{\overline{m}}_l=120$, 150, 155, 160 MeV for current quark mass ratios of $m_s^0/m_l^0=c_{{\sigma }_s}/c_{\sigma }\approx 14$, 27, 30, 34 in comparison to the lattice results in [187]. We observe quantitative agreement for all temperatures for ratios close to the physical one from $N_f=2+1$ flavor lattice simulations, $m_s^0/m_l^0\approx 27$ from [147].

Figure 24

Comparison between the $N_f=2$ propagators obtained from lattice gauge theory and fRG results at vanishing temperature and density. The solid red thin line shows the gluon propagator of the pure gauge theory we use as input. The open red circles are the corresponding lattice results [236]. The blue open squares are the lattice results for two-flavor QCD [237]. The solid blue thick line shows the fRG result using (d2) and the dotted black line shows the result based on (d4). Both results are taken from [28, 137]. Furthermore, results from the fRG calculations with the state-of-the-art truncations at vacuum for the pure Yang-Mills theory [31] and $N_f=2$ unquenched QCD [32] are shown for comparison, which are denoted by the green and cyan dashed lines, respectively.

Figure 25

$N_f=2$ flavor and $N_f=2+1$ gluon dressing functions $1/Z_A$ as function of momentum. The $N_f=2$ flavor gluon dressing computed with the fRG in [32] is the input in the present work. It is in quantitative agreement with the respective lattice results. The lattice data here are taken from [237]. The $N_f=2+1$ flavor gluon dressing shown here is a genuine prediction of the present computation. It is shown to be in quantitative agreement with the respective lattice results [134, 135, 136].

Figure 26

Two-flavor QCD quark mass functions from [32] and the present work for a fixed background $\sigma ={\sigma }_{\mathrm{EoM}}$ at vanishing cutoff scale. In [32] the full momentum-dependent mass $M_{q,k}(p)$ function in $N_f=2$ flavor QCD has been computed, and we have $m_{l.k}/m_{l.k=0}=M_{q,k}(0)/M_{q,k=0}(0)$.

Figure 27

$N_f=2$ and $N_f=2+1$ light quark mass functions as a function of the cutoff scale $k$, normalized by their value at vanishing cutoff scale $k=0$.

Figure 28

$N_f=2$ quark mass function from [32] and the $N_f=2+1$ light quark mass function from present work in the fixed expansion.

Figure 29

$N_f=2+1$ gluon propagator $G_A=1/(Z_{A,k}{k}^2)$ at vanishing frequency $p_0=0$ as a function of spatial momenta $|\mathit{p}|$ at different temperatures (left panel) and baryon chemical potentials (right panel).

Figure 30

Left panel: comparison between the flow of the $N_f=2+1$ four–light-quark coupling arising from the gluon exchange in (l1), ${\mathrm{\Phi }}_i=A$, and that from the meson exchange in (l2), ${\mathrm{\Phi }}_i=\varphi$, as a function of $k$ at $T={\mu }_B=0$. Right panel: the ratio of their respective flow over the total flow, i.e., the summation of (l1) and (l2), where the absolute value for the ratio is chosen. The blue vertical line in both plots shows the position where the flow resulting from the meson exchange changes sign.