Phenomenological Equation of State of Strongly Interacting Matter with First-Order Phase Transitions and Critical Points

An extension of the relativistic density functional approach to the equation of state for strongly interacting matter is suggested which generalizes a recently developed modified excluded-volume mechanism to the case of temperature and density dependent available-volume fractions. A parametrisation of this dependence is presented for which at low temperatures and suprasaturation densities a first-order phase transition is obtained. It changes for increasing temperatures to a crossover transition via a critical endpoint. This provides a benchmark case for studies of the role of such a point in hydrodynamic simulations of ultrarelativistic heavy-ion collisions. The approach is thermodynamically consistent and extendable to finite isospin asymmetries that are relevant for simulations of neutron stars, their mergers and core-collapse supernova explosions.


Introduction
The simulation of astrophysical phenomena, such as core-collapse supernovae (CCSN) or neutron-star (NS) mergers, requires a careful modeling of strongly interacting matter in a wide range of densities and temperatures. The same applies to the theoretical description of heavy-ion collisions (HIC) that study compressed baryonic matter in the laboratory from low to high beam energies. The properties of such matter are represented by the equation of state (EoS) that provides information on pressure, entropy, energies, and other thermodynamic variables of interest.
A particular feature of QCD matter is the supposed phase transition (PT) from hadronic matter to quark matter when density or temperature increase to sufficiently high values. A strong first-order PT could allow for the existence of a third branch of compact stars and the occurrence of the twin-star phenomenon [1][2][3]. Signals of the PT might also have direct consequences in dynamical processes when matter in the quark phase expands and cools down, e.g., the release of a second neutrino burst in CCSN [4][5][6]. For a recent review on the role of the EoS in CCSN simulations, see [7].
The theoretical description of the hadron-quark PT in strongly interacting matter often relies on a construction employing different models for the two phases. With such an approach, the coexistence line of the first-order PT will usually connect a point on the zero-temperature axis at finite baryon chemical potential µ B with a point at finite temperature T on the zero baryon chemical potential axis. By this construction, the QCD hadron-quark PT is of first order in the whole temperature-density plane, see for instance [8]. However, from lattice QCD studies, it is known that there is a smooth crossover at µ B = 0 with increasing T [9,10], so at least one critical point at finite µ B and T is expected to exist. Other possibilities are that the character of the transition is crossover all over the QCD phase diagram [11] or, as is advocated in studies of the BEC-BCS crossover transition in low-temperature QCD, that a second critical endpoint exists [12][13][14]. Since lattice QCD studies are presently incapable of exploring the EoS close to the presumed critical point with much confidence, unified models are needed that can account for the existence of these features, see, e.g., [15]. There are dedicated microscopic models available that incorporate the major expected features in the QCD phase diagram, e.g., chiral mean-field models [16] or parity-doublet quark-hadron models [17]. Simulations of CCSN or HIC that are based on a hydrodynamic description of matter during dynamical evolution use the thermodynamic properties of matter encoded in the EoS as an input. Such data can be provided by phenomenological models that do not need to incorporate all of the details of the underlying physics.
In this work, a novel approach is introduced to provide a phenomenological EoS of baryonic matter that exhibits a first-order PT and a critical point at densities and temperatures expected in QCD matter. The parameters of the model can be adjusted to place the coexistence line at arbitrary positions in the phase diagram. The description uses an extension of a relativistic energy density functional for hadronic matter assuming a medium-dependent change in the number of degrees of freedom. This approach employs a recently developed version of a modified excluded-volume (EV) mechanism that gives a thermodynamically consistent EoS with nuclear matter properties that are consistent with present constraints. Here, we concentrate on the hadron-quark transition but not on the liquid-gas PT, which is also contained in our model. The model allows us to study the PT for arbitrary isospin asymmetries; however, only isospin-symmetric matter is considered in this first exploratory study for simplicity. In the present work, no attempt was made to reproduce the EoS of QCD matter at vanishing baryon chemical potential obtained in lattice QCD studies. With appropriately chosen EV parameters, the crossover transition with increasing temperature can be well modeled, even for imaginary chemical potentials, e.g., in a hadron resonance gas model [18]. With improved parametrizations, the structure of the phase diagram in the full space of variables, i.e., temperature, baryon density/chemical potential, and isospin asymmetry, can be investigated in the future.
The theoretical formalism of the model is presented in Section 2, which includes the main equations that define the relevant thermodynamic quantities in Section 2.1. In Section 2.2, details of the parametrization of the interaction and of the effective degeneracy factors are given. They account for the change in the number of degrees of freedom with density and temperature. The phase transitions are explored in Section 3 for isospin-symmetric matter. Conclusions follow in Section 4.

