Of How Neutron Stars Observations Point (NotWithout Nuance) Towards Exotic Matter: Existing Explanations and a Novel Model

1. Introduction

In the last 15 years, several estimations of the masses and radii of neutron stars (NSs) have become available. Observations of NSs are being made through pulsar timing in binary systems, such as the canonical $2M_{\odot }$ pulsars [1,2,3,4,5]; NICER observations of isolated NS [6,7,8,9,10,11,12] and multimessenger astronomy after GW170817 [13,14,15]. There also exist spectroscopic and photometric observational techniques for pulsars in binary systems, such as the black widow PSR J0952-0607, with a very high mass estimation, $M=2.35\pm 0.17M_{\odot }$ [16]. In addition, two recent estimations of the mass and radius of compact objects have appeared: the very low mass object HESS J1731-347 ($M=0.{77}_{-0.17}^{+0.20}{M}_{\odot }$ and $R=10.4_{-0.78}^{+0.86}$ km) [17] and XTE J1814-338 which has a canonical mass ($M=1.21\pm 0.05M_{\odot }$) but a very small radius ($R=7.0\pm 0.4$ km) [18]. Besides these astronomical observations, other constraints on the equation of state (EOS) of dense matter arise from chiral effective field theory (cEFT) calculations (see, for example, Drischler et al. [19] and references therein) and perturbative QCD (pQCD) (see, for example, Ref. [20] and references therein). All these estimations, when considered simultaneously, produce tension with the theoretical models for cold dense matter used to describe matter in the inner core of NSs.

Regarding the uncertainties over recent observations, we can particularly mention that changes to the atmospheric model used to determine the restrictions for HESS J1731-347 produce a higher value for the mass, thus releasing some of the tension of such low-value. This issue has been already stated in Ref. [17] and also discussed in Ref. [21]. Something similar happens for object XTE J1814-338, where our poor understanding of the bursting mechanism might be introducing systematic errors leading to such an extreme object. In this case, though not preferred, an alternative model -also presented in Ref. [18]- produces more canonical estimations of the mass and radius of this compact star. In relation to the challenging estimations of the parameters of PSR J1231-1411, the lack of a reliable estimation for the distance and less restrictive priors on mass and observer inclination are certainly plotting against the determination of its mass and radius. The PSR J0952-0607 mass measurement also introduces a very high uncertainty with respect to other binary massive NSs mass detections [22]. As we will discuss, the mentioned tension does not necessarily relieve if these constraints upon which a cloak of distrust exists are rectified.

It is important to recall that there are no experimental facilities on Earth capable of studying matter in the density range between cEFT and pQCD calculations for cold dense matter. Moreover, theoretical schemes to construct its EOS from first principles systematically fail in this regime. For these reasons, there are large uncertainties regarding the behavior of matter in this range of densities. In this context, the most promising possibilities for understanding matter at these extreme densities are related to astronomical observations of NSs. In this sense, an important approach to addressing this issue is to produce model-independent studies that allow to obtain a broad set of results to test against modern astronomical data. Such results are also important in the context of future astronomical observations that will become available with, for example, third generation gravitational wave detectors -such as the Einstein Telescope or Cosmic Explorer- (see, for example, Refs. [23,24,25]).

