The possibility of twin star solutions in a model based on lattice QCD thermodynamics

The properties of compact stars and in particular the existence of twin star solutions are investigated within an effective model that is constrained by lattice QCD thermodynamics. The model is modified at large baryon densities to incorporate a large variety of scenarios of first order phase transitions to a phase of deconfined quarks. This is achieved by matching two different variants of the bag model equation of state, in order to estimate the role of the Bag model parameters on the appearance of a second family of neutron stars. The produced sequences of neutron stars are compared with modern constrains on stellar masses, radii, and tidal deformability from astrophysical observations and gravitational wave analyses. It is found that most of the possible scenarios disfavor a strong phase transition to quark matter and do not support the conjecture of a second family of neutron stars.


I. INTRODUCTION
It is presumed that neutron stars (NS) can contain deconfined quark matter due to the high densities achieved in their interiors and might therefore play a decisive role along with both low and high energy nuclear physics in the exploration of the strong interaction and the Quantum Chromodynamics (QCD) phase diagram. Next to particle accelerators [1], neutron stars open an alternative window into the structure of the densest matter in our universe [2,3]. The extremes of QCD matter manifest in similar form in both the stellar phenomena of merging neutron stars and in the laboratory through relativistic heavy ion collisions [4][5][6][7], implying that similar densities and temperatures are excited in rather different physical systems. One particular possible feature of the QCD phase diagram is of great interest: the possible existence of a first order phase transition from a hadronic phase to a system of deconfined quarks. Such a phase transition can be discovered astrophysically in the properties of compact stars and their corresponding mass-radius relations from the TolmanâĂŞOppenheimerâĂŞVolkoff equation (TOV) [8,9], tidal deformability [10], as well as dynamical observables from the binary mergers of neutron stars as in the GW170817 event [3,11]. New methods, such as machine learning [12,13] and Bayesian analysis [14][15][16] are being developed in order to directly extract the equation of state (EoS) from available data. In context of the QCD phase diagram, astrophysical searches for twin stars, stellar objects with identical masses but different radii, are of particular interest. The different radii are assumed to be a result of distinct particle compositions in stellar interiors. If twin stars are formed, this would give a direct hint towards a sharp phase transition in QCD matter, e.g. a phase with hadronic degrees of freedom and a second one composed of deconfined quarks. The most prominent scenario is a realization of three families of compact stars: white dwarfs, neutrons stars and their stable twins [17][18][19][20][21][22][23][24]. Recently the concept of a "delayed phase transition" has been proposed, in which a metastable hypermassive star, developed some time after the merger event, exhibits a quark core. Two distinct post-merger gravitational-wave frequencies, before and after the phase transition could as well be a promising signature for the existence of quark matter [25,26]. The EoS for compact stars is often calculated on a basis of nuclear interaction models [27,28], including additional hyperonic degrees of freedom [29][30][31][32][33] and models based on quark interactions [34][35][36], to name a few. In this work we will take a slightly different approach, by employing a quark hybrid EoS in which the model parameters are not only fixed to known nuclear matter properties, but also describe the smooth deconfinement transition and thermodynamics at high temperatures and vanishing densities, obtained from state-of-the-art lattice QCD calculations [37,38]. To account for a possibility of twin star solutions we modify the model, which usually only contains a crossover transition to a phase with hadrons and quarks, at high densities by constructing a transition to a deconfined phase of quark matter. This new construction allows to study various possible scenarios involving a first order phase transition to deconfined quarks.
We use the SU(3)-flavor parity-doublet Polyakov-loop quarkhadron mean-field model (CMF) to describe a hadronic and quark system in which interactions are driven by mean field meson exchange and repulsive excluded volume interactions. This model is in agreement with the mass-radius and tidal deformability constraints of astrophysical observations and has been used to investigate the properties of compact stars [38]. Because the high density region of the EoS cannot be constraint by methods of lattice QCD [39][40][41], we will modify the model to investigate various possibilities of 1st order phase transitions at large density. Perturbative QCD calculations in this regime suggest that the pressure of QCD matter is below the Stefan Boltzmann limit of massless non-interacting gas of three quark flavors [42].
For the high density regime of the quark phase, we test two different models: The MIT Bag model where we vary the Bag parameter B and the Bag model enhanced by vector meson repulsion (vMIT), in which we study different coupling constants g ω q between quarks and the ω vector field. The g ω q defines the strength of repulsive interactions of quarks by a vector ω meson exchange. Even though such a new EoS with explicit quark vector repulsion may violate lattice QCD constraints for QCD matter at vanishing density [43]. Density dependent couplings to a repulsive field may cure this problem at finite density region of the QCD phase diagram which is relevant for astrophysics. This paper is structured as follows: In Sections II and III, we present the CMF and Bag models used in this framework, as well as the Maxwell construction formalism to implement the transition. We present our results and discuss them in context of feasibility of twin star solutions in Section IV. We summarize our results in section V.

