#  Dark Matter or Regular Matter in Neutron Stars?
How to tell the difference from the coalescence of compact objects

Maurício Hippert  Illinois Center for Advanced Studies of the Universe,
Department of Physics, University of Illinois at Urbana-Champaign, Urbana, IL
61801, USA  Emily Dillingham  Department of Physics, Berea College, 101
Chestnut Street, Berea, KY 40404, USA  Hung Tan  Illinois Center for Advanced
Studies of the Universe, Department of Physics, University of Illinois at
Urbana-Champaign, Urbana, IL 61801, USA  David Curtin  Department of Physics,
University of Toronto, Toronto, ON M5S 1A7, Canada  Jacquelyn Noronha-Hostler
Illinois Center for Advanced Studies of the Universe, Department of Physics,
University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA  Nicolás
Yunes  Illinois Center for Advanced Studies of the Universe, Department of
Physics, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA

###### Abstract

The mirror twin Higgs model is a candidate for (strongly-interacting) complex
dark matter, which mirrors SM interactions with heavier quark masses. A
consequence of this model are mirror neutron stars – exotic stars made
entirely of mirror matter, which are significantly smaller than neutron stars
and electromagnetically dark. This makes mergers of two mirror neutron stars
detectable and distinguishable in gravitational wave observations, but can we
observationally distinguish between regular neutron stars and those that may
contain some mirror matter? This is the question we study in this paper,
focusing on two possible realizations of mirror matter coupled to standard
model matter within a compact object: (i) mirror matter captured by a neutron
star and (ii) mirror neutron star-neutron star coalescences. Regarding (i), we
find that (non-rotating) mirror-matter-admixed neutron stars no longer have a
single mass-radius sequence, but rather exist in a two-dimensional mass-radius
plane. Regarding (ii), we find that binary systems with mirror neutron stars
would span a much wider range of chirp masses and completely different binary
Love relations, allowing merger remnants to be very light black holes. The
implications of this are that gravitational wave observations with advanced
LIGO and Virgo, and X-ray observations with NICER, could detect or constrain
the existence of mirror matter through searches with wider model and parameter
priors.

##  I Introduction

Figure 1: Schematic representation of the coalescence between a NS and a MNS.
After inspiralling, radiating energy in GWs (top panel), the two stars collide
and merge (middle panel). Depending on the stellar masses, the system either
relaxes to an admixed mirror star (bottom left panel) or collapses to a small
black hole (lower right panel).

Compact objects present a tantalizing opportunity to study dark matter
candidates. Their strong gravitational fields present a unique environment,
and their properties are coming under increased scrutiny thanks to
gravitational wave (GW) and X-ray observations. Most dark matter candidates
consider one weakly self-interacting species due to simplicity Lisanti (2017);
Lin (2019). One of the consequences of this assumption is that neutron stars
(NSs) can then contain dark matter only through direct (gravitational) capture
or through non-gravitational interaction channels, leading to admixed stars
(i.e. those with a mixture of standard model (SM) matter and dark matter) Li
_et al._ (2012); Mukhopadhyay and Schaffner-Bielich (2016); Baym  _et al._
(2018); Ellis  _et al._ (2018); Dengler  _et al._ (2022); Kain (2021); Sen and
Guha (2021); Jiménez and Fraga (2022); Miao  _et al._ (2022). Admixed NSs can
have a small dark matter core Jiménez and Fraga (2022), or a dark matter halo
that would not affect the visible radius Miao  _et al._ (2022); sometimes,
such stars can be composed almost entirely of dark matter, leading to
extremely tiny dark compact objects (with masses of $\sim
10^{-8}M_{\odot}$)Horowitz and Reddy (2019). Alternatively, some studies have
looked into NSs composed of hidden sector nucleons from a dark-QCD sector
Maedan (2020), fundamental asymmetric dark matter fermions Narain  _et al._
(2006); Goldman  _et al._ (2013); Tolos  _et al._ (2015); Kouvaris and Nielsen
(2015); Mukhopadhyay  _et al._ (2017); Gresham and Zurek (2019); Jiménez and
Fraga (2022); Collier  _et al._ (2022); Das  _et al._ (2022); Leung  _et al._
(2022), or asymmetric bosons Rutherford  _et al._ (2022); Giangrandi  _et al._
(2022). These studies generally consider simplified interactions of only 1 or
2 particle species. An extensive review of dark matter capture within NSs can
be found in Deliyergiyev  _et al._ (2019).

In recent years, however, the possibility of complex, strongly self-
interacting dark matter candidates has been suggested Alexander  _et al._
(2016); Curtin and Setford (2020, 2021) and the consequence for compact
objects explored Hippert  _et al._ (2022); Ryan and Radice (2022).
Specifically, mirror matter Blinnikov and Khlopov (1982, 1983); Khlopov  _et
al._ (1991); Berezhiani  _et al._ (1996); Chacko  _et al._ (2006) within the
mirror Twin Higgs model Chacko  _et al._ (2006, 2017, 2018); Curtin and
Setford (2020, 2021) is a nearly identical copy of the SM in its matter
content and gauge interaction, except that the masses of fundamental particles
are scaled up by a factor $f/v$, where $f$ and $v$ are vacuum expectation
values of SM and mirror-sector Higgs fields, respectively. This model is
highly motivated as a solution to the hierarchy problem for $f/v\sim 3-7$. The
twin top quarks are neutral under all SM gauge interactions, but their
interaction with the SM Higgs nevertheless stabilizes its mass and solves the
little hierarchy problem without predicting large signals at the LHC Burdman
_et al._ (2015); Craig  _et al._ (2015). This is in contrast to other
solutions to the hierarchy problem like TeV-scale supersymmetry, which
predicts large LHC signals due to strong production cross sections for new
particles, like stops and gluinos Aad  _et al._ (2021). Additionally, because
mirror matter contains multiple species that are strongly interacting, it is
possible for mirror matter to clump together to form mirror neutron stars
(MNSs) Hippert  _et al._ (2022). MNSs were found to be very similar to SM NSs,
except that they are significantly smaller, with masses
$M\sim(0.5-1)M_{\odot}$ and radii $R\sim(4-8)$ km Hippert  _et al._ (2022).
These MNSs are entirely new hypothesized, electromagnetically dark, compact
objects with completely different mass-radius sequences from SM NSs.
Importantly, since MNSs are entirely electromagnetically dark GWs from MNS
mergers are the best method to detect them. Their distinct tidal deformability
and mass range makes them distinguishable from SM NSs or SM stellar-mass black
holes.

While mergers between MNSs are very distinct from standard astrophysical
gravitational wave signals and represent a spectacular discovery opportunity
for new physics, the precise abundance and distribution of MNSs, and hence the
rate of their mergers, is almost impossible to predict for a given
microphysical model of mirror matter. Depending on the mirror matter
distribution in our universe, other astrophysical events involving mirror
matter may be much more common. It is therefore vital to consider the full
range of possible phenomena involving mirror matter and neutron stars.

What would happen if a NS and a MNS were to coalesce? Schematically, the
evolution of such a binary is depicted in Fig. 1. First, the two compact
objects would spiral around one another. While their separation is large, they
could be effectively treated as point particles in the post-Newtonian
approximation Blanchet (2014). Because of GW emission, their orbit would
shrink and eventually their separation would become small enough that the
stars would tidally deform each other. As their separation continues to
decrease, the surface of the NS and that of the MNS would cross each other and
their stellar interiors would begin to overlap. Unlike mergers of two NSs or
two MNSs, this would not be a standard collision because mirror matter and
regular matter only interact gravitationally.111Mirror matter and SM matter do
interact via Higgs exchange, but the effect is negligible for astrophysical
processes. After the MNS sinks to the center of the NS, the final remnant
would either be a black hole or a stable NS with a mirror matter core,
depending on the mass of the progenitor stars.

The development of a quantitative description of the coalescence sketched
above requires the combination of various techniques. In the inspiral, post-
Newtonian theory can be used to describe the orbital motion of the stars and
their tidal deformations. Indeed, the previous study Hippert  _et al._ (2022)
has focused on this inspiral stage and the possibility of using measurements
of the tidal deformabilities to distinguish between NS binaries and MNS
binaries. In this work we extend the inspiral analysis to NS-MNS binaries.

The merger stage can only be described via a two-fluid, numerical relativity
simulation in which the fluids only interact gravitationally. However, at this
point in time, such a numerical relativity code does not exist. Even without
such a code, however, one can still study the hypothetical remnant object: a
“mirror-matter admixed NS” (MANS), partially composed of SM matter and mirror
matter. There are three possibilities to create MANS : (i) mirror matter
accreting into the cores of a SM NS, (ii) SM matter accreting into the cores
of a MNS, and (iii) MNSs and NSs merging into a stable admixed mirror star
remnant. In case (i) [(ii)], one would have very tiny mirror matter [SM] cores
that depend on the available mirror matter [SM matter] from their surroundings
and the precise nature of the matter-mirror and matter interactions (see e.g.
Curtin and Setford (2021) for a study of Twin Higgs mirror matter accretion in
white dwarfs). In case (iii), most NS and MNS mergers would produce light
black holes, but in certain rare cases, stable admixed mirror star remnants
would remain with large mirror matter cores.

In this paper, we study the properties of objects formed in both formation
scenarios of MANSs and address their observational prospects. We essentially
ask ourselves the following question: how can we distinguish between SM NSs,
MNSs, and MANSs from GW observations with advanced LIGO and Virgo, or X-ray
observations with NICER? To answer this question, we first construct two-fluid
solutions to the Tolman-Oppenheimer-Volkoff (TOV) equations, which describe
non-rotating stars (either SM NSs with $f/v=1$, or MNSs or MANSs with $f/v>1$)
in general relativity Goldman  _et al._ (2013); Kain (2020); Naidu  _et al._
(2021); Gleason  _et al._ (2022). Our primary findings are the following:

  * •

Given an equation of state (EoS), the mass-radius sequence of stable NSs is
continuously connected to the mass-radius sequence of stable MNSs through
mass-radius sequences of stable MANS. This means that MANSs have a two-
dimensional (2D) mass-radius plane, instead of a one-dimensional mass-radius
curve, given a SM EoS.

  * •

MANSs can have the ultimate _mass-radius twins_ , i.e. stars with the exact
same mass and radius, but different DM fraction, and thus, composition. These
ultimate twins, however, can be told apart by their tidal deformabilities,
$\Lambda$.

  * •

The mass-radius plane of stable MANSs is nearly independent of the dark
sector’s $\Lambda_{QCD}$ and extends to much lighter (smaller mass) and much
tinier (smaller radius) objects, relative to stable NSs above the
Chandresekhar limit.

  * •

The relation between the tidal deformability and the mass (or the compactness)
of MANS is also two-dimensional. This plane extends to a minimum $\Lambda\sim
10$, but not just for very massive objects, but rather also for light objects
with masses as low as $\sim 0.75M_{\odot}$.

  * •

The inspiral of a NS and a MNS and that of a MNS and another MNS fills a
unique phase space of symmetric mass ratio and chirp mass, unreachable by NS
binaries. The symmetric mass ratio can be as low as $0.16$ (corresponding to a
mass ratio of $0.3$), while the chirp mass can be as low as $0.5M_{\odot}$
(corresponding to a total mass of $1M_{\odot}$).

  * •

The remnant of the merger between a NS and a MNS primarily produces very light
black holes, with masses between $\sim(1.3,3.7)M_{\odot}$. Occasionally, such
a merger may also produce a stable MANS with a mass of $M_{\textrm{\mbox{\tiny
MANS}}}\sim(1.2,1.6)M_{\odot}$ and radius $R_{\rm{\rm MANS}}\sim(8,9.5)$ km.
The merger of a MNS binary can lead to black holes as small as $M\sim
0.8M_{\odot}$.

  * •

If a NS and MNS merge and both are above their respective Chandrasekhar
limits222The NS Chandrasekhar limit is about $M\sim 1.4M_{\odot}$. There have
already been observations of objects below that limit Strobel  _et al._
(1999); Doroshenko  _et al._ (2022) with X-Ray observation. From those X-Ray
observations, one can estimate the mass and sometimes the radius as well.
However, time-integrated measures have been found to be susceptible to
systematic error (especially due to assumptions about the atmospheres) such
that radii can differ by up to $\sim 50$% Catuneanu  _et al._ (2013);
Servillat  _et al._ (2012). , the admixed core of the remnant MANS can have
masses up to $M_{\textrm{\mbox{\tiny DM}}}\sim 0.4M_{\odot}$ and radii as
large as $4.5$ km. Thus, these cores can make up approximately one third of
the star’s mass and are distributed up to half of their total radius.

  * •

If MANS are formed by other means (e.g. accreting matter), then a much wider
2D mass-radius sequence is possible that could lead to either DM halos or
cores.

The results summarized above differ from other recent studies on mirror matter
admixture on NSs Sandin and Ciarcelluti (2009); Ciarcelluti and Sandin (2011);
Kain (2021) and binary NS mergers Emma  _et al._ (2022). These other studies
explored mirror matter candidates that, apart from discrete symmetries, are
identical to SM particles. Our work is unique in that we consider mirror Twin
Higgs matter, with fundamental masses scaled up by $f/v>1$, and we employ
first-principle results Hanhart  _et al._ (2008); McNeile  _et al._ (2009);
Walker-Loud  _et al._ (2009); Bratt  _et al._ (2010); Syritsyn  _et al._
(2010); Albaladejo and Oller (2012); Horsley  _et al._ (2014), from lattice
QCD and chiral perturbation theory, to scale the EoS with $f/v$ Hippert  _et
al._ (2022). For $f/v=1$, our model EoS is also tuned to meet current
observational constraints on the NS mass-radius relation, while mirror matter
is considered in the range $f/v\gtrsim 3$ favoured by LHC constraints Burdman
_et al._ (2015); Craig  _et al._ (2015). This makes our EoS models for the
dark matter sector directly relevant for the solution to the Hierarchy Problem
as well as highly realistic.

The remainder of this paper presents the details that led to the conclusions
summarized above, and it is organized as follows. Section II briefly reviews
our microscopic model for the EoS of NS matter developed in Hippert  _et al._
(2022), and how we extend it to the mirror sector. Section III outlines
general methods for the calculation of two-fluid stellar properties and their
stability. Section IV presents results regarding the properties of isolated
NSs, MNSs and MANS. Section V investigates and discusses the GW imprint of a
NS/MNS inspiral. Section VI describes the possible remnants of a NS/MNS
merger. Section VII summarizes the signatures of mirror from various merger
scenarios, and Sec. VIII concludes and presents an outlook of future work.

##  II Microscopic Model

In this short section, we briefly describe our microscopic model for the EoS
of NSs and MNSs Hippert  _et al._ (2022). All of the presented results are
derived under the simplifying approximations of vanishing temperature and slow
rotation, $T\simeq 0$ and $\Omega\simeq 0$.

