The phase structure of QCD

TL;DR

This paper reviews recent advances in understanding the phase structure of QCD, emphasizing lattice simulations, critical scaling, susceptibilities, and implications for heavy ion collision experiments.

Contribution

It synthesizes recent lattice QCD results on critical phenomena, susceptibilities, and phase diagram exploration, providing insights into QCD matter properties.

Findings

Universal critical scaling observed in lattice QCD simulations.

Susceptibilities characterize properties of QCD matter at various temperatures.

Lattice results inform understanding of matter in heavy ion collisions.

Abstract

We review recent results on the phase structure of QCD and bulk QCD thermodynamics. In particular we discuss how universal critical scaling related to spontaneous breaking of the chiral symmetry manifests itself in recent lattice QCD simulations and how the knowledge on non-universal scaling parameter can be utilized in the exploration of the QCD phase diagram. We also show how various (generalized) susceptibilities can be employed to characterize properties of QCD matter at low and hight temperatures, related to deconfining aspects of the QCD transition. Finally, we address the input these lattice QCD calculation can provide for our under standing of the matter created in heavy ion collisions and in particular on the freeze-out conditions met in the hydrodynamic evolution of this matter.

Equations

$$T=\frac{1}{a{N}_{\tau }}\text{mbox}andV=(a{N}_{\sigma }{)}^3.$$

$$U_{x,\nu }=P\exp \left\{ig{\int }_x^{x+a\hat{\nu }}\mathrm{d}{x}_{\nu }{A}_{\nu }(x)\right\},$$

$$Z(T,V,\vec{\mu })=\int \prod _{n,\nu }\mathrm{d}{U}_{n,\nu }\int \prod _{n,f}\mathrm{d}{\psi }_n^{(f)}\mathrm{d}{\bar{\psi }}_n^{(f)}{e}^{-S_E(T,V,\vec{\mu })},$$

$$S_G=\sum _n\sum _{0\le \nu <\rho \le 3}\beta \left(1-\frac{1}{3}\mathrm{Re}\mathrm{Tr}{U}_{n,\nu \rho }\right).$$

$$S_F=\sum _f\sum _{n,m}{\bar{\psi }}_n^{(f)}{M}_{n,m}^{(f)}{\psi }_m^{(f)},$$

$$Z(T,V,\vec{\mu })=\int \prod _{n,\nu }\mathrm{d}{U}_{n,\nu }\prod _f\det \left[M^{(f)}({\mu }_f)\right]e^{-S_G(T,V)},$$

$$⟨\mathcal{O}⟩=\frac{1}{Z}\int \prod _{n,\nu }\mathrm{d}{U}_{n,\nu }\mathcal{O}\prod _f\det \left[M^{(f)}({\mu }_f)\right]e^{-S_G(T,V)}.$$