II. CHIRAL MEAN FIELD MODEL
The Chiral SU(3)-flavor parity-doublet Polyakov-loop quark-hadron mean-field model, CMF, describes matter composed of hadrons and quarks. It incorporates several concepts of QCD phenomenology, meson exchange interactions in the baryon octet [44], excluded volume repulsive interactions amongst all hadrons [45,46], parity doubling1 amongst baryons [48] and quarks within a Polyakov loop extended Nambu Jona-Lasinio model [49]. The parity doubling assumes that the mass splitting of the baryon masses and their parity partners is generated by scalar mesonic fields, formulated within a mean field approach. As the energy density, (and therefore the scalar density), increases, the mass gap between baryons and their parity partners decreases until degeneracy between the states occurs. The CMF model includes the full PDG list of hadrons [50] which are attributed excluded volume parameters to mimic hadron finite size and their repulsive interactions. The coupling constants of the hadronic sector are chosen such that properties of nuclear matter are  s (shaded area, 2σ credible interval from the deep neural network) and CMF model lie within a similar range. A first order phase transition at 4ρ 0 is reflected in the squared speed of sound by the sudden drop to zero. At this density chiral symmetry is restored. The dashed grey line is the Stefan Boltzmann limit for an ultrarelativistic ideal fluid with squared sound speed v 2 s = 1/3. The CMF model converges to this value at large densities ≥ 30ρ B (not shown here). The black cross at ∼ 8ρ/ρ 0 is the maximum sound speed of the CMF EoS. This regime appears due to the strong repulsion among hadrons.
The speed of sound from the CMF model at T = 0 is shown in Fig. 1 and compared to estimates of a deep neural network analysis, which is based on training data from available massradius observations [12]. The CMF model neutron star EoS follows the trend of the 2σ confidence interval from the neural network at low-moderate densities. However, the uncertainties are still significant and the results of the neural network do not quantitatively constrain the high density regime. Note that this model does not naturally exhibit a first order phase transition from nuclear to quark matter, but a 1st order phase transition due to chiral symmetry restoration amongst baryon parity partners [51]. The abrupt decrease of the speed of sound to the zero value locates the phase transition of the CMF EoS. However this chiral phase transition is too narrow to be reflected in the structure of neutron stars [38] and can be easily hidden in the errorbands of the neutral network analysis. According to the results, the CMF model predicts hybrid neutron stars with masses up to 2 M . Such stars have a total quark mass fraction up to 30%. There is no sharp transition in the CMF EoS from hadrons to quarks, so that no second family of quark stars exists within this framework. More details and discussion on this model can be found in [38]. The CMF model includes scalar-isovector ρ-mesons which control isospin asymmetry and are thus relevant for NS matter, where the amount of neutrons is much larger then the amount of protons [29,37,44,52]. The baryon octet couples to the ω, ρ and (hidden strange) φ field [44]. The baryon masses are dynamically generated by their couplings to the scalar σ and strange ζ field. These two fields are order parameters for the chiral transition and directly affect the effective baryon masses, see Eq. 1 and [37]. With increasing baryon density ρ, the σfield decreases and causes the effective masses of the particles to restore chiral symmetry. The effective masses read where m 0 is an explicit mass term of the baryon octet m 0 = 759 MeV, n s is the number of strange quarks in baryons and m s = 130 MeV is the mass of the strange quark. The signs ± indicate the parity quantum number of the particle. Finally, g 1 σi , g 1 ζi , g 2 σi and g 2 ζi are coupling constants to scalar σ and ζ fields for i-th baryon of the octet. At high densities, quarks are expected to be dominant so a deconfinement mechanism should be incorporated in the model. This is done in analogy to the Polyakov-loop-extended Nambu Jona-Lasinio (PNJL) model [49] which is an effective chiral field model for describing quark matter. The Polyakov-loop Φ which effectively represents gluon degrees of freedom is controlled by the temperature dependent potential U(Φ) which is zero for the case of cold neutron star matter [38]. The quark masses m * i are dynamically generated and controlled by the σand ζ-field. The effective masses for up, down and strange quarks read The σ-meson controls masses for up and down quarks and the ζ-meson generates the strange quark mass. The light u and d quarks have the explicit ground state mass term δm u = δm d = 5 MeV and the heavier strange quark has a mass δ s = 150 MeV and m 0q = 235 MeV. An additional mass m 0q is introduced to take into account quarks sizable thermal masses which usually appear in EoS models for the quark gluon plasma [53][54][55][56], this term also prevents quark appearance in nuclear matter. An explicit volume term v B = v is added to the hadrons to suppress them in the quark phase [37]. Consequently, as soon as quarks contribute to the pressure P, they suppress hadrons by lowering their chemical potential. For neutron star matter, leptons are taken into account in order to obey charge neutrality and β equilibrium.

