Slow Stable Self-bound Hybrid Star Can Relieve All Tensions

Phase transitions within compact stars are a long-studied possibility Dexheimer et al. (2018); Orsaria et al. (2019); Heiselberg and Hjorth-Jensen (1999); Glendenning (2001); Yang and Shen (2008); Reddy et al. (2000); Glendenning (1992); Burgio et al. (2002), which would provide a unique opportunity to understand the fundamental strong interaction between quarks, in both the perturbative and the non-perturbative regimes. More and more researches indicate that quark matter (QM), a state comprised of deconfined free-flowing quarks, can possibly exist inside the neutron star core from hadron-quark phase transitions, forming so-called hybrid stars (e.g., Alford et al., 2005, 2013).

It is possible that QM is absolutely stable (i.e., energy per baryon number $E/A<930$ MeV) even at zero pressure, either in the form of strange quark matter (SQM) Bodmer (1971); Witten (1984); Terazawa (1979); Farhi and Jaffe (1984) or up-down quark matter ($ud$QM) Holdom et al. (2018); Wang et al. (2021); Xia et al. (2022); Bai and Korwar (2025); Bai and Chen (2024). They can constitute the entire star Haensel et al. (1986); Alcock et al. (1986); Xu et al. (1999); Yang et al. (2023); Zhang (2020); Zhao et al. (2019); Ren and Zhang (2020); Cao et al. (2022); Yuan et al. (2022); Restrepo et al. (2023); Xia et al. (2020); Zhang and Mann (2021); Pretel and Zhang (2024); Zhou et al. (2024); Wang et al. (2024a); Xie and Xia (2025); Chen et al. (2025), forming the so-called quark star that is self-bound featuring a high-density surface. Possible phase transitions can also occur inside the quark stars, leading to self-bound hybrid stars such as hybrid quark stars (HybQSs) from a QM-QM (different QM phases) transition Yang et al. (2025a), or inverted hybrid stars (IHSs) from a speculated quark-hadron phase transition Zhang and Ren (2023); Zhang et al. (2024); Negreiros et al. (2025). Clustered SQM, or in another name Strangeon Matter (SM), like strange quark matter, comprises nearly equal numbers of $u,d,s$ quarks Xu (2003); Lai and Xu (2009, 2017) as favoured by the approximate SU(3) flavour symmetry, and could also be absolutely stable, forming a solid state with localized quark clusters (strangeons) due to their large mass and strong interaction. Strangeon stars Xu (2003); Lai and Xu (2017, 2009); Miao et al. (2022); Gao et al. (2021); Zhang et al. (2023a), composed of SM and thus self-bound, have a stiff equation of state (EOS) and high compactness that helps address various astrophysical observations Zhang et al. (2023a); Li et al. (2023, 2024); Wang et al. (2024b); Yuan et al. (2025). A transition from SM to SQM is possible Miao et al. (2022), forming the so-called Hybrid Strangeon Stars (HybSSs) with a strangeon crust and a SQM core Zhang et al. (2023b).

The stability of compact stars is generally determined by the fundamental ($f$) mode frequency of radial oscillations, whose squared eigenfrequencies’ ($\omega_{0}^{2}$) positive signs indicate the star’s stability against small adiabatic perturbations. Studies of radial pulsations in hybrid stars are typically categorized into two types: slow and rapid transitions, depending on the phase-transition timescale relative to the oscillation timescale. For rapid transitions, the region where $\omega_{0}^{2}\geq 0$ coincides with the conventional stability condition $\partial M/\partial P_{c}\geq 0$ (i.e., star mass $M$ increases with increasing center pressure $P_{c}$), implying that the zero fundamental eigenfrequency point coincides with the maximum mass point Glendenning (1997); Pereira et al. (2018); Pretel et al. (2024). However, this coincidence does not hold for the slow-transition case, where an extended region of stability develops after a sharp phase transition even with $\partial M/\partial P_{c}<0$, corresponding to the so-called slow stable hybrid stars (SSHS) Pereira et al. (2018); Lugones and Grunfeld (2021); Lugones et al. (2023); Goncalves et al. (2022), which can be formed from early phase transitions triggered by quantum or thermal nucleations in realistic astrophysical scenarios such as protoneutron stars or the post-merger hot stars Lugones et al. (2023).

