Maximum mass and universal relations of rotating relativistic hybrid hadron-quark stars

TL;DR

This study models rotating hybrid hadron-quark stars with phase transitions, revealing that their maximum mass ratios vary with EOS and confirming universal relations under certain conditions, impacting gravitational wave interpretations.

Contribution

It introduces new equilibrium models of hybrid stars with phase transitions, analyzing their maximum mass limits and universal relations, which differ from purely hadronic models and influence gravitational wave constraints.

Findings

Maximum mass ratio for rotating hybrid stars varies between 1.15 and 1.31.

Universal relations hold for uniform rotation but depend on differential rotation.

GW170817 constraints do not fully apply to hybrid hadron-quark EOSs.

Abstract

We construct equilibrium models of uniformly and differentially rotating hybrid hadron-quark stars using equations of state (EOSs) with a first-order phase transition that gives rise to a third family of compact objects. We find that the ratio of the maximum possible mass of uniformly rotating configurations - the supramassive limit - to the Tolman-Oppenheimer-Volkoff (TOV) limit mass is not EOS-independent, and is between 1.15 and 1.31,in contrast with the value of 1.20 previously found for hadronic EOSs. Therefore, some of the constraints placed on the EOS from the observation of the gravitational wave event GW170817 do not apply to hadron-quark EOSs. However, the supramassive limit mass for the family of EOSs we treat is consistent with limits set by GW170817, strengthening the possibility of interpreting GW170817 with a hybrid hadron-quark EOSs. We also find that along constant angular momentum sequences of uniformly rotating stars, the third family maximum and minimum mass models satisfy approximate EOS-independent relations, and the supramassive limit of the third family is approximately 16.5 % larger than the third family TOV limit. For differentially rotating spheroidal stars, we find that a lower-limit on the maximum supportable rest mass is 123 % more than the TOV limit rest mass. Finally, we verify that the recently discovered universal relations relating angular momentum, rest mass and gravitational mass for turning-point models hold for hybrid hadron-quark EOSs when uniform rotation is considered, but have a clear dependence on the degree of differential rotation.

Figures

Tables

Table 1. Table 1: Labels of the equations of state employed in this paper and their corresponding labels in Paschalidis2018 . EOSs T 1 1 1 –TT belong to Set I and have a first-order phase transition with energy-density jump parametrized by ξ 𝜉 \xi as shown in Equation ( 2 ). Their quark phase is described by Equation ( 1 ). On the other hand, EOSs A 4 4 4 –A 7 7 7 belong to Set II, and have a quark phase parametrized as piecewise polytropes as shown in Equation ( 3 ). We also report the maximum mass that a TOV star can have on the hadronic ( M Max, hadronic TOV superscript subscript 𝑀 Max, hadronic TOV M_{\text{Max, hadronic}}^{\text{TOV}} ) and the hybrid ( M Max, hybrid TOV superscript subscript 𝑀 Max, hybrid TOV M_{\text{Max, hybrid}}^{\text{TOV}} ) branches, respectively.

EOS	EOS in Paschalidis2018	$M_{\text{Max, hadronic}}^{\text{TOV}}$	$M_{\text{Max, hybrid}}^{\text{TOV}}$
		$\left[M_{\odot }\right]$	$\left[M_{\odot }\right]$
T1	ACS-I $\xi =0.10$	$2.00$	$2.47$
T2	ACS-I $\xi =0.27$	$2.00$	$2.31$
T4	ACS-I $\xi =0.43$	$2.00$	$2.17$
T6	ACS-I $\xi =0.60$	$2.00$	$2.05$
T7	ACS-II $\xi =0.70$	$1.50$	$2.14$
T8	ACS-II $\xi =0.80$	$1.50$	$2.08$
T9	ACS-II $\xi =0.90$	$1.50$	$2.02$
TT	ACS-II $\xi =1.00$	$1.50$	$1.98$
A4	ACB4	$2.01$	$2.11$
A5	ACB5	$1.40$	$2.00$
A6	ACB6	$2.01$	$2.00$
A7	ACB7	$1.50$	$2.00$

Table 2. Table 2: Supramassive limit M Max supra subscript superscript 𝑀 supra Max M^{\text{supra}}_{\text{Max}} and maximum mass for uniformly rotating stars belonging to the third family M Max, third supra subscript superscript 𝑀 supra Max, third M^{\text{supra}}_{\text{Max, third}} for the EOSs considered in this paper (Table  1 ). With the exception of T6, A4 and A6, for all the EOSs the most massive star belongs to the third family. For the other three EOSs, the maximum mass occurs in the phase transition (see Figure  4 ).

EOS	$M_{\text{Max}}^{\text{supra}}$	$M_{\text{Max}}^{\text{supra}}$	$M_{\text{Max, third}}^{\text{supra}}$	$M_{\text{Max, third}}^{\text{supra}}$
	$\left[M_{\odot }\right]$	$M_{\text{Max}}^{\text{TOV}}$	$\left[M_{\odot }\right]$	$M_{\text{Max}}^{\text{TOV}}$
T1	$2.39$	$1.20$	$2.39$	$1.20$
T2	$2.37$	$1.19$	$2.37$	$1.19$
T4	$2.36$	$1.18$	$2.36$	$1.18$
T6	$2.46$	$1.23$	$2.34$	$1.17$
T7	$2.32$	$1.16$	$2.32$	$1.16$
T8	$2.31$	$1.16$	$2.31$	$1.16$
T9	$2.31$	$1.15$	$2.31$	$1.15$
TT	$2.32$	$1.15$	$2.32$	$1.15$
A4	$2.56$	$1.27$	$2.34$	$1.17$
A5	$2.30$	$1.15$	$2.30$	$1.15$
A6	$2.63$	$1.31$	$2.33$	$1.16$
A7	$2.32$	$1.16$	$2.32$	$1.16$

Table 3. Table 3: Maximum rest mass spheroidal models for the EOSs considered in this work. Shown are the values of the degree of differential rotation A ^ − 1 superscript ^ 𝐴 1 \hat{A}^{-1} , the maximum energy density ϵ max subscript italic-ϵ max \epsilon_{\rm max} in units of ​ 10 15   g   cm − 3 times E15 times gram centimeter 3 {10}^{15}\text{\,}\mathrm{g}\text{\,}{\mathrm{cm}}^{-3} , the ratio of polar to equatorial radius r p / r e subscript 𝑟 𝑝 subscript 𝑟 𝑒 {r_{p}}/{r_{e}} , the ratio of kinetic to gravitational potential energy T / | W | 𝑇 𝑊 {T}/{\lvert W\rvert} , the ratio of central to equatorial angular velocity Ω c / Ω e subscript Ω 𝑐 subscript Ω 𝑒 {\Omega_{c}}/{\Omega_{e}} , the circumferential radius at the equator R e subscript 𝑅 𝑒 R_{e} in units of km kilometer \mathrm{km} , the dimensionless spin J / M 2 𝐽 superscript 𝑀 2 {J}/{M^{2}} , the mass M ADM subscript 𝑀 ADM M_{\rm ADM} in units of M ⊙ sunmass \mathrm{\text{$M_{\odot}$}} , the ratio of the mass M ADM subscript 𝑀 ADM M_{\rm ADM} to the TOV limit mass M Max TOV superscript subscript 𝑀 Max TOV M_{\text{Max}}^{\text{TOV}} , the rest mass M 0 subscript 𝑀 0 M_{0} in units of M ⊙ sunmass \mathrm{\text{$M_{\odot}$}} and the ratio of the rest mass to TOV limit rest mass M 0 , Max TOV superscript subscript 𝑀 0 Max TOV M_{0,\text{Max}}^{\text{TOV}} . The increase in maximum mass compared to the non-rotating limit is between ∼ 40   % similar-to absent times 40 percent \sim$40\text{\,}\mathrm{\char 37\relax}$ and ∼ 123   % similar-to absent times 123 percent \sim$123\text{\,}\mathrm{\char 37\relax}$ . Quantities of particular interest are T / | W | 𝑇 𝑊 T/\lvert W\rvert that can be larger than ∼ 0.26 similar-to absent 0.26 \sim$0.26$ , which is the threshold for the dynamical bar mode instability (see e.g.  PaschalidisStergioulas2017 for a review), and the rest mass which can be larger than twice the TOV limit – such configurations are classified as ubermassive according to  Espino2019 .

EOS	${\hat{A}}^{-1}$	${ϵ}_{\max }$	$r_p$	$T$	${\Omega }_c$	$R_e$	$J$	$M$	$M$	$M_0$	$M_0$
		$\left[{10}^{15}g{\mathrm{cm}}^{-3}\right]$	$r_e$	$\left|W\right|$	${\Omega }_e$	$\left[\mathrm{km}\right]$	$M^2$	$\left[M_{\odot }\right]$	$M_{\text{Max}}^{\mathrm{TOV}}$	$\left[M_{\odot }\right]$	$M_{0,\text{Max}}^{\mathrm{TOV}}$