III. THE HIGH DENSITY TRANSITION
In order to allow a possible phase transition to a fully deconfined system of quarks, the CMF model is matched to different realizations of the Bag Model. In the following, we investigate a transition from the CMF model to two versions of the Bag model: the standard MIT Bag Model for the Quark Gluon plasma, with massless and non interacting quarks [35,57] and the vector MIT Bag model [58].

A. Bag Model
The first formulation of the Bag model came from Bogoliubov in 1968 who built a theory where three massless quarks inside a spherical volume with radius R are bound in an infinite square well potential [57]2. Then in 1974, an enhanced Bag model was rediscovered independently [35] which then was named the MIT Bag model. In this formalism, quarks are confined in the bag and the surface of the bag has boundary conditions such that the quark current perpendicular to the surface is zero. Vanishing quark current at the boundary of the bag effectively enforces quark confinement. Several modifications were proposed to the original MIT Bag model, such as the inclusion of medium effects [59], density [60] and temperature dependencies [61], and recently it was modified to take into account Hagedorn mass spectra of hadrons [56]. Among these modified versions, the ones which include vector interactions among quarks are particularly prominent in the description of massive neutron stars [58]. We relate the baryon and quark chemical potentials by a factor 1/3, using µ B = u,d,s µ q i ). For such a configuration, the pressure and energy density of the basic Bag Model at zero temperature read: For two flavors, the degeneracy factor ν u,d = 2 × 3 × 2 = 12 respectively for two spin states, 3 colours and 2 flavors. The three flavor version has a degeneracy factor of ν u,d,s = 2 × 3 × 3. The additional degree of freedom, the strange quark, in the three flavor Bag model may allow (depending on the bag constant) for a bound strange quark matter state [62] where up to 1/3 of matter is composed of strange quarks which are absent in ordinary matter. Resulting stable exotic nuclear states where conjectured in [63]. The creation of such a type of matter in relativistic nuclear collisions in the laboratory was proposed as signature of quark-gluon plasma formation [64][65][66], thus never confirmed experimentally.
In the bag model with vector interactions, quarks have nonzero masses and repulsive interaction is taken into account, represented by a coupling g V of the vector-isoscalar meson V to the quarks. The free leptonic e − and µ − degrees of freedom are included as well. The modified quark chemical potential of quarks at T = 0 reads with the Fermi momentum vector k F and the bare quark masses m q . The ω-field suppresses hadronic abundances.
The effective chemical potential m * q is reduced by the vector interactions. The phenomenological vMIT Bag model includes chiral symmetry breaking and repulsive vector repulsion [58] in its Lagrangian. It has been used in the literature to fulfill two solar mass constraint for neutron star masses from observational astrophysics [58,67,68]. There are precise Shapiro time delay measurements that observed high NS masses like pulsar PSR J0740+6620 (2.17 +0.11 −0.10 M ) [69], PSR J0348+0432 (2.01 ± 0.04 M ) [70], and, recently PSR J0348+0432 (2.27 +0. 17 −0.15 M ) [71]. The vMIT Bag model is an attempt to describe a stiff quark EoS that supports such high masses. It has an additional term that describes repulsive vector interactions coming from a non-vanishing mean field in the vector meson interaction channel. However, in the other regime of QCD at high temperatures and vanishing densities, analysis of lattice QCD data disfavors repulsion among quarks [43,72].