$$\frac{\frac{p}{T^{4}}}{(\frac{p}{T^{4}})_{\rm SB}}=1+\left\{\!\begin{array}[]{rcrcrl}\frac{248}{147}(\frac{\pi}{N_{\tau}})^{2}&\!\!+\!\!&\frac{635}{147}(\frac{\pi}{N_{\tau}})^{4}&\!\!+\!\!&\frac{13351}{8316}(\frac{\pi}{N_{\tau}})^{6}&\mbox{std. Wilson}\\ \frac{248}{147}(\frac{\pi}{N_{\tau}})^{2}&\!\!+\!\!&\frac{635}{147}(\frac{\pi}{N_{\tau}})^{4}&\!\!+\!\!&\frac{3796}{189}(\frac{\pi}{N_{\tau}})^{6}&\mbox{std. staggered}\\ &\!\!-\!\!&\frac{1143}{980}(\frac{\pi}{N_{\tau}})^{4}&\!\!-\!\!&\frac{365}{77}(\frac{\pi}{N_{\tau}})^{6}&\mbox{Naik}\\ &\!\!-\!\!&\frac{1143}{980}(\frac{\pi}{N_{\tau}})^{4}&\!\!+\!\!&\frac{73}{2079}(\frac{\pi}{N_{\tau}})^{6}&\mbox{p4}\end{array}\right..$$

$$\frac{{\Theta }^{\mu \mu }(T)}{T^4}=\frac{ϵ-3p}{T^4}=T\frac{d}{dT}\left(\frac{p}{T^4}\right)=-\frac{1}{T^{3}V}\frac{d\ln Z}{d\ln a},$$

$$\frac{p}{T^4}=\frac{1}{V{T}^3}\ln Z(T,V,m_l,\vec{\mu })=-f_s(T,V,m_l,\vec{\mu })-f_r(T,V,m_l,\vec{\mu }),$$

$$t\equiv \frac{1}{t_0}\left(\frac{T-T_c}{T_c}+{\kappa }_l{\hat{\mu }}_l^2+{\kappa }_{ls}{\hat{\mu }}_l{\hat{\mu }}_s+{\kappa }_s{\hat{\mu }}_s^2\right)\text{mbox}andh=\frac{1}{h_0}\frac{m_l}{m_s},$$

$$f_s(t,h)=h_0{h}^{1+1/\delta }{f}_f(z)=h_0{h}^{(2-\alpha )/\beta \delta }{f}_f(z)\text{mbox}withz=t/h^{1/\beta \delta }.$$

$$f_G(z)=-\left(1+{\delta }^{-1}\right)f_f(z)+z{\left(\beta \delta \right)}^{-1}{f}_f^{\prime }(z).$$

$$P({\hat{\mu }}_S,{\hat{\mu }}_B)$$

$${\chi }_n^S(\hat{\mu }=0)$$

$${\chi }_{kl}^{BS}(\hat{\mu }=0)$$

$$v_1={\chi }_{31}^{BS}-{\chi }_{11}^{BS},v_2=\frac{{\chi }_2^S-{\chi }_4^S}{3}-2{\chi }_{31}^{BS}-4{\chi }_{22}^{BS}-2{\chi }_{13}^{BS}.$$

$$P_C(T,{\mu }_B,{\mu }_C)=P_Q(T)\cosh \left(\frac{{\mu }_C+{\mu }_B/3}{T}\right)+$$

$$P_{B,C=1}(T)\cosh \left(\frac{{\mu }_C+{\mu }_B}{T}\right)+P_M(T)\cosh \left(\frac{{\mu }_C}{T}\right).$$

$$R_{12}^{QB}({\mu }_B,T)=R_{12}^{QB}(0,T)\left[R_{12}^{QB,2}-{\kappa }_2^f{T}_{f,0}\frac{d{R}_{12}^{QB}(0)}{dT}{|}_{T_{f,0}}\right]\frac{{\mu }_B^2}{T^2}$$

$$R_{12}^{QB}({\mu }_B)=R_{12}^{QB}(0)\left[1+c_{12}{\left(R_{12}^B\right)}^2\right]+\mathcal{O}{\left(R_{12}^B\right)}^4.$$

$$c_{12}(T_{f,0})=\frac{1}{(R_{12}^B{)}^2}\frac{R_{12}^{QB,2}(0)}{R_{12}^{QB}(0)},D_{12}(T)=\frac{1}{(R_{12}^B{)}^2}\frac{d{R}_{12}^{QB}(0)}{dT}{|}_{T_{f,0}}.$$

Full text

1 Introduction

The theory that describes the fundamental interaction between quarks and gluons is called Quantum Chromodynamics (QCD). Ever since its formulation people have speculated on the QCD phase diagram. QCD is constructed as a non abelian gauge theory and posses two most astonishing features: asymptotic freedom [1, 2] and a global symmetry (chiral symmetry) that is spontaneously broken at small temperatures and net baryon number densities [3]. These are the principles that not only guide our thinking of the spectrum of QCD bound states in vacuum, but give also rise to a rich and fascinating phase structure of QCD. Even before the fundamental properties of QCD had been fully revealed, it was realized that the spectrum of hadronic resonances that becomes exponentially dense at high energies, would naturally lead to the concept of a limiting temperature in a gas of hadrons [4]. Later on, with the analysis of the running QCD coupling that due to asymptotic freedom vanishes in the high energy limit, it was quickly understood that hadronic matter undergoes a transition to a new state of matter, the quark gluon plasma (QGP) [5, 6]. At low temperatures and net baryon number densities the coupling constant is so large, that quarks are tightly bound into hadrons. The static quark potential is in fact linear rising with distance, thus quarks are confined into hadrons. However, once the density of hadrons increases coursed by thermal fluctuations or a large baryon chemical potential, they start to overlap and the potential gets screened. Hence, the quarks get liberated and are allowed to move over macroscopic distances [7, 8].

It was also speculated that the deconfinement transition of QCD goes along with the global restoration of chiral symmetry, that is manifest in the QCD Lagrangian [9]. Based on the symmetries of QCD many effective models have been designed and analyzed at non zero temperature and baryon number density, which greatly enhanced our understanding of the QCD phase diagram. From these model calculations a complex phase structure with various homogeneous as well as in-homogeneous phases has emerged [10]. At very high net-quark number densities, several different pairing mechanisms among quarks have been proposed, that would lead to various color super-conducting phases, similar to electrical super-conductors in condensed matter systems [11, 12]. Unfortunately, details of those phases are highly model dependent and could so far not be addressed with first principle QCD calculations. For this reason we show in Fig. 1 a simplified version of the QCD phase diagram that mainly features a speculative critical point, that appears as a critical end-point (CEP) of a line of first order transitions. This CEP is currently subject to many lattice QCD calculations as well as experimental heavy ion programs.

The framework of lattice regularized QCD provides a non perturbative method, which has proven to be very successful since its inception [13, 14]. As the QCD coupling constant is rather large at temperatures close to the QCD transition [15] and since it is running logarithmically, it becomes small only at exponentially high energy scales. This is the reason why ordinary perturbation theory is expected to be reliable only at very high temperatures and net baryon number densities. Using appropriate re-summation schemes, as hard thermal loop (HTL) perturbation theory [16] or electric QCD (EQCD) [17] the range of validity can be extended down to a few times the transition temperature ($T_{c}$). However, close to $T_{c}$ QCD is deeply non perturbative. Here lattice QCD is the only known method that – after a proper continuum extrapolation – can address the thermodynamics of quarks and gluons. This is definitely true for the case of a vanishing net baryon number, or equivalently, a vanishing baryon chemical potential ($\mu_{B}$). At non zero chemical potential, lattice QCD suffers from the infamous sign problem. At small $\mu_{B}$, however, this problem can be overcome by a Taylor expansion in $\mu_{B}/T$ [18]. Other methods that have the potential to solve the sign problem completely are currently being developed [19, 20].

At the same time, huge experimental effort is undertaken to create and measure the thermal properties of QGP. Shortly after the formation of man-made QGP matter in a heavy ion experiment was first announced at the Super Proton Synchrotron (SPS) at CERN [21], the Relativistic Heavy Ion Collider (RHIC) at BNL went in operation. At RHIC the envisaged goal was exactly to crate and study QGP. Also at the Large Hadron Collider (LHC) a dedicated heavy ion program is conducted. Related to the dynamical expansion of the QGP medium, after its creation in the central collision region of two heavy nuclei, RHIC discovered a few very remarkable transport properties of QGP. Based on the observation of elliptic flow and jet quenching effects, it is believed that in a short temperature range directly above $T_{c}$, the QGP behaves as an almost perfect fluid, with a small shear viscosity and small transport coefficients such as electrical conductivity and diffusion coefficients. It is tempting to study those transport properties also on the lattice, which will provide crucial input to the modeling of the expansion of the QGP.

More recently, RHIC initiated a beam energy scan program (BES) aiming on the exploration of the QCD phase diagram and the discovery of the postulated CEP. Reducing the bombarding energy of the collision, leads to a variation of the initial conditions at which the QGP is created. Lower bombarding energies translate into colder but more dense systems at formation times. With the planed future facilities FAIR and NICA, one wants to increase the density even further. Approximate initial conditions of various heavy ion experiments are sketched in Fig. 1. Event-by-event fluctuation observables are supposed to be sensitive to the phase structure of QCD [22]. Among these the moments of conserved charge distributions are very interesting observables as they can also be obtained from lattice QCD calculations. From a comparison of results, as obtained on the lattice and from the experiment, one can infer state variables of the system at the time of freeze-out [23].

In this article we review recent selected (lattice) QCD results on QCD bulk thermodynamics that are relevant for heavy ion phenomenology. After a short presentation of the basic lattice methods, we start with the equation of state in Sec. 2. We will then discuss universal critical behavior, that can be observed in the lattice results close to second order phase transition points. We emphasize how non-universal constants that match the lattice data to the universal scaling functions can be utilized to extract information on the phase diagram in Sec. 3. One important class of observables in which critical behavior can be made manifest are the before mentioned cumulants of charge fluctuations. These can, however, also be used to demonstrate the liberation of quarks from the hadronic bound states for $T>T_{c}$, which is done in Sec. 4. And finally in Sec. 5, we discuss several phenomenological results that can be obtained from a comparison of experimental cumulants of conserved charge fluctuations with the corresponding lattice results.

2 QCD bulk thermodynamics on the lattice

2.1 The path integral on the lattice

The fundamental quantity to calculate thermal expectation values is the partition function. As we are dealing with a quantum field theory where particle pair creation and annihilation is at the heart of things, we will perform our simulations in the grand canonical ensemble. The partition function $Z$ is obtained in the path integral formalism as functional integral over all fundamental fields, weighted by the exponential of the negative Euclidean action functional. It is important to notice that thermal averages are calculated in Euclidean space, which is obtained after Wick rotation $t\to-it=\tau$. The temperature is introduced through a compactification of the $\tau$ direction, in analogy to classical statistical physics. The Euclidean space-time is discretized, thus coordinates can take values on a four dimensional hyper-cubic lattice, $x_{\nu}=an_{\nu}$, with $n_{\nu}\in\{1,2,\dots,N_{\nu}\}$, where $\nu\in\{\tau,x,y,z\}$ and $a$ denotes the dimensionful lattice constant. Here we restrict ourselves to isentropic lattices where $a$ is not depending on the direction $\nu$. Furthermore, we consider cubic volumes with $N_{x}=N_{y}=N_{z}\equiv N_{\sigma}$ and periodic boundary conditions. For a given lattice $N_{\tau}\times N_{\sigma}^{3}$, temperature and spatial volume are obtained as

[TABLE]

In lattice QCD, the fundamental fields are the Grassmannian quark $\psi^{(f)}_{n}$ and anti-quark fields $\bar{\psi}^{(f)}_{n}$ for each quark flavor ($f$), as well as the link fields $U_{n,\nu}$ that are associated with the links between neighboring sits. They describe a parallel transport of the continuum gluon field $A_{\nu}(x)$ from site $x$ to $x+a\hat{\nu}$

[TABLE]

where $\hat{\nu}$ is the unit vector in $\nu$-direction and $P$ denotes path ordering. The Link fields take values in the $SU(3)$ color group. Similarly, the Dirac 4-spinors $\psi^{(f)}_{n}$ and $\bar{\psi}^{(f)}_{n}$ are color 3-vectors. The spinor and color indices have been suppressed here. The number of the quark flavors $N_{f}$, as well as the corresponding quark masses $m_{f}$, with $f\in\{1,2,\dots,N_{f}\}$ are parameters of QCD. Here we will – if not stated differently – discuss (2+1)-flavor QCD, i.e. we will have two flavor of mass degenerate light quark ($u,d)$, with masses $m_{u}=m_{d}$ and one heavier strange quark $s$, with mass $m_{s}$. We also introduce one quark chemical potentials for each quark flavor $\vec{\mu}=(\mu_{u},\mu_{d},\mu_{s})$. Finally, we arrive at an expression for the partition function, which is given as

[TABLE]

where the measure ${\rm d}U_{n,\nu}$ denotes the Haar measure on the compact gauge group $SU(3)$. Note that on the lattice the path integral is not a functional integral any more but a mathematically well defined multi-dimensional integral.

The Euclidean action $S_{E}$ is composed of a fermionic and a gluonic part $S_{E}=S_{F}+S_{G}$. The gluon action is to leading order in $a$ given by Wilson’s plaquette action [13], where a plaquette variable $U_{n,\nu\rho}$ is defined as the product of link variables along an elementary plaquette $U_{n,\nu\rho}\equiv U_{n,\nu}U_{n+\hat{\nu},\rho}U^{\dagger}_{n+\hat{\rho},\nu}U^{\dagger}_{n,\rho}$. With this definition the gluon action is given as

[TABLE]

The fermionic part has the structure of a quadratic form

[TABLE]

with a fermion matrix $M^{(f)}_{n,m}$. Different forms of this matrix define different fermion discretizations on the lattice. Without giving a specific form of the fermion matrix, we will discuss properties of different lattice fermions in Sec. 2.2.

The Grassmannian quark fields $\psi^{(f)}_{n}$ and $\bar{\psi}^{(f)}_{n}$ may be integrated out analytically from Eq. (3), which is most desirable since handling Grassmannians on the computer is impractical. The partition function may then be written as

[TABLE]

featuring the fermion determinant ${\rm det}\left[M^{(f)}(\mu_{f})\right]$. Here we have made explicit that the chemical potential enters only in the fermion part of the action, i.e. in the fermion matrix. It is introduced by a modification of the temporal links as $U_{n,0}\to e^{a\mu_{f}}U_{n,0}$ and $U^{\dagger}_{n,0}\to e^{-a\mu_{f}}U^{\dagger}_{n,0}$ [24]. Other possibilities to introduce the chemical potential have also been considered [25, 26].

The connection to bulk thermodynamics is done by standard thermodynamic relations among the thermodynamic potentials. Here we employ the grand canonical potential $\Omega(T,V,\vec{\mu})\equiv-T{\rm ln}Z(T,V,\vec{\mu})$ from which the pressure can be obtained most easily as $p(T,V,\vec{\mu})=-\Omega(T,V,\vec{\mu})/V$. However, any other thermal expectation value $\left<\mathcal{O}\right>(T,V,\vec{\mu})$ of an observable $\mathcal{O}$ can be obtained by

[TABLE]

In fact, since the partition function itself is hardly accessible from MC simulations, it is common practice to study derivatives of $\ln Z$, taken with respect to various external parameter. Such derivatives are then obtained through the calculation of thermal expectation values, e.g. we have $\partial\ln Z/\partial m_{u}=<Tr[(M^{(u)})^{-1}]>$. This strategy is also followed for the calculation of the energy density and pressure (see Sec. 2.3). Equally important to infer on the properties of the QGP is the study of correlation-functions, resulting in the evaluation of observables of the type $<\mathcal{O}_{x}\mathcal{O}_{y}>$. Although seemingly very difficult, it is possible to reconstruct spectral-functions from those Euclidean correlators by either Bayesian or variational methods. The review of of correlation- and spectral-functions are, however, beyond the scope of this article.

2.2 Lattice discretization schemes and the continuum limit

The discretization that has been performed in order to perform the MC Integration Eq. (6) bears certain difficulties. The most obvious is certainly the introduction of discretization errors of the action functionals. The gauge action as given in Eq. (4) receives corrections which are of order $\mathcal{O}(a^{2})$. These discretization errors can be systematically removed by adding additional terms to the gauge action that are composed of Wilson loops of length greater than that of the plaquette [27, 28]. In Fact, it is common to add in a suitable linear combination the sum of all Wilson loops of rectangular shape. Such a Symanzik improved gauge action is correct up to correction of order $\mathcal{O}(a^{4})$.

With in the fermionic sector the problems are more severe. Here the naive nearest-neighbor discretization of the covariant derivative in the Dirac operator not only introduces discretization errors of order $\mathcal{O}(a)$ it also gives rise to unphysical degrees of freedom. This is due to the fact that a finite lattice spacing introduces a cut-off in momentum space and thus a periodicity in the lattice dispersion relation. The fermion modes that are associated with the roots of the dispersion relation in the $2^{d}$ corners of the Brillouin zone are known as fermion doublers. Moreover, the fermion doublers are intimately connected with the explicit breaking of the chiral symmetry on the lattice [29]. Several discretization schemes have been invented to deal with these problems. Among them are Wilson [13] and staggered [30] fermions which are most commonly used due to their moderate computational costs, but also Domain Wall [31] and Overlap [32, 33] fermions.

In lattice calculations at nonzero temperature, the consequences of the chosen fermion discretization scheme are three-fold. For a study of universal critical behavior that arises due to the spontaneous breaking of the chiral symmetry of QCD, it is most important to use an action that possess at least a remnant of this symmetry. The situation here is the following: The Wilson action breaks the chiral symmetry completely at any finite lattice spacing. The staggered action preserves a $U(1)_{\rm even}\times U(1)_{\rm odd}$ subgroup of the symmetry, i.e. the fermion fields on even and odd sides can be multiplied separately by a phase. Domain-Wall fermions are chiral in the limit of a large fifth dimension and Overlap fermions posses an exact (albeit modified) chiral symmetry at any finite lattice spacing. The latter two action are so costly that so far no scaling study close to the chiral limit has been possible. First studies of the transition temperature [34], thermodynamic quantities [35, 36] and the anomalously broken axial symmetry $U_{\rm A}(1)$ and the Dirac spectrum [37, 38, 39, 40, 41, 42] have, however, been performed recently.

For studies of bulk thermodynamic quantities such as the pressure $p$ and the energy density $\epsilon$ the high and the low temperature regimes suffer from two different types of discretization errors. At very high temperatures, deviations from the Stefan-Boltzmann (SB) limit stem from the discretization errors of the covariant derivative. Just as in the gluonic case, such errors can be eliminated systematically order by order using so called improved actions, which are widely used in case of Wilson and staggered fermions. For Wilson fermions, the clover-term [43] reduces the discretization errors of the action to $\mathcal{O}(a^{2})$, for staggered fermions, straight (Naik [44]) or bended (p4 [45]) three-link terms diminish the error to $\mathcal{O}(a^{4})$. It is interesting to note that even though the standard Wilson action features cutoff errors of $\mathcal{O}(a)$, corrections to bulk thermodynamic quantities are universal to standard Wilson and standard staggered actions and start at $\mathcal{O}(a^{2})$. To summarize, an expansion of the dimensionless quantity $p/T^{4}$ around the ideal gas limit yields [46]

[TABLE]

At low temperatures, bulk thermodynamic quantities are, as we will see, mainly determined by the hadronic spectrum of QCD. Due to the existence of fermion doublers on the lattice – also called fermion tastes – the hadronic states are duplicated as well. Moreover, high energy gluons that scatter fermions from one corner of the Brillouin zone to another, introduce a nontrivial interaction between different fermion tastes, which leads to the so called taste splitting of hadronic states. This splitting is most prominent in the light pion sector. Within the staggered fermions formulations, the doublers are reduced from 16 to 4 tastes. Consequently, with these 4 tastes many more (15) pions can be formed, which become degenerated only in the continuum limit. In order to suppress large momentum modes of the gluons and thus reduce the taste splitting, smeared gauge fields are used in the Dirac operator. Popular staggered type actions are the 2stout and 4stout actions, that use 2 and 4 levels of stout-smearing [47], respectively as well as the HISQ (highly improved staggered quark) action [48]. The latter combines two levels of 5- and 7-link smearing with the Naik-term that also improves the high-T limit of this action. In terms of taste splitting, which can be measures, e.g., by the root-mean-square mass $M_{\pi}^{\rm RMS}$ of the pions, it is found that the 4stout actions performs similarly well as the HISQ action.

Having understood the leading correction in the lattice spacing $a$, we eventually extrapolate to the continuum limit ($a=0$) by using calculations performed at several lattice spacings and fitting appropriate polynomials to it. This extrapolation is done at fixed volume $V$ and temperature $T$. The lattice spacing is controlled by the lattice coupling $\beta=6/g^{2}$. Given the relation Eq. (1), for a lattice with a fixed number of points in temporal direction $N_{\tau}$, a temperature scan is performed by varying $\beta$. At small values of the bare coupling $g$ (large $\beta$) the relation between $\beta$ and $a$ and thus $\beta$ and $T$ is given by the well known perturbative $\beta$-function of QCD. In practice, at temperatures close to the transition, the coupling is, however, not small and the scale has to be set by using an experimentally measured input observable, which is compared to the same quantity measured on the lattice, in units of the lattice spacing $a$. As this input observable is almost always a hadronic quantity, experimentally determined in vacuum, the scale setting is performed at $T=0$, i.e. on lattices with $N_{\tau}\leq N_{\sigma}$. In principle one is completely free to chose the input observable, as in the continuum limit all scale ambiguities should vanish. For recent thermodynamic calculations, quantities related to the static quark potential ($r_{0}$ [49], $r_{1}$ [50]), the kaon decay constant ($f_{K}$) or quantities related to the Wilson flow ($t_{0}^{1/2}$ [51], $w_{0}$ [52]) have been used for scale setting. It is worth noting that the continuum extrapolation is done with fixed renormalized quark masses. This requires a non trivial tuning of the bare quark masses as function of $\beta$. At first, the strange quark mass $m_{s}$ is fixed to its physical value by demanding that the mass of the $\eta_{s\bar{s}}$ meson, which is related to the pion and kaon by chiral perturbation theory, is given by $M_{\eta_{s\bar{s}}}=\sqrt{(2M_{K}^{2}-M_{\pi}^{2})}=686$ MeV. Having fixed $m_{s}$, we fix the light quark mass by holding the ration $m_{l}/m_{s}$ constant. Some of the results we will present have been obtained with $m_{l}=ms/20$, although the physical value is $m_{l}=ms/27$. The relation between the coupling $\beta$ and the bare quark masses is called line of constant physics (LCP).

2.3 The equation of state

The trace of the energy-momentum tensor, also called trace anomaly or the interaction measure ($\Theta^{\mu\mu}$), is related to the pressure $p$ by a temperature derivative (see, e.g. Ref. [53]) but can also be written as a total derivative of $\ln Z$ with respect to the lattice spacing

[TABLE]

where $\epsilon$ is the energy density. It is straight forward to evaluate the right-hand side of the last equality on the lattice, thus $\Theta^{\mu\mu}(T)$ is the fundamental observable. The pressure $p$ can be calculated by the integral method, i.e. by integrating $\Theta^{\mu\mu}(T)/T^{5}$, up to an integration constant. This also allows to reconstruct the energy density $\epsilon$ as well as the entropy density $s/T^{3}=(\epsilon+p)/T^{4}$. Recently continuum extrapolated results have been reported [54, 55], using stout and HISQ actions, respectively. The calculations are in good agreement, although the light quark masses are slightly different ($m_{l}=ms/27$ [54] vs. $m_{l}=ms/20$ [55]). In Fig. 2 we show the HISQ results [55].

The rise of the pressure, energy density, and entropy density close to the transition region of $T_{c}=159\pm 9$ MeV [56] is a consequence of the liberation of many new quark degrees of freedom. It is interesting to note that the critical energy density, associated with the transition region $\epsilon_{c}=0.34\pm 0.16$ $\mbox{GeV/fm}^{3}$ is only slightly larger than the critical energy density of ordinary nuclear matter $\epsilon^{\mbox{nuclear matter}}\approx 0.15$ $\mbox{GeV/fm}^{3}$. In fact, $\epsilon_{c}$ is close to the energy density reached in the dense packing limit of nucleons, assuming a nucleon radius of $R_{N}=0.8$ fm.

Also shown in Fig. 2 are dark solid lines that indicate the hadron resonance gas results at low temperature. The hadron resonance gas is a non interacting gas of all known hadronic bound states and resonances. It is eminent from Fig. 2 that at low temperatures the QCD results come very close the HRG. A more detailed comparison of the HRG and QCD results of the pressure, is given in section Sec. 4. Indeed it is found that a departure of the QCD pressure, and thus a fade of hadronic degrees of freedom starts just after the transition regions. On the other hand, the bulk thermodynamic quantities obtained from HISQ fermions [55] comes close to the hard thermal loop perturbative calculation [57] and the dimensional reduced QCD (EQCD) calculation [58] at temperatures $T\gtrsim 400$ MeV. At such temperatures all three calculations differ on the 10% level (for a comparison see Ref. [55]). We conclude that the temperature region $T_{c}\leq 2T_{c}$ is deconfined but highly nonperturbative, i.e. the QGP is strongly coupled in this temperature range.

Of high interest is also the specific heat, $C_{V}={\rm d}\epsilon/{\rm d}T$, and speed of sound square, $c_{s}^{2}={\rm d}p/{\rm d}\epsilon$, which are also shown in Fig, 2. In a hydrodynamic description of the fireball that is generated in an heavy ion collision, these quantities are, e.g., directly related to the life time of the fireball, and thus to the question of whether the plasma has enough time to equilibrate or not. We find that the specific head does not show a peak although in the chiral limit ($m_{l}\to 0$) it should develop a cusp, i.e. critical behavior dictates $C_{V}\sim\left(\frac{T-T_{c}}{T_{C}}\right)^{-\alpha}+C_{V}^{\rm reg}$, with $\alpha\approx-0.2$. The regular term is thus dominating.

3 Universal critical behavior and the QCD phase diagram

The chiral symmetry of the QCD Lagrangian plays an important role for the determination of the QCD phase structure. In the limit of vanishing quark masses this symmetry is $SU_{L}(N_{f})\times SU_{R}(N_{f})\times U_{A}(1)\times U_{V}(1)$, which implies that left handed and right handed spinors can be rotated separately. The latter $U_{V}(1)$ symmetry represents the baryon number conservation. The $U_{A}(1)$ is anomalously broken by quantum fluctuations. Nonzero quark masses break the $SU_{L}(N_{f})\times SU_{R}(N_{f})$ symmetry explicitly, however, more importantly it is also broken spontaneously at low temperatures. This gives, e.g. rise to three (almost) massless pseudo-scalar Goldstone bosons, that are identified with the pions.

The nature of the chiral transition crucially depends on the quark masses and the number of flavors. In Fig. 3 (left), a sketch of the nature of the QCD phase diagram with two mass degenerate light flavors (up and down) and one heavier strange quark at $\mu_{B}=0$ is show.

For the bulk part of the diagram, the QCD transition at $\mu_{B}=0$ is not a true phase transition, but a rapid continuous and analytic crossover. This includes also the physical mass point [59]. In the limit of infinite quark masses ($N_{f}=0$, pure gauge theory), which is realized in the upper right corner of the diagram we have a region of first order transitions. The same is true for the lower left corner, the massless limit of three degenerate quark masses ($N_{f}=3$). On the boundary of both of these regions the transition becomes second order, with an universality class belonging to the 3d-Z(2) (Ising) model. For the lower left corner this has been demonstrated with std. staggered [60, 61] and Wilson-clover fermions [62, 63]. Simulations with improved staggered fermions have so far not seen a first order region. Only an upper limit for the degenerate quark mass $m_{u}=m_{d}=m_{s}<m_{s}^{\rm pyhs.}/270$ could be given, where $m_{s}^{\rm pyhs.}$ [64, 65, 66] denotes the physical value of the strange quark mass.

The situation in the chiral limit of two degenerate flavor is yet to be clarified. It is dependent on the strength of the axial anomaly ($U_{A}(1)$ breaking) close to the transition temperature $T_{c}$. If the $U_{A}(1)$ breaking is still significant at $T_{c}$, the transition is of second order and since $SU(2)\times SU(2)$ is isomorphic to $O(4)$, it is in the universality class of the 3dim. $O(4)$ model [9]. However, if the anomaly is already sufficiently restored at $T_{c}$, the transition is either first or second order with a different symmetry breaking pattern [67, 68]. For recent lattice studies of the $U_{A}(1)$ anomaly see, e.g. [37, 38, 39, 40, 41, 42]. The situation is expected to carry over to the chiral limit of the (2+1)-flavor theory (holding $m_{s}$ fixed and sending $m_{u,d}\to 0$), as long as the strange quark mass is lager than a tri-critical strange quark mass $m_{s}^{\rm tri}$ where we enter the first order region in the lower left corner of the diagram. Unfortunately, the relative location of $m_{s}^{\rm tri}$ with respect to the physical strange quark mass $m_{s}^{\rm pyhs.}$ is unclear. It is thus not known whether $m_{s}^{\rm tri}<m_{s}^{\rm pyhs}$ or $m_{s}^{\rm tri}>m_{s}^{\rm pyhs}$ is realized in nature, which means that in the chiral limit of a (2+1)-flavor theory with a fixed strange quark mass we could either find a first or second order transition. The former would give rise to a critical $Z(2)$ point at finite quark masses. Even the possibility of $m_{s}^{\rm tri}\to\infty$ is not excluded. It was recently discussed and supported by lattice simulations at imaginary chemical potential with std. staggered fermions on coarse lattices [69]. In the following discussion we will, however, assume $m_{s}^{\rm tri}<m_{s}^{\rm pyhs.}$, which is consistent with the universal scaling observed with improved staggered fermions [70, 71, 72].

The extension of this diagram to a nonzero baryon chemical potential $\mu_{B}$ is also controversially discussed. Especially interesting is whether the first order region in the lower left corner is shrinking or growing with increasing $\mu_{B}$. One possible extension is shown in the right-hand side of Fig. 3, where we show the $(T-\mu_{B})$-phase diagram of the (2+1)-flavor theory with $m_{s}=m^{\rm phys}$, in dependence of the light quark mass $m_{u,d}$. The shown scenario assumes that $m_{s}^{tri}$ grows with $\mu_{B}$, implying that the second order transition in the chiral limit turns into a first order transition at some value of the chemical potential $\mu_{B}=\mu_{B}^{tri}$. A possible $Z(2)$-critical endpoint at physical quark masses would then be smoothly connected with the tri-critical point in the chiral limit by a line of second order transition. A different scenario of a shrinking first order region is supported by lattice simulations at imaginary chemical potential with std. staggered fermions on coarse lattices [73]. In this scenario, it is unlikely that a possible end-point at physical quark masses can be smoothly connected to the critical surface of $Z(2)$ end-points in the low mass limit of the three flavor theory.

Direct searches for the critical end-point on the lattice are so far not very conclusive. First studies that reported a positive result [74, 75] used a finite volume scaling analysis of the Lee-Yang zeros together we a reweighting approach. It is, however not clear whether this method could produce fake signals [76]. Later on, a number of positive results have been reported from calculations with the canonical partition function [77, 78, 80, 79], which suffer from small volumes, large cut-off effects and large quark masses. More recent calculations are based on estimating the radius of convergence of the Taylor expansion in $\mu_{B}$ of some suitable observables [81]. The problems here are twofold, first, only a few expansion coefficients are available which makes the estimate of the radius of convergence difficult. Second, the first few expansion coefficients are in agreement with calculations from the hadron resonance gas model for $T<Tc$. This models exhibits an infinite convergence radius and thus no critical point. Nevertheless two groups report estimates for the convergence radius [82, 83].

3.1 Universal scaling in (2+1)-flavor QCD

According the the scaling hypothesis, the free energy has a singular part that is responsible for the power laws that thermodynamic response functions exhibits near the critical point. For two degenerate light flavors of mass $m_{l}\equiv m_{u}=m_{d}$ and one heavier strange quark mass $m_{s}=m_{s}^{phys}$, we thus make the Ansatz

[TABLE]

where the singular part $f_{s}$ will become a generalized homogeneous function of its arguments once the correct scaling fields have been chosen as its natural variables. In general, all $O(N)$ and $Z(N)$ models exhibit two relevant scaling fields, the temperature-like ($t$) and external field-like scaling fields ($h$). Considering a critical point in the chiral limit ($m_{l}=0$) and at zero chemical potentials ($\vec{\mu}=0$) the leading order dependence of the scaling fields on $T,m_{l},\vec{\mu}$ are determined by symmetry arguments, i.e. to leading order $h$ depends only on couplings that break chiral symmetry in the light quark sector, while $t$ depends on all other couplings. We find

[TABLE]

where $T_{c}$ is the phase transition temperature and $t_{0},h_{0}$ are non-universal scale parameters. Here we have adopted the $\vec{\mu}=(\mu_{q},\mu_{I},\mu_{s})$ basis and assume $\mu_{I}=0$ as a finite isospin chemical potential $\mu_{I}$ would break the flavor symmetry in such a way that a different symmetry breaking pattern would arise. With the above scaling fields we can exploit the property of homogeneity of the function $f_{s}$ to obtain

[TABLE]

The function $f_{f}$ is a universal scaling function that depend on a single scaling variable $z$. In this way we have singled out the leading order singular behavior, but have neglected sub-leading terms that are produced by irrelevant scaling variables and are known as corrections to scaling. From here we can easily obtain the magnetic equation of state $M=h^{1/\delta}f_{G}(z)$ by taking a derivative with respect to $H=h_{0}h$, where the scaling function $f_{G}$ is connected to $f_{f}$ by

[TABLE]

Note that this scaling form is not the famous Widom-Griffiths from, as the scaling variables do not depend on the Magnetization. It is, however, well suited for a comparison with QCD where the calculated magnetization, which can, e.g., be defined as $M_{b}=m_{s}\left<\bar{\psi}\psi\right>/T^{4}$ has statistical and systematical errors. As an example of a scaling fit to the magnetization $M_{b}$ see the upper part of Fig. 4.

The lattice data was obtained with the HISQ action on $N_{\tau}=8$ lattices and quark mass ratios $m_{l}/m_{s}=1/20,1/27,1/40$ [72]. The fit parameter are $t_{0}$, $h_{0}$ and $T_{c}^{0}$ and three additional parameter for the regular part. The right panel shows the data collapse obtain when proper scaling variables are used. It can be seen that the scaling Ansatz works fine for quark masses as large as $m_{l}/m_{s}=1/20$. The relative contribution of the regular part has been found to be small.

A more detailed scaling analysis was performed in Ref. [72], using a joint fit to the chiral condensate and chiral susceptibility. As a result we show in the lower left panel of Fig. 4 the continuum extrapolation of the transition temperature at the physical mass point, performed jointly with two different lattice discretization schemes, HISQ and asqtad. The obtained value is $T_{c}=154\pm 9$ MeV. Other studies with improved staggered fermions, that use also other criteria for the definition of the crossover temperature, which are not directly related to universal critical behavior obtain similar values [84, 85, 86]. Also a study using chiral Domain Wall Fermions [34] yields a compatible value. For the value for the phase transition temperature in the chiral limit, we obtain values of the oder of $T_{c}^{0}\approx 145$ MeV.

Given the parameters $h_{0},t_{0}$, we are able to determine the non universal parameter $\kappa_{l}$ as well, which is the curvature of the phase transition temperature in the chiral limit. The scaled susceptibility $m_{s}t_{0}h^{(1-\beta)/\beta\delta}\chi_{m,\mu}/T^{2}$ can be fitted to the scaling function $f_{G}^{\prime}(z)$, with $\kappa_{l}$ being the only parameter left, which enters as a normalization factor. We obtain $\kappa_{q}=0.059(5)$ [71], which translates into a curvature of $\kappa_{B}=0.0066(5)$, associated with a change of the baryon chemical potential $\mu_{B}$. This result is not yet continuum extrapolated, but ultimately linked the critical behavior in the chiral limit. Other results based on imaginary chemical potential or a Taylor expansion of chiral and deconfining order parameter are in the range $0.007<\kappa_{B}<0.02$ [87, 88, 89].

4 The degrees of freedom in QCD near deconfinement

An important question in QCD is to understand what happens to the microscopic degrees of freedom at deconfinement. For the light quark sector in QCD, lattice studies point to the fact that colored degrees of freedom start to appear already when the chiral symmetry is approximately restored [56]. Since the strange quark mass is $\sim\Lambda_{QCD}$, it is not a priori obvious whether the same phenomena would be true even for the strangeness sector. Also of particular interest are the bound states of charm quark. It is now known from the first lattice studies [90, 91] that $c\bar{c}$ bound states like the $J/\psi$ may survive in the hot QCD medium upto $1.4~{}T_{c}$. The interesting problems to understand non-perturbatively are what is the fate of hadrons consisting of charm and light quarks near $T_{c}$, and what are the dominant degrees of freedom carrying charm in the deconfined medium.

4.1 Fate of strange hadrons

Below $T_{c}$, thermodynamic quantities like the pressure in QCD can be fairly well described as consisting of non-interacting hadrons and resonances. The total pressure of a thermal ensemble of a non-interacting gas of mesons, baryons and resonances containing strangeness quantum number $S$ at finite chemical potentials $\hat{\mu}_{X}$ $=$ $\mu_{X}/T$, $X\equiv B,S$ is given as,

[TABLE]

where $P_{M}$ and $P_{B,S=n}$ denote the partial pressure of open-strange mesons and respectively baryons containing strangeness quantum number $n$ at vanishing chemical potentials. The derivatives of $P(\hat{\mu}_{S},\hat{\mu}_{B})$ with respect to $\mu_{X}$, $\chi_{kl}^{BS}$ $=$ $\left.\frac{\partial^{(k+l)}[P(\hat{\mu}_{B},\hat{\mu}_{S})/T^{4}]}{\partial\hat{\mu}_{B}^{k}\partial\hat{\mu}_{S}^{l}}\right|_{\vec{\mu}=0}$, define different generalized susceptibilities which encode the fluctuations of $B,S$ and their correlations at vanishingly small chemical potentials. In a non-interacting gas of strange hadrons, these quantities can be written as,

[TABLE]

Inverting above relations would then give us the partial pressures in the meson and baryon sectors, $P_{M},P_{B}$ at $\mu_{X}=0$ [92]. The crucial input for this calculation are the susceptibilities $\chi_{kl}^{BS}$ measured on the lattice. If one measures susceptibilities upto fourth order, then in the $B$-$S$ sector one has $6$ independent measurements $\chi_{2}^{S},\chi_{4}^{S},\chi_{11}^{BS},\chi_{22}^{BS},\chi_{31}^{BS},\chi_{13}^{BS}$. If indeed the HRG model describes the data, then one could extract the $4$ independent parameters of the model which are $P_{M}$ and three $P_{B}$’s in terms of these susceptibilities, along with $2$ independent constraint equations. The lattice QCD data should satisfy the constraints imposed on the model, in order to validate the model. In terms of the susceptibilities, the $2$ independent constraints have been derived in Ref. [92],

[TABLE]

The values of $v_{1}$ and $v_{2}$ based on measurements on $N_{\tau}=6,8$ lattices with obtained with the HISQ action from Ref. [92] are summarized in Fig. 5. It can be concluded [92] that the HRG description for strange hadrons starts to break down at the chiral crossover transition region. In fact, the ratio of the partial pressure of strange baryons to the total HRG pressure given by the quantity $-\chi_{11}^{BS}/\chi_{2}^{S}$ shown in the right panel of Fig. 5 from [93] start to show deviation from the HRG description already at $T_{c}$, implying deconfinement of strangeness degrees of freedom. Similar consequences for strange hadrons have been observed independently from a different set of observables in Ref. [94] motivated from the above outlined analysis.

4.2 Melting of open charm hadrons

For the charm sector, the ground states of the baryons carrying charm quantum numbers $C=2,3$ respectively, are heavier than the $C=1$ baryons by about $1$ GeV, and hence their partial pressure contributions are suppressed near $T_{c}$. As a consequence the total pressure of the charm degrees of freedom in the hadronic medium is simply due to the sum of the partial pressures of charmed mesons $P_{M}$ and $C=1$ baryons $P_{B,C=1}$ respectively [95]. The $P_{M}$ can be expressed as two equivalent combinations of the susceptibilities, $P_{M}=\chi_{2}^{C}-\chi_{22}^{BC}=\chi_{4}^{C}-\chi_{13}^{BC}$. For the heavy quarks, additional complication arise due to strong mass dependent cut-off effects in the susceptibilities. Instead, it has been proposed in Ref. [95] that the ratio of the two independent measures of $P_{M}$ would be a good observable with significantly less contamination from the cut-off effects of order $\mathcal{O}(ma)$. The ratio should be unity in the hadron phase and any deviation from unity would be a signal for the melting of these charmed mesons. In the baryon sector, the relative contribution of the partial pressures of $C=2,3$ baryons given by $\chi_{4}^{C}-\chi_{2}^{C}$ has been indeed shown to be vanishingly small compared to $P_{B,C=1}$  [95] near $T_{c}$. This considerably simplifies the partial pressure due to the baryons carrying a unit charm to $P_{B=1}\simeq\chi_{mn}^{BC}~{},m,n>0$ and $m+n=$ even. Hence, ratios like $\chi_{13}^{BC}/\chi_{22}^{BC}$ would be unity in a phase consisting of charmed baryons. The partial pressures for the charmed meson and baryon sectors from Ref. [95] are summarized in Fig. 6. It is evident that both open charm meson and baryon partial pressure ratios depart from unity already in the crossover region signaling the deconfinement of charm at $T_{c}$. The conclusions based on these ratios are quite robust since these are independent of the details of the hadron spectrum and have a very mild lattice spacing dependence seen from the $N_{\tau}$ independence of the results in Fig. 6.

4.3 Evidence of additional charm hadrons from QCD thermodynamics

Since these newly proposed combinations of generalized susceptibility, proposed for the first time in Refs. [92, 95] and discussed in Sec. 4.1 and 4.2, allow for the disentaglement of the partial pressures of meson and baryon sectors with different quantum numbers, we can further understand the microscopic composition of individual sectors at the chiral crossover transition. Calculations in the Quark Model(QM) [96, 97] as well hadron spectrum studies on the lattice [98, 99] predict many more charmed baryons than experimentally detected and tabulated in the Particle Data Table [100], whereas for the charmed meson sector there is quite a good agreement [101, 102]. For instance, excited states for the $D,D_{s}$ mesons which are upto $700$ MeV above the ground state are fairly well measured experimentally. On the other hand, the spectrum of $\Lambda_{c}$ baryon is not well measured [95] as shown in the left panel of Fig. 7.

Only three states above the ground state have been measured in the experiments and assigned to a well defined spin and parity. However, the ground state spin is not yet measured with certainty. There are many more bound states in different spin-parity channels predicted from the QM. Similar results are obtained from charm baryon spectroscopy on the lattice [99]. The presence of these additional baryon states and resonances have been shown to have a measurable impact on the pressure of open-charm hadrons near $T_{c}$  [95]. The baryon states which have masses $1$ GeV or more than the ground state lying above the blue band in Fig. 7 may be unstable or would have negligible contribution to the partial pressures. For a systematic understanding of the relative contributions of these states by comparing to the lattice data, three different sets have been proposed [95]. These consist of a hadron resonance gas(HRG) with all predicted states in the Quark Model (QM-HRG)  [96, 97] denoted by brown lines, a QM-HRG-3 consisting of Quark Model states upto $3$ GeV mass and a PDG-HRG model consisting only of the experimentally known states. Using appropriate combinations of susceptibilities, the partial pressures of these channels with specific quantum numbers have been in derived in Ref. [95]. For the charged sector, charmed baryon and meson partial pressures are given as $P_{B}\simeq\chi_{112}^{BQC}$ and $P_{M}\simeq\chi_{13}^{QC}-\chi_{112}^{BQC}$ respectively whereas for the charm baryons additionally carrying strangeness, $P_{B}\simeq\chi_{112}^{BSC}$ and $P_{M}\simeq\chi_{13}^{SC}-\chi_{112}^{BSC}$. The comparative results between lattice data and predictions from different HRG models are summarized from Ref. [95] in the right panel of Fig. 7. In all these sectors a significant contribution of the additional baryon states over PDG-HRG at the crossover transition has been found. This is a remarkable example where finite temperature lattice QCD data provides an independent motivation for further experimental and theoretical studies on charm hadron spectroscopy [99].

4.4 Nature of quasi-particles beyond $T_{c}$

Evidence from the results for screening masses in the charm sector [103] show that the baryon and meson-like excitations do not immediately disappear at the chiral crossover transition. The screening masses only reaches the free fermion limit at $T>250$ MeV. Moreover comparison of the lattice data for off-diagonal susceptibilities like $\chi_{22}^{uc}$ normalized by $\chi_{2}^{c}$, to perturbative results within Hard thermal Loop resummation, show that these are in fairly good agreement only beyond $1.2~{}T_{c}$ [104]. This further hints to the fact that the dominant degrees of freedom carrying charm quantum number changes at around $1.2~{}T_{c}$ where quark degrees of freedom starts becoming dominant in the plasma. Owing to the large charm quark mass, the charm-carrying degrees of freedom could be well described in terms of quasiparticles. It allows one to further decompose the total pressure due to charm degrees of freedom into contributions from baryons, mesons and quark-like excitations carrying unit charm quantum number, following Ref. [104],

[TABLE]

Again using the lattice data for charm correlations and fluctuations upto fourth order, the partial pressures $P_{M}$, $P_{B,C=1}$ and $P_{Q}$ can be solved, as outlined in Ref.  [104] in terms of six observables $\chi_{2}^{C}$, $\chi_{11}^{BC}$, $\chi_{13}^{BC}$ $\chi_{22}^{BC},\chi_{31}^{BC},\chi_{4}^{C}$. The trivial constraints in the charm sector $\chi_{4}^{C}=\chi_{2}^{C}~{},~{}\chi_{11}^{BC}=\chi_{13}^{BC}~{}$, leaves only four independent observables. However, since three independent partial pressures need to be calculated within this model there would only be a unique non-trivial constraint for this model $c_{1}\equiv\chi_{13}^{BC}-4\chi_{22}^{BC}+3\chi_{31}^{BC}=0$. The contribution of the three different partial pressures within this model has been calculated [104]and summarized in Fig. 8. It is was concluded that the charm meson and baryon-like excitations contribute significantly to the total pressure till $T\sim 200$ MeV with the quark-like excitations dominating only for temperatures beyond. Moreover a specific sub-sector of this model was looked into detail [104] consisting of quasi-particles which carry both charm and strangeness quantum numbers, since these two quantum numbers are completely decorrelated when quark-like excitations play a dominant role. Indeed it was found that the meson and baryon-like excitations within this sector may survive well beyond $T_{c}$, with the meson-like contributions dominating [104]. The quark-like excitations accounts for more than $50\%$ of the total pressure only around $200$ MeV. These studies also looked for evidences for the presence of di-quarks in the deconfined phase of QCD. Atleast for the charm sector, the presence of di-quark like excitations has not been found above $T_{c}$ [104].

To check whether indeed a model can “mimic” QCD, the constraints that define the model should be exactly satisfied by the data from lattice QCD calculations. It was infact found in Ref. [104] that the best available lattice data satisfies constraints for all the above mentioned models, which is shown in the right panel of Fig. 8. So in summary, there is a very strong evidence atleast for the charm sector, that the quark-like degrees of freedom start appearing already near $T_{c}$. However, these remain strongly correlated as color neutral baryon and meson-like excitations till about $1.2~{}T_{c}$. For the first time it was shown [104] that colored quark-like excitations only start dominating thermodynamically beyond this temperature. This has important consequences for the phenomenological modeling of heavy-quark diffusion and quenching which needs the deconfinement temperature of the heavy-quarks and the nature and abundance of resonances that carry charm quantum number in QGP as crucial inputs [105].

5 Freeze-out and connection to heavy ion experiments

Many new insights on the nature of the quasi-particles in the hot QCD medium and their interactions have been recently obtained from the extensive lattice data on fluctuation measurements especially in the charm sector [104]. As motivated in the introductory sections, one of the outstanding issues is to understand when does the HRG picture breaks down for QCD with light quarks. This is not only important for our fundamental understanding of QCD, but it also allows us to relate to the wealth of experimental data coming from the heavy-ion experiments at Relativistic Heavy Ion Collider (RHIC), BNL, and at CERN. In the experiments, the fireball created due to the collision of two heavy nuclei rapidly expands and cools down until hadronization occur and many aspects of this dynamics are now fairly well understood in terms of the hydrodynamic modes. This hints to the existence of a local thermally equilibrated plasma and the lattice data can accurately describe the thermal QCD contribution inherent in the experimental data. The question is at what temperature and chemical potential can one compare the lattice data on fluctuations of conserved charges to those measured in the experiments. The traditional method is to determine the temperature and baryon chemical potential of the fireball at the time of chemical freeze-out i.e., when all the inelastic scatterings cease to exist between the hadrons. This involve statistical fits to the ratios of hadron abundances measured in the experiments [108], inspired from the HRG model, and extracting $T$ and $\mu_{B}$ such that the $\chi^{2}$ per degree of freedom for the fits is closest to unity. Once these parameters for the fireball are known, it was outlined in Ref. [106] how lattice QCD data on particular ratios of fluctuations, proposed in Ref. [109], can be used to predict the hadron multiplicities if indeed a thermalized medium is formed. One of the biggest source of uncertainties is in the determination of the freeze-out temperature [110]. It is thus important to calculate the freeze-out temperature from first principles in a model independent way since already several thermodynamic quantities start to show a departure from the HRG predictions for $T>140$ MeV. One other issue is to estimate the curvature of the freeze-out line at finite $\mu_{B}$ to understand how far apart in the QCD phase diagram are the chiral and freeze-out lines. If the outcome of this study would suggest a close proximity between these lines, then there is a hope to observe any effects of the QCD phase transition in the heavy-ion collisions assuming thermalization occurs.

5.1 Curvature of the freeze-out line: Input from the lattice

To estimate the curvature of the freeze-out line from first principles without resorting to hadron resonance gas model inspired fits, a new method has been proposed by combining lattice and experimental data in Ref. [111]. The freeze-out temperature to the leading order in baryon chemical potential $\mu_{B}$ can be parametrized as $T^{f}(\mu_{B})=T^{f}(0)\left[1-\kappa_{2}^{f}\mu^{2}_{B}/T^{2}\right]$ where $\kappa_{2}^{f}$ is the curvature of the freeze-out curve for small values of $\mu_{B}$. The inverse ratio of fluctuations of charges to the total number of charged species $R_{12}^{Q}(T,\mu_{B})=\chi_{1}^{Q}(T,\mu_{B})/\chi_{2}^{Q}(T,\mu_{B})$ was shown to act as a “baryometer” for the system i.e. at the leading order, give an estimate of $\mu_{B}$ of the fireball [23]. Equivalently $R_{12}^{B}$ is also another suitable candidate for a “baryometer” [23]. The ratios of these two simple quantities in the charge and the baryon sectors respectively, $R_{12}^{QB}=R_{12}^{Q}/R_{12}^{B}$, is recently proposed [111] as a suitable observable to extract the curvature since, in leading order in $\mu_{B}$

[TABLE]

In the above expression, the main systematic uncertainty comes from the uncertainty in estimating the value of $\mu_{B}$. In order to measure this ratio in a model independent way under the experimental conditions assuming thermalization, it could be simply expanded in terms of the “baryometer” $R_{12}^{B}$. A simple Taylor expansion of $R_{12}^{QB}(\mu_{B},T)$ thus gives,

[TABLE]

where the curvature of the freeze-out line is then extracted from the relation $c_{12}(T_{f,0},\kappa_{2}^{f})$ $=$ $c_{12}(T_{f,0})-\kappa_{2}^{f}D_{12}$ with quantities $c_{12},D_{12}$ at $\mu_{B}=0$ defined as,

[TABLE]

The superscript $2$ in all these expression denotes second derivatives of the quantities with respect to $\mu_{B}$ calculated at $\mu_{B}=0$. The observables at $\mu_{B}=0$ are very precisely measured on the lattice and even the continuum estimates are known [111]. In order that these ratios are sensitive to the experimental conditions, the quantities entering in these ratios at finite values of $\mu_{B}$ were extracted from the experimental data, to provide an estimate of $\kappa_{2}$ if indeed thermalization is achieved [111]. There are two important caveats to this procedure though. Experimentally, only the net-proton yields are measured and it needs to be ascertained how much of the net-baryon number is represented within this measurement. The second caveat one has to take into account is that the detectors have a finite momentum acceptance and can only measure hadrons which have transverse momenta beyond some threshold [112]. Surprisingly in the ratio $R_{12}^{QB}$, these two systematic uncertainties are argued to almost cancel one another [111]. Within HRG approximation atleast, one can show that the deviation between $R_{12}^{B}$ and for the protons denoted by $R_{12}^{P}$ is within $10\%$ which is equal and opposite in magnitude to the error introduced while neglecting the contribution of the undetected hadrons with momenta $p_{T}\leq 0.4$ GeV. The results for the ratio $R_{12}^{QB}(0)\equiv\Sigma$ at $\mu_{B}=0$ from Ref. [111] is summarized in Fig. 9. It is clearly visible that for this observable, the HRG approximation breaks down already at $T\sim 145$ MeV. Similar conclusions could be made for the quantities $c_{12}$ and $D_{12}(T)$.

A fit to the experimental data on $R_{12}^{QP}(\mu_{B})$ according to Eq. (20), yields an intercept $R_{12}^{QP}(0)$, which after comparing with the lattice data on $R_{12}^{QB}(0)$ gives a freeze-out temperature of $T_{f}(0)=147(2)$ MeV. Furthermore from the slope of the fit, the curvature comes out to be $\kappa_{2}^{f}(0)=-0.001(13)$ [111]. It agrees quite well with the curvature of the line of chiral transitions $\kappa_{2}^{B}=0.007$ using Taylor expansion method in $\mu_{B}$  [71, 113], but there is disagreement with the recent results of $\kappa_{2}^{B}\geq 0.015$ using imaginary chemical potential techniques [87, 88, 89].

The experimental measurements have effects due to the rapidly expanding medium and other systematic biases apart from the caveats mentioned earlier. Surprisingly the general systematic trends of the heavy ion data both from RHIC and LHC seem to suggest that these are consistent with predictions from finite temperature QCD using Lattice techniques. The current status of our understanding is outlined succinctly in Fig. 10 taken from Ref. [114]. The yellow band correspond to the current lattice estimate of the chiral curvature line and the blue band [115] and the black line [107] correspond to different statistical model fits to the experimental data. The blue boxes are predictions from a transport model fit to the LHC data for different center of mass energies of heavy-ion collisions [116]. It is evident that at vanishingly small $\mu_{B}$, the freeze-out temperatures predicted from lattice are lower than the ones obtained from phenomenological freeze-out curves. Though the experimental data points are currently consistent with both lattice and phenomenological curves within errors, their central values are systematically lower that the phenomenological freeze-out curves. A more precise determination of the freeze-out curve from the lattice would enable to clearly distinguish between all these different scenarios. With the new observables and methods like in Ref. [111] developed within Lattice QCD thermodynamics will allow us to elucidate the pattern within the complex and rich set of data from heavy-ion experiments with an ultimate goal to understand the phase diagram of QCD.

6 Summary and outlook

We have reviewed selected results on the QCD phase structure and bulk thermodynamics, obtained by recent lattice calculations. We have discussed results on derivatives of the logarithm of the QCD partition function ($\ln Z$) with respect to temperature and quark mass ($T,m$) and various chemical potentials $\vec{\mu}$. Results on correlation-functions have been not touched and are beyond the scope of this review.

We have argued that the QCD equation of state (EoS), i.e. the pressure, energy and entropy density as a function of temperature, obtained by the integral of the temperature derivative of $\ln Z$, has been reliably determined in the temperature interval of $(130-400)$ MeV. The EoS is the basis for every hydro dynamical model of the fireball created in heavy ion collisions. We find that its low and high temperature limits seem to be in good agreement with HRG model calculations and (re-summed) perturbative QCD calculations, respectively. However, more refined calculations in both of these temperature limits are needed.

We have shown that the behavior of the chiral condensate and (mixed) chiral susceptibility exhibit universal critical behavior in the chiral limit of ($2+1$)-flavor QCD, which is taken by sending the quark mass of the two light flavors to zero while keeping the strange quark mass constant. Matching QCD to the universal theory through the determination of a set of non-universal parameter, yields valuable informations on the QCD phase diagram. As an example we have discussed the pseudo-critical line as function of $m$ and $\mu_{B}$. Although we find a scaling behavior that is consistent with an $O(4)$ symmetric critical point in the chiral limit, the possibility of a first order transition in the chiral limit and correspondingly a critical $Z(2)$ end-point at very small but finite quark masses is not excluded yet. This issue is closely related with the possibility of a sufficiently restored axial anomaly ($U_{A}(1)$-symmetry) close to $T_{c}$, which is currently subject to many lattice QCD calculations with chiral fermions. Most studies do, however, support the opposite scenario that $U_{A}(1)$ is strongly broken around $T_{c}$.

Finally we have discussed derivatives of $\ln Z$ with respect to various chemical potentials (generalized susceptibilities). We have introduced combinations of these generalized susceptibilities which allow for the disentanglement of the partial pressures of different hadronic sectors and thus lead to a better understanding of the relevant microscopic degrees of freedom. We presented evidence for experimentally not yet measured strange and charmed hadrons and showed that also flavored hadrons start to melt at the chiral transition temperature. Vise versa, we used similar combinations of these generalized susceptibility to investigate the approach to the free gas behavior. We find that the temperature range of $1<T/T_{c}<2$ can be characterized by a strongly interacting plasma of quarks, gluons and some remnants of hadron like states, while above $T/T_{c}\gtrsim 2$ the quasi-particle spectrum resembles that of quarks and gluons alone.

The same generalized susceptibilities also encode the fluctuations of conserved charges, which are measured at RHIC and LHC. We discussed certain ratios of cumulants of conserved charges and showed that from a fit to the corresponding experimentally measured data, the freeze-out temperature, as well as the freeze-out curvature can be obtained. We find that the so determined freeze-out line is in good agreement with current determinations of the pseudo-critical line.

We emphasize that several of the $4^{\rm th}$ and $6^{\rm th}$ order susceptibilities that are used in the above summarized studies are not yet continuum extrapolated. Although we do not expect large cut-off effects for these observables, continuum extrapolations are necessary to consolidate our findings and are currently work in progress. This will also allow to extent the range of the QCD equation of state to $\mu_{B}/T\lesssim 2$, which corresponds to the full range of the energy scan program at RHIC [117]. Even though some results on $8^{\rm th}$ order susceptibilities are already available and will be refined in future, it is evident that there are many interesting aspects on the QCD phase diagram that will never be in reach by a Taylor expansion method. Therefore a reformulation of QCD in terms of other variables that allows for direct simulations at nonzero chemical potential will be needed.

Acknowledgment

We thank all members of the Bielefeld-BNL-CCNU collaboration for valuable discussions and comments. S.S. acknowledges support by the U.S. Department of Energy under Grant No. DE-SC0012704.