Different models have been proposed to explain the complex scenario of astronomical observations. One of the most explored alternatives is the inclusion of deconfined quark matter inside these compact objects. We must stress that quark degrees of freedom are allowed by Quantum Chromodynamics (QCD), where the theoretical QCD phase diagram predicts the occurrence of a hadron-quark phase transition for an unknown density at $T=0$; there also exist results that suggest that the existence of this kind of matter inside NSs is favored [20]. For this reason, it is fair to include deconfined quark matter in the unknown EOS that describes matter in the inner depths of compact objects. EOSs that allow a phase transition involving quark deconfinement are commonly known as hybrid EOSs, which produce hybrid stars (HSs). Two limiting cases to describe such phase transitions are available in the literature: hybrid EOSs in which a sharp transition, characterized by a gap in the energy density, $\Delta \epsilon$, occurs at a given transition pressure, $P_t$, or phase transitions in which a mixed phase, where hadrons and quarks coexist, is formed. In the first case, the sharp phase transitions is described using the Maxwell formalism (see, for example, Refs. [26,27] for the canonical treatment in which solutions to different Lagrangians are stitched together and Refs. [28,29,30] for unified models). In the second case, the soft phase transitions with a mixed phase are constructed using the (bulk) Gibbs formalism (see, for example, Refs. [31,32] and references therein). Moreover, in this last scenario, geometrical structures can appear as a consequence of the interplay of the Coulomb and nuclear surface energies (see, for example, Refs. [33,34] and references therein). The main physical quantity that determines the nature of the hadron-quark phase transition is the (highly) unknown value of the hadron-quark surface tension (see, for example, Refs. [35,36] and references therein). Moreover, if the sharp hadron-quark phase transition is favored, an additional key aspect needs to be taken into account: the conversion timescale between hadrons and quarks. This feature is of central relevance, as it has been shown that if this timescale is slow compared to the typical values of the fundamental mode of radial perturbations, a new branch of slow stable HSs (SSHS) can exist where the condition $\partial M/\partial {\epsilon }_c>0$ is not satisfied, as proven in Ref. [37]1. On the contrary, if the timescale is rapid, the traditional stability result remains valid, where the stability turning points lie at the critical points $\partial M/\partial {\epsilon }_c=0$.

Several theoretical frameworks -including, but not restricted to, hadron-quark HSs- have been proposed in the literature to account for the observational and theoretical constraints related to compact objects. We present them briefly in order to discuss them in detail in the following sections. The possibility of two coexisting families of compact objects with different nature has been proposed in the works by Drago et al. [44], Di Clemente et al. [45], Drago et al. [46], Drago and Pagliara [47]. This two family scenario explains the high-mass objects as strange quark stars (QSs) and the low-mass and small radii objects as hadronic NSs. Totally stable HSs (i.e. HSs that are stable in both slow and rapid conversion regimes) with sharp interfaces and HSs where a mixed phase is formed are used in the works by Mariani and Lugones [34], Abgaryan et al. [48], Ranea-Sandoval et al. [49], Carlomagno et al. [50], Alvarez-Castillo et al. [51], Pal and Chaudhuri [52]. Within the same theoretical approach, there are also works presenting the impact of the occurrence of sequential phase transitions (a hadron-quark one and a second one between different quark phases) in the inner core of compact stars, as presented in the works by Alford and Sedrakian [53], Rodríguez et al. [54], Li et al. [55]. Within the slow conversion regime, SSHS are good candidates to help explain every constraint (see, for example, the works of Gonçalves and Lazzari [56], Lugones et al. [57], Rau and Salaben [58], Rau and Sedrakian [59], Mariani et al. [60], Laskos-Patkos and Moustakidis [61]). Within this theoretical approach, slow QSs including phase transitions between different quark phases have also been recently proposed by Zhang et al. [62], Yang et al. [63]. In the work by Sagun et al. [64], a summary of different proposals, which include a dark matter (DM) contribution to NSs and HSs, is presented.

Given the current scenario of recent observations and constraints, and the different theoretical proposals, the general objective of this work is to present an updated SSHS approach to explain the current tension among modern constraints and to analyze it in the context of the other proposals existing in the literature. In this context, the first part of this work presents novel results obtained from a model-independent perspective, including a new phenomenological non-CSS parametrization that may be helpful for constructing hybrid EOSs consistent with microscopic results from pQCD. The second part, motivated by the aforementioned tension, reviews several theoretical approaches proposed to account for the full set of observations, discussing their strengths and limitations and comparing them with our own proposal.

The present work is organized in the following manner. In Section 2, we present the hybrid EOSs we use to describe dense matter, including a detailed presentation of the novel non-CSS parametrization for quark matter. The results obtained from implementing these EOSs are presented in Section 3, where we also provide an astrophysical interpretation and compare them with other proposals in the literature. Section 4 is devoted to summarizing our main findings and presenting the conclusions.

2. Hybrid EOS and HS Models