B. Constructing the combined model
We construct a first order phase transition at high baryonic densities to quark matter. The phase transition in our approach is modeled by a Maxwell construction which is well adopted for such a scenario [73,74]3. The transition occurs at a point where the two pressures of both EoSs intersect as functions of chemical potential. At the phase coexistence P cr low = i P cr q i and µ cr B = i µ cr q i , where P q i is the pressure contribution of the quarks and P low is the pressure of lower density phase, herethe CMF EoS. At the intersection point of both EoSs, baryon number density n B jumps as well as the energy density . In the following, we refer to the jump in energy density ∆ as latent heat, i.g. the discontinuity in energy density at the first order phase transition from the CMF to the Bag model. Visually, this can be seen in Fig. 2, where the combined EoS corresponds to the maximum at a given baryon chemical potential µ B of the dashed black curve (CMF model) and the colored curves (Bag model). The different colors in Fig. 2 correspond to different Bag constants in the model, see color code, which we will discuss more in detail in section IV. We will vary values of the Bag constant B and additionally, in the case of the vMIT Bag model, the coupling g ω q to study all possible scenarios of the phase transition. The instability of a star is proportional to the value of p trans and inverse proportional to the gap ∆ trans in energy density. Stable twin star branches in our model can only occur if the following condition is fulfilled This condition is called Seidov limit [78], it is a generic condition for stellar equilibrium of a star with a phase change. It Note, in the present setup there are no Bag matter stars so this relation can not be straightly employed in our calculations since the lower density matter is described with the CMF model.