For the description of matter in the NS core, we employ a relativistic mean-
field nuclear model, suplemented by a free gas of electrons and muons, so that
charge neutrality and chemical equilibrium under weak interactions are
satisfied. This model features scalar, vector, and vector-isovector
interactions and is tuned to reproduce properties of nuclear matter,333The
nuclear physics inputs of our model are the incompressibility, symmetry
energy, nucleon Dirac mass and binding energy per nucleon at saturation
density. The set of values that are employed here and the corresponding
parameters can be found in Hippert  _et al._ (2022). as well as current
constraints from NS observations Miller  _et al._ (2019); Riley  _et al._
(2019); Miller  _et al._ (2021); Riley  _et al._ (2021); Abbott  _et al._
(2017, 2018a, 2019). At lower densities, we match the model to a Baym-Pethick-
Sutherland EoS for the NS crust Baym  _et al._ (1971) via an interpolation
procedure. A detailed account of our model can be found in Hippert  _et al._
(2022) and references therein. Results for NS properties obtained with this
model are reviewed in Sec. IV.1.

For our description of matter inside a MNS, we adopt the same model, but
rescale parameters to account for the different microphysics of the Twin Higgs
mirror sector. Our dark-matter model is then fixed by the ratio $f/v$ between
the expectation values of the Higgs and Mirror Higgs. We assume this ratio to
lie in the range $f/v\simeq(3-7)$, as larger values would fail to naturally
solve the hierarchy problem in the minimal model Burdman  _et al._ (2015), and
smaller values are ruled out by Higgs coupling measurements at the LHC Aad
_et al._ (2020).

All the elementary masses in our model —that is, the ones of leptons and
quarks— are thus rescaled by the same factor $f/v$, defined such that the SM
is recovered when $f/v=1$. However, because of the non-perturbative scale
$\Lambda_{\textrm{QCD}}$ and the spontaneous breakdown of chiral symmetry of
QCD, this translates to non-trivial scalings for the nuclear-model couplings
and the baryon mass with $f/v$. These scalings are extracted from low-energy
nuclear-physics phenomenology and first-principle results from chiral
perturbation and lattice QCD (see Hippert  _et al._ (2022) for extensive
details and suggestions for new calculations from lattice QCD). In the absence
of these inputs, extra couplings are scaled according to dimensional analysis,
but for those couplings the details of this procedure have minimal impact on
the overall EoS scaling. Details of these scalings can be found in Hippert
_et al._ (2022).

##  III Stellar Structure Equations

In this section, we review the equations that we use to calculate the
structure of MNSs. We also review simple criteria for the stability of these
stars and the calculation of their Love number. These results are very well
known in the case of ordinary compact stars (see e.g. Yagi and Yunes (2013)
and references therein), but here we present them with an extra layer of
generality Li  _et al._ (2012); Mukhopadhyay and Schaffner-Bielich (2016);
Baym  _et al._ (2018); Ellis  _et al._ (2018); Nelson  _et al._ (2019);
Dengler  _et al._ (2022); Kain (2021); Sen and Guha (2021); Jiménez and Fraga
(2022); Miao  _et al._ (2022).

Besides the single component stars already discussed in Refs. Hippert  _et
al._ (2022), we address a third possibility, consisting of two-component
isolated stars, made of both SM and mirror matter. This choice is motivated
not only by completeness but by the possibility of such stars existing, as a
result of dark-matter admixture in NSs. We refer the reader to Sec. IV.2 for
details. In the case of such stars, because interactions between ordinary and
dark matter are assumed to be negligible, the structure equations, stability
analysis and calculation of Love number must be modified accordingly, as we
present below Dengler  _et al._ (2022).

###  III.1 Two-Fluid Configurations

The structure of (mirror) NSs follows from the equilibrium between the push of
pressure and the pull of gravity, summarized through the equation of
hydrostatic equilibrium. In the case of two-component stars, however, the lack
of significant interactions between dark and visible matter prevents their
equilibration with one another.444It is interesting to consider the impact of
small possible couplings, like a kinetic mixing between the SM photon and the
twin photon, which may result in non-negligible interactions between matter
and mirror matter in this case Curtin and Setford (2021, 2020, 2020), but we
leave this for future work. . Instead, each fluid component must independently
neutralize gravitational forces, resulting in two independent “equilibrium”
pressures, say, $p_{\textrm{SM}}$ and $p_{\textrm{\mbox{\tiny DM}}}$.
Nonetheless, the two fluids interact via gravitation, and the Einstein
equations must be solved in the presence of both fluids at the same time
Goldman  _et al._ (2013); Li  _et al._ (2012).

Two-component stars are odd creatures, and possess two zero-pressure surfaces,
which encompass the two independent fluids. In principle, these fluids can
rotate with different angular velocities without violating hydrostatic
equilibrium. Nonetheless, after sufficiently long times, tidal effects tend to
synchronize their rotations. In the following sections, for the sake
simplicity, we assume no differential rotation between the two fluid
components. Under this assumption, the two-fluid system we consider is no
different from a single-fluid system, except that there are two stress-energy
tensors on the right-hand side of Einstein’s field equation.

Assuming two independent fluids amounts then to writing the stress-energy
tensor as the linear combination

$T^{\mu\nu}=T_{\textrm{\mbox{\tiny SM}}}^{\mu\nu}+T_{\textrm{\mbox{\tiny
DM}}}^{\mu\nu}\,,$ (1)

where each component is independently conserved, i.e.

${T_{\textrm{\mbox{\tiny SM}}}^{\mu\nu}}_{;\nu}={T_{\textrm{\mbox{\tiny
DM}}}^{\mu\nu}}_{;\nu}=0$ (2)

Also, assuming global thermodynamic equilibrium, the stress-energy tensor for
each component can be written as

$\begin{split}T_{\textrm{\mbox{\tiny
SM}}}^{\mu\nu}&=(\epsilon_{\textrm{\mbox{\tiny SM}}}+p_{\textrm{\mbox{\tiny
SM}}})u_{\textrm{\mbox{\tiny SM}}}^{\mu}u_{\textrm{\mbox{\tiny
SM}}}^{\nu}+p_{\textrm{\mbox{\tiny SM}}}g^{\mu\nu}\\\ T_{\textrm{\mbox{\tiny
DM}}}^{\mu\nu}&=(\epsilon_{\textrm{\mbox{\tiny DM}}}+p_{\textrm{\mbox{\tiny
DM}}})u_{\textrm{\mbox{\tiny DM}}}^{\mu}u_{\textrm{\mbox{\tiny
DM}}}^{\nu}+p_{\textrm{\mbox{\tiny DM}}}g^{\mu\nu}.\end{split}$ (3)

From the above equations, it is clear that the standard case of a single fluid
that leads to a SM-only NS can be recovered by taking $p_{\textrm{\mbox{\tiny
DM}}}=0=\epsilon_{\textrm{\mbox{\tiny DM}}}$. Thus, from here on, we
concentrate on calculations for the more general case in which two fluids are
present.

###  III.2 Hydrostatic Equilibrium

Let us review the hydrostatic equilibrium equations in the case of two
independent fluids, relevant for two-component stars with admixed cores
Goldman  _et al._ (2013); Li  _et al._ (2012). These equations generalize the
well-known Tolman-Oppenheimer-Volkoff (TOV) equations for single-fluid stars,
which can be recovered as a particular case, by setting to zero the density
and pressure of one of the fluids.

We consider a static two-fluid hybrid star in spherical symmetry. The line
element squared can therefore be written as
$ds^{2}=-e^{2\tau}dt^{2}+e^{2\sigma}dr^{2}+r^{2}d\Omega^{2}$, with signature
$(-,+,+,+)$, where $\tau$ and $\sigma$ are functions of the (areal) radial
coordinate $r$ and $d\Omega^{2}$ is the line element squared of a two-sphere.
The only non-zero component of the four velocity of a static-fluid is the time
component. Therefore, $u_{\textrm{\mbox{\tiny
SM}}}^{\mu}=u_{\textrm{\mbox{\tiny DM}}}^{\mu}=\mathcal{N}(1,0,0,0)$, where
$\mathcal{N}$ is a normalization factor that ensures $u_{\mu}u^{\mu}=-1$.

With the four velocity known and a few algebraic manipulations, Eqs. (2) and
the Einstein field equation simplify to a set of TOV-like equations for two
fluids:

$\displaystyle\frac{dp_{i}}{dr}=-\left(\epsilon_{i}+p_{i}\right)\frac{m+4\pi
r^{3}\left(p_{\textrm{\mbox{\tiny SM}}}+p_{\textrm{\mbox{\tiny
DM}}}\right)}{r\left(r-2m\right)},$ (4a)
$\displaystyle\frac{dm}{dr}=4\pi\left(\epsilon_{\textrm{\mbox{\tiny
SM}}}+\epsilon_{\textrm{\mbox{\tiny DM}}}\right)r^{2},$ (4b)
$\displaystyle\frac{dN_{i}}{dr}=4\pi\frac{n_{i}}{\sqrt{1-{2m}/{r}}}r^{2},$
(4c)

with $i=\\{\textrm{SM},\textrm{DM}\\}$. Here, $p_{\textrm{\mbox{\tiny
SM}},\textrm{\mbox{\tiny DM}}}$ and $\epsilon_{\textrm{\mbox{\tiny
SM}},\textrm{\mbox{\tiny DM}}}$ are the pressure and energy density of the SM
and DM fluids, respectively, while $m=[1-\exp(-\lambda)]/2$ is the enclosed
mass function. Notice that Eq. (4a) governs the hydrostatic equilibrium of the
fluids separately. In Eq. (4c), $n_{\textrm{\mbox{\tiny
SM}},\textrm{\mbox{\tiny DM}}}$ and $N_{\textrm{\mbox{\tiny
SM}},\textrm{\mbox{\tiny DM}}}$ are the baryon number density and the total
number of baryons and mirror baryons, respectively.

We follow Yagi and Yunes (2013) closely for the standard numerical
implementation to solve the two-fluid TOV equation, except for the initial
conditions and stopping conditions. For the initial conditions, we need to
specify the energy densities of both the SM species and the DM species,
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ and
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$, at the center of the star. Let us
say then that we choose some value of $\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}$ and $\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ and begin to integrate
the equations out from the center. At some radial coordinate, the pressure of
one of the fluids (either the SM or the DM fluid) drops to $10^{-8}$ of the
corresponding central pressure. That value of the radial coordinate defines
the “radius of the admixed core $R_{\textrm{\mbox{in}}}$,” which contains both
fluids. Of course, that is not all of the star, since the density and pressure
of the other fluid can still be large at that radius. The integrator then
continues, with the energy density and pressure of the first fluid set to
zero, until the pressure of the second remaining fluid drops to $10^{-8}$ of
the corresponding central pressure. The radial coordinate value at that place
defines the outer radius of the entire star $R_{\textrm{\mbox{out}}}$ and the
integrator stops. The total mass $M_{\textrm{\mbox{\tiny MANS}}}$ of the star
is then simply the value of $m$ evaluated at the radius
$R_{\textrm{\mbox{out}}}$.

The description of this implementation presented above anticipates one of the
main results of this paper: stable two-fluid stars have (2D) mass-radius
planes instead of (1D) mass-radius curves. Mathematically, this tracks back to
the need to specify two initial conditions, $\epsilon_{c}^{\textrm{\mbox{\tiny
DM}}}$ and $\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$, to find a stellar
solution to the two-fluid TOV equations. This means that one can independently
change $\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ while keeping
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ constant (thus, increasing the
number of mirror baryons in the star), or alternatively independently change
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ while keeping
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ constant (thus, increasing the
number of SM baryons in the star). As we will show in the next section, the
sequence of stable stellar configurations is, indeed, represented by a plane
in mass and radius.

###  III.3 Tidal Deformability

Tidal deformation emerges under strong external tidal fields. If a MANS is in
a binary orbit with another compact object, at some point the MANS and the
other compact object will be a distance $r_{12}$ appart such that the MANS
will become tidally deformed due to the gravitational pull of the other
compact object. In this case, the tidal deformability can be calculated within
this two fluid approach as well, see Rafiei Karkevandi  _et al._ (2021);
Dengler  _et al._ (2022). Just as SM NSs, the two-fluid MANS can be tidally
deformed when the distance $r_{12}$ between the stars is small enough, i.e.
$R_{\textrm{\mbox{out}}}\ll r_{12}\ll R_{\rm ext}$, where $R_{\rm ext}$ is the
radius of curvature of the source of the external field (the other compact
object) and $R_{\textrm{\mbox{out}}}$ is the radius of the entire two-fluid
MANS. A multipolar expansion of the metric exterior to the star can be written
as Hinderer (2008):

$\begin{split}-\frac{1+g_{tt}}{2}=&-\frac{M_{\textrm{\mbox{\tiny
MANS}}}}{R_{\textrm{\mbox{out}}}}-\frac{3Q_{ij}^{(tid)}}{2R_{\textrm{\mbox{out}}}^{3}}\left(n^{i}n^{j}-\frac{1}{3}\delta^{ij}\right)\\\
&+{\cal{O}}\left(\frac{1}{R_{\textrm{\mbox{out}}}^{4}}\right)+\frac{1}{2}{\cal{E}}_{ij}R_{\textrm{\mbox{out}}}^{2}n^{i}n^{j}+{\cal{O}}\left(R_{\textrm{\mbox{out}}}^{3}\right)\,,\end{split}$
(5)

where ${\cal{E}}_{ij}$ is the (quadrupole) tidal tensor field,
$Q^{(tid)}_{ij}$ is the corresponding tidally-induced and traceless quadrupole
moment tensor, and $n^{i}=x^{i}/R_{\textrm{\mbox{out}}}$ is a field-point unit
vector. The tidal deformability $\lambda$ can be defined by
$Q^{(tid)}_{ij}=-\lambda\,{\cal{E}}_{ij}$. The dimensionless tidal
deformability can then be written as

$\Lambda=\lambda/M_{\textrm{\mbox{\tiny MANS}}}^{5},$ (6)

where $M_{\textrm{\mbox{\tiny MANS}}}$ is the mass of the entire two-fluid
star obtained from the two-fluid TOV equations.

The quantity $\lambda$ can be calculated by solving the Einstein field
equations interior to the star at first order in the tidal perturbations. The
version of this equation for a single fluid star is presented e.g. in Yagi and
Yunes (2013). To obtain the corresponding equation for a two-fluid star, one
may naively replace $\epsilon\to\epsilon_{\textrm{\mbox{\tiny
SM}}}+\epsilon_{\textrm{\mbox{\tiny DM}}}$, $p\to p_{\textrm{\mbox{\tiny
SM}}}+p_{\textrm{\mbox{\tiny DM}}}$, and $m\to m_{\textrm{\mbox{\tiny
SM}}}+m_{\textrm{\mbox{\tiny DM}}}$ in the relevant equations, where $m_{i}$
can be obtained using the equation $dm_{i}/dr=4\pi\epsilon_{i}r^{2}$ with
$i\in\\{\textrm{\mbox{\tiny SM}},\textrm{\mbox{\tiny DM}}\\}$. However,
because both fluids exert gravitational attraction on each other, one must
employ Eq. (4a) to write $d\epsilon_{i}/dr=(d\epsilon_{i}/dp_{i})(dp_{i}/dr)$,
instead of applying a naive replacement in the corresponding term for one-
fluid stars.

In the end, the Einstein field equations at first order in the tidal
perturbation become

$r\frac{dy_{*}}{dr}+y_{*}^{2}+F(r)y_{*}+r^{2}Q(r)=0$ (7)

with