In 2022, Doroshenko et al. Doroshenko et al. (2022) reported puzzlingly low mass-radius values ($M=0.77^{+0.20}_{-0.17}\,M_{\odot}$, $R=10.4^{+0.86}_{-0.78}$ km at $1\sigma$) for the central compact object in the supernova remnant HESS J1731-347, which have prompted extensive studies proposing several explanations such as hybrid star model with early phase transition Tsaloukidis et al. (2023); Sagun et al. (2023); Laskos-Patkos et al. (2024, 2025); Mariani et al. (2024); Gao et al. (2024); Pal et al. (2025); Hong and Ren (2024a), strange quark star model Di Clemente et al. (2024); Horvath et al. (2023); Oikonomou and Moustakidis (2023); Rather et al. (2023); Gholami et al. (2025); Issifu et al. (2025); Kourmpetis et al. (2025), a soft nuclear EOS neutron star Brodie and Haber (2023); Huang et al. (2024); Li and Sedrakian (2023); Kubis et al. (2023), or mixed dark matter models Hong and Ren (2024b); Routaray et al. (2024). Recently, Kini et al. Kini et al. (2024) reported similarly puzzling ultracompact parameters for XTE J1814-338 ($M=1.21^{+0.05}_{-0.05}\,M_{\odot}$, $R=7.0^{+0.4}_{-0.4}\,\rm km$ at $1\sigma$), with proposed explanations including a hybrid star branch Laskos-Patkos and Moustakidis (2025); Zhou and Huang (2025), a bosonic star with nuclear core Pitz and Schaffner-Bielich (2025), a dark matter admixed neutron star Lopes and Issifu (2025) or a strange star admixed with dark matter Yang et al. (2025b); Lopes (2025). Nevertheless, it is very difficult to reconcile both puzzling objects at the same time Veselsky et al. (2025); Laskos-Patkos and Moustakidis (2025); Kourmpetis et al. (2025). Besides, a binary merger event GW190814 was reported Abbott et al. (2020), featuring a primary black hole with mass $23.2^{+1.1}_{-1.0}\,M_{\odot}$, and a secondary companion of $2.59^{+0.08}_{-0.09}\,M_{\odot}$, which is much larger than the upper bound $M_{\rm TOV}\lesssim 2.3M_{\odot}$ of the maximum mass for a non-rotating neutron star, set by various analyses of GW170817 Margalit and Metzger (2017); Rezzolla et al. (2018); Ruiz et al. (2018); Shibata et al. (2019). To the best of our knowledge, no previous studies have managed to build models that could accommodate all these anomalous observations at the same time.

In the absence of a comprehensive model capable of describing all observational mass-radius and tidal deformability measurements, this work examines slow stable self-bound hybrid stars with quite promising results. The general finding that self-bound stars like strange stars can saturate HESS J1731-347 easily Doroshenko et al. (2022), and the slow stable hybrid stars benefit to meet XTE J1814-338 Laskos-Patkos and Moustakidis (2025), naturally motivates us to propose the slow stable self-bound hybrid stars that combine these two advantages. As we will explicitly show in the main text, the new branch that gives rise to the slow stable self-bound hybrid stars helps explain all these anomalous objects at the same time, while also fulfilling other astrophysical observations such as GW170817 Abbott et al. (2017), recent NICER analysis of PSR J0030+0451 Miller et al. (2019); Riley et al. (2019), PSR J0740+6620 Miller et al. (2021); Riley et al. (2021), and the more recent PSR J0437-4715 Choudhury et al. (2024); Rutherford et al. (2024).

This paper is organized as follows: first, we introduce models for the EOSs of three types of self-bound hybrid stars: HybQS, HybSS, and IHS. Then, we consider the mass-radius relations, radial oscillations, and tidal deformabilities with varying parameters, where we look for radially stable solutions that satisfy all the notable astrophysical constraints, particularly aiming to explain HESS J1731-347 and XTE J1814-338 at the same time. We employ the Maxwell construction for phase transitions and geometric units ($G=c=1$) throughout this paper.

To model the HybQS, we use the conventional MIT bag model EOS