In order to obtain results with a model independent approach, we adopt generic parametric models to construct both the hadron and quark phases of the proposed hybrid EOSs. In this sense, our approach is a continuation of the previous studies presented in Refs. [57,60]. However, in this work, we do not model, as in the previous ones, only hadron-quark (HQ) HSs; we also consider the possibility of quark-quark (QQ) HSs, in which a phase transition between two different quark phases is proposed. QSs with a QQ phase transition have been studied in recent works [62,63], motivated by several studies of the QCD phase diagram that suggest the existence of multiple distinct phases after deconfinement occurs, such as 2SC, 2SC$+s$, CFL, quark-gluon plasma, LOFF phase [65,66,67]. Even more, within the QQ-HS scenario, we explore the possibility of including (or not) a crust in the outermost layer, motivated by the potential existence of a floating crust above the strange matter core [68,69]. In the following paragraphs of this section, we describe in detail the models, prescriptions, and parametrizations implemented in each case.

For the hadron sector, we implement a new parametrization of the generalized piecewise polytropic (GPP) EOS originally proposed by O’Boyle et al. [70]. This hadronic sector, when applied, is used for modeling the outer core region with a number density above $0.5n_0$, where $n_0$ is the nuclear saturation density; for densities below $0.5n_0$, we adopt the BPS-BBP crust [71,72]. Within the GPP formalism, the junction density $0.5n_0$ corresponds to the value ${log}_{10}\left({\rho }_0\right)=14.127$, which represents the starting point of the first piece of the GPP EOS. The other details of the particular selected parametrization are presented in Table 1 and are discussed in the following section along with the results. In addition to the astrophysical motivations discussed later, the BPS-BBP crust and the adopted GPP parametrization have been calibrated to agree with the boundaries of the region permitted by the chiral Effective Field Theory (cEFT) constraint reported in Drischler et al. [19], as depicted in the low-density region of Figure 1 and Figure 3.

For the quark sector, depending on the type of compact object, we use three different models: the MIT Bag, the CSS, and the novel non-CSS one. To model the quark phase for the outer core in the QQ-HS cases, we use the MIT Bag model, which has only a single free parameter, the bag B that models confinement (see, for example, the review presented in Ref. [73] and references therein). When a crust is also added, as we did for the hadronic sector of the EOS, we use the BPS-BBP model up to $0.5n_0$; we call these objects Crust QQ (CQQ) stars. For the inner cores of both the HQ and QQ-HS cases, we implement the Constant Speed of Sound (CSS) parametrization [74]; finally, we also present a novel non-CSS parametric model to describe matter in the inner core of HSs.

Regarding the CSS model, it has been used to capture the generic behavior of HSs with a sharp first-order hadron-quark phase transition through three parameters: the phase transition pressure, $P_t$, the energy density gap between phases, $\Delta \epsilon$, and the squared speed of sound of the quark phase, ${(c_s/c)}^2$. However, when SSHS are taken into account, central energy densities might be as high as those typical of perturbative QCD (pQCD) calculations, $n_{\mathrm{pQCD}}\gtrsim 40n_0$ [20]. In particular, it is expected that the squared speed of sound for the quark sector ${(c_s/c)}^2$ approaches the conformal limit, ${(c_s/c)}^2=1/3$, for pQCD densities. Consequently, within the SSHS scenario, the CSS parameters should be carefully selected. For this reason, within the CSS model, there exists a limitation in exploring the parameter space when studying SSHSs, so as not to violate the pQCD constraint (the pQCD constraint in the P-$\epsilon$ constraint can be seen in the high density region of Figure 1 and Figure 3).

With this motivation, the novel non-CSS model proposal is to keep the $P_t$ and $\Delta \epsilon$ parameters, only modifying the ${(c_s/c)}^2\equiv constant$ condition. We propose a pressure dependent functional form for ${(c_s/c)}^2$ which, for high energy densities, behaves as expected from pQCD calculations:

$${(c_s\left(P\right)/c)}^2={\displaystyle \frac{1}{3}}+\left({(c_t/c)}^2-{\displaystyle \frac{1}{3}}\right)exp\left[-\beta {\displaystyle \frac{P-P_t}{P_t}}\right].$$

(1)

