aa r X i v : . [ nu c l - t h ] F e b Hadron-quark mixed phase in the quark-meson coupling model
Min Ju, Xuhao Wu, Fan Ji, Jinniu Hu, ∗ and Hong Shen † School of Physics, Nankai University, Tianjin 300071, China School of Physics, Peking University, Beijing 100871, China
We explore the possibility of a structured hadron-quark mixed phase forming in the interior ofneutron stars. The quark-meson coupling (QMC) model, which explicitly incorporates the internalquark structure of the nucleon, is employed to describe the hadronic phase, while the quark phase isdescribed by the same bag model as the one used in the QMC framework, so as to keep consistencybetween the two coexisting phases. We analyze the effect of the appearance of hadron-quark pastaphases on the neutron-star properties. We also discuss the influence of nuclear symmetry energyand the bag constant B in quark matter on the deconfinement phase transition. For the treatmentof the hadron-quark mixed phase, we use the energy minimization method and compare it with theGibbs construction. The finite-size effects like surface and Coulomb energies are taken into accountin the energy minimization method; they play crucial roles in determining the pasta configurationduring the hadron-quark phase transition. It is found that the finite-size effects can significantlyreduce the region of the mixed phase relative to that of the Gibbs construction. Using a consistentvalue of B in the QMC model and quark matter, we find that hadron-quark pasta phases are formedin the interior of massive stars, but no pure quark matter can exist. I. INTRODUCTION
It is generally believed that hadronic matter undergoesa phase transition to deconfined quark matter at highbaryon densities, which may occur in the interior of mas-sive neutron stars [1–3]. Recent advances in astrophysicalobservations have led to increasing interest in exploringvarious properties of neutron stars, such as their inter-nal structure, mass-radius relation, and tidal deformabil-ity. Several precise mass measurements of massive pul-sars, PSR J1614-2230 [4–6], PSR J0348+0432 [7], andPSR J0740+6620 [8], constrain the maximum neutron-star mass to be larger than about 2 M ⊙ , which chal-lenges our understanding of the equation of state (EOS)of superdense matter. The first detection of gravita-tional waves from a binary neutron-star merger, knownas GW170817, provides valuable constraints on the tidaldeformability [9–11], which also restricts the radii of neu-tron stars [12–16]. More recently, the gravitational-waveevent GW190425 was reported by LIGO and Virgo Col-laborations [17] with the total mass of the binary sys-tem as large as 3 . +0 . − . M ⊙ , which may offer importantinformation for the EOS at high densities. The latestgravitational-wave event GW190814 was detected fromthe merger of binary coalescence involving a 23 . +1 . − . M ⊙ black hole with a 2 . +0 . − . M ⊙ compact object [18], wherethe secondary component could be interpreted as eitherthe lightest black hole or the heaviest neutron star everobserved. Furthermore, the new observations by the Neu-tron Star Interior Composition Explorer ( NICER ) pro-vided a simultaneous measurement of the mass and ra- ∗ Electronic address: [email protected] † Electronic address: [email protected] dius for PSR J0030+0451 [19, 20]. It is interesting tonotice that the constraints on the neutron-star radiusfrom various observations are consistent with each other,and suggest relatively small radii of neutron stars. Allof these exciting developments in astrophysical observa-tions provide a wealth of information about neutron starinteriors, where the deconfinement phase transition is ex-pected to take place.Theoretically, a hadron-quark mixed phase is predictedto exist in the region between hadronic matter and quarkmatter based on various approaches [1, 21–31]. However,large uncertainties in the structure and density range ofthe mixed phase are present due to the different modelsand methods used. The Gibbs construction [21] is gener-ally adopted for the description of hadron-quark mixedphase, where the coexisting hadronic and quark phasesare allowed to be charged separately, but the finite-sizeeffects like surface and Coulomb contributions are ne-glected. When the surface tension of the hadron-quarkinterface is sufficiently large, the mixed phase shouldbe described by the Maxwell construction, where localcharge neutrality is imposed and coexisting hadronic andquark phases have equal pressures and baryon chemicalpotentials but different electron chemical potentials. Itis evident that Gibbs and Maxwell constructions corre-spond to two limits of vanishing and large values of thesurface tension, respectively, while both of them involveonly bulk contributions without finite-size effects [32, 33].For a moderate surface tension, the hadron-quark mixedphase with some geometric structures, known as pastaphases, is expected to appear as a consequence of thecompetition between the surface and Coulomb ener-gies [33–38]. A realistic description of the structuredmixed phase was proposed in Refs. [35, 39], where aconsistent treatment of the electric field in the Wigner–Seitz cell would lead to inhomogeneous distributions ofcharged particles in both hadronic and quark phases. Ingeneral, for simplicity, the particle densities in the twocoexisting phases are assumed to be spatially constant,and the charge screening effect is neglected. A simplecoexisting phases (CP) method [33] is used for the de-scription of the hadron-quark pasta phases, where thecoexisting phases satisfy the Gibbs conditions for phaseequilibrium, and the finite-size effects like the surfaceand Coulomb energies are perturbatively taken into ac-count. An improved energy minimization (EM) methodincorporates the finite-size effects in a more consistentmanner [33], where the equilibrium conditions for coex-isting phases are derived by minimization of the totalenergy including surface and Coulomb contributions. Asa result, the equilibrium conditions obtained in the EMmethod are significantly different from the Gibbs con-ditions. The EM method, known as the compressibleliquid-drop (CLD) model, has been widely used in thestudy of nuclear liquid-gas phase transition at subnu-clear densities [40–43]. In the present work, we intendto use the EM method for studying the hadron-quarkpasta phases in the interior of massive neutron stars.To describe the hadron-quark pasta phases, we employthe Wigner–Seitz approximation, in which the systemis divided into equivalent and charge-neutral cells witha given geometric symmetry. The hadronic and quarkphases inside the cell are assumed to have constant den-sities and are separated by a sharp interface. For the de-scription of hadronic matter, we employ the quark-mesoncoupling (QMC) model, which explicitly incorporates thequark degrees of freedom. The QMC model, initially pro-posed by Guichon [44], has been extensively developedand applied to various nuclear and hadronic phenomenain the past decades [45–51]. In the QMC model, the nu-cleons in the nuclear medium are described by static MITbags that interact through the self-consistent exchange ofscalar and vector mesons in the mean-field approxima-tion. The mesons couple directly to the confined quarksinside the bags, not to pointlike nucleons as in the rela-tivistic mean-field (RMF) model. In contrast to the RMFapproach, the internal structure of the nucleon is explic-itly included in the QMC model, so it could be used forinvestigating the medium modification of nucleon struc-ture [49]. At sufficiently high densities, the nucleons andother hadrons are expected to dissolve into quarks. Thequark matter is then treated as a degenerate Fermi gasof u , d , and s quarks within the MIT bag model. In thissense, both hadronic and quark phases are consistentlydescribed within the same bag model. Therefore, theQMC model that incorporates the internal quark struc-ture of the nucleon is a promising approach for exploringthe hadron-quark phase transition in neutron stars. InRefs. [52–55], Panda et al. extended the QMC model tostudy the properties of compact stars, including the pos-sibility of hyperon formation and quark deconfinement.They employed the Gibbs construction to describe thehadron-quark mixed phase without considering finite-sizeeffects and possible geometric structures. In the present work, we aim to study the hadron-quark pasta phasesusing the QMC model, where the finite-size effects areconsidered within the EM method, and the properties ofa structured mixed phase are investigated.This article is organized as follows. In Sec. II, wedescribe the QMC model for the hadronic matter. InSec. III, the MIT bag model for the quark matter isbriefly reviewed. The EM method for the descriptionof hadron-quark pasta phases is presented in Sec. IV,where the surface tension at the hadron-quark interfaceis discussed. Section V contains numerical results of thehadron-quark pasta phases and corresponding neutron-star properties. Section VI is devoted to the conclusions. II. HADRONIC PHASE
The hadronic matter is described by the QMC model,where the nucleon is treated as a static MIT bag contain-ing three confined quarks. The nucleon-nucleon interac-tion is realized by the exchange of isoscalar-scalar meson σ , isoscalar-vector meson ω , and isovector-vector meson ρ in the mean-field approximation. In the QMC model,the effective meson fields couple directly to the confinedquarks inside the nucleon, not to the pointlike nucleonas in the RMF model. Therefore, the internal structureof the nucleon can be self-consistently determined, andis influenced by the meson fields in nuclear matter. Thequark field inside the bag satisfies the Dirac equation (cid:2) iγ µ ∂ µ − ( m q + g qσ σ ) − γ (cid:0) g qω ω + g qρ τ ρ (cid:1)(cid:3) ψ q = 0 , (1)where m q is the current quark mass and g qσ , g qω , and g qρ denote the quark-meson coupling constants. The nor-malized ground state for a quark in the bag is given by ψ q ( r , t ) = N q e − iǫ q t/R (cid:18) j ( x q r/R ) iβ q σ · ˆ r j ( x q r/R ) (cid:19) χ q √ π , (2)where ǫ q = Ω q + R (cid:0) g qω ω + g qρ τ ρ (cid:1) , (3) β q = s Ω q − Rm ∗ q Ω q + Rm ∗ q , (4) N − q = 2 R j ( x q )[Ω q (Ω q −
1) + Rm ∗ q / /x q , (5)with Ω q = q x q + ( Rm ∗ q ) and m ∗ q = m q + g qσ σ . Here R is the bag radius and χ q is the quark spinor. The valueof x q is determined by the boundary condition at the bagsurface j ( x q ) = β q j ( x q ) . (6)The energy of a nucleon bag consisting of three groundstate quarks is then given by E bag = 3 Ω q R − ZR + 43 πR B, (7)where the model parameter Z accounts for the zero-pointmotion and center-of-mass corrections, and B denotes thebag constant. The effective nucleon mass in the mediumis taken to be M ∗ N = E bag . The bag radius R is de-termined by the equilibrium condition ∂M ∗ N /∂R = 0.In the present work, we use the current quark mass m q = 5 . u and d quarks. The model param-eters, B / = 210 .
854 MeV and Z = 4 . M N = 939 MeVand the bag radius R = 0 . p and n ) and leptons ( e and µ ), we start with theeffective Lagrangian density involving the internal struc-ture of the nucleon together with the meson fields in themean-field approximation, L QMC = X b = n,p ¯ ψ b (cid:2) iγ µ ∂ µ − M ∗ N − γ ( g ω ω + g ρ τ ρ ) (cid:3) ψ b − m σ σ + 12 m ω ω + 12 m ρ ρ + Λ v g ω g ρ ω ρ + X l = e,µ ¯ ψ l [ iγ µ ∂ µ − m l ] ψ l , (8)where σ = h σ i , ω = (cid:10) ω (cid:11) , and ρ = (cid:10) ρ (cid:11) are the nonvan-ishing expectation values of meson fields in homogeneousnuclear matter. The effective nucleon mass M ∗ N is calcu-lated from Eq. (7) at the quark level, and it depends onthe σ field. Generally, the effective nucleon mass can beexpanded in terms of σ , which is written in the practicalform: M ∗ N = M N + a ( g qσ σ ) + b ( g qσ σ ) + c ( g qσ σ ) , (9)where the parameters a = 1 . b = 7 . × − MeV − , and c = 1 . × − MeV − are de-termined by fitting to the results of M ∗ N from Eq. (7).In the QMC model, the fundamental couplings are thoseto the quarks, which are related to the correspondingnucleon-meson couplings as g ω = 3 g qω and g ρ = g qρ . Inpractical calculations, we use the quark-meson couplingconstants g qσ = 5 . g qω = 3 . g qρ = 5 . ω - ρ coupling Λ v = 0 . E/A = − . n = 0 .
15 fm − with the symmetryenergy E sym = 31 MeV and its slope L = 40 MeV. Itis well known that the ω - ρ coupling term plays a cru-cial role in determining the density dependence of thesymmetry energy [40]. We include this coupling term inEq. (8), so that a small value of the slope L could beobtained according to current constraints. The mesonmasses are taken to be m σ = 550 MeV, m ω = 783 MeV,and m ρ = 770 MeV.In uniform hadronic matter, the coupled equations ofmotion for meson fields can be easily solved. The total energy density of hadron matter is calculated by ε HP = X b = p,n π Z k bF q k + M ∗ N k dk + 12 m σ σ + 12 m ω ω + 12 m ρ ρ + 3Λ v g ω g ρ ω ρ + X l = e,µ π Z k lF q k + m l k dk, (10)and the pressure is given by P HP = X b = p,n π Z k bF k dk q k + M ∗ N − m σ σ + 12 m ω ω + 12 m ρ ρ + Λ v g ω g ρ ω ρ + X l = e,µ π Z k lF k dk p k + m l . (11)For hadronic matter consisting of nucleons and leptonsin β equilibrium, the chemical potentials satisfy the re-lations µ p = µ n − µ e , (12) µ µ = µ e . (13)The charge neutrality condition is expressed as n HP c = n p − n e − n µ = 0 . (14)It is known that muons may appear when the electronchemical potential exceeds the muon mass, which occursat a density of 0 .
11 fm − in the QMC model. III. QUARK PHASE
To describe the quark matter at high densities, we em-ploy the same bag model as the one used to describethe nucleon in the QMC model. In its simplest form, thequarks are treated as a noninteracting Fermi gas confinedin a large bag, where the zero-point energy characterizedby the parameter Z is neglected. We consider the quarkmatter consisting of three flavor quarks ( u , d , and s ) andleptons ( e and µ ) in β equilibrium, where the relationsbetween chemical potentials are expressed as µ s = µ d = µ u + µ e , (15) µ µ = µ e . (16)The charge neutrality condition is written as n QP c = 23 n u − n d − n s − n e − n µ = 0 . (17)We use the current quark masses m u = m d = 5 . m s = 150 MeV in the calculation of quark matter,which is consistent with that used in the QMC model. Asfor the bag constant, we mainly take the value B / =210 .
854 MeV determined in the QMC model. In fact,the value of B in quark matter may be different fromthe one used to describe the nucleon in the QMC modeldue to the large density difference between the two cases.Therefore, we will compare results with different choicesof B for quark matter, so as to examine the influence ofthe bag constant.For quark matter described in the MIT bag model, theenergy density and pressure are given by ε QP = X i = u,d,s π Z k iF q k + m i k dk + B + X l = e,µ π Z k lF q k + m l k dk, (18) P QP = X i = u,d,s π Z k iF k dk p k + m i − B + X l = e,µ π Z k lF k dk p k + m l , (19)where the bag constant B is added to the energy den-sity and subtracted from the pressure. It is well knownthat the bag constant could significantly affect the EOSof quark matter and consequently influence the hadron-quark phase transition in neutron stars [23]. IV. MIXED PHASE WITHIN THE ENERGYMINIMIZATION METHOD
For the description of a hadron-quark mixed phase,we take into account the finite-size effects by using theWigner–Seitz approximation in our calculations. In theWigner–Seitz approximation, the system is divided intoequivalent and charge-neutral cells, where the coexistinghadronic and quark phases are separated by a sharp inter-face with finite surface tension. The leptons are assumedto be uniformly distributed throughout the cell. In theEM method, the equilibrium conditions between coex-isting phases are determined by minimization of the to-tal energy including surface and Coulomb contributions.Due to the competition between surface and Coulombenergies, the geometric structure of the mixed phasemay change from droplet to rod, slab, tube, and bubblewith increasing baryon density; these are known as pastaphases. The total energy density of the mixed phase isgiven by ε MP = χε QP + (1 − χ ) ε HP + ε surf + ε Coul , (20)where χ = V QP / ( V QP + V HP ) denotes the volume fractionof the quark phase. The energy densities, ε HP and ε QP ,are calculated from Eqs. (10) and (18), respectively. The first two terms of Eq. (20) represent the bulk contribu-tions, while the last two terms come from the finite-sizeeffects. The surface and Coulomb energy densities aregiven by ε surf = Dσχ in r D , (21) ε Coul = e δn c ) r D χ in Φ ( χ in ) , (22)withΦ ( χ in ) = D +2 (cid:18) − Dχ − /D in D − + χ in (cid:19) , D = 1 , , χ in − − ln χ in D +2 , D = 2 . (23)Here D = 1 , , r D being the size of the inner phase. χ in representsthe volume fraction of the inner phase, i.e., χ in = χ fordroplet, rod, and slab configurations, and χ in = 1 − χ for tube and bubble configurations. e = p π/
137 is theelectromagnetic coupling constant. δn c = n HP c − n QP c isthe charge-density difference between hadronic and quarkphases. σ denotes the surface tension at the hadron-quark interface, which plays a key role in determiningthe structure of the mixed phase [33–38]. However, con-siderable uncertainties exist regarding the value of σ , soit is often treated as a free parameter in the literature.Estimates based on different models suggest relativelysmall values σ <
30 MeV/fm [56–60], but much largervalues σ >
100 MeV/fm have also been obtained inRefs. [61, 62].In the present work, we prefer to calculate the surfacetension σ in the MIT bag model by using the multiplereflection expansion (MRE) method [56]. The derivationof σ within the MRE framework is based on the modifieddensity of states for the quark species i in an MIT bag,which is approximately given by dN i dk i = 6 (cid:20) k i V π − k i S π (cid:18) − π arctan k i m i (cid:19)(cid:21) , (24)where V and S are the volume and surface area of thebag, respectively. According to the relation between thesurface tension and the thermodynamic potential at zerotemperature, one can calculate the surface tension con-tributed from the quark species i by σ i = Z k bF k i π (cid:18) − π arctan k i m i (cid:19) (cid:20) µ i − q k i + m i (cid:21) dk i = 34 π ( k iF µ i − m i ( µ i − m i )3 − π (cid:20) µ i arctan k iF m i − k iF m i µ i + m i ln (cid:18) k iF + µ i m i (cid:19)(cid:21)(cid:27) , (25) where µ i = q k iF + m i is the chemical potential ofquarks. In the MRE method, the surface tension of thebag is given by a sum over all flavors, σ = P i = u,d,s σ i .Note that the dominant contribution to the surface ten-sion comes from the s quark, since its mass is much largerthan those of u and d quarks. On the other hand, σ isalso dependent on the density of quark matter.At a given baryon density, the thermodynamically sta-ble state is the one with the lowest energy among allconfigurations considered, which could be determined inthe EM method by minimizing the total energy with re-spect to all variables. The energy density of the mixedphase given in Eq. (20) is considered as a function of thefollowing variables: n p , n n , n u , n d , n s , n e , n µ , χ , and r D . The minimization should be performed under theconstraints of global charge neutrality and baryon num-ber conservation, which are expressed as n e + n µ − χ n u − n d − n s ) − (1 − χ ) n p = 0 , (26) χ n u + n d + n s ) + (1 − χ ) ( n p + n n ) = n b . (27)We introduce the Lagrange multipliers µ e and µ n for theconstraints, and then perform the minimization for thefunction w = ε MP − µ e h n e + n µ − χ n u − n d − n s ) − (1 − χ ) n p i − µ n h χ n u + n d + n s ) + (1 − χ ) ( n p + n n ) i . (28) By minimizing w with respect to the particle densities,we obtain the following equilibrium conditions for chem-ical potentials: µ u − ε Coul χ δn c = 13 µ n − µ e , (29) µ d + 2 ε Coul χ δn c = 13 µ n + 13 µ e , (30) µ s + 2 ε Coul χ δn c = 13 µ n + 13 µ e , (31) µ p + 2 ε Coul (1 − χ ) δn c = µ n − µ e , (32) µ µ = µ e . (33)Minimizing w with respect to the volume fraction χ leadsto the equilibrium condition for the pressures, P HP = P QP − ε Coul δn c (cid:20) χ (2 n u − n d − n s ) + 11 − χ n p (cid:21) ∓ ε Coul χ in χ in Φ ′ Φ ! , (34)where the sign of the last term is “ − ” for droplet, rod,and slab configurations, while it is“+” for tube and bub-ble configurations. It is clear that equilibrium conditionsfor two-phase coexistence are altered due to the inclu-sion of surface and Coulomb terms in the minimizationprocedure, and, as a result, they are different from theGibbs equilibrium conditions. Some additional terms ap-pearing in the equilibrium equations are caused by thesurface and Coulomb contributions. If we neglect the finite-size effects by taking the limit σ →
0, these addi-tional terms vanish and the equilibrium equations reduceto the Gibbs conditions. Generally, the pressure of themixed phase can be calculated from the thermodynamicrelation, P MP = n b ∂ ( ε MP /n b ) ∂n b , which is somewhat differ-ent from P HP and P QP . This is similar to the case ofnuclear liquid-gas phase transition at subnuclear densi-ties [40–43]. Furthermore, minimizing w with respect tothe size r D yields the well-known relation ε surf = 2 ε Coul ,which leads to the formula for the size of the inner phase, r D = " σDe ( δn c ) Φ / . (35)The properties of hadron-quark pasta phases are ob-tained by solving the above equilibrium equations at agiven baryon density n b . We compare the energy den-sity of the mixed phase with different pasta configura-tions and then determine the most stable configurationwith the lowest energy density. The hadron-quark mixedphase exists only in the density range where its energydensity is lower than that of both hadronic matter andquark matter. V. RESULTS AND DISCUSSION
In this section, we present the numerical results ofthe hadron-quark pasta phases using the EM methoddescribed in the previous section. To examine the de-pendence of results on the bag constant of quark mat-ter, we calculate and compare the results with B / =210 .
854 MeV and B / = 180 MeV in the present work.Note that B / = 210 .
854 MeV is determined in theQMC model by fitting the free nucleon properties. As forthe surface tension at the hadron-quark interface, we self-consistently calculate its value by using the MRE methodwithin the MIT bag model. Since the quark degrees offreedom are explicitly involved in the QMC model, thedescriptions of hadronic and quark phases are consideredto be consistent with each other at the quark level. Theproperties of neutron stars are then calculated using theEOS with the inclusion of quarks at high densities.
A. Hadron-quark pasta phases
Due to the competition between the surface andCoulomb energies, the geometric structure of a hadron-quark mixed phase is expected to change from dropletto rod, slab, tube, and bubble with increasing baryondensity. In the present calculations, we employ the EMmethod to investigate the pasta phases, where the equi-librium conditions between coexisting phases are derivedby the minimization of the total energy including sur-face and Coulomb contributions. Since the surface andCoulomb energies are positive, the energy densities of n b (fm −3 ) ε − ε ( σ = ) ( M e V / f m ) = 180 MeV (b)HadronQuarkdropletrodslabtubebubbleMC0.6 0.8 1.0 1.2 1.4 1.6 1.8 n b (fm −3 ) ε − ε ( σ = ) ( M e V / f m ) = 210.854 MeV (a) FIG. 1: Energy densities of the mixed phase obtained using the EM method relative to those of the Gibbs construction ( σ = 0).The filled circles indicate the transition between different pasta phases. The results of the Maxwell construction (MC) areshown by green dotted lines. pasta phases are higher than that of the Gibbs construc-tion without finite-size effects. Generally, the energy den-sity difference between different pasta shapes is rathersmall compared to the total energy density, but it is cru-cial to determine the transition of pasta shapes. In Fig. 1,we show the energy densities of pasta phases obtainedusing the EM method relative to those of the Gibbs con-struction (i.e., σ = 0). The filled circles indicate thetransition between different pasta phases. For compari-son, the energy densities of pure hadronic and pure quarkphases are respectively plotted by black dot-dashed anddashed lines, whereas the results of the Maxwell con-struction are shown by green dotted lines. To exam-ine the influence of the bag constant in quark matter,we present the results with B / = 210 .
854 MeV and B / = 180 MeV in the left and right panels, respec-tively. We note that the onset of the mixed phase occurswhen its energy density becomes lower than that of purehadronic matter. In the case of B / = 210 .
854 MeV thatis consistently used in the QMC model and quark mat-ter, the formation of quark droplets occurs at 0 . − ,which is significantly higher the onset of the mixed phaseat 0 . − obtained in the Gibbs construction. This isbecause the finite-size effects like surface and Coulombcontributions increase the energy density and delay theappearance of the mixed phase. As n b increases, otherpasta shapes, such as rod, slab, tube, and bubble, mayappear when one has the lowest energy density amongall configurations. Furthermore, the mixed phase ends at1 .
47 fm − where the energy density of pure quark matterbecomes lower than that of pasta phases. In the rightpanel of Fig. 1, we can see that using a smaller bagconstant of quark matter B / = 180 MeV leads to anearly onset of the structured mixed phase at 0 .
43 fm − and pure quark phase at 0 . − , while the behavior is qualitatively similar to that displayed in the left panelof Fig. 1. The influence of the bag constant B on thehadron-quark phase transition has been extensively dis-cussed in the literature [1, 23]. It is seen that the energydensity in the Maxwell construction is clearly higher thanthat of pasta phases, so the Maxwell construction is dis-favored in the present calculations. n b (fm −3 ) σ ( M e V / f m ) = 210.854 MeV d r o p l e t r o d s l a b t u b e b u bb l e FIG. 2: Surface tension σ as a function of the baryon density n b obtained in pasta phases. The surface tension at the hadron-quark interface playsa crucial role in determining the structure of pastaphases. In Fig. 2, we display the surface tension σ asa function of the baryon density n b obtained in pastaphases with B / = 210 .
854 MeV. In the present work,we use the MRE method to calculate the surface tension
200 300 400 500 600 μ i (MeV) σ i ( M e V / f m ) u x10σ s FIG. 3: Surface tension of the quark flavor i as a function ofits chemical potential calculated by Eq. (25). n b (fm −3 ) μ i ( M e V ) = 210.854 MeVμ s μ u d r o p l e t r o d s l a b t u b e b u bb l e FIG. 4: Chemical potential of the quark flavor i as a functionof the baryon density n b obtained in pasta phases. that is the sum over quark flavors, σ = P i = u,d,s σ i , with σ i given by Eq. (25). One can see that σ slightly de-creases with increasing n b and its value lies in the rangeof 43 −
48 MeV/fm . This behavior can be understoodby analyzing Eq. (25), where σ i is a function of the quarkmass m i and chemical potential µ i . The dominant con-tribution to the surface tension comes from the s quark,since its mass is much larger than that of u and d quarks.In Fig. 3, we plot σ i as a function of µ i for i = s and i = u . It is found that σ s is about one order higher than σ u , and it increases significantly as the chemical poten-tial increases. However, in the pasta phases, µ s shownin Fig. 4 decreases with increasing n b , which leads to thedecline of the surface tension σ as shown in Fig. 2. n b (fm −3 ) χ = 210.854 MeV GCEM FIG. 5: Volume fraction of the quark phase χ as a functionof the baryon density n b obtained in the mixed phase. Theresults using the EM method are compared to those of theGibbs construction (GC). −3 −2 −1 npeμ d su EM0.0 0.5 1.0 1.5 2.0n b μ(fm −3 )10 −3 −2 −1 npeμ d su GC Y i FIG. 6: Particle fraction Y i as a function of the baryon den-sity n b obtained using the EM method (upper panel) andGibbs construction (lower panel). The shaded areas denotethe mixed phase regions. In Fig. 5, we show the volume fraction of the quarkphase χ as a function of the baryon density n b in themixed phase with B / = 210 .
854 MeV. The results ob-tained in the EM method are compared to those of theGibbs construction. It is shown that χ continuously in-creases with increasing n b in the mixed phase. The be-havior of χ in the EM method is very similar to that ofthe Gibbs construction, but the density range of pastaphases is reduced due to the inclusion of finite-size ef-fects. There is no large jump in χ at the transition pointof pasta configurations, so the hadron-quark phase tran-sition described in the EM method is relatively smooth.Because of the monotonic increase of χ , more and morehadronic matter is converted into quark matter duringthe phase transition. In Fig. 6, we display the rela-tive particle fractions Y i = n i /n b as a function of thebaryon density n b with B / = 210 .
854 MeV. The re-sults obtained in the EM method and Gibbs constructionare shown in the upper and lower panels, respectively.The shaded areas denote the mixed phase regions. Onecan see that the matter at low densities consists of neu-trons, protons, and electrons, whereas the muons appearat n b = 0 .
11 fm − playing the same role as electrons.When the quark matter is present in the mixed phase,the fractions of quarks, Y u , Y d , and Y s , increase rapidlytogether with a decrease of neutron fraction Y n . Mean-while, Y e and Y µ decrease significantly, since the quarkmatter is negatively charged that can take the role ofelectrons to satisfy the constraint of global charge neu-trality. On the other hand, the hadronic matter in themixed phase is positively charged, so the proton fraction Y p increases at the beginning of the mixed phase. As n b increases, Y p and Y n decrease to very low values dueto increasing χ as shown in Fig. 5. At sufficiently highdensities, the hadronic matter completely disappears andthe transition to a pure quark phase occurs. It is seenthat Y u ≈ Y d ≈ Y s ≈ / B. Nuclear symmetry energy effects
It is well known that the nuclear symmetry energy E sym and its slope L can significantly affect the proper-ties of neutron stars, especially the neutron-star radius,tidal deformability, and crust structure, which are par-ticularly sensitive to the slope parameter L [14, 63, 64].Many efforts have been devoted to constraining the val-ues of E sym and L at saturation density based on as-trophysical observations and terrestrial nuclear experi-ments (see Refs. [64–66] and references therein). Accord-ing to available constraints summarized in Ref. [64], itwas found that the most probable values for the sym-metry energy and its slope at saturation density are E sym = 31 . ± . L = 58 . ± . L than that for E sym . In order to explore the influence of nuclear symme- n b (fm −3 ) E s y m ( M e V ) QMC(L=40)QMC(L=100)
FIG. 7: Symmetry energy E sym as a function of the baryondensity n b obtained in the QMC( L =40) and QMC( L =100)models. n b (fm −3 ) r C , r D ( f m ) d r o p l e t r o d s l a b t u b e b u bb l e B = 210.854 MeV r D r C FIG. 8: Size of the Wigner-Seitz cell ( r C ) and that of theinner phase ( r D ) as a function of n b obtained using the EMmethod. The results with QMC( L =40) and QMC( L =100)are shown by thick and thin lines, respectively. try energy on the hadron-quark pasta phases, we adjustthe isovector couplings, g ρ and Λ v , so as to obtain a stiffEOS in the QMC model and compare with the soft oneused. By taking Λ v = 0 and g ρ = 4 . E sym = 35 . L = 100 MeV atsaturation density, whereas the isoscalar saturation prop-erties remain unchanged. In Fig. 7, we show the symme-try energy E sym as a function of the baryon density n b obtained in the QMC( L =40) and QMC( L =100) models.It is found that the two models have the same E sym at n b ≃ .
106 fm − , but they display rather different den-sity dependence due to different slope parameter L . Itis noteworthy that the binding energies of finite nucleiare directly related to the symmetry energy at subsatu-ration density of 0 . − .
11 fm − [40, 67], which corre-sponds to the average density in nuclei. Therefore, theQMC( L =40) and QMC( L =100) models are expected toprovide similar descriptions for finite nuclei, since theyhave the same isoscalar properties and equal values of E sym at n b ≃ .
106 fm − . On the other hand, the differ-ent slope parameter L can significantly alter the stiffnessof neutron-star matter EOS at high densities, which maylead to considerable differences in the hadron-quark pastaphases inside neutron stars.We calculate and compare the properties of thehadron-quark pasta phases by using the QMC( L =40)and QMC( L =100) models for the hadronic matter. InFig. 8, we display the size of the Wigner-Seitz cell ( r C )and that of the inner phase ( r D ) as a function of thebaryon density n b obtained using the EM method with B / = 210 .
854 MeV. The results with QMC( L =40)and QMC( L =100) are shown by thick and thin lines,respectively. It has been reported in Ref. [33] thata larger symmetry energy slope L in hadronic mat-ter corresponds to an earlier onset of the hadron-quarkphase transition. In the present work, a similar trendis seen, namely the formation of quark droplets withQMC( L =100) occurs at a lower density compared to thecase with QMC( L =40). Furthermore, there are also dif-ferences between QMC( L =40) and QMC( L =100) in thesize of the pasta structure, especially in r C . One can seethat as the density increases, r D in the droplet, rod, andslab phases increases, whereas r D in the tube and bub-ble phases decreases. This is consistent with the mono-tonic increase of χ shown in Fig. 5. There are obviousdiscontinuities in r D and r C at the transition betweendifferent pasta shapes, which exhibit the character of thefirst-order transition. Generally, a large surface tension σ leads to a large radius of quark droplets. In the presentwork, we use the MRE method to calculate σ that is inthe range of 43 −
48 MeV/fm , which yields the radiusof quark droplets as r D ≈ − L =40)and QMC( L =100), we present in Table I the onset den-sities of the hadron-quark pasta phases and pure quarkmatter obtained using the EM method. It is found thatthe onset densities with QMC( L =100) are noticeablysmaller than those with QMC( L =40), and the differencesgradually decrease from the droplet to pure quark phases.This is because the fraction of hadronic matter monotoni-cally decreases during the hadron-quark phase transition,so that the influence of nuclear symmetry energy becomesweaker and weaker. Due to the same reason, the influ-ence of the bag constant B gets stronger and strongerduring the phase transition. One can see that the on-set densities of pure quark matter with B / = 180 MeVare much smaller than those with B / = 210 .
854 MeV. We find that the results are very sensitive to the bagconstant adopted. By using a small bag constant of B / = 180 MeV, the hadron-quark pasta phases aresignificantly shifted toward the lower density region com-pared to the results with B / = 210 .
854 MeV.
C. Properties of neutron stars
To investigate the influence of quark matter onneutron-star properties, we provide the EOS includingthe hadron-quark pasta phases and pure quark matter athigh densities. In Fig. 9, we plot the pressures as a func-tion of the baryon density n b for hadronic, mixed, andquark phases. The results with B / = 210 .
854 MeVand B / = 180 MeV are displayed in the left and rightpanels, respectively.The results with QMC( L =40) shown in the upper pan-els are compared to those with QMC( L =100) in the lowerpanels.The pressures of pasta phases obtained using the EMmethod are compared to those of the Gibbs and Maxwellconstructions. It can be seen that the pressures of pastaphases lie between those of the Gibbs and Maxwell con-structions. The pressures with the Maxwell constructionremain constant during the phase transition, whereasthose with the Gibbs construction increase with n b overa relatively broad range. The influence of the bag con-stant B on the EOS can be observed by comparing theleft and right panels. It is clearly shown that a smaller B leads to early onset of the mixed phase and pure quarkphase, namely the mixed phase is shifted toward a lowerdensity region. As a result, the pressures of pasta phaseswith B / = 180 MeV are much lower than those with B / = 210 .
854 MeV.Meanwhile, the effect of nuclear symmetry energy onthe EOS is seen by comparing the upper and lower pan-els. In the case with the QMC( L =100) model (lowerpanels), the onset of the mixed phase occurs at lowerdensities than with the QMC( L =40) model (upper pan-els), but there is no visible difference at the end of themixed phase.One can see that the behavior of P in pasta phases issimilar to that of the Gibbs construction, and no abruptjump is observed at the transition between different pastaconfigurations.The properties of neutron stars, such as the mass-radius relations and tidal deformabilities, can beobtained by solving the Tolman-Oppenheimer-Volkoff(TOV) equation using the EOS over a wide range ofdensities. In the present work, we adopt the Baym–Pethick–Sutherland (BPS) EOS [68] for the outer crustbelow the neutron drip density, while the inner crustEOS is based on a Thomas–Fermi calculation using theTM1( L =40) effective interaction of the RMF model [69].The crust EOS is matched to the QMC EOS of uni-form neutron-star matter at the crossing point of thetwo segments. At high densities, the hadron-quark pasta0 TABLE I: Onset densities of the hadron-quark pasta phases and pure quark matter obtained from different models for describingthe hadronic and quark phases.QMC model Bag model Onset density (fm − )droplet rod slab tube bubble quark L = 40 MeV B / = 210 .
854 MeV 0.700 0.830 0.997 1.238 1.339 1.466 L = 100 MeV B / = 210 .
854 MeV 0.545 0.702 0.921 1.211 1.321 1.453 L = 40 MeV B / = 180 MeV 0.431 0.522 0.586 0.710 0.755 0.803 L = 100 MeV B / = 180 MeV 0.306 0.428 0.514 0.681 0.737 0.791 n b (fm −3 ) P ( M e V / f m ) (b) B = 180 MeVL = 40 MeVHadronQuarkGCMCdropletrodslabtubebubble0.0 0.4 0.8 1.2 1.6 2.0 n b (fm −3 ) P ( M e V / f m ) (a) B = 210.854 MeVL = 40 MeV0.0 0.4 0.8 1.2 1.6 2.0 n b (fm −3 ) P ( M e V / f m ) (c) B = 210.854 MeVL = 100 MeV 0.0 0.2 0.4 0.6 0.8 1.0 1.2 n b (fm −3 ) P ( M e V / f m ) (d) B = 180 MeVL = 100 MeV FIG. 9: Pressures as a function of the baryon density n b for hadronic, mixed, and quark phases. The results of pasta phasesobtained using the EM method are compared to those of the Gibbs and Maxwell constructions. phases are taken into account by using the QMC modelfor the hadronic phase and the MIT bag model forthe quark phase. In Fig. 10, we present the result-ing mass-radius relations of neutron stars with the in-clusion of quarks at high densities. The results with B / = 210 .
854 MeV and B / = 180 MeV are shown inthe left and right panels, respectively. The mass measure-ments of PSR J1614–2230 (1 . ± . M ⊙ ) [4–6], PSRJ0348+0432 (2 . ± . M ⊙ ) [7], and PSR J0740+6620(2 . +0 . − . M ⊙ ) [8] are indicated by the horizontal bars. The simultaneous measurement of the mass and radiusfor PSR J0030+0451 by NICER is also shown by the starwith 68% and 95% confidence intervals [20].The constraints on R . inferred from GW170817 [10]are indicated by the horizontal line with arrows at bothends.We find that the maximum mass of neutron stars usinga pure hadronic EOS in the QMC( L =40) model is about2 . M ⊙ , which is considerably reduced as the hadron-quark phase transition is included. The star masses ob-1
10 12 14 16 18
R (km) M / M ⊙ PSR J0740+6620
PSR J0348+0432PSR J1614-2230
PSR J0030+0451 B = 210.854 MeV L = 40 MeV L = 100 MeV
HPGCEM 10 12 14 16 18
R (km) M / M ⊙ PSR J0740+6620
PSR J0348+0432PSR J1614-2230 = 180 MeV
PSR J0030+0451
L = 40 MeV L = 100 MeV
HPGCEM
FIG. 10: Mass-radius relations of neutron stars for different EOS. The results of a pure hadronic phase (solid lines) are comparedto those including quarks in the EM method (dotted lines) and Gibbs construction (dashed lines) with B / = 210 .
854 MeV(left panel) and B / = 180 MeV (right panel). The filled circles indicate the onset of the star containing a hadron-quarkmixed phase. The results with QMC( L =40) and QMC( L =100) are shown by thick and thin lines, respectively. The horizontalbars indicate the observational constraints of PSR J1614–2230 [4–6], PSR J0348+0432 [7], and PSR J0740+6620 [8]. Thesimultaneous measurement of the mass and radius for PSR J0030+0451 by NICER is also shown by the star with 68% and95% confidence intervals [20]. The horizontal line with arrows at both ends represents the constraints on R . inferred fromGW170817 [10]. M/M ⊙ Λ L = M e V L = M e V GW170817 B = 210.854 MeVHPEMGC
FIG. 11: Dimensionless tidal deformability Λ as a functionof the neutron-star mass M . The vertical line with arrows atboth ends represents the constraints on Λ . from the analysisof GW170817 [10]. tained in the EM method are slightly larger than thoseof the Gibbs construction due to finite-size effects. Bycomparing the left and right panels, it is clear that theresults are sensitive to the bag constant B adopted. Inthe case of B / = 210 .
854 MeV, the EM method withQMC( L =40) leads to a maximum mass of 1 . M ⊙ , where the hadron-quark pasta phases could be formed in the in-terior of massive stars with M > . M ⊙ , but no purequark phase exists. In fact, a canonical 1 . M ⊙ neutronstar would not be affected by the quark phase, since itscentral density is lower than the onset of a hadron-quarkmixed phase. However, in the case of B / = 180 MeV,the maximum mass obtained in the EM method withQMC( L =40) is about 1 . M ⊙ , which is much lower thanthe constraint of 2 M ⊙ . Therefore, the early onset of themixed phase with B / = 180 MeV is disfavored since itleads to a large reduction of the maximum neutron-starmass.It is interesting to see the effect of nuclear symme-try energy on the neutron-star properties. When theQMC( L =100) model is adopted, the maximum massof neutron stars is slightly higher than that in theQMC( L =40) model, but the radius is significantly in-creased. This is because the QMC( L =100) model pro-vides a stiff EOS, which results in relatively large radiiof neutron stars. It is shown that the large radius of R . = 13 . L =100) model isdisfavored by the constraints from GW170817. In con-trast, we obtain a radius of R . = 12 .
43 km in theQMC( L =40) model, which is compatible with recent ob-servational constraints from NICER and GW170817.The dimensionless tidal deformability of a neutron starcan be calculated fromΛ = 23 k ( R/M ) , (36)where k is the tidal Love number which is computed to-2gether with the TOV equation as described in Refs. [69–71]. In Fig. 11, we display the dimensionless tidal de-formability Λ as a function of the gravitational mass M of the star. It is shown that Λ decreases rapidly withincreasing M .The results with the inclusion of quarks using B / =210 .
854 MeV are compared to those using a purehadronic EOS in the QMC( L =40) and QMC( L =100)models. It is seen that the effect of the hadron-quarkphase transition is almost invisible, since the mixed phaseis present only in massive stars. For a canonical 1 . M ⊙ neutron star, we obtain Λ . = 450 in the QMC( L =40)model, which is consistent with the constraints from theanalysis of GW170817 [10], while Λ . = 885 obtained inthe QMC( L =100) model is incompatible with the con-straints. One can see that a small discrepancy in Λ usingthe EM method is observed, and it becomes more obvi-ous in the Gibbs construction. The reduction of Λ is acombined result of changes in the Love number k andcompactness parameter M/R caused by the appearanceof the hadron-quark phase transition.
VI. CONCLUSIONS
In the present work, we have studied the propertiesof the hadron-quark pasta phases, which are expectedto occur in the interior of massive neutron stars. Wehave employed the QMC model to describe the hadronicphase, where the internal quark structure of the nucleonis explicitly taken into account based on the MIT bagmodel. For the description of the quark phase, we haveadopted the same bag model as the one used in the QMCframework, so that the coexisting hadronic and quarkphases are treated in a consistent way. We have usedthe Wigner–Seitz approximation to describe the hadron-quark pasta phases, where the system is divided intoequivalent cells with a given geometric symmetry. Thehadronic and quark phases inside the cell are assumedto have constant densities and are separated by a sharpinterface. We computed the surface tension σ consis-tently in the bag model by using the MRE method. Itwas found that σ in the pasta phases slightly decreaseswith increasing density, and its value lies in the rangeof 43 −
48 MeV/fm . The dominant contribution to thesurface tension comes from the s quark, which is aboutone order higher than those from u and d quarks.We have investigated the hadron-quark mixed phaseusing the EM method, where the equilibrium conditionsfor coexisting phases are derived by minimization of thetotal energy including surface and Coulomb contribu-tions. Due to these finite-size effects, some additionalterms appear in the equilibrium conditions for pressuresand chemical potentials, which are different from theGibbs conditions. It was found that including the finite-size effects could delay the onset of the hadron-quarkmixed phase and shrink its density range significantly.On the other hand, the results of pasta phases are very sensitive to the bag constant B of quark matter. Byusing a consistent value of B / = 210 .
854 MeV inthe QMC( L =40) model and quark matter, a structuredmixed phase could be formed in the density range of0 . − .
47 fm − . When B / = 180 MeV is adoptedfor quark matter, the density range of the mixed phaseis reduced to 0 . − .
80 fm − . It was observed that theresults obtained in the EM method lie between those ofthe Gibbs and Maxwell constructions.We examined the influence of nuclear symmetry energyon the hadron-quark phase transition by using the QMCmodels with different slope parameter L . It was foundthat the onset densities of the hadron-quark mixed phaseobtained with QMC( L =100) are obviously smaller thanthose with QMC( L =40), but there is no visible differenceat the end of the mixed phase.We applied the EOS including the hadron-quark phasetransition to study the properties of neutron stars. Inthe present calculations, a pure hadronic EOS of theQMC( L =40) model predicts a maximum neutron-starmass of 2 . M ⊙ , while the resulting radius and tidaldeformability of a canonical 1 . M ⊙ neutron star are R . = 12 .
43 km and Λ . = 450, respectively. Theseresults are compatible with current constraints from as-trophysical observations.However, the QMC( L =100) model results in ratherlarge radii and tidal deformabilities of neutron stars,which are disfavored by the constraints from GW170817.The recent gravitational-wave event GW190814 has trig-gered many efforts to explore the possibility of existinga super-massive neutron star of ≈ . M ⊙ . Within theQMC model used in the present work, it is unlikely thatsuch massive neutron star can be stable. It may be pos-sible to raise the maximum neutron-star mass by intro-ducing nonlinear self-couplings of the meson fields in theQMC model, as discussed in the RMF approach [72].This possibility will be explored in future studies.When the deconfinement phase transition is included,it is found that using the EM method with B / =210 .
854 MeV and QMC( L =40), the hadron-quark pastaphases could be formed in the interior of massive starswith M > . M ⊙ , but no pure quark matter exists.Meanwhile, the maximum neutron-star mass is reducedto 1 . M ⊙ in this case, while the properties of a canonical1 . M ⊙ neutron star remain unchanged. We emphasizethat although the hadron-quark mixed phases do not ap-preciably change the neutron-star bulk properties, theycould be important for studying cooling observations. Ifa small bag constant B / = 180 MeV is adopted forquark matter, the maximum mass is reduced to 1 . M ⊙ in the QMC( L =40) model, which is much lower than theconstraint of 2 M ⊙ . This implies that the early onset ofthe mixed phase with B / = 180 MeV for quark matteris disfavored due to its large reduction of the maximumneutron-star mass. In this work, we used a simple modelfor quark matter and neglected possible interactions be-tween quarks, which need to be investigated in futurestudies.3 Acknowledgment
This work was supported in part by the National Nat-ural Science Foundation of China (Grants No. 11675083 and No. 11775119). [1] G. Baym, T. Hatsuda, T. Kojo, P. D. Powell, Y. Song,and T. Takatsuka, Rept. Prog. Phys. , 056902 (2018).[2] N. K. Glendenning, Phys. Rep. , 393 (2001).[3] F. Weber, Prog. Part. Nucl. Phys. , 193 (2005).[4] P. B. Demorest, T. Pennucci, S. M. Ranson, M. S. E.Roberts, and J. W. T. Hessels, Nature (London) ,1081 (2010).[5] E. Fonseca et al. , Astrophys. J. , 167 (2016).[6] Z. Arzoumanian et al. , Astrophys. J. Suppl. , 37(2018).[7] J. Antoniadis et al. , Science , 6131 (2013).[8] H. T. Cromartie et al. , Nat. Astron. , 72 (2020).[9] B. P. Abbott et al. (LIGO Scientifc Collaboration andVirgo Collaboration), Phys. Rev. Lett. , 161101(2017).[10] B. P. Abbott et al. (LIGO Scientifc Collaboration andVirgo Collaboration), Phys. Rev. Lett. , 161101(2018).[11] B. P. Abbott et al. (LIGO Scientifc Collaboration andVirgo Collaboration), Phys. Rev. X , 011001 (2019).[12] I. Tews, J. Margueron, and S. Reddy, Phys. Rev. C ,045804 (2018).[13] S. De, D. Finstad, J. M. Lattimer, D. A. Brown, E.Berger, and C. M. Biwer, Phys. Rev. Lett. , 091102(2018).[14] F. J. Fattoyev, J. Piekarewicz, and C. J. Horowitz, Phys.Rev. Lett. , 172702 (2018).[15] T. Malik, N. Alam, M. Fortin, C. Providˆencia, B. K.Agrawal, T. K. Jha, B. Kumar, and S. K. Patra, Phys.Rev. C , 035804 (2018).[16] Z. Y. Zhu, E. P. Zhou, and A. Li, Astrophys. J. , 98(2018).[17] B. P. Abbott et al. (LIGO Scientifc Collaboration andVirgo Collaboration), Astrophys. J. Lett. , L3 (2020).[18] B. P. Abbott et al. (LIGO Scientifc Collaboration andVirgo Collaboration), Astrophys. J. Lett. , L44(2020).[19] T. E. Riley et al. , Astrophys. J. Lett. , L21 (2019).[20] M. C. Miller et al. , Astrophys. J. Lett. , L24 (2019).[21] N. K. Glendenning, Phys. Rev. D , 1274 (1992).[22] K. Schertler, S. Leupold, and J. Schaffner-Bielich, Phys.Rev. C , 025801 (1999).[23] K. Schertler, C. Greiner, J. Schaffner-Bielich, and M. H.Thoma, Nucl. Phys. A677 , 463 (2000).[24] A. W. Steiner, M. Prakash, and J. M. Lattimer, Phys.Lett. B , 239 (2000).[25] G. F. Burgio, M. Baldo, P. K. Sahu, and H.-J. Schulze,Phys. Rev. C , 025802 (2002).[26] D. P. Menezes and C. Providˆencia, Phys. Rev. C ,035804 (2003).[27] F. Yang and H. Shen, Phys. Rev. C , 025801 (2008).[28] J. Xu, L. W. Chen, C. M. Ko, and B. A. Li, Phys. Rev.C , 055803 (2010).[29] H. Chen, G. F. Burgio, H.-J. Schulze, and N. Yasutake,Astron. Astrophys. , A13 (2013). [30] M. Orsaria, H. Rodrigues, F. Weber, and G. A. Contrera,Phys. Rev. C , 015806 (2014).[31] X. H. Wu and H. Shen, Phys. Rev. C 96, 025802 (2017).[32] A. Bhattacharyya, I. N. Mishustin, and W. Greiner, J.Phys. G , 025201 (2010).[33] X. H. Wu and H. Shen, Phys. Rev. C 99, 065802 (2019).[34] H. Heiselberg, C. J. Pethick, and E. F. Staubo, Phys.Rev. Lett. , 1355 (1993).[35] T. Endo, T. Maruyama, S. Chiba, and T. Tatsumi, Prog.Theor. Phys. , 337 (2006).[36] T. Maruyama, S. Chiba, H.-J. Schulze, and T. Tatsumi,Phys. Rev. D , 123015 (2007).[37] N. Yasutake, R. Lastowiecki, S. Beni´c, D. Blaschke, T.Maruyama, and T. Tatsumi, Phys. Rev. C , 065803(2014).[38] W. M. Spinella, F. Weber, G. A. Contrera, and M. G.Orsaria, Eur. Phys. J. A , 61 (2016).[39] T. Tatsumi, M. Yasuhira, and D. Voskresensky, Nucl.Phys. A , 359 (2003).[40] S. S. Bao, J. N. Hu, Z. W. Zhang, and H. Shen, Phys.Rev. C , 045802 (2014).[41] G. Baym, H. A. Bethe, and C. J. Pethick, Nucl. Phys. A , 225 (1971).[42] J. M. Lattimer, C. J. Pethick, D. G. Ravenhall, and D.Q. Lamb, Nucl. Phys. A , 646 (1985).[43] J. M. Lattimer and F. D. Swesty, Nucl. Phys. A , 331(1991).[44] P. A. M. Guichon, Phys. Lett. B , 235 (1988).[45] K. Saito, K. Tsushima, and A. W. Thomas, Nucl. Phys.A , 339 (1996).[46] H. Muller and B. K. Jennings, Nucl. Phys. A , 55(1998).[47] P. Yue and H. Shen, Phys. Rev. C , 045807 (2006).[48] K. Saito, K. Tsushima, and A. W. Thomas, Prog. Part.Nucl. Phys. , 1 (2007).[49] P. A. M. Guichon, J. R. Stone, A. W. Thomas, Prog.Part. Nucl. Phys. , 262 (2018).[50] G. Grams, A. M. Santos, P. K. Panda, C. C. Providˆencia,D. P. Menezes, Phys. Rev. C , 055807 (2017).[51] T. F. Motta, A. M. Kalaitzis, S. Antic, P. A. M. Gui-chon, J. R. Stone, A. W. Thomas, Astrophys. J. ,159 (2019).[52] P. K. Panda, D. P. Menezes, and C. Providˆencia, Phys.Rev. C , 025207 (2004).[53] P. K. Panda, D. P. Menezes, and C. Providˆencia, Phys.Rev. C , 058801 (2004).[54] P. K. Panda, C. Providˆencia, and D. P. Menezes, Phys.Rev. C , 045801 (2010).[55] P. K. Panda, A. M. S. Santos, D. P. Menezes, and C.Providˆencia, Phys. Rev. C , 055802 (2012).[56] M. S. Berger and R. L. Jaffe, Phys. Rev. C , 213 (1987); , 125018(2010).[58] M. B. Pinto, V. Koch, and J. Randrup, Phys. Rev. C , , 015804(2017).[60] G. Lugones and A. G. Grunfeld, Phys. Rev. C , 035804(2019).[61] M. G. Alford, K. Rajagopal, S. Reddy, and F. Wilczek,Phys. Rev. D , 074017 (2001).[62] G. Lugones, A. G. Grunfeld, and M. A. Ajmi, Phys. Rev.C , 045803 (2013).[63] S. S. Bao and H. Shen, Phys. Rev. C , 015807 (2015).[64] M. Oertel, M. Hempel, T. Kl¨ahn, and S. Typel, Rev.Mod. Phys. , 015007 (2017).[65] B. A. Li, L. W. Chen, and C. M. Ko, Phys. Rep. ,113 (2008). [66] I. Tews, J. M. Lattimer, A. Ohnishi, and E. E. Kolomeit-sev, Astrophys. J. , 105 (2017).[67] Z. Zhang and L. W. Chen, Phys. Lett. B726 , 234 (2013).[68] G. Baym, C. Pethick, and P. Sutherland, Astrophys. J. , 299 (1971).[69] F. Ji, J. N. Hu, S. S. Bao, and H. Shen, Phys. Rev. C , 045801 (2019).[70] T. Hinderer, Astrophys. J. , 1216 (2008); , 964(E)(2009).[71] S. Postnikov, M. Prakash, and J. M. Lattimer, Phys. Rev.D , 024016 (2010).[72] F. J. Fattoyev, C. J. Horowitz, J. Piekarewicz, and B.Reed, Phys. Rev. C102