$F(r)=\frac{r+4\pi r^{3}(p_{\textrm{\mbox{\tiny SM}}}+p_{\textrm{\mbox{\tiny
DM}}}-\epsilon_{\textrm{\mbox{\tiny SM}}}-\epsilon_{\textrm{\mbox{\tiny
DM}}})}{r-2m_{\textrm{\mbox{\tiny SM}}}-2m_{\textrm{\mbox{\tiny DM}}}}$ (8)

and

$\begin{split}Q(r)&=\frac{4\pi r}{r-2m_{\textrm{\mbox{\tiny
SM}}}-2m_{\textrm{\mbox{\tiny DM}}}}\bigg{[}5(\epsilon_{\textrm{\mbox{\tiny
SM}}}+\epsilon_{\textrm{\mbox{\tiny DM}}})+9(p_{\textrm{\mbox{\tiny
SM}}}+p_{\textrm{\mbox{\tiny DM}}})+\\\ &\left.\frac{p_{\textrm{\mbox{\tiny
SM}}}+\epsilon_{\textrm{\mbox{\tiny SM}}}}{{c_{s}}_{\textrm{\mbox{\tiny
SM}}}^{2}}+\frac{p_{\textrm{\mbox{\tiny DM}}}+\epsilon_{\textrm{\mbox{\tiny
DM}}}}{{c_{s}}_{\textrm{\mbox{\tiny DM}}}^{2}}-\frac{3}{2\pi r^{2}}\right]+\\\
&-\left(2\frac{m_{\textrm{\mbox{\tiny SM}}}+m_{\textrm{\mbox{\tiny
DM}}}+4(p_{\textrm{\mbox{\tiny SM}}}+p_{\textrm{\mbox{\tiny DM}}})\pi
r^{3}}{r(r-2m_{\textrm{\mbox{\tiny SM}}}-2m_{\textrm{\mbox{\tiny
DM}}})}\right)^{2},\end{split}$ (9)

which agrees with Dengler  _et al._ (2022). Here the quantity $y_{*}$ is
defined by $y_{*}\equiv rh_{2}^{\prime}(r)/h_{2}(r)$, where $h_{2}$ is related
to time-time component of the metric by
$g_{tt}=-e^{2\tau}(1+h_{2}Y_{2m}(\theta,\phi))$, with $Y_{2m}$ the $\ell=2$
spherical harmonic. Notice that the quantity $h_{2}$ represents a first-order-
in-perturbation term due to tidal effects. The quantity ${c_{s}}_{i}^{2}$ in
Eq. (9) is the speed of sound squared of the $i_{th}$ fluid555Note that the
$Q(r)$ term depends on $1/c_{s1}^{2}$ and $1/c_{s2}^{2}$, which can diverge
during a first-order phase transition. In this radial regime, however, the
pressure and density drops violently to keep the ratio finite., which is
defined as ${c_{s}}_{i}^{2}\equiv dp_{i}/d\epsilon_{i}$. After solving the
interior equation, Eq. (7), one can match the interior solution to the
exterior metric of Eq. (5) to find the ratio of coefficients $Q^{tid}_{ij}$
and ${\cal{E}}_{ij}$, and thus, to obtain $\Lambda$. Doing so, one finds

$\begin{split}\Lambda&=\frac{16}{15}\left\\{\left(1-2{C}\right)^{2}\left[2-y+2{C}\left(y-1\right)\right]\right\\}\left\\{\right.\\\
&2{C}\left[6-3y+3{C}(5y-8)\right]\\\
&+4{C}^{3}\left[13-11y+{C}(3y-2)+2{C}^{2}(1+y)\right]\\\
&\left.+3(1-2{C})^{2}\left[2-y+2{C}(y-1)\right]\ln\left(1-2{C}\right)\right\\}^{-1}\,,\end{split}$
(10)

where $C=M_{\textrm{\mbox{\tiny MANS}}}/R_{\textrm{\mbox{out}}}$ is the
gravitational compactness of the entire star and
$y=y_{*}(R_{\textrm{\mbox{out}}})$.

###  III.4 Stability Analysis

The stability of stellar configurations may be determined by analyzing the
spectrum of radial density oscillations. In a complete analysis, one adds a
harmonic perturbation of the form $\xi(r)e^{-i\omega t}$ to the metric fields
and linearizes the Einstein field equations with respect to the perturbation.
The result of this procedure is a Sturm-Liouville problem, whose solution
determines the frequency eigenvalues $\\{\omega_{i}\\}$. In this approach,
unstable radial modes manifest as exponentially increasing solutions, with
$\operatorname{Im}\omega_{i}>0$ Comer  _et al._ (1999); Leung  _et al._
(2012); Kain (2020). The equation of radial pulsation for a single-fluid
compact star was first derived by ChandrasekharChandrasekhar (1964). Here, we
resort to a simpler, less rigorous, stability criterion, following Henriques
_et al._ (1990); Nyhan and Kain (2022). An analogous criterion, for single-
fluid stars, is explained in detail by Weinberg Weinberg (1972).

Assuming invariance under time reversal, each eigenmode $\xi_{i}(r)$ has a
pair of fundamental frequencies $\omega_{i}^{\pm}=\pm\sqrt{\omega_{i}^{2}}$,
where $\omega_{i}\in\mathds{R}$ is an eigenvalue. The onset of instability
then corresponds to the point at which the lowest-lying eigenvalue
$\omega_{0}^{2}$ flips sign and becomes negative. Precisely at this point,
$\omega_{0}^{\pm}=0$, and the field perturbations are simply $\xi_{0}(r)$,
which is a static configuration, and therefore, in equilibrium. Because
equilibrium configurations are uniquely specified by central densities, these
modes must correspond to shifts
$\epsilon_{c}^{i}\to\epsilon_{c}^{i}+\delta\epsilon_{c}^{i}$. Moreover, if we
approach $\omega_{0}^{2}\to 0^{+}$ from above, these configurations will be
joined by arbitrarily slow oscillations, which must leave the total (mirror)
baryon number $N_{i}$ unchanged.

All of this implies that at the onset of unstable radial modes, the total
mirror baryon and SM baryon numbers must be _stationary_ under variations of
$\epsilon_{c}^{i}$:

$\left(\begin{array}[]{c}\delta N_{\textrm{\mbox{\tiny SM}}}\\\ \delta
N_{\textrm{\mbox{\tiny
DM}}}\end{array}\right)=\left(\begin{array}[]{cc}{\partial
N_{\textrm{\mbox{\tiny SM}}}}/{\partial\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}}&{\partial N_{\textrm{\mbox{\tiny
SM}}}}/{\partial\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}}\\\ {\partial
N_{\textrm{\mbox{\tiny DM}}}}/{\partial\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}}&{\partial N_{\textrm{\mbox{\tiny
DM}}}}/{\partial\epsilon_{c}^{\textrm{\mbox{\tiny
DM}}}}\end{array}\right)\left(\begin{array}[]{c}\delta\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}\\\ \delta\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}\end{array}\right)=0.$
(11)

The existence of nontrivial solutions $\delta\epsilon_{c}^{i}\neq 0$ demands
that the matrix in Eq. (11) have zero determinant:

$\frac{\partial N_{\textrm{\mbox{\tiny
SM}}}}{\partial\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}}\frac{\partial
N_{\textrm{\mbox{\tiny DM}}}}{\partial\epsilon_{c}^{\textrm{\mbox{\tiny
DM}}}}-\frac{\partial N_{\textrm{\mbox{\tiny
SM}}}}{\partial\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}}\frac{\partial
N_{\textrm{\mbox{\tiny DM}}}}{\partial\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}}=0.$ (12)

This is the criteria for the onset of radial instability for two-fluid stars.

One might wonder if Eqs. (11) and (12) are consistent with the well-known
condition $\partial M/\partial\epsilon_{c}=0$ in the case of a single fluid.
For single-fluid stars, Eq. (11) will only have nontrivial solutions if
$\partial N/\partial\epsilon_{c}=0$. Under the assumption of uniform entropy
per baryon, the stellar mass is stationary under any transformation
$\epsilon(r)\to\epsilon(r)+\delta\epsilon(r)$ that leaves the total baryon
number $N$ unchanged Weinberg (1972). This holds if and only if the TOV
equation is satisfied, and thus, exclusively for equilibrium configurations.
In light of this result, demanding that $\partial N/\partial\epsilon_{c}=0$ is
equivalent to imposing that $\partial M/\partial\epsilon_{c}=0$, which is the
widely used criteria for the onset of instability Weinberg (1972).

The result mentioned above has recently been generalized to the case of two
fluid stars Goldman  _et al._ (2013). Given that the entropy per baryon in the
SM fluid and the entropy per mirror baryon in the DM fluid are kept uniform,
Eqs. (4a)-(4c) are satisfied if and only if any transformation
$\epsilon_{i}(r)\to\epsilon_{i}(r)+\delta\epsilon_{i}(r)$ satisfying $\delta
N_{i}=0$, with $i=\\{\textrm{SM},\textrm{DM}\\}$, also leaves the mass
unchanged, i.e. $\delta M=0$. Equation (11) is precisely the condition that
$N_{\textrm{\mbox{\tiny SM}}}$ and $N_{\textrm{\mbox{\tiny DM}}}$ are
stationary, under which case $M$ is stationary as well. That is, at the onset
of instability, $\delta N_{\textrm{\mbox{\tiny SM}}}=0$, $\delta
N_{\textrm{\mbox{\tiny DM}}}=0$ and $\delta M=0$ under variations of the
central energy densities, $(\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}},\epsilon_{c}^{\textrm{\mbox{\tiny
DM}}})\to(\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}}+\delta\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}},\epsilon_{c}^{\textrm{\mbox{\tiny
DM}}}+\delta\epsilon_{c}^{\textrm{\mbox{\tiny DM}}})$.

The physical intuition behind this simpler stability criterion is also left
unchanged. Suppose a single-fluid star is taken out of equilibrium, with
$\epsilon_{c}\to\epsilon_{c}+\delta\epsilon_{c}$, while its total baryon
number $N=N(\epsilon_{c})$ is kept fixed. The star is then at a baryon-number
difference $\Delta N=N(\epsilon_{c})-N(\epsilon_{c}+\delta\epsilon_{c})$ away
from equilibrium. Assuming $\partial M/\partial N>0$, if $\Delta N>0$, there
is a mass excess, indicating that the gravitational pull is too strong and the
star will contract. If $\Delta N<0$, gravitational attraction is not
sufficiently strong and the star expands. In the spirit of Le Chatelier’s
principle, equilibrium is stable if an increase in density
$\delta\epsilon_{c}>0$ leads to a restoring expansion, with $\Delta
N\approx-(\partial N/\partial\epsilon_{c})_{\textrm{eq}}\delta\epsilon_{c}<0$,
where the derivative is taken in hydrostatic equilibrium. Therefore, a
necessary condition for stable equilibrium is that $\delta
N/\delta\epsilon_{c}>0$, or equivalently $\delta M/\delta\epsilon_{c}>0$
Glendenning (2000).

For a two-fluid star, the two central densities and two conserved charges are
mutually intertwined, making the physical picture less clear. However, an
analogous heuristic argument can be found by diagonalizing the matrix
$\partial N_{i}/\partial\epsilon_{c\,j}$, with
$\\{i,j\\}=\\{\textrm{DM},\textrm{SM}\\}$, in Eq. (11). By doing so, one
obtains two independent sets of variables, $(\epsilon_{c}^{A},N_{A})$ and
$(\epsilon_{c}^{B},N_{B})$ corresponding to eigenvalues $\kappa_{A}$ and
$\kappa_{B}$, such that:

$\left(\begin{array}[]{c}\delta N_{A}\\\ \delta
N_{B}\end{array}\right)=\left(\begin{array}[]{cc}\kappa_{A}&0\\\
0&\kappa_{B}\end{array}\right)\left(\begin{array}[]{c}\delta\epsilon_{c}^{A}\\\
\delta\epsilon_{c}^{B}\end{array}\right).$ (13)

Because $N_{A}$ and $N_{B}$ are linear combinations of $N_{\textrm{\mbox{\tiny
SM}}}$ and $N_{\textrm{\mbox{\tiny DM}}}$, they are also conserved and are
kept fixed as the star is perturbed. Changes to $\epsilon_{c}^{A}$ and
$\epsilon_{c}^{B}$ are then performed independently to find that equilibrium
configurations can be stable only if both eigenvalues are positive

$\kappa_{A}>0,\qquad\kappa_{B}>0.$ (14)

This generalizes the widely used stability condition $\partial
M/\partial\epsilon_{c}>0$ to multi-fluid stars.

Finally, we observe that Eq. (12) provides one condition for two independent
variables $\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ and
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$. Therefore, as the mass-radius
relation for two-fluid stars is an area, the boundary of the stable region is
a set of curves. This is in contrast to the case of a single fluid, for which
the mass-radius relation is a curve and stable regions are delimited by
points.

Henceforth, we will determine stability by using the condition in Eq. (14),
which is in agreement with previous stability studies of two-fluid dark-matter
admixed stars Comer  _et al._ (1999); Leung  _et al._ (2012); Kain (2020). A
rigorous stability analysis of two-fluid stars, using e.g. our EoS, is left to
future work.

##  IV Properties of Individual Stars

Before we study coalescences involving NSs, MNSs and MANSs, we must first
address the properties of each of these stars in isolation. Masses and radii
follow directly from the EoSs discussed in Sec. II and the equations of
structure (including the two-fluid Tolman-Oppenheimer-Volkoff (TOV) equation),
which were reviewed in Sec. III.2. The calculation of tidal deformabilities,
can be computed as described in Sec. III.3. Results for isolated MNSs and NSs
were already presented in Hippert  _et al._ (2022), but we present them again
here for completeness. We also present results for a new type of object
(MANS), which we remind the reader are compact ojbects composed of an
admixture of mirror matter ($f/v>1$) and SM matter ($f/v=1$).

Figure 2: Upper panel: Mass-radius relations for NSs and MNSs with quark
masses scaled by $m_{q^{\prime}}/m_{q}=f/v$. Lower panel: tidal deformability-
mass for NSs and MNSs with quark masses scaled by $m_{q^{\prime}}/m_{q}=f/v$.
The blue and beige shaded regions correspond to $2\sigma$ confidence regions,
using X-ray observations of pulsars J0030+0451 and J0740+6620, respectively by
NICER Miller  _et al._ (2019); Riley  _et al._ (2019); Miller  _et al._
(2021); Riley  _et al._ (2021). Red shaded regions represent $2\sigma$
confidence intervals from GW event GW170817 observed by advanced LIGO/Virgo
Abbott  _et al._ (2017, 2018a, 2019).

###  IV.1 Neutron Stars and Mirror Neutron Stars

The NS EoS discussed in Sec. II Hippert  _et al._ (2022) can be inserted into
Eqs. (4a), (4b) and (7), with $p_{\textrm{\mbox{\tiny
DM}}}=\epsilon_{\textrm{\mbox{\tiny DM}}}=0$, to calculate the mass, radius
and tidal deformabilities of SM NSs. The resulting mass-radius and tidal
deformability-mass relation are shown in Fig. 2 and are consistent with the
current observational constraints Hippert  _et al._ (2022). The latest NICER
measurement of pulsar J0740+6620 Riley  _et al._ (2019); Miller  _et al._
(2019) was not used to constrain model parameters in Ref. Hippert  _et al._
(2022), but yet, the mass radius curve is consistent with these observations.
The aim of the present study is not to advocate for a particular model of NS
matter. Instead, we simply choose an EoS that is in agreement with known
astrophysical observations, and which we can extend to the mirror sector to
draw new predictions. An extensive discussion of uncertainties from the SM EoS
and mass-radius relation can be found in Ref. Hippert  _et al._ (2022), where
we show that our predictions for MNS matter are robust against these
uncertainties, as well as model details.

