# New temperature dependent hyperonic equation of state: Application to rotating neutron star models and $I-Q$ relations

Miguel Marques, Micaela Oertel, Matthias Hempel, Jérôme Novak

► To cite this version:

Miguel Marques, Micaela Oertel, Matthias Hempel, Jérôme Novak. New temperature dependent hyperonic equation of state: Application to rotating neutron star models and $I-Q$ relations. Physical Review C, 2017, 96 (4), pp.045806. 10.1103/PhysRevC.96.045806. hal-01645244

HAL Id: hal-01645244

https://hal.science/hal-01645244v1

Submitted on 5 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# New temperature dependent hyperonic equation of state: Application to rotating neutron star models and I-Q relations

Miguel Marques,¹,* Micaela Oertel,¹,† Matthias Hempel,²,‡ and Jérôme Novak¹,§

¹LUTH, Observatoire de Paris, PSL Research University, CNRS, Université Paris Diderot, Sorbonne Paris Cité,
5 place Jules Janssen, 92195 Meudon, France

²Universität Basel, Departement Physik, Klingelbergstrasse 82, 4056 Basel, Switzerland

(Received 12 June 2017; revised manuscript received 23 August 2017; published 10 October 2017)

In this work we present a newly constructed equation of state (EoS), applicable to stellar core collapse and neutron star mergers including the entire baryon octet. Our EoS is compatible with the main constraints from nuclear physics and, in particular, with a maximum mass for cold β-equilibrated neutron stars of $2M_⊙$ in agreement with recent observations. As an application of our new EoS, we compute numerical stationary models for rapidly (rigidly) rotating hot neutron stars. We consider maximum masses of hot stars, such as protoneutron stars or hypermassive neutron stars in the postmerger phase of binary neutron star coalescence. The universality of I-Q relations at nonzero temperature for fast rotating models, comparing a purely nuclear EoS with its counterparts containing Λ hyperons or the entire baryon octet, respectively, is discussed, too. We find that the I-Q universality is broken in our models when thermal effects become important, independent on the presence of entropy gradients. Thus, the use of I-Q relations for the analysis of protoneutron stars or merger remnant data, including gravitational wave signals from the last stages of binary neutron star mergers, should be regarded with care.

DOI: 10.1103/PhysRevC.96.045806

## I. INTRODUCTION

Neutron stars are among the most extreme objects in the universe. They represent unique laboratories for probing strongly interacting matter at ultrahigh densities, exceeding that in atomic nuclei, as well as gravity for strong fields. They are formed in a core-collapse supernovae (CCSN) and cool down mainly by neutrino emission to form a catalyzed cold neutron star on a time scale of several minutes. Thus, in the early post bounce phase, as protoneutron stars (PNSs), they do not contain only ultradense matter, but they are hot objects, too, reaching temperatures of the order ~50 MeV [1–3]. In addition, matter in a PNS is not transparent to neutrinos, being thus lepton rich. Temperature and lepton content are important ingredients to describe the physics of PNSs, be it matter composition and stability of the PNS against collapse to a black hole [2,4–13] or dynamical properties such as frequencies and damping times of quasinormal modes and consequently the emitted gravitational wave signal [14].

In the postmerger phase of a binary neutron star coalescence, a rapidly rotating neutron star could be formed which temporarily resists a black hole collapse [15], even if its mass exceeds the maximum mass of a cold nonrotating neutron star. Within these merger remnants, temperatures of the same order as for CCSN and PNSs are reached. Both, PNSs and merger remnants can rotate at rather high frequencies, with potentially a differential rotation profile.

The temperatures of 50–100 MeV reached in these astrophysical environments are such that thermal effects on

the equation of state (EoS) become important. They have in particular a non-negligible effect on the composition, favoring the production of non-nucleonic degrees of freedom such as hyperons, nuclear resonances, or mesons. Even a transition to the quark-gluon plasma could take place. The impact of these additional particles on the evolution of PNSs has a long history; see, e.g., [16] for an early review. Most models employ EoSs for homogeneous matter, neglecting inhomogeneous matter in the outer layers and the formation of a crust; see, e.g., [17–25].

Currently, only a few EoSs are available covering in a consistent way the whole necessary domain in temperature $T$, baryon number density $n_B$, and electron fraction, $Y_e = n_e/n_B$ where $n_e$ is the electron number density. We will call them general purpose EoSs. In the last years, a series of new EoS models has been developed, see, e.g., [10,12,26–34], focused mainly on the treatment of the inhomogeneous part and correct nuclear abundances, and/or nuclear interactions at high densities. Triggered by investigations of stellar black hole formation, some effort has recently been devoted to extend the existing purely nuclear models to include non-nucleonic degrees of freedom—hyperons, pions or quarks—at high densities and temperatures, too; see, e.g., [6,35–40]. The latter are very important for the description of PNSs and merger remnants in view of the high densities combined with high temperatures which are attained within these objects. However, up to now none of these extended models is really satisfactory, since they are either not compatible with constraints from nuclear physics or neutron star masses [41–43], or contain only a limited selection of additional degrees of freedom, typically Λ hyperons. Here, we will present for the first time an EoS taking into account the entire baryon octet and being well compatible with the main present constraints.

As an application of our new EoS, we will compute stationary models of (rotating) hot stars and study the influence of hyperons on PNS and merger remnant properties. Most studies

<footnote>*miguel.marques@obspm.fr</footnote>

<footnote>†micaela.oertel@obspm.fr</footnote>

<footnote>‡matthias.hempel@unibas.ch</footnote>

<footnote>§jerome.novak@obspm.fr</footnote>
of PNS evolution are based on sequences of quasiequilibrium models,<fnref n="1" /> an assumption which is well justified in view of the hydrodynamic time scale ($\lesssim 10^{-3}$ s) being much smaller than the time scale on which thermodynamic properties are modified considerably ($\approx 1$ s); see, e.g., [1,5,45]. Although merger remnants cannot be well approximated as being in a quasiequilibrium state, interest in stationary models of these stars arise in order to understand the physical mechanism stabilizing the hypermassive star without performing a complete numerical merger simulation.

Stationary models of (cold) relativistic stars have been extensively explored in the literature. The first models, the famous Tolman-Oppenheimer-Volkoff solutions, describing spherically symmetric (therefore nonrotating) stars, date from the late 1930s [46,47]. Hartle and Thorne proposed the first axisymmetric rotating solutions from a perturbative approach in the slow rotation approximation [48]. Nowadays several publicly available codes are able to obtain precise numerical solutions up to the mass shedding limit [49,50] at the Kepler frequency; see, e.g., the textbook by Friedman and Stergioulas [51]. However, all these solutions only treat cold $\beta$-equilibrated stars with a barotropic EoS. Goussard et al. [52,53] have introduced the first models including the effect of finite temperature, restricting their solution however to the isentropic or isothermal case in $\beta$ equilibrium with several fixed overall lepton fractions, where the EoS effectively reduces to a barotropic one. References [54,55] propose general solutions within a perturbative slow rotation approach. In this work, we follow Ref. [45] to consistently compute stationary rapidly rotating hot stars based on the publicly available numerical library LORENE [56].

The paper is organized as follows. In Sec. II we present the new EoS model and in Sec. III some of its properties and in particular its compatibility with available constraints. In Sec. IV we present the formalism to treat stationary rotating relativistic stars at nonzero temperature. Section V shows first applications of our models, discussing maximum masses of hot stars and I-Q relations, i.e., universal relations among the moment of inertia and the quadrupole moment. We conclude in Sec. VI. Throughout the paper we use natural units with $c = \hbar = k_B = 1$ where appropriate.

## II. EQUATION OF STATE

Although the transition to the quark-gluon plasma is very interesting, as it could facilitate the supernova explosion [36], explain some $\gamma$-ray bursts [57], or—within the scenario of “quark-novae”—some unusual supernova light curves [58–60], we will concentrate here on hyperonic degrees of freedom. Presently available general purpose EoS models including all hyperons and covering the entire range in baryon number density $n_B$, temperature $T$, and hadronic charge fraction $Y_Q = n_Q/n_B = Y_e$,<fnref n="2" /> necessary for applications in CCSN or binary mergers, are either not compatible with some constraints from nuclear physics and/or a neutron star

maximum mass of $2M_\odot$; see, e.g., [6,38] or consider only $\Lambda$ hyperons (e.g., [40]). Our new EoS, taking into account the entire baryon octet, is well compatible with the main present constraints; see Sec. III for details.

### A. Statistical model for inhomogeneous matter

At subsaturation densities and low temperatures, nucleonic matter is unstable with respect to variations in the particle densities and becomes inhomogeneous, i.e., nuclei or more generally nuclear clusters are formed. The critical temperature is of the order $\sim 15$ MeV just below saturation and decreases to about 1 MeV at lower densities. Below a density of roughly $n_B \sim 10^{-4} \text{ fm}^{-3}$, the cluster size is very small compared with its mean free path, such that matter can be described as a noninteracting gas of nuclei, nucleons, and leptons in thermodynamic equilibrium. This approach is generally called “nuclear statistical equilibrium” (NSE). In the last years several models have been developed to go beyond a pure NSE and take into account nucleon interactions and the interaction of clusters with the surrounding medium at higher densities (see, e.g., [27,28,31,61–64]). In stellar matter particular attention has to be paid to the interplay between the short-range nuclear interaction and the long-range Coulomb interaction, which determines sizes and shapes of the nuclear clusters and influences thus strongly the transition to homogeneous matter [65,66].

In the present EoS, clustered matter is described within the extended NSE model of Hempel and Schaffner-Bielich [10,28]. Nuclei are treated as classical Maxwell-Boltzmann particles. For the description of nucleons, a relativistic mean field (RMF) approach is employed (see Sec. II B for details) with the same parametrization as for the description of homogeneous matter. Several thousands of nuclei are considered, including light ones other than the $\alpha$ particle. If available, nuclear binding energies are taken from experimental measurements [67]. In particular for neutron rich nuclei, where no measurement exists, they are complemented with values from theoretical nuclear structure calculations [68]. Several corrections are considered to describe the modifications of cluster properties in medium: screening of the Coulomb energies by the surrounding gas of electrons, excited states, and excluded-volume effects.

### B. Homogeneous matter

Homogeneous matter is described within a phenomenological RMF. The basic idea of this type of model is that the interaction between baryons is mediated by meson fields inspired by the meson exchange models of the nucleon-nucleon interaction. Within RMF models, these are, however, not real mesons, but introduced on a phenomenological basis with their quantum numbers in different interaction channels. The coupling constants are adjusted to a chosen set of nuclear observables. Earlier models introduce nonlinear self-couplings of the meson fields in order to reproduce correctly nuclear matter saturation properties, whereas more recently density-dependent couplings between baryons and the meson fields have been widely used. The literature on those models is large and many different parametrizations exist (see, e.g., [69]).

<footnote n="1">¹However, see Ref. [44] for a first dynamical study.</footnote>