For the limiting cases, it values ${(c_s\left(P_t\right)/c)}^2={(c_t/c)}^2$ at the phase transition, and it satisfies the conformal limit for high pressures, ${lim}_{P\to \infty }{(c_s\left(P\right)/c)}^2=1/3$. The modified parameters are ${(c_t/c)}^2$, the quark speed of sound at the phase transition (which is analogous to ${(c_s/c)}^2$ in the original CSS model), and $\beta$, a dimensionless parameter that regulates the rapidity of the exponential approach to the conformal limit value. This proposal is integrable, having an analytical solution for the relationship $P\left(\epsilon \right)$. It is important to notice that when selecting $\beta =0$, this functional form reduces to the traditional CSS case.

Given all these specific models, we aim to obtain different kinds of HSs. In summary, we plan to construct:

• HQ-HSs, composed of a BPS-BBP crust, a hadronic GPP outer core detailed in Table 1, and an inner core made of quark matter, modeled through both the CCS and non-CSS parametrizations.
• QQ-HSs without a crust composed of an outer quark core modeled with the MIT bag model and an inner core of quark matter described using the CSS parametrization.
• CQQ-HSs, where the BPS-BBP crust outermost layer is added at $0.5n_0$ to the previous QQ-HSs configuration.

In all cases, we assume an abrupt first order phase transition in the inner-outer core interface and a slow conversion scenario between phases, which gives rise to the SSHS family. Using these different proposed configurations, we explore the possibility of satisfying the current astrophysical and microphysical constraints.

3. Results and Discussion

In this section, we present and analyze our results in comparison with other theoretical frameworks. As previously noted, constraints on neutron stars from $2M_{\odot }$ pulsars, gravitational-wave detections, and X-ray observations delineate a scenario that poses significant challenges to contemporary neutron star models.

Given the different EOSs considered and the stellar configurations proposed (HQ, QQ, and CQQ), we explored the corresponding parameter spaces to satisfy the existing observational constraints. Unlike our previous works [57,60], where various hadronic parametrizations were examined as qualitative representatives of different families or limiting cases, in this study we adopt a single, specific GPP parametrization. This choice aligns with the particular objectives of the present work; the discussion and implications of this decision will be addressed later, once the results are presented. The details of the hadronic GPP parametrization are listed in Table 1. Based on this single hadronic parametrization, we have selected several representative cases that serve as suitable scenarios for distinct stellar configurations. In this initial analysis, the non-CSS model was not considered, and only an inner quark core modeled through the CSS case ($\beta =0$) was included. The parameters of the five selected sets are presented in Table 2.

Figure 1 shows the P-$\epsilon$ relations for the five selected sets. In each curve, the constant-pressure plateau represents the energy density discontinuity associated with the first-order phase transition, while the small star symbol indicates the maximum energy density and pressure reached in the terminal mass configuration, the last stable SSHS. The figure also displays the microphysical constraints derived from cEFT and pQCD calculations. The pink dashed segments adjacent to the pQCD region roughly delineate the quadrant where neutron star matter should reside before entering the pQCD regime in order to not be incompatible with it. As mentioned previously, both the crust and the hadronic EOS satisfy the cEFT constraint, whereas the quark phases are not subject to it. Regarding the pQCD constraint, the SSHSs can reach several tens of the nuclear saturation density in their interiors, with the star symbols lying very close to the pQCD region, although the SSHS matter does not fully enter this regime. Thus, the modeled SSHS matter for these sets -though extreme- remains consistent with pQCD calculations, since the corresponding EOSs terminate within the allowed quadrant. Nevertheless, such extreme behavior highlights the need for a precise parameter adjustment to avoid inconsistency with pQCD, further motivating the consideration of the non-CSS model.