---
Figure 3:  Mass-radius plane of MANS (left panel) for different choices of
central SM and mirror-matter energy densities (right panel). The SM NS and the
MNS mass-radius sequences are shown with a cyan and a green dashed line
respectively. The cyan up triangle and the green down triangle correspond to
energy densities at which a SM NS and a MNS would have masses of $64\%$ their
Chandrasekhar limits respectively. Increasing the mirror matter and the SM
central densities respectively then leads to the red and orange boundaries.
Observe that the mapping from
$(\epsilon_{c}^{\mathrm{SM}},\epsilon_{c}^{\mathrm{DM}})$ to
$(R_{\mathrm{out}},M)$ is not bijective at low radii, allowing for the
possibility of ultimate twins. Two examples of ultimate twins are shown with
filled blue and open green circles and squares.

Similarly, the MNS EoS of Sec. II can be inserted into the same equations, but
this time with $p_{\textrm{\mbox{\tiny SM}}}=\epsilon_{\textrm{\mbox{\tiny
SM}}}=0$ and specific choices of $f/v$, to calculate the mass, radius and
tidal deformabilities of MNSs. The resulting mass-radius and $\Lambda$-mass
sequences of MNSs are also shown in Fig. 2, where Chandrasekhar-mass
configurations are marked with a filled circle in the mass-radius panel
Hippert  _et al._ (2022), and we calculate the sequences down to $64\%M_{Ch}$.
Figure 2 allows us to compare predictions for NSs and MNSs. Observe that, in
the phenomenologically relevant range $f/v\simeq 3-7$, the mass-radius and
$\Lambda$-mass sequences for MNS are shifted to lower masses and smaller
radii, and to smaller $\Lambda$ and smaller masses. This is perhaps most
clearly seen for the Chandrasekhar-mass configurations, which go from the
well-known $\sim 1.4M_{\odot}$ value for SM matter to sub-solar masses for
MNS. Moreover, the $\Lambda$-mass sequences for NSs and MNSs are very clearly
separated, indicating they could be distinguished by their GW signatures
Hippert  _et al._ (2022).

###  IV.2 Mirror-Matter-Admixed Neutron Stars

Besides isolated NSs and MNSs, our model enables us to make predictions for a
third class of objects, not addressed in Hippert  _et al._ (2022): $\rm MANS$,
compact objects akin to NSs and MNSs, but containing both SM ($f/v=1$) and
mirror matter ($f/v>1$) in their interior. We will here consider the structure
of MANS with an arbitrary fraction of mirror matter, which our two-fluid set
up is tailored-made to handle. In particular, we set $p_{\textrm{\mbox{\tiny
SM}}}(\epsilon_{\textrm{\mbox{\tiny SM}}})$ to the SM EoS and
$p_{\textrm{\mbox{\tiny DM}}}(\epsilon_{\textrm{\mbox{\tiny DM}}})$ to the
mirror matter EoS with a choice of $f/v$ in the 2-fluid structure equations
[Eqs. (4a) and (4b)] to solve for MANS. We vary the fraction of mirror matter
to SM through the choice of SM and mirror matter central energy densities. By
then varying $f/v$ and these central densities, we are able to study how the
astrophysical properties of $\rm MANS$ are affected by the Higgs VEV in the
mirror sector.

Before proceeding, let us comment on how the analysis we carry out here
differs from previous related work. Previous studies explored admixed stars
with mirror matter, but they set $f/v=1$ i.e. they assumed that the SM EoS was
identical to the mirror matter EoS Sandin and Ciarcelluti (2009); Ciarcelluti
and Sandin (2011); Kain (2021). Our study is quite different because we use
first-principle lattice QCD results to scale the SM EoS with $f/v$, thus using
a consistent EoS for both sectors. Most previous work fixed the fraction of
dark matter to SM because the specific dark matter capture model studied
predicted a certain fraction. The details of capture, however, depend
sensitively on the precise nature of the DM self-interactions, the
interactions between DM and SM matter, and unknown astrophysical details, such
as the precise distribution of mirror matter in our galaxy. Lifting this
restriction, the phase space of solutions of admixed stars is clearly two-
dimensional, and the properties of MANS vary depending on the choice of
fraction, as we will see below.

####  IV.2.1 Full mass-radius region

The full mass-radius plane of $\rm MANS$ is shown in the left panel of Fig. 3
for $f/v=5$ and various choices of central densities shown in the right panel
of this figure. Pure SM NSs are shown with a dashed cyan line, while pure MNSs
are shown with a green dashed line. Everything in between is a MANS, where we
now clearly see that what used to be a one-dimensional mass-radius sequence
for a NS or a MNS is now a two-dimensional mass-radius plane. In this figure,
we fixed $f/v=5$ as a representative example, but we have checked that we find
qualitatively similar results for other values of $f/v$.

How low of a mass can MANS have? SM NSs have a Chandrasekhar limit of $\sim
1.4M_{\odot}$, but lower NS mass measurements have been claimed using X-ray
observations Strobel  _et al._ (1999); Doroshenko  _et al._ (2022).666These
time-integrated measures have been found to be susceptible to systematic error
(especially due to assumptions about the atmospheres) such that radii can
differ by up to $\sim 50$% Catuneanu  _et al._ (2013); Servillat  _et al._
(2012), but we still adopt this measurement as an interesting benchmark. One
of the lowest SM NS observed is approximately $0.9M_{\odot}$ (i.e. $\sim 64\%$
of the Chandrasekhar limit), and this is shown with a cyan up triangle in Fig.
3. Fixing the SM central density to this value and then increasing the central
density of mirror matter leads to the solid red line. For illustration, let us
now assume that the smallest MNS possible is also $64\%$ of the Chandrasekhar
mass of that MNS sequence, shown with a green down triangle in Fig. 3. Fixing
the mirror matter central density, $\epsilon^{\textrm{\mbox{\tiny DM}}}_{c}$,
to $64\%$ of its Chandrasekhar mass and increasing the SM central density,
$\epsilon^{\textrm{\mbox{\tiny SM}}}_{c}$, leads to the orange line in the
figure. The region in the $(\epsilon_{c}^{\textrm{\mbox{\tiny
SM}}},\epsilon_{c}^{\textrm{\mbox{\tiny DM}}})$ below the red line and the
orange line is shown with a light gray color in both panels.

One of the most interesting features of the two-dimensional mass-radius plane
is the possibility of ultimate twins: stars with the same mass and radius but
different dark matter fraction. This occurs because the mapping from
$(\epsilon_{c}^{\mathrm{SM}},\epsilon_{c}^{\mathrm{DM}})$ to
$(R_{\mathrm{out}},M)$ is not bijective at low radii. Examples of such twins
are shown with a filled blue circle and square, and an open green circle and
square. As is clear from the right panel, these MANS have different ratios of
mirror matter to SM matter central densities; the green cases are actually
MNSs, while the blue cases are MANSs. Although not shown in the figure, there
are actually an infinite number of ultimate twins, as one varies the central
densities, with the pattern shown in the figure.

Figure 4:  Typical outer radius $R_{\textrm{\mbox{out}}}$ for stellar
configurations with different mirror baryon fractions $Y_{\textrm{\mbox{\tiny
DM}}}\equiv N_{\textrm{\mbox{\tiny DM}}}/(N_{\textrm{\mbox{\tiny
SM}}}+N_{\textrm{\mbox{\tiny DM}}})$. The inset illustrates the definition of
the inner and outer radii, $R_{\textrm{\mbox{in}}}$ and
$R_{\textrm{\mbox{out}}}$. The minimum value for the outer radius corresponds
to a configuration where both fluids occupy the entirety of the star, that is,
when $R_{\textrm{\mbox{in}}}=R_{\textrm{\mbox{out}}}$, which happens at
$Y_{\textrm{\mbox{\tiny DM}}}\sim 0.75$. For $0<Y_{\textrm{\mbox{\tiny
DM}}}\lesssim 0.75$, we find NSs develop a DM core, $Y_{\textrm{\mbox{\tiny
DM}}}\gtrsim 0.75$, we find a DM halo instead.

Let us pause at this point and discuss an unexpected feature of the results
presented above: the allowed region of $\rm MANSs$ is not bounded from the
left by the $\rm MNS$ sequence (thick long-dashed green curves in Fig. 3).
Consider then what happens to MANS as we start with a SM NS and we begin to
add mirror matter while keeping the total mass fixed, as shown schematically
in the cartoon of Fig. 4. As one increases the DM fraction
$Y_{\textrm{\mbox{\tiny DM}}}\equiv N_{\textrm{\mbox{\tiny
DM}}}/(N_{\textrm{\mbox{\tiny SM}}}+N_{\textrm{\mbox{\tiny DM}}})$ by
increasing the mirror matter energy density, one can identify the following
stages:

  * •

SM NSs. Initially, $Y_{\textrm{\mbox{\tiny DM}}}=0$ and the star is purely
made out of SM matter

  * •

MANS with a DM core and a SM halo. As we increase
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ and thus $Y_{\textrm{\mbox{\tiny
DM}}}$, the mirror matter settles at the center of the star (because of the
extra concentration of gravitational mass in the admixed region) and the star
now has two radii: one which contains some of the SM matter and all of the
mirror matter ($R_{\rm in}=R_{\textrm{\mbox{\tiny DM}}}$) and one that
contains all of the SM and mirror matter in the star ($R_{\rm
out}=R_{\textrm{\mbox{\tiny SM}}}$). The spherical region inside the inner
radius, $r<R_{\rm in}$, will be called the core, and for these MANSs, it is
where all of the mirror matter resides. The shell region outside the inner
radius but inside the outer radius, $R_{\rm out}>r>R_{\rm in}$, will be called
the halo, and for these MANS, it is composed entirely of SM matter.

  * •

MANS with equal DM and SM matter. As we increase
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ and $Y_{\textrm{\mbox{\tiny DM}}}$
further, the radius that contains the mirror matter grows, while the radius
that contains the SM matter shrinks. Eventually, one reaches a critical value
$Y_{\textrm{\mbox{\tiny DM}}}=Y_{\textrm{\mbox{\tiny DM}}}^{*}$ at which
$R_{\rm in}=R_{\textrm{\mbox{\tiny DM}}}=R_{\textrm{\mbox{\tiny SM}}}=R_{\rm
out}$ and the entire star is occupied by both fluids. This is when the MANS
sequence reaches its minimum radius.

  * •

MANS with a SM core and a DM halo. Increasing
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ and $Y_{\textrm{\mbox{\tiny DM}}}$
even further, the radius that contains SM matter is now inside the radius that
contains the DM, so $R_{\textrm{\mbox{\tiny SM}}}<R_{\textrm{\mbox{\tiny
DM}}}$, and therefore, the MANS now has a SM core ($R_{\rm
in}=R_{\textrm{\mbox{\tiny SM}}}$) and a DM halo ($R_{\rm
out}=R_{\textrm{\mbox{\tiny DM}}}$). This is the flip case to the MANS with a
DM core and a SM halo.

  * •

MNSs. Eventually, $\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ is increased
enough that $Y_{\textrm{\mbox{\tiny DM}}}$ tends to unity and the MANS becomes
a MNS. Since the radius that contains the DM grows as we increase
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$, $R_{\rm out}=R_{\textrm{\mbox{\tiny
DM}}}$ while $R_{\rm in}=R_{\textrm{\mbox{\tiny SM}}}=0$. In particular, note
that the radius of the MNS is here larger than the radius of the MANS with
equal DM and SM matter.

From this analysis, one can clearly see that the stellar radius
$R_{\textrm{\mbox{out}}}$ decreases for $Y_{\textrm{\mbox{\tiny
DM}}}<Y_{\textrm{\mbox{\tiny DM}}}^{*}$ and it increases for
$Y_{\textrm{\mbox{\tiny DM}}}>Y_{\textrm{\mbox{\tiny DM}}}^{*}$, until at
$Y_{\textrm{\mbox{\tiny DM}}}=Y_{\textrm{\mbox{\tiny DM}}}^{*}$ the MANS has
the minimum allowed value for $R_{\textrm{\mbox{out}}}$. This critical value
is the point at which the star is entirely filled with both fluids, so that
the maximum possible concentration of mass is achieved. The sequence of stars
with $Y_{\textrm{\mbox{\tiny DM}}}=Y_{\textrm{\mbox{\tiny DM}}}^{*}(M)$
therefore delimits the mass-radius region from the left in Fig. 3.

A consequence of the results presented above is that the mirror-baryon
fraction $Y_{\textrm{\mbox{\tiny DM}}}$ is indeed zero at the right most
boundary (at the SM NS sequence), but it is not unity at the left most
boundary, and therefore, $Y_{\textrm{\mbox{\tiny DM}}}$ must change non-
monotonically as the MANS radius changes near the MNS sequence. Because the
radius is non-monotonic in $Y_{\textrm{\mbox{\tiny DM}}}$, the region to the
left of the $\rm MNS$ sequence (to the left of the thick long-dashed green
curve in Fig. 3) is composed of two overlapping areas with different values of
$Y_{\textrm{\mbox{\tiny DM}}}$. Let us investigate how these overlapping areas
are connected . Since $Y_{\textrm{\mbox{\tiny DM}}}$ is one-to-one in
$\epsilon_{c}^{DM}$, one can focus on any constant
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ curve, such as the orange curve in
Fig. 3. Starting at the intersection of this sequence with the MNS sequence
(i.e. the intersection of the thick dashed green curve and the orange curve,
which also corresponds to the right-most star in Fig. 4)), as
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ increases, the mass-radius point
moves left along the orange curve, as shown in the top panel of Fig. 5, which
is just a zoomed version of the left panel of Fig. 3. After hitting the left
boundary, where $R_{\textrm{\mbox{out}}}=R_{\textrm{\mbox{in}}}$ (the central
star in Fig. 4), the mass-radius sequence turns around and begins to move
right. When this happens, the region left to the pure $\rm MNS$ sequence
becomes double-valued in $Y_{\textrm{\mbox{\tiny DM}}}$. One can see this more
clearly in the bottom panel of Fig. 5, which shows the mass-radius plane
again, with points of different $Y_{\textrm{\mbox{\tiny DM}}}$ color-coded.
This region between the MNS sequence and the left-most boundary is where the
ultimate mass-radius twins live.

---
Figure 5: Mass-radius regions, zoomed to the region inside which ultimate
twins exist (top) and color-coded to indicate stars with different DM fraction
(bottom), fixing $f/v=5$. The arrows in the top panel indicate the direction
in which the SM central energy density increases. Observe that that there is a
region (between the left-most boundary and the MNS sequence) inside which
different $Y_{\textrm{\mbox{\tiny DM}}}$ leads to stars with the same mass and
radius, which we have defined in this paper as ultimate twins.

####  IV.2.2 MANSs through Accretion of Dark Matter into Neutron Stars