for the quark matter composing the crust before the transition pressure $P_{\rm trans}$, with $B$ the effective bag constant characterizing the QCD vacuum difference inside and outside the quark matter bag. We then adopt the constant-sound-speed (CSS) parametrization Alford et al. (2013) for the quark matter core with EOS

for the core quark matter after $P_{\rm trans}$, with $c_{s}^{2}$ denoting the sound speed squared, $\rho_{\rm trans}$ denoting the energy density at $P_{\rm trans}$, and $\Delta\rho$ as the energy density discontinuity of the phase transition interface. Therefore, we have four independent parameters ($B$, $P_{\rm trans}$, $\Delta\rho$, $c_{s}^{2}$) for this HybQS construction.

To model HybSSs, the EOS of strangeon matter Xu (2003); Lai and Xu (2017) located at the outer layer ($P<P_{\rm trans}$) can be written into the following form Zhang et al. (2023a, b)

where $a=A_{6}^{2}/A_{12}=8.4^{2}/6.2\approx 11.38$, $m_{q}$ is the average constituent quark mass, $n$ is the number density of strangeons with $n_{s}$ being the surface baryon number density. The normalized quantities $\bar{n}=N_{q}\,n/n_{s}$ and $\tilde{\epsilon}=\epsilon/N_{q}$ with $\epsilon$ describing the depth of the interaction potential between strangeons, and $N_{\rm q}$ being the number of quarks in a strangeon. A larger $\epsilon$ will then indicate a larger repulsive force at short range and thus maps to a stiffer EOS. The CSS parameterization (Eq. 2) is then used to model the quark matter core ($P>P_{\rm trans})$. Thus, we have five independent parameters ($\tilde{\epsilon}$, $n_{s}$, $P_{\rm trans}$, $\Delta\rho$, $c_{s}^{2}$) for this HybSS construction.

To model IHSs Zhang and Ren (2023); Zhang et al. (2024); Negreiros et al. (2025), we use the simple MIT bag model Eq. (1) for the quark matter crust like the hybrid quark star model, and with the $e$-type polytrope EOS¹¹1Note that the polytope EOS can have different forms depending on whether the directly related variable is energy density ($e$-type) or rest mass density ($\mu$-type). For details, see discussions below Eq. (22) of Ref. Damour and Nagar (2009a).  Damour and Nagar (2009a)

for the hadronic core. The constant pressure at the interface of the two matter phases yields

Therefore, in this IHS model we have four independent parameters ($B$, $P_{\rm trans}$, $\Delta\rho$, $\gamma$). Due to the uncertainties associated with hadronic EOS, one can piecewise the polytrope EOS of different parameter sets at different density regimes Hebeler et al. (2013); Huang and Sourav (2025).

The static spherically symmetric background spacetime has the following line element

where $e^{-2\Psi(r)}=1-2m(r)/r$, with $m(r)$ being a mass function within the hybrid star. To obtain the configuration of relativistic stellar structure, we incorporate the combined EOS described above into the so-called Tolman-Oppenheimer-Volkov (TOV) equations Oppenheimer and Volkoff (1939); Tolman (1939)

with the boundary conditions $\rho(0)=\rho_{\mathrm{c}},\,\Phi(R)=-\Lambda(R)\,,$ where the star’s radius $R$ and mass $M$ are determined by the conditions $P(R)=0$ and $m(R)=M$, respectively.

To investigate the radial stability of our hybrid stellar models, we assume that a fluid element is displaced from its equilibrium position $r$ to the perturbed position $r+\Delta r$, where such a perturbation has a harmonic time dependence $e^{i\omega t}$. The equations for solving infinitesimal radial oscillations $\xi=\Delta r/r$ and the corresponding pressure perturbations $\Delta P$ of a spherical object are Chanmugam (1977); Väth and Chanmugam (1992); Pereira et al. (2018); Gondek et al. (1997); Vásquez Flores and Lugones (2010); Pretel and da Silva (2020)

with the coefficients given by