In Figure 2, we present the M-R relationships for the five selected sets. Each curve is presented up to the terminal configuration, which is marked with a small star; in this sense, for all sets, the phase transition occurs at the maximum mass configuration, and the long branches after the maximum mass configurations are SSHS branches. With different color regions, we also present the astrophysical constraints; of particular relevance for this work are the four extreme constraints from HESS J1731-347, PSR J0952-0607, PSR J1231-1411, and XTE J1814-338. As previously discussed, although these measurements are complex and, in some cases, subject to debate, they remain viable and represent a significant challenge for theoretical NS models. In line with the objectives of this work, we selected five representative sets to explore different approaches to satisfying the aforementioned observational constraints, or at least a significant subset of them. The HQ curve, in particular, satisfies all constraints, reaching the XTE J1814-338 region only marginally. In this case, the hadronic GPP parameters were tuned to produce a purely hadronic curve that meets all constraints except for XTE J1814-338, while the SSHS configuration accounts for this final measurement. The shape and location of the most recent constraints in this plane is the reason why we are presenting a single new GPP parametrization rather than the multiple enveloping hadron EOS previously proposed in Ref. [60]. The QQ1 case has a low-valued B parameter in the outer core and a very high speed of sound, ${(c_s/c)}^2=0.9$ in the inner core. This very stiff EOS set yields large radii and masses for the stellar configurations along the traditional branch, as well as an extended SSHS branch that satisfies the J0614-3329 and XTE J1814-338 measurements. In contrast, the QQ2 set represents a similar case to QQ1 but with a softer EOS, achieved through a larger B value and a lower speed of sound. This less extreme QQ configuration would be more appropriate if the J1231-1411 measurement gets revised in the future. The inclusion of a crust in the CQQ1 and CQQ2 sets allows these models to satisfy the latter constraint without requiring very small B values or excessively stiff EOSs.

Considering both the P-$\epsilon$ and M-R analyzes, and despite the limitations in evaluating the CSS parameters within the SSHS and pQCD frameworks, all the selected cases indicate that the SSHS scenario is a promising candidate for simultaneously satisfying the current astrophysical and microphysical constraints. Although a precise tuning of the model parameters is required, our results show that more than one family of hybrid stars can yield consistent outcomes. Within the slow-conversion regime, we obtain models of HQ, QQ, and CQQ-HSs capable of addressing the present observational tensions. A common feature among all studied cases is the necessity of adopting a very high speed of sound, ${(c_s/c)}^2\sim 0.8$-$0.9$. Such values are difficult to achieve in many effective microphysical quark EOS models, though repulsive interactions -particularly vector interactions- can, in some cases, provide the necessary stiffness. Nevertheless, as previously discussed, the conformal limit predicted by pQCD at extremely high densities imposes a stringent upper bound on this quantity.

In addition to the CSS cases analyzed, we also investigate the novel non-CSS scenario. For this purpose, we consider a specific HQ-HS configuration that, within the CSS framework ($\beta =0$), satisfies the constraints in the M-R plane but fails to meet the pQCD requirement. For this scenario, we only present HQ-HS configurations, and not the QQ and CQQ alternatives, in order to avoid redundancy since the three alternatives share the same qualitative results regarding the non-CSS effects. In this case, we take, by construction, a particularly extreme set to effectively test the non-CSS proposal:

$$\begin{array}{cc}\hfill P_t& =400{\mathrm{MeV}/\mathrm{fm}}^3,\hfill \\ \hfill \Delta \epsilon & =2500{\mathrm{MeV}/\mathrm{fm}}^3,\hfill \\ \hfill {(c_t/c)}^2& =1.0.\hfill \end{array}$$

Given this set, we explore different EOSs by varying the $\beta$ parameter. In Figure 3 and Figure 4, we present the corresponding results in the P-$\epsilon$ and M-R planes, respectively. In the P-$\epsilon$ diagram, one can observe that increasing the $\beta$ value gradually softens the EOS, and eventually brings it into agreement with the pQCD constraints. In the particular case presented, with $\beta =0.1$, the pQCD constraint is already (marginally) satisfied, while the $\beta =0.4,0.7$ satisfies it much more naturally. Beyond the specific details, it can be observed that, in general, the implementation of the non-CSS approach helps to alleviate the tension with the pQCD regime in the SSHS scenario. However, the M-R plane evidences the shortcoming of this novel method: the results in that plane show that increasing $\beta$ strongly shortens the SSHS branch, making it difficult or even impossible to reach the most extreme XTE J1814-338 measurement. Therefore, although this novel model provides a promising framework for satisfying the conformal limit of the speed of sound while maintaining stiffness in the lower-density regions of the quark EoS, its inherent limitations appear unavoidable as long as the current extreme observational constraints hold.