---
Figure 6: Left panel: Mass-radius region for MANSs that contain a DM core
(instead of a DM halo). The black dashed line is the isolated NS mass-radius
sequence. Green shaded regions are constraints from LIGO and NICER Miller  _et
al._ (2019); Riley  _et al._ (2019); Miller  _et al._ (2021); Riley  _et al._
(2021); Abbott  _et al._ (2017, 2018a, 2019). The overlapping shaded regions
connected to the isolated NS sequence are the mass-radius regions for MANSs
with a DM halo for different values of $f/v$. Right panel: DM core mass as a
function of its radius $R_{in}$ for different values of $f/v$. Observe that
the mass-radius planes are mostly insensitive to the value of $f/v$, although
the relation between the mass and radius of the DM core is not.

Having explored the entire mass-radius plane that connects NSs and MNSs, we
now explore predictions for MANSs formed from accretion of dark matter into SM
NS. Since SM NSs cannot have masses smaller than some percentage (64% for this
paper) of their respective Chandrasekhar masses ($0.9M_{\odot}$ for this
paper), the mass of the MANS after dark matter accretion is also constrained
from below. Such MANS have a DM core instead of a DM halo, with $R_{\rm
in}=R_{\textrm{\mbox{\tiny DM}}}<R_{\textrm{\mbox{\tiny SM}}}=R_{\rm out}$.

The left panel of Fig. 6 shows the mass-radius plane for these objects for
various values of $f/v$. Observe that the mass-radius region for these MANS
depends very weakly on $f/v$, and thus, on the mirror quark mass
$m_{q}^{\prime}=(f/v)\,m_{q}$. The right panel of Fig. 6 shows the mirror-
matter mass $M_{\textrm{\mbox{\tiny DM}}}$ as a function of the inner radius
$R_{\textrm{\mbox{in}}}$. While the inner radius $R_{\textrm{\mbox{in}}}$ can
change by a factor of $\sim 2$ with varying $f/v$, the maximum value of
$M_{\textrm{\mbox{\tiny DM}}}$ changes at most by $\sim 15\%$.

Why is the mass-radius plane approximately universal with $f/v$ but the
$M_{\textrm{\mbox{\tiny DM}}}$-$R_{\textrm{\mbox{in}}}$ plane is not? Let us
consider the mass-radius plane first. Near the SM NS sequence, for small
$M_{\textrm{\mbox{\tiny DM}}}$, the DM core is very small (compare to the
total radius of the star). The influence of such a core is then to only
increase the total mass (by adding DM mass), while keeping the size of the
inner radius small (as one can see on the right panel). As the
$M_{\textrm{\mbox{\tiny DM}}}$ increases further, the size of the inner core
also increases and the $f/v$ universality begins to be lost because different
$f/v$ lead to different sizes of DM cores (as seen again on the right panel).

A more detailed view of the maximum DM mass, $M_{\textrm{\mbox{\tiny
DM}}}^{\textrm{max}}$, as a function of the initial SM mass
$M_{\textrm{\mbox{\tiny SM}}}$, is presented in Fig. 7. Indeed, one observes
that the curves for different $f/v$ lie very close to each other for small
$M_{\textrm{\mbox{\tiny DM}}}^{\textrm{max}}$, partly explaining the behavior
observed in Fig. 6. The dependence of $M_{\textrm{\mbox{\tiny
DM}}}^{\textrm{max}}$ on $f/v$ becomes weaker for decreasing
$M_{\textrm{\mbox{\tiny DM}}}$ and for increasing $f/v$, as
$R_{\textrm{\mbox{in}}}$ becomes larger.

Figure 7: Maximum amount of mirror matter that can be added to a SM NS before
it becomes unstable, as a function of its total SM mass. Different values of
f/v are considered for the mirror matter contribution.

One might naïvely expect that adding DM to a NS would correspond to increasing
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ at a fixed
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$, but this, however, is not the case.
As DM is added to a NS, the increase in stellar mass leads to an increase in
gravitational pull, which leads to an increase in the SM central density
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$. In fact, this change in
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ can be used to investigate the
effect of DM-admixture on the matter distribution and hydrostatic balance of
the SM fluid. Results for the fractional change of
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$ as a function of
$M_{\textrm{\mbox{\tiny DM}}}$ are shown in Fig. 8, where we fix SM mass at
$M_{\textrm{\mbox{\tiny SM}}}=0.9M_{\odot}$. The point at which the star
becomes radially unstable is marked by a star and coincides with the maximum
of $M_{\textrm{\mbox{\tiny DM}}}$, as one would expect. Observe that the SM
central density can increase by a factor of $\sim 5-7$ due to the admixture of
DM. Surprisingly, this change tends to be smaller for smaller values of $f/v$,
even though $M_{\textrm{\mbox{\tiny DM}}}^{\textrm{max}}$ tends to be larger.

Figure 8: Fractional change in the central energy density of SM matter
$\epsilon_{c}^{\textrm{\mbox{\tiny SM}}}$, normalized to its value in the
absence of dark matter $\epsilon_{c}^{\textrm{\tiny NS}}$, as a function of
the admixed DM mass $M_{\textrm{\mbox{\tiny DM}}}$. Here, we increase
$M_{\textrm{\mbox{\tiny DM}}}$ at a fixed SM mass of
$M_{\textrm{SM}}=0.9\,M_{\odot}$. The point at which the star becomes radially
unstable is marked by a star, with faint lines corresponding to unstable
configurations. Different colors correspond to different values of $f/v$.

###  IV.3 Compactness and Tidal deformability of MANSs

Figure 9:  Compactness of MANSs as a function of total mass (top) and tidal
deformability as a function of total mass (bottom) for different values of
$f/v$ (shown with different colors). The black solid line corresponds to SM
NSs, while the colored solid lines correspond to MNSs with different values of
$f/v$. Observe that the compactness of MANSs is at most $0.3$, which is
comparable to that of SM NSs. Such high-compactness stars also correspond to
the maximum mass configurations, which occurs at lower masses for MANSs than
for SM NSs. Similarly, the tidal deformability of MANSs can be extremely low,
even at low masses.

As we have seen above, MANSs can have rather small radii, so one may expect
their compactness $\mathcal{C}=M/R_{out}$ to also be large and their tidal
deformability to be small. The top panel of Fig. 9 shows the compactness as a
function of the total mass $M$ for MANSs with different values of $f/v$.
Observe that, even with the addition of DM, MANS never reach a compactness
above $C=0.3$, which is comparable to the maximum compactness of SM NSs.
However, for MNSs and DM-rich MANS, this large value of compactness is
achieved at much lighter masses. As the mass of a star increases, the star
acquires a higher compactness, and therefore, it should become harder to
deform and thus it possess a smaller tidal deformability (i.e. a smaller
$\Lambda$). The relation between $\Lambda$ and the total mass $M$ is shown on
the middle panel of Fig. 9, which corroborates this expectation.

---
Figure 10: Compactness versus total mass (left) and tidal deformability versus
total mass (right) for the $f/v=5$ case. The color coding of curves is the
same as in Fig. 3, except that here we also introduce a new sequence (colored
magenta), with arrows indicating the direction in which
$\epsilon^{\textrm{\mbox{\tiny DM}}}_{c}$ increases, which we discuss in the
main text. Observe that the behavior of the compactness and the tidal
deformability is not monotonic in the total mass. Observe also that ultimate
twins have the same mass, radius and compactness, but different tidal
deformabilities.

The tidal deformability of MANSs, however, is not monotonic in the mass of the
star, contrary to the intuition presented above. Figure 10 shows the
compactness as a function of mass (left) and the tidal deformability as a
function of mass (right) for the $f/v=5$ case. For concreteness, let us focus
on a particular stellar sequence, depicted through the dotted magenta line.
When the mirror-matter central density is zero
($\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}=0$), the MANS is simply a SM NS,
shown at the point where the magenta line connects to the cyan line. As one
increases $\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$, the sequence moves to the
left (as shown by the arrows on the magenta line), terminating in the ultimate
twins shown with the green empty circle and blue filled circle. Observe that
for $M\gtrsim 0.85M_{\odot}$, as $\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$
increases, the mass and compactness decrease, while the tidal deformability
increases. The compactness decreases because, although the radius decreases,
it does so slower than the mass. The tidal deformability increases because, as
the star becomes lighter, it also becomes easier to deform, as expected.

This behavior, however, changes drastically when
$\epsilon_{c}^{\textrm{\mbox{\tiny DM}}}$ has increased enough that $M\lesssim
0.85M_{\odot}$. Now, as the mass continues to decrease, the compactness
increases, while the tidal deformability decreases. The compactness increases
because the radius begins to decrease faster than the mass. In turn, the tidal
deformability decreases because, as the compactness increases, the star
becomes more difficult to deform. We conclude from this then that the tidal
deformability is indeed tied to how easy or hard it is to deform a star, but
the latter is connected to the compactness of the star, and not just its total
mass, when one changes $\epsilon_{\textrm{\mbox{\tiny DM}}}^{c}$ while fixing
$\epsilon_{\textrm{\mbox{\tiny SM}}}^{c}$.

Another interesting feature of Fig. 10 refers to the ultimate twins, which we
recall are depicted with a green empty circle and a green filled square. As
explained earlier, ultimate twins are MANSs with the same mass and radius, but
different internal DM composition and, thus, different
$\epsilon_{\textrm{\mbox{\tiny DM}}}^{c}$. As such, ultimate twins must also
have the same compactness, which is corroborated by the left panel of Fig. 10.
The tidal deformabilities of these twins, however, is not the same, as shown
in the right panel of this figure. Indeed, the twin with a combination of
mirror matter and SM matter (the blue filled square and filled circle in Figs.
10 and 3) has a lower tidal deformability than the twin that only possesses
mirror matter (the green empty circle and empty square in Figs. 10 and 3).

Why is this? The reason that for a two-fluid star the tidal deformability does
not just scale with the compactness can be gleaned from the definition of
$\Lambda$ in Eq. (7) and (9). As one can see from these equations, $\Lambda$
is a function of both the compactness of the star and the variable $y_{*}$,
evaluated at the surface of the star. From the differential equation that
$y_{*}$ must satisfy, one can schematically argue that $y_{*}(R_{\rm out})$
will depend on the pressure and energy density inside the star, the averaged
value of which can be related to the compactness of a star in a single-fluid
model. For a two-fluid model, however, the average pressure and density do not
just scale with the compactness. In fact, for a two-fluid model, the averaged
value of the pressure is much larger for a two-fluid star than for a single-
fluid star. A larger amount of interior pressure translates into a star that
is more difficult to deform, and thus one with a lower tidal deformability.
This is indeed what we find with ultimate twins, which always have one member
that is very near the MNS sequence (essentially a single fluid star with
$\epsilon^{\textrm{\mbox{\tiny SM}}}_{c}\approx 0$) and one member in the MANS
plane (a two fluid star with $\epsilon^{\textrm{\mbox{\tiny SM}}}_{c}\neq
0\neq\epsilon^{\textrm{\mbox{\tiny DM}}}_{c}$). Therefore, the MANS member of
the twin has a larger averaged pressure, is more difficult to deform, and thus
has a smaller $\Lambda$.

The different tidal deformabilities suggests the possibility that GW
observations of the inspiral of stellar-mass compact objects could be used to
detect MANSs. We have already seen that the $\Lambda$-$M$ relation of MANSs
lies on a plane instead of a line. Therefore, a set of measurements of
$\Lambda$, each with sufficient accuracy, could allow us to reconstruct the
$\Lambda$-$M$ relation and determine its dimensionality. To further determine
whether there are ultimate twins in this plane, one would have to measure the
mass, radius and tidal deformability of the same star, which may be achievable
if one can infer not just $\Lambda$, which is the quadrupolar (electric-type)
tidal deformability, but also the octopolar one. This could be achievable with
third-generation detectors.

Inferences on the compactness, and thus the radius, of MANS through the
$\Lambda$-C relation, however, must be done with great care. This is because,
although SM NSs share a nearly EoS insensitive relation between the tidal
deformability and the compactness, this universality is lost in MANS, as shown
in Fig. 11. Indeed, observe how SM NS sequences and MNS sequences all share a
nearly identical $\Lambda$-$C$ relation, MANSs do not, and in fact, they
possess a two-dimensional $\Lambda$-$C$ relation. Therefore, the measurement
of the tidal deformability of MANSs cannot be associated with a single
compactness (and thus a single radius).

---
Figure 11: Tidal deformability as a function of compactness for different
values of $f/v$. Observe that while SM NSs and MNSs (solid black and color
lines) all share an approximately insensitive $\Lambda$-$C$ relation, but MANS
do not. Instead, the $\Lambda$-$C$ relation for MANS is not a one-dimensional
curve, but rather a two-dimensional place, which can differ significantly from
the relation for SM NSs and MNSs.  Figure 12: Mass ratio $q$ vs. total mass
$M$ of binary systems for three different values of $f/v$. The shaded region
enclosed by solid lines corresponds to binaries composed of one NS and one
MNSs, while the region delimited by dashed lines corresponds to binaries of
two MNSs. For reference, we also show the region corresponding to SM NS
binaries, delimited by the black dotted lines.

##  V Inspiralling mixed binaries

In this section, we present predictions for the observable signatures of
inspiraling binaries of NSs and MNSs. We present results for the chirp mass
and mass ratio between these stars and find that MNS-MNS and NS-NS inspirals
produce unique GW signatures.

###  V.1 Mass Ratio Properties

Having determined the structure of isolated NSs, MNSs and MANSs, we now
investigate binary systems composed of NSs and MNS, specifically NS-$\rm MNS$,
and $\rm MNS$-$\rm MNS$ binaries.

The properties of these binaries are limited by the formation mechanisms for
these objects. For instance, there are no known mechanisms to create NSs below
a certain threshold $M^{(\textrm{\mbox{\tiny SM}})}_{\textrm{min}}$. To limit
ourselves to viable NS and MNS masses, we implement minimum mass thresholds
according to the estimated constraint $M\gtrsim
M^{(\textrm{SM})}_{\textrm{min}}\simeq 0.9-1\,M_{\odot}$ for NSs Strobel  _et
al._ (1999). For definiteness, and to avoid the risk of being overly
restrictive, we take $M^{(\textrm{SM})}_{\textrm{min}}=0.9\,M_{\odot}$. To
extrapolate this estimate to the mirror sector, we fix the ratio between
$M_{\textrm{min}}$ and the Chandrasekhar mass of MNSs at
$M^{(\textrm{SM})}_{\textrm{min}}/M_{\textrm{Ch.}}=0.64$. In practice, this
yields a scaling of $M_{\textrm{min}}$ with the (mirror) baryon mass
$m_{B}^{\prime}$:

$\frac{M^{(\textrm{\mbox{\tiny
DM}})}_{\textrm{min}}}{M^{(\textrm{SM})}_{\textrm{min}}}=\frac{M^{(\textrm{\mbox{\tiny
DM}})}_{\textrm{Ch.}}}{M^{(SM)}_{\textrm{Ch.}}}=\left(\frac{m_{B}}{m_{B}^{\prime}}\right)^{2}\,,$
(15)

