MUSES Collaboration
Phase Stability in the 3-Dimensional Open-source Code for the Chiral mean-field Model
Abstract
In this paper we explore independently for the first time three chemical potentials (baryon , charged , and strange ) in the Chiral mean-field (CMF) model. We designed and implemented CMF++, a new version of the CMF model rewritten in C++ that is optimized, modular, and well-documented. CMF++ has been integrated into the MUSES Calculation Engine as a free and open source software module. The runtime improved in more than 4 orders of magnitude across all 3 chemical potentials, when compared to the legacy code. Here we focus on the zero temperature case and study stable, as well as metastable and unstable, vacuum, hadronic, and quark phases, showing how phase boundaries vary with the different chemical potentials. Due to the significant numerical improvements in CMF++, we can calculate for the first time high-order susceptibilities within the CMF framework to study the properties of the quark deconfinement phase transition. We found phases of matter that include a light hadronic phase, strangeness-dominated hadronic phase, and quark deconfinement within our , , phase space. The phase transitions are of first, second (quantum critical point), and third order between these phases and we even identified a tricritical point.
I INTRODUCTION
In the past decades, the increase of colliding energy in particle accelerators and the unprecedented accuracy in astronomical observations allowed us to grasp a better understanding of the building blocks of matter, the quarks, and gluons. In a way, this allows us to glimpse at the matter that existed in the first s after the Big Bang. In the laboratory, it was shown that at extremely high temperatures the quark-gluon plasma created presents very low viscosity, behaving like an ideal fluid [1]. On the other hand, neutron stars reach extremely large baryon densities, the value being model dependent but attaining more than 14 times nuclear saturation density, in extreme cases [2], 10 for the heaviest neutron stars. Around these densities, several microphysical models have predicted deconfined quark matter within the core of neutron stars (starting with Ivanenko et. al in the 60’s [3]), while being consistent with astrophysical data, see e.g. [4]. Finally, mergers of neutron stars provide both hot and dense environments, where deconfined quarks may be observed not only electromagnetically, but also gravitationally [5, 6].
At lower energies, due to asymptotic freedom, quarks and gluons are confined within hadrons. At even lower energies baryons form atomic nuclei. These different “phases” of matter, which can be produced both in the laboratory and in the cosmos are usually depicted in a phase diagram, the Quantum chromodynamics (QCD) phase diagram, referring to the theory that describes quarks, gluons, and their interactions. The phase transition from nuclei to hadronic matter (composed of baryons with 3 quarks and mesons with one quark and one antiquark) is referred to as the Liquid-Gas phase transition, while the one from hadronic to deconfined quark matter is referred to as deconfinement. Both are expected to be first-order phase transitions in the low-temperature and high baryon density () regime and present a crossover region beyond a critical point [7, 8] (see Figure 1).
The thermodynamical description of equilibrated matter is done through the equation of state (EoS), usually given as the relation between pressure and energy density. The dimensionality of the complete EoS depends on the characteristics of the system being described, such as temperature, number of conserved charges, and other effects (e.g. magnetic fields, spin, etc). In the case of QCD, the conserved charges typically considered are baryon number (), electric charge (), and strangeness (). In principle, all quark flavors should be conserved on QCD time scales but not all quarks are produced in enough abundance to be considered in chemical equilibrium for the EoS (although some studies have considered charm [9]). In this work, the dimensionality of our equation state is 3D , where is the respective (number) density associated with the conserved charge . We plan to add finite temperature in future works.
Furthermore, different quantities can be conserved globally or locally. Changing this assumption can reduce the dimensionality of the EoS, but this is not always a completely accurate assumption. For example, an electrically neutral system could contain 10 protons and 10 electrons. That system could be distributed such that the protons and electrons are paired close enough in phase space that locally it appears that there is no net-electric charge. However, it is also possible that all the protons are clumped together and the electrons are clumped together. In the case of a clumped-up charge, the system may look more like a dipole and locally net-electric charge is not zero, even though globally the system is electrically neutral. When describing multiple phases, conservation of specific quantities can be applied either to each phase separately, or allowing mixtures of phases, see Ref. [10] and references therein. In this work, we do not impose conservation of any quantity, but instead freely vary the chemical potentials associated with the conservation of such that we can vary in the phase space of .
To describe fully evolved (beyond the protoneutron-star stage) cold neutron-star matter, one assumes charged neutrality and -equilibrium with leptons. At -equilibrium, the charge chemical potential is related to the electron and muon chemical potentials via . Electric charge neutrality is enforced i.e. , where stands for all particles involved, for electric charge, and for (number) density of each particle. The time scales associated with neutron stars and their mergers allow for the creation of net strangeness through weak interactions, meaning that there is no strange chemical potential, . Applying -equilibrium and , reduces the dimensionality of the system from 3 dimensions of (4 if one includes ) into 1 dimension (or 2), meaning that it only depends on baryon chemical potential, (and ). That being said, numerical relativity simulations of, e.g., neutron-star mergers require a 3D EoS of typical temperature, baryon density, and electric charge fraction, i.e., . Additionally, out-of-equilibrium effects can be important (e.g. bulk viscosity) when there is a delay in the system to reach -equilibrium, such that information about the EoS out of equilibrium is also required [11].
On the other hand, the conserved charges in heavy-ion collisions are dictated by the choice of nuclei that are collided and the energy of the beam that collides. At high beam energies, the nuclei are extremely Lorentz contracted, so they appear as nearly 2D objects that pass through each other nearly instantaneously, dumping energy but not stopping, so in the initial state reproduces . As the beam energies are lowered, the nuclei are less Lorentz contracted, such that they become 3D objects that take a finite amount of time to pass through each other. Due to this longer timescale, there is enough time for baryons to be stopped in the initial state, such that collisions present a finite baryon number i.e. . The ions themselves have a specific number of protons and nucleons , such that one can define the initial charge fraction . Since both electric charge and baryon number are exactly conserved globally within heavy-ion collisions, then is determined by the choice of ions collided and how many (and which type) of baryons are stopped in an individual collision. The colliding nuclei do not contain net strangeness, however, due to gluon splittings into quark anti-quark pairs, strangeness is produced over time in heavy-ion collision while preserving strangeness neutrality, i.e., , where is the particle strangeness. This results in a non-zero when due to strange baryons and antibaryons. The time scales are short enough that the system cannot undergo weak decays, so the flavor of the strange quarks is preserved in the system (although quark-antiquark pairs can be annihilated or produced from gluons pairs). Experiments can measure strange mesons and baryons and find that approximately of hadrons produced have non-zero strangeness (predominately kaons) at mid to high energy collisions, see e.g. [12, 13]. Thus, one typically reduces the dimensionality of heavy-ion collisions from a 4D phase space into 2D () because the strange and electric charge chemical potentials become functions of i.e. and .
Note that the temperature of heavy-ion collisions is always non-negligible, even at some of the lowest beam energies (estimates from HADES suggest an average temperature of MeV [14]). However, the EoS limit is still interesting to study as an input for theoretical models (see e.g. [15, 16] for transport simulations where the temperature is introduced through kinetic contributions and how they connect to neutron stars [17]). Heavy-ion collisions are close to the limit of symmetric nuclear matter (SNM) where (or ), although data exist for nuclei with a range of . In the limit of symmetric nuclear matter, there is exactly the same number of protons and neutrons in the colliding nuclei. At the limit of SNM, there are no antiparticles, meaning that in this special case, there cannot be strange particles as well and becomes irrelevant.
In this work, we do not impose neutron star, protoneutron star (almost isospin symmetric, but charge neutral and with ), neutron-star merger, nor heavy-ion collision conditions. Rather, we study the much more general full 3-dimensional (, , ) space assuming that the temperature is low enough compared to the chemical potentials that we can approximate . While the conditions we discuss above are relevant for equilibrium physics they are not the only type of physics that plays a role in these systems. For instance, in neutron star mergers the system may have some delay to return to -equilibrium, such that the EoS out-of--equilibrium is relevant to describe bulk viscosity effects [11]. In heavy-ion collisions, local fluctuations of baryon, strangeness, and electric charge can play a role at the Large Hadron Collider (LHC) due to gluons splitting into quark anti-quark pairs and also at the Relativistic Heavy-Ion Collider Beam Energy Scan (RHIC BES) due to fluctuations in the position of the protons and neutrons in the initial state [18]. In these examples one cannot simply reduce the EoS down to 2-dimensions because information about the full 3D () or 4D () is required to understand local fluctuations of charges, see e.g. [19]. Thus, our work is an important first step in the direction of eventually developing 3D, 4D, and 5D (when including additional magnetic field, ) equations of state needed for these simulations.
While the Lagrangian of QCD is well-known, solving QCD is far from easy. The most common approach that has been extremely successful is lattice QCD, which represents space-time as a crystalline lattice with quarks at vertices connected by lines where the gluons travel [20]. In the limit of small vertex separation, this approach reaches the true continuum theory of QCD. However, lattice QCD calculations can only be performed at vanishing densities due to the fermion sign problem (otherwise, it exhibits the sign problem when trying to integrate highly oscillatory functions [21, 22]). In order to circumvent the fermion sign problem, it is possible to perform calculations at either imaginary chemical potentials or perturb around vanishing chemical potentials in order to obtain susceptibilities, , of the EoS. Then these susceptibilities can be used to recreate the EoS in up to 4D () through a Taylor series expanded in , where [23, 24]. Given that lattice QCD results are only available at temperatures of MeV and the expansion is currently only valid up to with the furthest reaching expansion scheme [25], there is only a limited regime where lattice QCD results can be applied. For this reason, lattice QCD cannot be used to describe neutron stars, where MeV at vanishing temperatures (in the MeV scale).
Due to the limitations of lattice QCD, several effective approaches have been developed in the regime of physics relevant to heavy-ion collisions. One such example is based on a bottom-up approach from holography [26, 27, 28]. Since the initial papers, the holographic approach has been significantly improved, such that one can tune its description to the latest lattice QCD results and predict the location of the QCD critical point [29]. Other approaches have found the QCD critical point in a similar location as well [30, 31, 32, 33, 34]. Thus, the EoS from heavy-ion collisions is beginning to converge at finite (concerning certain important features), going beyond the current regime of validity for lattice QCD. Still, the EoS at finite is still not known precisely at this time, especially when one considers effects that are off the strangeness neutral trajectory (see [19] for missing regimes in the EoS given the current lattice QCD results).
Outside the region covered by lattice QCD, perturbative QCD (pQCD) calculations are possible at extremely large and/or large . They are performed using a perturbative expansion in the QCD coupling constant, which is small in this regime [35, 36, 37]. However, near the deconfinement phase transition, the QCD coupling constant becomes large and the truncated sum from perturbation theory no longer approximates the infinite sum. On the other hand, chiral effective theory (EFT) calculations are possible at nearly vanishing temperatures and baryon densities around nuclear saturation density [38, 39]. They include every allowed particle interaction and order them by the number of powers of mass and momentum [40]. However, even with the combination of lattice QCD, pQCD, and EFT the vast majority of the QCD phase diagram is not yet possible to map out from first principle calculations (see Fig.1 of [41]).
As a result, one must turn to effective models, utilizing phenomenological constraints to construct Lagrangians that contain the right degrees of freedom (quarks at high , hadrons at intermediate , and nuclei at very low ) to describe the entire space of 4D or 5D phase diagrams. These effective models should smoothly connect to first principle QCD calculations in their regime of validity and should also include known particles, their masses, and their known interactions.
At low , the smooth (crossover) deconfinement to quark matter approximately coincides with the restoration of chiral symmetry. The spontaneous breaking of chiral symmetry (related to spin handedness) is what gives baryons approximately of their masses, with a small bare mass produced via the Higgs mechanism [42]. Spontaneous refers to the fact that the physical state of the system may be asymmetric even though the underlying physical laws are symmetric. This can be achieved by using a description in which hadronic masses are generated by interactions with the medium, and depend on density and/or temperature. Additional explicit symmetry-breaking terms ensure that pseudo-scalar mesons such as pions (the Goldstone bosons of the theory) have small but finite masses. Chiral symmetry breaking is also related to the formation of scalar condensates, which can for this reason be used as order parameters for this symmetry. These condensates are associated with scalar mesons that mediate the attraction between hadrons, while vector mesons mediate the repulsion between hadrons. Effective (chiral) models include the Nambu-Jona-Lasinio (NJL) model, the linear sigma model, and the parity doublet model, all of which account for chiral symmetry but have no mechanism to describe confinement [43, 44, 45, 46].
In particular, the Chiral mean-field (CMF) model is a very successful relativistic approach based on a non-linear realization of chiral symmetry [47, 48], which allows for a very good agreement with experimental nuclear data. In addition, only the mean values of the mesons are used in the CMF model, as the mesonic field fluctuations are expected to be small at high densities. As an effective model, once it is calibrated to work on a certain regime of energies, it can produce reliable results for the matter EoS and particle compositions, which can then be used in e.g, hydrodynamical simulations of heavy-ion collisions [49, 50] or astrophysics [51, 52, 53, 54, 55, 56], including core-collapse supernova explosions [57], stellar cooling [58, 59], and compact star mergers [5, 60, 61]. See Ref. [62, 63, 64] for a recent reviews that compare the CMF with other multidimensional models available in the CompOSE repository [65, 66].
The CMF model can be applied at as well as intermediate, and larger temperatures. It has also been extended to include the effects of strong magnetic fields [67, 68, 69, 70]. It includes degrees of freedom (d.o.f) expected to appear in different laboratory and astrophysical scenarios (leptons, baryons, mesons, and quarks) within a single framework. Both isospin asymmetry and strangeness (from hyperons and strange quarks) are included in the formalism, in order to describe the different environments. Inspired by unified approaches for the nuclear liquid-gas phase transition (between a phase with nuclei and a bulk hadronic one) [71], a unified approach for quark deconfinement was implemented in the CMF model [72], as explained below. All degrees of freedom are a priori included in the description, allowing deconfinement to appear as a first-order phase transition or crossover (in this case referring to higher than first-order phase transition), as expected at large temperatures [73]). The transition from baryons to quarks as the density and temperature increase is done utilizing an order parameter, a scalar field named in analogy with the Polyakov loop [74, 75]. This order parameter is introduced in the effective mass of baryons and quarks. Within this approach, full QCD phase diagrams can be built, showing both the liquid gas and deconfinement phase transitions [72, 76, 10, 77, 78, 79].
The CMF model in its current form has been used for over 2 decades. However, the software developed for those calculations written in Fortran 77 was inefficient, not well-documented, proprietary, and had most variables hard-coded. Thus, the legacy CMF software placed many numerical limitations on the CMF model. For instance, while the theory allows for 4D or 5D calculations of the EoS, the legacy version of the CMF model was only calculated in maximum 3D [78] due to the very long run time.
In this paper, we report on a brand-new version of the CMF model in C++ in collaboration with computer scientists through the MUSES collaboration (Modular Unified Solver of the Equation of State [80]) that increases the efficiency of the code by several orders of magnitude. It also allows for more accurate solutions, such that high-order derivatives of the EoS are now possible for the first time, and captures not just the stable region of the phase diagram but also the metastable and unstable regions across first-order phase transitions. For this work, we consider the vanishing temperature limit of the CMF model () and no magnetic fields (). However, future work is ongoing to extend the C++ version of the CMF model both to finite and .
The paper is organized as follows: in Sec. II the theoretical development of the CMF model is outlined. First, the full chiral Lagrangian is built in Sec. II.1, then the mean-field Lagrangian is obtained in Sec. II.2, followed by the equations of motion in Sec. II.3, the thermodynamical observables in Sec. II.4, and the coupling constants used within CMF in Sec. II.5. The numerical implementation of the CMF model in C++ is outlined including a discussion of thermodynamical stability in Sec. III, and the corresponding benchmark tests from this new code are presented. Finally, the results from the upgraded version of CMF are shown in Sec. IV for different chemical potential combinations using different couplings. High-order derivatives of the EoS known as the susceptibilities are calculated and first-order phase transitions are explored, including the discussion of metastable and unstable regions. In Sec. V concluding remarks and future plans are discussed.
II CMF Formalism
This section outlines the equations used to calculate the thermodynamical properties of bulk hadronic and quark matter within the CMF model. For the first time, we show in detail in this paper the derivation of the entire mean-field Lagrangian, equations of motion, and the thermodynamic properties. Due to the large densities found in compact stars, the particles in their interior become relativistic, each with their momentum comparable to their mass. And so, a relativistic space-time metric must be adopted to describe them. The CMF model is relativistic, therefore, it respects causality, provided that the repulsive vector interactions are not too strong (see the footnote in [81]). Following the literature of relativistic models, we make use of natural units throughout the paper.
The CMF model is based on a non-linear realization of the SU(3) sigma model [48] in which hadrons interact by mesonic exchange: , , , , , and . The scalar-isoscalar field corresponds loosely to the light quark composed meson and is the main driver for chiral symmetry restoration. A strange scalar-isoscalar field () with hidden strangeness (which is assumed to couple with strange particles) is also introduced to provide necessary attraction and is associated with the meson [82]. The scalar-isovector field couples differently to particles with different isospin, introducing a mass splitting between isospin multiples and making the EoS sensitive to isospin. It is associated with the meson [83, 84]. These fields mediate interactions among nucleons, hyperons, and quarks, contributing to attractive medium-range forces (scalar fields) and short-range repulsion (vector fields: vector-isoscalar , strange vector-isoscalar , and vector-isovector ). The scalar dilaton field, , representing the hypothetical glue ball field, is introduced to replicate QCD’s trace anomaly [48].
While in reality the strong force is propagated by gluons, the CMF model approximates this interaction as an exchange of mesons. They are not the typical particle physics mesons, such as pions or kaons (bound states of a quark and an antiquark), instead they are virtual particles that serve as force carriers for the strong force, much like how the photon is the force carrier of the electromagnetic force. Unlike electromagnetism, the strong force changes sign (and therefore whether it is attractive or repulsive) based on the separation between particles. At low , mesons do not contribute kinetically.
The CMF model is built in a chirally invariant way, as the masses of the particles are built from interactions with the medium and, as a result, the masses decrease with the energy. Note that the commonly used linear sigma model with linear realization approach in meson-baryon coupling leads to imbalanced hyperon potentials due to symmetric spin-0 and antisymmetric spin-1 meson interactions, and additional attraction from the field without counterbalancing repulsion. Moreover, explicit symmetry-breaking terms cannot correct these potentials without disrupting partially conserved axial current relations. The non-linear realization, incorporating pseudoscalar mesons as angular parameters of chiral transformation, allows explicit symmetry-breaking terms to be added without affecting partially conserved axial current relations and decouples strange and non-strange condensates, ensuring a balanced interaction that gives correct hyperon potentials [48]. While in the linear sigma model, the different left- and right-handed chirality wave functions transform differently within the SU(3) SU(3)R group, in the nonlinear realization we apply a transformation to the left- and right-handed chirality wave functions that allow them to transform in the same way, see Refs. [85, 86, 87] for more details.
The strength of the (confining) strong force changes with momentum transfer between particles, where the strong force becomes weaker with increased momentum transfer. This means that for high energies, temperature, or chemical potential, quarks become effectively deconfined [88, 89]. For this reason, quarks are also included in the CMF model (within the same description) but with different couplings [72]. The different phases, hadronic and quark, are characterized and distinguished from each other by the values of the condensates (such as ) and the order parameter for deconfinement, .
Although there are six known “flavors” of quarks, effectively, only up, down, and strange quarks are present in the energy regime we are discussing in this work (given in Table 1). The gluons serve to carry both the attractive and repulsive attributes of the strong force, but in the CMF model, these attributes are split between scalar (spin-) mesons mediating attractive interactions and vector (spin-) mesons mediating repulsive interactions. The baryons included in the CMF model are the baryonic octet (Table 2) and the decuplet (Table 3). An alternative version of the CMF model also exists that includes the chiral partners of the baryons, see Ref. [76, 56] but this approach is not included in CMF++.
| Quark | Symbol | Mass | Electric Charge | Isospinz | Strangeness |
|---|---|---|---|---|---|
| (MeV) | () | () | |||
| up | |||||
| down | |||||
| strange |
| Symbol | Valence | Mass | Electric Charge | Isospinz | Strangeness |
|---|---|---|---|---|---|
| Quarks | (MeV) | () | () | ||
| Symbol | Valence | Mass | Electric Charge | Isospinz | Strangeness |
|---|---|---|---|---|---|
| Quarks | (MeV) | () | () | ||
| 0 | |||||
| -1 | |||||
| 0 |
The construction of the CMF model is described in detail in the following subsections, however, the general procedure is shown in Figure 2. One develops an effective theory by determining the relevant symmetry group and then constructing the appropriate Lagrangian that contains all the particles and their interactions. Once the Lagrangian is established, then the mean-field approximation is made to simplify the Lagrangian to a form that can be solved straightforwardly. After the deconfinement mechanism is implemented, the model is named CMF. Next, the equations of motion are obtained from the Euler-Lagrange equation and the ideal fluid approximation is assumed, such that one can diagonalize the energy-momentum tensor. At this point, input from experimental and theoretical constraints for the model parameters are applied (e.g. mass and couplings of the baryons, etc.). Then, for a given set of chemical potentials , , (and ), the equations of motion for the mesons can be used to determine the particle population and to calculate the thermodynamic observables that allow one to obtain a multidimensional EoS.
II.1 Full chiral Lagrangian
For the non-linear realization of the sigma model, the full hadronic Lagrangian density reads [48]
| (1) |
where is the kinetic energy term, is the baryon-meson interaction term, is the scalar meson self-interaction term, is the vector meson self-interaction term, and is the term for explicit breaking of chiral symmetry. We now cover each of these five terms in more detail.
II.1.1 the kinetic-energy term
The kinetic energy term expands as
| (2) |
where is a covariant derivative defined by
| (3) |
with being any of the following particle matrix: stands for the baryon octet matrix, for the scalar meson nonet and for the pseudoscalar singlet. They are shown in the mean-field approximation in Appendix A. The represents the operator commutator and is a vector-type field that assures chiral invariance and is defined by
| (4) |
The kinetic energy term of the pseudoscalar mesons is introduced (in analogy to Eq. 4) by defining the axial vector,
| (5) |
with
| (6) |
where is the pseudoscalar octet matrix defined in LABEL:Pmatrix. are the Gell-mann matrices, the fifth Dirac gamma matrix, which is the chirality operator, and are the components of the pseudoscalar meson octet.
The vector and axial-vector field tensors are = and =, with the associated vector and axial field vectors and . The vector meson nonet is shown in Appendix A and represents the dilaton field, a.k.a. glueball field.
The first term in Eq. 2 is a Dirac term for the baryons, the second, fourth, and fifth terms are Klein-Gordon terms for their respective scalar, pseudo-scalar singlet, and dilaton fields, the sixth and seventh terms are Proca terms for the vector and axial-vector mesons, whereas the third term contains interaction between the scalar mesons and the pseudo-scalar meson nonet, including the pseudo-scalar kinetic term, .
II.1.2 the scalar-meson self-interaction term
The scalar meson self-interaction couplings are governed solely by SU(3)V symmetry, resulting in three lowest independent invariants,
| (7) |
For integer , are invariant but not independent, as they can be written in terms of , , and ; for example, it can be shown using the matrix from LABEL:Pmatrix that
| (8) |
where
| (9) |
Using these invariants (excluding linear terms in the scalar mesonic fields that would generate a non-zero scalar density in vacuum), we define the scalar Lagrangian density up to order as
| (10) |
where each term has been multiplied by an appropriate power of the dilaton field to allow the coupling constants to be dimensionless and thus make the Lagrangian scale invariant [91, 87]. The parameter is related to the nuclei-scalar meson interaction in the chiral model [48]. It is not considered in the CMF, as it currently does not include nuclei as degrees of freedom, .
Whenever there are remaining dimensionful terms in the Lagrangian, the dilaton field must be multiplied with appropriate power to keep the coupling constants dimensionless and the Lagrangian scale invariant. Additionally, to mimic the broken scale invariance property of QCD , a scale breaking term is added to the effective Lagrangian (with as a gluon field tensor) [92]
| (11) |
In analogy to the scale-breaking term discussed in Ref. [91], the first term is added to the effective Lagrangian at tree level, where represents a field associated with a spin glueball. This term disrupts scale invariance, resulting in the proportionality , which follows from the definition of scale transformations [93]. However, this form of the glueball potential is strictly applicable only to the effective low-energy theory of pure, quarkless QCD. To generalize the glueball potential for the case of massless quarks, a second term is introduced. Moreover, a third term is introduced to generate a phenomenological consistent finite vacuum expectation value [94]. The second and third terms extend the logarithmic term introduced in Ref. [95] within the context of SU(3), ensuring that holds. In Eq. 11, is the vacuum expectation value of the scalar matrix, is the vacuum expectation value of the dilaton field, and we set , which is related to the quark contribution to the QCD beta function [91]. Adding these two pieces together gives us the full scalar mesonic self-interaction term
| (12) |
II.1.3 the vector-meson interaction term
The vector-meson interaction term is
| (13) |
where the first term is the mass term of each vector meson . The second one presents different possibilities for the self-interaction terms of the vector mesons that are chiral invariant [53]
| (14) |
where C2 is a combination of other terms. Two more chiral invariant combinations can be used, but they were never studied in detail because they did not seem to produce physical results.
II.1.4 the explicit symmetry-breaking term
As previously discussed in Sec. I, when chiral symmetry is spontaneously broken, Goldstone bosons emerge, which leads to large fluctuations that can lead to divergences or instabilities in the model. To remove the effects of these fluctuations, we add explicit symmetry-breaking terms to the Lagrangian density which also give rise to pseudoscalar-meson mass terms,
| (16) |
where is the matrix of explicit symmetry breaking parameters [96], with and being the decay constants of pions and kaons. This term gives rise to a pion mass and leads to partially conserved axial current relations for and mesons. The choice of power for the dilaton field matches the dimension of the fields of the chiral condensates [92].
In contrast to linear realization, a symmetry-breaking term can be explicitly introduced in the nonlinear realization to accurately reproduce the hyperon potentials without impacting the partially conserved axial current relations [48]. We introduce the following term with a free parameter , which contributes to the bare mass of hyperons (with typical values given in Table 9)
| (17) |
where . The net Lagrangian for the explicit symmetry-breaking contribution then reads
| (18) |
II.1.5 the baryon-meson interaction term
This interaction term is similar for all mesons, with the only difference being the Lorentz space occupied by the mesons. Therefore we can write all the interactions with two compact terms, and any baryon () - meson () interaction expands as
| (19) |
where depends on the specific interaction, with values listed in Table 4, and are the coupling constants related to the singlet and octet (discussed in the following), and controls the mixing between the -type (symmetric) and -type (anti-symmetric) terms, that read
| (20) |
and
| (21) |
The last term in the -type interaction term is added to cancel out the singlet contribution to the octet interaction when a nonet meson matrix is utilized.
| Interaction with baryons | ||
|---|---|---|
| Scalar | ||
| Pseudo-scalar | ||
| Vector (vector interaction) | ||
| Vector (tensor interaction) | ||
| Axial-vector |
II.2 The mean-field Lagrangian
The full quantum operator fields in the Lagrangian (Eq. 1) lead to nonlinear quantum field equations with large couplings, making perturbative approaches infeasible and challenging to solve. Hence, reliable non-perturbative approximations are essential for solving these complex many-body interactions and achieving accurate comparisons between theory and experiment [97]. To describe dense matter, we apply the mean-field approximation, as first proposed in Ref. [97]. Within the mean-field approximation, we assume homogeneous and isotropic infinite baryonic matter with defined parity (+) and charge (0). Thus, only mean-field mesons with positive parity (scalar mesons and time-like component of vector mesons) and zero third component of isospin (mesons along the diagonal of the matrices (see Eq. 93) and (see Eq. 94)) are non-vanishing. The mean-field mesons with negative parity (space-like component of vector mesons, time-like component of axial-vector mesons, and pseudoscalars) do not follow parity conservation, and there is no source term for them in mean-field infinite baryonic matter111The ground state expectation value of space-like components of axial-vector mesons, despite their positive parity is zero because of the homogeneous and isotropic medium assumption..
Furthermore, in this approximation, fluctuations around the constant ground state expectation values of the scalar and vector field operators are neglected, for example,
| (22) |
As a consequence, , , , , , and are all reduced to time-and space-independent quantities. For simplicity, we omit the time index of vector mesons (, , and ) and also omit the third component of the isospin index of the isovector mesons ( and ). We refer to the resulting Lagrangian, also including quark degrees of freedom [72], as the Chiral Mean-Field (CMF) Lagrangian density
| (23) |
In the regimes we examine, the field has a weak coupling to the baryons, resulting in little overall contribution to the baryon thermodynamic quantities regardless of the value of the field. Thus, for the remainder of this work, we set ( remains “frozen” at its vacuum value, ) and apply further simplifications. For details regarding , see Ref. [91, 98, 99]. As a result, no equation of motion for the field is shown.
II.2.1 the kinetic-energy term
In the mean-field approximation, , thus the commutator meaning that the covariant derivative reduces to the partial derivative, . The mesons are taken as static, and thus no longer have kinetic terms (, where is the matrix from Table 4), such that all of Eq. 2 reduces to
| (24) |
Quarks are discussed in the following.
II.2.2 the scalar-meson self-interaction term
Applying the mean-field approximation to the scalar-meson self-interaction term and calculating the terms explicitly (see Section B.1 for a detailed calculation), we obtain
| (25) |
where and are the vacuum values of the and fields, respectively.
II.2.3 the vector-meson interaction term
After applying the mean-field approximation, the total vector-meson Lagrangian can be expressed as a sum of the mass and self-interaction terms
| (26) |
or explicitly
| (27) |
The coupling scheme is a linear combination of and and which exhibits no mixing. The coupling scheme denoted as for the self-interaction of vector mesons is quite different from the other schemes, as it includes a term that exhibits a linear dependence on the strange vector meson . Because of this linear dependence on , the scheme requires a different parametrization that includes a bare mass term, Eq. 15 to ensure that the compressibility of nucleons is in a better agreement with nuclear physics data [100, 101]. See Ref. [102] for combinations of the couplings C1-C4 that allow one to separate each coupling term (in the non-strange case, ).
II.2.4 the explicit symmetry-breaking term
In the mean-field approximation, the first explicit symmetry-breaking term Eq. 16 (together with Eq. 93) simplifies to
| (28) |
From Eq. 17, the expanded form of symmetry-breaking Lagrangian related to the hyperon () potential is given by
| (29) |
It leads to an additional contribution to the coupling between the hyperons and the mesons and through the parameter , and to a constant (bare) mass term .
Note that also receives a contribution from the bare mass term (in the case of C4 vector-coupling), Eq. 15 as follows
| (30) | ||||
| (31) |
Similarly, a mass correction due to an explicit breaking term with parameter for the baryon decuplet () is written as
| (32) | ||||
| (33) | ||||
| (34) |
II.2.5 the baryon-meson interaction term
Applying the mean-field approximation to the baryon-meson interaction term, we only get non-zero values for the cases where and . This is due to the -matrix for pseudovector mesons having vanishing expectation values and the pseudoscalar mesons only coupling to the baryons with a pseudovector coupling. By doing the explicit calculation of the interaction term (see Section B.2), we can rewrite the interaction Lagrangian, Eq. 19, as
| (35) |
where the couplings are written in terms of (from Eq. 19), , , and , as shown in Table 5 for the scalar-mesons (). We can identify the effective mass terms for the baryons in terms of these as
| Particle | |||
| 0 | |||
| 0 | |||
| 0 | |||
| 0 | |||
| 0 | |||
| 0 |
| (36) | ||||
which can be written compactly as
| (37) |
It must be noted that an additional contribution to the effective mass must be accounted for when the deconfinement order parameter is introduced in the next section. If we disregard the -meson contribution, the baryons masses of the nucleon doublet and hyperon triplets are degenerate. The inclusion of the isovector meson breaks this multiplet mass equality.
For the baryon decuplet , we follow [94] and assume they are described by Dirac spinors such that, from the interactions between the baryon resonances and scalar mesons, we may extract the effective mass terms for the isospin degenerate baryon decuplet
| (38) |
Similarly to what has been done in Section B.2 for the baryon-scalar meson coupling constants, we can calculate the baryon-vector meson coupling constants . Based on the vector dominance model (VDM) and the universality principle, it can be inferred that the -type coupling is likely to be minimal [103]. Therefore, in our analysis, we employ only -type coupling by choosing for all fits. Additionally, we can decouple the nucleons from the strange vector meson by setting such that . Following a similar pattern, we assign =0, resulting in the absence of coupling between the and the baryons. The remaining couplings to the strange baryons are subsequently determined by symmetry relations (the quark model) [104] in terms of (the only free parameter for the baryon-vector mesonic coupling), such that the and -meson couplings are given in Table 6. Note that the -meson couplings follow the sign convention of the -meson. The scheme described is known as -type or SU(6) as it includes SU(3) flavor symmetry and SU(2) spin symmetry [105, 106]. Nevertheless, in the CMF model we break this scheme and use, e.g. for C4 (instead of ) to slightly modify allowing a better fit of experimental data for the symmetry energy (as small differences matter). Moreover, a parameter called is introduced in the decuplet baryons’ vector coupling () allowing a better fit of experimental data for the -nucleon potential. More general couplings will be explored in the future.
II.2.6 adding quarks to the model
To reproduce quark deconfinement, we include up, down, and strange quarks in the model. We assume the same Lagrangian as the baryonic one, with kinetic, mass, and interaction terms given by
| (39) |
with and masses MeV for up and down quarks and MeV for the strange quark. We write the effective quark mass like the baryonic one, Eq. 37, by defining
| (40) |
CMF parameters associated with the quark sector have scalar couplings that are set to be roughly one-third of the nucleon scalar couplings, while the vector couplings are set to zero, in agreement with the findings of Ref. [107]. The coupling values are discussed in Sec. II.5.
II.2.7 the deconfinement order-parameter potential
To obtain a unified quark-hadron EoS, we implement a Polyakov-inspired potential term (referred to as the deconfinement potential) of the form [72]
| (41) |
which at reduces to
| (42) |
where the ’s and are constants. Here we introduced a scalar field , which serves as an order parameter for the quark-hadron phase transition. was modified from its original form in the PNJL model [108, 109] to also contain baryon chemical potential dependent terms (of even order), to be used to study low-temperature and high-density environments, such as neutron stars. It has been shown that a term in (instead of ) would significantly weaken the deconfinement phase transition at [81, 110, 111]. The form of dictates the shape and location of the quark-hadron phase transition in the QCD phase diagram. If future information from the RHIC Beam Energy Scan and theoretical developments further constrain the QCD critical point, one could redefine to reproduce these new constraints.
This bosonic scalar field also appears in an additional contribution to the effective masses of the baryons.
| (43) |
Similarly, the quark effective masses have the form
| (44) |
Considering to be large, quark masses are large and baryon masses are small when (and vice-versa when ) The larger the mass of a particle, the more energy is required to create it. Therefore, when is large it causes the baryon masses to be so large that it suppresses their influence and one is in a quark-dominated phase. On the other hand, when is small then the quark masses are large such that they are suppressed and one is in a hadron-dominated phase. Putting this all together, corresponds to having only hadrons and corresponds to having only quarks, with intermediate values corresponding to having both hadrons and deconfined quarks (only reproduced at large temperatures).
II.3 Equations of motion
To derive the equations of motion for the seven bosons (the six mean-field mesons and the deconfinement order parameter, ), we apply the Euler-Lagrange equation to the CMF Lagrangian density
| (45) |
with , and , resulting in
| (46) |
where is the scalar number density and is the baryon (vector) number density. The index always indicates a summation of the baryon octet, decuplet, and quark flavors.
Note that we do not derive equations of motion for fermions, as the expected value of their fields does not come from their equations of motion in our formalism, but instead directly from their effective chemical potentials and effective masses, which come from the chemical potentials looped over, . This is discussed in detail in the following.
II.4 Thermodynamical observables
The CMF Lagrangian density can be alternatively seen as consisting of a fermion part, a boson part, and a vector interaction term . For fermions, the Lagrangian density reads
| (47) |
where the scalar-meson interactions are hiding within the effective mass . This is the relativistic free Fermi gas Lagrangian but with effective masses (see Appendix C). The vector-meson interaction term is
| (48) |
and it leads to an effective chemical potential for the fermions, given by
| (49) |
The individual particle chemical potentials are given by
| (50) |
where , the particle baryon number, is for baryons and for quarks.
Within our formalism, bosons do not acquire effective masses. There is also no contribution from the kinetic term, resulting in the bosonic Lagrangian as
| (51) |
where . The total energy density, pressure, vector or baryon (number) density, and scalar density include the sum of contributions from individual fermions
| (52) |
with the individual particle contributions being calculated from the energy momentum-tensor as shown in Appendix C. Here, is the number density of the fermions without the scalar field contribution and, therefore, is different from the baryon density defined as (containing a contribution) and discussed in the following. ensures that each quark counts as of a baryon.
At vanishing temperature, the entropy density is identically zero in our framework. Then the thermodynamic variables can be calculated directly using:
| (53) |
| (54) |
| (55) |
| (56) |
where is the total degeneracy (spin and color), is the Fermi momentum of particle , and at we can also write the effective chemical potential as the effective energy level
| (57) |
The asterisks represent the influence of the strong interaction. For a given set of , and , once one determines the effective particle chemical potentials Eq. 49 and masses Eqs. 43 and 44 (solving for the mean fields), at the Fermi momenta and thermodynamical properties easily follow.
The baryon-vector meson interactions modify the solution of the Dirac equation, Eq. 135, by modifying its energy in the plane wave exponential as . Then, the derivative of the Dirac spinor in Eq. 139, applied to Eq. 135 leads to a contribution to the energy density as
| (58) |
We note here that the terms of the form only contribute to the energy density due to being a simplification of . The pressure, on the other hand, does not receive extra contributions from the mean-field mesons, since their spatial components are taken to be zero (see Eq. 139).
Unlike the mesons, has explicit temperature and chemical potential dependence (see Eq. 41). This means that to satisfy thermodynamic consistency, we must have (with having no electric charge or strangeness), with and For the mesonic contribution to the thermodynamic quantities, we start from (Eq. 136) and acknowledge that for all mesons and due to the mean-field approximation. This gives the energy density and pressure as
| (59) | ||||
| (60) |
Furthermore, in a vacuum, all thermodynamic quantities should be zero. This is not the case for the scalar meson contribution, which acts as self-energy. However, their vacuum values are constant, and we can always add a constant to the Lagrangian density. Accordingly, we make a final alteration to the CMF Lagrangian density by subtracting the constant vacuum state. We do not have the same issue with fermions because their thermodynamic variables are already zero in the vacuum. The net result is
| (61) |
where
| (62) |
is the aforementioned constant vacuum state value, achieved by taking and with all other meson fields vanishing, keeping in mind that we already have and .
Thus, we can write the expressions for the thermodynamic quantities (including mesonic and contributions)
| (63) |
where the last term is calculated analytically and would represent some sort of gluonic interaction contributing to the baryon density. Finally, we get the total contribution to the thermodynamic quantities by adding the fermion and boson contributions together
| (64) |
Throughout this paper, we often refer to fractions instead of densities. They are defined as ratios of sums of quantum numbers (weighted by the number density of particles)
| (65) |
for the charge fraction and
| (66) |
for the strangeness fraction. Since the field possesses no quantum number, it does not contribute to these fractions.
II.5 Coupling constants
| Parameter | Interaction | Constraint |
|---|---|---|
| , , | ||
| 222For the decuplet the correct approach is to use the Rarita-Schwinger equation, which can be written as a Dirac equation with extra constraints [112]. Here, we simply follow the results presented in Ref. [94]. | , , | |
| splitting | ||
| one-loop function | ||
| , | ||
| , | ||
| , | ||
| , , | , | , , |
| , , , - mixing (VDM) | ||
| 11footnotemark: 1 | ||
| , , | ||
| 11footnotemark: 1 | ||
| 11footnotemark: 1 |
In Table 7, we list the free parameters of the CMF model and the corresponding constraints used to fix them. Note that the CMF couplings are constant, e.g., are not dependent on quantities like density. In the first set of rows, we present the scalar coupling constants concerning the interaction between scalar mesons and baryon octet (decuplet), which are determined based on the vacuum masses of baryons. For the octet, the following coupling relations are obtained through Eq. 36 in the vacuum
| (67) |
| (68) |
| (69) |
and for the decuplet, the coupling constants are obtained through Eq. 38 in the vacuum
| (70) |
| (71) |
Additionally, parameters like , , and governing scalar self-interactions are adjusted to the Lagrangian minima for , , and in the vacuum, while and are tuned to match the vacuum masses of (which is uncertain) and the , splitting, respectively. The parameter , linked to scalar scale breaking Lagrangian, is calibrated to the one-loop QCD beta function. Moreover, the vacuum value of is set to reproduce zero pressure at saturation. The vacuum value of the scalar field is fitted to the decay constant of the meson, while is fitted to the decay constants of the and mesons. See Table 8 for a complete list.
The vector-baryon coupling constants have been fitted to reproduce nuclear saturation properties for isospin-symmetric matter and asymmetric matter, together with neutron-star observations. This includes , , and the bare mass of baryons fitting simultaneously saturation density fm-3 and binding energy per nucleon MeV (which results in compressibility of MeV) and the asymmetry energy at saturation MeV (by using ) producing a slope MeV) separately for all for the vector couplings (C1-C4).
There is also a requirement to reproduce M⊙ stars with radii consistent with observations. Reproducing these values requires a setting of vector coupling constants given in Table 9. The remaining baryon-vector-meson coupling constants relate to the value of associated with . Non-strange particles do not couple to and . Finally, parameter relates to experimental vector meson vacuum masses.
We fit to reproduce reasonable hyperon potentials () [113] for symmetric matter at saturation, in particular MeV (reproducing MeV and MeV). We fit to reproduce a reasonable baryon potential for symmetric matter at saturation, MeV (similar to the nucleon one MeV). This procedure is done separately for each of the couplings C1-C4. We use a fixed value for , since there is little data available for the strange members of the baryon decuplet. Additionally, a full list of constants shared among coupling schemes is provided in Table 8, and a list of constants that are different in different coupling schemes is provided in Table 9.
| MeV | MeV | |
| MeV | MeV | MeV |
| MeV | ||
| MeV | ||
| MeV | MeV | MeV |
| Coupling | |||||||
|---|---|---|---|---|---|---|---|
| C1 | 58.40 | 13.66 | 11.06 | -11.06 | 0 | 1.24 | 1.07 |
| C2 | 58.40 | 13.66 | 3.51 | -3.51 | 0 | 1.24 | 1.07 |
| C3 | 58.40 | 13.66 | 3.82 | -3.82 | 0 | 1.24 | 1.07 |
| C4 | 38.90 | 11.90 | 4.03 | -4.03 | 150 | 0.86 | 1.2 |
| MeV | MeV | MeV |
|---|---|---|
| =200 MeV | MeV |
Following this, we detail the parameters related to the deconfinement potential (not including the decuplet) and quarks. For the C4 coupling scheme, the quark and coupling constants (listed in Table 10) have been fitted to reproduce lattice results at zero and small chemical potential and known physics of the phase diagram. Lattice QCD predicts the first-order deconfinement phase transition (for pure glue Yang-Mills) observed at a temperature of MeV [109]. At , we fit the parameter and together to as well as the pressure function which mirrors patterns seen in previous works ([108, 109]) for pure glue Yang-Mills theories. At vanishing chemical potential, when including fermions, the hadron to quark phase change is a crossover rather than a sharp transition. The mid value of the crossover band is known as the pseudo-critical temperature of chiral symmetry restoration marked by a transition temperature . In the CMF model, this temperature is identified through the peak change in the condensate and field . The parameters and (coupling between quarks and ) are fitted together to reproduce =171 MeV in agreement with results from 2001 [114].
Furthermore, is fitted to the critical number density ( ) at the onset of deconfinement transition at for neutron stars, and is constrained by the critical temperature (=167 MeV) and critical baryon chemical potential (=354 MeV) for isospin symmetric matter, aligned with findings from 2004 [115]. Additionally, is tuned to maintain value within 0 and 1. It is noteworthy that parameters from the Polyakov-inspired potential (Eq. 42) and quark couplings (Eq. 44) have been fitted solely for C4, determining the location of the deconfinement phase transition at specific and EoS behavior in the quark regime. Since the paper aims to compare C++ and Fortran solutions while also analyzing stability, we employ the quark sector parameters of C4 (refer to Table 10) for all other coupling schemes. Adjusting the parameters to C1-C3 coupling schemes would lead to shifts in the location of the deconfinement phase transition as well as the behavior of EoS post-deconfinement transition. See Ref. [79] for a recent work in which we broke the mass degeneracy of vector mesons in the CMF model using their field redefinition. This required us to fit the C1-C4 coupling schemes to the up-to-date constraints coming from lattice QCD, low-energy nuclear, and astrophysics.
III Code Implementation
III.1 Code Overview
Figure 3 shows the CMF++ code flowchart with the Python and C++ layers, where the shaded gray region highlights the C++ main driver routine. The code can be divided into three sections: input preprocessing, main algorithm, and output postprocessing. In the first section, yaml_preprocess.py validates the YAML input configuration file required for the main execution. Details about the validation procedure are illustrated in Sec. III.2. In the second section, the main routine is responsible for finding solutions for Eq. 46 and computing derived thermodynamic quantities (see Eq. 64) for the valid solutions found. More details about this section are described in Sec. III.3. The last section covers postprocessing and output. In the postprocessing section, the solutions found in the main algorithm are cleaned and classified as stable, metastable, or unstable. Additionally, the output is divided by the underlying degrees of freedom, i.e., quarks vs. baryons. The criteria and procedure for separating the solutions are detailed in Sec. III.4. Finally, the output adapters are called in postprocess.py, which transforms the final output files into either CSV or HDF5 format via the MUSES Porter library for the consumption of other MUSES modules. Details on how to run CMF++ can be found in Appendix D.
III.2 Input preprocessing
The only input required to execute the code is a YAML-formatted configuration file to ensure human and machine readability, which we named config.yaml. The YAML file contains all the computational options and the physical parameters required to run. The computational options detailed in Table 11 encompass the model hyperparameters like the name of the run, see Table 12. The file structure is detailed in the OpenAPI specifications for the model version.
The config.yaml file is processed by yaml_preprocess.py which validates it via the openapi-core library and flattens it for the ingestion of the main algorithm.
| Category | Variable | Value | Description |
| computational_parameters | run_name | default | name of the run |
| solution_resolution | 1.e-8 | resolution for mean-field solutions | |
| maximum_for_residues | 1.e-4 | threshold for solution residues | |
| production_run | true | Is this a production run? | |
| options | baryon_mass_coupling | 1 | baryon-meson coupling scheme |
| use_ideal_gas | false | use ideal gas? | |
| use_quarks | true | use quarks? | |
| use_octet | true | use baryon octet? | |
| use_decuplet | true | use baryon decuplet? | |
| use_pure_glue | false | use gluons only (no baryons nor quarks)? | |
| use_hyperons | true | are hyperons included? | |
| use_constant_sigma_mean_field | false | fix sigma mean-field to chosen value | |
| use_delta_mean_field | true | is delta mean-field included? | |
| use_Phi_order | true | use Polyakov-inspired potential? | |
| use_constant_Phi_order | false | fix Phi field value to chosen value | |
| vector_potential | 4 | vector coupling scheme C1-C4 | |
| use_default_vector_couplings | true | use default vector couplings? | |
| output_files | output_Lepton | true | create output file for Lepton module |
| output_debug | false | create output file for debugging | |
| output_flavor_equilibration | true | create output file for Flavor equilibration module | |
| output_format | CSV | create output files either in CSV or HDF5 format | |
| output_particle_properties | true | create output file for particle populations and properties | |
| chemical_optical_potentials | muB_begin | 900.0 | initial baryon chemical potential (MeV) |
| muB_end | 1800.0 | final baryon chemical potential (MeV) | |
| muB_step | 1.0 | step for baryon chemical potential (MeV) | |
| muS_begin | 0.0 | initial strange chemical potential (MeV) | |
| muS_end | 1.0 | final strange chemical potential (MeV) | |
| muS_step | 5.0 | step for strange chemical potential (MeV) | |
| muQ_begin | 0.0 | initial charge chemical potential (MeV) | |
| muQ_end | 1.0 | final charge chemical potential (MeV) | |
| muQ_step | 5.0 | step for charge chemical potential (MeV) | |
| mean_fields_and_Phi_field | sigma0_begin | -100.0 | initial mean-field (MeV) |
| sigma0_end | -10.0 | final mean-field (MeV) | |
| sigma0_step | 30.0 | step for mean-field (MeV) | |
| zeta0_begin | -110.0 | initial mean-field (MeV) | |
| zeta0_end | -40.0 | final mean-field (MeV) | |
| zeta0_step | 23.333 | step for mean-field (MeV) | |
| delta0_begin | 0.0 | initial mean-field (MeV) | |
| delta0_end | 1.0 | final mean-field (MeV) | |
| delta0_step | 10.0 | step for mean-field (MeV) | |
| omega0_begin | 0.0 | initial mean-field (MeV) | |
| omega0_end | 100.0 | final mean-field (MeV) | |
| omega0_step | 33.333 | step for mean-field (MeV) | |
| phi0_begin | -40.0 | initial mean-field (MeV) | |
| phi0_end | 0.0 | final mean-field (MeV) | |
| phi0_step | 13.333 | step for mean-field (MeV) | |
| rho0_begin | 0.0 | initial mean-field (MeV) | |
| rho0_end | 1.0 | final mean-field (MeV) | |
| rho0_step | 10.0 | step for mean-field (MeV) | |
| Phi0_begin | 0.0 | initial mean-field (MeV) | |
| Phi0_end | 0.9999 | final mean-field (MeV) | |
| Phi0_step | 0.333 | step for mean-field (MeV) |
| Category | Variable (Symbol) | Value | Description |
| physical_parameters | d_betaQCD () | 0.0606060606 | fit parameter for beta QCD function |
| f_K () | 122.0 | decay constant (MeV) | |
| f_pi () | 93.3000031 | decay constant (MeV) | |
| hbarc () | 197.3269804 | (MeV) | |
| chi_field_vacuum_value () | 401.933763 | vacuum value (MeV) | |
| Phi_order_optical_potential | a_1 () | -0.001443 | fit parameter for deconfinement phase transition |
| a_3 () | -0.396 | fit parameter to keep between 0 and 1 | |
| T0 (crossover) () | 200 | fit parameter for pseudo critical transition temperature (MeV) | |
| T0 (pureglue) () | 270 | fit parameter for deconfinement critical temperature (MeV) | |
| scalar_mean_field_equation | k_0 () | 2.37321880 | fit parameter to minimize scalar Lagrangian with respect to |
| k_1 () | 1.39999998 | fit parameter for mass of meson | |
| k_2 () | -5.54911336 | fit parameter to minimize scalar Lagrangian with respect to | |
| k_3 () | -2.65241888 | fit parameter to account splitting | |
| explicit_symmetry_breaking | m_3H () | 0.85914584 | fit parameter for potential of strange octet baryons |
| m_3D () | 1.25 | fit parameter for potential of strange decuplet baryons | |
| V_Delta () | 1.2 | fit parameter for potential of decuplet particles | |
| vector_nucleon_couplings | gN_omega () | 11.90 | Nucleon coupling to field |
| gN_rho () | 4.03 | Nucleon coupling to field | |
| g_4 () | 38.90 | Self-coupling of the vector mesons | |
| mean_field_vacuum_masses | omega_mean_field_vacuum_mass () | 780.562988 | mean-field vacuum mass (MeV) |
| phi_mean_field_vacuum_mass () | 1019. | mean-field vacuum mass (MeV) | |
| rho_mean_field_vacuum_mass () | 761.062988 | mean-field vacuum mass (MeV) | |
| quark_bare_masses | up_quark_bare_mass () | 5.0 | up quark bare mass (MeV) |
| down_quark_bare_mass () | 5.0 | down quark bare mass (MeV) | |
| strange_quark_bare_mass () | 150.0 | strange quark bare mass (MeV) | |
| vacuum_masses | Delta_vacuum_mass () | 1232. | vacuum mass (MeV) |
| Lambda_vacuum_mass () | 1115. | vacuum mass (MeV) | |
| Sigma_vacuum_mass () | 1202. | vacuum mass (MeV) | |
| Sigma_star_vacuum_mass () | 1385. | vacuum mass (MeV) | |
| Omega_vacuum_mass () | 1691. | vacuum mass (MeV) | |
| Kaon_vacuum_mass () | 498. | vacuum mass (MeV) | |
| Nucleon_vacuum_mass () | 937.242981 | Nucleon vacuum mass (MeV) | |
| Pion_vacuum_mass () | 139. | vacuum mass (MeV) | |
| mass0 () | 150. | Bare vacuum mass (MeV) | |
| quark_to_fields_couplings | gu_sigma () | -3.0 | up quark coupling for mean-field |
| gd_sigma () | -3.0 | down quark coupling for mean-field | |
| gs_sigma () | 0 | strange quark coupling for mean-field | |
| gu_zeta () | 0 | up quark coupling for mean-field | |
| gd_zeta () | 0 | down quark coupling for mean-field | |
| gs_zeta () | -3.0 | strange quark coupling for mean-field | |
| gu_delta () | 0.0 | up quark coupling for mean-field | |
| gd_delta () | 0.0 | down quark coupling for mean-field | |
| gs_delta () | 0.0 | strange quark coupling for mean-field | |
| gu_omega () | 0.0 | up quark coupling for mean-field | |
| gd_omega () | 0.0 | down quark coupling for mean-field | |
| gs_omega () | 0.0 | strange quark coupling for mean-field | |
| gu_phi () | 0.0 | up quark coupling for mean-field | |
| gd_phi () | 0.0 | down quark coupling for mean-field | |
| gs_phi () | 0.0 | strange quark coupling for mean-field | |
| gu_rho () | 0.0 | up quark coupling for mean-field | |
| gd_rho () | 0.0 | down quark coupling for mean-field | |
| gs_rho () | 0.0 | strange quark coupling for mean-field | |
| gq_Phi () | 500.0 | quark coupling for field (MeV) | |
| baryon_to_Phi_field_coupling | gbar_Phi () | 1500.0 | baryon coupling to field (MeV) |
III.3 Algorithm
In computational terms, the CMF model is a coupled system of nonlinear algebraic equations for the mean-field mesons , and field (see Eq. 46), therefore, a root solver algorithm is required. In our implementation, we adopted the numerical root solver fsolve [116], which is inspired by the fsolve function from MATLAB and is based on MINPACK [117, 118]. MINPACK is a Fortran library designed to solve systems of nonlinear equations by residual’s least-squares minimization employing a pseudo-Gauss-Newton algorithm in conjunction with gradient descent.
This validated config.yaml file is read by the C++ layer via an input class and stored within a structure. The coupling constants for each particle respective to every mean-field (Tables 5 and 6) are computed. The different particle classes (quarks, baryons from octet, and/or decuplet) are initialized and filled with their quantum numbers read from the PDG table 2021+ [119] and the couplings just computed.
The code loops over desired , , , so the chemical potential per particle is computed via Eq. 50, then loops over every mean-field initial guesses (, , , , , , ) follow. The fsolve routine is called where these initial values, in conjunction with the input parameters provided by the user, are used to compute the right hand side (RHS) of Eq. 46. To compute the left hand side (LHS) of Eq. 46, the scalar Eq. 56 and vector Eq. 55 densities must be obtained for each particle involved, which implies the calculation of the effective chemical potential via Eq. 49, the effective masses (see Eq. 43 for hadrons and see Eq. 44 for quarks), and the Fermi momentum.
The fsolve routine then computes the gradient for every field equation involved and updates the mean fields and field to an improved guess that minimizes the difference in RHS and LHS of Eq. 46. If the new solution for the fields lies outside of the domains, then the code skips to the next initial conditions guess. If the solution found is not self-consistent (LHS not equal to RHS), recompute the effective masses, chemical potentials, and scalar and vector densities using the new guesses and evolve the field solutions along the gradient.
The previous procedure is performed until self-consistency is achieved, which means that LHS is equal to RHS within a certain threshold and that the solution has not been achieved before. Given that a valid solution has been found, the code now proceeds to compute a collection of thermodynamic observables like pressure, energy density, density (see Eq. 64), and other relevant quantities like strangeness density, charge density, density without , and densities per particle sector (quarks, baryon octet, baryon decuplet). This data is written to an intermediate file, and the C++ layer continues its execution into the next field’s initial condition guess.
Once all the domains of interest have been exhausted, the main algorithm execution finishes.
III.4 Stability and phase transition criteria
Let us begin the discussion by defining susceptibilities of the pressure:
| (72) |
where whatever chemical potential is not being varied is kept constant as well. Due to the symmetries in QCD, the ordering of the derivatives does not matter, i.e.
| (73) |
where and are any combinations. The first susceptibilities relate to the respective density of each conserved charge i.e.
| (74) |
Additionally, the second-order susceptibilities are then equivalent to
| (75) | ||||
| (76) | ||||
| (77) |
which have been shown to have interesting connections to the speed of sound in neutron stars [120]. The susceptibilities are also important to provide connections to the search for the QCD critical point at finite and understanding the deconfinement phase transition [121, 122, 123, 124, 125, 126, 127, 128, 129, 130].
A first-order phase transition is defined as a jump in at a specific . Higher-order phase transitions appear as jumps in the higher-order susceptibilities. Thus, an -order phase transition occurs at the point if diverges. When an -order phase transition occurs then all higher order derivatives also diverge i.e. diverges where at . However, we only determine the order of the phase transition by the first derivative where either a jump or divergence occurs.
In the grand canonical ensemble, in the infinite volume limit, stability corresponds to minimizing the grand potential density or maximizing the pressure (see Appendix E). For this case, and assuming BSQ conserved charges, the 4-dimensional Hessian matrix is shown in Section E.5. In the following, we show results only for . In this case, the Hessian matrix is 3D:
| (78) |
where the matrix is symmetric due to Eq. (73. Then, the determinant of each submatrix must be zero or positive. Thus, for the matrix
| (79) |
and for the matrix
| (80) | ||||
| (81) |
Using Eq. (81), then it also implies that because is real. Finally, the matrix gives the condition that we show later on in Eq. (86).
The matrix defined in Eq. (78) was somewhat arbitrarily built in that one could also have ordered it as SQB or QBS (or any other ordering). Thus, when considering all perturbations of the matrix we then arrive at the following independent conditions:
| (82) | ||||
| (83) | ||||
| (84) | ||||
| (85) |
| (86) |
The susceptibilities can be thought of as moments of the net-BSQ distributions (again recalling that the first moment implies the respective BSQ charge densities). Then, Eq. (82) implies that the variance of each net-BSQ distribution is positive and Eqs. (83-85) imply that the covariances must also be semi-negative definite. Finally, we note that the matrix in Eq. (78) becomes more complicated at finite temperatures. However, we leave finite studies to future work.
For multiple solutions of the EoS, if more than one solution obeys the stability criteria, then the one with the highest pressure (or lowest grand potential density) at fixed a value is denoted as the stable EoS. The other EoS’ that obey Eqs. (82-85) are called metastable 333If we have two mixtures, I and II, I is stable and II is metastable if and both satisfy Eqs. (82-85).. Additionally, if the pressure of an EoS is negative, but it obeys the stability criteria, it is also considered metastable. If there is a unique EoS with that obeys Eqs. (82-85), then vacuum solutions are considered stable. We summarize our stability criteria in Table 13.
| stability label | thermodynamic criteria | multiple solution criteria for phase |
|---|---|---|
| stable | Eqs. (82-86) | single solution |
| metastable | Eqs. (82-86) | |
| unstable | At least 1 fails: Eqs. (82-86) | |
| stable vacuum |
The variables related to the stability of the system also dictate the occurrence of a phase transition. For example, a first-order phase transition, such as the quark deconfinement transition at low temperatures, occurs when
| (87) |
where the superindex indicates the hadronic phase and the superindex indicates the quark phase. The presence of strangeness and charge chemical potentials offers different possibilities for the phase transition, such as making the phase transition at fixed charge fraction or strangeness fraction [10, 78]. In this paper, we assume all charges are conserved during the phase transition (a non-congruent transition), where
| (88) |
In Figure 4 different phases relevant, e.g., for the description of the core of neutron stars, are shown: vacuum, hadronic matter, and quark matter. The top panel shows a first-order phase transition in vs space and the bottom shows the first-order phase transition in vs space. Because vs must be continuous, we see a clear maximum solution at each point in . In contrast, the stable solution for the first-order phase transition demonstrated here has a jump in such that across a range of we see only metastable and unstable solutions.
As the chiral model in its current version does not include nuclei, we reproduce the liquid-gas phase transition as being from vacuum to bulk baryonic matter. Depending on how they are connected and which ones are present, these phases can be unstable, metastable, or stable (see Table 13). If a system is in equilibrium, then a Maxwell construction can be performed across the metastable/unstable regime, such that the EoS remains stable even across the phase transition. However, dynamical simulations often require metastable/unstable regimes in order to accurately describe the time spent in each phase of matter (see e.g. [15, 16]). Thus, in CMF++ we build Maxwell constructions, but also preserve the metastable/unstable regimes.
Given an EoS with a metastable/unstable regime across a phase transition, one can obtain the Maxwell construction in one of two ways:
-
•
using the equal area method in space, in which one finds the line (dashed line) such that the two areas in Figure 4 d) are equal i.e. . See [131] for examples and discussion in a van der Waals model for the liquid gas phase transition. This method is more typical for models within the canonical ensemble;
-
•
choosing the maximum pressure (minimizing the grand potential density) at a specific point in (one can also do the same at a specific point in ), which is demonstrated in Figure 4 (a-b). This method is more typical for models within the grand canonical ensemble.
In this work, we follow the procedure depicted in Table 13 and apply the second method to find the (most) stable phase such that at each point in our phase diagram, we choose the given multiple solutions where . The second method is significantly easier in CMF++ because the metastable/unstable regime in CMF++ can become significantly more complicated than the sinusoidal appearance shown in Figure 4 d). Rather, depending on the degrees of freedom one may find more than 2 solutions or even solutions that cross each other. Thus, it is not always obvious what the definition of the areas is with so many solutions present, such that the equal area method would be impractical. To differentiate between unstable and metastable phases, we also follow the procedure depicted in Table 13.
III.5 Benchmark
In Figure 5 a benchmark of the time it takes to run CMF++ compared to the legacy CMF in Fortran is shown. Along the x-axis, we demonstrate the typical total number of EoM systems solved (the set of equations in Eq. 46) for 1D, 2D, and 3D EoSs. For the 1D case, and were kept at zero and was varied from 10 points (10k EoM systems) to 100 points (100k) and finally 1000 points (1M). For the 2D case, was kept at zero, 1000 points were used in and 10 points in . Finally, for the 3D case, 1000 points were used in , 10 points in , and 10 points in . Along the y-axis, we find the average runtime for 16 different particle configuration combinations (4 vector potentials, decuplet on/off, quarks on/off). The dashed line represents an extrapolation given Fortran’s extreme runtime, where the star at the end is an educated guess based on the behavior at 20M. It is important to note that the time complexity for CMF++ is whereas Fortran has an one. We find that CMF++ significantly improves the performance of the calculations of the EoS by at least an order of magnitude (for the simplest calculations) and up to 4 orders of magnitude for the 3-Dimensional case.
IV Results
In this section, we present our numerical findings, exploring all four vector couplings across various combinations of degrees of freedom (d.o.f.) within the CMF model, using different combinations of independent chemical potentials.
Our general approach in the following sections is to demonstrate our results for the mean fields, , and certain thermodynamic variables. Then, we show population plots for the individual species of hadrons and/or quarks. Finally, the charge fractions and susceptibilities are shown. Initially, we demonstrate that the new CMF++ can both reproduce the legacy Fortran version of CMF and also obtain more precise results in 1D. Later, new results across the 3D phase space of are only shown for CMF++ due to the extremely long run times that they would take in the legacy Fortran code.
IV.1
We begin by examining various sets of d.o.f, considering the simplest case where .
IV.1.1 C3 and C4 with baryon octet + quarks
We start by exploring the behavior of mean-field mesons, the deconfinement phase transition order parameter, , and thermodynamical properties, including the baryon octet plus quarks as d.o.f under the influence of C3 and C4 vector couplings (see Tables 8, 9 and 10 for related parameters). Displaying all coupling schemes would involve an excessive amount of quantitative detail; therefore, we only present the results for C3 and C4 couplings, as C3 behaves similarly to C2 and C1 (see Sec. F for details on the other couplings).
We begin with the C3 coupling in Figure 6 in the limit of . The mean-field mesons (, , , and ) as a function of are shown in panels a)-d), respectively. Additionally, Figure 6 contains in e), the pressure vs the energy density in f), the speed of sound squared vs the baryon density in g), and the baryon density vs the baryon chemical potential in h). Also shown in Figure 6 are the different types of solutions obtained in CMF++ that we classify as stable, metastable, and unstable solutions. The legacy Fortran code only provides stable solutions, which are also shown for comparison. The criteria for stability are comprehensively discussed in Sec. III.4, as well as Figure 4.
We observe a consistent correspondence between the results obtained from C++ and Fortran across all mean fields for stable solutions in Figure 6 panels a)-d). Additionally, and the thermodynamic quantities are all precisely reproduced in the CMF++ version of the code, as one can see in the comparison of the black solid vs red long-dashed lines in panels e)-h). In all panels, we see identical results from stable solutions, except for the region with very low , where Fortran has trouble finding solutions.
Within a given phase of matter (either hadron or quark), the stable solutions for the non-strange scalar (panel a)) and strange scalar (panel b)) ratios, hereon simply referred to as and , monotonically decrease with increasing . The decreasing trend in scalar condensates is indicative of chiral symmetry restoration. Panel b) resembles panel a) but we can see that the field has a kink just below MeV, coinciding with the emergence of hyperons (refer to the discussion of hyperon population Figure 7). Also in the quark phase, we see a kink in at MeV marking the appearance of strange quarks. The non-strange vector field (panel c)) exhibits the opposite behavior than (panel a)), while the strange vector field (panel d)) remains zero until the emergence of strange particles (in this case hyperon) and then presents a very similar behavior to panel b). Post deconfinement phase-transition, both vector fields become zero, as they do not couple to the quarks. Since we are discussing isospin symmetric matter (), there is no finite value for isovector mesons.
As discussed previously in Sec. III.4, CMF reproduces three distinct phases of matter at : vacuum, hadronic phase, and quark phase. In Figure 6 the Maxwell constructions of the first-order phase transitions generate vertical lines at a fixed for all the mean fields in panels a)-d). The liquid-gas phase transition (vacuum to hadrons) occurs at low and the deconfinement phase transition (hadronic to quark) occurs at intermediate to high . Within the phase transition itself, unstable phases may appear.
Additionally, within both the hadronic and quark phases, there are differences between phases that only/mainly contain light hadrons or light quarks vs those that contain strange hadrons or strange quarks. Separating the light vs strange dominant regimes within a given hadronic or quark phase may also be a phase transition of various order. In fact, there can appear sometimes even first-order phase transitions leading into strangeness dominant phases (light hadrons to strange hadron dominated) that occur before the deconfinement phase transition for specific parameter sets of CMF. In the following, we discuss the appearance of all possible phases in Figure 6 and in subsequent CMF parametrizations across different combinations of .
At low , CMF reproduces a first-order liquid-gas phase transition. In Figure 6 the liquid-gas phase transition occurs at MeV. One can see the telltale vertical line in , , and in Figure 6. However, the meson does not experience the liquid-gas phase transition because strangeness is not relevant at such a low . Similarly, , the order parameter for the deconfinement phase transition remains zero throughout the liquid-gas phase transition.
The phase transition at higher is the one related to quark deconfinement. There is a strong relation between and the meson mean fields. The value of is shown in panel e) of Figure 6. The change in its value from to at MeV signals the change in values for the effective masses of baryons (which become too large for them to be present) and quarks (which become light enough to appear).
Looking at the meson, we find that strangeness begins to play a role at MeV when the meson begins to deviate from 0. At the same time, we see that both the field deviates further from the vacuum at that point as well. However, there is not a first-order phase transition as the strangeness begins to play a role because there is no clear vertical line in any of the mean fields between stable phases when the strange mean fields begin to become non-zero. We later analyze the susceptibilities in order to determine the order of the strangeness phase transition.
Even within the first-order phase transition, most of the solutions fulfill the stability criteria (although they do not have the maximum pressure) such that they are labeled metastable (shown in green). However, one can see small regions of unstable phases that appear (shown in pink). One surprising outcome can be seen fairly clearly in the plot in panel b) and in the plot in panel d). There is a small unstable region around MeV that connects two metastable regions, which is indicative of a first-order phase transition that connects two metastable phases. The two metastable phases contain strange hadrons (indicated because the and mean fields mediate the strange interactions). Thus, the phase transition (between metastable phases) goes from a hadronic phase with some strangeness into a strangeness-dominated hadronic phase. However, this strangeness-dominated hadronic phase has a lower overall pressure at a fixed than the quark phase, such that it is not considered a stable solution.
The EoS is a relationship between the pressure and energy density , which is shown in panel f) of Figure 6. The stable branches with the Maxwell construction show a monotonically increasing EoS with a slight kink around MeV/fm3 where the hyperons switch on. Then, one can see the first-order phase transition at the plateau where remains constant while increases. In the quark phase, the pressure monotonically increases once again, with a very small kink appearing again when strange quarks become relevant (although that is quite hard to see).
The metastable/unstable solutions are shown in green and pink, respectively, which include both quark and hadronic solutions. As discussed previously for the strange mean fields, we can see in panel f) a phase transition within metastable phases (going from a light-dominated hadronic phase into a strangeness-dominated hadronic phase) that appears at higher pressure for a fixed energy density). The lower branches that appear at low pressure below the first-order phase transition correspond to light quark phases. Note that if quarks had not been considered, the hadronic phases would be the stable ones and a Maxwell construction would have to be built across the unstable hadronic phase.
Some of the features of the EoS may be difficult to pinpoint in the plot (specifically where strange hadrons or quarks become relevant). However, in panel g) showing vs these features become clearer (being the derivative ). For instance, looking at and starting first at low , we can see several interesting features. First, strange hadrons appear, creating a small peak/kink in (denoting a higher-order phase transition) when the ’s appear. Then the next feature is a drop in that leads to a plateau. The plateau corresponds to the first-order phase transition (when you have used a Maxwell construction). The metastable and unstable regimes are also shown across the phase transition. They provide a different structure instead of the plateau from the Maxwell construction. Finally, after the phase transition, we see that the rises again because of the quark phase. We can see that there is a small peak in the quark phase that arises from strange quarks appearing. Afterwards, remains near the QCD conformal limit of and is consistent with the pQCD results [132] as well as their constraints following stability and causality [133].
Finally, the vs plot is shown in panel h). Note that the number density contains an additional non-fermionic contribution from (Eqs. 63 and 64). At the liquid-gas phase transition, there is a small vertical jump at low , followed by a monotonically increasing until the deconfinement phase transition is reached just after MeV. Then, in the quark phase, there is a steeper rise in the compared to what occurred in the hadronic phase (due to the contribution). The metastable/unstable branches are shown in green and pink, respectively. The branch that is continuous with the hadronic phase at lower is for the hadronic into strangeness-dominated hadronic phases, whereas the upper branch that connects to the quark phase is the metastable light quark phase. We can see that for this specific parametrization, the deconfinement phase transition happens at relatively large (for ).
To analyze the individual contributions of particles to the EoS, we plot the population of particles against in Figure 7. The populations are defined at the baryonic number density of specific species such that if they were all added together we would recover the baryon number density (minus the contribution in the quark phase) i.e.
| (89) | ||||
| (90) |
where Eq. (89) provides the baryon density in the hadronic phase and Eq. (90) provides the baryon density in the quark phase (given the species that appear for this specific C3 coupling at this specific choice in the space). Given that we have , we are dealing with the symmetric nuclear matter or in other words have isospin symmetry. Due to isospin symmetry, the in-medium mass and density of protons and neutrons exhibit degeneracy i.e. and at all .
The nucleons begin to populate around
| (91) |
where is the average nucleon vacuum mass and is the binding energy per nucleon. Additionally, at MeV, hyperons emerge, opening the Fermi sea and thereby softening the EoS (this coincides in the first peak in in Figure 6). Following the deconfinement phase transition around MeV, baryons are replaced by quarks, resulting in a growth in the number density of quarks with . The density of and quarks are equal across all . Around MeV, the quarks switch on, which corresponds to the second peak in during the quark phase in Figure 6. One significant difference in the hadronic vs quark phase is that baryon shares a significantly larger fraction of the total baryon number density (with respect to protons and neutrons) than the strange quark does (with respect to up and down quarks).
In Figure 8 we have our last set of figures for C3 and . In panel a) we show the strangeness fraction vs the baryon chemical potential and in panel b) we show the electric charge fraction vs the baryon chemical potential. The baryon susceptibilities vs the baryon chemical potential are shown for second-order in panel c) and third-order in panel d).
Naively, one may expect that for , we have and . However, as we can see in panel a), it is possible to obtain due to the switching on of the baryon in the hadronic phase and later strange quarks in the quark phase (recall that a strange quark carries strangeness, such that ). After the emergence of the hyperon (refer to Figure 7), the magnitude of steadily increases until MeV, when it jumps to zero at the deconfinement, at which point the strange quarks slowly appear increasing the strangeness magnitude (see again Figure 7). In the limit of isospin symmetric matter, the baryon contributes no net-isospin since its quark content is and, therefore, is more preferred compared to the (uus) or (dds), which would create an isospin imbalance. Note the increase in the magnitude of in the metastable and stable phases at large in panel a). This increase in the magnitude of in the metastable regime arises due to the hadronic phase with many hyperons (i.e. a strangeness-dominated hadronic phase) that would have appeared were quarks not included in the calculation.
In the charge fraction plot in panel b) we find that remains at 0.5 (representing an equal amount of positive and neutral nucleons) until the appearance of the neutral hyperon. After the hyperons appear, continues to decrease because they create a further imbalance between charged and neutron hadrons. At densities above the deconfinement transition, begins to decrease once more with increasing due to the quark’s increasing relevance. Note that the isospin fraction remains zero for all . The relation between the different fractions is related to quantum numbers and can be derived from the Gell-Mann–Nishijima formula, resulting in (for negative strangeness). This is discussed in detail in ref. [78]. We note that the metastable and unstable branches are also shown in and mirror the same qualitative behavior as was already discussed for .
To gain a clearer insight into the different phase transitions, second and third-order susceptibilities ( with ) are displayed against in panels c) and d), respectively. The susceptibilities are normalized by different orders of to ensure they are dimensionless. Recall that the first-order susceptibility is just , which was already shown in Figure 6 where we already saw the first-order discontinuities for the liquid-gas phase transition ( MeV) and deconfinement ( MeV). Then, the two first-order phase transitions are propagated into as two divergences and into also as two divergences (although they are significantly larger).
In we can see that two kinks appear at MeV and MeV that correspond to the transition where strange baryons and strange quarks appear, respectively. Then, in the kinks turn into divergences at these same locations, which indicates that the onset of strangeness is a third-order phase transition. Thus, we can draw an interesting connection here to that displays kinks at precisely these locations as well. Given that one can show that for a single conserved charge that [120] and is smooth and continuous, then this kink in also appears in as an inverted kink.
Figure 9 depicts the same set of plots that include mean fields, , and thermodynamic properties as in Figure 6, but for the C4 coupling scheme. We once again include the information on both the stable, metastable, and unstable phases. We also compared the legacy CMF to CMF++. We found a strong agreement between the results obtained from C++ and Fortran solutions.
In Figure 10 the corresponding population plot is shown for the C4 coupling for the combination of octet+quarks. Finally, in Figure 11 the charge and strange fractions are shown as well as the susceptibilities of the pressure. In the following, we discuss these three plots and compare and contrast them with the previous C3 coupling that we showed as well.
The key difference between the C3 and C4 couplings lies in the presence of mixed couplings and in the C4 coupling (Eq. 27), leading to a different sensitivity to strangeness. Due to these changes, there is no stable hadronic phase that includes strange baryons for the C4 couplings at (see Figure 10). However, strange quarks do appear in the deconfined quark phase.
The mean fields associated with strangeness ( and ) include self-interactions for the meson and the meson that couples with the non-strange scalar field (see Eq. 46), resulting in additional attractive interactions. Thus, due to these self-interactions and couplings, even when no strange baryons are present in the hadronic phase, we find that and (see Figure 9). More specifically due to the and coupling terms in C4, has a much larger absolute value in the hadronic phase, see panel d) in Figure 9, compared to C3 in panel d) in Figure 6. On the other hand, comparing in panel b) to C3 in panel b) in Figure 6, we find that is closer to unity for C4.
We only focus on the deconfinement phase transition for C4 because the liquid-gas phase transition has the same properties as C3. For the C4 coupling, across panels a) to b) of Figure 9, we observe a jump in the ratio of the scalar mean fields, after which the scalar mean fields continue to decrease. The vector fields become zero at this point as they do not couple to quarks. Based on the behavior of the mean fields, we conclude that the phase transition for the C4 coupling shifts to a somewhat lower MeV. The shift in the phase transition to lower is confirmed by the behavior shown in panel e) where we see the vertical line at MeV and the shift is also confirmed later on in the population plot in Figure 10. The first-order phase transition also appears in panels f) for and g) for of Figure 9. The signatures of the first-order phase transition are similar to C3, even though they appear at a different location. C4 has a slightly larger phase transition than C3 (larger jump in ) because it is a sharper change from protons and neutrons into quarks than if other hadrons had appeared before deconfinement.
In order to understand strangeness in the coupling C4, let us begin with the population plot in Figure 10. We see that for the hadronic regime, we only have protons and neutrons (in precisely equal amounts since this is for symmetric nuclear matter) such that there are no strange baryons. Thus, in the thermodynamic properties of C4 the lack of strange baryons implies that there is no kink in in the stable hadron phase (panel g) of Figure 9).
However, in the quark regime, we see in the population plot in Figure 10 that strange quarks are present and they play a more significant role compared to the C3 coupling. Because of the presence of these strange quarks, we see there is a kink in (panel b) of Figure 9) when the strange quark switches on and then becomes significantly larger in magnitude in the quark phase as the amount of strangeness increases with . In terms of thermodynamics, we see a kink in at high when the strange quarks appear (panel g) of Figure 9). Additionally, we find that is significantly steeper in the quark phase (panel h) of Figure 9). The influence on and are both consequences of the fact that the strange quarks play a much larger role in the C4 coupling than in the C3 coupling, even though they appear at roughly the same .
As discussed above, the C4 coupling does not produce strange hadrons in the hadron phase. Since strange hadrons soften the EoS, then the C4 coupling has a stiffer EoS at low compared to the C3 coupling. In fact, the larger value in C4 leads to a stronger repulsive force, that gives us a stiffer at low . The stiffer EoS at low then results in the modeling of more massive neutron stars when finite isospin is included [51, 72, 134, 53, 54, 81, 68, 110, 67, 69, 135, 10, 77, 78, 136].
In the quark phase, the looks very similar between C3 and C4. However, because of the stiffer low EoS for C4 (and the fact that is a derivative of ), then that leads to also a stiffer high EoS for C4 as well (even though both of their look very similar at high ).
One of the most significant differences between the C4 and C3 (see panels g) of Figure 9 and Figure 6) couplings is that the metastable/unstable regimes are significantly different for C3 compared to C4. We find that the order parameter has the same qualitative shape, although the phase transition does take place at a different spot (see panels g) of Figure 9 and Figure 6). However, when we look at we see that the metastable region is significantly more simplistic for C4, such that there is no longer a phase transition within hadronic metastable phases that existed previously for C3 (see panels f) of Figure 9 and Figure 6). Thus, C4 does not present a strangeness-dominated metastable hadronic phase (that would be stable if it wasn’t for the quarks) as C3 did. However, one can see that for the metastable regime in around , there is a small kink, which implies that the strange hadrons would become non-zero if quarks were not considered (see panel g) of Figure 9).
Due to the absence of hyperons in the baryonic sector (for C4 in the case with ), both and in panels a) and b) of Figure 11 remain and until the strange quarks appear (compare with Figure 8). However, once the strange quarks appear after the phase transition, we see that in the stable phase is significantly larger in overall magnitude for C4 than for C3, demonstrating once again the larger role that strange quarks play in C4. That effect also leads to a much larger deviation of from 0.5 in the quark phase. The metastable hadronic phase shown for and also makes it clear what we discussed previously when we looked at the metastable region of that for the C4 coupling strange baryons appear in the metastable regime.
Finally, we return to the susceptibilities to better understand the order of the phase transitions (panels c) and d) of Figure 11). To remind the reader, we already saw jumps in the first-order baryon susceptibility, i.e in Figure 9 at MeV for the liquid-gas phase transition and at MeV for the first-order deconfinement phase transition. Then in the second-order susceptibility, we see significantly larger divergences at these points. As expected, for , those divergences are amplified further (higher-order susceptibilities scale with high-orders of the correlation length such that they are more sensitive to phase transitions, see e.g. [121]). We can also see the telltale kink in when the strange quarks switch on at high . The kink in then leads to a divergence in such that we also see a third-order phase transition for the C4 coupling within the quark phase. This third-order phase transition also has the kink in , providing another example that a third-order phase transition leads to a kink in .
IV.1.2 C1-C4 with baryon octet + decuplet
Up until now, we have considered a combination of the baryon octet and quark phases. Here we explore the alternative where no quark phase is present but a larger, more complex hadronic phase is possible wherein both the baryon octet and decuplet are possible. Thus, in this scenario, as increases we anticipate that a wealth of new baryonic states are switched on. Since the decuplet includes ’s, the new baryons do not necessarily carry strangeness, but maybe light states as well. Note that the appearance of first ’s vs hyperons depends on the couplings and parameters in the CMF module.
At this time, the quark couplings have not yet been fitted while also including the decuplet baryons. Thus, we have not yet incorporated the combination of the octet+decuplet+quarks in this analysis. We have performed some initial testing with the possibility of octet+decuplet+quarks and found that in that case (within the current parametrization) the quarks appear at very large ’s, outside the regime of neutron stars.
In Figure 12 we show the mean fields, , and thermodynamic properties for C1-C4 for the case of the octet+decuplet at . While we no longer show the comparison between C++ and legacy Fortran versions of CMF, we note that we have checked their consistency for all following plots. However, we do not show the comparisons to improve the readability of the plots. Finally, we note that the C4 coupling reproduces for the chosen d.o.f. a much lower density and energy density than the other couplings (for a given , so it only appears at the beginning of the plots that use those quantities as the x-axis. In the presence of both the baryon octet and decuplet, the differences between C1-C4 are overall more significant than before. Previous work [53] has shown that C1-C3 produce nearly identical results in the case of just the baryon octet+quarks (although there are small subtle differences in and close to the phase transition from vs ). Now with the decuplet present, we find that all lead to distinct solutions that can be seen quite clearly in both the mean-field mesons and the thermodynamics.
In panels a)-d) of Figure 12 we find general trends in the mean-field mesons that are the same as what we saw for the baryon octet+quark configuration i.e. , , and always decrease with while increases. However, their exact behavior within these general trends can differ quite significantly and present jumps/kinks as different particles switch on. The C1 and C4 couplings are simpler because they do not contain first-order phase transitions from a light hadronic phase into a strangeness or -dominated hadronic phase. Thus, we find that the C1 and C4 coupling schemes change more smoothly as one increases . Specifically, the C4 coupling scheme even with the baryon octet+decuplet, has significantly fewer contributions from hadrons that are not protons and neutrons (see Figure 13). Because C4 is dominated by protons and neutrons, it leads to a significantly stiff low EoS, and the means fields have a smooth behavior across all . In contrast, has a relatively smooth behavior in and , but there is a sudden change (not quite first order) in and , where we can see sharp drops in their values around MeV (again, in Figure 12). The increase in the magnitude of and still occurs across a range of , such that we would not classify it as a first-order phase transition, but it appears to be a sharper higher-order phase transition (we discuss the exact order in the susceptibilities later on).
We now come to the C2 and C3 couplings, which demonstrate quite interesting behavior. For every panel a)-d) in Figure 12, we see a clear vertical line around MeV for C2 and MeV for C3. For the mesons, these vertical lines translate to a small jump in and , and a very significant jump in and . Thus, what we find is that we have a first-order phase transition from a light baryon-dominated regime to a strange and baryon-dominated regime. We later explore the implications of what baryons switch on across this phase transition when we explore the population plots. This first-order phase transition in the hadronic phase is the same phase transition that we saw previously for C3 couplings when quarks were present in Figures 6 and 9 that fell in the metastable regime.
Panel e) in Figure 12 confirms the absence of quarks, as indicated by the field remaining zero across the entire range. Thus, the first-order hadronic phase transition that we see in Figure 12 arises entirely from the couplings, effective masses, etc. and is not driven by the explicit order parameter that we have built into our model.
Panel f) in Figure 12 presents the EoS for all coupling cases. We see significant differences from C4 compared to the other couplings. The C4 couplings lead to a very stiff EoS that steadily rises with . If we then compare this to the speed of sound in panel g), we see that we also reach the largest with C4, such that . One can see a small kink in for C4, which is our first hint that other particles beyond protons and neutrons switch on for C4, but it is a rather subtle effect. It is clear then, looking at in panel h), that the lack of other hadrons that switch on for C4 leads to the lowest corresponding for a given . In other words, the fewer hadronic species possible at a given implies a lower .
Let us now return to the EoS in panel f) in Figure 12 and explore the C1 coupling that appears to fall in between what is seen for C4 and C1-C2. Recall that the C1 coupling did not have a first-order phase transition into a strangeness-dominated hadronic phase. However, it did demonstrate steep changes in its mean-field mesons, potentially indicating a higher-order phase transition. In we can see signs of this higher-order phase transition by the slight flattening in the pressure, but never reaches a true plateau like one would anticipate for a first-order phase transition. If we then look at the speed of sound in panel g), we can see the quite interesting non-monotonic behavior in C1. The kinks in correspond to the appearance of new baryons at a given . Altogether, we find 5 distinct peaks in vs for C1, which likely indicates the presence of 5 new hadronic species switching on at certain ’s. Additionally, we see a softening where drops to but does not hit exactly , as one would expect for a first-order phase transition. Finally, this behavior leads to a steeper rise in in panel h) because of the new hadronic degrees of freedom that switch on.
Finally, we come to C2-C3 which has shown indications that a first-order phase transition occurred within the mean-field mesons, even though the deconfinement order parameter is always . At low for C2 and C3 we see a tiny peak in in panel g) of Figure 12 around that occurs at the same location as what is seen for C1. Thus, we likely have at least one new particle switching on before the first-order phase transition occurs at higher . Looking at the relationship in panel f) we see that there is a plateau consistent with a first-order phase transition that then translates into a region of in panel g) and a jump in in panel h) for C2-C3. Following the first-order phase transition, we find a very steep increase in in panel h), which is consistent with what we understood before - new degrees of freedom leads to a larger at a fixed .
The evidence is quite clear that we have a first-order phase transition, but that the phases of matter are always hadronic on both sides of the phase transition. It is also interesting to compare this hadronic phase transition to what we saw previously for the baryon octet into quarks in Figure 6 and Figure 9. The general behavior of the C2-C3 EoS is very similar for the deconfinement phase transition vs the hadronic phase transition. In fact, the hadronic phase transitions for C2-C3 take place at nearly the same location as the deconfinement phase transition, which is why in Figure 6 we saw a phase transition within the metastable regime.
Here we show the population plots for C1-C4 in Figure 13. The C4 coupling has the simplest population, so let us begin with that. We find that for the C4 coupling the system is heavily dominated by just protons and neutrons, although baryons appear at high . However, their contribution is only a very small fraction of the baryon number.
Let us now discuss the C1 coupling, for which previously we suspected 5 new particles to switch on due to the peaks in . We also previously found that there was no first-order phase transition within the hadronic phase. When we look at the population plot, we find confirmation of both of these facts. The particles slowly turn on (there is no distinct jump in their population numbers) while the influence of the proton/neutron begins to decrease around MeV. Additionally, we can confirm that new particles (all degeneracies at once) switch on in the (increasing ) order of: first ’s, second ’s, third and , fourth and , fifth and . Thus, it is true that 5 different species of particles turn on at those peaks, but each peak sometimes includes multiple versions of that particle with different electric charges. We find that once the baryons switch on we reach the peak in protons/neutrons such that they decrease as the contribution of ’s increases. We also observe that as the baryons switch on, we see the influence of ’s wane as well.
In the C2 and C3 coupling, we observe that in the low region the hadron phase consists of predominately protons, neutrons, and ’s that appear a bit before MeV. For the C3 coupling, we even see a small contribution of ’s just before the first-order phase transition appears into the strangeness and -dominated phase. We then see a first-order phase transition at MeV for C2 and at MeV for C3 to a different mixture of hyperons and baryon decuplet (roughly going from top to bottom in density): first , nucleons, , , , and , then , then . C2 and C3 have slightly different mixtures of all these baryons after the first-order phase transition, however, the key result is that are the most populous baryon in this strangeness-dominated phase and even baryons are allowed.
Panels a) and b) of Figure 14 depicts the evolution of and with respect to . Because the C4 coupling has the least populous strange baryons, and remain at 0 and 0.5, respectively, until the emergence of the hyperon at large . Conversely, for C1, and remain constant until the first strange particle switches on () where smoothly because non-zero and smoothly decreases. As more new hadrons switch on, and both begin to drop more rapidly. The largest drop occurs around [MeV] where the becomes relevant and rapidly increases in importance. Finally, for C2-C3 we see a vertical line in and where the first-order hadronic phase transition occurs. The overall magnitude of is the largest for C2-C3 because even hadronic states are switched on.
The susceptibilities in panels c) and d) illustrate the location of the first-order phase transition in C2 and C3 at MeV and MeV, respectively (in addition to the liquid-gas phase transition). For C1 we see an especially surprising effect. In we have two second-order phase transitions (or maybe nearly second-order) that appear in the range of [MeV] that relates to ’s and ’s rapidly switching on. These (nearly) second-order phase transitions lead to the softening in that was seen earlier in panel g) of Figure 12. We specify that these are “nearly” second-order phase transitions because we do not see and it does not appear that is quite diverging, although it looks very close to that behavior. in panel d) confirms that the other phase transitions that we saw in C1 are all of third-order. We do not find any phase transitions in C4 (beyond the liquid-gas phase transition).
IV.2
In the following discussion, we analyze the role of the strangeness chemical potential while keeping the charge potential fixed, . This scenario is similar to what would be seen in low-energy heavy-ion collisions (although they are typically at finite ) because there are local fluctuations of both baryon charge and strangeness such that fluctuations in also appear, but the system is nearly isospin symmetric such that (the amount of isospin asymmetry depends on the choice of initial colliding nuclei) is close to zero.
Since the strangeness number is negative (see Table 2), a positive at large leads to a smaller magnitude of , and a negative leads to a larger magnitude of strangeness density. Nevertheless, a more positive could have no effect on the system in the case that there was not enough energy for strange particles to appear in the first place. In principle, an extraordinarily large would lead to the preference for anti-strange particles, however, this is not a limit that is relevant to neither heavy-ion collisions nor neutron stars, so we do not explore it here. In astrophysics, . In heavy-ion collisions, the system is globally strangeness neutral and if one is calculating quantities related to the entire system, can easily reach up to a third of to ensure strangeness neutrality. That being said, local fluctuations of strangeness absolutely exist and have been measured experimentally (in fact, strange particles account for roughly of the final particles produced at high energies [12], it is just that they are exactly balanced by particles that carry anti-strangeness). Thus, it is possible that one subsection of the fluid experiences whereas another section experiences (see e.g. [19]).
IV.2.1 C3 and C4 with baryon octet + decuplet + quarks
We now focus on the C3 and C4 coupling schemes, which are shown using 3D density plots in Figure 15 and Figure 16, respectively, and allow the presence of the baryon octet, decuplet, and quarks. In these figures, only the stable solutions of C++ results are shown. Let us now explain how we can interpret the 3D plots. The x-axis corresponds to the baryon chemical potential and the y-axis corresponds to the strange chemical potential. Then, each plot uses the color scheme to depict a given variable’s value, with its range indicated on the right side of the image. Panel a) shows the ratio of the mean-field to its vacuum value, b) shows the deconfinement order parameter , panels c) and d) show the strangeness and charge fractions, panel e) depicts the ratio of baryon density to the saturation density and panel f) the speed of sound squared in units of the speed of light.
For both C3-C4 in the light hadronic phase (not strangeness or dominated), we find that the field is always maximum at low and then continuously decreases with (regardless of the value in , as discussed in Sec. IV.1. At some point, a discontinuity appears that we will discuss later, but then a new phase of matter appears at large . The nature of this new phase of matter at large depends on both the coupling and the values of .
For the C4 coupling, the new phase of matter is a quark phase wherein within the quark phase the (panel a) of Figure 16) has a maximum immediately after the phase transition, followed by a monotonic decrease as increases. Referring to the order parameter in panel b), we can see clearly that this phase transition corresponds to the deconfinement phase transition, since in the quark phase. We find that the deconfinement phase transition for C4 does not depend on the value, which is because ’s appear significantly later in the C4 coupling such that they do not significantly influence the mean fields/thermodynamics in the hadronic phase.
However, for the C3 coupling, the intermediate/high phase of matter strongly depends on . More positive leads to a quark phase for C3 coupling, however, more negative leads to a strangeness-dominated hadronic phase. We can distinguish between the quark phase by using the order parameter in panel b) of Figure 15 that only approaches 1 for most positive values of (although some small negative values of at large do also lead to a quark phase). Then, it is clear that the other phase of matter must be a strangeness-dominated phase by referring to in panel c) in Figure 15 where in this phase. In the quark phase for C3, is significantly larger than for the strangeness-dominated phase wherein approaches zero, which is indicative of the chiral restoration at very high densities. This highlights the anti-correlation of and at the deconfinement phase transition (see Figure A21 in Sec. F) related through hadronic and quark effective masses (Eq. 43 and Eq. 44).
For C3, one very interesting consequence of these two different discontinuities is that it gives rise to a triple (critical) point. In fact, this is a quantum triple critical point because these calculations are all performed at . Given that heavy-ion collisions can have fluctuations of , it may mean that it is possible to experience this quantum critical point at very low beam energy heavy-ion collisions. We have yet to study how the location/existence of this quantum triple critical point varies with temperature, but we leave that for the next stage of CMF++ when the temperature is included.
In the C3 coupling, as decreases and the magnitude of the strangeness fraction increases, the phase transition to the strange-dominated phase rapidly moves to lower values of (one can see this clearly in panel c) as the green/yellow colors switch abruptly to red). The shift in the phase transition with happens because the pressure increases with decreasing (for a given ) due to the increasing amount of strange particles, which also increases the baryon density (panel e), filling up the Fermi sea. The shift in the phase transition is not related to the stiffening of the EoS (by stiffening we mean an increase in with respect to because instead hyperons typically soften the EoS). Meanwhile, in the light hadronic phase, the pressure is not as strongly affected, since only the hyperon appears.
We find that for C3 there is already some strangeness (likely the ’s) that switches on during the light hadronic phase. One can see that in panel c) in Figure 15 since the blue (consistent with ) shifts to green/yellow (consistent with to ). Then for one would go to a light quark phase (strangeness appears to play almost no role in the C3 quark phase for positive ) or into the strangeness-dominated hadronic phase for .
In the strangeness-dominated phase, both the octet hyperons and their corresponding spin 3/2 excitations are present, as well as the s and the baryons. The light-to-strangeness-dominated hadronic phase transition was already present in the case. However, at this hadronic phase transition was only observed in the metastable region when quarks were present (see Figure 6) but in the absence of quarks, the hadronic phase transition appeared in the stable region (see Figure 12).
In contrast, for C4 in Figure 16, we see that in panel c) for the entire regime where in panel b). However, in the quarks phase for C4, a large, negative is strongly correlated with a large, negative , such that large switches on a large number of strange quarks. In the regime in the quark phase for C4, we find almost no strangeness, such that light quarks dominate.
We have already discussed quite a bit, but here we briefly discuss the relationship between and . For these results, we have made the assumption that , which corresponds to isospin symmetric matter, where (as already discussed) . This clearly holds in panels c) and d) in Figs. 15-16. In fact, this relationship is easy to calculate at the limit of that corresponds to , which corresponds exactly to what can be seen in the bottom right-hand corner of panels c) and d) in Figs. 15-16.
We do not show the relationship directly but rather in panel f) because it is easier to see changes in the degrees of freedom and the influence of phase transition in a derivative plot. Due to the complexities that appear in across the plane, we discuss the C3 and C4 couplings separately.
The C3 coupling has a number of new particles that switch on as increases (this was previously apparent through the population plots at the limit of ). Here in the plane we can see new particles switch on through discontinuities in the color spectra, which correspond to kinks in . For C3, the particle is the first one to appear, at MeV at . The discontinuity in color (or rather kink in ) is correlated with as we vary such that large, negative values of see the kink at low whereas large positive values of see the kink at large . Additionally, the sharpness of the kink changes with , negative values of have a smoother kink vs positive values of that sharpen the kink. We can understand the correlated behavior in the location of the kink in because decreases the energy required to produce strange baryons, allowing for them to appear at lower . For positive there is a very small kink near the deconfinement phase transition, which is a result of the onset of baryons.
Following the hadronic phase for C3, there is a first-order phase transition into quarks. Because allows for the possibility of new strange states (at lower ), another kink appears, followed very closely behind by a jump to (not shown here because the Maxwell construction is not shown) at the first or second-order phase transition marking the onset of the strangeness-dominated phase. In the negative region, the kink near the hadronic phase transition is due to the onset of the particles, which now appear before s (the one at lower is due to the ’s). The lack of points in the large , positive region indicates that no convergence is found when solving the field equations Eq. 46.
The C4 coupling results for in panel f) are significantly easier to understand. Generally, is stiffer than the other couplings and we can see that steady increase up until the quark phase appears. At the quark phase, there is a sharp drop in that stays close to the conformal limit of . For we do see that there is the possibility of hyperons switching on in the hadronic phase when a small discontinuity corresponding to a kink in can be seen (although the kink appears rather smooth). Even in the quark phase, there appears a small kink as well, but it is quite difficult to see.
IV.3
Now we analyze the interplay of the charge chemical potential, , and , while holding . The charge chemical potential is directly related to the appearance of electrically charged particles and, unlike the strangeness potential, it breaks the degeneracy between particles of the same family. From Eq. 50, we can expect a significantly negative to reduce the overall net positive charge of the system, bringing it farther away from symmetric nuclear matter. How this shift away from symmetric nuclear matter depends on the degrees of freedom of the system, their couplings, and interactions. For instance, for a system of just neutrons and protons, suppresses the of the system such that either fewer protons appear or more neutrons appear. For a system that allows for richer hadron chemistry, a negative suppresses the appearance of positively charged particles and/or enhances the appearance of negatively charged particles. We explore only the region of negative , since this is the physically relevant region both for heavy-ion collisions where (leading to small ) and for asymmetric nuclear matter found in neutron stars (leading to large ), where equilibrium is achieved.
IV.3.1 C3 and C4 with baryon octet + decuplet + quarks
The figures follow the same structure as in Sec. IV.2, but where we have changed the y-axis from to . We show the mean fields a), the deconfinement order parameter b), charge d) and strange c) fractions, baryon density e) and speed of sound squared f) for C3 in Figure 17 and for C4 in Figure 18.
At finite , the mean-field (panel a)) still generally decreases with , within a given phase of matter. However, we can see that a non-monotonic behavior in can appear across . In the case of the C4 coupling, we see a non-monotonic behavior in where there is a minimum in the hadronic phase right before the deconfinement phase transition (see panel b)) wherein then increases in the quark phase, before decreasing again at high . In the case of the C3 coupling, this behavior is more complex and depends on because the existence of the deconfinement phase transition only appears for MeV. In the C3 coupling, the strangeness-dominated hadronic phase further complicates this non-monotonic behavior of . For instance, if one fixes MeV, one can see that steadily decreases until just before MeV where a sharp drop appears (going into the strangeness-dominated hadronic phase). Then, going towards higher , we see a sharp increase in at the onset of quark deconfinement at around MeV.
At negative MeV for both the C3 and C4 couplings, a small region with is seen (dark blue seen at very low ), which indicates the presence of the liquid-gas transition. When we see , this is an indication that the vacuum solution is the stable solution. Of course, the actual liquid-gas phase transition should not go to a vacuum phase but rather switch to nuclei (which are not included in our model).
As previously discussed, we use the baryon density (panel e)) to determine the existence of a first-order phase transition (a jump in seen in ) and/or the order of any given phase transition through the susceptibilities. In the C3 coupling scheme, we can immediately see in Figure 17 two first-order lines that appear when . We can then use a combination of the deconfinement order parameter in panel b) (recall that if we are still in a hadronic phase) and the strangeness fraction in panel c) (recall that in the strangeness-dominated phase) to disentangle the deconfinement phase transition vs a first-order phase transition into the strangeness-dominated phase. At exactly these two first-order lines converge into a triple point such that one goes directly from the light hadronic phase into the deconfined quark phase (see that the point where the lines converge also corresponds to an order parameter of ). However, at it is clear that as one increases , then one first reaches a strangeness-dominated phase, and then only at even higher is the quark deconfinement phase reached. If MeV, then the quark deconfinement phase transition disappears entirely because there is a region where the code finds no solution. The critical of the hadronic phase transition moves to lower values as decreases since the pressure in the strangeness-dominated phase rises while in light hadronic phase it decreases (both at a given ) due to the larger amount of negatively charged particles in the light hadronic phase, which soften the EoS (pressure vs. energy density).
For C4, there is no hadronic phase transition. Thus, we never reach a strangeness-dominated hadronic phase in the C4 coupling. However, we can still use in panel f) to identify new particles that have switched on and understand the role that values have in the possibility of opening up these new particles. We find that lowering the charge chemical potential to more negative values moves the onset of the to lower , as marked by the strangeness fraction in panel c) and the discontinuity in the speed of sound, shown in panel f), starting around MeV and MeV. The is affected by the charge chemical potential, even though it is not charged, due to the coupling between the and mesons in C4, Eq. 46, which increases the effective chemical of the , moving its onset to a lower as decreases.
In the C3 coupling, the location of the deconfinement transition in is anti-correlated with (in other words more negative leads to a phase transition at a higher ). The anti-correlation of the location of the critical for deconfinement with occurs because the quark phase presents a very large density and has a lower pressure for lower (for fixed ). Within the deconfined quark phase, we find that at low absolute value of there are very few strange quarks. However, as the absolute value of increases, it opens up more strange quarks (because they carry electric charge ). Eventually, at MeV, the code no longer finds solutions consistent with a quark phase such that solutions end with the strangeness-dominated hadronic phase.
For the deconfinement phase transition, we find an opposite effect when it comes to the correlations between the critical and for C4 coupling. At we find that the deconfinement phase transition occurs at the maximum in panel b) for the order parameter . As becomes increasingly negative, then there is a slow shift into the critical for deconfinement to lower values. Because C4 has a much larger quark phase than C3, we can see more interesting effects that appear in the quark phase at large . Both the strange and down quarks are preferred for low because they carry a negative charge, whereas the positive up quarks are suppressed. The consequence of a dominant down/strange quark phase is that becomes very negative and as well in the quark phase. This isospin asymmetry causes the pressure to rise as a function of decreasing (for a given ).
While we have discussed how we can use in panel c) and in panel d) in Figure 17 to interpret the underlying properties of the phase of matter, here we discuss their general behavior across the phase space. Starting with the C3 coupling, we find that is mostly 0 in the light hadronic phase. Only after MeV do we begin to see a slightly negative due to the appearance of the . However, opens up the strangeness-dominated hadronic phase in the range of MeV (depending on ) wherein where almost all hadrons carry strangeness (and some carry multiple strangeness). Then the deconfined quark phase is significantly less strange and only has some significant contribution of strange quarks when is very negative.
The behavior is quite different than . The plot for C3 in panel d) has different regimes of interest. Unsurprisingly, close to one is close to the symmetric nuclear matter and (the one exception is in the strangeness-dominated phase where, even at low , ). Then at low and very negative , we find that because pure neutron matter is reached. For and intermediate (still in the light hadronic phase) we find small, positive values of . In this regime protons have switched on and eventually some other particles like ’s (that are neutral but still suppress because they increase ). At very largely negative and large but still, in the light hadronic phase, we can even see . The strangeness-dominated phase generally has mostly negative values of because protons are heavily suppressed while strange baryons with a negative charge are preferred. Finally, in the quark phase, the value of strongly depends on the value of . At low , the quark phase is primarily an even mix of up and down quarks such that , but as becomes more negative then strange and down quarks are preferred, which decreases .
For the C4 coupling in Figure 18 we find that is essentially zero for the light hadronic phase and only because significantly negative once the deconfined quark phase is reached. Then at large and very negative , there is a region where strange quarks are preferred such that .
The for C4 reaches 0 for the pure neutron matter region at low and very negative . Then at large (or smaller for large ) it is close to 0.5 in the light hadronic phase. For the quark phase, only at do we find . As becomes more negative in the quark phase, we find slowly becomes smaller until it eventually becomes negative (and then more and more negative).
Next, let us discuss the properties of in panel e) for both the C3 and C4 couplings. We find that in the light hadronic phase, both have very small values of that slowly increase with (regardless of ). For C3 in the strangeness-dominated phase, the new hadronic states open up new degrees of freedom, leading to much larger . For the deconfined quark phase in C3, we find very large values of such that they are much likely well beyond the reach of neutron stars. In contrast, C4 goes directly to the quark phase at large . In the quark phase, we find a large jump in from hadrons into quarks across the phase transition. Then steadily increases with in the quark phase, independently of .
Next, we can use to determine when new hadronic species are switching on in the C3 coupling due to the kinks/discontinuities that appear. From the speed of sound plot (in panel f)), we can identify how the appearance of different baryons changes with .
For C3, at MeV the first discontinuity in as one increases is due to the ’s. This is clear because changes at this point, but is only mildly affected. Then we see a smaller kink that appears that causes a to become more negative and to slightly decrease. Then, the strangeness-dominated phase appears around MeV, which switches on the , , baryons, such that and . There is then one final discontinuity in that indicates the quark deconfinement phase transition.
For MeV an even richer hadronic phase appears and many new states open up. At low we have a discontinuity, which occurs as a transition from pure neutron matter into one that includes protons as well (as often seen in neutron star calculations, since the proton chemical potential is less than the neutron one, ). Then, the baryon switches on at MeV, with an associated discontinuity that is only slightly modified by due to the changes in the meson fields. Three more discontinuities appear at sufficiently low , associated with the onset of the , , and particles, respectively. Then a sharp transition in appears for the first-order phase transition into the strangeness-dominated regime. In summary, the order or appearance is , , , , in the hadronic phase, and then, in the strangeness-dominated, the ’s become dominant, but all the octet particles appear, as well as the ’s.
Overall, we find that is pretty dependent on for both C3 and C4. At for C3 we see one kink at around MeV and then the triple critical point at around MeV that leads to some kinks with a very brief dip in at the phase transition(s). As becomes more negative, then has a large number of kinks that appear at may see 1-2 regions (depending on the exactly ) where . For C4 we already discussed previously that symmetric nuclear matter has a very stiff EoS with no kinks in the hadronic phase (and only a very tiny one in the quark phase). However, as one goes to a dip in appears when the proton switches one and a kink appears when the appears as well, followed by the jump across the first-order deconfinement phase transition.
IV.4
Up until this point, we have always set one or two chemical potentials to zero. Due to the limitations of 2-dimensional plots, it is significantly more difficult to vary all 3 chemical potentials at once in a meaningful manner. Thus, in the following, we hold either while varying or hold while varying .
IV.4.1 C3 with baryon octet + decuplet + quarks
In this section, we consider the baryon octet, decuplet, and quarks to allow for the widest possible range of degrees of freedom while we explore the interplay of . We focus on the C3 coupling, which has the most distinct features and has a larger variety of particle species that regularly appear in the EoS compared to C4. Note that, because heavy particles tend to soften the EoS, the lower amount of heavy particles in the C4 coupling generally has an easier time reproducing astrophysical constraints of neutron star masses and radii [53].
In Figure 19, the charge chemical potential is kept fixed at MeV, which is a typical value in hadronic neutron star matter, and the strangeness chemical potential is varied. Then in Figure 20 we study the opposite scenario wherein the strangeness chemical potential is fixed at MeV, and the is varied. In each of the panels, we show the mean-field meson a), the order parameter b), charge d) and strange c) fractions, baryon density e) and speed of sound squared f).
Let us begin with the fixed MeV in Figure 19. To explain the , we begin at the low end and work our way up to larger . At very low and largely negative we find a very tiny phase of matter that forms a blue triangle, i.e. . While one might be tempted to assume that this is a vacuum state since it occurs at low , we later see that it is dominated by strange baryons (although at very low density). Then at larger the behavior of has the same qualitative appearance as what was seen in panel a) from Figure 15 and even has similar values for the hadronic phase, the strangeness-dominated phase, and the quark deconfined phase.
However, the exact location of these phases of matter and the shape of their first-order phase transitions are different at finite . The quark deconfinement in Figure 19 occurs in a region of positive and larger and the strangeness-dominated hadronic phase occurs at lower values of the compared to what we saw previously in Figure 15. The main role of the negative is to push the triple point from to MeV. A consequence of the change in the location of the triple point at finite is that the deconfinement phase transition shifts to lower . At, e.g., MeV, only the hadronic phase transition occurs, as confirmed by the deconfinement field in panel b).
For MeV the deconfinement transition shifts to lower in the positive region (when compared with Figure 15), due to the higher amount of negatively charged particles at the same in the hadronic phase. These extra negatively charged strange particles increase the pressure for a given (see panels c), d) and e) in Figure 19) of the hadronic phase (while softening the EoS). In contrast, in the quark phase, finite (Figure 19 vs. Figure 15) only changes the ratio of up and down quarks, since the phase transition occurs before the onset of strangeness (panel c)).
One conclusion that we can draw from these plots is if C3 correctly describes the matter within neutron stars and heavy-ion collisions, then we could expect that neutron stars have deconfinement phase transition at lower than what we expect in heavy-ion collisions (where is negative but normally smaller than what is seen in neutron stars). Additionally, because heavy-ion collisions can experience fluctuations in due to variations in the local strangeness content, it could be that some regions of the fluid see a first-order phase transition into a strangeness-dominated regime, but other parts of the fluids experience a first-order phase transition into quarks, and yet other parts may fluctuate from quarks into strange baryons (or vice versa).
When we hold MeV, fixed in Figure 20, we find that we remove the deconfinement phase transition entirely for C3 (see that in panel b) across all ). However, we can see that there is significant variation in across in that we see a small region where at low (panel a)), very negative , then a somewhat steadily decreasing values of followed by a sharp drop for . Since there is no deconfined phase transition, then this regime only contains hadronic states. This sharp drop in corresponds to the strangeness-dominated hadronic phase (see panel C)).
In the case of fixed MeV and varying in Figure 20, we find that the phase transition into the strangeness-dominated hadronic phase is at approximately the same location as what was previously seen for in Figure 17, it is just that the strangeness-dominated hadronic phase is now always the stable solution regardless of at finite (to the exclusion of the quark phase). Also, the strangeness-dominated phase is able to find solutions out to large as well when as compared to the case (e.g., Figure 17).
In both Figure 19 and Figure 20, the appears before the proton. For MeV when MeV and for MeV when MeV. The strangeness and charge fraction, displayed in panels c) and d) of Figure 19 and Figure 20, show the same tendency as discussed previously: in each phase, both present an overall decrease with after the onset of the strangeness and with decreasing or at the lowest ’s analyzed. This feature is mostly due to the substantial contribution of the and hyperons, which become more relevant than protons at sufficiently low . We can tell the strong roles of these hyperons by the very negative in panel c) and the also negative or small in panel d). In fact, for C3 we find an interesting phenomena at MeV and MeV, where there is a phase where only the particle appears when MeV (Figure 19).
The baryon density, shown in panel e) in both Figure 19 and Figure 20 increases with and with decreasing or , in all phases. Additionally, the value of the density in the quark and the strangeness-dominated phases are very similar immediately after the transition, with .
For MeV, the vertical discontinuities in the speed of sound, displayed in the leftmost part of panel f) in Figure 19, indicate the appearance of the proton and the . The is then followed by the ’s, , and as one increases . However, , and do not appear above MeV, but as decreases, they eventually overcome protons and neutrons. While we do not show the population plots in 3D due to their immense complexity, in the hyperonic phase the most abundant particles are in order of abundance , , , n, , , , , , , p, , , , and .
For MeV (Figure 20), the appearance of protons and are identified as the first two kinks from the figure. The change on the onset of the with is due to the coupling of the vector fields Eq. 46. After the ’s appear, the and appear as one increases in the MeV. Below MeV, the eventually becomes more abundant than the and the proton. Across the first-order phase transition, into the strangeness-dominated phase, we see that , and the population follows the one in Figure 20.
V Final remarks
We presented in this paper new results from the Chiral Mean Field (CMF) model at vanishing temperatures calculated using the new CMF++ code that will be integrated as a module in the MUSES cyberinfrastructure and be available to the public soon as open-source [137]. The runtime improved more than 4 orders of magnitude in 3 dimensions (), when compared to the legacy Fortran code, while showing good agreement in a wide variety of configurations, couplings, thermodynamic quantities, etc. Numerical improvements also allowed us to calculate higher-order susceptibilities for the first time, which allowed us to identify first, second (quantum critical points), and third-order phase transitions for the first time within CMF.
For the sake of clarity, we have presented a thorough review of the CMF model, including for the first time a thorough derivation of the formalism focusing on . The derivation includes the Lagrangian density, equations of motion, and thermodynamic properties. We are now able to calculate the unstable and metastable regimes between first-order phase transitions within CMF++ and have outlined the complex stability criteria for 3 chemical potentials (BSQ) in the grand canonical ensemble for both (7 stability conditions) and finite (15 stability conditions).
While exploring different CMF vector couplings (C1-C4 defined in the paper) and a large 3-Dimensional parameter space in chemical potential, , , , we were able to identify different phases of matter that may appear at : pure neutron matter, light hadronic matter (protons, neutrons, and sometimes ’s), strangeness and -dominated hadronic phase ( baryons dominate but other baryons may appear such as ’s, ’s, and even ’s), and deconfined quark matter (although there is a strong flavor dependence of the quark phase such that we have seen quarks phases dominated by light quarks, others by strange quarks, and yet others by down and strange quarks). The strangeness and dominated hadronic phase is especially interesting because often a first-order phase transition appears into this phase of matter such that the EoS may look very similar to one that goes from a light hadronic phase into a deconfined quark phase (although in our model we find that the phase transition into the strangeness-dominated phase tends to be smaller than into quarks). This strangeness-dominated hadronic phase is triggered by large amounts of hyperons or heavier (spin 3/2) baryons appearing in the system. We found that the strangeness-dominated phase may be a stable solution of the EoS or even hidden within the metastable regime.
Within one of the vector couplings we studied, we found that a tricritical point appeared in the phase space such that at high (depending on the ) one could either find a first-order phase transition into deconfined quarks or into a strangeness dominated hadronic phase (or even a first-order phase transition separating quarks from a strangeness dominated hadronic phase). Given that heavy-ion collisions can have local fluctuations of due to gluons splitting into quark anti-quark pairs, then it is not unreasonable to think that effects of this tricritical point could be potentially observed in low-energy heavy-ion collisions.
Our next step is to develop the finite version of the CMF++ code, which is already underway. The finite version of CMF++ will allow direct comparisons to lattice QCD as well as the potential to couple to other finite codes that reach large such as the holography EoS [29].
Acknowledgements.
The authors acknowledge useful discussions with Jorge Noronha and Johannes Jahan on the stability constraints. All authors of this paper are a part of the MUSES collaboration, which is supported by the NSF under OAC-2103680. Additional support for the collaboration members includes: RH is supported under OAC-2005572, and OAC-2004879. J.N.H. acknowledges support from the US-DOE Nuclear Science Grant No. DE-SC0023861, and within the framework of the Saturated Glue (SURGE) Topical Theory Collaboration. R.K., J.P., and V.D were funded by the National Science Foundation, grant number PHY1748621. The work is also supported by the Illinois Campus Cluster, a computing resource that is operated by the Illinois Campus Cluster Program (ICCP) in conjunction with the National Center for Supercomputing Applications (NCSA), which is supported by funds from the University of Illinois Urbana Champaign. Figures in the manuscript were produced using matplotlib [138, 139].References
- Heinz and Snellings [2013] U. Heinz and R. Snellings, Collective flow and viscosity in relativistic heavy-ion collisions, Ann. Rev. Nucl. Part. Sci. 63, 123 (2013), arXiv:1301.2826 [nucl-th] .
- Lattimer and Prakash [2005] J. M. Lattimer and M. Prakash, The Ultimate energy density of observable cold matter, Phys. Rev. Lett. 94, 111101 (2005), arXiv:astro-ph/0411280 .
- Ivanenko and Kurdgelaidze [1970] D. D. Ivanenko and D. F. Kurdgelaidze, Quark stars, Sov. Phys. J. 13, 1015 (1970).
- Annala et al. [2020] E. Annala, T. Gorda, A. Kurkela, J. Nättilä, and A. Vuorinen, Evidence for quark-matter cores in massive neutron stars, Nature Phys. 16, 907 (2020), arXiv:1903.09121 [astro-ph.HE] .
- Most et al. [2019] E. R. Most, L. J. Papenfort, V. Dexheimer, M. Hanauske, S. Schramm, H. Stöcker, and L. Rezzolla, Signatures of quark-hadron phase transitions in general-relativistic neutron-star mergers, Phys. Rev. Lett. 122, 061101 (2019), arXiv:1807.03684 [astro-ph.HE] .
- Bauswein et al. [2019] A. Bauswein, N.-U. F. Bastian, D. B. Blaschke, K. Chatziioannou, J. A. Clark, T. Fischer, and M. Oertel, Identifying a first-order phase transition in neutron star mergers through gravitational waves, Phys. Rev. Lett. 122, 061102 (2019), arXiv:1809.01116 [astro-ph.HE] .
- Stephanov et al. [1998] M. A. Stephanov, K. Rajagopal, and E. V. Shuryak, Signatures of the tricritical point in QCD, Phys. Rev. Lett. 81, 4816 (1998), arXiv:hep-ph/9806219 .
- Elliott et al. [2013] J. B. Elliott, P. T. Lake, L. G. Moretto, and L. Phair, Determination of the coexistence curve, critical temperature, density, and pressure of bulk nuclear matter from fragment emission data, Phys. Rev. C 87, 054622 (2013).
- Capellino et al. [2022] F. Capellino, A. Beraudo, A. Dubla, S. Floerchinger, S. Masciocchi, J. Pawlowski, and I. Selyuzhenkov, Fluid-dynamic approach to heavy-quark diffusion in the quark-gluon plasma, Phys. Rev. D 106, 034021 (2022), arXiv:2205.07692 [nucl-th] .
- Hempel et al. [2013] M. Hempel, V. Dexheimer, S. Schramm, and I. Iosilevskiy, Noncongruence of the nuclear liquid-gas and deconfinement phase transitions, Phys. Rev. C 88, 014906 (2013), arXiv:1302.2835 [nucl-th] .
- Alford and Harris [2018] M. G. Alford and S. P. Harris, Beta equilibrium in neutron star mergers, Phys. Rev. C 98, 065806 (2018), arXiv:1803.00662 [nucl-th] .
- Adam et al. [2017] J. Adam et al. (ALICE), Enhanced production of multi-strange hadrons in high-multiplicity proton-proton collisions, Nature Phys. 13, 535 (2017), arXiv:1606.07424 [nucl-ex] .
- Adamczyk et al. [2017] L. Adamczyk et al. (STAR), Bulk Properties of the Medium Produced in Relativistic Heavy-Ion Collisions from the Beam Energy Scan Program, Phys. Rev. C 96, 044904 (2017), arXiv:1701.07065 [nucl-ex] .
- Adamczewski-Musch et al. [2019] J. Adamczewski-Musch et al. (HADES), Probing dense baryon-rich matter with virtual photons, Nature Phys. 15, 1040 (2019).
- Sorensen and Koch [2021] A. Sorensen and V. Koch, Phase transitions and critical behavior in hadronic transport with a relativistic density functional equation of state, Phys. Rev. C 104, 034904 (2021), arXiv:2011.06635 [nucl-th] .
- Oliinychenko et al. [2023] D. Oliinychenko, A. Sorensen, V. Koch, and L. McLerran, Sensitivity of Au+Au collisions to the symmetric nuclear matter equation of state at 2–5 nuclear saturation densities, Phys. Rev. C 108, 034908 (2023), arXiv:2208.11996 [nucl-th] .
- Yao et al. [2024] N. Yao, A. Sorensen, V. Dexheimer, and J. Noronha-Hostler, Structure in the speed of sound: From neutron stars to heavy-ion collisions, Phys. Rev. C 109, 065803 (2024), arXiv:2311.18819 [nucl-th] .
- Carzon et al. [2022] P. Carzon, M. Martinez, M. D. Sievert, D. E. Wertepny, and J. Noronha-Hostler, Monte Carlo event generator for initial conditions of conserved charges in nuclear geometry, Phys. Rev. C 105, 034908 (2022), arXiv:1911.12454 [nucl-th] .
- Plumberg et al. [2024] C. Plumberg et al., BSQ Conserved Charges in Relativistic Viscous Hydrodynamics solved with Smoothed Particle Hydrodynamics (2024), arXiv:2405.09648 [nucl-th] .
- Troyer and Wiese [2005] M. Troyer and U.-J. Wiese, Computational complexity and fundamental limitations to fermionic quantum Monte Carlo simulations, Phys. Rev. Lett. 94, 170201 (2005), arXiv:cond-mat/0408370 .
- Muroya et al. [2003] S. Muroya, A. Nakamura, C. Nonaka, and T. Takaishi, Lattice QCD at finite density: An Introductory review, Prog. Theor. Phys. 110, 615 (2003), arXiv:hep-lat/0306031 .
- de Forcrand [2009] P. de Forcrand, Simulating QCD at finite density, PoS LAT2009, 010 (2009), arXiv:1005.0539 [hep-lat] .
- Monnai et al. [2019] A. Monnai, B. Schenke, and C. Shen, Equation of state at finite densities for QCD matter in nuclear collisions, Phys. Rev. C 100, 024907 (2019), arXiv:1902.05095 [nucl-th] .
- Noronha-Hostler et al. [2019] J. Noronha-Hostler, P. Parotto, C. Ratti, and J. M. Stafford, Lattice-based equation of state at finite baryon number, electric charge and strangeness chemical potentials, Phys. Rev. C 100, 064910 (2019), arXiv:1902.06723 [hep-ph] .
- Borsányi et al. [2021] S. Borsányi, Z. Fodor, J. N. Guenther, R. Kara, S. D. Katz, P. Parotto, A. Pásztor, C. Ratti, and K. K. Szabó, Lattice QCD equation of state at finite chemical potential from an alternative expansion scheme, Phys. Rev. Lett. 126, 232001 (2021), arXiv:2102.06660 [hep-lat] .
- Gubser and Nellore [2008] S. S. Gubser and A. Nellore, Mimicking the QCD equation of state with a dual black hole, Phys. Rev. D 78, 086007 (2008), arXiv:0804.0434 [hep-th] .
- DeWolfe et al. [2011] O. DeWolfe, S. S. Gubser, and C. Rosen, A holographic critical point, Phys. Rev. D 83, 086005 (2011), arXiv:1012.1864 [hep-th] .
- Rougemont et al. [2016] R. Rougemont, A. Ficnar, S. Finazzo, and J. Noronha, Energy loss, equilibration, and thermodynamics of a baryon rich strongly coupled quark-gluon plasma, JHEP 04, 102, arXiv:1507.06556 [hep-th] .
- Hippert et al. [2023] M. Hippert, J. Grefa, T. A. Manning, J. Noronha, J. Noronha-Hostler, I. Portillo Vazquez, C. Ratti, R. Rougemont, and M. Trujillo, Bayesian location of the QCD critical point from a holographic perspective (2023), arXiv:2309.00579 [nucl-th] .
- Gunkel and Fischer [2021] P. J. Gunkel and C. S. Fischer, Locating the critical endpoint of QCD: Mesonic backcoupling effects, Phys. Rev. D 104, 054022 (2021), arXiv:2106.08356 [hep-ph] .
- Fu et al. [2020] W.-j. Fu, J. M. Pawlowski, and F. Rennecke, QCD phase structure at finite temperature and density, Phys. Rev. D 101, 054032 (2020), arXiv:1909.02991 [hep-ph] .
- Basar [2021] G. Basar, Universality, Lee-Yang Singularities, and Series Expansions, Phys. Rev. Lett. 127, 171603 (2021), arXiv:2105.08080 [hep-th] .
- Dimopoulos et al. [2022] P. Dimopoulos, L. Dini, F. Di Renzo, J. Goswami, G. Nicotra, C. Schmidt, S. Singh, K. Zambello, and F. Ziesché, Contribution to understanding the phase structure of strong interaction matter: Lee-Yang edge singularities from lattice QCD, Phys. Rev. D 105, 034513 (2022), arXiv:2110.15933 [hep-lat] .
- Sorensen and Sorensen [2024] A. Sorensen and P. Sorensen, Locating the critical point for the hadron to quark-gluon plasma phase transition from finite-size scaling of proton cumulants in heavy-ion collisions (2024), arXiv:2405.10278 [nucl-th] .
- Andersen and Strickland [2002] J. O. Andersen and M. Strickland, The Equation of state for dense QCD and quark stars, Phys. Rev. D 66, 105001 (2002), arXiv:hep-ph/0206196 .
- Fraga et al. [2014] E. S. Fraga, A. Kurkela, and A. Vuorinen, Interacting quark matter equation of state for compact stars, Astrophys. J. Lett. 781, L25 (2014), arXiv:1311.5154 [nucl-th] .
- Kurkela and Vuorinen [2016] A. Kurkela and A. Vuorinen, Cool quark matter, Phys. Rev. Lett. 117, 042501 (2016), arXiv:1603.00750 [hep-ph] .
- Holt et al. [2010] J. W. Holt, N. Kaiser, and W. Weise, Density-dependent effective nucleon-nucleon interaction from chiral three-nucleon forces, Phys. Rev. C 81, 024002 (2010), arXiv:0910.1249 [nucl-th] .
- Tews et al. [2013] I. Tews, T. Krüger, K. Hebeler, and A. Schwenk, Neutron matter at next-to-next-to-next-to-leading order in chiral effective field theory, Phys. Rev. Lett. 110, 032504 (2013), arXiv:1206.0025 [nucl-th] .
- Weinberg [1979] S. Weinberg, Phenomenological Lagrangians, Physica A 96, 327 (1979).
- Kumar et al. [2024a] R. Kumar et al. (MUSES), Theoretical and experimental constraints for the equation of state of dense and hot matter, Living Rev. Rel. 27, 3 (2024a), arXiv:2303.17021 [nucl-th] .
- Higgs [1964] P. W. Higgs, Broken Symmetries and the Masses of Gauge Bosons, Phys. Rev. Lett. 13, 508 (1964).
- Nambu and Jona-Lasinio [1961] Y. Nambu and G. Jona-Lasinio, Dynamical model of elementary particles based on an analogy with superconductivity. II., Phys. Rev. 124, 246 (1961).
- Gell-Mann and Levy [1960] M. Gell-Mann and M. Levy, The axial vector current in beta decay, Nuovo Cim. 16, 705 (1960).
- Hatsuda and Prakash [1989] T. Hatsuda and M. Prakash, Parity Doubling of the Nucleon and First Order Chiral Transition in Dense Matter, Phys. Lett. B 224, 11 (1989).
- Detar and Kunihiro [1989] C. E. Detar and T. Kunihiro, Linear Model With Parity Doubling, Phys. Rev. D 39, 2805 (1989).
- Weinberg [1968] S. Weinberg, Nonlinear realizations of chiral symmetry, Phys. Rev. 166, 1568 (1968).
- Papazoglou et al. [1999] P. Papazoglou, D. Zschiesche, S. Schramm, J. Schaffner-Bielich, H. Stoecker, and W. Greiner, Nuclei in a chiral SU(3) model, Phys. Rev. C 59, 411 (1999), arXiv:nucl-th/9806087 .
- Steinheimer et al. [2010] J. Steinheimer, V. Dexheimer, H. Petersen, M. Bleicher, S. Schramm, and H. Stoecker, Hydrodynamics with a chiral hadronic equation of state including quark degrees of freedom, Phys. Rev. C 81, 044913 (2010), arXiv:0905.3099 [hep-ph] .
- Steinheimer et al. [2008] J. Steinheimer, M. Bleicher, H. Petersen, S. Schramm, H. Stocker, and D. Zschiesche, (3+1)-dimensional hydrodynamic expansion with a critical point from realistic initial conditions, Phys. Rev. C 77, 034901 (2008), arXiv:0710.0332 [nucl-th] .
- Dexheimer and Schramm [2008] V. Dexheimer and S. Schramm, Proto-Neutron and Neutron Stars in a Chiral SU(3) Model, Astrophys. J. 683, 943 (2008), arXiv:0802.1999 [astro-ph] .
- Schürhoff et al. [2010] T. Schürhoff, S. Schramm, and V. Dexheimer, Neutron stars with small radii – the role of delta resonances, Astrophys. J. Lett. 724, L74 (2010), arXiv:1008.0957 [astro-ph.SR] .
- Dexheimer et al. [2015] V. Dexheimer, R. Negreiros, and S. Schramm, Reconciling Nuclear and Astrophysical Constraints, Phys. Rev. C 92, 012801 (2015), arXiv:1503.07785 [astro-ph.HE] .
- Dexheimer et al. [2019] V. Dexheimer, R. de Oliveira Gomes, S. Schramm, and H. Pais, What do we learn about vector interactions from GW170817?, J. Phys. G 46, 034002 (2019), arXiv:1810.06109 [nucl-th] .
- Roark et al. [2019] J. Roark, X. Du, C. Constantinou, V. Dexheimer, A. W. Steiner, and J. R. Stone, Hyperons and quarks in proto-neutron stars, Mon. Not. Roy. Astron. Soc. 486, 5441 (2019), arXiv:1812.08157 [astro-ph.HE] .
- Motornenko et al. [2020] A. Motornenko, J. Steinheimer, V. Vovchenko, S. Schramm, and H. Stoecker, Equation of state for hot QCD and compact stars from a mean field approach, Phys. Rev. C 101, 034904 (2020), arXiv:1905.00866 [hep-ph] .
- Jakobus et al. [2023] P. Jakobus, B. Müller, A. Heger, S. Zha, J. Powell, A. Motornenko, J. Steinheimer, and H. Stoecker, Gravitational Waves from a Core g Mode in Supernovae as Probes of the High-Density Equation of State, Phys. Rev. Lett. 131, 191201 (2023), arXiv:2301.06515 [astro-ph.HE] .
- Negreiros et al. [2010] R. Negreiros, V. A. Dexheimer, and S. Schramm, Modeling Hybrid Stars with an SU(3) non-linear sigma model, Phys. Rev. C 82, 035803 (2010), arXiv:1006.0380 [astro-ph.SR] .
- Dexheimer et al. [2013] V. Dexheimer, J. Steinheimer, R. Negreiros, and S. Schramm, Hybrid Stars in an SU(3) parity doublet model, Phys. Rev. C 87, 015804 (2013), arXiv:1206.3086 [astro-ph.HE] .
- Most et al. [2020] E. R. Most, L. Jens Papenfort, V. Dexheimer, M. Hanauske, H. Stoecker, and L. Rezzolla, On the deconfinement phase transition in neutron-star mergers, Eur. Phys. J. A 56, 59 (2020), arXiv:1910.13893 [astro-ph.HE] .
- Most et al. [2023] E. R. Most, A. Motornenko, J. Steinheimer, V. Dexheimer, M. Hanauske, L. Rezzolla, and H. Stoecker, Probing neutron-star matter in the lab: Similarities and differences between binary mergers and heavy-ion collisions, Phys. Rev. D 107, 043034 (2023), arXiv:2201.13150 [nucl-th] .
- Raduta et al. [2021] A. R. Raduta, F. Nacu, and M. Oertel, Equations of state for hot neutron stars, Eur. Phys. J. A 57, 329 (2021), arXiv:2109.00251 [nucl-th] .
- Raduta [2022] A. R. Raduta, Equations of state for hot neutron stars-II. The role of exotic particle degrees of freedom, Eur. Phys. J. A 58, 115 (2022), arXiv:2205.03177 [nucl-th] .
- Tsiopelas et al. [2024] S. Tsiopelas, A. Sedrakian, and M. Oertel, Finite-temperature equations of state of compact stars with hyperons: three-dimensional tables, Eur. Phys. J. A 60, 127 (2024), arXiv:2406.00484 [nucl-th] .
- Typel et al. [2022] S. Typel et al. (CompOSE Core Team), CompOSE Reference Manual, Eur. Phys. J. A 58, 221 (2022), arXiv:2203.03209 [astro-ph.HE] .
- Dexheimer et al. [2022] V. Dexheimer, M. Mancini, M. Oertel, C. Providência, L. Tolos, and S. Typel, Quick Guides for Use of the CompOSE Data Base, Particles 5, 346 (2022), arXiv:2311.04715 [nucl-th] .
- Franzon et al. [2016] B. Franzon, V. Dexheimer, and S. Schramm, A self-consistent study of magnetic field effects on hybrid stars, Mon. Not. Roy. Astron. Soc. 456, 2937 (2016), arXiv:1508.04431 [astro-ph.HE] .
- Dexheimer et al. [2021a] V. Dexheimer, K. D. Marquez, and D. P. Menezes, Delta baryons in neutron-star matter under strong magnetic fields, Eur. Phys. J. A 57, 216 (2021a), arXiv:2103.09855 [nucl-th] .
- Marquez et al. [2022] K. D. Marquez, M. R. Pelicer, S. Ghosh, J. Peterson, D. Chatterjee, V. Dexheimer, and D. P. Menezes, Exploring the effects of baryons in magnetars, Phys. Rev. C 106, 035801 (2022), arXiv:2205.09827 [astro-ph.HE] .
- Peterson et al. [2021] J. Peterson, V. Dexheimer, R. Negreiros, and B. G. Castanheira, Effects of Magnetic Fields in Hot White Dwarfs, Astrophys. J. 921, 1 (2021), arXiv:2105.03387 [astro-ph.SR] .
- Oertel et al. [2017] M. Oertel, M. Hempel, T. Klähn, and S. Typel, Equations of state for supernovae and compact stars, Rev. Mod. Phys. 89, 015007 (2017), arXiv:1610.03361 [astro-ph.HE] .
- Dexheimer and Schramm [2010] V. A. Dexheimer and S. Schramm, A Novel Approach to Model Hybrid Stars, Phys. Rev. C 81, 045201 (2010), arXiv:0901.1748 [astro-ph.SR] .
- Aoki et al. [2006] Y. Aoki, G. Endrodi, Z. Fodor, S. D. Katz, and K. K. Szabo, The Order of the quantum chromodynamics transition predicted by the standard model of particle physics, Nature 443, 675 (2006), arXiv:hep-lat/0611014 .
- Pisarski [2000] R. D. Pisarski, Quark gluon plasma as a condensate of SU(3) Wilson lines, Phys. Rev. D 62, 111501 (2000), arXiv:hep-ph/0006205 .
- Fukushima [2004] K. Fukushima, Chiral effective model with the Polyakov loop, Phys. Lett. B 591, 277 (2004), arXiv:hep-ph/0310121 .
- Steinheimer et al. [2011] J. Steinheimer, S. Schramm, and H. Stocker, The hadronic SU(3) Parity Doublet Model for Dense Matter, its extension to quarks and the strange equation of state, Phys. Rev. C 84, 045208 (2011), arXiv:1108.2596 [hep-ph] .
- Roark and Dexheimer [2018] J. Roark and V. Dexheimer, Deconfinement phase transition in proto-neutron-star matter, Phys. Rev. C 98, 055805 (2018), arXiv:1803.02411 [nucl-th] .
- Aryal et al. [2020] K. Aryal, C. Constantinou, R. L. S. Farias, and V. Dexheimer, High-Energy Phase Diagrams with Charge and Isospin Axes under Heavy-Ion Collision and Stellar Conditions, Phys. Rev. D 102, 076016 (2020), arXiv:2004.03039 [nucl-th] .
- Kumar et al. [2024b] R. Kumar, Y. Wang, N. C. Camacho, A. Kumar, J. Noronha-Hostler, and V. Dexheimer, Modern nuclear and astrophysical constraints of dense matter in a redefined chiral approach, Phys. Rev. D 109, 074008 (2024b), arXiv:2401.12944 [nucl-th] .
- [80] MUSES, https://musesframework.io/.
- Dexheimer et al. [2021b] V. Dexheimer, R. O. Gomes, T. Klähn, S. Han, and M. Salinas, GW190814 as a massive rapidly rotating neutron star with exotic degrees of freedom, Phys. Rev. C 103, 025808 (2021b), arXiv:2007.08493 [astro-ph.HE] .
- Schaffner et al. [1993] J. Schaffner, C. B. Dover, A. Gal, C. Greiner, and H. Stoecker, Strange hadronic matter, Phys. Rev. Lett. 71, 1328 (1993).
- Kubis and Kutschera [1997] S. Kubis and M. Kutschera, Nuclear matter in relativistic mean field theory with isovector scalar meson, Phys. Lett. B 399, 191 (1997), arXiv:astro-ph/9703049 .
- Hofmann et al. [2001] F. Hofmann, C. M. Keil, and H. Lenske, Density dependent hadron field theory for asymmetric nuclear matter and exotic nuclei, Phys. Rev. C 64, 034314 (2001), arXiv:nucl-th/0007050 .
- Gerstein et al. [1971] I. S. Gerstein, R. Jackiw, S. Weinberg, and B. W. Lee, Chiral loops, Phys. Rev. D 3, 2486 (1971).
- Coleman et al. [1969] S. R. Coleman, J. Wess, and B. Zumino, Structure of phenomenological Lagrangians. 1., Phys. Rev. 177, 2239 (1969).
- Callan et al. [1969] C. G. Callan, Jr., S. R. Coleman, J. Wess, and B. Zumino, Structure of phenomenological Lagrangians. 2., Phys. Rev. 177, 2247 (1969).
- Gross and Wilczek [1973] D. J. Gross and F. Wilczek, Ultraviolet Behavior of Nonabelian Gauge Theories, Phys. Rev. Lett. 30, 1343 (1973).
- Politzer [1973] H. D. Politzer, Reliable Perturbative Results for Strong Interactions?, Phys. Rev. Lett. 30, 1346 (1973).
- Zyla et al. [2020] P. A. Zyla et al. (Particle Data Group), Review of Particle Physics, PTEP 2020, 083C01 (2020).
- Schechter [1980] J. Schechter, Effective Lagrangian with Two Color Singlet Gluon Fields, Phys. Rev. D 21, 3393 (1980).
- Papazoglou et al. [1998] P. Papazoglou, S. Schramm, J. Schaffner-Bielich, H. Stoecker, and W. Greiner, Chiral Lagrangian for strange hadronic matter, Phys. Rev. C 57, 2576 (1998), arXiv:nucl-th/9706024 .
- Schechter and Ueda [1971] J. Schechter and Y. Ueda, Symmetry breaking and spin-zero mass spectrum, Phys. Rev. D 3, 168 (1971).
- Zschiehcse [2003] D. Zschiehcse, Excited Hadronic Matter in Chiral SU(3) Model, Thesis (2003).
- Heide et al. [1992] E. K. Heide, S. Rudaz, and P. J. Ellis, Implications of a modified glueball potential for nuclear matter, Phys. Lett. B 293, 259 (1992).
- Koch [1997] V. Koch, Aspects of chiral symmetry, Int. J. Mod. Phys. E 6, 203 (1997), arXiv:nucl-th/9706075 .
- Serot and Walecka [1997] B. D. Serot and J. D. Walecka, Recent progress in quantum hadrodynamics, Int. J. Mod. Phys. E 6, 515 (1997), arXiv:nucl-th/9701058 .
- Bonanno and Drago [2009] L. Bonanno and A. Drago, A Chiral lagrangian with Broken Scale: Testing the restoration of symmetries in astrophysics and in the laboratory, Phys. Rev. C 79, 045801 (2009), arXiv:0805.4188 [nucl-th] .
- Sasaki and Mishustin [2012] C. Sasaki and I. Mishustin, The Phase structure of a chiral model with dilatons in hot and dense matter, Phys. Rev. C 85, 025202 (2012), arXiv:1110.3498 [hep-ph] .
- Dexheimer et al. [2008a] V. Dexheimer, S. Schramm, and D. Zschiesche, Nuclear matter and neutron stars in a parity doublet model, Phys. Rev. C 77, 025803 (2008a), arXiv:0710.4192 [nucl-th] .
- Dexheimer et al. [2008b] V. Dexheimer, G. Pagliara, L. Tolos, J. Schaffner-Bielich, and S. Schramm, Neutron stars within the SU(2) parity doublet model, Eur. Phys. J. A 38, 105 (2008b), arXiv:0805.3301 [nucl-th] .
- Malik et al. [2024] T. Malik, V. Dexheimer, and C. Providência, Astrophysics and Nuclear Physics Informed Interactions in Dense Matter: Insights from PSR J0437-4715 (2024), arXiv:2404.07936 [nucl-th] .
- Sakurai [1969] J. J. Sakurai, Currents and Mesons (University of Chicago Press, Chicago, 1969).
- Dover and Gal [1985] C. B. Dover and A. Gal, HYPERON NUCLEUS POTENTIALS, Prog. Part. Nucl. Phys. 12, 171 (1985).
- Weissenborn et al. [2012] S. Weissenborn, D. Chatterjee, and J. Schaffner-Bielich, Hyperons and massive neutron stars: vector repulsion and SU(3) symmetry, Phys. Rev. C 85, 065802 (2012), [Erratum: Phys.Rev.C 90, 019904 (2014)], arXiv:1112.0234 [astro-ph.HE] .
- Lopes et al. [2023] L. L. Lopes, K. D. Marquez, and D. P. Menezes, Baryon coupling scheme in a unified SU(3) and SU(6) symmetry formalism, Phys. Rev. D 107, 036011 (2023), arXiv:2211.17153 [hep-ph] .
- Steinheimer and Schramm [2014] J. Steinheimer and S. Schramm, Do lattice data constrain the vector interaction strength of QCD?, Phys. Lett. B 736, 241 (2014), arXiv:1401.4051 [nucl-th] .
- Ratti et al. [2006] C. Ratti, M. A. Thaler, and W. Weise, Phases of QCD: Lattice thermodynamics and a field theoretical model, Phys. Rev. D 73, 014019 (2006), arXiv:hep-ph/0506234 .
- Roessner et al. [2007] S. Roessner, C. Ratti, and W. Weise, Polyakov loop, diquarks and the two-flavour phase diagram, Phys. Rev. D 75, 034007 (2007), arXiv:hep-ph/0609281 .
- Clevinger et al. [2022] A. Clevinger, J. Corkish, K. Aryal, and V. Dexheimer, Hybrid equations of state for neutron stars with hyperons and deltas, Eur. Phys. J. A 58, 96 (2022), arXiv:2205.00559 [astro-ph.HE] .
- Kumar et al. [2024c] R. Kumar, K. Aryal, A. Clevinger, and V. Dexheimer, Effects of hyperon potentials and symmetry energy in quark deconfinement, Phys. Lett. B 849, 138475 (2024c), arXiv:2311.15968 [nucl-th] .
- de Paoli et al. [2013] M. G. de Paoli, D. P. Menezes, L. B. Castro, and C. C. Barros, Jr, The Rarita-Schwinger Particles Under de Influence of Strong Magnetic Fields, J. Phys. G 40, 055007 (2013), arXiv:1207.4063 [math-ph] .
- Wang et al. [2001] P. Wang, Z.-Y. Zhang, and Y.-W. Yu, Nuclear matter in a chiral SU(3) quark mean field model, Commun. Theor. Phys. 36, 71 (2001).
- Fodor and Katz [2002] Z. Fodor and S. D. Katz, Lattice determination of the critical point of QCD at finite T and mu, JHEP 03, 014, arXiv:hep-lat/0106002 .
- Fodor and Katz [2004] Z. Fodor and S. D. Katz, Critical point of QCD at finite T and mu, lattice results for physical quark masses, JHEP 04, 050, arXiv:hep-lat/0402006 .
- Burkardt [2023] J. Burkardt, Fsolve: Nonlinear equation solver (2023), accessed on 9/29/23.
- More et al. [1980] J. More, B. Garbow, and K. Hillstrom, User Guide for MINPACK-1, Technical Report ANL-80-74 (Argonne National Laboratory, 1980).
- More et al. [1984] J. More, D. Sorenson, B. Garbow, and K. Hillstrom, The minpack project, in Sources and Development of Mathematical Software, edited by W. Cowell (Prentice-Hall, 1984).
- San Martin et al. [2023] J. S. San Martin, R. Hirayama, J. Hammelmann, J. M. Karthein, P. Parotto, J. Noronha-Hostler, C. Ratti, and H. Elfner, Thermodynamics of an updated hadronic resonance list and influence on hadronic transport (2023), arXiv:2309.01737 [nucl-th] .
- McLerran and Reddy [2019] L. McLerran and S. Reddy, Quarkyonic Matter and Neutron Stars, Phys. Rev. Lett. 122, 122701 (2019), arXiv:1811.12503 [nucl-th] .
- Stephanov [2009] M. A. Stephanov, Non-Gaussian fluctuations near the QCD critical point, Phys. Rev. Lett. 102, 032301 (2009), arXiv:0809.3450 [hep-ph] .
- Stephanov [2011] M. A. Stephanov, On the sign of kurtosis near the QCD critical point, Phys. Rev. Lett. 107, 052301 (2011), arXiv:1104.1627 [hep-ph] .
- Karsch [2012] F. Karsch, Determination of Freeze-out Conditions from Lattice QCD Calculations, Central Eur. J. Phys. 10, 1234 (2012), arXiv:1202.4173 [hep-lat] .
- Critelli et al. [2017] R. Critelli, J. Noronha, J. Noronha-Hostler, I. Portillo, C. Ratti, and R. Rougemont, Critical point in the phase diagram of primordial quark-gluon matter from black hole physics, Phys. Rev. D 96, 096026 (2017), arXiv:1706.00455 [nucl-th] .
- Borsanyi et al. [2018] S. Borsanyi, Z. Fodor, J. N. Guenther, S. K. Katz, K. K. Szabo, A. Pasztor, I. Portillo, and C. Ratti, Higher order fluctuations and correlations of conserved charges from lattice QCD, JHEP 10, 205, arXiv:1805.04445 [hep-lat] .
- Pradeep and Stephanov [2019] M. S. Pradeep and M. Stephanov, Universality of the critical point mapping between Ising model and QCD at small quark mass, Phys. Rev. D 100, 056003 (2019), arXiv:1905.13247 [hep-ph] .
- Bazavov et al. [2020] A. Bazavov et al., Skewness, kurtosis, and the fifth and sixth order cumulants of net baryon-number distributions from lattice QCD confront high-statistics STAR data, Phys. Rev. D 101, 074502 (2020), arXiv:2001.08530 [hep-lat] .
- Mroczek et al. [2021] D. Mroczek, A. R. Nava Acuna, J. Noronha-Hostler, P. Parotto, C. Ratti, and M. A. Stephanov, Quartic cumulant of baryon number in the presence of a QCD critical point, Phys. Rev. C 103, 034901 (2021), arXiv:2008.04022 [nucl-th] .
- Dore et al. [2022] T. Dore, J. M. Karthein, I. Long, D. Mroczek, J. Noronha-Hostler, P. Parotto, C. Ratti, and Y. Yamauchi, Critical lensing and kurtosis near a critical point in the QCD phase diagram in and out of equilibrium, Phys. Rev. D 106, 094024 (2022), arXiv:2207.04086 [nucl-th] .
- Borsanyi et al. [2024] S. Borsanyi, Z. Fodor, J. N. Guenther, S. D. Katz, P. Parotto, A. Pasztor, D. Pesznyak, K. K. Szabo, and C. H. Wong, Continuum-extrapolated high-order baryon fluctuations, Phys. Rev. D 110, L011501 (2024), arXiv:2312.07528 [hep-lat] .
- Vovchenko et al. [2015] V. Vovchenko, D. V. Anchishkin, and M. I. Gorenstein, Van der Waals Equation of State with Fermi Statistics for Nuclear Matter, Phys. Rev. C 91, 064314 (2015), arXiv:1504.01363 [nucl-th] .
- Gorda et al. [2021] T. Gorda, A. Kurkela, R. Paatelainen, S. Säppi, and A. Vuorinen, Cold quark matter at N3LO: Soft contributions, Phys. Rev. D 104, 074015 (2021), arXiv:2103.07427 [hep-ph] .
- Komoltsev and Kurkela [2022] O. Komoltsev and A. Kurkela, How Perturbative QCD Constrains the Equation of State at Neutron-Star Densities, Phys. Rev. Lett. 128, 202701 (2022), arXiv:2111.05350 [nucl-th] .
- Dexheimer et al. [2012] V. Dexheimer, R. Negreiros, and S. Schramm, Hybrid Stars in a Strong Magnetic Field, Eur. Phys. J. A 48, 189 (2012), arXiv:1108.4479 [astro-ph.HE] .
- Peterson et al. [2023] J. Peterson, P. Costa, R. Kumar, V. Dexheimer, R. Negreiros, and C. Providencia, Temperature and strong magnetic field effects in dense matter, Phys. Rev. D 108, 063011 (2023), arXiv:2304.02454 [nucl-th] .
- Negreiros et al. [2018] R. Negreiros, L. Tolos, M. Centelles, A. Ramos, and V. Dexheimer, Cooling of Small and Massive Hyperonic Stars, Astrophys. J. 863, 104 (2018), arXiv:1804.00334 [astro-ph.HE] .
- [137] Coming soon in october 2024.
- Hunter [2007] J. D. Hunter, Matplotlib: A 2d graphics environment, Computing in Science & Engineering 9, 90 (2007).
- Caswell et al. [2023] T. A. Caswell, A. Lee, E. S. de Andrade, M. Droettboom, T. Hoffmann, J. Klymak, J. Hunter, E. Firing, D. Stansby, N. Varoquaux, J. H. Nielsen, B. Root, R. May, O. Gustafsson, P. Elson, J. K. Seppänen, J.-J. Lee, D. Dale, hannah, D. McDougall, A. Straw, P. Hobson, K. Sunden, G. Lucas, C. Gohlke, A. F. Vincent, T. S. Yu, E. Ma, S. Silvester, and C. Moad, matplotlib/matplotlib: Rel: v3.7.1 (2023).
- Veltman et al. [2007] M. Veltman, B. de Wit, and G. ’t Hooft, Lie groups in physics (2007).
Appendix A Particle multiplets
The following baryon and meson matrices are constructed from the triplet tensor products and , respectively, where [140].
-
•
Baryon Matrix
(92) -
•
Scalar-Meson Matrix
(93) -
•
Vector-Meson Matrix
(94) -
•
Pseudoscalar-Meson Matrix
where .
-
•
Pseudoscalar-Meson Singlet Matrix
(96)
Appendix B Lagrangian density calculations
B.1 The self-interaction term for scalar mesons
In this study, we consider a term with the form
| (97) |
In the mean-field approximation, Eq. 93 becomes
| (98) |
and from Eq. 7,
| (99) |
For the vacuum expectation value of , , the isovector meson does not contribute
| (100) |
For the other terms in the Lagrangian, we use
| (101) |
such that from Eqs. 7, 8, and 9 becomes
| (102) |
The scalar Lagrangian is then given by
| (103) | ||||
| (104) |
B.2 The baryon-meson interaction term
having for a diagonal
| (106) |
| (107) |
| (108) |
with trace
| (109) |
On the other hand,
and again, explicitly for a diagonal ,
| (111) |
| (112) |
| (113) |
with trace
| (114) |
Also,
| (115) |
and
| (116) |
Combining those into Eqs. 20 and 21, we obtain
| (117) |
| (118) |
where the terms are not considered as they do not contribute to the mean-field approximation.
Case 1) For (the scalar-meson matrix in the mean-field approximation)
To illustrate how to find the coupling for the baryons, we take the particular case of the proton. Its couplings to the scalar mesons , and are found by replacing Eq. 116, Eq. 117, and Eq. 118 in Eq. 19, such that
| (119) |
From the mass term of Dirac Lagrangian for fermions, we find
| (120) |
which is a term appearing in the total effective mass of the proton in Eq. 36. By rearranging the terms for a particular scalar meson, we get
| (121) |
The neutron has the same coupling to the and mesons, but it couples to the with opposite sign, so we can write
| (122) | ||||
Indeed, this holds for all baryons in the octet: the and couplings are equal for the baryon families, while the coupling differentiates them due to isospin. Additionally, for the hyperons, the addition of the symmetry-breaking term to fix their potentials, Eq. 29, adds an additional contribution proportional to in the and couplings. For the hyperon:
| (123) |
For the ’s:
| (124) |
And for the ’s:
| (125) |
From this discussion, we can identify the effective masses, written explicitly in terms of the original couplings. In Eq. 120, the singlet term is identical for all baryons. The second term exists for both nucleons , and the third term differentiates the nucleons due to isospin. We can repeat this for all hyperons as well, the terms that are identical for the hyperon multiplets are
| (126) |
The full effective mass expressions, including the constant mass term , are
| (127) |
Case 2) For (the vector-meson matrix in the mean-field approximation)
| (128) |
| (129) |
For the vector mesons, an additional complication arises: the and fields are defined as a combination of the singlet and octet mesons ( and ) as
| (130) |
which adds a dependence on in the baryon-meson couplings. In principle, this is also the case for the and mesons, but for them, the mixing angle is , such that is purely singlet and is purely octet, and there is no angular contribution. For the vectors on the other hand, it is customary to take the ideal mixing angle, , with the condition (see Sec. II.2.5 for more). See [104] for a more complete discussion.
Appendix C Equations of motion and thermodynamics for a free Fermi gas
For a given fermion (or antifermion) of mass , which must obey Fermi-Dirac statistics, the distribution function as a function of energy level, temperature, and chemical potential reads
| (131) |
The Dirac Lagrangian density for spin 1/2 fermions is given by
| (132) |
from which applying the Euler-Lagrange equations
| (133) |
for each particle (or antiparticle) resulting in the equations of motion
| (134) |
These are both linear and first-order, indicating a plane-wave solution of the form
| (135) |
where is a four-vector spinor for Fermi momentum and spin . If satisfies the equations of motion, then the Lagrangian is zero. We can use this in the energy-momentum tensor together with the dimensional Minkowski metric to obtain the energy-momentum tensor
| (136) |
yielding
| (137) |
The energy density and pressure are then obtained in the ideal fluid approximation, where there is no dissipation and all non-diagonal terms vanish, by
| (138) |
giving
| (139) |
Additionally, we assume there is rotational symmetry, which is broken by the presence of magnetic fields, that would to different pressures in the directions longitudinal and perpendicular to the local direction of the field. Applying plane-wave and periodic boundary conditions it can be shown that
| (140) |
and
| (141) |
where for baryons and leptons and for quarks is the particle degeneracy. are particle energy levels. Likewise, the number density and entropy density can be calculated as
| (142) |
and
| (143) |
Additionally, in the presence of interactions, the scalar (number) density (or the source for scalar fields) is
| (144) |
Appendix D How to use the software
There are multiple ways to run the CMF solver software to calculate equations of state:
-
•
To run a calculation using a standalone script, download the source code package from the associated software publication [137], where you can also find instructions detailing how to compile and execute the code along with directions to the MUSES support community.
-
•
You may also run the CMF solver as a MUSES module in a processing workflow executed by the MUSES Calculation Engine. This method allows you to optionally include other MUSES modules in your workflows to perform more complex data processing. The Calculation Engine is also free and open-source software, available both for download and as an online service offered by the MUSES collaboration. Although a local installation requires Docker Compose, the use of containerization means you can run the software without installing the complex set of specific dependencies required by the CMF module. See the MUSES project website to learn more [80].
Appendix E Ensembles
Given a multi-variable function dependent on parameters F(a,b,c), we can ensure it possesses a minimum by showing that it has an extremum
| (145) |
and that it has positive concavity. The latter follows from the determinant of the Hessian matrix
| (146) |
and its submatrices being . In the case of, e.g., , it is implied that this means . The order of derivatives does not matter, as this matrix is symmetric, but permutations of the variables must be included, as there is no physical justification for any particular order.
E.1 Microcanonical ensemble
In this ensemble, based on the conservation of energy , the fixed variables are a number of particles , volume , and energy . Minimization of energy implies that the differential
| (147) |
and det with (, , and )
| (148) |
where we used that , , and .
In the zero-temperature limit, reduces to
| (149) |
and stability requires
| (150) | ||||
| (151) |
and
| (152) |
Using the Maxwell relation , we obtain
| (153) |
Using Eq. 151 and the definition of , one can also write
| (154) | ||||
| (155) |
E.2 Canonical ensemble
In this ensemble, based on the conservation of (Helmholtz) free energy , the fixed variables are a number of particles , volume , and energy . Minimization of free energy implies that the differential
| (156) |
and det with (, , and )
| (157) |
where we used that , , and .
In the zero-temperature limit, reduces once more to Eq. 149.
E.3 Grand canonical ensemble
In this ensemble, based on the conservation of grand potential , the fixed variables are chemical potentials , volume , and energy . Minimization of grand potential implies that the differential
| (158) |
and det with (, , and )
| (159) |
where we used that , , and .
In the zero-temperature limit, reduces to
| (160) |
and stability requires
| (161) | ||||
| (162) |
and
| (163) |
Using the Maxwell relation , we obtain
| (164) |
E.4 Infinite volume limit
For bulk matter, it is convenient to divide our grand potential with respect to the volume to get , the fixed variables now being only and . Maximization of the pressure implies that the differential
| (165) |
and det with (, , and )
| (166) |
where we used that , and .
In the zero-temperature limit, reduces to
| (167) |
and stability requires
| (168) |
We can also write for our particular case
| (169) | ||||
| (170) |
which using Eq. 168 means
| (171) |
This is the case for baryon number, , at . Note that this is similar to Eq. 155, one of the stability conditions in the microcanonical or canonical ensembles.
E.5 Multiple chemical potentials
Expanding Eq. 166 to the case of chemical potentials, , , and
| (172) |
| (173) |
Then, to find the stability constraints, one must ensure that the determinants of all submatrices are positive. Thus, one can show that at finite and infinite the constraints are:
| (174) | ||||
| (175) | ||||
| (176) | ||||
| (177) |
| (178) |
| (179) | ||||
| (180) | ||||
| (181) | ||||
| (182) |
| (183) | |||||
| (184) | |||||
| (185) |
Appendix F C1-C4 with baryon octet + quarks
Below we make a direct comparison between all four couplings C1-C4 for stable solutions only for the baryon octet+quarks at . Within panels a)-d) of Figure A21, we depict mean-field mesons against . As the previous discussions have covered in detail the C3 and C4 coupling schemes, we proceed to examine C1 and C2 in detail here. We note that the stable solutions of C++ and legacy Fortran solutions match for these couplings, except for the liquid-gas first-order phase transition, where Fortran presents no points. The and mean fields (panel a) and b)) follow a trend similar to C3 and C4. The C1 and C2 coupling schemes showcase a pronounced deconfinement first-order phase transition (panel e)) around MeV, and MeV, respectively.
Concerning the strange mean-field mesons, for C1 and C2, they decrease with increasing , particularly upon the emergence of hyperons, until the deconfinement phase transition. Above the deconfinement phase transition, the strange field decreases in value, while, the field drops to zero in the quark phase (due to its lack of coupling with quarks, see Table 10), exactly as C3 and C4.
For C1 and C2, (panel e)) behaves just like C3. For the EoS (panel f)), all couplings exhibit first-order deconfinement phase transition. In the hadronic phase, the wiggle observed in C1 and C2, just like C3, is due to the emergence of hyperons, indicating a higher-order phase transition as confirmed by the speed of sound plot in panel g). In the quark phase, all coupling schemes overlap because vector fields do not couple with quarks. Panel h) shows that the density of C1 and C2 aligns with C3.
In panels a) and b) of Figure A22, we observe and , with C1 and C2 resembling C3. Lastly, panels c) and d) depict susceptibilities for all coupling cases, with the discontinuities marking the first-order phase transitions (marked by discontinuities in ). The discontinuities in indicate that the onset of strangeness in the stable phases (both hadronic and quark) is a third-order phase transition.