where $\Gamma=\frac{n}{P}\frac{dP}{dn}=\frac{\rho+P}{P}\frac{dP}{d\rho}$ is the adiabatic index. The initial conditions are $(\Delta P)_{r=0}=-3(\xi\Gamma P)_{r=0}$ and $\xi(0)=1$. We obtain the squared eigenfrequencies $\omega^{2}$ by solving Eq. (8) and Eq. (9) using the shooting method, so that the boundary condition $(\Delta P)_{r=R}=0$ is satisfied at the surface. The fundamental zero mode that determines the radial stability is the lowest-lying one, the frequency of which is commonly denoted as $\omega_{0}$. For rapid transitions, the matching conditions at the interface are $[\Delta P]^{+}_{-}=0$ and $\left[\xi-\frac{\Delta P}{rP^{\prime}}\right]^{+}_{-}=0$, where $[x]^{+}_{-}\equiv x^{+}-x^{-}$ and primes on the variables denote the partial derivatives with respect to $r$. In this context, the stable region determined from these radial oscillation calculations generally coincides with the simple criterion $\partial M/\partial P_{c}>0$ Glendenning (1997); Pereira et al. (2018). For slow transitions, in which we are most interested, the matching conditions across the interface of the two matter phases are $[\xi]^{+}_{-}=0$ and $[\Delta P]^{+}_{-}=0$. The classical stability criterion $\partial M/\partial P_{c}>0$ generally fails in this context, and it is necessary to use the signs of squared eigenfrequencies to determine the stability Pereira et al. (2018); Lugones and Grunfeld (2021); Pretel et al. (2024); Goncalves et al. (2022).

For comparison with gravitational wave observations, the dimensionless tidal deformability $\Lambda=2k_{2}/(3C^{5})$ can be computed, where $C=M/R$ is the compactness and $k_{2}$ is the Love number characterizing the stellar response to external perturbations Hinderer (2008); Hinderer et al. (2010); Postnikov et al. (2010). The Love number $k_{2}$ is determined by solving a differential equation for $y(r)$ Postnikov et al. (2010) concurrently with the TOV equation Eq. (7), using the boundary condition $y(0)=2$. For self-bound hybrid stars, the matching condition $y(r_{d}^{+})-y(r_{d}^{-})=-4\pi r_{d}^{3}\Delta\rho_{d}/(m(r_{d})+4\pi r_{d}^{3}P(r_{d}))$ must be applied at $r_{d}$ (i.e., at both the core radius where $r_{\rm core}=r(P=P_{\rm trans})$ and the star surface $R$), where an energy density jump $\Delta\rho_{d}$ occurs Damour and Nagar (2009b); Takátsy and Kovács (2020).

We aim to explain the puzzling XTE J1814-338 and HESS J1731-347 together with other astrophysical observations utilizing the slow stable branch of self-bound hybrid stars. There are two possible scenarios:

• •

  Scenario 1: GW170817’s constraint needs to be met by the non-hybrid branch (i.e., $\Lambda_{1.4M_{\odot}}<800$ for bare quark stars of HybQS and IHS cases, and for bare strangeon stars of HybSS case).

• •

  Scenario 2: GW170817’s constraint is only met by the hybrid branch (i.e., after $P_{\rm trans}$).

In all scenarios, the HESS J1731-347 constraint can be easily met by all the non-hybrid self-bound branches, while the XTE J1814-338 is conveniently met by the (slow stable) hybrid branches in a viable EOS parameter space.

Studies Abbott et al. (2019a); Radice et al. (2018a); Most et al. (2018); De et al. (2018); Burgio et al. (2018); Shibata et al. (2019) have shown that GW170817’s tidal deformability constraints $\Lambda_{1.4M_{\odot}}<800$ exclude very stiff EOSs and large radius around $1.4\,M_{\odot}$ ($R_{1.4M_{\odot}}$), constraining possible maximum mass $M_{\rm TOV}$ well below that of the massive object ($2.6\,M_{\odot}$) observed in GW190814. Interestingly, Scenario 2 can relax $R_{1.4M_{\odot}}$ and $M_{\rm TOV}$ to large values, and reconcile the tension between the gravitational wave events GW170817 and GW190814 by leaving GW170817’s constraints satisfied only by the hybrid branch, as we will explicitly show later.