Theoretical Model
The theoretical description of strongly interacting matter in the present work is adapted from the model introduced in [19]. It combines a relativistic mean-field (RMF) approach for hadronic matter with density-dependent nucleon-meson couplings and a modified EV mechanism. Here, it is sufficient to provide only the main equations without a detailed derivation. The essential quantities that determine the position of the PT and the critical point in the phase diagram are the effective degeneracy factors that depend on the number densities of the particles and the temperature.

Relativistic Energy Density Functional with Modified Excluded-Volume Mechanism
The present model assumes neutrons and protons as well as their antiparticles as the basic degrees of freedom. These particles interact by the exchange of mesons, and the model effectively describes the short-range repulsion (ω meson), the intermediate-range attraction (σ meson), and the isospin dependence of the nuclear interaction (ρ and δ mesons), as is common of RMF models. The contribution of leptons or other degrees of freedom like nuclei, hyperons or photons, as required for multi-purpose EoS for astrophysical applications [20], is not considered here.
The nucleons i = n, p,n,p with rest masses m i are treated as quasi-particles of energy, which depends on the particle momentum k and the scalar (S i ) and vector (V i ) potentials. Denoting the particle chemical potentials with µ i , the contribution of the quasi-particles to the total pressure of the system can be written as where the medium-dependent effective degeneracy factors are a product of the usual degeneracy factor g i = 2 for nucleons and the available-volume fraction Φ i , which is defined in Section 2.2. The meson contribution to the total pressure in Equation (2) contains the coupling factors of the mesons given as a ratio of the density-dependent coupling functions Γ j and the meson masses m j . The source densities for vector mesons (j = ω, ρ) and for scalar mesons (j = σ, δ) in Equation (5) are obtained from the quasi-particle vector densities and scalar densities n (s) with the Fermi-Dirac distribution function The scaling factors in Equations (7) and (8) determine the coupling between mesons and nucleons. They also appear in the vector potential and the scalar potential in the quasi-particle energy (Equation (1)). The rearrangement potential contributes to the vector potential (Equation (16)) because the couplings Γ j in Equation (6) are assumed to depend on the baryon density is the baryon number of particle i, and the quantities C ′ j = dC j /dn B are the derivatives of the coupling factors. The dependence of the available-volume fractions Φ i in the effective degeneracy factor (Equation (4)) on the vector or scalar quasi-particle densities (9) and (10) also generates rearrangement contributions and S (r) in the potentials (Equations (16) and (17)), respectively. Furthermore, there is a rearrangement term in the total pressure (Equation (2)) with two contributions from the density dependence of the couplings p (r) and p (r) from the EV effects. The free energy density of the system is obtained with the total pressure (Equation (2)) and the chemical potentials of the particles µ i . They are not independent since, for nucleons with baryon number B i and charge number Q i , they are given by with the baryon chemical potential µ B and the charge chemical potential µ Q . Only the latter two are independent quantities. For the internal energy density the entropy density is needed. Besides the standard contribution depending on the distribution functions f i , there is a term from the possible temperature dependence of the available-volume fractions Φ i . In order to guarantee the third law of thermodynamics, i.e., lim T→0 s = 0, the temperature derivative of the available-volume fractions has to vanish for T → 0. After solving the equations above for a given T, µ B , and µ Q , a fully consistent thermodynamic EoS is obtained. For practical purposes, however, the baryon density n B and the hadronic charge fraction are used as independent variables instead of µ B and µ Q . A possible shortcoming of models that consider EV effects is the potential appearance of a superluminal speed of sound in certain regions of the space of thermodynamic variables, see, e.g., [21,22]. This causality constraint has to be checked case by case depending on the specific implementation of the EV mechanism.