T1	$0.390$	$1.03$	$0.35$	$0.265$	$1.930$	$18.25$	$0.922$	$3.860$	$1.563$	$4.526$	$1.542$
T2	$0.410$	$1.09$	$0.32$	$0.276$	$1.965$	$19.35$	$0.943$	$3.974$	$1.720$	$4.662$	$1.699$
T4	$0.415$	$1.23$	$0.31$	$0.280$	$1.990$	$19.53$	$0.948$	$4.008$	$1.847$	$4.699$	$1.842$
T6	$0.425$	$1.38$	$0.30$	$0.283$	$2.038$	$19.58$	$0.952$	$4.018$	$2.004$	$4.705$	$1.971$
T7	$0.415$	$1.71$	$0.36$	$0.247$	$2.064$	$14.70$	$0.892$	$3.097$	$1.447$	$3.633$	$1.393$
T8	$0.395$	$2.07$	$0.38$	$0.238$	$2.040$	$13.55$	$0.878$	$2.944$	$1.415$	$3.450$	$1.368$
T9	$0.420$	$2.08$	$0.38$	$0.235$	$2.079$	$13.49$	$0.875$	$2.845$	$1.408$	$3.325$	$1.361$
TT	$0.430$	$2.22$	$0.38$	$0.234$	$2.116$	$13.11$	$0.872$	$2.758$	$1.393$	$3.216$	$1.356$
A4	$0.405$	$1.11$	$0.26$	$0.304$	$1.985$	$22.69$	$0.985$	$4.688$	$2.222$	$5.593$	$2.275$
A5	$0.510$	$0.89$	$0.28$	$0.293$	$1.990$	$23.07$	$1.013$	$3.951$	$1.976$	$4.577$	$1.937$
A6	$0.410$	$1.31$	$0.25$	$0.305$	$2.036$	$22.00$	$0.983$	$4.557$	$2.279$	$5.431$	$2.334$
A7	$0.490$	$0.96$	$0.28$	$0.292$	$1.992$	$22.44$	$1.000$	$4.013$	$2.007$	$4.675$	$1.987$

Table 4. Table 4: Maximum rest mass differentially rotating models for the T1, T8, and A6 EOSs. Shown are the solution types (A, B, or C), the value of the degree of differential rotation A ^ − 1 superscript ^ 𝐴 1 \hat{A}^{-1} , the maximum energy density ϵ max subscript italic-ϵ max \epsilon_{\rm max} in units of ​ 10 15   g   cm − 3 times E15 times gram centimeter 3 {10}^{15}\text{\,}\mathrm{g}\text{\,}{\mathrm{cm}}^{-3} , the ratio of polar to equatorial radius r p / r e subscript 𝑟 𝑝 subscript 𝑟 𝑒 {r_{p}}/{r_{e}} , the ratio of kinetic to gravitational potential energy T / | W | 𝑇 𝑊 {T}/{\lvert W\rvert} , the ratio of central to equatorial angular velocity Ω c / Ω e subscript Ω 𝑐 subscript Ω 𝑒 {\Omega_{c}}/{\Omega_{e}} , the circumferential radius at the equator R e subscript 𝑅 𝑒 R_{e} in units of km kilometer \mathrm{km} , the dimensionless spin J / M 2 𝐽 superscript 𝑀 2 {J}/{M^{2}} , the ADM mass M 𝑀 M in units of M ⊙ sunmass \mathrm{\text{$M_{\odot}$}} , the ratio of the mass to the TOV limit mass, the rest mass in units of M ⊙ sunmass \mathrm{\text{$M_{\odot}$}} and the ratio of the rest mass to TOV limit rest mass M 0 , Max TOV superscript subscript 𝑀 0 Max TOV M_{0,\text{Max}}^{\text{TOV}} .

EOS	Type	${\hat{A}}^{-1}$	${ϵ}_{\max }$	$r_p$	$T$	${\Omega }_c$	$R_e$	$J$	$M$	$M$	$M_0$	$M_0$
			$\left[{10}^{15}g{\mathrm{cm}}^{-3}\right]$	$r_e$	$\left|W\right|$	${\Omega }_e$	$\left[\mathrm{km}\right]$	$M^2$	$\left[M_{\odot }\right]$	$M_{\text{Max}}^{\mathrm{TOV}}$	$\left[M_{\odot }\right]$	$M_{0,\text{Max}}^{\mathrm{TOV}}$
T1	A	$0.35$	$1.22$	$0.38$	$0.247$	$1.787$	$17.50$	$0.900$	$3.700$	$1.498$	$4.331$	$1.476$
	B	$0.35$	$0.56$	$0.01$	$0.332$	$1.756$	$27.18$	$1.047$	$5.764$	$2.333$	$6.744$	$2.299$
	C	$0.50$	$0.69$	$0.01$	$0.314$	$2.589$	$21.44$	$0.979$	$4.897$	$1.983$	$5.757$	$1.962$
T8	A	$0.39$	$2.14$	$0.39$	$0.236$	$2.037$	$13.27$	$0.874$	$2.920$	$1.404$	$3.421$	$1.357$
	B	$0.35$	$0.56$	$0.01$	$0.332$	$1.756$	$27.18$	$1.047$	$5.764$	$2.771$	$6.744$	$2.675$
	C	$0.55$	$0.66$	$0.02$	$0.309$	$2.646$	$21.70$	$0.983$	$4.728$	$2.273$	$5.519$	$2.189$
A6	A	$0.37$	$1.21$	$0.38$	$0.243$	$1.584$	$20.97$	$0.923$	$3.638$	$1.819$	$4.281$	$1.840$
	B	$0.35$	$0.50$	$0.01$	$0.332$	$1.759$	$27.35$	$1.045$	$5.849$	$2.925$	$6.943$	$2.984$
	C	$0.45$	$0.58$	$0.01$	$0.318$	$2.380$	$22.99$	$0.987$	$5.327$	$2.664$	$6.392$	$2.748$

Equations