First, for Scenario 1, we choose the EOS parameters of the out-layer matter (i.e., bag constant $B$ of quark matter for HybQSs and IHSs, ($\epsilon/N_{q},n_{s}$) of strangeon matter for HybSSs) so that GW170817’s $\Lambda_{1.4M_{\odot}}<800$ bound²²2Note that other constraints from analyses of GW170817 commonly assumes neutron star EOS Abbott et al. (2018, 2019b); Radice et al. (2018b); Kiuchi et al. (2019), so that they do not apply to our study of self-bound stars. For details, see related discussions in Ref. Zhang (2020). is satisfied by the non-hybrid branch Zhang (2020); Zhang et al. (2023a) before the phase transition (i.e., $P_{c}<P_{\rm trans}$). We display the results of Scenario 1 in Fig. 1. The non-hybrid branch is naturally identified as the compact object associated with HESS J1731-347, GW170817, and all the NICER-related PSR J0030+0451, PSR J0740+6620 and the very recent PSR J0437-4715. At the more massive region above $2\,M_{\odot}$, a hybrid branch is developed from a sharp 1st-order phase transition, where the large densities discontinuities at the phase transition interface will render $\partial M/\partial P_{c}<0$, however, still be radially stable in the context of slow phase transition even until the XTE J1814-338’s constrained $M$-$R$ region. Note that for the case of inverted hybrid quark stars (second row of Fig. 1), some of the $\omega_{0}^{2}(P_{c})$ curves that are not terminated at zero frequency because they violate causality ($c_{s}^{2}>1$) earlier before reaching zero $\omega_{0}^{2}$, due to a combined effect of their large polytrope index for the core EOS and the large center pressure extended.

We can see that all three types of self-bound stars can have viable parameter space for their slow stable branch (colored curves) to satisfy the constraints of XTE J1914-338, while the non-hybrid star branch (dot-dashed black curves) before the phase transition to meet all other constraints without the need for much fine-tuning. Besides, it can be seen that the slow stable branch of all hybrid star types manifests the following common features:

• •

  $\omega_{0}^{2}$ increases after the phase transition, with a larger slope and eigenfrequency peak for a smaller $P_{\rm trans}$ (dashed) or a larger $\Delta\rho/{\rho_{\rm trans}}$ (darker) or stiffer core EOS (orange).

• •

  a larger $\Delta\rho/{\rho_{\rm trans}}$ or a stiffer core EOS (i.e., larger $c_{s}^{2}$ for HybQS and HybSS, or larger $\gamma$ for IHS) maps to a more extended slow stable branch.

We expect these features to also apply to conventional slow stable hybrid stars.

Besides, Fig. 1 also shows that HybSSs have some peculiar behaviour compared to other types: for cases of transition near $2\,M_{\odot}$ (solid curves), they develop hybrid branches that cross the XTE J1914-338’s $M$-$R$ region with $\partial M/\partial P_{c}>0$, meaning that they are even stable in the rapid transition context. Furthermore, without invoking Scenario 2, they can already satisfy GW190814’s $M_{\rm TOV}>2.59^{+0.08}_{-0.09}\,M_{\odot}$ constraint, provided a late phase transition (dashed curves).

For Scenario 2 which addresses GW170817’s tidal deformability constraints only with the hybrid branch, we directly choose large $\Delta\rho/\rho_{\rm trans}$, stiff EOS and large $P_{\rm trans}$ to get a more extended slow stable branch, saving us efforts of tuning parameters, benefited from the general features observed in Scenario 1 as discussed above. We display the results in Fig. 2, skipping showing the drawing of $\omega_{0}^{2}(P)$ since they share a similar feature as Scenario 1, but instead put the results of tidal deformability to show the main difference. We can easily identify the viable solutions for HybQSs and HybSSs that satisfy all the constraints, even including the GW190814, benefited from a very stiff non-hybrid branch with a late sharp phase transition, leaving the extended hybrid branch satisfying the tidal deformability bounds of GW170817. Note that the HybSS examples cross the XTE J1914-338’s $M$-$R$ region with $\partial M/\partial P_{c}>0$, thus are even stable in the rapid transition context like Scenario 1. The very recent PSR J0437-4715 constraint is easily met by the slow stable branches of all star types. However, by varying the simple polytrope model and even its piecewised form (i.e., piecewise another one or two polytrope EOSs of different parameter set at higher densities) for the hadronic EOS, we find that the slow stable branches of IHSs that overlaps with GW190814’s large mass region can hardly reach the XTE J1814-338’s $M$-$R$ region, mainly because the polytrope model of their core hadronic matter would easily violate causality ($c_{s}^{2}>1$) in order to make its EOS stiff enough to meet GW190814 while at even much higher center pressure have a long extended stable slow branch to reach XTE J1814-338’s $M$-$R$ region. So we neglect this case in this figure.