A. Transition to MIT Bag model
In the following we present different mass-radius relations, their corresponding EoS, and the dimensionless tidal deformability Λ for our combined CMF-MIT Bag model. Two different scenarios are tested: a combination of the CMF model with a two flavor MIT bag model and the combination of the CMF with a three flavor MIT Bag model, which only differ in the amount of quarks in the Bag model regime. The Bag constant influences the onset of the phase transition and the number of flavors changes the latent heat of the Maxwell construction. Both of these quantities influence possible twin star solutions. Smaller values for the Bag parameter B 1/4 than 145 MeV are excluded. Otherwise two-flavor quark matter would have a lower energy than 56 Fe and will form a ground state for ordinary matter different from the one we observe [80]. Starting with the 2-flavor Bag-CMF model model, the pressure in units of µ 4 B as function of chemical potential µ B is presented in Fig. 3 (left). The combined EoS are constructed such that they, for all B, converge to the Stefan Boltzmann limit, see Eq. 5. The SB limit for the 3-flavor Bag-CMF model is shifted upwards, see Fig. 3 (right). The increase of the number of flavors from two to three increases the Stefan Boltzmann limit by 50%. We define the lower limit for the Bag parameter B by requiring that the CMF and Bag models EoS still intersect. For three flavors, the Bag parameter is shifted upwards to B 1/4 ≈ 170 MeV. We observe stable TOV solutions with a significant contribution of Bag matter for Bag constants B 1/4 175 MeV in Fig. 4 (left). These solutions correspond to the 2-flavor Bag-CMF EoS in Fig. 3 (left). In Fig. 4  [MeV] FIG. 5. Dimensionless tidal deformability parameter Λ for the 2-flavor (left) and 3-flavor (right) Bag-CMF model. The blue shaded area is the prediction of a neutral network trained on observed neutron star masses and radii in [12]. The vertical black line is a 2σ confidence bound 375 < Λ 1.4M < 800 from Bayesian analysis based on the merger event GW170817 for a 1.4M NS [11]. The color indicates the value of the Bag constant, see colorbar on the right-hand side. The 3-flavor Bag-CMF model predicts significantly lower values of Λ for masses ≥ 1.2 M NS which fall below the constraints on Λ from the merger event, as well as below the neural network analysis of NS data [12].
in Fig. 1.The tidal deformability for the two flavor Bag-CMF EoS in Fig. 5 (left) lies within the blue area constraint as well as in the GW170817 merger constraint for Λ assuming a mass ratio of 1 with M 1 = M 2 = 1.4 M [3,11,12]. For three flavors in Fig. 5 (right) the values for Λ lie below both 2σ confidence intervals, assuming neutron star masses ≥ 1.2M . NSs with lower tidal deformability are more compact. The investigation above demonstrates that the latent heat ∆ of the 2-flavor Bag-CMF phase transition is not sufficient to obtain twin star solutions. For the 2-flavor Bag-CMF model, the combination of latent heat ∆ , transition pressure p trans , and transition energy density trans do not fulfill the condition in Eq. 8. Thus, either transition pressure/energy density need to be larger/smaller, and/or the latent heat is of considerable extent. Stars with Bag constants B ≥ 180 MeV become unstable, see Fig. 4 (left). A reason for that is the upper mass limit for pure Bag model stars, see Eq. 9. On the other hand, for a 3-flavor Bag-CMF model twin star solutions occur for lower values of B with masses below 2M . A reason for these low twin masses can arise from the nature of the soft MIT Bag EoS with a constant speed of sound v s = 1/3 for any MIT Bag parametrization. The first order phase transition has a larger latent heat, compared to our 2-flavor Bag-CMF model. The reason for that is the larger deviation in slopes, that both models each have at the Maxwell construction. This discontinuity leads to a jump in the baryon number density.

B. vMIT Bag model
We use a 3-flavor vector enhanced Bag model to construct a stiffer EoS for the quark phase in order to obtain sufficient high masses of 2M and consider repulsive vector interactions in the quark phase. This stiffening could possibly lead to higher mass twin star solutions. The masses of up, down and strange quarks are 1 MeV for up and down and 100 MeV for the strange quark. We vary the coupling constant g ω q in the range g ω q /m ω =[0, 1.75] fm where we use a ω-meson mass m ω = 728 MeV, see the color code in Fig. 6 and Fig. 7 [81]. Fig. 7 and Fig. 6 show M-R relations and corresponding EoS of 3 different g ω q -parameterizations where B 1/4 ∈ {166, 171, 180} MeV. The different color lines with corresponding color code show different quark coupling parameters g ω q /m ω in each plot. g ω q is the coupling parameter to the ω field and controls the strength of repulsive force amongst quarks. Increasing the vector repulsion decreases the pressure at a fixed chemical potential µ B , hence, an increase of g ω q will lead to a stiffer EoS [2]. The black star in each EoS figure is the maximum mass as obtained from the CMF model. A higher Bag constant leads to a later phase transition, this is the same behaviour as we observed for the MIT Bag model. For B 1/4 = 180 MeV we find twins with masses below 1.5 M . The motivation to choose a Bag model with repulsion was to increase the maximum masses of twins, however it seems that a stiff enough quark EoS alone is not sufficient for stable stars. In analogy to the reasoning for the flavor variation of the combined CMF MIT Bag model, the quantitative change in slope of CMF and vMIT EoS is a direct measure for the latent heat ∆ . A higher repulsion eventually leads to a too small jump in the energy density such that the Seidov limit in Eq. 8, which defines if a star with a specific central pressure can be stable, is not fulfilled and the star is not destabilized by the Maxwell construction. We summarize the interplay of the coupling strength a 0 and the Bag constant B in the framework of the vMIT Bag model as follows: 1. The repulsive coupling g ω q /m ω influences the onset of the transition. A smaller coupling constant leads to a smaller discontinuity in the baryon number density n B and thus the latent heat ∆ . Bag constant decreases the latent heat and the onset of the phase transition is shifted towards lower chemical potential. Possible twins only occur if the Bag constant is above ≈ 180 MeV, this softens the vMIT EoS because it shifts the curve parallel along the x-axis to 8  higher chemical potentials. The star is then immediately unstable after the transition.
Following these points, problems arise when the latent heat ∆ is too small. This is the case if the intersection of both EoS lie nearly parallel. A stiff EoS in the deconfined phase could help to increase the latent heat. One way for that to happen could be a larger repulsive coupling constant g ω q so that the quark EoS is stiffer than the hadronic EoS.
However, within this framework this is not feasible since g ω q has an upper limit that arises by requirement of EoS curves to intersect. If the quark EoS is too stiff, then deconfined quark matter always has a higher pressure then nuclear matter, that is not the case in nature since at lower densities hadrons dominate. Having a stiff hadronic EoS at the transition, followed by a soft quark EoS that stiffens quickly after the transition could possibly lead to twin star solutions. This could be formulated within a density dependent repulsive quark coupling framework, as it has been investigated in [82]. Instead of stiffening the quark EoS one could consider to soften the hadronic EoS at intermediate densities.
A new analysis of the NICER data gives hint that an extremely soft nuclear EoS and a strong phase transition are mutually exclusive [83]. A softening of the hadronic phase is possible through the appearance of additional baryonic degrees of freedom. In the CMF model, all hadronic species are included, but at T = 0 only nucleons and their parity partners appear while other hadronic species are suppressed by their EV interactions. However, additional softening could result from the appearance of ∆-Baryons or hyperons in the NS EoS due to a decrease of their repulsion or increase of their attractive interactions. The analysis of both lattice QCD data and heavy ion collisions indeed suggest that strange hadrons are subject to smaller EV repulsion due to their smaller size [84,85]. We leave the investigation of these systematics for future studies.