where $m_{B}$ is the SM baryon mass and $m_{B}^{\prime}/m_{B}$ as a function
of $f/v$ is taken from Ref. Hippert  _et al._ (2022). These are the same
restrictions we placed on the minimum mass of MANS in Sec. IV.2.1.

With this in mind, let us now look at the allowed mass ratio $q\equiv
M_{1}/M_{2}$ and total mass $M_{bin}=M_{1}+M_{2}$ regions for NS-NS, NS-MNS
and MNS-MNS binaries, using the convention $M_{1}\leq M_{2}$. Results are
shown in Fig. 12, where the shaded areas correspond to the predicted region
for NS-MNS inspirals. The areas corresponding to NS-NS and MNS-MNS binaries
are enclosed by dotted and dashed lines, respectively. Each panel corresponds
to a different value of $f/v=3,5$ and $7$. For all values of $f/v$, we find
significant regions in the $q$-$M_{bin}$ space that are exclusive of NS-MNS
and MNS-MNS inspiral and mergers, indicating that these systems can be
effectively distinguished from the binary masses alone. For $f/v\gtrsim 5$, we
find no overlap between the three different regions. For example, one could
have a $(0.5,0.5)M_{\odot}$ binary MNS inspiral, which would lead to a mass
ratio of 1 and a total mass of $1M_{\odot}$, which is not possible for a
binary NS inspiral.

Figure 13: Symmetric mass ratio $\eta$ vs. chirp mass $\mathcal{M}$ of binary
systems for three different values of $f/v$. The shaded region enclosed by
solid lines corresponds to binaries composed of one NS and one MNSs, while the
region delimited by dashed lines corresponds to binaries of two MNSs. For
reference, we also show the region corresponding to SM NS binaries, delimited
by the black dotted lines.

###  V.2 Gravitational Wave Signatures

In the event of a NS-MNS binary coalescence, GW radiation from the inspiral
phase would provide the best candidate for signatures of its exotic nature. In
a coalescence, the evolution of the orbital and, therefore, GW, frequency is
determined by the chirp mass $\mathcal{M}\equiv\eta^{3/5}(M_{1}+M_{2})$ of the
binary, where $\eta\equiv M_{1}M_{2}/(M_{1}+M_{2})^{2}$ is the symmetric mass
ratio. The chirp mass and symmetric mass ratio regions relevant to NS-NS, NS-
MNS and MNS-MNS inspirals are shown in Fig. 13. The low values of chirp mass
for MNS-MNS and NS-MNS binaries indicate that the inspiral phase would be
responsible for most of the detectable GW radiation emitted by such binaries.
The merger of such compact objects would be outside the sensitivity band of
second-generation ground-based detectors.

The inspiral of binaries containing MNSs would also be less loud than those
composed of SM NSs. The chirp mass of the binary system determines the
magnitude of the GW signal, since the latter scales with
$\sim\mathcal{M}^{5/6}/D_{L}$, where $D_{L}$ is the luminosity distance to the
source. As can be seen from the figure, the chirp mass of a NS-MNS binary is
$\mathcal{M}\sim 0.9\,M_{\odot}$ for $f/v=3$, and $\mathcal{M}\sim
0.7\,M_{\odot}$ for $f/v=7$, which is smaller than $\mathcal{M}\sim
1.2M_{\odot}$ for NS-NS binaries. This suggests that GW signals from NS-MNS
inspirals should be $\sim 60-80\%$ that of NS-NS binaries, and therefore,
although weaker, they should still be detectable. Such a detection would then
allow us to distinguish a NS-MNS binary from a NS-NS binary just from the
chirp mass measurement.

GW emission in the late inspiral phase is also characterized by the tidal
deformabilities of the inspiraling stars. The tidal deformabilities of the two
stars, $\Lambda_{1}$ and $\Lambda_{2}$, can be calculated by solving Eq. (7)
with a given EoS Flanagan and Hinderer (2008). Using current GW detectors,
only a certain combination of the individual tidal deformabilities
$\Lambda_{1}$ and $\Lambda_{2}$, the so-called “chirp deformability”Favata
(2014); Yagi and Yunes (2014), can be measured. This implies that without any
additional information, the individual tidal deformabilities are degenerate
and cannot be solved for from a GW measurement. One way to break this
degeneracy is to use EoS insensitive relations, including the binary Love
relation Yagi and Yunes (2016), or some functional form for the EoS Abbott
_et al._ (2018b). However, both approaches assume that the merging compact
objects have a single EoS that can describe both stars. In the case of an
inspiraling NS-MNS, this is a poor assumption because each compact object has
a different EoS, therefore completely breaking the universality of the binary
Love relation.

Figure 14: Binary Love relations for binaries of two MNSs (top), a MNS with a
heavier NS companion (middle), and a MNS with a lighter NS companion (bottom).
The asymmetric binary Love number
$\Lambda_{a}\equiv(\Lambda_{1}-\Lambda_{2})/2$ is shown as a function of the
symmetric Love number $\Lambda_{s}\equiv(\Lambda_{1}+\Lambda_{2})/2$.
Different colors correspond to different values of $f/v$, while the distinct
line styles represent different mass ratios $q\equiv M_{1}/M_{2}$. Observe
that, in the top panel, curves for fixed values of $q$, but different values
of $f/v$, lie on top of each other, and therefore present EoS-insensitivity.
In the middle panel, two of the curves are very close to each other, but they
correspond to different mass ratios. In the bottom panel, the EoS-
insensitivity is completely lost.

To explore the potential breakage of the binary Love relations, we plot the
relations for MNS-MNS and NS-MNS binaries in Fig. 14. The top panel shows the
binary Love relations of MNS-MNS binary systems with both stars with the same
$f/v$. One can clearly see that the EoS insensitivity is preserved for
different values of $f/v$, keeping $q$ fixed, which is unsurprising for MNS-
MNS mergers because they both come from the same EOS. In that panel, NS-NS
binaries are shown in black, and they lie on top of the other curves. The
middle panel shows the binary Love relations of a lighter NS with a heavier
MNS. We cannot calculate the relations for lower mass ratios and high $f/v$
because, from Fig. 2, MNSs become lighter for higher $f/v$, and for some $q$
and $f/v$, there are no MNS that are heavy enough. In this case, we see that
the EoS insensitivity is lost777The $f/v=2$ and $q=0.9$ curve is almost on top
of the $f/v=3$ and $q=0.60$ curve, by pure coincidence., with the $f/v=2$
curve lying far from the $f/v=3$ curve for fixed $q=0.90$. The bottom panel
shows the binary Love relations for a higher NS with a lighter MNS. Because
MNSs tend to be lighter than NS, this is the most common scenario of a MNS-NS
binary. Observe, once more, that the EoS insensitivity is completely lost
because, for the same mass ratio $q$ (dashed curves for example), there is no
overlap for curves with different values of $f/v$. Observe also that, since
MNSs are typically less deformable than NSs (due to their higher compactness),
the asymmetric Love number $\Lambda_{a}$ can be negative. This occurs for both
$f/v>4$ and $q>0.60$, and tends toward $\Lambda_{a}=-\Lambda_{s}$ as $f/v$ and
$q$ increase, making $\Lambda_{1}$ much smaller than $\Lambda_{2}$.

Before continuing with a discussion of the post-merger phase, let us end with
a warning. Many GW searches today use a prior on the chirp mass and symmetric
mass ratio that may exclude MNS stars all together. Indeed, only a fraction of
the allowed $\eta$-$\mathcal{M}$ region overlaps with the region expected to
be occupied by NS-NS binaries. Because of this, GW events from binary systems
including one or two MNSs run the risk of being overlooked in a template-based
coherent analysis. To detect such events, the priors on the chirp mass and the
symmetric mass ratio should be updated to include the relevant mass ranges
(i.e. much lower masses). Alternatively, one could rely on wavelet methods or
excess power methods, both of which do not rely on templates, to detect such
events, although such approaches may be less efficient than template-based
searches.

##  VI Post-Merger Remnant

In this section, we discuss the possible remnants of a MNS-NS collision. In
particular, we present results for possible stable MANS remnants. We also
comment on the formation of black holes when a stable stellar configuration is
not achieved.

Figure 15: Left panel: Mass-radius region for stable remnants of a NS-MNS
merger. The black dashed line is the isolated NS mass-radius sequence. Green-
shaded regions are constraints from LIGO and NICER Miller  _et al._ (2019);
Riley  _et al._ (2019); Miller  _et al._ (2021); Riley  _et al._ (2021);
Abbott  _et al._ (2017, 2018a, 2019). The overlapping shaded regions are the
allowable mass-radius region for stable remnants, for different values of
$f/v$. Right panel: The shaded regions correspond to the allowed dark-matter
mass and radius for the mirror-matter core of the stable remnants of NS-MNS
mergers. Stable remnants are not expected for $f/v\lesssim 4$.  Figure 16:
Mirror matter fraction $Y_{\textrm{\mbox{\tiny DM}}}$ versus total mass
$M_{\textrm{\mbox{\tiny MANS}}}$ for stable remnants of NS-MNS coalescences.
Allowed regions for different $f/v$ are represented as shaded areas of
different color. Stable remnants are not expected for $f/v\lesssim 4$.

###  VI.1 NS-MNS Collision

A NS-MNS binary merger would be a unique event. Because there are no
interactions via QCD, QED, or weak interactions that exist between mirror
matter and the SM, one cannot directly change the SM EoS of NSs to add in the
contribution of mirror matter. In other words, collisions of NSs and MNSs
would not be like NS-NS collisions at all because they would be like “ghosts
passing through each other.” Rather, they are likely to continue to inspiral
well passed the regime where their surfaces are overlapping. Once inside each
other, the nested stars are likely to keep on spinning at different rotational
angular frequencies, even as differential rotation within each fluid subsides,
because the only mechanism for angular-momentum transfer between the two
fluids in the remnant is gravitational, including tidal effects. While the
observational signature of NS-MNS mergers can be interesting due to the above
differences to NS-NS mergers, conclusions would have to rely on input from
computer simulations, which are not yet available. We thus leave the
discussion of the merger phase of NS-MNS binaries for future work.

###  VI.2 MANS Remnants

After a sufficiently long time, transient effects should die out and a NS-MNS
merger would lead to a black hole or to a MANS, like the ones discussed in
Sec. IV.2. In both cases, the resulting configurations are restricted by the
relevant mass thresholds for each of the colliding stars. Let us then consider
these end-states and discuss whether they would lead to interesting
observational signatures.

Let us once again consider the viable star requirement ($M_{\rm NS}\geq
0.9\,M_{\odot}$) but this time we will also apply this constraint to MNS
($M_{\rm MNS}\geq 0.9(m_{B}/m_{B}^{\prime})^{2}\,M_{\odot}$). To accomplish
this, we follow the procedure in Sec. V.1 and use Eq. (15). For simplicity, we
neglect the effects of rotation on the stellar structure. We also neglect
possible (SM and DM) particle number losses during the merger. The resulting
range for MANS remnants is shown in the left panel of Fig. 15. The right panel
of this figure shows the DM mass $M_{\textrm{\mbox{\tiny DM}}}$ versus inner
radius $R_{\textrm{\mbox{in}}}$ for values of $f/v$ of the stable remnants.
The effect of the lower-mass threshold $M_{\textrm{min}}$ on possible stable
post-merger remnants is noticeable. In combination with radial stability
constraints (see Sec. III.4), the minimal mass thresholds lead to a tight
region in $M_{\textrm{\mbox{\tiny MANS}}}$–$R_{\textrm{\mbox{out}}}$ inside
which MANS remnants can exist and only if $f/v\gtrsim 5$. When $f/v\lesssim
4$, on the other hand, all remnants either violate the minimum threshold
condition or are radially unstable. Consequently, NS-MNS binaries in the left
panels of Fig. 12 and Fig. 13, corresponding to a dark sector with f/v=3, will
all collapse to black holes after coalescence.

Figure 16 shows the corresponding ranges in mirror baryon fraction
$Y_{\textrm{\mbox{\tiny DM}}}=N_{\textrm{\mbox{\tiny
DM}}}/(N_{\textrm{\mbox{\tiny SM}}}+N_{\textrm{\mbox{\tiny DM}}})$ and total
mass $M_{\textrm{\mbox{\tiny MANS}}}$ for MANS remnants of NS/MNS mergers,
when $f/v=3,5$ and $7$. The mass of these objects is tightly limited to the
regime $(1.2-1.55)\,M_{\odot}$. The DM fraction $Y_{\textrm{\mbox{\tiny DM}}}$
is also on a tight range, from $\sim 15\%$ to $\sim 25\%$. Observe that the
mirror-matter fraction $Y_{\textrm{\mbox{\tiny DM}}}$ is anti-correlated with
the total mass of the star, which we can understand intuitively as follows.
Increasing the amount of DM at the center of the star leads to stronger
gravitational forces. That additional force must be balanced by larger
pressure gradients, which means the pressure in the $\rm MANS$ will drop
faster in the radial direction. Because the energy density is a monotonic
function of pressure, the energy density will also drop faster, leading to a
smaller total mass. Another consequence of Figure 16 is that all MANS produced
from a MNS-NS merger leading to a stable remnant are of the form of a DM core
with a SM halo.

###  VI.3 Black-Hole Remnants

If a stable MANS cannot form, the system will eventually collapse to a black
hole. Because of the lower mass thresholds on the colliding stars, this turns
out to be the most likely scenario. Black-hole remnants from MNS-NS
coalescences, however, would be very light, lying in the mass range
$\sim(1.3,3.7)\,M_{\odot}$. Black holes from MNS-MNS coalescences would be
even lighter, in the interval $\sim(0.8,3.2)\,M_{\odot}$. To explain how we
come to these numbers, we look at the minimum mass $M_{min}$ from $f/v=4$ that
produces no stable remnants. Then, $M_{min}^{f/v=4}\sim 0.4M_{\odot}$ and the
$M_{min}^{\textrm{\mbox{\tiny SM}}}=0.9M_{\odot}$ so our total minimum black
hole mass from a coalescence of a SM NS and a MNS with $f/v=4$ is then
$M_{min}^{BH}\sim 1.3M_{\odot}$. Of course, for larger f/v then smaller
$M_{min}$ are possible. For instance, $M_{min}^{f/v=7}=0.3M_{\odot}$. However,
the coalescence of this MNS with a $M_{min}^{\textrm{\mbox{\tiny
SM}}}=0.9M_{\odot}$ NS leads to a stable configuration (as shown in Fig. 15)
so it would not produce a black hole. The maximum values of a black hole mass
for a MNS-NS coalescence comes from taking the maximum mass for $f/v=2$ and
adding that to the maximum mass of the standard model sequence such that
$M_{max}^{BH}=1.6+2.1=3.7M_{\odot}$. In Sec. VII we will focus specifically on
the difference between $f/v<5$ and $f/v\geq 5$ for two reasons: i.) $f/v\geq
5$ obtain stable configurations whereas $f/v<5$ only produce black hole
remnants ii.) $f/v=5$ also is the smallest $f/v$ that provides no overlap in
the $\eta$ vs $\mathcal{M}$ relation shown in Fig. 13. Thus, for much of our
following calculations, it is important to know that $M_{min}^{f/v=5}\sim 0.4$
and $M_{max}^{f/v=5}\sim 0.9$.

