QCD equation of state at finite baryon density with fugacity expansion

TL;DR

This paper develops a model for the QCD equation of state at finite baryon density using lattice data and fugacity expansion, enabling improved simulations of heavy-ion collisions.

Contribution

It introduces a crossover model based on the cluster expansion model fitted to lattice QCD data for finite baryon density.

Findings

Accurate description of lattice data at imaginary and zero chemical potential.

A parameterized pressure function for finite baryon density.

Application potential in fluid dynamical simulations of heavy-ion collisions.

Abstract

We explore the QCD equation of state at finite baryon density through an expansion in baryon number fugacity, making use of the recent lattice data on the four leading Fourier coefficients of the expansion. A state-of-the-art description of the lattice data at both imaginary $\mu_B$ and $\mu_B = 0$ is provided by the cluster expansion model (CEM). The CEM pressure function has three temperature dependent input parameters, which we fix by parameterizing the available lattice QCD data. This results in a crossover model of the QCD equation of state at finite baryon density, which can be used in fluid dynamical simulations of heavy-ion collisions.

Tables

Table 1. Table 1: Values of the switching function parameters T 0 subscript 𝑇 0 T_{0} and r 𝑟 r for the parameterizations of b 1 ​ ( T ) subscript 𝑏 1 𝑇 b_{1}(T) , b ​ ( T ) 𝑏 𝑇 b(T) , and p 0 ​ ( T ) subscript 𝑝 0 𝑇 p_{0}(T) .

	$b_1\left(T\right)$	$b\left(T\right)$	$p_0\left(T\right)$
$T_0$ [MeV]	170	175.5	100
$r$	7.95	8	7

Equations

$$\frac{p(T,{\mu }_B)}{T^4}=\frac{1}{2}\sum _{k=-\infty }^{\infty }{p}_{|k|}(T)e^{k{\mu }_B/T},$$

$$\frac{{\rho }_B(T,{\mu }_B)}{T^3}$$

$$=$$

$${\frac{{\rho }_B(T,{\mu }_B)}{T^3}|}_{{\mu }_B=i{\theta }_{B}T}$$

$$b_k(T)=\frac{2}{\pi }{\int }_0^{\pi }\mathrm{Im}\left[\frac{{\rho }_B(T,i{\theta }_{B}T)}{T^3}\right]\sin (k{\theta }_B)d{\theta }_B.$$

$${\alpha }_k=\frac{[b_1(T){]}^{k-2}}{[b_2(T){]}^{k-1}}{b}_k(T),k=3,4,\dots$$

$$b_k^{\mathrm{SB}}$$

$${\alpha }_k^{\mathrm{SB}}$$

$$b_k(T)={\alpha }_k^{\mathrm{SB}}\frac{[b_2(T){]}^{k-1}}{[b_1(T){]}^{k-2}},k=3,4,\dots$$

$$\frac{{\rho }_B(T,{\mu }_B)}{T^3}=-\frac{2}{27{\pi }^2}\frac{{\hat{b}}_1^2}{{\hat{b}}_2}\left\{4{\pi }^2\left[{\text{Li}}_1\left(x_{+}\right)-{\text{Li}}_1\left(x_{-}\right)\right]+3\left[{\text{Li}}_3\left(x_{+}\right)-{\text{Li}}_3\left(x_{-}\right)\right]\right\}.$$

$${\hat{b}}_k^{\mathrm{RFM}}(T)=\frac{c(T)}{1+k/k_0(T)},k=3,4,\dots$$

$$k_0(T)$$

$$c(T)$$

$${\chi }_{2n}^B(T)\equiv {\frac{{\partial }^{2n-1}({\rho }_B/T^3)}{\partial ({\mu }_B/T{)}^{2n-1}}|}_{{\mu }_B=0}=\sum _{k=1}^{\infty }{k}^{2n-1}{b}_k(T).$$

$$\frac{p(T,{\mu }_B)}{T^4}$$

$$b_1(T)$$

$$S(T)$$

$$b_1^{\mathrm{lat}}(T)=\frac{b_1^{\mathrm{sb}}+a_n/t_b+b_n/t_b^2}{1+a_d/t_b+b_d/t_b^2},t_b=T/T_0,$$