Available-Volume Fractions and Model Parameters
For a quantitative evaluation of the EoS in the present approach, the functional forms of the meson-nucleon couplings Γ j and the available-volume fractions Φ i have to be specified as well as all parameters. Here, we use the masses of nucleons and mesons and the coupling functions of the DD2 parametrization presented in [23]. It only considers ω, ρ, and σ mesons for the effective description of the nuclear interaction but not the δ mesons. The parameters were obtained by fitting observables (binding energies, radii, etc.) of selected nuclei. With this set, the EoS of nuclear matter at zero temperature exhibits characteristic nuclear matter parameters, e.g., the saturation density (n sat = 0.149065 fm −3 ), binding energy at saturation (B = 16.02 MeV), incompressibility (K = 242.7 MeV), symmetry energy (J = 32.73 MeV), and slope (L = 57.94 MeV), that are consistent with modern constraints from experiment and theory.
EV effects are frequently employed to describe an effective repulsive interaction between particles, in particular in calculations of the EoS in the framework of hadron resonance gas models, see, e.g., [24]. For zero baryon density, a comparison of the resulting EoS with results from lattice QCD calculations can be used to fix the volume parameters. If a finite volume v i is attributed to each particle i, the available volume for the motion of the particle is reduced from the total system volume V to VΦ (cl) with the classical available-volume fraction Clearly, there is a limiting density above which a compression of the system becomes impossible. In general, the volumes and available-volume fractions can depend on the particle species, and the EV mechanism can be used to suppress particles, e.g., nuclei, in a mixture, causing them to disappear above a certain density, see, e.g., [25] for applications to the low-temperature and low-density EoS in astrophysical simulations.
In [19], the interpretation of the EV mechanism was changed by moving the factor Φ i from the system volume V to the degeneracy factor g i as in Equation (4) and allowing the available-volume fractions to be arbitrary functions of the particle densities and temperature. The medium dependence of the effective degeneracy factors is interpreted as a change in the effective number of degrees of freedom. A decrease in g (eff) i has the effect of a repulsive interaction between the particles, whereas an increase can be seen as the action of an attractive interaction, c.f., the softening of the nuclear EoS when hyperons are included, see, e.g., [26,27]. This freedom leaves the room to modify the properties of an EoS in a favored way.
In the present application of the modified EV mechanism, the available-volume fraction is defined to be identical for all particles i as depending on the temperature T and an auxiliary quantity depending on T and the quasi-particle vector densities n (v) i with parameters s, v, and a cutoff density of n cut . The θ function in Equation (30) guarantees that Φ i = 1 for x ≤ 0 and the EV mechanism has no effect on the EoS. The functions g 1 in Equation (30) and g 8 in Equation (31) are defined as with parameters t and T 0 . In the limit T → 0, the function g t approaches one and it decreases with increasing temperature. Furthermore, the derivatives ∂Φ i /∂T approach zero for T → 0, as required for the thermodynamic consistency, because of the choice of the function g t (T). For T → T 0 , the function g t vanishes very smoothly and there are no effects at higher temperatures because Φ i = 1. In order to reproduce the correct high-temperature limit, given by a Stefan-Boltzmann-type behavior, a modification of the available-volume fractions for temperatures well above T 0 is required. This is left to future extensions of the model. According to Equation (31), the quantity x is only positive for baryon densities larger than g 8 n cut , a quantity that decreases with increasing temperature. There are no artificial singularities due to the presence of the θ functions in Equations (30) and (32) because all derivatives of the exponential functions are zero when the arguments of the θ functions vanish. The actual values of the parameters for the modified EV mechanism used in the present study are s = 3, v = 2 fm 3 , T 0 = 270 MeV, and n cut = n sat of the DD2 parametrization.