In what follows, we would like to discuss and compare the proposals and results of this work (along with the previous SSHS results presented in the preceding works by Lugones et al. [57] and Mariani et al. [60]) in the context of other recent theoretical proposals existing in the literature (already mentioned in Section 1). The proposals discussed in this section are not meant to form an exhaustive survey, but instead to provide a qualitative representation of the existing theoretical models. Firstly, there exists the so called two family scenario, proposed originally by Drago et al. [44] and also developed in more recent works [45,46,47]. In this proposal, a hadron NS family and a self-bound QS family would coexist, where the NSs would explain the low mass-low radius objects and the QSs the high radius-high mass objects. The recent work by Shirke et al. [76] also points in the direction of a strange QS family favored by the recent PSR J0614-3329 observation. A comprehensive list of alternatives is presented in the work of Sagun et al. [64], which explores the potential explanations for the measurement of HESS J1731-347. In this work, the authors mentioned the possibility of including different exotic degrees of freedom inside NSs, such as an early HQ phase transition or a DM contribution. On the other side, the recent work by Pal and Chaudhuri [52] studied totally stable HQ-HSs, including both early and late HQ phase transition, considering only NICER observations. Finally, along with the most recent NS detections from XTE J1814-338 and HESS J1731-347, there appear a couple of works that also aim to explain current observations in different SSHS scenarios: Laskos-Patkos and Moustakidis [61] implemented HQ-HSs and Zhang et al. [62] proposed the existence of self-bound QQ-HSs.

In Figure 5, we present the M-R relationships for some of the existing proposals available in the mentioned literature. All the curves presented in the figure are obtained from the articles mentioned in the last paragraph. In this sense, it is important to note that the specific curves chosen to illustrate each proposal are not intended to be fully representative of the complete set of results reported in the corresponding references. Rather, they are selected to provide a graphical illustration of each proposal and to facilitate a qualitative assessment of its overall shape, as well as its potential strengths and limitations. In this framework, we seek to compare and discuss these various proposals in light of our results and the current astronomical constraints.

Firstly, the two family scenario, presented with the blue and cyan pair of curves, appears to be promising since the hadron NS family (blue curve) could be adjusted to also meet the XTE J1814-338 constraint, and the QS family (cyan curve) is able to reach most of the other measurements; however, this proposal does not seem to be fully capable of satisfying the J1231-1411 NICER observation without missing other constraints. The Sagun et al. [64] proposals -comprising an HQ-HS configuration with an early phase transition (orange curve) and an admixed DM-NS model (red curve)- can, through suitable parameter adjustments, reproduce HESS J1731-347 and J1231-1411, but fail to simultaneously explain XTE J1814-338. The Pal and Chaudhuri [52] model of fully stable HQ-HSs with both early (olive curve) and late (dark blue curve) phase transition displays similar strengths and limitations to those of the Sagun et al. [64] framework. The Laskos-Patkos and Moustakidis [61] (green curve when considering the APR hadronic EOS and purple curve for DD2 one) and Zhang et al. [62] models (brown curve) serve as precursors to the results presented here; however, they are not tuned to reproduce the extreme measurements of J1231-1411 and J0952-0607. Furthermore, both of these earlier works exhibit the need for extreme ${(c_s/c)}^2$ values; this previously discussed potential difficulty of the SSHS models should be further verified and analyzed in these works to assess the viability of their proposals. In all the cases presented, although we do not show the corresponding P-$\epsilon$ relations for these models, the fact that none of the hadronic proposals we considered satisfy the cEFT M-R region suggests that their respective EOSs likely fail to meet this constraint. While different cEFT calculations exist -and thus the particular cEFT constraint shown here is not the only one available- we would like to emphasize that a suitable compact star proposal should, at the very least, take one of them into account. Finally, the constraint by Shibata et al. [75], $M_{\max }\le 2.3M_{\odot }$, becomes relevant in this context: there exists a narrow viable mass range between the J0952-0607 measurement and this maximum mass constraint that may warrant some consideration when proposing models.

4. Conclusions and Future Perspectives

In this work we aim to explain the current astrophysical and microphysical constraints on NSs invoking the SSHS scenario in a model-independent manner. We construct, within this hypothesis, different kind of compact objects, as the HQ, QQ and CQQ-HSs and, exploring the parameter space, select some representative sets capable of satisfying the mentioned constraints. We compare our results with other existing NS theoretical proposals of the literature, aiming to review the current state-of-the-art in this area.