<footnote n="2">²$n_Q$ represents the total hadronic charge density.</footnote>
In the present paper, we will use models with density dependent couplings. The Lagrangian density can be written in the following form:<fnref n="3" />

$$ \mathcal{L} = \sum_{j \in B} -\bar{\psi}_j (\gamma_\mu \partial^\mu + m_j - g_{\sigma j} \sigma - g_{\sigma^* j} \sigma^* \\ - i g_{\omega j} \gamma_\mu \omega^\mu - i g_{\phi j} \gamma_\mu \phi^\mu - i g_{\rho j} \gamma_\mu \bar{\rho}^\mu \cdot \vec{I}_j) \psi_j \\ - \frac{1}{2} (\partial_\mu \sigma \partial^\mu \sigma + m_\sigma^2 \sigma^2) \\ - \frac{1}{2} (\partial_\mu \sigma^* \partial^\mu \sigma^* + m_{\sigma^*}^2 \sigma^{*2}) \\ - \frac{1}{4} W_{\mu\nu}^\dagger W^{\mu\nu} - \frac{1}{4} P_{\mu\nu}^\dagger P^{\mu\nu} - \frac{1}{4} \vec{R}_{\mu\nu}^\dagger \cdot \vec{R}^{\mu\nu} \\ - \frac{1}{2} m_\omega^2 \omega_\mu \omega^\mu \\ - \frac{1}{2} m_\phi^2 \phi_\mu \phi^\mu - \frac{1}{2} m_\rho^2 \bar{\rho}_\mu \cdot \bar{\rho}^\mu, \quad (1) $$

where $\psi_j$ denotes the field of baryon j, and $W_{\mu\nu}, P_{\mu\nu}, \vec{R}_{\mu\nu}$ are the field tensors of the vector mesons, $\omega$ (isoscalar), $\phi$ (isoscalar), and $\rho$ (isovector), of the form

$$ V^{\mu\nu} = \partial^{\mu} V^{\nu} - \partial^{\nu} V^{\mu}. \quad (2) $$

$\sigma, \sigma^*$ are scalar-isoscalar meson fields, coupling to all baryons ($\sigma$) and to strange baryons ($\sigma^*$), respectively. Some models introduce an additional scalar-isovector coupling via a $\delta$ meson, which we do not consider here. The values of the baryon masses $m_j$ are chosen as follows: $m_n = 939.565\,346, m_p = 938.272\,013, m_\Lambda = 1115.683, m_\Sigma = 1190, m_{\Xi^-} = 1321.68, m_{\Xi^0} = 1314.83$ MeV.

In mean field approximation, the meson fields are replaced by their respective mean-field expectation values, which are given in uniform matter as

$$ m_{\sigma}^{2} \bar{\sigma} = \sum_{j \in B} g_{\sigma j} n_{j}^{s}, \quad (3) $$

$$ m_{\sigma^*}^2 \bar{\sigma}^* = \sum_{j \in B} g_{\sigma^* j} n_j^s, \quad (4) $$

$$ m_{\omega}^{2} \bar{\omega} = \sum_{j \in B} g_{\omega j} n_{j}, \quad (5) $$

$$ m_{\phi}^{2} \bar{\phi} = \sum_{j \in B} g_{\phi j} n_{j}, \quad (6) $$

$$ m_{\rho}^{2} \bar{\rho} = \sum_{j \in B} g_{\rho i} t_{3j} n_{j}, \quad (7) $$

where $\bar{\rho} = \langle \rho_3^0 \rangle, \bar{\omega} = \langle \omega^0 \rangle, \bar{\phi} = \langle \phi^0 \rangle$, and $t_{3j}$ represents the third component of isospin of baryon j with the convention

that $t_{3p} = 1/2$. The scalar density of baryon j is given by

$$ \begin{aligned} n_j^s &= \langle \bar{\psi}_j \psi_j \rangle \\ &= \frac{1}{\pi^2} \int k^2 \frac{M_j^*}{\sqrt{k^2 + M_j^{*2}}} \{ f[\epsilon_j(k)] + \bar{f}[\epsilon_j(k)] \} dk, \end{aligned} \quad (8) $$

and the number density by

$$ n_j = i \langle \bar{\psi}_j \gamma^0 \psi_j \rangle = \frac{1}{\pi^2} \int k^2 (f[\epsilon_j(k)] - \bar{f}[\epsilon_j(k)]) dk. \quad (9) $$

$f$ and $\bar{f}$ represent here the occupation numbers of the respective particle and antiparticle states with $\epsilon_j(k) = \sqrt{k^2 + M_j^{*2}}$, and effective chemical potentials $\mu_j^*$. They reduce to a step function at zero temperature. The effective baryon mass $M_j^*$ depends on the scalar mean fields as

$$ M_j^* = M_j - g_{\sigma j} \bar{\sigma} - g_{\sigma^* j} \bar{\sigma}^*, \quad (10) $$

and the effective chemical potentials are related to the chemical potentials via

$$ \mu_j^* = \mu_j - g_{\omega j} \bar{\omega} - g_{\rho j} t_{3j} \bar{\rho} - g_{\phi j} \bar{\phi} - \Sigma_0^R. \quad (11) $$

The rearrangement term $\Sigma_0^R$ is present in models with density-dependent couplings of meson $M$ to baryon j,

$$ g_{Mj}(n_B) = g_{Mj}(n_0)h_M(x), \quad x = n_B/n_0, \quad (12) $$

to ensure thermodynamic consistency. It is given by

$$ \Sigma_0^R = \sum_{j \in B} \left( \frac{\partial g_{\omega j}}{\partial n_j} \bar{\omega} n_j + t_{3j} \frac{\partial g_{\rho j}}{\partial n_j} \bar{\rho} n_j + \frac{\partial g_{\phi j}}{\partial n_j} \bar{\phi} n_j \right. \\ \left. - \frac{\partial g_{\sigma j}}{\partial n_j} \bar{\sigma} n_j^s - \frac{\partial g_{\sigma^* j}}{\partial n_j} \bar{\sigma}^* n_j^s \right). \quad (13) $$

The density $n_0$ is a normalization constant, usually taken to be the saturation density $n_0 = n_{\text{sat}}$ of symmetric nuclear matter.

In the present paper we will consider the DD2 parametrization [27], where the following form for the density dependence of the isoscalar couplings is assumed [27]:

$$ h_M(x) = a_M \frac{1 + b_M(x + d_M)^2}{1 + c_M(x + d_M)^2} \quad (14) $$

and

$$ h_M(x) = a_M \exp[-b_M(x - 1)] - c_M(x - d_M) \quad (15) $$

for the isovector ones. The values of the parameters $a_M, b_M, c_M$, and $d_M$ are listed in Ref. [27].

Similar to many recent works [40, 70, 71], for the hyperonic coupling constants, we will follow a symmetry inspired procedure. The individual isoscalar vector meson-baryon couplings are expressed in terms of $g_{\omega N}$ and a few additional parameters, $\alpha, \theta, z = g_1/g_8$, see e.g., [72], as follows:

$$ \begin{aligned} \frac{g_{\omega \Lambda}}{g_{\omega N}} &= \frac{1 - \frac{2z}{\sqrt{3}}(1 - \alpha) \tan \theta}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha) \tan \theta}, \\ \frac{g_{\phi \Lambda}}{g_{\omega N}} &= - \frac{\tan \theta + \frac{2z}{\sqrt{3}}(1 - \alpha)}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha) \tan \theta}, \end{aligned} $$

<footnote n="3">³Note that we work here with a locally flat Minkowski metric $\eta^{\mu\nu}$. For the $\gamma$ matrices, we use the anticommutation relation $\{\gamma^{\mu}, \gamma^{\nu}\} = 2\eta^{\mu\nu}$.</footnote>
$$
\begin{aligned}
\frac{g_{\omega\Xi}}{g_{\omega N}} &= \frac{1 - \frac{z}{\sqrt{3}}(1 + 2\alpha)\tan\theta}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha)\tan\theta}, \\
\frac{g_{\phi\Xi}}{g_{\omega N}} &= -\frac{\tan\theta + \frac{z}{\sqrt{3}}(1 + 2\alpha)}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha)\tan\theta}, \\
\frac{g_{\omega\Sigma}}{g_{\omega N}} &= \frac{1 + \frac{2z}{\sqrt{3}}(1 - \alpha)\tan\theta}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha)\tan\theta}, \\
\frac{g_{\phi\Sigma}}{g_{\omega N}} &= \frac{-\tan\theta + \frac{2z}{\sqrt{3}}(1 - \alpha)}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha)\tan\theta}, \\
\frac{g_{\phi N}}{g_{\omega N}} &= -\frac{\tan\theta + \frac{z}{\sqrt{3}}(1 - 4\alpha)}{1 - \frac{z}{\sqrt{3}}(1 - 4\alpha)\tan\theta}.
\end{aligned}
\qquad (16)
$$

Assuming an underlying SU(6) symmetry, we will take $\tan\theta = 1/\sqrt{2}$, corresponding to ideal $\omega-\phi$ mixing, $\alpha = 1$, and $z = 1/\sqrt{6}$. Extending the above procedure to the isovector sector would lead to contradictions with the observed nuclear symmetry energy. $g_{\rho N}$ is therefore left as a free parameter and the remaining hyperonic isovector couplings are fixed by isospin symmetry.

The information from hypernuclear data on hyperonic single-particle mean field potentials is then used to constrain the scalar coupling constants. The potential for particle j in k-particle matter is given by

$$U_j^{(k)}(n_k) = M_j^* - M_j + \mu_j - \mu_j^* \quad (17)$$

We will assume here standard values [40,70,73] in symmetric nuclear matter at saturation density $n_{\text{sat}}$: $U_{\Lambda}^{(N)}(n_{\text{sat}}) = -30$ MeV, $U_{\Xi}^{(N)}(n_{\text{sat}}) = -18$ MeV, and $U_{\Sigma}^{(N)}(n_{\text{sat}}) = +30$ MeV. The resulting values are in the range obtained by calculating directly properties of single $\Lambda$ hypernuclei; see Refs. [73,74].