$$a_n=-0.940,b_n=0.345,a_d=-1.336,b_d=\mathrm{0.502.}$$

$$b_2(T)=-[b_1(T){]}^2{T}^{3}b(T).$$

$$b(T)=\left[1-S(T)\right]\cdot b^{\mathrm{hrg}}+S(T)\left[-\frac{b_2^{\mathrm{sb}}}{(b_1^{\mathrm{sb}}{)}^2{T}^3}\right].$$

$$p_0(T)=\left[1-S(T)\right]p_0^{\mathrm{hrg}}+S(T)2\left[\frac{p^{\mathrm{lat}}(T,{\mu }_B=0)-\Delta p(T,{\mu }_B=0)}{T^4}\right].$$

Full text

1 Introduction

The equation of state of hot QCD matter is available for $\mu_{B}=0$ from first-principle lattice QCD simulations [1, 2]. Finite $\mu_{B}$ simulations, on the other hand, are hindered by the sign problem. All the available lattice methods here are indirect and restricted to small values of $\mu_{B}/T$. These methods include the analytic continuation from imaginary chemical potential [3, 4], the Taylor expansion around $\mu_{B}=0$ [5, 6], the reweighing technique [7, 8], and the canonical approach [9, 10]. At larger baryon density the equation of state is usually described by effective models.

Recent advances in lattice QCD simulations at zero and imaginary chemical potential include the evaluation of higher-order baryon number susceptibilities at $\mu_{B}=0$ [11, 12, 13] and the Fourier coefficients of net baryon density at imaginary $\mu_{B}$ [14], which put stringent constraints on effective models for the QCD equation of state at finite baryon density. Here we construct a model for the QCD equation of state at finite baryon density which incorporates all of the above-mentioned constraints. Our considerations are based on the recently developed cluster expansion model (CEM) formalism [15, 16]. It uses the cluster expansion in baryonic fugacity,

[TABLE]

which represents the pressure as a Laurent series in $\lambda_{B}\equiv e^{\mu_{B}/T}$. This expansion incorporates two important QCD symmetries: the CP-symmetry ($\mu_{B}\to-\mu_{B}$) and the Roberge-Weiss periodicity ($\mu_{B}\to\mu_{B}+i\,2\pi T$) of the partition function. The net baryon density reads

[TABLE]

At purely imaginary baryochemical potential the net baryon density has the form of a trigonometric Fourier series expansion

[TABLE]

The expansion coefficients $b_{k}$ become Fourier expansion coefficients which can be extracted through the standard Fourier analysis:

[TABLE]

These Fourier coefficients can be evaluated in lattice QCD simulations at imaginary chemical potential [17, 18]. The four leading coefficients were evaluated at the physical point in Ref. [14] in the temperature range $135<T<230$ MeV.

2 Cluster expansion model

2.1 Formulation

The CEM [15, 16] is a model for the QCD equation of state at finite baryon density which makes use of the expansion in fugacities (2). It is based on an empirical observation of temperature independent Fourier coefficient ratios

[TABLE]

as observed in the lattice data of Ref. [14] for $k=3,4$. The available lattice data are consistent with a temperature independent $\alpha_{3,4}$ values, as evaluated in the Stefan-Boltzmann limit of massless quarks:

[TABLE]

The comparison of the lattice data for $\alpha_{3}$ and $\alpha_{4}$ with Eq. (8) is shown in Fig. 1.

The CEM assumes $\alpha_{k}=\alpha_{k}^{\rm SB}$ for all $k\geq 3$. This implies that all higher-order Fourier coefficients are given in the CEM in terms of the two leading ones:

[TABLE]

The expansion in Eq. (2) can be analytically summed in the CEM (9). The result is

[TABLE]

Here $\hat{b}_{k}=\displaystyle\frac{b_{k}(T)}{b_{k}^{\rm SB}}$, $x_{\pm}=\displaystyle-\frac{\hat{b}_{2}}{\hat{b}_{1}}\,e^{\pm\mu_{B}/T}$, and $\textrm{Li}_{s}(z)=\displaystyle\sum_{k=1}^{\infty}\frac{z^{k}}{k^{s}}$ is the polylogarithm.

2.2 Comparison with the rational function model