V. SUMMARY
The viability of twin star solutions due to a sharp phase transition to deconfined quark matter was studied. The transition was implemented by a Maxwell construction between the CMF model, a realistic hadron-quark EoS model for QCD matter which supports lattice QCD data and NS constraints, and the Bag model, EoS model for deconfined quark matter, at lower and higher densities respectively. To investigate different scenarios of the transition, the original MIT Bag model and vector-enhanced Bag model were considered. In the study the parameters of the CMF model remain fixed, but parameters of the Bag model were varied, namely, Bag constant, number of quark flavors and strength of the quark vector repulsion. The variation of the parameters leads the phase transition to occur at different densities and with different latent heat. In the present framework the mass-radius relations and tidal deformabilities were analysed. The mass-radius relations suggest that 2-flavor Bag-CMF model is stiff enough to produce stable configurations with significant fraction of deconfined quark matter for values of bag constant B 1/4 175 MeV. However, no twin solutions appear for 2-flavor case. For the 3-flavor Bag-CMF model, stable configurations with large quark content appear for bag values B 1/4 175 MeV, when the latent heat of the phase transition is large enough to destabilize the M-R branch. However these solutions only support NS masses up to 1.5M . The 3-flavor Bag-CMF model was further investigated by an inclusion of the vector repulsion amongst quarks where repulsion coupling was varied as well. The stiffening of the EoS allowed to produce twin star solutions with masses M ∼ 1.3M , however higher masses for twin solutions are not supported by this EoS. The analysis of tidal deformabilities Λ and comparison with available constraints also disfavors the suggested scenario of a sharp phase transition to quark matter. The transition to quark matter in NS must respect the absence of the vector repulsion among quarks at low density and high temperature regime of QCD [43,72]. The results suggest that stable high mass twin-star solutions may appear with a density dependent repulsive interaction scheme which incorporates a soft behaviour of quark matter at the density of the Maxwell construction, followed by a stiff quark phase at higher densities. These two characteristics seem necessary in order to obtain a sufficient latent heat at the 1st order phase transition and a second stable branch respectively. It may be worthwhile to investigate whether such behaviour can be brought in agreement with measurements of susceptibilities from lattice QCD simulations which are sensitive probes of density dependent interactions.