Apart from a few light double-$\Lambda$ hypernuclei, that constrain only the low density behavior, almost no information is available on the hyperon-hyperon (YY) interaction and the corresponding couplings, in particular $\sigma^*$ and $\phi$ are only very poorly constrained. As mentioned above, we fix the $\phi$ couplings via the relations in Eqs. (16) and neglect $\sigma^*$ for simplicity in the main version of our EoS, named “DD2Y” hereafter. Without the coupling to $\sigma^*$, the YY interaction is very repulsive already at low densities. We obtain $U_{\Lambda}^{(\Lambda)}(n_{\text{sat}}/5) = 7$ MeV, $U_{\Xi}^{(\Xi)}(n_{\text{sat}}/5) = 47$ MeV, and $U_{\Sigma}^{(\Sigma)}(n_{\text{sat}}/5) = 26$ MeV, whereas the data on double-$\Lambda$-hypernuclei suggest a weakly attractive potential at least for $\Lambda$ hyperons, $U_{\Lambda}^{(\Lambda)}(n_{\text{sat}}/5) \approx -1$ to $-5$ MeV [73,75,76]. Although, as shown, e.g., in Refs. [77,78], the $\sigma^*$ has only a weak influence on the EoS and (proto)neutron star properties, we include a second version (named “DD2Y$\sigma^*$” hereafter) of the EoS with a $\sigma^*$ coupling adjusted to have $U_{\Lambda}^{(\Lambda)}(n_{\text{sat}}/5) = -0.4$ MeV, $U_{\Xi}^{(\Xi)}(n_{\text{sat}}/5) = -0.4$ MeV, and $U_{\Sigma}^{(\Sigma)}(n_{\text{sat}}/5) = -0.4$ MeV. Table I summarizes the values of the scalar meson hyperon couplings in both models obtained from the above described procedure. Note that the couplings to $\Lambda$ in model DD2Y are the same as in the BHB$\Lambda\phi$ EoS [40], where exactly the same procedure has been followed.

<table><caption>TABLE I. Coupling constants of the scalar mesons to different hyperons within the two models presented here, normalized to the $\sigma N$ coupling from the DD2 parameter set, i.e., $R_{Mj} = g_{Mj}/g_{\sigma N}$.</caption><thead><tr><th>Model</th><th>R<sub>σΛ</sub></th><th>R<sub>σ*Λ</sub></th><th>R<sub>σΣ</sub></th><th>R<sub>σ*Σ</sub></th><th>R<sub>σΞ</sub></th><th>R<sub>σ*Ξ</sub></th></tr></thead><tbody><tr><td>DD2Y</td><td>0.62</td><td>0</td><td>0.48</td><td>0</td><td>0.32</td><td>0</td></tr><tr><td>DD2Yσ*</td><td>0.62</td><td>0.46</td><td>0.48</td><td>0.84</td><td>0.32</td><td>1.11</td></tr></tbody></table>

### C. Combining different parts of the EoS

The HS(DD2) EOS contains the transition from inhomogeneous or clusterized matter to uniform nucleonic matter. This is done via the excluded volume mechanism, which suppresses nuclei around and above nuclear saturation density. On top of that, for some thermodynamic conditions a Maxwell construction over a small range in density is necessary; for details see Ref. [10].

Here the situation is slightly more complicated, since homogeneous matter might contain hyperons. In the simplest case, hyperons appear within homogeneous (nucleonic) matter and it is sufficient to minimize the free energy of the homogeneous system to decide upon the particle content of matter. Such a situation occurs at low temperatures and high densities.

In some parts of the $T-n_B$ diagram, however, a transition from inhomogeneous matter directly to hyperonic homogeneous matter is observed. This is the case at low densities and high temperatures, i.e., the density regions up to the bumps in Fig. 1. There, light clusters compete with hyperonic degrees of freedom with only very small differences in free energy which are of the order of the numerical accuracy of the EoS calculation. To technically construct the transition in this region, we follow a similar prescription as in Ref. [40] and introduce a threshold value for the total hyperon fraction, $Y_{\text{hyperons}} = \sum_{j \in B_y} n_j / n_B$. We let hyperonic matter appear only if $Y_{\text{hyperons}} > 10^{-6}$. Note that the hyperon fraction is not the same as the strangeness fraction $Y_S$, defined as the sum of all particle fractions multiplied by their respective strangeness quantum numbers, $Y_S = \sum_{j \in B} S_j n_j / n_B$.

Although the above described procedure allows us to construct a smooth transition between the different parts of the EoS, it is of course not completely consistent. In principle, whenever hyperons compete with light nuclear clusters, the free energy of the system should be minimized allowing simultaneously for all different possibilities, e.g., a coexistence of light clusters with hyperons. In view of the tiny differences in free energy and the small fractions of particles other than nucleons, electrons, and photons in the transition region, a completely consistent treatment is left for future work.

## III. EQUATION OF STATE PROPERTIES

### A. Compatibility with constraints

The interaction between nucleons can be constrained by data of finite nuclei and nuclear matter properties. The latter are chosen in general as the coefficients of a Taylor expansion of the energy per baryon of isospin symmetric nuclear matter around saturation. Values with a reasonable
<figure><img src="image_2.png" /><figcaption>FIG. 1. The lines delimit the regions in temperature and baryon number density for which the overall hyperon fraction exceeds $10^{-4}$, which are situated above the lines. The dark thick purple line corresponds to the $BHB\Lambda\phi$ model and light thin red line to the DD2Y model. Different charge fractions are shown as indicated within the panels.</figcaption></figure>

precision can be obtained for the saturation density ($n_{sat}$), binding energy ($E_B$), incompressibility ($K$), symmetry energy ($E_{sym}$), and its slope ($L$). In addition, much effort has been recently devoted to theoretical *ab initio* calculations of pure neutron matter in order to constrain the equation of state. This is particularly interesting for the EoS of compact stars, completing the information about symmetric matter. The only robust constraint on the interactions at supersaturation density arises from the recent observation of two massive neutron stars, indicating that the maximum mass of a cold, non-, or slowly rotating (therefore spherically symmetric) neutron star should be above $2M_{\odot}$. A summary and discussion of some of the most important available constraints can be found, e.g., in Ref. [79].

The present parametrization, DD2, has been chosen since it agrees well with most of the established constraints. The values for $n_{sat} = 0.149 \text{ fm}^{-3}$, $E_B = 16.0$ MeV, and $K = 243$ MeV are within standard ranges [79]. The compatibility of $E_{sym}$ and $L$ with ranges derived in Ref. [80] (light gray rectangle) and in Ref. [79] (dark gray rectangle), respectively, are shown in Fig. 2. For comparison we show the values for two other interactions, that of the Lattimer and Swesty EoS (LS) [81] and that for the TM1 parametrization [82], too. These two interactions have been employed in other recently developed general purpose EoSs, including non-nucleonic degrees of freedom, e.g., [6, 37, 38].

In Fig. 3 pressure and energy per baryon for pure neutron matter are shown below saturation density. The blue band represents the results from the *ab initio* calculations from Ref. [83] including an estimate of the corresponding uncertainties. In contrast to LS and TM1, the interaction DD2 employed here is in reasonable agreement with the *ab initio* calculations.

The mass-radius relation of cold<fnref n="4" /> spherically symmetric neutron stars within different general purpose EoS models is displayed in Fig. 4. Purely nucleonic versions are shown with solid lines, models including $\Lambda$ hyperons with dotted lines,

and those including the entire baryon octet with dashed-dotted lines. These are the LS EoS [81], its extension with $\Lambda$ hyperons ("LS220$\Lambda$") [11], the EoS by Shen *et al.* ("STOS") employing the TM1 interaction [84], its extension with $\Lambda$ hyperons ("STOS$\Lambda$") [37] and all hyperons ("STOSY") [6], as well as the two models including $\Lambda$ hyperons within the same nuclear model as the present one from Ref. [40], ("BHB$\Lambda$") and ("BHB$\Lambda\phi$"). It is evident from the figure that there are only two EoSs including hyperons compatible with the $2M_{\odot}$ constraint: $BHB\Lambda\phi$ containing only $\Lambda$ hyperons and the present DD2Y. Both models are the same, except for the particle content. The additional hyperonic degrees of freedom in DD2Y slightly reduce the maximum mass with respect to $BHB\Lambda\phi$, but it remains above $2M_{\odot}$. The additional attractive $YY$ interaction in $DD2Y\sigma^*$ reduces the maximum mass to $1.87M_{\odot}$, thus slightly below the observational limit. A summary of cold neutron star properties for the different EoSs is given in Table II.

<figure><img src="image_3.png" /><figcaption>FIG. 2. Values of $E_{sym}$ and $L$ in different nuclear interaction models. The two gray rectangles correspond to the range for $E_{sym}$ and $L$ derived in Ref. [80] (light gray) and Ref. [79] (dark gray) from nuclear experiments and some neutron star observations.</figcaption></figure>

<footnote n="4">⁴For convenience we have chosen a temperature of $T = 0.1$ MeV for producing this figure. In the following discussion of our results we always refer to this temperature upon speaking about “cold” stars.</footnote>
<figure><img src="image_4.png" /><figcaption>FIG. 3. Pressure (left panel) and energy per baryon (right panel) of pure neutron matter as functions of baryon number density within different nuclear interaction models compared with the *ab initio* calculations of Ref. [83], indicated by the blue band.</figcaption></figure>

## B. Hyperon content and thermodynamic properties

As already mentioned in Ref. [78], the overall hyperon content within the EoS remains similar between the models containing only $\Lambda$ hyperons and the corresponding ones with the full baryonic octet. For cold NSs, this can be seen from Table II. In Fig. 1, the regions where the overall hyperon fraction exceeds $10^{-4}$ are compared for $BHB\Lambda$ and DD2Y. Although, as expected, hyperons are slightly more abundant in the full model, the shape of the regions remains the same and only small quantitative differences are observed. The bump in the curves, i.e., the part of the lines above approximately 20 MeV, where the abundance of hyperons is still below $10^{-4}$, arises from the competition between light nuclear clusters and hyperons in this particular temperature and density domain and does not exist in the EoSs built on nuclear models without light clusters; see Ref. [79].

<figure><img src="image_5.png" /><figcaption>FIG. 4. Gravitational mass vs circumferential equatorial radius for cold spherically symmetric neutron stars within different EoS models. The two horizontal bars indicate the two recent precise NS mass determinations, PSR J1614-2230 [41,43] (hatched gray) and PSR J0348+0432 [42] (green).</figcaption></figure>

In Fig. 5, the overall strangeness fraction is shown as function of baryon number density for different values of fixed temperature and electron fraction for both models, $BHB\Lambda\phi$ and DD2Y. As mentioned before, the hyperon onset density remains similar in both models and the decrease in $\Lambda$ fraction in DD2Y with respect to $BHB\Lambda\phi$ at high densities is compensated by the presence of other hyperons such that the overall strangeness fraction is larger in DD2Y. Note that here the strangeness fraction $Y_S$ has been taken and not the hyperon fraction. Naturally, the difference between both models increases with increasing temperature. With increasing $Y_e$, as expected, in both models the overall strangeness fraction decreases. The effect is, however, less pronounced in DD2Y since the population of neutral cascades and $\Sigma^+$ compensates partially the suppression of other hyperonic degrees of freedom.

Pressure and free energy per baryon are considerably reduced above roughly 2-3 times nuclear saturation density in the models with hyperons compared with the purely nucleonic HS(DD2) EoS; see Figs. 6 and 7. It is not surprising that the