$$P(ϵ)=\{\begin{array}{ll}{P}_{\text{tr}}& {ϵ}_1\le ϵ\le {ϵ}_2,\\ P_{\text{tr}}+c_s^2(ϵ-{ϵ}_2)& ϵ\ge {ϵ}_2,\end{array}$$

$${ϵ}_2=(1+\xi ){ϵ}_1.$$

$$P(n)={\kappa }_i{\left(\frac{n}{n_0}\right)}^{{\Gamma }_i}\text{for }{n}_i<n<n_{i+1},$$

$$u=(u^t,0,0,u^{\varphi }),$$

$$\begin{gathered} {ds}^2=-{\mathrm{e}}^{\gamma +\rho }{dt}^2+{\mathrm{e}}^{\gamma -\rho }{r}^2{\sin }^2\theta (d\varphi -\omega dt{)}^2+ \\ {\mathrm{e}}^{2\alpha }({dr}^2+r^2{d\theta }^2), \end{gathered}$$

$$u^t{u}_{\varphi }=A^2({\Omega }_c-\Omega ),$$

$${\hat{A}}^{-1}=\frac{r_e}{A},$$

$$\beta=-\left(\frac{r_{e}}{r_{p}}\right)^{2}\frac{\operatorname{d}\!{(}z^{2})}{\operatorname{d}\!{(}\varpi^{2})}\biggr{\rvert}_{\varpi=r_{e}}$$

$$\hat{\beta }=\frac{\beta }{1+\beta }.$$

$$\frac{\partial M(\epsilon_{c})}{\partial\epsilon_{c}}\biggr{\rvert}_{J}=0\,.$$

$$\frac{M_{\uparrow}}{M_{\text{Max}}^{\text{TOV}}}=1+$0.1316$\left(\frac{J}{J_{\text{Max}}}\right)^{2}+$0.0711$\left(\frac{J}{J_{\text{Max}}}\right)^{4}\,.$$

$$M_{\text{Max}}^{\text{supra}}=\alpha M_{\text{Max}}^{\text{TOV}},$$

$$M_{\text{tot}}\approx$2.74\text{\,}\mathrm{\text{$M_{\odot}$}}$>M^{\text{supra}}_{\text{Max}}=\alpha M^{\text{TOV}}_{\text{Max}}\,.$$

$$M^{\text{TOV}}_{\text{Max}}<\frac{M_{\text{tot}}}{\alpha}\approx$2.28\text{\,}\mathrm{\text{$M_{\odot}$}}$\,.$$

$$$2.07\text{\,}\mathrm{\text{$M_{\odot}$}}$\lesssim M^{\text{TOV}}_{\text{Max}}\lesssim$2.38\text{\,}\mathrm{\text{$M_{\odot}$}}$\,.$$

$$\frac{M_{\uparrow}}{M^{\text{TOV}}_{\text{Max}}}=1+$0.215$\left(\frac{J}{J_{\text{Kep}}}\right)^{2}-$0.050$\left(\frac{J}{J_{\text{Kep}}}\right)^{4}\,.$$

$$M_{\text{Max, third}}^{\text{supra}}=($1.165$\pm$0.035$)\,M^{\text{TOV}}_{\text{Max}}\,,$$

$$\frac{M_{\uparrow }}{M_{\text{Max}}^{\text{TOV}}}$$

$$\frac{M_{0,\uparrow }}{M_{0,\text{Max}}^{\text{TOV}}}$$

$$\frac{M_{\downarrow}}{M^{\text{TOV}}_{\downarrow}}=1+$0.33$\left(\frac{J}{J_{\downarrow,\text{Kep}}}\right)^{2}-$0.10$\left(\frac{J}{J_{\downarrow,\text{Kep}}}\right)^{4}\,.$$

$$\frac{M_{\downarrow }}{M_{\downarrow }^{\text{TOV}}}$$

$$\frac{M_{0,\downarrow }}{M_{0,\downarrow }^{\text{TOV}}}$$

$$M_{\text{Max, third}}^{\text{supra}}=($1.165$\pm$0.035$)\,M^{\text{TOV}}_{\text{Max}}\,.$$

$$\frac{\partial\hat{A}^{-1}}{\partial(r_{p}/r_{e})}\biggr{\rvert}_{\epsilon_{\text{max}}}=0=\frac{\partial\hat{A}^{-1}}{\partial\hat{\beta}}\biggr{\rvert}_{\epsilon_{\text{max}}}\,.$$

$$ϵ(P)={ϵ}_0,$$

Full text

1 Introduction

The properties of the equation of state (EOS) of dense nuclear matter constitute an important open problem that brings together several different fields ranging from astrophysics and gravitation to particle and nuclear physics. While the EOS is understood up to nuclear saturation density $\rho_{0}\approx$2.7\text{\times}{10}^{15}\text{,}\mathrm{g}\text{,}{\mathrm{cm}}^{-3}$$ Haensel1994 ; Haensel2007 ; Lattimer2012 , the situation is less clear at super-nuclear densities, with several different models being possible (see, e.g. Heiselberg2000 ; Lattimer2012 ). Regimes of such high density are expected to occur in Nature in the central regions of compact stars, where the density can reach up to several times $\rho_{0}$.

The EOS determines several of the macroscopic properties of compact stars, such as mass and radius, thus affecting observable quantities (for reviews, see Lattimer2001 ; Hebeler2013 ; Lattimer2016 ; Ozel2016 ). For this reason, astronomical observation is one of the primary tools to probe the high-density end of the EOS (see, e.g., Lattimer2007 ). Recently, the detection of gravitational waves from a binary neutron star Abbott2017a (event GW170817), accompanied by electromagnetic counterparts (events GW170817-GBM, and AT2017gfo Abbott2017a ; Abbott2017b ; Abbott2017c ; Fermi-LAT-Collaboration2017 ; Valenti2017 ; Savchenko2017 ) provided a wealth of information that, under certain assumptions, helped to place constraints on the EOS Bauswein2017b ; Shibata2017 ; Margalit2017 ; Annala2018 ; Ruiz2018 ; Radice2018b ; Abbott2018b ; Raithel2018 ; Most2018 ; Rezzolla2018 ; Burgio2018 ; Koppel2019 (see also Raithel:2019uzi for a review). Among the EOSs that are compatible with GW170817, there are those with a first-order phase transition from hadronic to quark matter at densities larger than $\rho_{0}$, which can lead to the so-called third family of stable stars. The third family consists of objects that are more compact than white dwarfs and neutron stars. These stars are hybrid hadron-quark stars (or just hybrid stars), i.e., configurations with a quark core surrounded by a hadronic shell.111It is also possible for a fourth family to arise Alford2017 if the EOS has an additional phase transition. However, this scenario will not be considered here. The third family has been studied for a long time Gerlach1968 ; Kampfer1981 ; Schertler1997 ; Schertler1998 ; Schertler2000 ; Glendenning2000 ; Ayvazyan2013 ; Blaschke2013 ; Benic2015 ; Haensel2016 ; Bejger2017 ; Kaltenborn2017 ; Alvarez-Castillo2017 ; Abgaryan2018 ; Maslov2018 ; Blaschke2018 ; Alvarez-Castillo2019 , with several of these earlier works focusing on “twin stars”, which are hybrid stars in the third family whose mass is in the same range as the mass of neutron stars, but have smaller radii. Interest in hybrid stars has recently increased due to GW170817, because GW170817 is also compatible with at least one component being a hybrid star as first pointed out in Paschalidis2018 . Subsequent studies also focused on prospects of future multi-messenger detections of the merger of binary compact stars with hybrid hadron-quark EOSs (e.g., Drago2018 ; Burgio2018 ; Han2018 ; Chen2018 ; Annala2019 ; Most2019 ; Bauswein2019 ; DePietri:2019khb ).

Motivated by the fact that such hybrid stars are compatible with GW170817, in this paper we focus on improving our understanding of uniformly and differentially rotating hybrid stars. These have received only little consideration Ayvazyan2013 ; Haensel2016 and many of their properties are largely unexplored. We adopt the same EOSs with first-order phase transitions to quark matter introduced in Paschalidis2018 , and focus on answering the following set of questions: What is the maximum mass for uniformly and differentially rotating hybrid stars? Are these compatible with constraints from GW170817? Do rotating hybrid stars satisfy the same existing universal relations that hadronic EOSs satisfy?

The maximum mass that a non-rotating relativistic star222We will interchangeably also use the terms “static” and “TOV” to indicate non-rotating stars Tolman1934 ; Oppenheimer1939 . can support – the TOV limit – is determined by the EOS. However, rotation can increase the maximum supportable mass Friedman1987 ; Cook1994 ; Eriguchi1994 ; PaschalidisStergioulas2017 . In the case of hadronic EOSs, uniform rotation increases the maximum mass by about $20\text{\,}\mathrm{\char 37\relax}$ Cook1994 ; Cook1994b ; Lasota1996 ; Morrison2004 ; Breu2016 .333However, it has been shown that strange quark-matter stars can support up to $40\text{\,}\mathrm{\char 37\relax}$ more than the non-rotating limit Gourgoulhon1999 ; Rosinska2000 ; Bhattacharyya2016 . Equilibrium configurations whose mass exceeds the TOV limit are commonly referred to as “supramassive” Cook1992 . The maximum mass supported when allowing for maximal uniform rotation is known as the “supramassive limit”. Given the approximate EOS-independence (for hadronic EOSs), the supramassive limit has been used to place constraints on the TOV limit mass using observations (e.g. Most2018 ; Ruiz2018 ). Whether this universality for the ratio of the supramassive limit to the TOV limit mass holds for hybrid stars has not been investigated before and we test this here. Using the code of Cook1992 , we construct mass shedding sequences and find that for the hybrid EOSs we treat, this ratio is between $1.15$ and $1.30$. As a result, the universality does not hold for hybrid stars or, at least, the spread of the universality is significantly enhanced. However, using the same arguments as in Most2018 ; Ruiz2018 , we find that the supramassive limit of our hybrid EOSs is consistent with GW170817, lending further support to the results of Paschalidis2018 that GW170817 is compatible with hybrid hadron-quark EOSs.

Differential rotation can greatly boost the maximum mass that a star can have Baumgarte2000 ; Morrison2004 . The term “hypermassive” Baumgarte2000 is used to describe those stars that are more massive than the supramassive limit, and are supported by the additional centrifugal support provided by differential rotation. There also exist stars that can support more than two times the TOV limit mass (these were termed “ubermassive” in Espino2019 ), and which can arise only when differential rotation is present. Differential rotation plays an important role in temporarily stabilizing a binary neutron star merger remnant, and likely many merger remnants go trough a hypermassive phase (see, e.g., Baiotti2017 ; PaschalidisStergioulas2017 ; Paschalidis2017 for recent reviews on binary neutron star mergers). Differentially rotating stars have been shown to exhibit a rich and interesting solution space (see Ansorg2009 ; Studzinska2016 ; Rosinska2016 for polytropes, Espino2019 for hadronic EOSs, and Zhou:2019hyy ; Szkudlared2019 for strange quark stars). For instance, as the degree of differential rotation varies, the topology of the star can change from spheroidal to quasi-toroidal, where the maximum energy density does not occur at the geometric center of the star, but in a ring around the stellar center of mass. Studies of the dynamical stability of such quasi-toroidal configuration are under way and it is found that these configurations are dynamically unstable Espino2019b . For this reason, in this paper we focus primarily on the maximum mass that differentially rotating spheroidal hybrid star configurations can support using the differential rotation law of Komatsu1989a ; Komatsu1989b . The motivation for this study is that quark deconfinement can happen not only for isolated stars, but also following a binary neutron star merger Most2019 , that can result in a differentially rotating hybrid star. Using the code of Cook1992 , we find that a lower-limit on the maximum supportable rest mass is $123\text{\,}\mathrm{\char 37\relax}$ more than the TOV limit rest mass. This means that ubermassive spheroidal hybrid hadron-quark configurations can exist.

In the context of binary mergers, EOS-independent relations have proven useful as tools to devise tests for general relativity and break degeneracies in gravitational wave signals. For instance, the universal I-Love-Q relations Yagi2013 ; Yagi2013b ; Yagi2014 ; Yagi2017 were used in Abbott2018b to tighten the constraints placed on the tidal deformability of the two components of the GW170817 binary. In Paschalidis2018 it was verified that, for the EOSs we adopt here, the I-Love-Q relations are still satisfied with a spread of at most $3\text{\,}\mathrm{\char 37\relax}$ in both slowly and rapidly rotating hybrid stars.444Note that this result is in tension with what was found in Bandyopadhyay2018 where different equations of state were adopted. In Bozzola2018 , relations that do not depend on the EOS or the degree of differential rotation were found for turning-point models, which are the critical points of the gravitational mass along equilibrium sequences with fixed angular momentum or rest mass. For uniformly rotating configurations, turning points locate the onset of the unstable branch, and for this reason are useful to study the stability properties of rotating stars Friedman1988 ; Takami2012 ; Bauswein2017 ; Weih2018 . We verify that the universal relations relating angular momentum, rest mass and gravitational mass for turning-point models found in Bozzola2018 hold for hybrid hadron-quark EOSs when uniform rotation is considered, but have a clear, albeit weak, dependence on the degree of differential rotation.

The rest of the paper is structured as follows. In Section 2 we review the EOSs we consider in this work. In Section 3 we briefly describe the code we use to perform our study. In Section 4 we compute the supramassive limit for hybrid stars and discuss its implications for GW170817. We also test the previously found universal relations of Breu2016 and Bozzola2018 for the maximum mass along constant angular momentum sequences, and point out that the minimum mass twin stars also satisfies a universal relation. In Section 5 we compute the maximum rest mass of differentially rotating spheroidal hybrid stars, and test the dependence on the degree of differential rotation of the aforementioned universal relations. We conclude in Section 6 with a summary of our main findings. Unless otherwise specified, throughout this paper we adopt geometrized units where $c=G=1$ (where $c$ is the speed of light in vacuum and $G$ the gravitational constant). Square brackets are used to designate the units of quantities. For instance, $[$\mathrm{\text{$M_{\odot}$}}$]$ indicates the solar mass. In addition, when we use the word “mass” we refer to the gravitational (Arnowitt-Deser-Misner – ADM) mass for which we use the symbol $M$ (or $M_{\rm ADM}$) with other sub- or superscripts except for the subscript [math]. We will specify “rest” or “baryon” mass otherwise, and to designate rest mass we use the subscript [math] with other sub- or superscripts, e.g., $M_{0}$.

2 Equations of State

The hybrid hadron-quark EOSs we employ are the same as those developed in Paschalidis2018 for the EOSs introduced in Alford2017 (Set I) and Alvarez-Castillo2019 (Set II). Here we briefly describe these EOSs, and refer the interested reader to Paschalidis2018 for more details. Note that we use a different naming convention compared to Paschalidis2018 , and the map between the two is listed in Table 1, where we also report the maximum mass for static stars in both the hadronic and third family branches.

The EOSs correspond to zero-temperature matter in beta-equilibrium with a low-density hadronic phase and a high-density quark phase that are matched through a first-order phase transition. High tension at the quark-hadron interface is assumed, leading to a sharp transition boundary between the two phases. A low-density crust component based on the model of Baym1971 ; Negele1973 is also added to both sets of EOSs. In all EOS models the pressure matching between the two phases is performed via a Maxwell construction.

The Set I EOSs have a hadronic part that follows Colucci2013 in applying a covariant density functional theory Ring2010 with density-dependent couplings Lalazissis2005 . The quark phase in Set I is described by a MIT bag model Witten1984 ; Haensel1986 ; Alcock1986 ; Zdunik2000 , adopting the constant sound speed ($c_{s}$) parametrization Seidov1971 ; Zdunik2013 ; Alford2014 ; Alford2017 . Hence, the Set I quark phase pressure $P$ as a function of the energy density $\epsilon$ is given by

[TABLE]

with $P_{\text{tr}}$ value of the pressure at the phase transition (which occurs in the energy density range $\epsilon_{1}\leq\epsilon\leq\epsilon_{2}$). The energy density jump is parametrized by the value of $\xi$ defined through555In Paschalidis2018 the letter $j$ is used instead of $\xi$.

[TABLE]

The values of $\epsilon_{1},\xi$, $P_{\text{tr}}$ and $c_{s}$ that produce the Set I EOSs (labeled here T$1$–T$9$ and TT) are listed in Table I of Paschalidis2018 .

Set II Alvarez-Castillo2019 consists of a piecewise polytropic representation Read2009 ; Hebeler2013 ; Raithel2016 of the quark phase in which the pressure is

[TABLE]

where $n$ is the baryon number density and $n_{0}$ its value at the nuclear saturation density. The second polytrope has $\Gamma_{2}=0$ so that the pressure is constant, corresponding to the phase transition. The various coefficients that specify these EOSs (labeled here as A$4$–A$7$) are listed in Table II of Paschalidis2018 . Note that EOSs with similar properties as the EOSs in Set I and II have been obtained recently within a relativistic density functional approach to quark matter Bonanno2012 ; Ayvazyan2013 ; Kaltenborn2017 .

Figure 1 shows the pressure $P$ as a function of the energy density $\epsilon$ for the different EOSs, and Figure 2 depicts the resulting mass-areal radius relationship for non-rotating compact stars. The EOSs are compatible with the observational constraint of a pulsar mass of $($1.97$\,\pm\,$0.04$)\,$\mathrm{\text{$M_{\odot}$}}[Demorest2010 ](#bib.bib106); [Antoniadis2013 ](#bib.bib107); [Fonseca2016 ](#bib.bib108); [Arzoumanian2018 ](#bib.bib109), and half of them are also compatible with the newly discovered2.17$^{+0.11}_{-0.10}~{}$\mathrm{\text{$M_{\odot}$}} pulsar J0740+6620 [Cromartie2019 ](#bib.bib110). The EOSs span a relatively large range of compactness and radii (10$-$13.5\text{,}\mathrm{km}$$), and as such allow us to test for universal relations. Despite the fact that some of the EOSs in our sample are incompatible with J0740+6620, we include them in our study to check whether existing universal relations are respected independently of the EOS.

3 Methods

We assume that neutron star matter can be described as a perfect fluid with no meridional currents. An equilibrium model for a rotating general-relativistic star is a stationary and axisymmetric solution of Einstein’s equations coupled with the equation of hydrostationary equilibrium. Under these assumptions, the fluid four-velocity is given by

[TABLE]

and the line element ${\operatorname{d}\!{s}}^{2}$ of the spacetime is

[TABLE]

with $(t,r,\theta,\phi)$ quasi-isotropic coordinates, and $\gamma,\rho,\alpha,\omega$ spacetime potentials that depend on the coordinates $r,\theta$ only Bardeen1971 ; Friedman2013 ; PaschalidisStergioulas2017 .

To solve the equilibrium equations, an EOS and a rotation law have to be supplied (see PaschalidisStergioulas2017 for a recent review on theoretical and numerical approaches to the construction of rotating relativistic stars). For uniform rotation, we fix $\Omega=u^{\phi}/u^{t}$ (the local angular velocity of the fluid as seen by an observer at rest at infinity) to be constant. For differential rotation, we adopt the Komatsu-Eriguchi-Hachisu (KEH) $j$-constant law Komatsu1989a ; Komatsu1989b , which is described by

[TABLE]

where $\Omega_{c}$ it is the angular velocity evaluated on the rotation axis, and $A$ is a parameter with units of length that determines the lengthscale of variation of angular velocity within the star. This rotation law is not the only possible choice, but it is most commonly studied due to its simplicity (for a summary of other differential rotation laws, see PaschalidisStergioulas2017 ). In the following, instead of using $A$, we work with the dimensionless parameter $\hat{A}^{-1}$ defined by

[TABLE]

where $r_{e}$ is the coordinate radius of the star at the equator. Stars with $\hat{A}^{-1}=0$ are uniformly rotating.

We use the code developed in Cook1994 ; Cook1994b (the “Cook code”) to solve the structure equations and construct equilibrium models for rotating stars. For a given EOS and degree of differential rotation $\hat{A}^{-1}$, a rotating equilibrium model is built by providing the maximum energy density $\epsilon_{\rm max}$ and the ratio of the polar ($r_{p}$) to equatorial coordinate radius $r_{e}$. In the case of differential rotation with the KEH law, the set of parameters $\{\epsilon_{\rm max},r_{p}/r_{e},\hat{A}^{-1}\}$ does not uniquely specify an equilibrium model, in spite of the fact that it is sufficient to find a solution Ansorg2009 . To distinguish between physically distinct models that have the same $\{\epsilon_{\rm max},r_{p}/r_{e},\hat{A}^{-1}\}$, the parameter

[TABLE]

was introduced in Ansorg2009 . Here $\varpi=r\sin(\theta)$ and $z=r\cos(\theta)$ are cylindrical coordinates, and the derivative is evaluated on the surface of the star at the equator. In practice, it is more convenient to work with a new quantity $\hat{\beta}$ defined in terms of $\beta$ as

[TABLE]

The parameter $\hat{\beta}$ describes how close to the mass-shedding limit a configuration is, and approaches three limiting values according to the stellar shape in the solution space: in the case of a spherical solution, $\beta=1$ and so $\hat{\beta}$ is ${1}/{2}$; for a model close to the mass-shedding limit, the stellar surface is highly pinched near the equator, which implies that the derivative appearing in Equation (8) vanishes, and hence $\hat{\beta}$ approaches [math]; finally, for a quasi-toroidal shape the ratio ${r_{p}}/{r_{e}}$ becomes very small, and hence $\beta$ becomes large. In this limit $\hat{\beta}$ approaches $1$.

In the case of the KEH rotation law, the solution space of differentially rotating stars is fully specifiable by the quadruplet of parameters $\{\epsilon_{\rm max},r_{p}/r_{e},\hat{A}^{-1},\hat{\beta}\}$ Ansorg2009 . Our code does not have the capability of fixing $\hat{\beta}$ directly, but, as in Espino2019 , we produce models with different values of $\hat{\beta}$ given the other parameters $\{\epsilon_{\rm max},r_{p}/r_{e},\hat{A}^{-1}\}$ by following specific trajectories in the parameter space. This is especially relevant for our discussion in Appendix A, where we study the different solution types for hydrid stars.

4 Uniform Rotation

In this section we focus on the case of uniform rotation. We compute the supramassive limit mass and the maximum mass on sequences of constant angular momentum, and test the applicability of the universal relations found in Breu2016 and Bozzola2018 to the hydrid hadron-quark EOSs listed in Table 1. Differential rotation is treated in Section 5.

One of the most important macroscopic parameters related to an equation of state is the maximum mass that an equilibrium rotating star can support.666The most massive models are typically dynamically unstable Friedman1988 ; Takami2012 ; Weih2018 . In this paper we ignore the issue of stability and focus on constructing equilibrium models. To study supramassive hybrid stars we construct equilibrium sequences with constant angular momentum $J$. For these computations, we set the angular resolution of the Cook code to 300 points, the radial to 500 and used a basis of 20 Legendre polynomials. Then, each sequence is constructed with at least 500 models logarithmically equispaced in the central energy density. Examples of such sequences are depicted in Figure 3 for EOS T9, where boxes indicate the stars with maximum mass along each $J$-constant sequence. In Figure 3, we also highlight other notable models: the turning points indicated with filled diamonds and filled circles. In Figure 3, most of the circles lie inside a box, since those turning points are also the most massive model of the sequence.

Given a sequence of equilibrium models with constant angular momentum $J$, and varying central energy density $\epsilon_{c}$, the turning points are the stationary points of the function $M(\epsilon_{c})$, where $M(\epsilon_{c})$ is the gravitational mass along the sequence. For a constant $J$ sequence, a turning point is defined by the condition

[TABLE]

For constant angular momentum sequences and hadronic EOSs, the turning points are local maxima, whereas for constant $M_{0}$ sequences the turning points are local minima. By virtue of the turning point theorem Sorkin1982 ; Friedman1988 that applies to uniformly rotating configurations, the turning-point lines for $J$-constant sequences coincide with those of $M_{0}$-constant sequences. In this work we located the turning points along $J$-constant sequences by interpolating with a cubic spline around the local maxima and minima model of the sequence, and then finding the roots of the derivative of the cubic spline interpolant.

In the case of EOSs exhibiting a hadron-quark phase transition with a large jump in energy density, typically each $J$-constant sequence has multiple turning points. In this paper we focus on the two turning points of $J$-constant sequences that occur at central energy densities higher than the phase-transition. We call them “top” (labeled with an upward pointing arrow $\uparrow$) and “bottom” (labeled with a downward pointing $\downarrow$), indicating the higher and lower mass turning point, respectively. In all figures diamonds represent bottom turning points, filled circles indicate top turning points, and boxes designate the most massive models. Top turning-point models are interesting because for uniformly rotating stars they detect the onset of the instability to collapse Friedman1988 ; Takami2012 . On the other hand, bottom turning points indicate the least massive twin star which must be unstable, too. Not all EOSs in our set have bottom turning points: T1, T2 and T4 have too small a jump $\xi$ to produce such configurations (as can be seen in Figure 2, where these three EOSs do not undergo a dip above the phase transition energy density). These EOSs do not exhibit twin star solutions.

The TOV sequences are constant angular momentum ($J=0$) sequences. Their top turning points indicate the maximum mass allowed for static configurations,777This is not true for EOS A6, since, for this EOS, non-rotating stars with central energy density within the phase transition have mass $M=$2.01\text{,}\mathrm{\text{$M_{\odot}$}}, whereas the top turning point has mass $M=$2.00\text{\,}\mathrm{\text{$M_{\odot}$}}. Since our focus is on the stable branch of the third family of stars, the top turning point is the maximum mass of the third family. and their bottom turning points indicate the minimum mass of non-rotating twin stars. These maximum and minimum masses provide fundamental mass scales associated with each EOS, and will be used to normalize other quantities. We will append the superscript “TOV” to designate these mass, so that $M^{\text{TOV}}_{\text{Max}}$ and $M_{0,{\text{Max}}}^{\text{TOV}}$ are the maximum gravitational and baryonic mass of a non-rotating hybrid star in the third family, respectively, whereas $M^{\text{TOV}}_{\downarrow}$ and $M_{0,\downarrow}^{\text{TOV}}$ are the minimum gravitational and rest mass for a static twin star.888Using the maximum possible TOV mass $M=$2.01\text{,}\mathrm{\text{$M_{\odot}$}} (instead of the maximum TOV mass in the third family $M=$2.00\text{\,}\mathrm{\text{$M_{\odot}$}}) to normalize mass scales for EOS A6 has negligible impact on our results.

4.1 The supramassive limit

An EOS-insensitive relation for hadronic EOSs that relates the maximum mass ($M_{\text{Max}}$) and angular momentum along $J$-constant sequences was recently found in Breu2016 . The resulting universal relation (Equation (12) in Breu2016 ) is expressed in terms of $J$ normalized to the maximum possible angular momentum $J_{\text{Max}}$ that can be achieved with uniform rotation, and is the following

[TABLE]

This last equation can be used to estimate the supramassive limit for hadronic EOSs, and yields the approximately EOS-independent value of $20\text{\,}\mathrm{\char 37\relax}$ larger than the TOV limit. Here we are interested in testing whether Equation (11) applies to hybrid hadron-quark EOSs. Moreover, we explore whether a similar relation exists for hybrid EOSs not only for the top turning points, but also for the bottom ones. As a reminder, we note that for hybrid EOSs, the top turning point is not necessarily the most massive model along a $J$-constant sequence, as shown for example in Figure 4.

4.1.1 Maximum mass

In Figure 5 we plot the maximum possible gravitational mass normalized to $M^{\text{TOV}}_{\text{Max}}$ on $J$-constant sequences as a function of $J$ normalized to $J_{\text{Max}}^{\text{supra}}$, the angular momentum of the most massive uniformly rotating star i.e., the supramassive limit configuration. This choice ensures that the maximum mass is reached when $J/J_{\text{Max}}^{\text{supra}}=1$. As is clear from Figure 5, some of the EOSs in our set produce uniformly rotating stars that can support more mass than 1.20$\,M^{\text{TOV}}_{\text{Max}}$. On the other hand, some EOSs cannot even reach such a large enhancement in mass compared to the TOV limit mass. To be more specific, the supramassive limit mass varies from 1.15$\,M^{\text{TOV}}_{\text{Max}}$ for TT to $$1.31$\,M^{\text{TOV}}_{\text{Max}}$ for A6.

4.1.2 Implications for GW170817

In Most2018 ; Ruiz2018 the result that for hadronic EOSs the ratio of the supramassive limit to the TOV limit is EOS-independent was combined with the observation of gravitational waves from event GW170817 to set an upper bound on the TOV limit mass. We repeat here the same argument to highlight the differences when EOSs with a hadron-quark phase transition are considered.

As shown by simulations, the observed properties of the gamma-ray burst and the ejecta following GW170817 suggest that the event left a rotating remnant supported by differential rotation. The remnant then underwent delayed collapsed to a black hole (see arguments in Most2018 ; Ruiz2018 ), meaning that its mass was above the supramassive limit $M^{\text{supra}}_{\text{Max}}$. As we already mentioned, for hadronic EOSs it is possible to write

[TABLE]

with $\alpha\approx$1.2, independently of the EOS. Then, assuming the low-spin prior LIGO analysis for GW170817, the total gravitational mass of the system as inferred by the gravitational wave signal is $M_{\text{tot}}\approx{2.74^{+0.04}_{-0.01}}~{}$\mathrm{\text{$M_{\odot}$}} Abbott2017a , and since this quantity has to be greater than the supramassive limit, the following holds

[TABLE]

With this inequality, the authors of Most2018 ; Ruiz2018 placed a constraint on the value of the maximum mass of a non-rotating configuration

[TABLE]

Figure 5 shows the loss of EOS-independence of the supramassive limit mass, which can be between $15\text{\,}\mathrm{\char 37\relax}$ to $31\text{\,}\mathrm{\char 37\relax}$ more than the TOV limit mass. Assuming that we can extend Equation (12) to include hybrid stars with 1.15$\leq\alpha\leq$1.31, using the same steps as in the case of hadronic EOSs the upper bound when hybrid EOSs are considered is revised to

[TABLE]

The lower bound in Equation (15) comes from most recent massive pulsar J0740+6620 Cromartie2019 . The upper bound is relaxed from 2.28\text{\,}\mathrm{\text{M_{\odot}}} to 2.38\text{\,}\mathrm{\text{M_{\odot}}} to encompass the EOSs treated in this work.

The assumption that we can extend Equation (12) (with 1.15$\leq\alpha\leq$1.31) to include hybrid stars, should be further tested using hybrid hadron-quark EOS with different baseline hadronic EOSs. But, this goes beyond the scope of the current work. However, that GW170817 is compatible with hybrid EOSs based on the fact that the supramassive limit mass should be smaller than the inferred mass from GW170817 is demonstrated in Table 2, where we show the supramassive limit mass and the supramassive limit in the third family with each of the EOSs we treat in this work. As is evident from the table, the largest supramassive limit corresponds to EOS A6 which is 2.63\text{\,}\mathrm{\text{M_{\odot}}}, and hence safely smaller than even 2.73\text{\,}\mathrm{\text{M_{\odot}}} (the lower limit on the mass from GW170817). This result lends further support to the finding of Paschalidis2018 that GW170817 can be interpreted as the inspiral of a binary hybrid star–neutron star.

4.2 Universal relations for turning points

Turning points mark the onset of a radial instability of the star, and provide a sufficient condition for secular axisymmetric instability in uniformly rotating (isentropic) stars Sorkin1982 ; Friedman1988 ; Takami2012 ; Weih2018 (see also Schiffrin:2013zta ; Prabhu:2016pei ). Top turning points identify an instability to collapse and in the case of our EOSs designate the most massive hybrid stars. On the other hand, the bottom turning points single out the least massive twin stars. For this reason, universal relations involving turning point models are useful when studying the stability of rotating stars. In this section, we focus on the turning points in the third family of compact objects and discover a new relation, similar to Equation (11), for both top and bottom turning points. We also test the universal relations reported in Bozzola2018 , where it was shown that the angular momentum and masses of differentially rotating turning-point models satisfy relations that are approximately independent of either the degree of differential rotation (and hence apply to uniformly rotating stars) or the EOS.

4.2.1 Top turning points

To normalize the top turning point mass and angular momentum we first define the “Kepler turning point” configuration in the $(\epsilon_{c},M)$ plane as the model that lies at the intersection of the mass-shedding sequence with the turning-point line, i.e., the line connecting the top turning points along $J$-constant sequences. Note that the Kepler turning point model is neither the most massive nor more massive than the supramassive limit in the third family (see, e.g, Figure 4), but is very close to the supramassive limit of the third family. To study universal relations of top turning points, we normalize their angular momentum with that of the Kepler turning point, which we denote $J_{\text{Kep}}$. In practice, we find the Kepler turning point by least-squares fitting the turning-point line with a fifth-order polynomial, and extending the resulting function until it intersects the mass-shedding sequence. Given that for each EOS we construct between 15 and 40 $J$-constant sequences, we never need to extrapolate more than $1\text{\,}\mathrm{\char 37\relax}$ out of the range of energy densities where we have data.

In Figure 6 we plot the top turning point mass normalized by the TOV limit mass versus the top turning point normalized angular momentum. The plot displays a stronger degree of universality compared to the one for the most massive models (Figure 5). The best-fitting function is

[TABLE]

This expression approximates all the EOSs we study here with a spread of at most $3\text{\,}\mathrm{\char 37\relax}$, which is only a little larger than the $2\text{\,}\mathrm{\char 37\relax}$ spread found for hadronic EOSs Breu2016 . As in the case of Breu2016 , the spread increases with the angular momentum, and for $J/J_{\text{Kep}}\lesssim$0.5$$ all the EOSs agree with Equation (16) to within approximately $1\text{\,}\mathrm{\char 37\relax}$. The spread at higher angular momenta is mostly due to EOSs T1, T2 and T4.

By use of Equation (16), we can also find the supramassive limit mass of the third family $M_{\text{Max, third}}^{\text{supra}}$. This is because the flatness of the mass-shedding sequence near the supramassive limit (see for example Figures 3 and 4) ensures that the Kepler turning point mass agrees with the supramassive limit mass of the third family to better than $0.1\text{\,}\mathrm{\char 37\relax}$ for all EOSs we consider. Moreover, the error in estimating the supramassive limit mass due to considering the Kepler turning point angular momentum instead of the supramassive limit angular momentum is smaller than the dependence from the EOS. Plugging $J=J_{\text{Kep}}$ in Equation (16), we obtain

[TABLE]

where the spread corresponds to the largest discrepancy from the mean value in the set of EOSs we study (shown as the gray shaded region in Figure 6).

In Figure 6 we also compare our finding with the universal relation of Breu2016 , since turning points were considered there, too. We find that if we allow for a $3\text{\,}\mathrm{\char 37\relax}$ spread, Equation (12) in Breu2016 holds up to $J/J_{\text{Kep}}\approx$0.85$$, but overestimates the mass for larger angular momenta. Equation (12) of Breu2016 can describe the EOSs studied here with larger angular momenta with an increased spread of $4.2\text{\,}\mathrm{\char 37\relax}$ (cyan shaded region in Figure 6). Thus, our expression (16) provides a better fit to the data, and as is clear from Figure 6 Equation (16) better captures the trend of the data.

Finally, we note that for the EOSs of Set I, there is a correlation between the jump size $\xi$ and the value of $M_{\uparrow}$: for a given value of $J/J_{\text{Kep}}$, $M_{\uparrow}/M^{\rm TOV}_{\rm max}$ is larger for EOSs with smaller $\xi$ (as is clear from the top three curves in Figure 6). As a result, T1 has the largest normalized maximum mass (1.19$\,M^{\text{TOV}}_{\text{Max}}$) and TT the smallest (1.15$\,M^{\text{TOV}}_{\text{Max}}$). This does not happen for the EOSs in Set II which are all comparable (regardless of the width of the phase transition) and have no emerging trend.

Next, we consider the EOS-independent relations found in Bozzola2018 between angular momentum, gravitational and rest mass of top turning points, which we report here for convenience

[TABLE]

Figure 7 depicts these quantities and shows that the top turning points are aligned with the previously known equation. In fact, we find that these models satisfy all the three relations (18) with the same spread as hadronic EOSs ($1.6\text{\,}\mathrm{\char 37\relax}$ for Equation (18a)999This relation carries similar information compared to Equation (16), and exhibits better universality. and $1.2\text{\,}\mathrm{\char 37\relax}$ for Equations (18b) and (18c)). This is a non-trivial result because in Bozzola2018 it was noted that purely quark stars described with the MIT bag model do not follow the same universal relations, even in the uniform rotation case. Nonetheless, here we find that for the EOSs of Set I, which employs a variant of the MIT bag model for the high-density phase, the relations are still satisfied. This suggests that the hadronic low-density part plays an important role in determining the universality.

4.2.2 Bottom turning points

We discovered that bottom turning points satisfy a relation similar to Equation (16), if the normalizing angular momentum is that of the Kepler bottom turning point (found at the intersection of the bottom turning-point line with the mass shedding sequence)

[TABLE]

The spread in this equation is at most $2\text{\,}\mathrm{\char 37\relax}$.101010Note that T1, T2, T4 do not have any bottom turning point.

Moreover, we find that the bottom turning points can be described with the same Equations (18), but with a spread of $3\text{\,}\mathrm{\char 37\relax}$. The universality becomes tighter if we consider top and bottom turning points separately. For the bottom turning points, the best-fitting functions are

[TABLE]

with largest spread of $1\text{\,}\mathrm{\char 37\relax}$ for the first two and $0.2\text{\,}\mathrm{\char 37\relax}$ for the last. Figure 8 shows Equation (20a), and in the inset we compare it to Equation (18a), which is the corresponding universal relation for top turning points. Among the three universal relations, Equations (20a) and (18a) are the ones in which top and bottom turning point differ the most.

5 Differential Rotation

In this section we consider differentially rotating stars with hybrid EOSs and compute the maximum mass that can be supported. In this context, it is more common to consider the rest mass instead of the gravitational mass Morrison2004 . Nevertheless, the results that will follow are qualitatively the same also for the gravitational mass.

The KEH rotation law of Equation (6) produces stars with two different topologies: spheroidal and quasi-toroidal. In quasi-toroidal configurations the maximum energy density occurs in a ring around the stellar geometric center, whereas in spheroidal stars the location of the maximum energy density coincides with geometric center of the star. Examples of such configurations are reported in Figure 9, where we show meridional energy density contours for a spheroidal and a quasi-toroidal solution. In this paper, we focus primarily on spheroidal models because quasi-toroidal configurations are likely to be unstable Espino2019b .

5.1 Maximum mass

To find the maximum rest mass for each EOS, we fix the degree of differential rotation $\hat{A}^{-1}$, and scan the maximum energy density in the range $\epsilon_{\text{max}}\in($0.5$,$3$)~{}${10}^{15}\text{,}\mathrm{g}\text{,}{\mathrm{cm}}^{-3}$$, which is where the most massive models always lie for the EOSs considered here. For each value of $\epsilon_{\text{max}}$ we construct equilibrium models with decreasing value of $r_{p}/r_{e}$ until we find a quasi-toroidal solution or we reach the mass-shedding limit. The last configuration built before the quasi-toroidal or mass-shedding limit is reached is the maximum mass for that value of $\epsilon_{\text{max}}$ in our search.

In Figure 10 we show the maximum rest mass normalized to the rest-mass of the TOV limit as a function of $\hat{A}^{-1}$. The behavior is similar to that seen for polytropes Studzinska2016 and hadronic EOSs Espino2019 and it is followed very closely by the gravitational mass: at relatively low values of $\hat{A}^{-1}$ ($\lesssim 0.25$), corresponding to rotational profiles closer to uniform rotation, there are only modest increases in the rest mass, not exceeding $40\text{\,}\mathrm{\char 37\relax}$. For higher values of $\hat{A}^{-1}$, the increase of $M_{0,\text{Max}}$ compared to the TOV limit grows until it reaches a maximum and then begins to decrease again as $\hat{A}^{-1}$ is increased further. The largest gain in rest mass is seen for EOS A6 ($123\text{\,}\mathrm{\char 37\relax}$ more than the TOV limit), which is also the EOS with the highest gain in the supramassive limit (see Section 4).

In Table 3 we report relevant properties of the most massive spheroidal models found in our search of the parameter space. Larger increases in the rest mass are possible if we consider quasi-toroidal stars and possibly other differential rotation laws. Thus, these masses should be considered as lower-limits on the maximum mass of even spheroidal stars. Moreover, finding the “absolute” maximum rest mass model for hybrid EOSs with the KEH law depends sensitively on the resolution in $\epsilon_{\rm max}$, $r_{p}/r_{e}$, and $\hat{A}^{-1}$ used to scan the parameter space. For the models presented in Table 3 we use a step in $\epsilon_{\rm max}$ and $r_{p}/r_{e}$ of 0.01, and a step in $\hat{A}^{-1}$ of 0.005. Even with this high resolution in $\hat{A}^{-1}$, we found large increases in the rest mass over small ranges of $\hat{A}^{-1}$. Such increases in mass correspond to a transition to a different solution type, as shown in the inset of Figure 10, where we show the maximum rest mass type A models along with the maximum rest mass spheroidal solutions (which also include type C models) for the A4 EOS. We explore the maximum rest mass of each solution type in Appendix A, where more details on the different types of solutions of differentially rotating stars are also presented.

5.2 Universal relations

In Bozzola2018 it was found that for hadronic EOSs and restricting to spheroidal models, Equation (18) applies not only to uniformly rotating stars, but also to differentially rotating stars, and, what is more, the relation is approximately insensitive to the degree of differential rotation (with maximum spread of about $1\text{\,}\mathrm{\char 37\relax}$ even for high degrees of differential rotation). Moreover, recently, in Zhou:2019hyy it was noted that for purely strange quark-matter stars the relations between angular momentum and masses of turning points follow equations approximately independent of $\hat{A}^{-1}$ (with spread below $2\text{\,}\mathrm{\char 37\relax}$), although the relation is not the same as for hadronic EOSs Bozzola2018 . Here, we construct sequences of differentially rotating equilibrium models with constant angular momentum and various degrees of differential rotation in the range $\hat{A}^{-1}\in($0$,$2$)$ to test whether this property holds for the hybrid hadron-quark EOSs we treat in this work.

We find that relations (18) apply approximately to hybrid EOSs for degrees of differential rotation $\hat{A}^{-1}\in($0$,$2$)$, but with a larger spread of about $3\text{\,}\mathrm{\char 37\relax}$ and with an evident (albeit weak) dependence on $\hat{A}^{-1}$. Furthermore, the larger the degree of differential rotation, the larger the spread becomes. For this reason, the universality found in Bozzola2018 appears to be broken in the case of hybrid stars.

For top turning points, we show an example of the above conclusions in Figure 11, where we report the results from $J$-constant sequences constructed with EOS TT. We found that there is a correlation between the jump size $\xi$ and the loss of the universality: the larger the jump, the stronger the deviation from universality. Hence, the example in Figure 11 is one of the cases that violates the universality the strongest. For most of the EOSs we treat, the bottom turning point universal relations exhibit a smaller spread with increasing $\hat{A}^{-1}$ (below $2\text{\,}\mathrm{\char 37\relax}$), but they still violate Equations (20) for some EOSs, as reported in Figure 12 for EOS A7. For this EOS, the deviation from Equations (20) can be up to $4\text{\,}\mathrm{\char 37\relax}$.

Both Figure 11 and 12 show that there is a clear dependence on $\hat{A}^{-1}$. Hence, we conclude that the universal relations found in Bozzola2018 are EOS-independent for uniform rotation, but not for differential rotation, as they acquire a dependence on $\hat{A}^{-1}$ which is not the same for all the EOSs. The regularity and smoothness of sequences in the figures suggests that it is possible to parametrize the dependence on $\hat{A}^{-1}$ using some scale characteristic of the EOS. However, exploring this goes beyond the scope of the current work.

Finally, we point out that for sufficiently high values of $\hat{A}^{-1}$ and $J$, in some EOSs the $J$-constant sequences exhibit no bottom turning points. This is because stars whose central energy density is in the phase transition region become more massive when spun up compared to the ones with higher central energy density. This means that $J$-constant sequences in the $(\epsilon_{c},M)$ plane tend to become flatter and have shallower dips (see, for example, Figure 4). In some other cases with high degree of differential rotation and angular momentum we could not even identify the top turning point, because the change of concavity is not resolved by the accuracy of our code.

6 Conclusions

The equation of state of nuclear matter is uncertain at densities larger than nuclear saturation density $\rho_{0}$. Some of the proposed EOS models undergo a phase transition from hadronic to quark matter which can result in hadron-quark hybrid stars. For EOSs with a sufficiently large jump in energy density over which the pressure remains constant, a third family of compact objects emerges. In this work we analyzed the properties of rotating relativistic stars using EOSs with first-order hadron-quark phase transitions for many of which a third family of stars emerges. To be more specific, we employed the quark matter parametrizations introduced in Paschalidis2018 for the EOSs developed in Alford2017 and Alvarez-Castillo2019 .

We employed the Cook code Cook1994 ; Cook1994b to build rotating relativistic stars and we studied equilibrium sequences with constant angular momentum and varying central energy density. We found that the maximum mass a uniformly rotating hybrid star can support (the so-called supramassive limit) is not an EOS-independent quantity, but it varies from $15\text{\,}\mathrm{\char 37\relax}$ to $31\text{\,}\mathrm{\char 37\relax}$ more than the TOV limit mass, in contrast to previous work Lasota1996 ; Breu2016 which highlighted an approximately EOS-independent increase of about $20\text{\,}\mathrm{\char 37\relax}$ for hadronic EOSs. This implies that some constraints placed on the equation of state based on the universality of the supramassive limit and GW170817 Most2018 ; Ruiz2018 do not apply to EOSs like the ones considered here. However, the supramassive limit mass of the EOSs we adopt is consistent with GW170817, providing further support that hybrid hadron-quark EOSs can describe GW170817.

We located the turning points on constant angular momentum sequences, which are identified by the stationarity condition of Equation (10). We defined the bottom and top turning points along a $J$-constant sequence as those that correspond to the least massive twin star and the most massive star in the third family, respectively. Top turning points satisfy a universal relation [Equation (16)] with spread of about $3\text{\,}\mathrm{\char 37\relax}$. With this relation, we found that the supramassive limit of the third family is given by

[TABLE]

Enabling differential rotation and focusing on spheroidal models, we found that the maximum mass (both baryonic and gravitational) can increase by as much as $123\text{\,}\mathrm{\char 37\relax}$ compared to the non-rotating case. This enhancement compared to the TOV limit mass, represents a lower limit of what can be achieved with differential rotation.

Finally, we investigated the applicability of universal relations for turning points found in Bozzola2018 , who reported relations between the angular momentum, and the masses of turning points that are independent of the EOS and the degree of differential rotation. We found that the relations in Bozzola2018 hold for the top turning points for uniformly rotating stars. We also discovered similar universal relations for the bottom turning points for uniformly rotating stars. However, both universal relations are violated when differential rotation is enabled.

In future work we will consider constant rest-mass sequences, and explore additional hybrid equations of state with varying baseline hadronic EOSs, to further test some of the results we reported in this work.

Acknowledgments

We are grateful to S. L. Shapiro for access to the code that we used to build equilibrium models for rotating relativistic stars, and to D. Alvarez-Castillo, D. Blanschke and A. Sedrakian for giving us permission to use the equations of state developed in Paschalidis2018 . We also thank N. Stergioulas for access to the RNS code Stergioulas1995 ; Nozawa1998 that was used to cross-check some of the results presented here. Computations were in part performed on the Ocelote cluster at The University of Arizona. This research used also resources provided by the Open Science Grid Pordes2007 , which is supported by the National Science Foundation and the U.S. Department of Energy’s Office of Science.

Appendix A Appendix: Solution Space of Differentially Rotating Stars

In this appendix we discuss the solution space of equilibrium models for stars rotating with the KEH law of Equation (6). We find features in the solution space which are consistent with both hadronic and strange quark star EOSs, while also finding others which are unique to hybrid EOSs. In the case of differential rotation with the KEH law, solutions for a given maximum energy density can be separated on the $(r_{p}/r_{e},\hat{\beta})$ plane into four classes depending on the degree of differential rotation. These classes are known as types A, B, C and D Ansorg2009 . The first, type A, consists of spheroidal configurations with a relatively low degree of differential rotation. All uniformly rotating stars, if seen as differentially rotating with $\hat{A}^{-1}=0$, belong to this family. Type B solutions are characterized by a high central angular velocity but relatively low degree of differential rotation. These models are quasi-toroidal, since the maximum energy density does not occur in the geometric center of the star. In Figure 9, we show a comparison between the energy density profiles of a type A and a type B solution. The third category (type C) contains spheroids and quasi-toroids with higher degrees of differential rotation than types A and B. Finally, type D stars are typically highly pinched at the equator and quasi-toroidal. We were unable to construct type D solutions with the Cook code, but this is not a severe limitation on our study since these configurations are not likely to be found in Nature because their sequences start from a mass-shedding limit and end at a mass-shedding limit Studzinska2016 .

The classification of different solution types works as follows: for a fixed value of $\epsilon_{\rm max}$, there exists a critical value of the degree of differential rotational $\hat{A}^{-1}_{\rm crit}$, defined by the condition

[TABLE]

The function $\hat{\beta}(r_{p}/r_{e})$ that corresponds to $\hat{A}^{-1}_{\rm crit}$ for a given energy density is referred to as the “separatrix”. This quantity partitions the solution space into four regions, each corresponding to a different solution type. Models are classified directly using the separatrix in the $(r_{p}/r_{e},\hat{\beta})$ plane, as shown in Figure 13, where we show the separatrix for EOS A5 at energy density $0.5\text{\times}{10}^{15}\text{\,}\mathrm{g}\text{\,}{\mathrm{cm}}^{-3}$. Figure 13 is consistent with the results observed for polytropic and hadronic EOSs Ansorg2009 ; Espino2019 . However, since we were not able to build type D stars, there are three regions instead of four in Figure 13.

In practice, instead of working directly with the separatrix, it is more convenient to consider $\hat{A}^{-1}_{\rm crit}$ as a function of maximum energy density in order to assign a type to a configuration. For a fixed EOS and $\epsilon_{\rm max}$, all equilibrium models with $\hat{A}^{-1}<\hat{A}^{-1}_{\rm crit}$ are of either type A or B, whereas those with $\hat{A}^{-1}>\hat{A}^{-1}_{\rm crit}$ are of type C (or D). Then, a solution is categorized as spheroidal if the maximum energy density occurs in the center of mass. Since we do not have control over $\hat{\beta}$, we cannot directly solve Equation (22) to find the critical value for $\hat{A}^{-1}$. Hence, as in Espino2019 , we estimate $\hat{A}^{-1}_{\rm crit}$ with the maximum of the function $\hat{A}^{-1}_{\rm min}(r_{p}/r_{e})$, where $\hat{A}^{-1}_{\rm min}$ is the smallest degree of differential rotation for which an equilibrium model exists with the given maximum energy density. This method allows the calculation of $\hat{A}^{-1}_{\rm crit}$ to within $1\text{\,}\mathrm{\char 37\relax}$ accuracy Espino2019 .

Figure 14 shows the critical degree of differential rotation as a function of the maximum energy density for all the hybrid EOSs we consider in this work. The figure exhibits features which are unique to hybrid EOSs but which may be reconciled with the results of hadronic and strange-quark-matter EOSs. Crucially, $\hat{A}^{-1}_{\rm crit}$ is not a monotonically decreasing function of $\epsilon_{\rm max}$,111111We refer to the increase and subsequent decrease of $\hat{A}^{-1}_{\rm crit}(\epsilon_{\rm max})$ at values of $\epsilon_{\rm max}$ above the phase transition (as shown in Figure 14) as a “bump” feature. and, as a result, for some values of $\hat{A}^{-1}$, an equilibrium sequence with varying $\epsilon_{\rm max}$ (such as the ones considered in Section 5.2) can jump back and forth between being class A or B and C. The significant increases in mass that we find in Section 5 correspond to a jump in models from type A to type C.

In Ansorg2009 , homogeneous stars described by an incompressible EOS in which

[TABLE]

where $\epsilon_{0}$ is a constant, were considered. For such an EOS, the pressure and energy density do not depend on one another, similar to the constant pressure regions of the EOSs studied in this work. In Ansorg2009 it was found that the function $\hat{A}^{-1}_{\rm crit}(\epsilon_{\rm max})$ for incompressible EOSs first increases to a maximum and subsequently decreases as the energy density increases. This same feature was observed for strange-quark-matter EOSs, which may be suitably approximated as homogeneous bodies Szkudlared2019 . We note that similar features appear in the curves presented in Figure 14 above the values of energy density corresponding to the constant pressure regions of the EOSs. It is possible that the “bump” features seen in the curves in Figure 14 arise for the same reasons as those found in Ansorg2009 and Szkudlared2019 for homogeneous and strange-quark-matter EOSs, respectively.

In Table 4 we list the maximum rest mass models found in our search for each solution type we were able to construct with the Cook code and for three representative EOSs. We focus on the T1, T8, and A6 EOSs to highlight particular features of the solution space. As can be seen from Figure 14, the “bump” feature of the function $\hat{A}^{-1}_{\rm crit}(\epsilon_{\rm max})$ results in the appearance of A and B solutions at values of $\epsilon_{\rm max}$ and $\hat{A}^{-1}$ which are high compared to hadronic EOSs. Despite the presence of type A and B stars in larger regions of the parameter space compared to hadronic EOSs, we find that the maximum rest mass configurations still correspond to modest degrees of differential rotation and not for the largest energy densities considered. For all three EOSs, we find that the type B models are the most massive, which agrees with the findings from studies with hadronic EOSs Espino2019 . The largest gain in mass is for EOS A6, with which it is possible to construct stars almost three times more massive than the corresponding TOV limit. This increase in rest mass is larger than the ones seen for the hadronic EOSs studied in Espino2019 . Note that the type B stars tend to have the lowest value of $\epsilon_{\rm max}$ among the maximum rest mass solutions, suggesting that these configurations predominantly sample from the low density part of the corresponding EOS. Note also that for this reason the T1 and T8 type B models presented in Table 4 are the same. The T1 and T8 EOSs are identical in the low density regime and begin to diverge at $\epsilon\approx$0.66\text{\times}{10}^{15}\text{,}\mathrm{g}\text{,}{\mathrm{cm}}^{-3}, whereas the maximum rest mass stars for each of these EOS have $\epsilon_{\rm max}=$0.56\text{\times}{10}^{15}\text{\,}\mathrm{g}\text{\,}{\mathrm{cm}}^{-3}. Because the maximum energy density of these configuration is the same, and, because the EOSs are identical below this energy density, the corresponding solutions are identical. Any EOS that is the same in the low density regime as the T1 EOS will have an identical maximum rest mass type B model (including T2, T4, T6, T7, T9, and TT). However, note that because of their different TOV masses, each of these EOSs exhibit a different maximal increase in the rest mass compared to the TOV limit rest-mass. The maximum rest mass type C solutions presented in Table 4 are different between T1 and T8, because the corresponding values of $\epsilon_{\rm max}$ are greater than the energy density at which the two EOS begin to diverge from one another. The models presented in Table 4 are not the definitively most massive, but only the largest we were able to locate in the solution space. We were unable to scan the entire solution space using the Cook code (in particular we could not construct type D stars or highly pinched type B ones). Moreover, because different solution types can emerge as the maximum energy density changes, locating the absolute maximum rest mass models for the solution types we were able to compute depends on the resolution in the parameters ${\epsilon_{\rm max},r_{p}/r_{e},\hat{A}^{-1}}$ used to scan the solution space. We used a step in $\epsilon_{\rm max}$ and $r_{p}/r_{e}$ of $0.01$ and a step in $\hat{A}^{-1}$ of $0.01$ for the class A, but a step of $0.05$ for families B and C.