In general, a determination of the higher-order Fourier coefficients from the two leading ones is not unique. Functional forms different from (9) can be considered. This point has been raised in Ref. [19], where a rational function model (RFM) has been introduced. The Fourier coefficients in the RFM are given as follows:

[TABLE]

with

[TABLE]

The coefficients in the RFM exhibit a different asymptotic behavior in the temperature range considered: a power-law damping at large $k$ instead of an exponential damping in the CEM. The RFM gives a similarly good description of the lattice data for the coefficients $b_{3}(T)$ and $b_{4}(T)$ as the CEM when these are viewed on the linear scale (see Ref. [19]). However, significant differences can be seen by considering the $\alpha_{3}$ and $\alpha_{4}$ ratios [Eq. (6)], which are plotted for the CEM and the RFM in Fig. 1. While the two models describe similarly well the lattice data at $T\gtrsim 200$ MeV, the RFM shows increased values of $\alpha_{3}$ and $\alpha_{4}$ at $T\lesssim 190$ MeV, which is not supported by the available data.

The subtle differences in the higher-order Fourier coefficients between the two models may show up in other observables, in particular those sensitive to the equation of state at finite baryon density. Particularly interesting ones are the baryon number susceptibilities at $\mu_{B}=0$, which can formally be given in terms of a series of the Fourier coefficients:

[TABLE]

The sum (14) can be carried out analytically in both the CEM and the RFM (see Refs. [16] and [19], respectively, for the details).

Figure 2 depicts the temperature dependence of the susceptibility ratios $\chi_{4}^{B}/\chi_{2}^{B}$ and $\chi_{6}^{B}/\chi_{2}^{B}$, as calculated within the CEM (red stars) and the RFM (dashed red lines) along with the lattice data of the Wuppertal-Budapest [13] and HotQCD collaborations [12]. The two models yield similar results in the low and high temperature limits while a notable difference is seen in the range $150\lesssim T\lesssim 180$ MeV. This difference can already be distinguished within the precision of the available lattice data, which appears to favor the CEM among the two models considered. The CEM also describes fairly well the available lattice estimates for $\chi_{8}^{B}$ [13] (see Ref. [16] for the agreement level).

While a determination of the higher-order Fourier coefficients from the lower order ones is not unique, we conclude that the CEM appears to be the best available tool for constructing a model of the QCD equation of state which is consistent with both the $\mu_{B}=0$ and imaginary $\mu_{B}$ lattice data for physical quark masses. We, therefore, proceed with constructing the full pressure function $p(T,\mu_{B})$ within the CEM framework.

3 The CEM equation of state

The CEM pressure is obtained by integrating the net baryon density (10) over $\mu_{B}$:

[TABLE]

The CEM input parameters are the temperature-dependent coefficients $p_{0}(T)$, $b_{1}(T)$, and $b_{2}(T)$. At high temperatures ($T\gtrsim 130$ MeV) these can be extracted from the available lattice data. At low temperatures they can be matched to the equation of state of the hadron resonance gas (HRG) model. This matching is achieved through a smooth switching function [20].

The first Fourier coefficient, $b_{1}(T)$, is parameterized as follows:

[TABLE]

Here $b_{1}^{\rm hrg}$ is the partial pressure of baryons and antibaryons at $\mu_{B}=0$ in the ideal hadron resonance gas (HRG) model. We evaluate it using the HRG model of Ref. [14]. $b_{1}^{\rm lat}$ corresponds to the lattice data of the Wuppertal-Budapest collaboration [14] at sufficiently high temperatures. We parametrize $b_{1}^{\rm lat}(T)$ as

[TABLE]

with the following parameter values:

[TABLE]

The values of the parameters $T_{0}$ and $r$ are listed in the first row of Table 1.

The second Fourier coefficient, $b_{2}(T)$, is represented in a form

[TABLE]

The parameter $b(T)$ can be interpreted as a baryonic eigenvolume parameter in a hadron gas with repulsive baryonic interactions [14]. We observe that the lattice data for $b(T)$ [14] can be well described at high temperatures by a power-law dependence, $b(T)\sim T^{-3}$. At low temperatures we assume that the matter corresponds to a HRG with residual repulsive baryonic interactions, which are characterized by the eigenvolume parameter value of $b^{\rm hrg}=1$ fm3. The complete temperature dependence of $b$ is, therefore, parametrized as