Such a light range of remnant black hole masses would have important
implications for the GWs emitted post-merger. In particular, the GW ringdown,
generated as the remnant settles to its final stable configuration, would
occur at very high frequencies. The GW frequency of the dominant fundamental
ringdown mode scales as $f_{\rm GW}\sim 0.44/(2\pi M_{\rm R})$, where $M_{\rm
R}$ is the mass of the remnant Berti  _et al._ (2009); Yunes  _et al._ (2022).
For a normal black-hole binary merger with component masses
$m_{1}=m_{2}=5M_{\odot}$, one then expects $f_{\rm GW}\sim 1400$ Hz. However,
for a MNS-NS coalescence (with $m_{1}=m_{2}=3M_{\odot}$), the black hole
remnant would ring at $f_{\rm GW}\sim 2300$ Hz, while for a MNS-MNS
coalescence (with $m_{1}=m_{2}=1M_{\odot}$), the black hole remnant would ring
at $f_{\rm GW}\sim 14,000$ Hz. This means that while the ringdown of a MNS-NS
coalescence may be observable, it will be incredibly hard to observe the
ringdown of a MNS-MNS collision. Moreover, if somehow only the ringdown were
observed, then it would be impossible to distinguish this small-mass black
hole remnant from a primordial black hole.

##  VII Observational prospects

Because of the unique signatures for binary MNSs mergers and NS-MNS mergers,
we find that there are multiple ways that one could detect the presence of
mirror matter within NSs. Here we will list all the unique signatures that
underlie various scenarios, summarizing several discussions presented earlier
in this paper.

  * •

Scenario 1 - binary MNS inspiral and merger with $f/v\geq 5$: Both objects are
dark (no electromagnetic signatures) and the symmetric mass ratio and chirp
mass is entirely distinguishable from that of SM NS binaries. The tidal
deformabilities are orders of magnitude smaller than that of SM NSs, and the
remnant is also dark, but may produce a more massive MNS or a black hole. The
black hole remnant would be very light and certainly lower than the maximum
mass of a neutron star (black hole masses starting at $M\sim 0.8M_{\odot}$ for
$f/v=7$).

  * •

Scenario 2 - binary MNS inspiral and merger with $f/v<5$: Since both objects
are dark (no electromagnetic signatures), we need to rely on GW measurements
and need to look into the parameter space of such measurements. For the
parameter space of mass ratio and $M_{bin}$ (or $\eta$ and $\mathcal{M}$),
there is a small overlap between SM NS binaries and MNS binaries. Therefore,
for most binary systems, a MNS binary can be distinguished from a SM NS
binary, and evidence of DM can be found once we detect a binary system in the
MNS-MNS region. Tidal deformabilities are approximately one order of magnitude
smaller than that of SM NSs, and the remnant is dark, but may produce a more
massive MNS or a black hole. In this case, however, the black holes produced
may be smaller than the maximum mass of SM NSs, i.e. the black hole remnants
can be as small as $M\gtrsim 1.3M_{\odot}$.

  * •

Scenario 3 - NS-MNS inspiral and merger with $f/v\geq 5$: One object is dark,
but the NS may produce electromagnetic signatures from tidal disruption. The
symmetric mass ratio and chirp mass are entirely distinguishable from those of
SM NS binaries. The tidal deformabilities of the binary components are orders
of magnitude different from each other. The remnant may produce a MANS or a
black hole. The resulting admixed star has a significantly smaller radius than
that of a SM NS. The black hole remnant may be lighter than the maximum mass
of the SM NS sequence (black hole masses starting at $M\sim 1.4M_{\odot}$).

  * •

Scenario 4 - NS-MNS inspiral and merger with $f/v<5$: One object is dark but
the NS may produce electromagnetic signatures from tidal disruption. The
symmetric mass ratio and chirp mass are mostly distinguishable from those of
SM NS binaries, but some small overlap region exists. The tidal
deformabilities of the binary components are an order of magnitude different
from each other. The remnant will be a black hole. The black hole remnant may
be lighter than the maximum mass of the SM sequence (masses starting at $M\sim
1.3M_{\odot}$).

  * •

Scenario 5 - MANS-MANS inspiral and mergers or NS-MANS inspiral and mergers:
Both objects may produce electromagnetic signatures from tidal disruption. The
symmetric mass ratio and chirp mass may overlap with that of SM NS binaries.
Tidal deformabilities would break universal relations and may appear similar
to twin stars Tan  _et al._ (2022a). The remnant will be a black hole in most
scenarios, unless the dark matter core is extremely small.

  * •

Scenario 6 - Radius measurement of $\rm MANS$ with X-rays: Because a $\rm
MANS$ can have an outer layer (or halo) entirely composed of SM matter, it can
have a hot spot on its surface, just like a SM NS. If the hot spot emits
X-rays, NICER could potentially observe them and obtain the
$M_{\textrm{\mbox{\tiny MANS}}}$ and $R_{\textrm{\mbox{out}}}$. Notice that it
is possible to observe two $\rm MANSs$ with the same mass, but with very
different radii, because $\rm MANSs$ exist in a 2D area in the mass-radius
plane. Therefore, from an observational perspective, MANSs may be confused
with 1st-order phase transition (1st-order PT) twin stars, which also allows
two SM NSs with the same mass but very different radii. The difference between
$\rm MANSs$ and 1st-order PT twin stars (known as mass twins Gerlach (1968);
Kampfer (1981)) is that mass twins are produced by a single SM EoS (one fluid
only), and they lie on a single non-monotonic mass-radius curve, instead of a
non-bijective mapping from
$(\epsilon_{c}^{\mathrm{SM}},\epsilon_{c}^{\mathrm{DM}})$ to the
$(R_{\mathrm{out}},M)$ plane that occurs for MANS. However, it may still be
possible to distinguish the 2D mirror matter plane from mass twins with just
the measurement of two compact objects with the same mass, if they have
drastically different radii. For instance, applying the heavy maximum mass
constraint of $M_{\rm max}\sim 2M_{\odot}$ Cromartie  _et al._ (2019), it is
difficult to produce mass twins; only extremely large first-order phase
transitions can produce a large radius difference between the different
branches (e.g. see Fig. 12 and Fig. 14 from Tan  _et al._ (2022b)). Thus,
admixed mirror stars may initially mask as mass twins, but once a 2D mass-
radius plane is measured, they would be clearly identified as mirror matter.

  * •

Scenario 7 - Black hole-NS Twins: Binary inspirals and mergers involving MNSs
have the potential to create very light black holes. A striking possibility is
that this would lead to black hole-NS mass twins—that is, a black hole and a
NS with the same gravitational mass. In order to clearly distinguish these
objects from each other, one would need a sufficiently good measurement of the
their tidal deformability, since black hole are expected to have none. This
could prove very difficult for very massive NSs, because their deformability,
$\Lambda\sim 10$, would not be sufficiently different from that of black
holes. Therefore, it would be more feasible to distinguish between lighter
black hole-NS twins, say with $M\approx 1.4\,M_{\odot}$.

We have outlined some of the specific scenarios for the inspiral and merger of
compact objects containing mirror matter. However, some signatures will only
appear if a large population of these objects is measured. The primary
signature with large statistics is that if admixed mirror stars exist, then
one could find a 2D plane in the mass-radius or the tidal deformability-mass
relation. The size and shape of the 2D mass-radius plane depend on a number of
factors, such as the value of $f/v$, the minimum mass of a MNSs (i.e. if they
are subject to the Chandrasekhar limit or not), if the admixed stars are
produced due to accreted mirror matter or as a remnant of a MNS-NS merger, and
details of the EOS. However, the existence of the 2D plane is robust, only
possible in a two-fluid model, and would not occur from any known SM EOS.

Thus, if it were to become clear from observations (say with GW detectors such
as advanced LIGO/Virgo/KAGRA, or the X-ray telescope NICER) that the mass-
radius posteriors of compact objects cannot be described through a one-
dimensional sequence alone (even considering disconnected sequences), then
such a measurement would be a clean-cut signature for admixed mirror stars.
Such a measurement would require the confidence ellipses on the mass-radius
(or mass-tidal deformability) plane of two independent observations of at
least two different NSs to be non-overlapping. This could happen, for example,
if one were to observe a compact object with $M=1\,M_{\odot}$ and $R\lesssim
6$ km, with an accuracy of $\delta R\lesssim 4$ km. Black holes of such a low
mass could potentially be distinguishable from admixed mirror star remnants
through the quasi-normal GWs they would emit as they settle down to their
final stationary configuration.

##  VIII Conclusions and Outlook

In this paper, we have investigated the effects of Twin Higgs mirror matter in
isolated MANS, and in NS-MNS and MNS-MNS binaries. We have found that, instead
of forming a single curve, MANS span over an 2D region on the mass-radius
plan. This region ranges from the SM NS mass-radius curve on the left, past
the MNS sequence, to a curve of minimal radius, where
$R_{\textrm{\mbox{out}}}\sim 4$ km, and the star is filled with both DM and SM
matter throughout its extension. Because the outer radius changes non-
monotonically with the DM fraction $Y_{\textrm{\mbox{\tiny DM}}}$ in the star,
different $Y_{\textrm{\mbox{\tiny DM}}}$ intervals can overlap in the mass-
radius diagram, leading to “ _ultimate twins_ ,” i.e. stars with the same
exact mass and radius, but different internal composition. While these twins
have exactly the same mass and radius, their tidal deformability is different
and would distinguish the two compact objects.

We have also explored the properties of MANS that would have formed via DM
capture in NSs above a generous assumed NS formation threshold
$M_{\textrm{\mbox{\tiny SM}}}\gtrsim 0.9\,M_{\odot}$. Surprisingly, we find
that the corresponding mass-radius region is nearly independent of the mirror
Higgs scale, yielding approximately the same region for $f/v=3-7$. If we allow
for the possibility of mirror matter capture, then a DM-admixed NS can support
masses in the same range as a regular NS, but with a smaller radius (as low as
$R\sim 8$ km). Conversely, if we begin with a mirror neutron star and allow
for SM capture, then masses on the small end of NS are possible but the radii
are significantly smaller $R\sim 4$ km. A smoking-gun signature of MANS is the
mass-radius plane that they inhabit, which presents a 2D structure, unlike any
other object studied so far. Indeed, even objects with distinct stable
branches due to first-order phase transition (mass twins) still exist in a
one-dimensional sequence.

We then moved to a discussion of the effects of mirror matter in NS-MNS and
MNS-MNS binary systems. We found that the amplitude of the GWs emitted by
these binary systems (at fixed distance) would be smaller than, but still
comparable to, the one from NS-NS binaries due to the smaller total mass of
MNSs. We also emphasize that, by measuring the chirp mass $\mathcal{M}$ and
the symmetric mass ratio $\eta$ alone, one can distinguish a $\rm MNS$-NS or
$\rm MNS$-$\rm MNS$ system from a NS-NS one. This distinction is clear-cut
only when $f/v\gtrsim 5$. We have also shown that binary Love relations cannot
be naïvely applied to break the degeneracy of a GW measurement from a NS-MNS
system, because that relation is no longer EoS insensitive. These relations
can, however, be applied without change to MNS-MNS binaries.

NS-MNS mergers are of particular interest because they provide a promising way
to form $\rm MANSs$. However, stability against radial oscillations and NS and
MNS mass thresholds impose very stringent constraints on the resulting
remnant. For $f/v<5$, we have found no possibility of stable $\rm MANSs$
remnants. For $f/v\geq 5$, we find a small stable mass-radius region. In cases
where a stable MANS remnant is not formed (the most likely scenario), collapse
to a black hole will follow. In that case, one could potentially measure the
post-merger ringdown.

One very intriguing possibility that arises from the coalescence of MNS-MNS
and NS-MNS is that they can produce very light black-hole remnants, with
masses as low as $0.8\,M_{\odot}$. Thus, depending on the masses of the MNS-
MNS and NS-MNS binaries, it is possible to produce black hole remnants in the
range $(0.8-3.7)\,M_{\odot}$, MNSs remnants in the range
$(0.8-1.5)\,M_{\odot}$, and $\rm MANS$ remnants in the range
$(1.2-1.6)\,M_{\odot}$. Thus, mirror matter provides a mechanism that allows
for black holes, MNSs and $\rm MANS$s remnants - all in the same mass range,
but with very different radii.

Our work opens up a number of new possibilities and questions. It is not clear
what the precise merger and post-merger signals of such a two fluid model
would be like. Would the two-fluid nature of the compact objects affect the
peak frequency of the GWs emitted at merger? What properties would we expect
if spin is considered? Is there a possibility of a Kilonova if $\rm MANS$
merge? In order to answer these questions, a full two-fluid numerical
relativity simulation is required, which we leave to future work.

###### Acknowledgements.

H.T. and N.Y. acknowledge financial support through NSF grants No.
PHY-1759615, PHY-1949838 and PHY-2207650. M.H. was supported in part by the
National Science Foundation (NSF) within the framework of the MUSES
collaboration, under grant number OAC-2103680. The research of DC was
supported in part by a Discovery Grant from the Natural Sciences and
Engineering Research Council of Canada, the Canada Research Chair program, the
Alfred P. Sloan Foundation, the Ontario Early Researcher Award, and the
University of Toronto McLean Award. J.N.H. acknowledges financial support by
the US-DOE Nuclear Science Grant No. DESC0020633. E.D. acknowledges financial
support from Internship Grant Funding – Berea College Internship Program. This
material is based upon work supported by the National Science Foundation and
the Department of Defense under Grant PHY-1950744. Any opinions, findings, and
conclusions or recommendations expressed in this material are those of the
author(s) and do not necessarily reflect the views of the National Science
Foundation or the Department of Defense.