Results
In order to illustrate the characteristic effects of the modified EV mechanism on the EoS, we limit the presentation to the case of symmetric matter, i.e., Y q = 0.5. The pressure p and baryon chemical potential µ B are calculated as a function of the baryon density n B . Due to the increase in the available-volume fractions Φ i or the effective degeneracy factors g (eff) i , a considerable softening of the EoS is observed in a certain range of densities. Below the critical temperature T crit , the pressure is not a monotonous function of the baryon chemical potential for an isotherm. A Maxwell construction is used to determine the two densities of the coexisting phases where p and µ B are identical. The pressure p as a function of the baryon density n B is depicted in Figure 1 for selected temperatures. Note that an integral part of the underlying RMF model is a detailed description of the liquid-gas PT in nuclear matter and the formation and dissociation of nuclear clusters in compact-star matter. For details, see, e.g., [28].  In the coexistence region of the pseudo hadron-quark PT between the low and high density phases, the pressure is constant as typical for a first-order PT. The area of coexistence is enclosed by the binodal, and, above the critical temperature T crit ≈ 155.5 MeV, there is no PT anymore. The peculiar shape of the binodal is a result of the specific form (30) of the available-volume fractions. It can be adjusted with appropriate changes in the functional form and parameters.
The binodals of the liquid-gas and pseudo hadron-quark PT in the temperature-baryon density plane are shown in panel (a) of Figure 2. At vanishing temperatures, the coexistence region of the pseudo hadron-quark PT covers a density range from 0.270 to 0.408 fm −3 , well above the nuclear saturation density n sat . At higher temperatures, it moves to lower densities with an almost constant extension in baryon density except for temperatures close to T crit . Here, the critical density is found as 0.201 fm −3 , still above n sat . The dashed line in panel (a) marks the boundary between regions without (lower left) and with (upper right) effects of the modified EV mechanism in the present parametrization. It corresponds to the condition x = 0. There is another region in the phase diagram without modified EV effects at temperatures above T 0 , outside the figure.   Figure 2 depicts the lines of the first-order PT in the temperature-baryon chemical potential diagram ending in critical points. With increasing temperature, the baryon chemical potential at the pseudo hadron-quark PT reduces from 979.1 to 591.8 MeV at the critical point. By crossing the transition line, an abrupt change in the density occurs that becomes continuous at the critical point.

Conclusions
The extension of the modified EV approach to a density-and temperature-dependent parametrization of the available-volume fractions as introduced in this work was successful in achieving the main goal of this study: As a generic structure of the QCD phase diagram, a first-order pseudo hadron-quark phase transition at low temperatures and a crossover for low baryon densities could be modeled that also includes a critical endpoint at T crit = 155.5 MeV and µ B,crit = 591.8 MeV. Other patterns of the QCD phase diagram that have been theoretically motivated could also be modeled within the present approach. Further extensions of the model are straightforward. They include the extension to a larger number of components of the hadron resonance gas in the underlying RMF model and an isospin dependence of the EV model. It would be worthwhile to study further thermodynamic quantities such as the speed of sound, heat capacities, or the susceptibilities in such an enlarged model. It would also be interesting to elaborate on a parametrization that would result in a second endpoint at low temperatures, as suggested by Hatsuda et al. [14]. The so-generalized parametrization of the QCD EoS can be used in Bayesian analysis studies for astrophysical applications pertaining to compact stars [29,30], their mergers, and core-collapse supernova explosions, as well as heavy-ion collisions analogous to those studied in [31].