[TABLE]

The values of the switching function parameters $T_{0}$ and $r$ for $b(T)$ are listed in the second row of Table 1.

In order to parameterize $p_{0}(T)$ we use the lattice QCD data for the pressure at zero chemical potential, as well the already obtained parameterizations for $b_{1,2}(T)$. At small temperatures $p_{0}(T)$ is matched with the doubled partial pressure of mesons in an ideal HRG model. A complete parameterization is thus given by the following:

[TABLE]

Here $p^{\rm lat}(T,\mu_{B}=0)$ is evaluated using the $p/T^{4}$ parametrization of the HotQCD collaboration from Ref. [2] whereas $\Delta p(T,\mu_{B}=0)$ depends on $b_{1}(T)$ and $b_{2}(T)$ only, both of which are already parametrized. The values of the switching function parameters $T_{0}$ and $r$ for $b(T)$ are listed in the third row of Table 1.

Note that within the parametrization considered, $p_{0}$, $b_{1}$, and $b_{2}$ all smoothly approach the Stefan-Boltzmann limit of massless quarks and gluons in the high-temperature limit.

The temperature dependence of the resulting pressure $p/T^{4}$ at $\mu_{B}=0$, the pressure contribution $p_{0}/2$ of the 0th Fourier harmonic, and the Fourier coefficients $b_{1}$ and $b_{2}$ is depicted in Fig. 3 for a broad temperature range $50<T<800$ MeV. It is seen that $p_{0}/2$ corresponds to the bulk of the total pressure at $\mu_{B}=0$ across all temperatures. Note that $p_{0}/2$ can be associated with the partial pressure of mesons at sufficiently small temperatures, where a HRG picture can be applied.

Different isentropic trajectories in the $\mu_{B}$-$T$ plane are depicted in Fig. 4 to illustrate the behavior of the CEM equation of state at finite baryon densities for conditions which can be expected to be created in relativistic heavy-ion collisions at various energies. Seven entropy-per-baryon values are considered, $S/N_{B}=10$, 15, 20, 30, 50, 150, and 400, which cover collision energies ranging from low SPS to top RHIC. The overall behavior of the isentropes is reasonable, smoothly interpolating between the isentropes of an ideal HRG model (dashed lines in Fig. 4) and those of an ideal massless gas of quarks and gluons (dot-dashed lines in Fig. 4).

The CEM equation of state yields a broad crossover transition across the whole $\mu_{B}$-$T$ plane. The model, therefore, corresponds to the no-critical-point scenario. It should be noted that the model is not suited for the cold and dense region of nuclear matter at low temperatures, as noted in Ref. [15]. A considerably more involved approach might be necessary to incorporate contraints from both the nuclear matter properties and lattice QCD data at $\mu_{B}=0$ simultaneously [21].

4 Summary

The QCD equation of state at finite baryon density, based on the parameterized form of the cluster expansion model, has been presented. This provides a first example of an equation-of-state model, which is originally formulated at purely imaginary baryochemical potential, and then analytically continued to real values of $\mu_{B}$. A particular feature of the CEM equation of state is that it is able to incorporate constraints from both the lattice simulations at $\mu_{B}=0$ (net baryon number susceptibilities) and imaginary $\mu_{B}$ (Fourier coefficients of net baryon density) simultaneously. The model provides a reasonable behavior of the isentropes in the $\mu_{B}$-$T$ plane and can be used in fluid dynamical simulations of heavy-ion collisions from the lowest SPS energies up to the highest LHC energies. The CEM equation of state presented here is available in a tabulated format [22].

The present CEM equation of state exhibits a regular behavior at real values of $\mu_{B}$ and thus corresponds to the no-critical-point scenario for the QCD phase diagram, exhibiting a crossover-type transition all over the $\mu_{B}$-$T$ plane. A generalization of the present model to incorporate a critical point at a certain location at finite baryon density is a possible future endeavor. It should also be noted that the CEM is currently restricted to the baryochemical potential only. A generalization of the approach to include the strangeness and electric charge chemical potentials is another possible future avenue to explore.

Acknowledgments

H.St. acknowledges the support through the Judah M. Eisenberg Laureatus Chair at Goethe University, and the Walter Greiner Gesellschaft, Frankfurt.