## References

  * Lisanti (2017) M. Lisanti, in  _Theoretical Advanced Study Institute in Elementary Particle Physics: New Frontiers in Fields and Strings_ (2017) pp. 399–446, arXiv:1603.03797 [hep-ph] .
  * Lin (2019) T. Lin, PoS 333, 009 (2019), arXiv:1904.07915 [hep-ph] .
  * Li  _et al._ (2012) X. Li, F. Wang, and K. S. Cheng, JCAP 10, 031 (2012), arXiv:1210.1748 [astro-ph.CO] .
  * Mukhopadhyay and Schaffner-Bielich (2016) P. Mukhopadhyay and J. Schaffner-Bielich, Phys. Rev. D 93, 083009 (2016), arXiv:1511.00238 [astro-ph.HE] .
  * Baym  _et al._ (2018) G. Baym, D. H. Beck, P. Geltenbort, and J. Shelton, Phys. Rev. Lett. 121, 061801 (2018), arXiv:1802.08282 [hep-ph] .
  * Ellis  _et al._ (2018) J. Ellis, G. Hütsi, K. Kannike, L. Marzola, M. Raidal, and V. Vaskonen, Phys. Rev. D 97, 123007 (2018), arXiv:1804.01418 [astro-ph.CO] .
  * Dengler  _et al._ (2022) Y. Dengler, J. Schaffner-Bielich, and L. Tolos, Phys. Rev. D 105, 043013 (2022), arXiv:2111.06197 [astro-ph.HE] .
  * Kain (2021) B. Kain, Phys. Rev. D 103, 043009 (2021), arXiv:2102.08257 [gr-qc] .
  * Sen and Guha (2021) D. Sen and A. Guha, Mon. Not. Roy. Astron. Soc. 504, 3354 (2021), arXiv:2104.06141 [hep-ph] .
  * Jiménez and Fraga (2022) J. C. Jiménez and E. S. Fraga, Universe 8, 34 (2022), arXiv:2111.00091 [hep-ph] .
  * Miao  _et al._ (2022) Z. Miao, Y. Zhu, A. Li, and F. Huang, Astrophys. J. 936, 69 (2022), arXiv:2204.05560 [astro-ph.HE] .
  * Horowitz and Reddy (2019) C. J. Horowitz and S. Reddy, Phys. Rev. Lett. 122, 071102 (2019), arXiv:1902.04597 [astro-ph.HE] .
  * Maedan (2020) S. Maedan, PTEP 2020, 033B07 (2020), arXiv:1908.00711 [hep-ph] .
  * Narain  _et al._ (2006) G. Narain, J. Schaffner-Bielich, and I. N. Mishustin, Phys. Rev. D 74, 063003 (2006), arXiv:astro-ph/0605724 .
  * Goldman  _et al._ (2013) I. Goldman, R. N. Mohapatra, S. Nussinov, D. Rosenbaum, and V. Teplitz, Phys. Lett. B 725, 200 (2013), arXiv:1305.6908 [astro-ph.CO] .
  * Tolos  _et al._ (2015) L. Tolos, J. Schaffner-Bielich, and Y. Dengler, Phys. Rev. D 92, 123002 (2015), [Erratum: Phys.Rev.D 103, 109901 (2021)], arXiv:1507.08197 [astro-ph.HE] .
  * Kouvaris and Nielsen (2015) C. Kouvaris and N. G. Nielsen, Phys. Rev. D 92, 063526 (2015), arXiv:1507.00959 [hep-ph] .
  * Mukhopadhyay  _et al._ (2017) S. Mukhopadhyay, D. Atta, K. Imam, D. N. Basu, and C. Samanta, Eur. Phys. J. C 77, 440 (2017), [Erratum: Eur.Phys.J.C 77, 553 (2017)], arXiv:1612.07093 [nucl-th] .
  * Gresham and Zurek (2019) M. I. Gresham and K. M. Zurek, Phys. Rev. D 99, 083008 (2019), arXiv:1809.08254 [astro-ph.CO] .
  * Collier  _et al._ (2022) M. Collier, D. Croon, and R. K. Leane, (2022), arXiv:2205.15337 [gr-qc] .
  * Das  _et al._ (2022) A. Das, T. Malik, and A. C. Nayak, Phys. Rev. D 105, 123034 (2022), arXiv:2011.01318 [nucl-th] .
  * Leung  _et al._ (2022) K.-L. Leung, M.-c. Chu, and L.-M. Lin, Phys. Rev. D 105, 123010 (2022), arXiv:2207.02433 [astro-ph.HE] .
  * Rutherford  _et al._ (2022) N. Rutherford, G. Raaijmakers, C. Prescod-Weinstein, and A. Watts, (2022), arXiv:2208.03282 [astro-ph.HE] .
  * Giangrandi  _et al._ (2022) E. Giangrandi, V. Sagun, O. Ivanytskyi, C. Providência, and T. Dietrich, (2022), arXiv:2209.10905 [astro-ph.HE] .
  * Deliyergiyev  _et al._ (2019) M. Deliyergiyev, A. Del Popolo, L. Tolos, M. Le Delliou, X. Lee, and F. Burgio, Phys. Rev. D 99, 063015 (2019), arXiv:1903.01183 [gr-qc] .
  * Alexander  _et al._ (2016) J. Alexander _et al._ (2016) arXiv:1608.08632 [hep-ph] .
  * Curtin and Setford (2020) D. Curtin and J. Setford, JHEP 03, 041 (2020), arXiv:1909.04072 [hep-ph] .
  * Curtin and Setford (2021) D. Curtin and J. Setford, JHEP 03, 166 (2021), arXiv:2010.00601 [hep-ph] .
  * Hippert  _et al._ (2022) M. Hippert, J. Setford, H. Tan, D. Curtin, J. Noronha-Hostler, and N. Yunes, Phys. Rev. D 106, 035025 (2022), arXiv:2103.01965 [astro-ph.HE] .
  * Ryan and Radice (2022) M. Ryan and D. Radice, Phys. Rev. D 105, 115034 (2022), arXiv:2201.05626 [astro-ph.HE] .
  * Blinnikov and Khlopov (1982) S. I. Blinnikov and M. Y. Khlopov, Sov. J. Nucl. Phys. 36, 472 (1982).
  * Blinnikov and Khlopov (1983) S. I. Blinnikov and M. Khlopov, Sov. Astron. 27, 371 (1983).
  * Khlopov  _et al._ (1991) M. Y. Khlopov, G. M. Beskin, N. E. Bochkarev, L. A. Pustylnik, and S. A. Pustylnik, Sov. Astron. 35, 21 (1991).
  * Berezhiani  _et al._ (1996) Z. G. Berezhiani, A. D. Dolgov, and R. N. Mohapatra, Phys. Lett. B 375, 26 (1996), arXiv:hep-ph/9511221 .
  * Chacko  _et al._ (2006) Z. Chacko, H.-S. Goh, and R. Harnik, Phys. Rev. Lett. 96, 231802 (2006), arXiv:hep-ph/0506256 .
  * Chacko  _et al._ (2017) Z. Chacko, N. Craig, P. J. Fox, and R. Harnik, JHEP 07, 023 (2017), arXiv:1611.07975 [hep-ph] .
  * Chacko  _et al._ (2018) Z. Chacko, D. Curtin, M. Geller, and Y. Tsai, JHEP 09, 163 (2018), arXiv:1803.03263 [hep-ph] .
  * Burdman  _et al._ (2015) G. Burdman, Z. Chacko, R. Harnik, L. de Lima, and C. B. Verhaaren, Phys. Rev. D 91, 055007 (2015), arXiv:1411.3310 [hep-ph] .
  * Craig  _et al._ (2015) N. Craig, A. Katz, M. Strassler, and R. Sundrum, JHEP 07, 105 (2015), arXiv:1501.05310 [hep-ph] .
  * Aad  _et al._ (2021) G. Aad _et al._ (ATLAS), JHEP 04, 165 (2021), arXiv:2102.01444 [hep-ex] .
  * Blanchet (2014) L. Blanchet, Living Rev. Rel. 17, 2 (2014), arXiv:1310.1528 [gr-qc] .
  * Kain (2020) B. Kain, Phys. Rev. D 102, 023001 (2020), arXiv:2007.04311 [gr-qc] .
  * Naidu  _et al._ (2021) N. F. Naidu, S. Carloni, and P. Dunsby, Phys. Rev. D 104, 044014 (2021), arXiv:2102.05693 [gr-qc] .
  * Gleason  _et al._ (2022) T. Gleason, B. Brown, and B. Kain, Phys. Rev. D 105, 023010 (2022), arXiv:2201.02274 [gr-qc] .
  * Strobel  _et al._ (1999) K. Strobel, C. Schaab, and M. K. Weigel, Astron. Astrophys. 350, 497 (1999), arXiv:astro-ph/9908132 .
  * Doroshenko  _et al._ (2022) V. Doroshenko, V. Suleimanov, G. Pühlhofer, and A. Santangelo, Nature Astronomy (2022), 10.1038/s41550-022-01800-1.
  * Catuneanu  _et al._ (2013) A. Catuneanu, C. O. Heinke, G. R. Sivakoff, W. C. G. Ho, and M. Servillat, Astrophys. J. 764, 145 (2013), arXiv:1301.3768 [astro-ph.HE] .
  * Servillat  _et al._ (2012) M. Servillat, C. O. Heinke, W. C. G. Ho, J. E. Grindlay, J. Hong, M. v. d. Berg, and S. Bogdanov, Mon. Not. Roy. Astron. Soc. 423, 1556 (2012), arXiv:1203.5807 [astro-ph.HE] .
  * Sandin and Ciarcelluti (2009) F. Sandin and P. Ciarcelluti, Astropart. Phys. 32, 278 (2009), arXiv:0809.2942 [astro-ph] .
  * Ciarcelluti and Sandin (2011) P. Ciarcelluti and F. Sandin, Phys. Lett. B 695, 19 (2011), arXiv:1005.0857 [astro-ph.HE] .
  * Emma  _et al._ (2022) M. Emma, F. Schianchi, F. Pannarale, V. Sagun, and T. Dietrich, Particles 5, 273 (2022), arXiv:2206.10887 [gr-qc] .
  * Hanhart  _et al._ (2008) C. Hanhart, J. R. Pelaez, and G. Rios, Phys. Rev. Lett. 100, 152001 (2008), arXiv:0801.2871 [hep-ph] .
  * McNeile  _et al._ (2009) C. McNeile, C. Michael, and C. Urbach (ETM), Phys. Lett. B 674, 286 (2009), arXiv:0902.3897 [hep-lat] .
  * Walker-Loud  _et al._ (2009) A. Walker-Loud _et al._ , Phys. Rev. D 79, 054502 (2009), arXiv:0806.4549 [hep-lat] .
  * Bratt  _et al._ (2010) J. D. Bratt _et al._ (LHPC), Phys. Rev. D 82, 094502 (2010), arXiv:1001.3620 [hep-lat] .
  * Syritsyn  _et al._ (2010) S. N. Syritsyn _et al._ , Phys. Rev. D 81, 034507 (2010), arXiv:0907.4194 [hep-lat] .
  * Albaladejo and Oller (2012) M. Albaladejo and J. A. Oller, Phys. Rev. D 86, 034003 (2012), arXiv:1205.6606 [hep-ph] .
  * Horsley  _et al._ (2014) R. Horsley, Y. Nakamura, A. Nobile, P. E. L. Rakow, G. Schierholz, and J. M. Zanotti, Phys. Lett. B 732, 41 (2014), arXiv:1302.2233 [hep-lat] .
  * Miller  _et al._ (2019) M. C. Miller _et al._ , Astrophys. J. Lett. 887, L24 (2019), arXiv:1912.05705 [astro-ph.HE] .
  * Riley  _et al._ (2019) T. E. Riley _et al._ , Astrophys. J. Lett. 887, L21 (2019), arXiv:1912.05702 [astro-ph.HE] .
  * Miller  _et al._ (2021) M. C. Miller _et al._ , Astrophys. J. Lett. 918, L28 (2021), arXiv:2105.06979 [astro-ph.HE] .
  * Riley  _et al._ (2021) T. E. Riley _et al._ , Astrophys. J. Lett. 918, L27 (2021), arXiv:2105.06980 [astro-ph.HE] .
  * Abbott  _et al._ (2017) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 119, 161101 (2017), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ (2018a) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 121, 161101 (2018a), arXiv:1805.11581 [gr-qc] .
  * Abbott  _et al._ (2019) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. X 9, 011001 (2019), arXiv:1805.11579 [gr-qc] .
  * Baym  _et al._ (1971) G. Baym, C. Pethick, and P. Sutherland, Astrophys. J. 170, 299 (1971).
  * Aad  _et al._ (2020) G. Aad _et al._ (ATLAS), Phys. Rev. D 101, 012002 (2020), arXiv:1909.02845 [hep-ex] .
  * Yagi and Yunes (2013) K. Yagi and N. Yunes, Phys. Rev. D 88, 023009 (2013), arXiv:1303.1528 [gr-qc] .
  * Nelson  _et al._ (2019) A. Nelson, S. Reddy, and D. Zhou, JCAP 07, 012 (2019), arXiv:1803.03266 [hep-ph] .
  * Rafiei Karkevandi  _et al._ (2021) D. Rafiei Karkevandi, S. Shakeri, V. Sagun, and O. Ivanytskyi, in  _16th Marcel Grossmann Meeting on Recent Developments in Theoretical and Experimental General Relativity, Astrophysics and Relativistic Field Theories_ (2021) arXiv:2112.14231 [astro-ph.HE] .
  * Hinderer (2008) T. Hinderer, Astrophys. J. 677, 1216 (2008), arXiv:0711.2420 [astro-ph] .
  * Comer  _et al._ (1999) G. L. Comer, D. Langlois, and L. M. Lin, Phys. Rev. D 60, 104025 (1999), arXiv:gr-qc/9908040 .
  * Leung  _et al._ (2012) S. C. Leung, M. C. Chu, and L. M. Lin, Phys. Rev. D 85, 103528 (2012), arXiv:1205.1909 [astro-ph.CO] .
  * Chandrasekhar (1964) S. Chandrasekhar, Astrophys. J. 140, 417 (1964), [Erratum: Astrophys.J. 140, 1342 (1964)].
  * Henriques  _et al._ (1990) A. B. Henriques, A. R. Liddle, and R. G. Moorhouse, Phys. Lett. B 251, 511 (1990).
  * Nyhan and Kain (2022) J. E. Nyhan and B. Kain, Phys. Rev. D 105, 123016 (2022), arXiv:2206.07715 [gr-qc] .
  * Weinberg (1972) S. Weinberg,  _Gravitation and Cosmology: Principles and Applications of the General Theory of Relativity_ (John Wiley and Sons, New York, 1972).
  * Glendenning (2000) N. K. Glendenning,  _Compact Stars_ (Springer, New York, 2000).
  * Flanagan and Hinderer (2008) E. E. Flanagan and T. Hinderer, Phys. Rev. D 77, 021502 (2008), arXiv:0709.1915 [astro-ph] .
  * Favata (2014) M. Favata, Phys. Rev. Lett. 112, 101101 (2014), arXiv:1310.8288 [gr-qc] .
  * Yagi and Yunes (2014) K. Yagi and N. Yunes, Phys. Rev. D 89, 021303 (2014), arXiv:1310.8358 [gr-qc] .
  * Yagi and Yunes (2016) K. Yagi and N. Yunes, Class. Quant. Grav. 33, 13LT01 (2016), arXiv:1512.02639 [gr-qc] .
  * Abbott  _et al._ (2018b) B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 121, 161101 (2018b), arXiv:1805.11581 [gr-qc] .
  * Berti  _et al._ (2009) E. Berti, V. Cardoso, and A. O. Starinets, Class. Quant. Grav. 26, 163001 (2009), arXiv:0905.2975 [gr-qc] .
  * Yunes  _et al._ (2022) N. Yunes, M. C. Miller, and K. Yagi, Nature Rev. Phys. 4, 237 (2022), arXiv:2202.04117 [gr-qc] .
  * Tan  _et al._ (2022a) H. Tan, V. Dexheimer, J. Noronha-Hostler, and N. Yunes, Phys. Rev. Lett. 128, 161101 (2022a), arXiv:2111.10260 [astro-ph.HE] .
  * Gerlach (1968) U. H. Gerlach, Phys. Rev. 172, 1325 (1968).
  * Kampfer (1981) B. Kampfer, J. Phys. A 14, L471 (1981).
  * Cromartie  _et al._ (2019) H. T. Cromartie _et al._ (NANOGrav), Nature Astron. 4, 72 (2019), arXiv:1904.06759 [astro-ph.HE] .
  * Tan  _et al._ (2022b) H. Tan, T. Dore, V. Dexheimer, J. Noronha-Hostler, and N. Yunes, Phys. Rev. D 105, 023018 (2022b), arXiv:2106.03890 [astro-ph.HE] .