Despite obtaining satisfactory results within our models, they show the tension and zero-sum game that the modern astrophysical and microphysical constraints configure currently. In this context, the SSHS proposal, even considering its shortcomings, remains as one of the promising strategies to explain and describe the measured compact objects. However, it should also be reminded that some of the extreme detections, HESS J1731-347, PSR J0952-0607, PSR J1231-1411, and XTE J1814-338, generate some skepticism due to the potential impact of the hypothesis used in each estimation process, and any rectification over its measured values should imply a revision of these conclusions.

To summarize our findings, we present a takeaway of the major conclusions obtained in this work:

• The modern picture of astronomical constraints related to compact objects produce strong tensions, and models with some type of exotic matter seem to be favored. In particular, if the current estimations of either XTE J1814-338 or PSR J1231-1411 are confirmed (or not strongly rectified) by future analysis, the need to include some type of exotic matter in the inner core of compact stars might be strongly favored.
• In accordance to previous proposals presented by Lugones et al. [57] and Mariani et al. [60] (and also in line with alternative scenarios recently presented in Refs. [61,62]), we show that SSHSs could lead to an appropriate description of modern astronomical observations of compact objects (even considering the extreme ones). Despite this being true, large values of the speed of sound are needed, and potential issues with the conformal limit of pQCD might arise for long SSHSs branches.
• Contrary to the traditional CSS model, the novel non-CSS parametrization proposed in this work is useful for avoiding potential problems with pQCD calculations, but it introduces issues, particularly when explaining the challenging XTE J1814-338 observation.
• The analysis of other recent proposals from the literature -including regular hadronic NSs, the two-family scenario, admixed DM HSs and NSs, and QQ-HSs- shows that, while all leave some room for further refinement or updating, none of them is entirely suitable. Whether considered jointly or separately, XTE J1814-338 and PSR J1231-1411 place stringent constraints on these models.
• Despite these limitations, the other proposed hadronic and hybrid models are in tension with cEFT calculations. If these ideas are used in the future, the low pressure region needs to be adjusted.

In order to give closure to this work, we would like to emphasize some potentially relevant observational aspects where attention must be focused in order to shed some light into the behavior of matter under extreme conditions of pressure.

Future observational signaling of a g-mode associated with a sharp first-order phase transition might be a strong evidence to favor the existence of SSHSs as those non-radial oscillation modes are only excited in the slow conversion scenario [77]. For this reason, this might be the most promising observational aspect to discern between purely hadronic objects and hybrid ones; moreover, it could shed some light into the nature of the hadron-quark phase transition [54,78,79,80]. Another quantity to which attention must be put in is the dimensionless tidal deformability, as it has been proven to show extremely different behavior for NSs and SSHSs of a given gravitational mass. This is particularly evident for values $M\lesssim 1.6M_{\odot }$ (see, for example, Refs. [52,57,81,82,83]). Upcoming detections of GWs emitted by isolated (proto-)NSs are expected to provide an even more direct probe of the internal composition of these compact objects. Gravitational-wave asteroseismology can reveal characteristic oscillation modes whose frequencies and damping times depend very sensitively on the internal structure and composition of NSs. Observing these signals could therefore yield decisive evidence regarding the occurrence, nature, and properties of a hadron-quark phase transition in the depths of these compact objects.

The next generation of observational facilities promises to revolutionize our understanding of compact stars and the behavior of ultra-dense matter. The continuing operation and planned upgrades of gravitational-wave detectors, along with the forthcoming third-generation observatories, will enable the detection of a larger population of binary NS mergers and the follow up of post-merger remnants with unprecedented sensitivity (see, for example, Ref. [25] and references therein). Simultaneously, high-precision X-ray timing observations will refine measurements of NSs masses and radii (and potentially restrict also their moment of inertia), providing tighter constraints on the dense matter EOS (see, for example, Refs. [84,85] and references therein). Together, these multi-messenger observations are expected to deliver decisive empirical evidence regarding the highest densities known in the Universe and the existence and nature of a hadron-quark phase transition.