Recent astrophysical observations challenge conventional neutron star models: Compact stars identified in HESS J1731-347 ($M\approx 0.77\,M_{\odot}$, $R\approx 10.4\,\rm km$) and XTE J1814-338 ($M\approx 1.21\,M_{\odot}$, $R\approx 7.0\,\rm km$) exhibit anomalously small radii, while GW190814’s secondary component ($M\approx 2.6\,M_{\odot}$) significantly exceeds typical neutron star masses without violating GW170817 constraints from various analyses. In this work, we resolve these tensions by proposing a novel class of self-bound hybrid stars featuring a sharp but slow first-order phase transition, i.e., a large $\Delta\rho/\rho_{\rm trans}$ but with phase transition timescales that exceed oscillation periods. A larger energy density jump ($\Delta\rho/\rho_{\rm trans}$) or stiffer core EOSs extend this slow stable branch to smaller radii, enabling explanations of ultracompact objects.

Our key findings demonstrate that self-bound hybrid stars such as hybrid quark stars (HybQSs), inverted hybrid stars (IHSs), and hybrid strangeon stars (HybSSs) all accommodate HESS J1731-347 and XTE J1814-338 together with other observational constraints via their slow stable branches, and even via the rapid stable branches for HybSSs. In general, the self-bound feature (i.e., mass grows with radius as center pressure increases at the small-mass regime) benefits saturating HESS J1731-347, while the extended stable hybrid branch benefits meeting XTE J1814-338. We also showed the advantage of slow stable self-bound hybrid stars in explaining GW190814’s anomalously massive object by allowing a very stiff EOS of the non-hybrid branch if leaving the slow stable hybrid branches to meet the GW170817 constraint, which turns out to be very feasible for HybQSs and HybSSs. Therefore, our work has demonstrated the first working example that reconciles all the astrophysical tensions related to compact stars’ masses, radii, and tidal deformabilities, revealing new intriguing possibilities for the behavior of matter and stars under extreme conditions.

While the IHS case turns out to have difficulties in meeting GW190814 and XTE J1814-338 at the same time using the simple $e$-type polytrope model and its piecewised version for its core hadronic matter EOS, the chance is still open for other EOS constructions of hadronic matter, such as $\mu$-type polytrope Read et al. (2009), more generalized polytope constructions O’Boyle et al. (2020), spectral representation Lindblom (2010), or the speed-of-sound parameterization Greif et al. (2019). It is also more exhaustive to explore a general Bayesian analysis of all parameter space on this new subject. We leave these for future studies.

There are potential distinct signatures of slow stable self-bound hybrid stars compared to other types of compact stars. The gravitational-wave asteroseismology of this new object may also exhibit interesting features Zhang et al. (2024), particularly their $g$-mode of nonradial oscillations Sun et al. (2025); Guha et al. (2025); Zhao et al. (2025), which may have significantly large frequencies Tonetto and Lugones (2020). In addition, the unusually large density discontinuities at their phase transition interface will likely cause a large violation of neutron star universal relations Han and Steiner (2019); Hoyos et al. (2022); Annala et al. (2018); Ranea-Sandoval et al. (2022); Ranea et al. (2023). Moreover, the strong phase transition is accompanied by a release of energy from the changes of internal energy and gravitational potential energy, potentially sourcing supergiant glitches Ma and Xie (1996), core-quakes Bejger et al. (2005) and electromagnetic wave radiations such as gamma-ray burst Loeb (1998); Menezes et al. (2006) and fast radio burst Shen et al. (2023); Wang et al. (2024b), which may also originate from the self-bound star’s surface Geng et al. (2021); Xu et al. (2025).