<table><caption>TABLE II. Properties of cold spherically symmetric neutron stars in neutrinoless $\beta$ equilibrium: Maximum gravitational and baryonic masses, respectively, radius at a fiducial mass of $M_g = 1.4M_\odot$, the total strangeness fraction $f_S$, representing the integral of the strangeness fraction $Y_S/3$ over the whole star, defined as in Ref. [70], and the central baryon number density. The latter two quantities are given for the maximum mass configuration. In addition to the EoSs presented here, for comparison the values for the purely nucleonic version HS(DD2) [85] and the two versions including only $\Lambda$ hyperons from Ref. [40] are listed.</caption><thead><tr><th>Model</th><th><i>M</i><sub>g</sub><sup>max</sup><br>(<i>M</i><sub>☉</sub>)</th><th><i>M</i><sub>B</sub><sup>max</sup><br>(<i>M</i><sub>☉</sub>)</th><th><i>R</i><sub>1.4</sub><br>(km)</th><th><i>f</i><sub>S</sub></th><th><i>n</i><sub>B</sub><sup>(c)</sup><br>(fm<sup>-3</sup>)</th></tr></thead><tbody><tr><td>HS(DD2)</td><td>2.43</td><td>2.90</td><td>13.27</td><td></td><td>0.84</td></tr><tr><td>BHB<sub>Λ</sub></td><td>1.96</td><td>2.26</td><td>13.27</td><td>0.05</td><td>0.95</td></tr><tr><td>BHB<sub>Λ</sub><i>φ</i></td><td>2.11</td><td>2.47</td><td>13.27</td><td>0.05</td><td>0.96</td></tr><tr><td>DD2Y</td><td>2.04</td><td>2.36</td><td>13.27</td><td>0.04</td><td>1.00</td></tr><tr><td>DD2Y<i>σ</i><sup>*</sup></td><td>1.87</td><td>2.15</td><td>13.27</td><td>0.04</td><td>0.98</td></tr></tbody></table>
<figure><img src="image_6.png" /><figcaption>FIG. 5. Total strangeness fraction as function of baryon number density for different values of fixed temperature and electron fraction $Y_e$ within $BHB\Lambda\phi$ (dotted lines) and DD2Y (solid lines) EoS. In (a)–(c), $Y_e = 0.3$ and in (d)–(f), $T = 30$ MeV. For information, the $\Lambda$ fraction in model DD2Y is indicated, too (dash-dotted lines).</figcaption></figure>

reduction is most important for high temperatures and low electron fractions. The presence of the full baryon octet in DD2Y leads only to a small further reduction with respect to the model $BHB\Lambda\phi$, containing only $\Lambda$ hyperons. This is due to the fact that the overall hyperon fraction is very similar in both models; see the discussion above.

## IV. TEMPERATURE DEPENDENT STELLAR STRUCTURE

In this section, we describe our strategy to solve for the star’s structure, following Ref. [45]. Equilibrium equations will be solved together with Einstein equations, assuming stationarity and axisymmetry. In addition, the matter content (represented by the energy momentum tensor) should fulfill the circularity condition, i.e., the absence of meridional convective currents. An EoS will close the system of equations. In full generality the EoS depends on temperature and on the different particle number densities or thermodynamically equivalent variables. Conditions for electromagnetic and strong equilibrium reduce the number of degrees of freedom in the EoS to 3, related to baryon number density $n_B$, electron number density $n_e$, and temperature $T$. In neutron stars older than several minutes, the temperature can be considered as vanishing and neutrinoless weak $\beta$ equilibrium is achieved, such that the EoS becomes effectively barotropic, i.e., depends only on baryon number density or a thermodynamically equivalent variable. Neither in PNSs nor in merger remnants are these conditions fulfilled and in particular a nonzero temperature has to be considered.

Here, we will allow for an EoS with an explicit temperature dependence. Under the current assumptions, in particular stationarity, the most general solution for the star’s structure becomes again barotropic, and a relation $T(n_B)$ (or thermodynamically equivalent) has to be provided [45,52,53]. For simplicity we will restrict the results within the present work either to neutrinoless $\beta$ equilibrium or to a constant lepton fraction $Y_L = (n_e + n_\nu)/n_B = n_L/n_B$ ($n_\nu$ and $n_L$ being, respectively, neutrino and lepton number densities). Following the standard presentation of the formalism (see, e.g., [45,86]), Latin letters $i, j, \dots$ are used for spatial indices only, whereas Greek letters $\alpha, \beta, \dots$ denote space-time indices.

### A. Einstein equations

General-relativistic models shall be described within the 3+1 formulation, where space-time is foliated by a family of spacelike hypersurfaces $\Sigma_t$, labeled by the time coordinate $t$. Introducing coordinates $(x^i)$ on each hypersurface, the line element can be written as

$$ds^2 = -N^2 dt^2 + \gamma_{ij}(dx^i + \beta^i dt)(dx^j + \beta^j dt). \quad (18)$$

$N$ represents the lapse function, $\beta^i$ the shift vector, and $\gamma_{ij}$ the three-metric on each hypersurface $\Sigma_t$, thus defining the space-time metric $g_{\alpha\beta}$. More details can be found, e.g., in [87].

The assumptions of stationarity, axisymmetry, and asymptotic flatness imply the existence of two commuting Killing vector fields, given as $\zeta = \partial/\partial t$ and $\vec{\chi} = \partial/\partial\varphi$ in an adapted coordinate system $(t, x^1, x^2, \varphi)$. The two remaining coordinates
<figure><img src="image_7.png" /><figcaption>FIG. 6. Pressure (d)–(f) and normalized free energy per baryon (a)–(c) as function of baryon number density for different values of fixed temperature and $Y_e = 0.3$ within the three different EoS, HS(DD2), BHB$\Lambda\phi$, and DD2Y. For information, the pressure in the classical models LS220 and STOS is displayed, too.</figcaption></figure>

are chosen to be spherical, i.e., $x^1 = r, x^2 = \theta$. These adapted coordinates simplify the expression of the metric: $\beta^r = \beta^\theta = 0$ and $\gamma_{r\varphi} = \gamma_{\theta\varphi} = 0$. Finally, following [45,86] we use a quasi-isotropic gauge, which additionally gives $\gamma_{r\theta} = 0$, such that the line element (18) becomes

$$ds^2 = -N^2 dt^2 + A^2(dr^2 + r^2 d\theta^2) + B^2 r^2 \sin^2\theta(d\varphi^2 + \beta^\varphi dt)^2, \quad (19)$$

with the notations $A^2 = \gamma_{rr} = \gamma_{\theta\theta}/r^2$ and $B^2 = \gamma_{\varphi\varphi}/(r^2 \sin^2\theta)$. All the metric potentials ($N, \beta^\varphi, A, B$) are functions of the coordinates $(r, \theta)$ only. Einstein equations for these four gravitational potentials, under our symmetry assumptions, reduce to a set of four elliptic (Poisson-like) partial differential equations, in which source terms contain both contributions from the energy-momentum tensor (matter) and nonlinear terms with noncompact support, involving the gravitational field itself. Explicit expressions and discussion of these equations can be found in Ref. [86].

## B. Equilibrium equations

Matter is described as a perfect fluid with an energy-momentum tensor of the form

$$T^{\alpha\beta} = (\varepsilon + p)u^{\alpha}u^{\beta} + p g^{\alpha\beta}. \quad (20)$$

$\varepsilon$ denotes here the total energy density (including rest mass), $p$ the pressure, and $u^\alpha$ is the fluid four-velocity; the fluid angular velocity is then defined as $\Omega := u^\varphi/u^t$. We also introduce the

pseudolog enthalpy

$$H = \ln \left( \frac{\varepsilon + p}{m_B n_B} \right), \quad (21)$$

with $m_B$ a constant mass, where we chose the value $m_B = 939.565$ MeV. Conservation of the energy-momentum tensor$^5 \nabla_\alpha T^{\alpha\beta} = 0$ yields the equation for the fluid equilibrium [52,53],

$$\partial_i (H + \ln N - \ln \Gamma) = \frac{\mathrm{Te}^{-H}}{m_B} \partial_i s_B - u_\varphi u^t \partial_i \Omega. \quad (22)$$

$\Gamma = \mathrm{Nu}^t$ represents the Lorentz factor of the fluid with respect to the Eulerian observer and $s_B$ the entropy per baryon in units of the Boltzmann constant.

In this work, we will restrict ourselves only to the case where matter is rigidly rotating ($\Omega = \text{const}$), which means that the last term in Eq. (22) is zero. This equation is then integrable in three cases. First, for a constant $s_B$, which is in particular the case at zero temperature. The second case is the isothermal one (constant $T^* = TN/\Gamma$) defined in Ref. [52]. Finally, the most general solution in rigid rotation is found introducing the heat function [45],

$$\hat{H}(n_B) = \int_0^{n_B} \frac{dp}{dn} \frac{1}{\varepsilon(n) + p(n)} dn, \quad (23)$$

<footnote>$^5 \nabla_\alpha$ denotes here the covariant derivative associated to the four-metric $g_{\alpha\beta}$.</footnote>
<figure><img src="image_8.png" /><figcaption>FIG. 7. Same as Fig. 6, but for $T = 30$ MeV and different fixed values of $Y_e$.</figcaption></figure>

where a parametrization $T(n_B)$ has been assumed such that the EoS effectively is again barotropic. Using $\hat{H}$, the equilibrium condition reduces to

$$ \hat{H} + \ln N - \ln \Gamma = \text{const.}, \quad (24) $$

which is pretty similar to the zero-temperature case [86]. It is obvious that the heat function $\hat{H}$ reduces to the pseudolog enthalpy $H$ at zero temperature, up to a constant factor of $\ln[m_B/\mu_B(n_B = 0)]$, which can be absorbed in the right-hand side of Eq. (24). For differentially rotating stars, allowing the rotation law to depend on the entropy profile, in principle, the condition of the EoS being barotropic could be relaxed. Such a scheme is, however, beyond the purpose of the present paper.

We have implemented the above described scheme, with a temperature dependent EoS within the numerical library LORENE [56]; see also Refs. [45,86]. The resolution of elliptic-type partial differential equations (Einstein equations in our case) is based on multidomain spectral methods [88] and is widely used for the computation of stationary rotating compact objects. The equilibrium condition (24) is integrated in a straightforward way and the heat function (23) computed using the trapezoidal rule. Finally, we can use either an analytic (polytropic type) EoS or a tabulated realistic one, which is interpolated in a thermodynamically consistent way using the scheme by Swesty [89].

Input parameters for a rigidly rotating neutron star model are a temperature vs density profile, a prescription for the lepton fraction (either $\beta$ equilibrium or constant $Y_L$), an EoS, a central value for the heat function $\hat{H}(r = 0)$, and a value for the rotation frequency $\Omega$. We can then compute the numerical solution of all field equations described above

and deduce global quantities such as gravitational mass $M_g$ (from the asymptotic behavior of the gravitational potential $N$), angular momentum $J$ (from the asymptotic behavior of the gravitational potential $\beta^\psi$), or circumferential equatorial radius [from the integration of the line element (19) along the star's equator]. More details about these calculations can be found in Ref. [86].

## V. MODELS OF HOT STARS

Within this section we will discuss results for both nonrotating (maximal masses, Sec. V A) and rotating ($I-Q$ relations, Sec. V B) stars with nonzero temperatures, employing different microscopic EoSs exposed in the preceding sections. For the study of their properties, $Y_Q$ will be fixed either by the condition of $\beta$ equilibrium and assuming that neutrinos freely leave the system, i.e., a vanishing electron lepton number chemical potential

$$ \mu_L = 0, \quad (25) $$

or by fixing the electron lepton fraction $Y_L = 0.4$. This value lies slightly above typical values obtained from simulations; see, e.g., [17]. We have chosen it in order to maximize the differences to the $\beta$-equilibrated case and thus show the maximal effect we would expect from composition. Muons will not be considered, although they might have a non-negligible influence on the EoS at the very center of the PNS [38]. Neither of these conditions might be very realistic, since the hydrodynamic evolution should be coupled to neutrino transport, fixing the corresponding evolution of $Y_e = Y_Q$ inside the star. A more complete study of PNS evolution, combining our models with neutrino transport, is left for future work.
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /><figcaption>FIG. 8. Left: $T(n_B)$ relations used for the computation of hot star models. For comparison, the profiles obtained for different values of constant $s_B$ within the DD2Y EoS are shown, too. Right: $s_B(n_B)$ resulting from the two chosen $T(n_B)$ relations within the DD2Y EoS.</figcaption></figure>

Hence, we certainly do not pretend to give a completely realistic picture of a PNS or a merger remnant. This simplified setup is nevertheless sufficient for the purpose of the present study, namely to demonstrate the usability of the newly developed EoS within a numerical code, and to get some ideas about the influence of hyperons on the properties of hot stars. Results with different temperature profiles will be presented: either yielding constant values of entropy per baryon $s_B$ or profiles shown in Fig. 8 (left panel), inspired by realistic calculations of PNS evolution. Profile $T_1$ is within the range of values from the “canonical” simulations of cooling PNSs [3]. The maximum temperatures of profile $T_2$ is slightly above typical values from the aforementioned simulations and corresponds to values reached for a PNS formed in the collapse of a very massive progenitor star, close to the eventual collapse to a black hole; see, e.g., Fig. 16 of Ref. [10]. Analytic expressions for both temperature profiles are given in Appendix B. The corresponding entropy profiles are displayed in the right panel of Fig. 8.

A. Maximal masses of nonrotating hot stars

The maximum baryonic mass a hot star can support is interesting both for the merger remnant of a binary coalescence, and for the PNS after the bounce occurs in a core collapse event, in order to determine the conditions for the formation of a black hole. Different mechanisms were evoked for stabilizing these objects against collapse to a black hole. First, these objects are supposed to be rotating and, as rotational effects on the maximum mass have been examined elsewhere, see, e.g., Refs. [52,53,90], we do not discuss them here. In addition, for the merger remnant and the PNS in the case of collapse of fast spinning progenitor stars, the rotation profile is strongly differential. Although it is not clear what are the time scales driving toward rigid rotation, strong differential rotation can help in supporting very massive configurations [90–92].

Next, in PNSs, the lepton rich environment certainly contributes to support a higher mass [3,16] and it is not only the cooling, but also the deleptonization via neutrino emission of the star which causes a potential collapse to a

black hole. In a merger remnant, which is supposed to be close to $\beta$ equilibrium, this mechanism cannot play the same role. Finally, canonical calculations suggested that thermal pressure is unlikely to be able to stabilize the star [3,16,93]; it might even slightly reduce the maximum mass due to the population of additional degrees of freedom at finite temperature. However, these studies were restricted to rather low entropy values. For PNSs formed in core collapse of massive progenitors, which eventually are expected to collapse to a black hole, it was found in Refs. [10,12] that thermal effects can increase the maximal gravitational mass by up to $0.6M_\odot$, where neutrinoless $\beta$ equilibrium and a constant entropy per baryon of $s_B = 4$ was considered.

When studying the maximum mass, previous works were considering cold stars [92,94], or a very restricted set of EoS, containing only homogeneous matter [3,16] or only nucleonic matter [10,12,90,93]. Our new EoS including hyperonic degrees of freedom allows us to check the influence of these new degrees of freedom on the mass, treating consistently nuclear clustering at low densities and temperatures. A recent study of PNSs with EoSs containing antikaons can be found in Ref. [95].

In order to discuss maximum masses, we have to consider the stability of the computed stellar configurations. At zero temperature for nonrotating stars, stable configurations verify simply $dM_g/dn_B^{(c)} \ge 0$. This criterion is a special case of the result by Friedman et al. [96], who have established a turning point criterion for determining whether a rotating star becomes secularly unstable with respect to axisymmetric perturbations. Here, we have to use an extended version for hot (nonrotating) stars.

For the following considerations, we will consider a more general case for hot rotating stars and we will denote by $J$ the total angular momentum of the star and $S$ the total entropy. $J$ is defined as [86]

$$J = \int A^2 B^2 (E + p) Ur^3 \sin^2 \theta dr d\theta d\phi. \quad (26)$$

$E$ denotes the energy density as measured by a locally nonrotating observer, $E = \Gamma^2(\epsilon + p) - p$, and $U$ the fluid
<table><caption>TABLE III. Maximum gravitational and baryonic masses in units of solar mass for nonrotating stars and different values of constant total entropy. The central temperature of the maximum mass configuration is given, too. Since the entropy per baryon $s_B$ is constant for each configuration, its value for the respective maximum mass configurations can be obtained simply by dividing $S$ by $M_B^{\text{max}}$; see Eq. (27). The upper part assumes neutrinoless $\beta$ equilibrium and in the lower part $Y_L = 0.4$. For sake of an easier comparison the maximum masses in the cold $\beta$-equilibrated case are recalled in the first two columns. No values are given for DD2Y and the $\beta$-equilibrated case at $S = 3M_\odot$ and $S = 5M_\odot$ since at high central densities the electron fraction lies below the limiting value of the table ($Y_e < 0.01$). The corresponding curves in Fig. 9 do not show a maximum; we could thus not determine the maximum masses.</caption><thead><tr><th rowspan="2">Model</th><th colspan="2">S = 0M<sub>☉</sub></th><th colspan="3">S = 3M<sub>☉</sub></th><th colspan="3">S = 5M<sub>☉</sub></th><th colspan="3">S = 7M<sub>☉</sub></th><th colspan="3">S = 9M<sub>☉</sub></th></tr><tr><th>M<sub>g</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>B</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>g</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>B</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>T<sup>(c)</sup><br>(MeV)</th><th>M<sub>g</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>B</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>T<sup>(c)</sup><br>(MeV)</th><th>M<sub>g</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>B</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>T<sup>(c)</sup><br>(MeV)</th><th>M<sub>g</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>M<sub>B</sub><sup>max</sup><br>(M<sub>☉</sub>)</th><th>T<sup>(c)</sup><br>(MeV)</th></tr></thead><tbody><tr><td>HS(DD2)</td><td>2.42</td><td>2.90</td><td>2.43</td><td>2.88</td><td>41</td><td>2.43</td><td>2.83</td><td>68</td><td>2.45</td><td>2.77</td><td>90</td><td>2.50</td><td>2.73</td><td>108</td></tr><tr><td>BHB&Lambda;&phi;</td><td>2.11</td><td>2.47</td><td>2.11</td><td>2.42</td><td>41</td><td>2.13</td><td>2.39</td><td>67</td><td>2.18</td><td>2.37</td><td>91</td><td>2.27</td><td>2.39</td><td>107</td></tr><tr><td>DD2Y</td><td>2.04</td><td>2.35</td><td></td><td></td><td></td><td></td><td></td><td></td><td>2.02</td><td>2.15</td><td>81</td><td>2.11</td><td>2.17</td><td>96</td></tr><tr><td>HS(DD2)</td><td></td><td></td><td>2.37</td><td>2.70</td><td>32</td><td>2.38</td><td>2.67</td><td>53</td><td>2.40</td><td>2.64</td><td>71</td><td>2.44</td><td>2.61</td><td>89</td></tr><tr><td>BHB&Lambda;&phi;</td><td></td><td></td><td>2.17</td><td>2.42</td><td>27</td><td>2.18</td><td>2.39</td><td>49</td><td>2.21</td><td>2.37</td><td>65</td><td>2.27</td><td>2.36</td><td>86</td></tr><tr><td>DD2Y</td><td></td><td></td><td>2.17</td><td>2.43</td><td>22</td><td>2.16</td><td>2.36</td><td>39</td><td>2.16</td><td>2.30</td><td>55</td><td>2.20</td><td>2.27</td><td>70</td></tr></tbody></table>

velocity as measured by the same observer. The latter is related to the factor $\Gamma$ as $\Gamma = (1 - U^2)^{-1/2}$. $S$ can be expressed in a similar way by

$$S = \int A^2 B \Gamma n_B s_B m_B r^2 \sin\theta dr d\theta d\phi. \quad (27)$$

As shown in Ref. [52], based on the work by Sorkin [97], a meaningful criterion for a configuration being secularly stable can be obtained for rigidly rotating stars with a constant $s_B$ (or $T^*$) throughout the star. In the former case, i.e., for constant $s_B$, the total entropy is simply given by $S = s_B M_B$. Following Ref. [52], a star becomes unstable at the extremal points,

$$\left(\frac{\partial J}{\partial n_B^{(c)}}\right)_{M_B, S} = 0, \quad \left(\frac{\partial M_B}{\partial n_B^{(c)}}\right)_{J, S} = 0, \quad \left(\frac{\partial S}{\partial n_B^{(c)}}\right)_{M_B, J} = 0. \qquad (28)$$

Obviously, upon varying the central baryon number density (or equivalently the central heat function) at constant angular momentum, the rotation frequency changes; see, e.g., the textbook [51], i.e., sequences at constant rotation frequency do not allow us to distinguish stable from unstable solutions. Equivalently, for sequences at constant total entropy, the entropy per baryon $s_B$ is not constant, and sequences at constant $s_B$ do not allow us to identify stable and unstable configurations.

We are mainly interested here in thermal effects on the star's mass, corresponding to the maximum mass a cooling star can support, i.e., the second criterion of Eq. (28) is the most interesting one. It determines the maximum mass at different given values of constant total entropy and angular momentum. In the following, we restrict the discussion to nonrotating stars, for which $J = 0$ and is therefore constant. Since the criterion for distinguishing secularly stable from unstable configurations is meaningful only for constant $s_B$ (or $T^*$), we will restrict our investigations of maximum masses to models with constant $s_B$, too.

The different values of the maximum mass are summarized in Table III. In the upper part, neutrinoless $\beta$ equilibrium is assumed, in the lower part a constant $Y_L = 0.4$.

1. Thermal effects

In Fig. 9 we thus display the gravitational mass versus central baryon number density obtained for nonrotating stars and different values of total entropy $S$. Results for three different EoSs are shown: the purely nucleonic one, HS(DD2), the one containing $\Lambda$ hyperons, BHB$\Lambda\phi$, and the new EoS considering the entire baryon octet, DD2Y. We do not show results for DD2Y$\sigma^*$ here since it does not respect the cold neutron star maximum mass constraint. $\beta$ equilibrium is assumed for all calculations.

It is obvious that for $S = 3M_\odot$, corresponding to configurations with $s_B$ roughly between 1 and 2, thermal effects on the maximum mass are small, and almost no difference can be observed with respect to the result for cold stars. A slight reduction of the gravitational mass for BHB$\Lambda\phi$ and, in particular, DD2Y, is due to the population of additional degrees of freedom at finite temperature for these two EoSs allowing for non-nucleonic particles. These findings confirm previous investigations, see, e.g., Refs. [3,16,93]. At low central densities, thermal effects are more important. This can be understood since for total entropy constant, with decreasing gravitational mass, the entropy per baryon $s_B$ of the configurations increases, reaching almost $s_B = 3$ at the lower end of the curves, modifying considerably the EoS.

In contrast, at $S = 9$, thermal effects on the gravitational masses are clearly non-negligible for all three EoSs. The maximum mass is increased by 4% for HS(DD2), 8% for BHB$\Lambda\phi$, and 7% for DD2Y, respectively. These values are of the same order as those expected for rigid rotation [86]. The temperatures and entropies of these configurations are reached typically for PNSs in the postbounce phase of core-collapse events with massive progenitors. The importance of thermal effects can be seen also from the shift in central density of the maximum mass configurations compared with the cold result. The central density is reduced with increasing value of $S$ since the hot star becomes less compact due to thermal excitations.

It should be pointed out that for a given entropy per baryon the temperature is significantly lower within an EoS including hyperons than in a purely nuclear one; see, e.g., [78]. This
<figure><img src="image_11.png" /><figcaption>FIG. 9. Gravitational mass vs central baryon number density for nonrotating stars for three different EoSs, without hyperons (left), with Λ hyperons (middle), and the complete baryon octet (right). Different values of constant total entropy $S$ have been used, indicated in units of solar masses. $\beta$ equilibrium has been assumed. For comparison, the cold result is shown, too. For DD2Y, at $S = 3M_{\odot}$ and $S = 5M_{\odot}$, the curves end at some central density above which no longer any $\beta$-equilibrated solution is found. The reason is that the electron fraction becomes lower than the limiting value of the EoS table ($Y_e = 0.01$). No maximum could be determined in this case.</figcaption></figure>

is a trivial thermodynamic effect: the appearance of hyperons implies that the energy is shared among an increased number of degrees of freedom, with consequently reduced thermal excitations for each of them. Therefore, although the value of $s_B$ for the maximum mass configurations is higher for the EoS with hyperons than for the purely nuclear one, the central temperature with DD2Y is only 96 MeV, whereas it is 108 MeV for HS(DD2).

## 2. Composition

It is known that a lepton rich environment disfavors hyperonic degrees of freedom and that generally with increasing hadronic charge fraction $Y_Q$, the EoS becomes stiffer due to the reduced number of degrees of freedom present [3, 16, 78, 79, 98]. Therefore neutrino trapping has been evoked for a long time already as one of the main mechanisms to stabilize a PNS with hyperons (or pions/kaons) against collapse to a black hole. From Fig. 10 it is evident that our results confirm previous findings. We display the gravitational mass of nonrotating stars as a function of central baryon number density for the three previously considered EoSs. A fixed lepton fraction of $Y_L = 0.4$ and high temperatures ($S = 9$) or low temperatures ($S = 3$) is compared with the respective $\beta$-equilibrated results for cold stars.

As expected, for DD2Y, the lepton rich environment with $Y_L = 0.4$ clearly contributes to increasing considerably the gravitational mass supported by the star. To a lesser extent,

this is true for $BHB\Lambda\phi$, too. The difference between DD2Y and $BHB\Lambda\phi$ becomes small since Λ hyperons, being charge neutral, are less affected by the higher electrons fraction than charged hyperons, essentially $\Sigma^-$. In contrast, for the purely nucleonic EoS HS(DD2) almost no difference between the lepton rich and the $\beta$-equilibrated case is observed at $S = 3$ and only a moderate increase for $S = 9$. The combination of thermal and composition effects leads to a maximum mass of $M_g = 2.2M_{\odot}$ for DD2Y with $S = 9$ and $Y_L = 0.4$, $0.16M_{\odot}$ ($\approx 8\%$) above the cold $\beta$-equilibrated maximum mass.

It should be noted, however, that with increasing temperature the hyperon suppression in a lepton rich environment becomes less pronounced. Therefore, for DD2Y—and to a lesser extent for $BHB\Lambda\phi$, too—the increase in maximum gravitational mass with increasing total entropy is moderate at $Y_L = 0.4$. This can be seen from Fig. 11, too, where the different particle fractions for the maximum mass configurations are shown. For $S = 9$, corresponding to $s_B = 3.96$ with DD2Y, all different hyperonic species have non-negligible fractions at the center of the star due to the high temperatures reached.

## B. I-Q relation

It has been shown [99,100] that there exist relations between the moment of inertia ($I$), the tidal deformability ($\lambda$), and the quadrupole moment ($Q$) of neutron stars which are approximately independent of the internal composition and the EoS. Originally proposed for slowly rotating cold neutron

<figure><img src="image_12.png" /><figcaption>FIG. 10. Same as Fig. 9, comparing stars with $Y_L = 0.4$ at different constant total entropy values with the $\beta$-equilibrated result for cold stars.</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>FIG. 11. Particle fractions in hot neutron star matter for the three different EoSs discussed here. The lepton fraction has been fixed to $Y_L = 0.4$ and the entropy per baryon corresponds to the value of the respective maximum mass configuration with the total entropy indicated in each panel. The vertical lines show the central density of the respective maximum mass configurations.</figcaption></figure>

stars, they remain EoS independent for fast rotation, too, and universal fits with a functional form,

$$ \ln y = a + b \ln x + c (\ln x)^2 + d (\ln x)^3 + e (\ln x)^4, \quad (29) $$

can be established [99,100]. The coefficients a,b,c,d,e are frequency dependent [101] but do not depend on the EoS. x and y represent any couple of the normalized quantities $\bar{I}, \bar{Q}, \bar{\lambda}$,

$$ \bar{I} = \frac{I}{M_g^3}, \quad \bar{Q} = \frac{Q}{M_g^3(J/M_g^2)^2}, \quad \bar{\lambda} = \frac{\lambda}{M_g^5}, \quad (30) $$

with $M_g$ being the star's gravitational mass and J its angular momentum. We will employ here the numerical values for the fit coefficients obtained from a fit to the results for cold stars with the APR EoS [102], the reference EoS in most papers in the literature. They are listed in Table IV.

<table><caption>TABLE IV. Values of the fit parameters in Eq. (29) relating the normalized moment of inertia and quadrupole moment obtained from the results for cold slowly rotating neutron stars with the APR EoS [102].</caption><thead><tr><th>a</th><th>b</th><th>c</th><th>d</th><th>e</th></tr></thead><tbody><tr><td>1.5196</td><td>0.4372</td><td>0.0687</td><td>0.013</td><td>0.000 897</td></tr></tbody></table>

Considering the difficulty of defining Love numbers for the case of a rapidly spinning object (see, e.g., Pani et al. [103]), we will focus here on the $\bar{I}-\bar{Q}$ relation. Nevertheless, a loss of universality in this relation would imply a loss of universality in the more general $\bar{I}-\bar{\lambda}-\bar{Q}$, too. The results for different EoSs are shown in Fig. 12.

Results for cold stars are shown in Fig. 12(a), for slow and fast rotating stars; see the symbols for $s_B = 0$. Different colors represent different EoSs. In addition to the classical nuclear LS and STOS EoSs, we include other general purpose models, not only purely nucleonic but, respectively, with $\Lambda$ hyperons and the entire baryon octet, too, always assuming neutrinoless $\beta$ equilibrium. The present results, considering in addition hyperonic EoSs, clearly confirm previous findings that $I-Q$ relations are independent of the EoS with frequency dependent fit coefficients [101].

In Ref. [54] a study of this relation has been performed, employing purely nuclear EoSs from Refs. [3,17], this time assuming different realistic entropy per baryon and electron fraction profiles for the PNS evolution during the minute following bounce. The main result the authors found was that universality of the so-called $I-Love-Q$ relations is violated in the early phases of PNS evolution and recovered as soon as the entropy gradients smooth out and the star becomes more or less isentropic. It should then be independent of the exact value of $s_B$.
<figure><img src="image_14.png" /><figcaption>FIG. 12. Normalized moment of inertia vs quadrupole moment for different EoSs. The color coding corresponds to different EoS models whereas different symbols indicate rotation frequencies and entropy per baryon of the stars. In panel (a) different nucleonic and hyperonic EoSs are shown for cold stars, in β equilibrium, for the slow as well as fast rotating case, respectively. The LS220 EoS and the counterpart with Λ hyperons are thereby shown by red symbols, STOS EoS with and without hyperons by green symbols, and the three EoSs based on DD2, HS(DD2), BHBAϕ, and DD2Y by violet symbols. For the latter three EoSs in addition the results for slowly rotating configurations with s<sub>B</sub> = 4 are displayed. In the three other panels the cold reference case with the HS(DD2) EoS is displayed by black symbols and different situations are considered with the three EoSs, HS(DD2), BHBAϕ, and DD2Y, indicated by violet symbols: constant Y<sub>L</sub> = 0.4 (b), profile T<sub>2</sub> (c), profile T<sub>1</sub> (d).</figcaption></figure>

Our results including hyperonic EoS confirm that indeed,
the $\bar{I}-\bar{Q}$ relation for an isentropic star with $s_B = 1$ or $s_B = 2$
agrees with the result for cold stars. The same is true for
fast rotation, and assuming $\beta$ equilibrium or a constant
lepton fraction $Y_L$ does only induce a small scatter in the
results. The results with constant $s_B = 4$—see Figs. 12(a) and
12(b)—although they remain universal in the sense that there
is only a small difference between different EoSs, deviate,
however, clearly from the results for cold stars. This can be
seen from Fig. 13, too, where we have plotted for the DD2Y
EoS the relative difference between the numerical results and
the fit function of Eq. (29), $\Delta\bar{I}/\bar{I}_{\text{fit}} = (\bar{I} - \bar{I}_{\text{fit}})/\bar{I}_{\text{fit}}$ at different
values of constant entropy per baryon. For $s_B = 1$ or $s_B = 2$,
the deviations remain below 2%, whereas at $s_B = 4$, they can
exceed 10%.

Both temperature profiles with entropy gradients—see
Figs. 12(c) and 12(d)—display obvious deviations from the
results for cold stars, too. With increasing temperature, the
differences induced by the lepton fraction increase, too.

In Ref. [54] the observed deviations from universality in the
early stages of PNS evolution were attributed to the presence
of entropy gradients. Our results suggest a slightly modified
picture, in the sense that universality is not a question of
entropy gradients, but of thermal effects. As we have seen
also during the preceding discussion on maximum masses,
at $s_B = 1$ or $s_B = 2$, which are typical values in the late
stages of PNS evolution probed in Ref. [54], thermal effects
on the EoS and thus on the star's structure remain small. At
higher entropies, thermal effects start to influence the EoS,
thus the star's structure and universality of $\bar{I}-\bar{Q}$ relations
are modified. Such entropy values can be reached in PNSs
or merger remnants, depending on many factors such as the
progenitors, rotation, or metallicity.

Since the $\bar{I}-\bar{Q}$ relation still seems independent of the
employed EoS, it might be tempting to try to obtain another
“universal” fit, depending this time on rotation frequency
and entropy/temperature. In contrast to the former, neither
temperature nor entropy of the star are quantities which are
<figure><img src="image_15.png" /><figcaption>FIG. 13. Relative difference between Eq. (29)—the fitted results for cold slowly rotating stars—and the results at different constant $s_B$ values for the DD2Y EoS, assuming neutrinoless $\beta$ equilibrium.</figcaption></figure>

observationally accessible. Therefore such a law would not help for data analysis and we refrain from giving one here. Anyway, in view of the present results doubts are allowed concerning the relevance of $\bar{I}-\bar{Q}$ relations for analysis of PNS or merger remnant data, including gravitational wave signals from the last stages of binary neutron star mergers. Let us stress here that entropy values of the order of $s_B = 4$ are quite realistic in such cases; see, e.g., the simulations in Refs. [10,12,104].

## VI. SUMMARY AND CONCLUSIONS

In this paper we have presented a new consistent general purpose EoS, including in particular thermal effects. The new EoS, including the entire baryon octet, is compatible with present constraints from nuclear physics and neutron star observations. The complete new EoS as function of $T, n_B, Y_e$ will be made publicly available in tabulated form on the COMPOSE database [105]; see Appendix A for details.

We have demonstrated the applicability of the new EoS, investigating maximum masses of hot stars, comparing a purely nuclear EoS with one including $\Lambda$ hyperons and the new one with all hyperons. To that end we have applied a numerical code able to provide stationary models of relativistic rotating stars, including the effect of nonzero temperature. The main motivation for studying hot (rotating) stars is the birth of neutron stars, i.e., the evolution of PNSs, and the neutron star created in the aftermath of a binary neutron star merger. In order to correctly identify the configurations which are secularly stable, we have constructed sequences at different values of constant total entropy $S$ in contrast to many previous works considering constant entropy per baryon $s_B$.

As we have seen, thermal effects and a lepton rich environment can considerably increase the maximally supported mass to a degree depending on the EoS. The lepton rich environment is important in particular if hyperons are present. If the entropy per baryon exceeds roughly $s_B = 2$, thermal effects become important in the EoS, too. Thus for a total entropy roughly above $5M_\odot$ thermal effects on the maximum mass

become noticeable. These high temperatures can be reached in both merger remnants and PNSs depending on the particular conditions. Let us recall again that previous works [90–92] suggest that the main effect stabilizing a merger remnant or a PNS above the maximum mass of its nonrotating cold neutrinoless $\beta$-equilibrated counterpart is differential rotation, which we did not consider here.

Following the work by Martinon *et al.* [54], the universality of $I-Q$ relations has been tested for fast rotating hot stars, retrieving their results that a low constant nonzero entropy does not modify the relations. Universality, tested before only for purely nuclear models, is maintained in the presence of hyperons, too. This is, however, no longer true if thermal effects in the EoS become non-negligible, independently of the presence of entropy gradients, i.e., it also occurs for high, but constant entropies.

## ACKNOWLEDGMENTS

We would like to thank D. Gondek and L. Villain for instructive discussions. This work has been partially funded by the SN2NS project ANR-10-BLAN-0503, the “Gravitation et physique fondamentale” action of the Observatoire de Paris, and the COST action MP1304 “NewCompstar.”

## APPENDIX A: TECHNICAL ISSUES OF THE NEW EOS TABLE

The new EoS in its version DD2Y is provided in a tabular form in the COMPOSE data base, http://compose.obspm.fr as a function of $T, n_B, Y_e$. The contribution from electrons is included. Note that the COMPOSE software allows us to calculate additional quantities, such as, e.g., sound speed, from those provided in the tables. Please see the COMPOSE manual [105] and the data sheet on the web site for more details about the definition of the different quantities:

(i) The grid is specified in Table V.

(ii) Thermodynamic quantities provided:

(1) pressure divided by baryon number density $p/n_B$ (MeV),

(2) entropy per baryon $s/n_B$,

(3) scaled baryon chemical potential $\mu_B/m_n - 1$,

(4) scaled charge chemical potential $\mu_Q/m_n$,

(5) scaled (electron) lepton chemical potential $\mu_L/m_n$,

(6) scaled free energy per baryon $f/(n_B m_n) - 1$,

(7) scaled energy per baryon $e/(n_B m_n) - 1$;

(iii) Compositional data provided:

(1) particle fractions of baryons and electrons, $Y_i = n_i/n_B$,

(2) particle fractions of deutons ($^2$H), tritons ($^3$H), $^3$He, and $\alpha$ particles ($^4$He),

(3) fraction of a representative (average) heavy nucleus, together with its average mass number and average charge.

Please note that only nonzero particle fractions are listed.
<table><caption>TABLE V. EoS tables are provided using the above grid for the thermodynamic parameters.</caption><thead><tr><th></th><th><i>T</i></th><th><i>n</i><sub>B</sub></th><th><i>Y</i><sub>e</sub></th></tr></thead><tbody><tr><td>No. of points</td><td>80</td><td>302</td><td>59</td></tr><tr><td>Minimum value</td><td>0.1 MeV</td><td>10<sup>-12</sup> fm<sup>-3</sup></td><td>0.01</td></tr><tr><td>Maximum value</td><td>158.5 MeV</td><td>1.202 fm<sup>-3</sup></td><td>0.6</td></tr><tr><td>Scaling</td><td>logarithmic</td><td>logarithmic</td><td>linear</td></tr></tbody></table>

(iv) Effective Dirac masses $M^*$ of all baryons with nonzero density are provided within homogeneous matter.

## APPENDIX B: EXPRESSIONS FOR THE TEMPERATURE PROFILES

Although they are inspired by results from simulations, for computational simplicity, analytic parametrizations for the

<table><caption>TABLE VI. Parameter values for the two analytic temperature profiles, see Eq. (B1).</caption><thead><tr><th></th><th><i>a</i><br>(MeV fm<sup>3</sup>)</th><th><i>b</i><br>fm<sup>6</sup></th><th><i>c</i><br>(MeV fm<sup>3α</sup>)</th><th><i>d</i><br>(MeV fm<sup>3</sup>)</th><th><i>α</i></th></tr></thead><tbody><tr><td><i>T</i><sub>1</sub></td><td>10.01</td><td>26.21</td><td>77.39</td><td>-65.15</td><td>0.35</td></tr><tr><td><i>T</i><sub>2</sub></td><td>470.0</td><td>26.21</td><td>77.39</td><td>-65.15</td><td><sup>1</sup>&frasl;<sub>3</sub></td></tr></tbody></table>

temperature profiles, $T_1$ and $T_2$, are employed of the form

$$T(n_B) = c n_B^\alpha + d n_B + \frac{a n_B}{1 + \exp[b(n_B - n_0)^2]} \quad (B1)$$

$n_0$ indicates here the saturation density, $n_0 = 0.155 \text{ fm}^{-3}$, and the values of the other parameters are listed in Table VI.

[1] A. Burrows and J. M. Lattimer, Astrophys. J. 307, 178 (1986).

[2] W. Keil and H. T. Janka, Astron. Astrophys. 296, 145 (1995).

[3] J. A. Pons, S. Reddy, M. Prakash, J. M. Lattimer, and J. A. Miralles, Astrophys. J. 513, 780 (1999).

[4] T. W. Baumgarte, S. A. Teukolsky, S. L. Shapiro, H. T. Janka, and W. Keil, Astrophys. J. 468, 823 (1996).

[5] J. A. Pons, A. W. Steiner, M. Prakash, and J. M. Lattimer, Phys. Rev. Lett. 86, 5223 (2001).

[6] C. Ishizuka, A. Ohnishi, K. Tsubakihara, K. Sumiyoshi, and S. Yamada, J. Phys. G 35, 085201 (2008).

[7] K. Sumiyoshi, C. Ishizuka, A. Ohnishi, S. Yamada, and H. Suzuki, Astrophys. J. Lett. 690, L43 (2009).

[8] K. Nakazato, K. Sumiyoshi, and S. Yamada, Astrophys. J. 721, 1284 (2010).

[9] K. Nakazato, S. Furusawa, K. Sumiyoshi, A. Ohnishi, S. Yamada, and H. Suzuki, Astrophys. J. 745, 197 (2012).

[10] M. Hempel, T. Fischer, J. Schaffner-Bielich, and M. Liebendörfer, Astrophys. J. 748, 70 (2012).

[11] B. Peres, M. Oertel, and J. Novak, Phys. Rev. D 87, 043006 (2013).

[12] A. W. Steiner, M. Hempel, and T. Fischer, Astrophys. J. 774, 17 (2013).

[13] P. Char, S. Banik, and D. Bandyopadhyay, Astrophys. J. 809, 116 (2015).

[14] V. Ferrari, G. Miniutti, and J. Pons, Mon. Not. R. Astron. Soc. 342, 629 (2003).

[15] Y. Sekiguchi, K. Kiuchi, K. Kyutoku, and M. Shibata, Phys. Rev. Lett. 107, 051102 (2011).

[16] M. Prakash, I. Bombaci, M. Prakash, P. J. Ellis, J. M. Lattimer, and R. Knorren, Phys. Rep. 280, 1 (1997).

[17] J. A. Pons, J. A. Miralles, M. Prakash, and J. M. Lattimer, Astrophys. J. 553, 382 (2001).

[18] J. A. Pons, S. Reddy, P. J. Ellis, M. Prakash, and J. M. Lattimer, Phys. Rev. C 62, 035803 (2000).

[19] O. E. Nicotra, M. Baldo, G. F. Burgio, and H.-J. Schulze, Phys. Rev. D 74, 123001 (2006).

[20] D. P. Menezes and C. Providencia, arXiv:astro-ph/0703649.

[21] V. Dexheimer and S. Schramm, Astrophys. J. 683, 943 (2008).

[22] N. Yasutake and K. Kashiwa, Phys. Rev. D 79, 043012 (2009).

[23] I. Bombaci, D. Logoteta, C. Providencia, and I. Vidana, Astron. Astrophys. 462, 1017 (2007).

[24] G. F. Burgio, H.-J. Schulze, and A. Li, Phys. Rev. C 83, 025804 (2011).

[25] H. Chen, M. Baldo, G. F. Burgio, and H.-J. Schulze, Phys. Rev. D 86, 045006 (2012).

[26] A. S. Botvina and I. N. Mishustin, Nucl. Phys. A 843, 98 (2010).

[27] S. Typel, G. Röpke, T. Klähn, D. Blaschke, and H. H. Wolter, Phys. Rev. C 81, 015803 (2010).

[28] M. Hempel and J. Schaffner-Bielich, Nucl. Phys. A 837, 210 (2010).

[29] G. Shen, C. J. Horowitz, and S. Teige, Phys. Rev. C 83, 035802 (2011).

[30] G. Shen, C. J. Horowitz, and E. O'Connor, Phys. Rev. C 83, 065808 (2011).

[31] A. R. Raduta and F. Gulminelli, Phys. Rev. C 82, 065801 (2010).

[32] S. Furusawa, K. Sumiyoshi, S. Yamada, and H. Suzuki, Astrophys. J. 772, 95 (2013).

[33] H. Togashi, K. Nakazato, Y. Takehara, S. Yamamuro, H. Suzuki, and M. Takano, Nucl. Phys. A 961, 78 (2017).

[34] S. Furusawa, H. Togashi, H. Nagakura, K. Sumiyoshi, S. Yamada, H. Suzuki, and M. Takano, J. Phys. G 44, 094001 (2017).

[35] K. Nakazato, K. Sumiyoshi, and S. Yamada, Phys. Rev. D 77, 103006 (2008).

[36] I. Sagert, M. Hempel, G. Pagliara, J. Schaffner-Bielich, T. Fischer, A. Mezzacappa, F.-K. Thielemann, and M. Liebendörfer, Phys. Rev. Lett. 102, 081101 (2009).

[37] H. Shen, H. Toki, K. Oyamatsu, and K. Sumiyoshi, Astrophys. J. Suppl. 197, 20 (2011).

[38] M. Oertel, A. F. Fantina, and J. Novak, Phys. Rev. C 85, 055806 (2012).

[39] F. Gulminelli, A. R. Raduta, M. Oertel, and J. Margueron, Phys. Rev. C 87, 055809 (2013).

[40] S. Banik, M. Hempel, and D. Bandyopadhyay, Astrophys. J. Suppl. 214, 22 (2014).

[41] P. Demorest, T. Pennucci, S. Ransom, M. Roberts, and J. Hessels, Nature (London) 467, 1081 (2010).
[42] J. Antoniadis, P. C. C. Freire, N. Wex, T. M. Tauris, R. S. Lynch, M. H. van Kerkwijk, M. Kramer, C. Bassa, V. S. Dhillon, T. Driebe, J. W. T. Hessels, V. M. Kaspi, V. I. Kondratiev, N. Langer, T. R. Marsh, M. A. McLaughlin, T. T. Pennucci, S. M. Ransom, I. H. Stairs, J. van Leeuwen, J. P. W. Verbiest, and D. G. Whelan, *Science* **340**, 448 (2013).

[43] E. Fonseca et al., Astrophys. J. **832**, 167 (2016).

[44] Y. Suwa, Publ. Astron. Soc. Jpn. **66**, L1 (2014).

[45] L. Villain, J. A. Pons, P. Cerda-Duran, and E. Gourgoulhon, Astron. Astrophys. **418**, 283 (2004).

[46] R. C. Tolman, Phys. Rev. **55**, 364 (1939).

[47] J. R. Oppenheimer and G. M. Volkoff, Phys. Rev. **55**, 374 (1939).

[48] J. B. Hartle and K. S. Thorne, Astrophys. J. **153**, 807 (1968).

[49] T. Nozawa, N. Stergioulas, E. Gourgoulhon, and Y. Eriguchi, Astron. Astrophys. Suppl. Ser. **132**, 431 (1998).

[50] M. Ansorg, D. Gondek-Rosinska, and L. Villain, Mon. Not. R. Astron. Soc. **396**, 2359 (2009).

[51] J. L. Friedman and N. Stergioulas, *Rotating Relativistic Stars* (Cambridge University Press, Cambridge, UK, 2013).

[52] J. O. Goussard, P. Haensel, and J. L. Zdunik, Astron. Astrophys. **330**, 1005 (1998).

[53] J.-O. Goussard, P. Haensel, and J. L. Zdunik, Astron. Astrophys. **321**, 822 (1997).

[54] G. Martinon, A. Maselli, L. Gualtieri, and V. Ferrari, Phys. Rev. D **90**, 064026 (2014).

[55] G. Camelio, L. Gualtieri, J. A. Pons, and V. Ferrari, Phys. Rev. D **94**, 024008 (2016).

[56] E. Gourgoulhon, P. Grandclément, J.-A. Marck, J. Novak, and K. Taniguchi, LORENE: Spectral methods differential equations solver, Astrophysics Source Code Library (2016), http://ascl.net/1608.018.

[57] A. G. Pili, N. Bucciantini, A. Drago, G. Pagliara, and L. Del Zanna, Mon. Not. R. Astron. Soc. **462**, L26 (2016).

[58] R. Ouyed, J. Dey, and M. Dey, Astron. Astrophys. **390**, L39 (2002).

[59] R. Ouyed, B. Niebergal, and P. Jaikumar, in *Proceedings of the Compact Stars in the QCD Phase Diagram III (CSQCD III) Conference, December 12–15, 2012* (Guarujá, SP, Brazil, 2012).

[60] M. Buballa, V. Dexheimer, A. Drago, E. Fraga, P. Haensel et al., J. Phys. G **41**, 123001 (2014).

[61] F. Gulminelli and A. R. Raduta, Phys. Rev. C **92**, 055803 (2015).

[62] K. Sumiyoshi and G. Röpke, Phys. Rev. C **77**, 055804 (2008).

[63] N. Buyukcizmeci, A. S. Botvina, and I. N. Mishustin, Astrophys. J. **789**, 33 (2014).

[64] S. Heckel, P. P. Schneider, and A. Sedrakian, Phys. Rev. C **80**, 015805 (2009).

[65] F. Gulminelli, P. Chomaz, A. H. Raduta, and A. R. Raduta, Phys. Rev. Lett. **91**, 202701 (2003).

[66] F. Gulminelli and A. R. Raduta, Phys. Rev. C **85**, 025803 (2012).

[67] G. Audi, A. H. Wapstra, and C. Thibault, Nucl. Phys. A **729**, 337 (2003).

[68] P. Möller, J. R. Nix, W. D. Myers, and W. J. Swiatecki, At. Data Nucl. Data **59**, 185 (1995).

[69] M. Dutra, O. Lourenço, S. S. Avancini, B. V. Carlson, A. Delfino, D. P. Menezes, C. Providência, S. Typel, and S. J. R. Stone, Phys. Rev. C **90**, 055203 (2014).

[70] S. Weissenborn, D. Chatterjee, and J. Schaffner-Bielich, Phys. Rev. C **85**, 065802 (2012).

[71] T. Miyatsu, S. Yamamuro, and K. Nakazato, Astrophys. J. **777**, 4 (2013).

[72] J. Schaffner and I. N. Mishustin, Phys. Rev. C **53**, 1416 (1996).

[73] M. Fortin, S. S. Avancini, C. Providência, and I. Vidaña, Phys. Rev. C **95**, 065803 (2017).

[74] E. N. E. van Dalen, G. Colucci, and A. Sedrakian, Phys. Lett. B **734**, 383 (2014).

[75] I. Vidaña, A. Polls, A. Ramos, and H. J. Schulze, Phys. Rev. C **64**, 044301 (2001).

[76] E. Khan, J. Margueron, F. Gulminelli, and A. R. Raduta, Phys. Rev. C **92**, 044313 (2015).

[77] M. Oertel, C. Providência, F. Gulminelli, and A. R. Raduta, J. Phys. G **42**, 075202 (2015).

[78] M. Oertel, F. Gulminelli, C. Providência, and A. R. Raduta, Eur. Phys. J. A **52**, 50 (2016).

[79] M. Oertel, M. Hempel, T. Klähn, and S. Typel, Rev. Mod. Phys. **89**, 015007 (2017).

[80] J. M. Lattimer and Y. Lim, Astrophys. J. **771**, 51 (2013).

[81] J. M. Lattimer and F. D. Swesty, Nucl. Phys. A **535**, 331 (1991).

[82] Y. Sugahara and H. Toki, Nucl. Phys. A **579**, 557 (1994).

[83] K. Hebeler, J. Lattimer, C. Pethick, and A. Schwenk, Astrophys. J. **773**, 11 (2013).

[84] H. Shen, H. Toki, K. Oyamatsu, and K. Sumiyoshi, Prog. Theor. Phys. **100**, 1013 (1998).

[85] T. Fischer, M. Hempel, I. Sagert, Y. Suwa, and J. Schaffner-Bielich, Eur. Phys. J. A **50**, 46 (2014).

[86] S. Bonazzola, E. Gourgoulhon, M. Salgado, and J. Marck, Astron. Astrophys. **278**, 421 (1993).

[87] E. Gourgoulhon, *3+1 Formalism in General Relativity* (Springer-Verlag, Berlin, 2012).

[88] P. Grandclément and J. Novak, Living Rev. Relativ. **12**, 1 (2009).

[89] F. D. Swesty, J. Comput. Phys. **127**, 118 (1996).

[90] W. Kastaun and F. Galeazzi, Phys. Rev. D **91**, 064027 (2015).

[91] T. W. Baumgarte, S. L. Shapiro, and M. Shibata, Astrophys. J. **528**, L29 (2000).

[92] I. A. Morrison, T. W. Baumgarte, and S. L. Shapiro, Astrophys. J. **610**, 941 (2004).

[93] J. D. Kaplan, C. D. Ott, E. P. O'Connor, K. Kiuchi, L. Roberts, and M. Duez, Astrophys. J. **790**, 19 (2014).

[94] M. Salgado, S. Bonazzola, E. Gourgoulhon, and P. Haensel, Astron. Astrophys. **291**, 155 (1994).

[95] N. D. Batra, K. P. Nunna, and S. Banik, arXiv:1705.06435.

[96] J. L. Friedman, J. R. Ipser, and R. D. Sorkin, Astrophys. J. **325**, 722 (1988).

[97] R. D. Sorkin, Astrophys. J. **257**, 847 (1982).

[98] G. Colucci and A. Sedrakian, Phys. Rev. C **87**, 055806 (2013).

[99] K. Yagi and N. Yunes, Science **341**, 365 (2013).

[100] K. Yagi and N. Yunes, Phys. Rev. D **88**, 023009 (2013).

[101] D. D. Doneva, S. S. Yazadjiev, N. Stergioulas, and K. D. Kokkotas, Astrophys. J. **781**, L6 (2013).

[102] A. Akmal, V. R. Pandharipande, and D. G. Ravenhall, Phys. Rev. C **58**, 1804 (1998).

[103] P. Pani, L. Gualtieri, and V. Ferrari, Phys. Rev. D **92**, 124003 (2015).

[104] A. Perego, S. Rosswog, R. M. Cabezón, O. Korobkin, R. Käppeli, A. Arcones, and M. Liebendörfer, Mon. Not. R. Astron. Soc. **443**, 3134 (2014).

[105] S. Typel, M. Oertel, and T. Klähn, Phys. Part. Nucl. **46**, 633 (2015).