# Dynamical Instability of Collapsed Dark Matter Halos

Wei-Xiang Feng  Hai-Bo Yu  and Yi-Ming Zhong

###### Abstract

A self-interacting dark matter halo can experience gravothermal collapse,
resulting in a central core with an ultrahigh density. It can further contract
and collapse into a black hole, a mechanism proposed to explain the origin of
supermassive black holes. We study dynamical instability of the core in
general relativity. We use a truncated Maxwell-Boltzmann distribution to model
the dark matter distribution and solve the Tolman-Oppenheimer-Volkoff
equation. For given model parameters, we obtain a series of equilibrium
configurations and examine their dynamical instability based on considerations
of total energy, binding energy, fractional binding energy, and adiabatic
index. Our numerical results indicate that the core can collapse into a black
hole when the fractional binding energy reaches $0.035$ with a central
gravitational redshift of $0.5$. We further show for the instability to occur
in the classical regime, the boundary temperature of the core should be at
least $10\%$ of the mass of dark matter particles; for a $10^{9}~{}{\rm
M_{\odot}}$ seed black hole, the particle mass needs to be larger than a few
keV. These results can be used to constrain different collapse models, in
particular, those with dissipative dark matter interactions. .

##  1 Introduction

The study of dynamical instability of a self-gravitating system and its
collapse to a black hole has a long history [1, 2]. Early work analyzed the
evolution of stellar clusters in general relativity and examined conditions
for their relativistic instability with linear perturbation theory [3, 4, 5].
The techniques and tools of numerical relativity and N-body simulations were
further developed in [6, 7, 8], which can be used to trace full evolution in
the nonlinear regime.

Recently, we proposed a scenario to explain the origin of supermassive black
holes in the early universe [9]; see also [10, 11, 12]. This is based on the
mechanism that a self-interacting dark matter halo can experience gravothermal
collapse. Dark matter self-interactions can thermalize the inner halo over
cosmological timescales [13, 14, 15, 16, 17, 18, 19, 20]; see [21] for a
review. As a self-gravitating system with a finite size, the halo has negative
heat capacity, and the self-interactions transport heat from the central
region at late stages of the evolution, resulting in a core with an ultrahigh
density [22, 23, 24, 25]. The core can further contact and collapse into a
seed black hole [26], which would grow into a supermassive one by accreting
baryonic matter. We used a semi-analytical method and derived the condition
for triggering dynamical instability of the core. Following Chandrasekhar’s
criterion [1], i.e., requiring the pressure averaged adiabatic index of the
gravothermal system to be less than its critical adiabatic index, we found the
instability occurs when the 3D central velocity dispersion of dark matter
particles reaches $\sim 0.57c$ at which the adiabatic index is $1.62$.

In this work, we systematically study the dynamical instability of a collapsed
halo. We use a truncated Maxwell-Boltzmann distribution to model the dark
matter distribution near the relativistic limit. This is well motivated, as
the self-interactions thermalize dark matter particles. In addition, the core
is gravitationally bound and particles with a sufficiently high velocity can
evaporate and escape from the gravitational pull of the core. We then
implement the distribution with the Tolman-Oppenheimer-Volkoff equation [27,
28] and find a series of equilibrium solutions. For each of them, we evaluate
its thermal dynamical properties and test its instability. Besides the
Chandrasekhar’s criterion, we will use the turning-point method [29, 30, 31]
to examine instability conditions based on considerations of total energy,
binding energy, and fractional binding energy, as illustrated in Figure 1
schematically.

We will compare our numerical results to those from relativistic N-body
simulations [7] and show that the agreement is excellent, i.e., they all
indicate that the system can collapse into a black hole when the fractional
binding energy reaches $0.035$ with a central gravitational redshift of $0.5$.
Thus the method developed in this work may have broad applications as it is
computationally inexpensive. We will further study conditions for the
classical Maxwell-Boltzmann distribution to be valid, and discuss their
implications for constraining models proposed to explain the origin of
supermassive black holes via the gravothermal collapse of dark matter halos.
In particular, we show that although the presence of dissipative interactions
could help speed up the gravothermal evolution of a halo, they may make it
difficult for the core to eventually collapse into a black hole because of
energy loss.

The paper is organized as follows: We present the classical truncated Maxwell-
Boltzmann distribution and its Tolman-Oppenheimer-Volkoff equation in Sec. 2.
We discuss instability conditions and numerical results in Sec. 3. We study
conditions for the classical distribution to be valid and constraints on dark
matter models in Sec. 4, discuss connections with the nonrelativistic fluid
model in Sec. 5, and conclude in Sec. 6.

Figure 1: Schematic illustration of the formation of a seed black hole via the
gravothermal collapse of a self-interacting dark matter halo. At late stages
of gravothermal evolution, the halo can be divided into two regimes, i.e., a
collapsed central core with an ultrahigh density (orange) and a cuspy outer
envelope (gray). As it further contracts, the total mass of the collapsed core
remains almost constant. The elliptical circles denote the sequence of
dynamical instability conditions when the core collapses into a seed black
hole.

##  2 The truncated Maxwell-Boltzmann model

We treat the high-density central region of a collapsed SIDM halo as a
gravitationally bound system. Since dark matter particles with sufficiently
high energies will evaporate and move to the outer envelope, it is natural to
introduce a distribution function with an energy cutoff to model the system.
In this work, we take a truncated Maxwell-Boltzmann distribution, based on
Michie-King models [32, 33]. Consider the following general form [34, 35]

$f(\epsilon\leq\epsilon_{c})=\frac{1-e^{(\epsilon-\epsilon_{c})/k_{B}T}}{e^{(\epsilon-\mu)/k_{B}T}-\eta},\quad
f(\epsilon>\epsilon_{c})=0,$ (2.1)

where $\epsilon$ is the kinetic energy, $\epsilon_{c}$ the cutoff energy, $T$
the temperature and $\mu$ the chemical potential. And they are a function of
radius. The number factor $\eta$ is $+1$ and $-1$ for bosons and fermions,
respectively, and $k_{B}$ is the Boltzmann constant. Note
$\epsilon=\sqrt{|{\bf p}|^{2}c^{2}+m^{2}c^{4}}-mc^{2}$, where $|{\bf p}|$ is
the momentum and $m$ the mass of dark matter particles, and we have subtracted
the rest mass in defining $\mu$. For a dilute gas of classical particles,
$\epsilon-\mu\gg k_{B}T$, the distribution function reduces to the truncated
Maxwell-Boltzmann form

$f(\epsilon\leq\epsilon_{c})=e^{\mu/k_{B}T}(e^{-\epsilon/k_{B}T}-e^{-\epsilon_{c}/k_{B}T}),\quad
f(\epsilon>\epsilon_{c})=0.$ (2.2)

We introduce the following dimensionless variables [36],
$w(r)\equiv\epsilon_{c}(r)/k_{B}T(r)$, $\alpha(r)\equiv\mu(r)/k_{B}T(r)$, and
$b\equiv k_{B}T(R)/mc^{2}$, where $R$ is the boundary radius of the system.
Following the Tolman-Klein law [37, 38] for a gravothermal system, we have the
relation $w(r)=\alpha(r)-\alpha(R)$. The temperature at a given radius is
related to the one at $r=R$ as $T(r)=T(R)/[1-bw(r)]$. It indicates that the
system does not follow an isothermal distribution globally in general
relativity, although it can be achieved locally.

Given the distribution function, one can readily derive the equation of state
and express the number density $n$, energy density $\rho$, thermal energy
density $u$, and pressure $p$ as

$\displaystyle n(r)=$ $\displaystyle 4\sqrt{2}\pi
gm^{3}(c^{3}/h^{3})e^{\alpha(R)}I_{n}(b,w),$ (2.3) $\displaystyle\rho(r)=$
$\displaystyle 4\sqrt{2}\pi gm^{4}(c^{3}/h^{3})e^{\alpha(R)}I_{\rho}(b,w),$
$\displaystyle u(r)=$ $\displaystyle 4\sqrt{2}\pi
gm^{4}(c^{5}/h^{3})e^{\alpha(R)}I_{u}(b,w),$ $\displaystyle p(r)=$
$\displaystyle(8/3)\sqrt{2}\pi gm^{4}(c^{5}/h^{3})e^{\alpha(R)}I_{p}(b,w),$

respectively, where $g=2s+1$ is the spin multiplicity of dark matter
particles, $h$ the Planck constant as a normalization factor, and $c$ the
speed of light; the $I(b,w)$ functions stand for integrals of [35]

$\displaystyle I_{n}(b,w)\equiv{}$
$\displaystyle\left(\frac{b}{1-bw}\right)^{3/2}\int^{w}_{0}(e^{w-x}-1)\left(1+\frac{bx/2}{1-bw}\right)^{1/2}\left(1+\frac{bx}{1-bw}\right)x^{1/2}\text{d}x,$
(2.4) $\displaystyle I_{\rho}(b,w)\equiv{}$
$\displaystyle\left(\frac{b}{1-bw}\right)^{3/2}\int^{w}_{0}(e^{w-x}-1)\left(1+\frac{bx/2}{1-bw}\right)^{1/2}\left(1+\frac{bx}{1-bw}\right)^{2}x^{1/2}\text{d}x,$
$\displaystyle I_{u}(b,w)\equiv{}$
$\displaystyle\left(\frac{b}{1-bw}\right)^{5/2}\int^{w}_{0}(e^{w-x}-1)\left(1+\frac{bx/2}{1-bw}\right)^{1/2}\left(1+\frac{bx}{1-bw}\right)x^{3/2}\text{d}x,$
$\displaystyle I_{p}(b,w)\equiv{}$
$\displaystyle\left(\frac{b}{1-bw}\right)^{5/2}\int^{w}_{0}(e^{w-x}-1)\left(1+\frac{bx/2}{1-bw}\right)^{3/2}x^{3/2}\text{d}x.$

For the model we consider, the Tolman-Oppenheimer-Volkoff equation can be
written as

$\frac{\text{d}{M}}{\text{d}r}=4\pi
r^{2}\rho,\quad\frac{\text{d}w}{\text{d}r}=-\frac{G}{rc^{2}}\left(\frac{1-bw}{b}\right)\frac{4\pi
pr^{3}+{Mc^{2}}}{rc^{2}-2GM},$ (2.5)

where $G$ is the Newton constant, $M(r)$ is the enclosed mass at radius $r$,
and $\rho(r)$ the density. We impose the following boundary conditions: $M=0$
and $w=w(0)$ at $r=0$; $M=M(R)$ and $w=0$ at $r=R$. To further simplify the
calculation, we introduce a fiducial length scale defined as [35]

$\zeta=\lambda_{\rm
C}\left(\frac{m_{\text{Pl}}}{m}\right)\left(\frac{8\pi^{3}}{ge^{\alpha(R)}}\right)^{1/2}~{}{\rm
with}\quad r=\zeta\hat{r},$ (2.6)

where $m_{\text{Pl}}=(\hbar c/G)^{1/2}$ is the Planck mass and $\lambda_{\rm
C}=\hbar/mc$ the Compton wavelength of the particle. With the fiducial length,
we can express thermal dynamical quantities of the system using their
corresponding dimensionless counterpart denoted with a “hat” as
$n=(c^{2}/Gm\zeta^{2})\hat{n}$, $\rho=(c^{2}/G\zeta^{2})\hat{\rho}$,
$u=(c^{4}/G\zeta^{2})\hat{u}$, $p=(c^{4}/G\zeta^{2})\hat{p}$ and
$M=(c^{2}\zeta/G)\hat{M}$, and uniquely determine their profiles for a given
set of $b$ and $w(0)$.

Figure 2: 3D velocity dispersion (left) and adiabatic index $\gamma$ (right)
vs. normalized cutoff energy $bw=\epsilon_{c}/(mc^{2}+\epsilon_{c})$ for
$b=k_{B}T(R)/mc^{2}=(0.5,~{}0.1,~{}0.03,~{}0.01,~{}0.001$). The dashed
horizontal lines denote $v/c=0.57$ (left) and $\gamma=1.59$ (right), at which
the system approaches the relativistic regime and dynamical instability may
occur.

To trigger the onset of dynamical instability, the system needs to be in the
relativistic limit. This requirement puts constraints on $b$ and $w(0)$. It is
useful to consider the product of $b$ and $w(r)$, the normalized cutoff energy

$bw(r)=\frac{\epsilon_{c}(r)/mc^{2}}{1+\epsilon_{c}(r)/mc^{2}},$ (2.7)

where we have used the relation
$T(r)=\left[1+\epsilon_{c}(r)/mc^{2}\right]T(R)$. In the ultrarelativistic
limit $\epsilon_{c}\gg mc^{2}$, $bw\rightarrow 1$. In the opposite limit,
$bw\rightarrow 0$. Since $b=k_{B}T(R)/mc^{2}$ determines the temperature at
the core boundary, a higher $b$ value indicates a hotter thermal bath. $w(r)$
is related to the cutoff energy as $w(r)=\epsilon_{c}(r)/k_{B}T(r)$. For $w\gg
1$, the distribution reduces to the usual Maxwell-Boltzmann form without a
truncation, as indicated in equation (2.2).

The 3D velocity dispersion
$v(r)\equiv\sqrt{{3p}/{\rho}}=c\sqrt{{2I_{p}(b,w)}/{I_{\rho}(b,w)}}$ also
characterizes the relativistic extent of the system. In [9], we showed that
when $v(0)$ approaches $0.57c$, dynamical instability can be triggered. In
addition, the adiabatic index of the system

$\gamma(r)=1+\frac{p}{u}=1+\frac{2}{3}\frac{I_{p}(b,w)}{I_{u}(b,w)}$ (2.8)

approaches $4/3$ and $5/3$ in ultra- and nonrelativistic limits, respectively;
see Appendix A for detailed derivation.

To gain insight into instability conditions, we first use equations (2.3) and
(2.4), and find the radially independent equation of state without imposing
constraints from equation (2.5). We choose $b=0.001,\;0.01,\;0.03,\;0.1$, and
$0.5$ and evaluate $v$ and $\gamma$ as a function of $bw$, as shown in the
left and right panels of Figure 2, respectively. We see that it is easier to
reach the relativistic limit if $b$ is higher, i.e., a hotter thermal bath.
For $b=0.5$, $bw\sim 0.3$ for achieving $v\sim 0.57c$ at which $\gamma\sim
1.59$. But when $b$ decreases to $10^{-3}$, $bw$ needs to be close to $1$. For
the latter case, both $v$ and $\gamma$ are hardly changed over the wide range
of $bw$. For low $b$, the equation of state is stiff, indicating that the
system is hard to compress and reach instability. The results shown in Figure
2 provide guidance in choosing boundary conditions as we discuss further in
the next section.

##  3 Dynamical instability

The instability of a self-gravitating spherical system sets in when the
gravitational energy becomes comparable to its mass energy, $GM^{2}/R\sim
Mc^{2}$, where $M$ and $R$ are the total mass and characteristic radius,
respectively. To be more quantitative, one often considers the compactness of
a system calculated as $C={GM}/{c^{2}R}$ [39]. For a typical neutron star,
$M\sim 2~{}{\rm M_{\odot}}$ and $R\sim 12~{}{\rm km}$, we have $C\sim 0.25$.
For a $M\sim 10^{9}~{}{\rm M_{\odot}}$ black hole, $R\sim 2GM/c^{2}\sim
3\times 10^{9}~{}{\rm km}$ and $C\sim 0.5$. Its average density is $\sim
6\times 10^{2}~{}{\rm kg/m^{3}}$, lower than the water density. As we will
discuss later, for a gaseous sphere with $M\sim 10^{9}~{}{\rm M_{\odot}}$,
$C\sim 0.04$ when the instability is triggered. Overall, the heavier the
system, the easier for it to collapse.

###  3.1 The adiabatic index

A more concrete way to determine the system’s instability is by checking its
adiabatic index. Chandrasekhar first derived the instability conditions for a
spherical system in the context of general relativity [1]. Consider the
background metric
$\text{d}s^{2}=g_{\alpha\beta}\text{d}x^{\alpha}\text{d}x^{\beta}=-e^{2\Phi(r)}c^{2}\text{d}t^{2}+e^{2\Lambda(r)}\text{d}r^{2}+r^{2}\text{d}\Omega$,
where $\Phi(r)$ and $\Lambda(r)$ satisfy the following conditions

$e^{2\Phi(r)}=\exp\left(\int^{\infty}_{r}\frac{4p(r^{\prime})r^{\prime
3}+M(r^{\prime})c^{2}}{r^{\prime}(r^{\prime}c^{2})-2GM(r^{\prime})}\text{d}r^{\prime}\right),~{}e^{2\Lambda(r)}=\left(1-\frac{2GM(r)}{c^{2}r}\right)^{-1}.$
(3.1)

The pulsation equation of a perfect fluid is given by [1]

$\displaystyle\omega^{2}e^{2(\Lambda-\Phi)}\left(\rho+\frac{p}{c^{2}}\right){\xi}={}$
$\displaystyle\frac{4}{r}\frac{\text{d}p}{\text{d}r}\xi-e^{-(2\Phi+\Lambda)}\frac{\text{d}}{\text{d}r}\bigg{[}e^{3\Phi+\Lambda}\frac{\gamma
p}{r^{2}}\frac{\text{d}(r^{2}e^{-\Phi}\xi)}{\text{d}r}\bigg{]}$
$\displaystyle+\frac{8\pi
G}{c^{2}}e^{2\Lambda}p\left(\rho+\frac{p}{c^{2}}\right)\xi-\frac{1}{\rho
c^{2}+p}\bigg{(}\frac{\text{d}p}{\text{d}r}\bigg{)}^{2}\xi,$ (3.2)

where $\xi$ is the Lagrangian displacement and $\omega$ is its corresponding
oscillation frequency. If $\omega^{2}<0$, the Lagrangian displacement $\xi$
receives unbounded growth and the self-gravitating system becomes unstable.
The boundary conditions for equation (3.1) are $\xi(0)=0$ and $p(R)=0$.

Choosing the Lagrangian displacement $\xi=re^{\Phi}$ that satisfies the
boundary conditions, one can show that the system becomes unstable if the
pressure-averaged adiabatic index $\left<\gamma\right>$ is less than the
critical adiabatic index $\gamma_{\rm cr}$

$\langle\gamma\rangle\equiv\frac{\int^{R}_{0}\gamma
e^{3\Phi+\Lambda}p\text{d}^{3}r}{\int^{R}_{0}e^{3\Phi+\Lambda}p\text{d}^{3}r}<\gamma_{\rm
cr},$ (3.3)

where

$\displaystyle\gamma_{{\rm cr}}\equiv$
$\displaystyle\frac{4}{3}+\frac{1}{36}\frac{\int^{R}_{0}e^{3\Phi+\Lambda}[16p+(e^{2\Lambda}-1)(\rho
c^{2}+p)](e^{2\Lambda}-1)r^{2}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}$
(3.4) $\displaystyle+\frac{4\pi
G}{9c^{2}}\frac{\int^{R}_{0}e^{3(\Phi+\Lambda)}[8p+(e^{2\Lambda}+1)(\rho
c^{2}+p)]pr^{4}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}+\frac{16\pi^{2}G^{2}}{9c^{4}}\frac{\int^{R}_{0}e^{3\Phi+5\Lambda}(\rho
c^{2}+p)p^{2}r^{6}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r};$

see Appendix B. The choice of $\xi$ is not unique, but the result is not
sensitive to the particular form of $\xi$ as long as the boundary conditions
are satisfied [1]. In the limit $\rho\gg p/c^{2}$ and $\Phi,\Lambda\rightarrow
0$, the pulsation equation (3.1) reduces to its Newtonian form and
$\gamma_{\rm cr}=4/3$; see Appendix C for a heuristic derivation in the
Newtonian limit. For a monatomic ideal gas, $4/3<\left<\gamma\right><5/3$ [40,
41, 39]. Thus the instability could hardly occur in the context of Newtonian
gravity. On the other hand, in general relativity $\gamma_{\rm cr}$ increases
due to relativistic corrections, which are ${\cal O}(p/\rho c^{2})$. As a
result, the spherical system can reach the dynamical instability condition
$\langle\gamma\rangle<\gamma_{{\rm cr}}$ before the particles become
ultrarelativistic, i.e., $\langle\gamma\rangle\rightarrow 4/3$ as
$p\rightarrow\rho c^{2}/3$.

###  3.2 The turning-point method

Aside from the instability condition based on the adiabatic index, we will
also use the turning-point method [29, 30, 31] and show the former could be
conservative, i.e., dynamical instability could occur before the condition
shown in equation (3.3) is satisfied. Once the boundary temperature parameter
$b$ is fixed, the equation of state only depends on one parameter, i.e., the
central energy cutoff $w(0)$. We can define an energy functional $S$ as a
function of those two variables,

$S=S\left[b,w(0)\right].$ (3.5)

The turning-point ansatz [31] states, for a fixed $b$ value, the marginally
stable configuration reaches at

$\frac{\partial S}{\partial w(0)}\bigg{|}_{b}=0$ (3.6)

This tuning point separates the stable and unstable branches along the one
parameter sequence of $w(0)$. The turning-point method has been applied to
study various stellar systems [42, 43, 44, 45] and gaseous spheres [46, 47].
Our application is similar to those to a rotating relativistic star [44] with
its angular momentum being replaced by $b$ in our model.

We consider $S=\\{E,~{}B,~{}\varepsilon\\}$, where $E$, $B$, $\varepsilon$ are
total energy, binding energy, fractional binding energy, respectively. The
total energy $E=Mc^{2}$ is associated with the Schwarzschild mass $M$ of the
sphere, the binding energy $B=E_{\text{rest}}-E$, and $E_{\rm rest}$ is the
total rest energy calculated as [41]

$E_{\text{rest}}=\int^{R}_{0}mn(r)c^{2}\left(1-\frac{2{GM}(r)}{rc^{2}}\right)^{-1/2}\text{d}^{3}r.$
(3.7)

The fractional binding energy is $\varepsilon=B/E_{\text{rest}}$. It is easy
to see that the internal energy of the system $E-E_{\text{rest}}=-B$ can be
written as the sum of kinetic and potential energies,

$E-E_{\rm
rest}=\int^{R}_{0}\left[1-\frac{2{GM}(r)}{rc^{2}}\right]^{-1/2}u(r)\text{d}^{3}r+\int^{R}_{0}\left[1-\left(1-\frac{2{GM}(r)}{rc^{2}}\right)^{-1/2}\right]\rho(r)c^{2}\text{d}^{3}r,$
(3.8)

respectively, where $u=(\rho-mn)c^{2}$.

We will find the parameter regions that satisfy ${\partial S}/{\partial
w(0)}|_{b}=0$ for $S=E,~{}B$ and $\varepsilon$, respectively. The
corresponding turning points separate the stable and unstable branches of the
sequence, and delineate various extents of instabilities. We again convert
$E$, $E_{\text{rest}}$, and $B$ into dimensionless quantities
$\hat{E}=(G/c^{4}\zeta)E=\hat{M}$,
$\hat{E}_{\text{rest}}=(G/c^{4}\zeta)E_{\text{rest}}$ and
$\hat{B}=(G/c^{4}\zeta)B$. We will also use the interior redshift

$Z(r)=e^{-\Phi(r)}-1=\left(1+\frac{\epsilon_{c}(r)}{mc^{2}}\right)\left(1-\frac{2GM(R)}{c^{2}R}\right)^{-1/2}-1$
(3.9)

to indicate the relativistic extent of the system. Either high $\epsilon_{c}$
or high compactness $GM(R)/c^{2}R$ will lead to high interior redshift, though
they are interrelated.

###  3.3 Numerical results

Figure 3: Dynamical variables vs. central 3D velocity dispersion $v(0)$ and
redshift $Z(0)$ of a gravothermal system near the onset of general
relativistic instability. From top left to bottom right panels, the blue
curves denote the pressured-averaged adiabatic index $\left<\gamma\right>$,
fractional binding energy $\varepsilon$, binding energy $\hat{B}$, and total
energy $\hat{E}$, for stable (solid) and unstable (dashed) configurations,
respectively. The vertical line indicates where the instability condition is
reached (dotted). In the top left panel, the orange curve denotes the critical
adiabatic index $\gamma_{\text{cr}}$ (solid), horizontal lines denote
$\langle\gamma\rangle=5/3$ in the Newtonian limit and $4/3$ in the
ultrarelativistic limit (dotted). In the other panels, the horizontal line
indicates the maximal value of the corresponding dynamical variable (dotted).
The boundary temperature is fixed to be $b=k_{B}T(R)/mc^{2}=0.1$.

We use the fourth-order Runge-Kutta algorithm [48] to solve the Tolman-
Oppenheimer-Volkoff equation (2.5), together with equation (2.4), assuming the
two input parameters $b=k_{B}T(R)/mc^{2}$ and
$w(0)=\epsilon_{c}(0)/k_{B}T(0)$. The algorithm is robust and well-suited to
solve implicit differential equations, especially when they are stiff as in
our case. Given the results shown in Figure 2, we choose
$b=\\{0.1,0.2,0.3,0.5\\}$. For each fixed $b$ value, we scan over the central
energy cutoff $w(0)$ and find corresponding equilibrium configurations. We
then evaluate their thermal quantities and examine their stability conditions.
We collect our numerical results in Table 0(d), Appendix D, and highlight the
main findings in what follows for $b=0.1$.

Figure 3 (top left) shows pressure-averaged adiabatic index
$\left<\gamma\right>$ (blue) and critical index $\gamma_{\rm cr}$ (orange) vs.
3D velocity dispersion $v(0)$, and gravitational redshift $Z(0)$. As $v(0)$
increases, $\left<\gamma\right>$ gradually decreases from its value in the
nonrelativistic limit $5/3$, while the $\gamma_{\rm cr}$ increases from $\sim
4/3$ due to corrections in general relativity. It reaches the critical value
$1.62$ when $v(0)=0.566c$ [9], which corresponds to $Z(0)=0.750$. The results
are largely insensitive to a specific value of $b$. For $b=0.1\textup{--}0.5$
we consider, according to the adiabatic index criterion, dynamical instability
occurs when $\left<\gamma\right>$ reaches $1.62$ at
$v(0)=(0.566\textup{--}0.564)c$, corresponding to
$Z(0)=0.750\textup{--}0.622$; see Table 0(d).

The “insensitivity” reflects the degeneracy between the boundary temperature
$b$ and the central cutoff energy $w(0)$ in determining the equation of state
as indicated in Figure 2. In our numerical study, for given $b$, we scan
$w(0)$ to find the configuration that satisfies the instability condition. If
a system has a low boundary temperature, the instability can be triggered only
if the central cutoff is high enough so that the evaporation effect is
suppressed and more particles are retained on the high-energy tail. For
example, to satisfy the adiabatic index criterion, $w(0)=0.662$ for $b=0.5$,
while $w(0)=4.05$ for $b=0.1$, as shown in Table 0(d).

Figure 3 (top right) shows the fractional binding energy
$\varepsilon=(\hat{E}_{\rm rest}-\hat{E})/\hat{E}_{\rm rest}$ (blue) vs.
$v(0)$ and $Z(0)$. As $v(0)$ increases, $\varepsilon$ first increases and
reaches its maximum $\varepsilon_{\rm max}=0.0352$ at $v(0)=0.506c$,
corresponding to $Z(0)=0.548$, then decreases. From the turning-point method,
$\varepsilon_{\rm max}$ separates the equilibrium configuration into two
branches, i.e., stable (solid) and unstable (dashed). The pattern is
universal, i.e., $\varepsilon_{\rm max}=0.0352\textup{--}0.0356$ at
$Z(0)=0.548\textup{--}0.522$ is the turning point for $b=0.1\textup{--}0.5$.
We also find that $\langle\gamma\rangle=1.63$ at the turning point of the
fractional binding energy, which is slightly higher than $1.62$ from the
adiabatic index criterion.

Earlier studies [3, 4, 5] suggest that a system becomes dynamically unstable
when its fractional binding energy reaches maximum. Fully relativistic N-body
simulations [7] show that the system can collapse to a black hole when
$\varepsilon\approx 0.035$ at $Z(0)\approx 0.5$, in excellent agreement with
what we find based on the semi-analytical method. In several unstable cases
found in [7], the oscillation frequency of radial linear perturbations is
still positive, i.e., the adiabatic index condition is not satisfied. Thus the
criterion $\gamma_{\rm cr}>\langle\gamma\rangle$ is a sufficient, but may not
be necessary condition for the dynamical stability.

Figure 3 (bottom left) shows binding energy $\hat{B}=(\hat{E}_{\rm
rest}-\hat{E})$ vs. $v(0)$ and $Z(0)$. $\hat{B}$ reaches its maximum
$\hat{B}_{\rm max}=8.49\times 10^{-4}$ at $v(0)=0.395c$, corresponding to
$Z(0)=0.297$, then decreases. $\hat{B}_{\rm max}$ separates the configuration
into stable (solid) and unstable (dashed) branches. Similarly, Figure 3
(bottom right) shows the total energy $\hat{E}$ vs. $v(0)$ and $Z(0)$. The
maximum value of the total energy is $\hat{E}_{\rm max}=3.47\times 10^{-2}$ at
$v(0)=0.225c$ and $Z(0)=0.087$. According to the turning-point method, the
system becomes unstable when $\hat{E}_{\rm max}$ or $\hat{B}_{\rm max}$ is
reached. However, it is unlikely that the system could collapse into a black
hole at this stage. Instead, it would further evolve until the instability
condition based on fractional binding energy or adiabatic index is met. Figure
1 summarizes our numerical results schematically and illustrates the sequence
of dynamical instability conditions when a self-interacting dark matter halo
collapses to a black hole.

To see whether the four collapsing stages denoted in Figure 1 occur
chronologically, we need to trace the time evolution of a collapsing system.
In Sec. 5, we will estimate the dynamical timescale for collapsing into a seed
black hole for the configurations satisfying the instability conditions shown
in Figure 3. It turns out that the timescale associated with the total energy
criterion is a factor of $\sim 3$ longer than the other three ones, which are
comparable. In this work, we search for quasi-equilibrium, static solutions to
the Tolman-Oppenheimer-Volkoff equation. It is interesting to see if the
system deviates from a quasi-equilibrium state after passing the stage of the
total energy criterion, and we will leave it for future work.

In Figure 4, we show radial profiles for normalized cutoff energy $bw$ (top
left), cutoff energy $\epsilon_{c}/mc^{2}$ (top right), density $\hat{\rho}$
(middle left), 3D velocity dispersion $v/c$ (middle right), temperature
$k_{B}T/mc^{2}$ (bottom left), and adiabatic index $\gamma$ (bottom right),
for marginally stable configurations with criteria based on the adiabatic
index (dash-dotted orange), fractional binding energy (dotted purple), binding
energy (dashed blue) and total energy (solid magenta). We fix the boundary
temperature to be $b=k_{B}T(R)/mc^{2}=0.1$ and adjust the central cutoff
function $w(0)$ to find the corresponding marginal configurations. It is clear
that $\epsilon_{c}/mc^{2}$, so as for $bw(0)$, becomes higher for a stronger
instability condition, which is expected. For a given configuration, the
cutoff energy drops significantly towards outer regions $\hat{r}\rightarrow
1$. The $v$ and $\hat{\rho}$ profiles follow a similar behavior. The
temperature becomes higher towards the center due to the gravitational
redshift effect. The adiabatic index decreases towards inner regions as the
pressure increases and the equation of state becomes softer accordingly. And
it is much softer for a stronger instability criterion, in particular, the
configuration with $\gamma(0)\approx 1.59$ has
$\langle\gamma\rangle=\gamma_{\rm cr}\approx 1.62$.

We have also performed a finer scan of $b$ values for marginally stable
configurations under the adiabatic index criterion, and the results are
summarized in Table D, Appendix D. For $b=0.09\textup{--}5.0$, the
configurations that satisfy $\langle\gamma\rangle\approx\gamma_{\rm cr}=1.62$
have a central velocity dispersion of $v(0)=(0.588\textup{--}0.566)c$,
compactness of $C\simeq 0.0236\textup{--}0.0793$, and the central redshift of
$Z(0)\simeq 0.889\textup{--}0.613$. Thus the instability condition exhibits a
universal pattern. For $b\gtrsim 1$, pair production of dark matter particles
could be relevant, and we will leave it for future work.

Figure 4: Radial profiles of $bw=\epsilon_{c}(r)/[mc^{2}+\epsilon_{c}(r)]$,
cutoff energy $\epsilon_{c}/mc^{2}$, density $\hat{\rho}$, 3D velocity
dispersion $v/c$, temperature $k_{B}T/mc^{2}$, and adiabatic index $\gamma$,
for marginally stable configurations with the criteria based on the adiabatic
index (dash-dotted orange), fractional binding energy (dotted purple), binding
energy (dashed blue) and total energy (solid magenta). We fix
$b=k_{B}T(R)/mc^{2}=0.1$.

##  4 Constraining dark matter models

We have demonstrated the conditions of dynamical instability for a
gravothermal system. Our study assumes a classical truncated Maxwell-Boltzmann
distribution, which neglects quantum statistics. To examine the validity of
this assumption, it is useful to calculate the de Broglie thermal wavelength
in relativistic thermodynamics $\lambda_{\rm
dB}=\hbar/\absolutevalue{\mathbf{p}}$. Setting
$\epsilon=[\sqrt{1+({\mathbf{p}c}/{mc^{2}})^{2}}-1]mc^{2}=3k_{B}T/2$, we have

$\lambda_{\text{dB}}=\lambda_{\text{C}}\left[\left(1+\frac{3k_{B}T}{2mc^{2}}\right)^{2}-1\right]^{-1/2}.$
(4.1)

In the nonrelativistic limit $k_{B}T/mc^{2}\ll 1$, this reduces to the
familiar expression
$\lambda_{\text{dB}}=\lambda_{\text{C}}\sqrt{mc^{2}/3k_{B}T}$, where
$\lambda_{\rm C}$ is the Compton wavelength, and typically
$\lambda_{\text{dB}}\gg\lambda_{\text{C}}$. They become compatible,
$\lambda_{\text{dB}}\sim\lambda_{\text{C}}$, when the temperature is
comparable to the particle rest mass, i.e., in the relativistic regime. When
the ultrahigh density core collapses, we demand the thermal de Broglie
wavelength much smaller than the average separation distance, i.e.,
$n\lambda_{\rm dB}^{3}\ll 1$, where $n$ is the number density of dark matter
particles. In the relativistic regime
$\lambda_{\text{dB}}\simeq\lambda_{\text{C}}$, so we have $n\lambda_{\rm
dB}^{3}\sim(M_{\rm seed}/mR^{3})\lambda_{\rm C}^{3}\ll 1$, where $M_{\rm
seed}$ is the mass of seed black holes. Using the compactness relation
$C=GM_{\rm seed}/Rc^{2}$, we can write $n\lambda_{\rm dB}^{3}=C^{3}(m_{\rm
Pl}/M_{\rm seed})^{2}(m_{\rm Pl}/m)^{4}$, and hence

$n\lambda_{\rm dB}^{3}\approx
0.27\left(\frac{C}{0.04}\right)^{3}\left(\frac{10^{9}~{}{\rm
M_{\odot}}}{M_{\rm seed}}\right)^{2}\left(\frac{5~{}{\rm
keV}}{mc^{2}}\right)^{4}\ll 1.$ (4.2)

Thus for given $M_{\rm seed}$, we can derive a lower limit on the particle
mass $m$.

Consider a benchmark case discussed in [9], where the core of a $6.8\times
10^{11}~{}{\rm M_{\odot}}$ halo collapses to a seed black hole with $M_{\rm
seed}=1.9\times 10^{9}~{}{\rm M_{\odot}}$. Such a seed could further grow into
a supermassive black hole with a mass of $2.2\times 10^{9}~{}{\rm M_{\odot}}$
through accreting baryonic matter, to be consistent with observations of the
J1205-0000 quasar at redshift $6.7$ [49]. J1205-0000 has a low accretion
efficiency, and hence a massive $\sim 10^{9}~{}{\rm M_{\odot}}$ seed is needed
if one assumes an Eddington accretion history. Taking these into account, the
particle mass needs to be larger than a few ${\rm keV}$ such that the
classical truncated Maxwell-Boltzmann distribution is valid.

We further check constraints on the boundary temperature $b=k_{B}T(R)/mc^{2}$
and the central cutoff energy $w(0)=\epsilon_{c}(0)/k_{B}T(0)$. Consider a
$10^{9}~{}{\rm M_{\odot}}$ core, when the onset of dynamical instability
occurs, the boundary radius is $R=GM_{\rm seed}/Cc^{2}\approx 10^{-3}~{}{\rm
pc}$, where we take $C=0.04$. Setting the characteristic length scale
$\zeta=\lambda_{\text{C}}(m_{\rm Pl}/m)(8\pi^{3}/ge^{\alpha(R)})^{1/2}$ to be
$R/\hat{R}$ and taking $m=1~{}{\rm GeV}/c^{2}$ and $\hat{R}=0.37$, see Table
0(d) ($b=0.1$), we find $\alpha(R)\approx-44$ for fermionic dark matter,
$g=2$. The central degeneracy $\alpha(0)=\alpha(R)+w(0)$ must be much less
than $-1$ for classical distribution to be valid, i.e., $-\mu(0)/k_{B}T(0)\gg
1$; see equation (2.2). Thus there is an upper limit on $w(0)$, i.e., $w(0)\ll
43$. As shown in Figure 4 (top left), the instability condition requires
$bw(0)\gtrsim 0.3$. For $b\sim 0.1$, we have $w(0)\gtrsim 3$. In this case,
the core can collapse to a seed black hole when the system is still in the
classical regime where quantum effects are negligible. On the other hand, for
$b\sim 0.001$, the required $w(0)$ would be larger than $300$. Thus for the
core to collapse to a black hole in the classical regime, the boundary
temperature needs to be $b\gtrsim 0.1$ at the onset of the instability.

This can be used to test the collapse models based on self-interacting dark
matter. In [9], the self-interactions are purely elastic and gravothermal
collapse could occur early enough to explain the origin of supermassive black
holes at high redshifts after taking into account the effects of baryons.
Recent studies [11, 12] consider dissipative interactions, as they could speed
up the onset of gravothermal collapse [24, 25, 50]. However, it is not clear
whether dynamical instability can be triggered in this case, as the
temperature of the central core may not reach the quasirelativistic limit,
i.e., $b\sim 0.1$, due to the energy loss induced by the dissipative
interactions. It is possible that the energy release could be confined within
the collapsing core as the radiation particle has a mean free path much less
than the core radius. A detailed study is needed to further assess those
models.

In this study, we have assumed a spherical symmetry. Dark matter self-
interactions lead to a spherical shape of the inner halo [14, 51]. They also
induce viscosity that can dissipate away net angular momentum of the core
inherited from the main halo in a short timescale [9]. Thus our spherical
assumption is self-consistent and well justified. It is interesting to see
whether our work can be generalized to axisymmetric cases. Fully relativistic
numerical simulations show that a marginally unstable system can collapse into
a black hole containing $90\%$ of the total mass even if it rotates at the
mass shedding limit (Keplerian speed) [52, 53, 54]. Given this encouraging
result, we expect our overall findings in this work could be valid for
axisymmetric systems as well.

##  5 Connecting to the conducting fluid model

Another direction is to develop a formalism that bridges the relativistic
truncated Maxwell-Boltzmann model in this work, the analytical SIDM halo model
[18], as well as the nonrelativistic conducting fluid model [22]. The latter
has been widely used to study the gravothermal collapse of SIDM halos. The two
essential parameters in our model, i.e., boundary temperature
$b=k_{B}T(R)/mc^{2}$ and central energy cutoff
$w(0)=\epsilon_{c}(0)/k_{B}T(0)$, would be ultimately related to the halo
parameters and the self-interacting cross section. To fully establish such
relations, one would need simulations with a relativistic fluid model, which
is beyond the scope of this work. Here, we highlight a few useful comparisons
between our results and those in the relevant literature.

In the nonrelativistic limit $bw(r)\ll 1$, the Tolman-Oppenheimer-Volkoff
equation (2.5) for $w(r)$ becomes ${\rm d}w/{\rm d}r=-GM(r)/(br^{2}c^{2})$.
For the truncated model, one can show [35]

$\frac{{\rm d}p}{{\rm d}r}=\frac{b}{1-bw}\left(p+\rho c^{2}\right)\frac{{\rm
d}w}{{\rm d}r}=b\rho c^{2}\frac{{\rm d}w}{{\rm d}r},$ (5.1)

where the last equality assumes the nonrelativistic limit. Thus we obtain the
Newtonian hydrostatic equation ${\rm d}p/{\rm d}r=-GM(r)\rho/r^{2}$. Ref. [18]
uses the hydrostatic equation and proposes an analytical model to describe an
SIDM halo when it reaches the maximal gravothermal expansion. It further
assumes an isothermal equation of state $p=\rho v^{2}/3$ over the inner halo,
where the 3D velocity dispersion $v$ is a constant, and hence the adiabatic
index is fixed to $\gamma=5/3$. In this case, $k_{B}T/m=v^{2}/3$ is a
constant, and we have a simple relation of $b=k_{B}T/mc^{2}=(v/c)^{2}/3$.
Using the relation $\epsilon_{c}(r)=w(r)k_{B}T$, we find ${\rm
d}\epsilon_{c}/{\rm d}r=-GM(r)m/r^{2}$, thus $\epsilon_{c}$ is the escape
energy as expected.

When the SIDM halo reaches its maximal expansion, $v\approx V_{\rm
max}=1.64r_{s}\sqrt{G\rho_{s}}$ [18], where $V_{\rm max}$ is the maximal
circular velocity of the initial halo, $r_{s}$ and $\rho_{s}$ are its scale
radius and density, respectively. Since $r_{s}\propto M^{1/3}_{200}/c_{200}$
and $\rho_{s}\propto c^{3}_{200}/f_{200}$, where $M_{200}$ is the halo mass,
$c_{200}$ is the concentration and
$f_{200}=\ln(c_{200}+1)-c_{200}/(c_{200}+1)$, we obtain a scaling relation of
$b\propto M^{2/3}_{200}c_{200}/f_{200}$ in the context of the analytical SIDM
halo model [18]. For a collapsed halo, $v$ increases with the evolution time
and it becomes much larger than $V_{\rm max}$, as we will show an example
later. It is interesting to see whether the scaling relation $b\propto
M^{2/3}_{200}c_{200}/f_{200}$ still holds in the collapse phase.
Qualitatively, if the initial halo has a high mass and concentration, the
collapsed inner halo has a large energy reservoir, leading to a high boundary
temperature.

We further compare our numerical results with simulations in [9] based on the
nonrelativistic conducting fluid model. Consider the J1205-0000 benchmark case
again [9], where $\rho_{s}\approx 8.1\times 10^{7}~{}{\rm M}_{\odot}/{\rm
kpc}^{3}$ and $r_{s}\approx 10~{}{\rm kpc}$ for the halo with a mass of
$6.8\times 10^{11}~{}{\rm M_{\odot}}$ and $M_{\rm seed}=1.9\times
10^{9}~{}{\rm M_{\odot}}$ [9]. Taking $b=0.1$, the instability occurs when the
compactness reaches $C=3.9\times 10^{-2}$, see Table 0(d), and we can obtain
the boundary radius as $R=GM_{\rm seed}/Cc^{2}\approx 2.3\times 10^{-3}~{}{\rm
pc}$. The average density within $R$ is $\left<\rho_{\rm seed}\right>=3M_{\rm
seed}/(4\pi R^{3})\approx 3.5\times 10^{25}~{}{\rm M_{\odot}/kpc^{3}}$. For
comparison, the simulations based on the nonrelativistic fluid model find that
the average density of the collapsed halo is $\left<\rho_{\rm in}\right>\sim
10^{11}\rho_{s}=8.1\times 10^{18}~{}{\rm M_{\odot}/kpc^{3}}$ at the last
snapshot shown in Figure 1 of [9]. We see that $\left<\rho_{\rm
seed}\right>/\left<\rho_{\rm in}\right>\sim 4.4\times 10^{6}$, a significant
difference. Thus the simulated halo needs to evolve further for matching
$\left<\rho_{\rm seed}\right>$, which requires a relativistic version of the
conducting fluid model towards the end.

When the instability occurs, the central 3D velocity dispersion reaches the
relativistic limit, $v(0)/c=\sqrt{3p/\rho c^{2}}\approx 0.57$, but on the
boundary $v(R)/c\approx 8.8\times 10^{-3}$ ($b=0.1$), see Table D.
Interestingly, $v(R)/c$ is comparable to the velocity dispersion found in the
fluid simulations $v/c\approx
10\times\sqrt{12\pi\rho_{s}r^{2}_{s}/c^{2}}=37.4V_{\rm max}/c\approx 3.8\times
10^{-2}$, derived from the last snapshot in Figure 1 of [9] assuming the
hydrostatic condition. During the gravothermal collapse the change in the
velocity dispersion is much milder than that in the density, as the pressure
increases as well and $v\propto\sqrt{p/\rho}$. One may consider exactly
mapping the thermal quantities, such as $\rho$ and $p$, from the truncated
Maxwell-Boltzmann model, to those from the fluid model. However, this is
challenging because the adiabatic index $\gamma$ varies dynamically towards
the onset of the instability, as illustrated in Figure 3, but the original
nonrelativistic fluid model assumes a specific value of $\gamma=5/3$ [22]. To
resolve this, we could consider a relativistic fluid model and find a class of
solutions by varying $\gamma$ from $\gamma=5/3$ to $4/3$.

Lastly, we examine timescales for the ultrahigh density core to collapse into
a black hole. For the gravothermal collapse of an SIDM halo, $t_{\rm
halo}\sim{\cal O}(100)/[r_{s}\rho_{s}(\sigma/m)\sqrt{4\pi G\rho_{s}}]$ [22],
where $\sigma/m$ is the self-scattering cross section. The presence of the
baryonic potential could shorten the collapse timescale by a factor of $\sim
100$ [9] and $r_{s}\rho_{s}(\sigma/m)$ is typically $0.1\textup{--}1$ [24], we
take $t_{\rm halo}\sim 1/\sqrt{4\pi G\rho_{s}}$ after neglecting ${\cal O}(1)$
numerical factors. We estimate the timescale for collapsing into a seed as
$t_{\rm seed}\sim 1/\sqrt{4\pi G\left<\rho_{\rm seed}\right>}$. For the
J1205-0000 benchmark, $\left<\rho_{\rm seed}\right>\approx 3.5\times
10^{25}~{}{\rm M_{\odot}/kpc^{3}}$ based on the adiabatic index criterion and
$\rho_{s}\approx 8.1\times 10^{7}~{}{\rm M_{\odot}/kpc^{3}}$, and hence
$t_{\rm seed}\sim\sqrt{\rho_{s}/\left<\rho_{\rm seed}\right>}=1.5\times
10^{-9}t_{\rm halo}\sim 69~{}{\rm days}$, where we take $t_{\rm
halo}=124~{}{\rm Myr}$ for J1205-0000 [9]. Thus the timescale for collapsing
into a seed black hole is extremely short compared to that of the gravothermal
collapse of an SIDM halo.

For the J1205-0000 benchmark, we have further checked that the $t_{\rm seed}$
values associated with the adiabatic index, fractional binding energy, and
binding energy criteria are comparable, but a factor of $\sim 3$ shorter than
the one with the total energy criterion. This indicates that the total energy
criterion may not be a sufficient condition for collapsing into a black hole.
On the other hand, the system may deviate from a hydrostatic equilibrium after
the total energy criterion is satisfied, which we will leave for future work.
In addition, we have a scaling relation of $\left<\rho_{\rm
seed}\right>\propto M_{\rm seed}/R^{3}\propto 1/M^{2}_{\rm seed}$, where
$R=GM_{\rm seed}/Cc^{2}$ is used. For the SIDM model, $M_{\rm seed}\sim
10^{-3}M_{200}$ [22, 9], and we expect that the compactness $C$, a
dimensionless quantity, is largely independent of specific halo parameters.
Thus the collapse time increases with the initial halo mass as $t_{\rm
seed}\propto 1/\sqrt{\left<\rho_{\rm seed}\right>}\propto M_{\rm 200}$.

##  6 Conclusions

The origin of supermassive black holes remains unknown and the gravothermal
collapse of dark matter halos is a promising mechanism to explain the puzzle.
In this work, we have investigated a key aspect of this mechanism, i.e.,
dynamical instability of the ultrahigh density core produced at late stages of
gravothermal evolution. We used a truncated Maxwell-Boltzmann distribution to
model the dark matter distribution in the core, solved the Tolman-Oppenheimer-
Volkoff equation in a self-consistent way, and obtained a series of
equilibrium configurations. We examined four instability conditions based on
considerations of total energy, binding energy, fractional binding energy, and
adiabatic index. As the core contracts, these conditions would be satisfied in
sequential order. The adiabatic index criterion is the strongest among the
four. We have also compared our results from the semi-analytical method to
those from fully relativistic N-body simulations and found a good agreement.
In particular, both show the instability can occur when the fractional binding
energy reaches $0.035$ with a central gravitational redshift of $0.5$.

We further found that to meet the instability condition in the classical
regime, the boundary temperature of the core should be at least $10\%$ of the
mass of dark matter particles. In addition, the classical Maxwell-Boltzmann
distribution is valid only if the particle mass is larger than a few ${\rm
keV}$ for a $10^{9}~{}{\rm M_{\odot}}$ seed black hole. We have also shown
that the timescale for collapsing into a seed black hole is extremely short
compared to that of the gravothermal collapse of an SIDM halo. In the future,
we could extend our work to study dynamical instability of a self-gravitating
quantum sphere, and whether the presence of a baryonic potential would help
trigger the instability. In addition, signatures of the gravothermal collapse
could be tested using observations of satellite dwarf galaxies of the Milky
Way [55, 56, 57, 58, 59, 60, 61, 62, 63] and substructures of galaxy clusters
[64], as the interplay between self- and tidal interactions could seed up the
process. It would be interesting to explore formation of seed black holes in
those systems.

## Acknowledgments

WXF acknowledges the Institute of Physics, Academia Sinica, for the
hospitality during the completion of this work. HBY was supported by the U.S.
Department of Energy under Grant No. de-sc0008541 and the John Templeton
Foundation under Grant ID #61884. YZ was supported by the Kavli Institute for
Cosmological Physics at the University of Chicago through an endowment from
the Kavli Foundation and its founder Fred Kavli. The opinions expressed in
this publication are those of the authors and do not necessarily reflect the
views of the John Templeton Foundation.

##  Appendix A Adiabatic index for an ideal fluid

The adiabatic index of a fluid is defined as $\gamma\equiv\left({\partial\ln
p}/{\partial\ln n}\right)_{s}$ locally in spacetime. The solution for an ideal
fluid is often parametrized as $p=K(mn)^{\gamma}$, where $K$ and $\gamma$ are
not explicit functions of $n$ in an adiabatic process. We show the derivation
of the adiabatic index of an ideal fluid. For an adiabatic process, the first
law of thermodynamics tells ${\rm d}U=-p\,{\rm d}V$, where $U$ is the total
internal energy and $V$ is the volume. Suppose $N$ is the total number of
particles, $u$ is the internal energy density and $n$ is the number density,
$U=Nu/n$ and $V=N/n$. Since $N$ is a constant, we have

${\rm d}\left(\frac{u}{n}\right)=-p\,{\rm
d}\left(\frac{1}{n}\right)=Km^{\gamma}n^{\gamma-2}{\rm d}n,$ (A.1)

where the ansatz $p=K(mn)^{\gamma}$ is used for the last equality. For ideal
gas, $K$ and $\gamma$ are independent of $n$. Integrating both sides of
equation (A.1) gives
$u=Km^{\gamma}(\gamma-1)^{-1}n^{\gamma}=(\gamma-1)^{-1}p$. Since $u=(\rho-
mn)c^{2}$, there is a general relation between $\rho$ and $p$, i.e., $(\rho-
mn)c^{2}=p(\gamma-1)^{-1}$ [65, 41, 39]. Thus we have $\gamma=1+p/u$, which
can be further expressed in terms of $b$ and $w$,

$\gamma(b,w)=1+\frac{p}{u}=1+\frac{2}{3}\frac{I_{p}(b,w)}{I_{u}(b,w)},$ (A.2)

where $I_{p}$ and $I_{u}$ are given in equation (2.4). In the nonrelativistic
limit $bw\rightarrow 0$ $(I_{p}\simeq I_{u})$, $\gamma\rightarrow 5/3$; in the
ultrarelativistic limit $bw\rightarrow 1$ $(I_{p}\simeq I_{u}/2)$,
$\gamma\rightarrow 4/3$.

##  Appendix B Chandrasekhar’s instability condition

The pulsation equation (3.1) is derived by perturbing the equilibrium solution
to the Einstein equation with a Lagrangian displacement $\xi$ [1]. Here we
take a series of steps and convert it into an integral form. Multiplying its
both sides by a factor of $r^{2}e^{\Phi+\Lambda}\xi$ and integrating it over
$r$, we get ($G=c=1$)

$\displaystyle\omega^{2}\int^{R}_{0}e^{3\Lambda-\Phi}(\rho+p)r^{2}{\xi}^{2}$
$\displaystyle\text{d}r=-\int^{R}_{0}(r^{2}e^{-\Phi}\xi)\left[e^{3\Phi+\Lambda}\frac{\gamma
p}{r^{2}}(r^{2}e^{-\Phi}\xi)^{\prime}\right]^{\prime}\text{d}r+4\int^{R}_{0}e^{\Phi+\Lambda}r\frac{\text{d}p}{\text{d}r}{\xi}^{2}\text{d}r$
$\displaystyle-\int^{R}_{0}e^{\Phi+\Lambda}\left(\frac{\text{d}p}{\text{d}r}\right)^{2}\frac{r^{2}\xi^{2}}{\rho+p}\text{d}r+8\pi\int^{R}_{0}e^{3\Lambda+\Phi}p(\rho+p)r^{2}\xi^{2}\text{d}r,$
(B.1)

where “′” denotes “$\text{d}/\text{d}r$” for simplicity. Taking the first term
on the right hand side of (B), and integrating it by parts, we have

$\displaystyle-\int^{R}_{0}(r^{2}e^{-\Phi}\xi)\left[e^{3\Phi+\Lambda}\frac{\gamma
p}{r^{2}}(r^{2}e^{-\Phi}\xi)^{\prime}\right]^{\prime}\text{d}r=$
$\displaystyle\int^{R}_{0}e^{3\Phi+\Lambda}\frac{\gamma
p}{r^{2}}\left[(r^{2}e^{-\Phi}\xi)^{\prime}\right]^{2}\text{d}r$
$\displaystyle-\left.\xi e^{2\Phi+\Lambda}\gamma
p(r^{2}e^{-\Phi}\xi)^{\prime}\right|_{0}^{R},$ (B.2)

where the total derivative term vanishes after imposing the boundary condition
$\xi(0)=0$ and $p(R)=0$. Integrating the second term by parts gives rise to

$\displaystyle
4\int^{R}_{0}e^{\Phi+\Lambda}r\frac{\text{d}p}{\text{d}r}{\xi}^{2}\text{d}r$
$\displaystyle=-4\int^{R}_{0}e^{\Phi+\Lambda}[\xi^{2}+2r\xi\xi^{\prime}+r\xi^{2}(\Phi^{\prime}+\Lambda^{\prime})]p\text{d}r$
$\displaystyle=-4\int^{R}_{0}e^{\Phi+\Lambda}(\xi^{2}+2r\xi\xi^{\prime})p\text{d}r-16\pi\int^{R}_{0}e^{3\Lambda+\Phi}p(\rho+p)r^{2}\xi^{2}\text{d}r,$
(B.3)

where we have used
$2e^{-2\Lambda}(\Phi^{\prime}+\Lambda^{\prime})/r=8\pi(\rho+p)$ from the
Einstein equation.

For the third term on the right hand side of (B), we substitute
$\text{d}p/\text{d}r$ with

$\frac{\text{d}p}{\text{d}r}=-(\rho+p)\left[\frac{{M}+4\pi
pr^{3}}{r(r-2{M})}\right]=-(\rho+p)\left[\frac{1}{2r}(e^{2\Lambda}-1)+4\pi
pre^{2\Lambda}\right]$ (B.4)

and find

$\displaystyle-\int^{R}_{0}$ $\displaystyle
e^{\Phi+\Lambda}\left(\frac{\text{d}p}{\text{d}r}\right)^{2}\frac{r^{2}\xi^{2}}{\rho+p}\text{d}r=$
$\displaystyle-\int^{R}_{0}e^{\Phi+\Lambda}(\rho+p)\bigg{[}\frac{1}{4}(e^{2\Lambda}-1)^{2}+4\pi
pr^{2}(e^{2\Lambda}-1)e^{2\Lambda}+16\pi^{2}p^{2}r^{4}e^{4\Lambda}\bigg{]}\xi^{2}\text{d}r.$
(B.5)

We take $\xi(r)=re^{\Phi}$ as the trial function, which satisfies the boundary
condition $\xi(0)=0$. From the Einstein equation, we have
$2\Phi^{\prime}e^{-2\Lambda}/r-(1-e^{-2\Lambda})/r^{2}=8\pi p$, thus

$\xi^{2}+2r\xi\xi^{\prime}=r^{2}e^{2\Phi}+2r^{2}(1+r\Phi^{\prime})e^{2\Phi}=\left[3r^{2}+8\pi
pr^{4}e^{2\Lambda}+r^{2}(e^{2\Lambda}-1)\right]e^{2\Phi}.$ (B.6)

Putting all the relevant terms together, we have

$\displaystyle\omega^{2}\int^{R}_{0}e^{3\Lambda+\Phi}(\rho+p)r^{4}\text{d}r=9\int^{R}_{0}e^{3\Phi+\Lambda}\gamma
pr^{2}\text{d}r$ $\displaystyle-4\int^{R}_{0}e^{3\Phi+\Lambda}[3r^{2}+8\pi
pr^{4}e^{2\Lambda}+r^{2}(e^{2\Lambda}-1)]p\text{d}r-8\pi\int^{R}_{0}e^{3(\Phi+\Lambda)}p(\rho+p)r^{4}\text{d}r$
$\displaystyle-\int^{R}_{0}e^{3\Phi+\Lambda}\left[\frac{r^{2}}{4}(e^{2\Lambda}-1)^{2}+4\pi
pr^{4}(e^{2\Lambda}-1)e^{2\Lambda}+16\pi^{2}p^{2}r^{6}e^{4\Lambda}\right](\rho+p)\text{d}r$
$\displaystyle={}$
$\displaystyle\int^{R}_{0}e^{3\Phi+\Lambda}(9\gamma-12)pr^{2}\text{d}r-\frac{1}{4}\int^{R}_{0}e^{3\Phi+\Lambda}[16p+(e^{2\Lambda}-1)(\rho+p)](e^{2\Lambda}-1)r^{2}\text{d}r$
$\displaystyle-4\pi\int^{R}_{0}e^{3(\Phi+\Lambda)}[8p+(e^{2\Lambda}+1)(\rho+p)]pr^{4}\text{d}r-16\pi^{2}\int^{R}_{0}e^{3\Phi+5\Lambda}(\rho+p)p^{2}r^{6}\text{d}r.$
(B.7)

We determine the critical stability condition by setting the right hand side
of equation (B) to $0$ and rewrite it as $\langle\gamma\rangle-\gamma_{\rm
cr}=0$, where

$\langle\gamma\rangle\equiv\frac{\int^{R}_{0}e^{3\Phi+\Lambda}\gamma
pr^{2}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}$ (B.8)

is the pressure-averaged adiabatic index of the system, and

$\displaystyle\gamma_{{\rm cr}}$
$\displaystyle\equiv\frac{4}{3}+\frac{1}{36}\frac{\int^{R}_{0}e^{3\Phi+\Lambda}[16p+(e^{2\Lambda}-1)(\rho+p)](e^{2\Lambda}-1)r^{2}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}$
$\displaystyle+\frac{4\pi}{9}$
$\displaystyle\frac{\int^{R}_{0}e^{3(\Phi+\Lambda)}[8p+(e^{2\Lambda}+1)(\rho+p)]pr^{4}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}+\frac{16\pi^{2}}{9}\frac{\int^{R}_{0}e^{3\Phi+5\Lambda}(\rho+p)p^{2}r^{6}\text{d}r}{\int^{R}_{0}e^{3\Phi+\Lambda}pr^{2}\text{d}r}$
(B.9)

is the critical adiabatic index. A similar derivation can also be found in
[66].

##  Appendix C Critical adiabatic index in the Newtonian limit: a heuristic
derivation

For the illustration purpose, we follow [40] and show a heuristic derivation
of the instability condition in the Newtonian limit. The idea is to obtain the
pulsation equation of a Newtonian star of mass $M$ and radius $R$ with
spherical symmetry, $\delta\ddot{R}+(k/M)\delta R=0$ and determine the
effective “spring constant” $k$ of the star. A tachyonic instability will
develop if $k/M<0$.

Consider a particle on the surface, it is pulled by an inward gravitational
force $\bar{f_{g}}={GM}/{R^{2}}\approx G{\bar{\rho}}^{2}R$ and an outward
force due to pressure $\bar{f_{p}}\approx{\bar{p}}/{R}$, with the boundary
condition $p(r=R)=0$. The system is in equilibrium when
$\bar{f_{g}}-\bar{f_{p}}=0$. Let’s perturb the system radius $R\rightarrow
R+\delta R$ while keep its total mass $M$ fixed. This leads to perturbations
in the density and pressure $\delta\bar{\rho}=-3{\bar{\rho}\delta R}/{R}$ and
$\delta\bar{p}=-3{\bar{\gamma}\bar{p}}\delta R/{R}$, respectively, where have
used
$\bar{\gamma}\approx\left({\partial\ln\bar{p}}/{\partial\ln\bar{\rho}}\right)_{s}$
and $\bar{\rho}=m\bar{n}$. The resulting changes in the force are
$\delta\bar{f}_{p}=-\left(3\bar{\gamma}+1\right)\bar{f}_{p}{\delta R}/{R}$ and
$\delta\bar{f}_{g}=-5\bar{f}_{g}{\delta R}/{R}$. The acceleration related to
the net force is
$\delta\ddot{R}=({\bar{\delta}f_{p}-\bar{\delta}f_{g}})/{\bar{\rho}}=-3\left(\bar{\gamma}-4/3\right)G\bar{\rho}\delta
R.$ We can identify $3\left(\bar{\gamma}-4/3\right)G\bar{\rho}$ as $k/M$. The
spherical system will undergo an exponential growth or decay under small
radical perturbation if $\bar{\gamma}-4/3<0$.

##  Appendix D Summary of numerical results

Properties of equilibrium configurations, where we fix the boundary
temperature parameter as $b=k_{B}T(R)/mc^{2}=0.1,~{}0.2,~{}0.3,$ and $0.5$ and
scan over the central energy cutoff $w(0)=\epsilon_{c}(0)/k_{B}T(0)$ for each
$b$. From the 2nd to 14th columns, we show their total energy
$\hat{E}=\hat{M}$, total rest energy $\hat{E}_{\rm rest}$, binding energy
$\hat{B}$, fractional binding energy $\varepsilon$, system radius $\hat{R}$,
compactness $C=GM(R)/c^{2}R=\hat{M}/\hat{R}$, central interior redshift
$Z(0)$, central energy cut off $\epsilon_{c}(0)$, central energy density
$\hat{\rho}(0)$, central pressure $\hat{p}(0)$, central velocity dispersion
$v(0)$, pressure averaged adiabatic index $\left\langle{\gamma}\right\rangle$,
and critical adiabatic index $\gamma_{\rm cr}$, respectively. For each case,
we underscore marginally stable configurations following instability criteria
based on total energy, binding energy, fractional binding energy, and
adiabatic index by underscoring $w(0)$ and the corresponding critical values.
The electronic version of the data can be found at
https://github.com/michaelwxfeng/truncated-Maxwell-Boltzmann.

(a) $b=k_{B}T(R)/mc^{2}=0.1$

$w(0)$ | $\hat{E}=\hat{M}$ | $\hat{E}_{\text{rest}}$ | $\hat{B}$ | $\varepsilon$ | $\hat{R}$ | $C=\hat{M}/\hat{R}$ | $Z(0)$ | $\epsilon_{c}(0)/mc^{2}$ | $\hat{\rho}(0)$ | $\hat{p}(0)$ | $v(0)/c$ | $\langle\gamma\rangle$ | $\gamma_{{\rm cr}}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
0.1 | $2.72074\times 10^{-2}$ | $2.72729\times 10^{-2}$ | $6.54443\times 10^{-5}$ | $2.39961\times 10^{-3}$ | $6.82000\times 10^{0}$ | $3.98936\times 10^{-3}$ | $1.41544\times 10^{-2}$ | $1.01010\times 10^{-2}$ | $4.99973\times 10^{-4}$ | $1.42585\times 10^{-6}$ | $9.24965\times 10^{-2}$ | 1.66545 | 1.34070
0.3 | $3.29993\times 10^{-2}$ | $3.32285\times 10^{-2}$ | $2.29264\times 10^{-4}$ | $6.89961\times 10^{-3}$ | $2.89200\times 10^{0}$ | $1.41053\times 10^{-2}$ | $4.28931\times 10^{-2}$ | $3.09277\times 10^{-2}$ | $8.69198\times 10^{-3}$ | $7.40805\times 10^{-5}$ | $1.59902\times 10^{-1}$ | 1.66301 | 1.35554
0.6 | $\underline{3.47240\times 10^{-2}}$ | $3.51790\times 10^{-2}$ | $4.54900\times 10^{-4}$ | $1.29311\times 10^{-2}$ | $1.63800\times 10^{0}$ | $2.11990\times 10^{-2}$ | $8.71214\times 10^{-2}$ | $6.38296\times 10^{-2}$ | $5.82063\times 10^{-2}$ | $9.86165\times 10^{-4}$ | $2.25450\times 10^{-1}$ | 1.65936 | 1.37799
1.1 | $3.29600\times 10^{-2}$ | $3.36735\times 10^{-2}$ | $7.13463\times 10^{-4}$ | $2.11877\times 10^{-2}$ | $9.64001\times 10^{-1}$ | $3.41909\times 10^{-2}$ | $1.64064\times 10^{-1}$ | $1.23591\times 10^{-1}$ | $3.55962\times 10^{-1}$ | $1.09378\times 10^{-3}$ | $3.03616\times 10^{-1}$ | 1.65331 | 1.41580
1.5 | $3.02662\times 10^{-2}$ | $3.10835\times 10^{-2}$ | $8.17290\times 10^{-4}$ | $2.62934\times 10^{-2}$ | $7.25001\times 10^{-1}$ | $4.17464\times 10^{-2}$ | $2.28846\times 10^{-1}$ | $1.76464\times 10^{-1}$ | $9.92442\times 10^{-1}$ | $4.12090\times 10^{-2}$ | $3.52943\times 10^{-1}$ | 1.64851 | 1.44620
1.9 | $2.72897\times 10^{-2}$ | $2.81385\times 10^{-2}$ | $\underline{8.48824\times 10^{-4}}$ | $3.01659\times 10^{-2}$ | $5.81001\times 10^{-1}$ | $4.69701\times 10^{-2}$ | $2.96903\times 10^{-1}$ | $2.34552\times 10^{-1}$ | $2.33058\times 10^{0}$ | $1.21454\times 10^{-1}$ | $3.95398\times 10^{-1}$ | 1.64377 | 1.47651
2.3 | $2.43394\times 10^{-2}$ | $2.51671\times 10^{-2}$ | $8.27673\times 10^{-4}$ | $3.28872\times 10^{-2}$ | $4.88001\times 10^{-1}$ | $4.98757\times 10^{-2}$ | $3.68726\times 10^{-1}$ | $2.98693\times 10^{-1}$ | $4.95282\times 10^{0}$ | $3.09617\times 10^{-1}$ | $4.33059\times 10^{-1}$ | 1.63912 | 1.50646
2.7 | $2.15635\times 10^{-2}$ | $2.23348\times 10^{-2}$ | $7.71245\times 10^{-4}$ | $3.45311\times 10^{-2}$ | $4.25001\times 10^{-1}$ | $5.07375\times 10^{-2}$ | $4.44974\times 10^{-1}$ | $3.69819\times 10^{-1}$ | $9.88713\times 10^{0}$ | $7.19304\times 10^{-1}$ | $4.67177\times 10^{-1}$ | 1.63459 | 1.53571
3.2 | $1.84498\times 10^{-2}$ | $1.91224\times 10^{-2}$ | $6.72639\times 10^{-4}$ | $\underline{3.51754\times 10^{-2}}$ | $3.78001\times 10^{-1}$ | $4.88088\times 10^{-2}$ | $5.47883\times 10^{-1}$ | $4.70526\times 10^{-1}$ | $2.22039\times 10^{1}$ | $1.89599\times 10^{0}$ | $5.06132\times 10^{-1}$ | 1.62919 | 1.57054
3.4 | $1.73316\times 10^{-2}$ | $1.79606\times 10^{-2}$ | $6.28996\times 10^{-4}$ | $3.50209\times 10^{-2}$ | $3.68001\times 10^{-1}$ | $4.70966\times 10^{-2}$ | $5.91899\times 10^{-1}$ | $5.15125\times 10^{-1}$ | $3.03421\times 10^{1}$ | $2.74353\times 10^{0}$ | $5.20825\times 10^{-1}$ | 1.62714 | 1.58366
3.7 | $1.58026\times 10^{-2}$ | $1.63652\times 10^{-2}$ | $5.62545\times 10^{-4}$ | $3.43745\times 10^{-2}$ | $3.61001\times 10^{-1}$ | $4.37745\times 10^{-2}$ | $6.61622\times 10^{-1}$ | $5.87271\times 10^{-1}$ | $4.80814\times 10^{1}$ | $4.71035\times 10^{0}$ | $5.42124\times 10^{-1}$ | 1.62423 | 1.60212
4.05150 | $1.42556\times 10^{-2}$ | $1.47424\times 10^{-2}$ | $4.86859\times 10^{-4}$ | $3.30243\times 10^{-2}$ | $3.68001\times 10^{-1}$ | $3.87378\times 10^{-2}$ | $7.50107\times 10^{-1}$ | $6.81035\times 10^{-1}$ | $8.17923\times 10^{1}$ | $8.74006\times 10^{0}$ | $5.66189\times 10^{-1}$ | $\underline{1.62117}$ | 1.62117
4.2 | $1.36891\times 10^{-2}$ | $1.41457\times 10^{-2}$ | $4.56567\times 10^{-4}$ | $3.22761\times 10^{-2}$ | $3.78001\times 10^{-1}$ | $3.62145\times 10^{-2}$ | $7.90144\times 10^{-1}$ | $7.24120\times 10^{-1}$ | $1.02219\times 10^{2}$ | $1.13099\times 10^{1}$ | $5.76136\times 10^{-1}$ | 1.62002 | 1.62810
4.5 | $1.27271\times 10^{-2}$ | $1.31271\times 10^{-2}$ | $4.00052\times 10^{-4}$ | $3.04753\times 10^{-2}$ | $4.13001\times 10^{-1}$ | $3.08160\times 10^{-2}$ | $8.76941\times 10^{-1}$ | $8.18182\times 10^{-1}$ | $1.60208\times 10^{2}$ | $1.89661\times 10^{1}$ | $5.95948\times 10^{-1}$ | 1.61808 | 1.63933

(b) $b=k_{B}T(R)/mc^{2}=0.2$

$w(0)$ | $\hat{E}=\hat{M}$ | $\hat{E}_{\text{rest}}$ | $\hat{B}$ | $\varepsilon$ | $\hat{R}$ | $C=\hat{M}/\hat{R}$ | $Z(0)$ | $\epsilon_{c}(0)/mc^{2}$ | $\hat{\rho}(0)$ | $\hat{p}(0)$ | $v(0)/c$ | $\langle\gamma\rangle$ | $\gamma_{{\rm cr}}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
0.05 | $3.87808\times 10^{-2}$ | $3.88743\times 10^{-2}$ | $9.34786\times 10^{-5}$ | $2.40464\times 10^{-3}$ | $9.67200\times 10^{0}$ | $4.00960\times 10^{-4}$ | $1.41753\times 10^{-2}$ | $1.01010\times 10^{-2}$ | $2.46406\times 10^{-4}$ | $7.04972\times 10^{-7}$ | $9.26448\times 10^{-2}$ | 1.66545 | 1.34072
0.15 | $4.77866\times 10^{-2}$ | $4.81207\times 10^{-2}$ | $3.34059\times 10^{-4}$ | $6.94210\times 10^{-3}$ | $4.12300\times 10^{0}$ | $1.15903\times 10^{-2}$ | $4.30859\times 10^{-2}$ | $3.09278\times 10^{-2}$ | $4.15840\times 10^{-3}$ | $3.57900\times 10^{-5}$ | $1.60686\times 10^{-1}$ | 1.66299 | 1.35569
0.4 | $\underline{5.18443\times 10^{-2}}$ | $5.27274\times 10^{-2}$ | $8.83084\times 10^{-4}$ | $1.67481\times 10^{-2}$ | $1.84400\times 10^{0}$ | $2.81151\times 10^{-2}$ | $1.18857\times 10^{-1}$ | $8.69556\times 10^{-2}$ | $5.93008\times 10^{-2}$ | $1.37078\times 10^{-3}$ | $2.63338\times 10^{-1}$ | 1.65680 | 1.39423
0.65 | $4.97407\times 10^{-2}$ | $5.09868\times 10^{-2}$ | $1.24612\times 10^{-3}$ | $2.44400\times 10^{-2}$ | $1.20100\times 10^{0}$ | $4.14160\times 10^{-2}$ | $2.00179\times 10^{-1}$ | $1.49422\times 10^{-1}$ | $2.47787\times 10^{-1}$ | $9.37848\times 10^{-3}$ | $3.36967\times 10^{-1}$ | 1.65052 | 1.43435
0.85 | $4.67140\times 10^{-2}$ | $4.81152\times 10^{-2}$ | $1.40117\times 10^{-3}$ | $2.91212\times 10^{-2}$ | $9.36001\times 10^{-1}$ | $4.99081\times 10^{-2}$ | $2.69793\times 10^{-1}$ | $2.04808\times 10^{-1}$ | $5.80946\times 10^{-1}$ | $2.89401\times 10^{-2}$ | $3.86583\times 10^{-1}$ | 1.64544 | 1.46752
1.05 | $4.32590\times 10^{-2}$ | $4.47129\times 10^{-2}$ | $\underline{1.45394\times 10^{-3}}$ | $3.25173\times 10^{-2}$ | $7.65001\times 10^{-1}$ | $5.65476\times 10^{-2}$ | $3.44027\times 10^{-1}$ | $2.65807\times 10^{-1}$ | $1.19109\times 10^{0}$ | $7.38002\times 10^{-2}$ | $4.31138\times 10^{-1}$ | 1.64034 | 1.50157
1.2 | $4.05844\times 10^{-2}$ | $4.20226\times 10^{-2}$ | $1.43828\times 10^{-3}$ | $3.42264\times 10^{-2}$ | $6.73001\times 10^{-1}$ | $6.03036\times 10^{-2}$ | $4.03116\times 10^{-1}$ | $3.15789\times 10^{-1}$ | $1.92873\times 10^{0}$ | $1.37321\times 10^{-1}$ | $4.62161\times 10^{-1}$ | 1.63649 | 1.52760
1.35 | $3.79259\times 10^{-2}$ | $3.93103\times 10^{-2}$ | $1.38436\times 10^{-3}$ | $3.52163\times 10^{-2}$ | $6.00001\times 10^{-1}$ | $6.32097\times 10^{-2}$ | $4.65490\times 10^{-1}$ | $3.69827\times 10^{-1}$ | $3.01906\times 10^{0}$ | $2.43204\times 10^{-1}$ | $4.91598\times 10^{-1}$ | 1.63264 | 1.55399
1.5 | $3.53324\times 10^{-2}$ | $3.66322\times 10^{-2}$ | $1.29981\times 10^{-3}$ | $\underline{3.54828\times 10^{-2}}$ | $5.43001\times 10^{-1}$ | $6.50687\times 10^{-2}$ | $5.31533\times 10^{-1}$ | $4.28522\times 10^{-1}$ | $4.60834\times 10^{0}$ | $4.14967\times 10^{-1}$ | $5.19751\times 10^{-1}$ | 1.62879 | 1.58064
1.65 | $3.28360\times 10^{-2}$ | $3.40276\times 10^{-2}$ | $1.19160\times 10^{-3}$ | $3.50187\times 10^{-2}$ | $4.98001\times 10^{-1}$ | $6.59355\times 10^{-2}$ | $6.01703\times 10^{-1}$ | $4.92480\times 10^{-1}$ | $6.90365\times 10^{0}$ | $6.88165\times 10^{-1}$ | $5.46849\times 10^{-1}$ | 1.62497 | 1.60740
1.7 | $3.20294\times 10^{-2}$ | $3.31808\times 10^{-2}$ | $1.15138\times 10^{-3}$ | $3.47003\times 10^{-2}$ | $4.86001\times 10^{-1}$ | $6.59040\times 10^{-2}$ | $6.26103\times 10^{-1}$ | $5.15151\times 10^{-1}$ | $7.87346\times 10^{0}$ | $8.10392\times 10^{-1}$ | $5.55681\times 10^{-1}$ | 1.62370 | 1.61631
1.73635 | $3.14517\times 10^{-2}$ | $3.25728\times 10^{-2}$ | $1.12106\times 10^{-3}$ | $3.44172\times 10^{-2}$ | $4.77001\times 10^{-1}$ | $6.59364\times 10^{-2}$ | $6.44181\times 10^{-1}$ | $5.31996\times 10^{-1}$ | $8.65548\times 10^{0}$ | $9.11403\times 10^{-1}$ | $5.62044\times 10^{-1}$ | 1.62278 | 1.62278
1.85 | $2.96954\times 10^{-2}$ | $3.07167\times 10^{-2}$ | $1.02133\times 10^{-3}$ | $3.32499\times 10^{-2}$ | $4.53001\times 10^{-1}$ | $6.55526\times 10^{-2}$ | $7.02657\times 10^{-1}$ | $5.87236\times 10^{-1}$ | $1.15907\times 10^{1}$ | $1.30712\times 10^{0}$ | $5.81653\times 10^{-1}$ | 1.61994 | 1.64293
2.0 | $2.75004\times 10^{-2}$ | $2.83817\times 10^{-2}$ | $8.81323\times 10^{-4}$ | $3.10525\times 10^{-2}$ | $4.29001\times 10^{-1}$ | $6.41034\times 10^{-2}$ | $7.84825\times 10^{-1}$ | $6.66595\times 10^{-1}$ | $1.69159\times 10^{1}$ | $2.07718\times 10^{0}$ | $6.06946\times 10^{-1}$ | 1.61625 | 1.66912

(c) $b=k_{B}T(R)/mc^{2}=0.3$

$w(0)$ | $\hat{E}=\hat{M}$ | $\hat{E}_{\text{rest}}$ | $\hat{B}$ | $\varepsilon$ | $\hat{R}$ | $C=\hat{M}/\hat{R}$ | $Z(0)$ | $\epsilon_{c}(0)/mc^{2}$ | $\hat{\rho}(0)$ | $\hat{p}(0)$ | $v(0)/c$ | $\langle\gamma\rangle$ | $\gamma_{{\rm cr}}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
0.033 | $4.75157\times 10^{-2}$ | $4.76291\times 10^{-2}$ | $1.13485\times 10^{-4}$ | $2.38268\times 10^{-3}$ | $1.19480\times 10^{1}$ | $3.97687\times 10^{-3}$ | $1.40394\times 10^{-2}$ | $9.99899\times 10^{-3}$ | $1.59375\times 10^{-4}$ | $4.51884\times 10^{-7}$ | $9.22282\times 10^{-2}$ | 1.66546 | 1.34065
0.133 | $6.14737\times 10^{-2}$ | $6.20370\times 10^{-2}$ | $5.63289\times 10^{-4}$ | $9.07989\times 10^{-3}$ | $4.03500\times 10^{0}$ | $1.52351\times 10^{-2}$ | $5.77953\times 10^{-2}$ | $4.15580\times 10^{-2}$ | $5.77271\times 10^{-3}$ | $6.64531\times 10^{-5}$ | $1.85835\times 10^{-1}$ | 1.66176 | 1.36328
0.267 | $\underline{6.48158\times 10^{-2}}$ | $6.59264\times 10^{-2}$ | $1.11061\times 10^{-3}$ | $1.68462\times 10^{-2}$ | $2.27100\times 10^{0}$ | $2.85406\times 10^{-2}$ | $1.19486\times 10^{-1}$ | $8.70740\times 10^{-2}$ | $3.81380\times 10^{-2}$ | $8.90310\times 10^{-4}$ | $2.64638\times 10^{-1}$ | 1.65674 | 1.39468
0.467 | $6.22397\times 10^{-2}$ | $6.38962\times 10^{-2}$ | $1.65647\times 10^{-3}$ | $2.59245\times 10^{-2}$ | $1.38700\times 10^{0}$ | $4.48736\times 10^{-2}$ | $2.18873\times 10^{-1}$ | $1.62921\times 10^{-1}$ | $1.94114\times 10^{-1}$ | $8.05085\times 10^{-3}$ | $3.52738\times 10^{-1}$ | 1.64910 | 1.44387
0.633 | $5.78329\times 10^{-2}$ | $5.96955\times 10^{-2}$ | $1.86259\times 10^{-3}$ | $3.12016\times 10^{-2}$ | $1.04400\times 10^{0}$ | $5.53954\times 10^{-2}$ | $3.09024\times 10^{-1}$ | $2.34409\times 10^{-1}$ | $5.08491\times 10^{-1}$ | $2.89754\times 10^{-2}$ | $4.13460\times 10^{-1}$ | 1.64265 | 1.48675
0.733 | $5.48232\times 10^{-2}$ | $5.67169\times 10^{-2}$ | $\underline{1.89374\times 10^{-3}}$ | $3.33894\times 10^{-2}$ | $9.07001\times 10^{-1}$ | $6.04445\times 10^{-2}$ | $3.67187\times 10^{-1}$ | $2.81887\times 10^{-1}$ | $8.33940\times 10^{-1}$ | $5.54914\times 10^{-2}$ | $4.46793\times 10^{-1}$ | 1.63871 | 1.51344
0.8 | $5.27484\times 10^{-2}$ | $5.46295\times 10^{-2}$ | $1.88110\times 10^{-3}$ | $3.44338\times 10^{-2}$ | $8.33001\times 10^{-1}$ | $6.33234\times 10^{-2}$ | $4.07959\times 10^{-1}$ | $3.15788\times 10^{-1}$ | $1.13380\times 10^{0}$ | $8.28133\times 10^{-2}$ | $4.68104\times 10^{-1}$ | 1.63606 | 1.53165
0.867 | $5.06592\times 10^{-2}$ | $5.25040\times 10^{-2}$ | $1.84479\times 10^{-3}$ | $3.51361\times 10^{-2}$ | $7.70001\times 10^{-1}$ | $6.57911\times 10^{-2}$ | $4.50296\times 10^{-1}$ | $3.51529\times 10^{-1}$ | $1.51812\times 10^{0}$ | $1.20871\times 10^{-1}$ | $4.88730\times 10^{-1}$ | 1.63340 | 1.55012
0.967 | $4.75566\times 10^{-2}$ | $4.93089\times 10^{-2}$ | $1.75236\times 10^{-3}$ | $\underline{3.55384\times 10^{-2}}$ | $6.92001\times 10^{-1}$ | $6.87233\times 10^{-2}$ | $5.16649\times 10^{-1}$ | $4.08626\times 10^{-1}$ | $2.29463\times 10^{0}$ | $2.05574\times 10^{-1}$ | $5.18428\times 10^{-1}$ | 1.62942 | 1.57809
1.033 | $4.55399\times 10^{-2}$ | $4.72098\times 10^{-2}$ | $1.66994\times 10^{-3}$ | $3.53727\times 10^{-2}$ | $6.49001\times 10^{-1}$ | $7.01692\times 10^{-2}$ | $5.62702\times 10^{-1}$ | $4.49016\times 10^{-1}$ | $2.97846\times 10^{0}$ | $2.86744\times 10^{-1}$ | $5.37417\times 10^{-1}$ | 1.62678 | 1.59677
1.067 | $4.45149\times 10^{-2}$ | $4.61367\times 10^{-2}$ | $1.62176\times 10^{-3}$ | $3.51513\times 10^{-2}$ | $6.30001\times 10^{-1}$ | $7.06585\times 10^{-2}$ | $5.87183\times 10^{-1}$ | $4.70792\times 10^{-1}$ | $3.39660\times 10^{0}$ | $3.38804\times 10^{-1}$ | $5.47032\times 10^{-1}$ | 1.62543 | 1.60645
1.12526 | $4.27849\times 10^{-2}$ | $4.43162\times 10^{-2}$ | $1.53130\times 10^{-3}$ | $3.45539\times 10^{-2}$ | $5.99001\times 10^{-1}$ | $7.14271\times 10^{-2}$ | $6.30406\times 10^{-1}$ | $5.09562\times 10^{-1}$ | $4.23693\times 10^{0}$ | $4.48081\times 10^{-1}$ | $5.63266\times 10^{-1}$ | 1.62311 | 1.62311
1.167 | $4.15683\times 10^{-2}$ | $4.30293\times 10^{-2}$ | $1.46104\times 10^{-3}$ | $3.39546\times 10^{-2}$ | $5.80001\times 10^{-1}$ | $7.16693\times 10^{-2}$ | $6.62419\times 10^{-1}$ | $5.38688\times 10^{-1}$ | $4.95051\times 10^{0}$ | $5.45061\times 10^{-1}$ | $5.74723\times 10^{-1}$ | 1.62145 | 1.63507
1.267 | $3.87411\times 10^{-2}$ | $4.00187\times 10^{-2}$ | $1.27768\times 10^{-3}$ | $3.19270\times 10^{-2}$ | $5.40001\times 10^{-1}$ | $7.17426\times 10^{-2}$ | $7.42982\times 10^{-1}$ | $6.13138\times 10^{-1}$ | $7.13493\times 10^{0}$ | $8.60887\times 10^{-1}$ | $6.01643\times 10^{-1}$ | 1.61750 | 1.66373

(d) $b=k_{B}T_{R}/mc^{2}=0.5$

$w(0)$ | $\hat{E}=\hat{M}$ | $\hat{E}_{\text{rest}}$ | $\hat{B}$ | $\varepsilon$ | $\hat{R}$ | $C=\hat{M}/\hat{R}$ | $Z(0)$ | $\epsilon_{c}(0)/mc^{2}$ | $\hat{\rho}(0)$ | $\hat{p}(0)$ | $v(0)/c$ | $\langle\gamma\rangle$ | $\gamma_{{\rm cr}}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
0.02 | $6.16063\times 10^{-2}$ | $6.17550\times 10^{-2}$ | $1.48683\times 10^{-4}$ | $2.40763\times 10^{-3}$ | $1.53190\times 10^{1}$ | $4.02156\times 10^{-3}$ | $1.41877\times 10^{-2}$ | $1.01010\times 10^{-2}$ | $9.77183\times 10^{-5}$ | $2.80107\times 10^{-7}$ | $9.27331\times 10^{-2}$ | 1.66544 | 1.34073
0.1 | $8.22824\times 10^{-2}$ | $8.32134\times 10^{-2}$ | $9.31015\times 10^{-4}$ | $1.11883\times 10^{-2}$ | $4.35600\times 10^{0}$ | $1.88894\times 10^{-2}$ | $7.30960\times 10^{-2}$ | $5.26315\times 10^{-2}$ | $6.19205\times 10^{-3}$ | $8.99255\times 10^{-5}$ | $2.08730\times 10^{-1}$ | 1.66049 | 1.37116
0.18 | $\underline{8.51525\times 10^{-2}}$ | $8.67687\times 10^{-2}$ | $1.61614\times 10^{-3}$ | $1.86259\times 10^{-2}$ | $2.66900\times 10^{0}$ | $3.19043\times 10^{-2}$ | $1.35727\times 10^{-1}$ | $9.89006\times 10^{-2}$ | $3.06823\times 10^{-2}$ | $8.12905\times 10^{-4}$ | $2.81927\times 10^{-1}$ | 1.65545 | 1.40298
0.28 | $8.26154\times 10^{-2}$ | $8.48251\times 10^{-2}$ | $2.20966\times 10^{-3}$ | $2.60497\times 10^{-2}$ | $1.80700\times 10^{0}$ | $4.57196\times 10^{-2}$ | $2.19883\times 10^{-1}$ | $1.62788\times 10^{-1}$ | $1.10025\times 10^{-1}$ | $4.61268\times 10^{-3}$ | $3.54643\times 10^{-1}$ | 1.64901 | 1.44473
0.36 | $7.85954\times 10^{-2}$ | $8.10666\times 10^{-2}$ | $2.47127\times 10^{-3}$ | $3.04845\times 10^{-2}$ | $1.43200\times 10^{0}$ | $5.48850\times 10^{-2}$ | $2.92500\times 10^{-1}$ | $2.19510\times 10^{-1}$ | $2.38590\times 10^{-1}$ | $1.30406\times 10^{-2}$ | $4.04932\times 10^{-1}$ | 1.64376 | 1.47970
0.44 | $7.38054\times 10^{-2}$ | $7.63684\times 10^{-2}$ | $\underline{2.56300\times 10^{-3}}$ | $3.35610\times 10^{-2}$ | $1.18200\times 10^{0}$ | $6.24411\times 10^{-2}$ | $3.70462\times 10^{-1}$ | $2.82048\times 10^{-1}$ | $4.59406\times 10^{-1}$ | $3.11249\times 10^{-2}$ | $4.50833\times 10^{-1}$ | 1.63843 | 1.51604
0.48 | $7.12704\times 10^{-2}$ | $7.38231\times 10^{-2}$ | $2.55277\times 10^{-3}$ | $3.45795\times 10^{-2}$ | $1.08600\times 10^{0}$ | $6.56264\times 10^{-2}$ | $4.11664\times 10^{-1}$ | $3.15783\times 10^{-1}$ | $6.18568\times 10^{-1}$ | $4.60437\times 10^{-2}$ | $4.72555\times 10^{-1}$ | 1.63574 | 1.53470
0.54 | $6.73963\times 10^{-2}$ | $6.98728\times 10^{-2}$ | $2.47644\times 10^{-3}$ | $3.54421\times 10^{-2}$ | $9.67001\times 10^{-1}$ | $6.96962\times 10^{-2}$ | $4.76532\times 10^{-1}$ | $3.69836\times 10^{-1}$ | $9.40223\times 10^{-1}$ | $7.95828\times 10^{-2}$ | $5.03912\times 10^{-1}$ | 1.63168 | 1.56326
0.58 | $6.48033\times 10^{-2}$ | $6.71929\times 10^{-2}$ | $2.38964\times 10^{-3}$ | $\underline{3.55639\times 10^{-2}}$ | $9.02001\times 10^{-1}$ | $7.18439\times 10^{-2}$ | $5.21991\times 10^{-1}$ | $4.08438\times 10^{-1}$ | $1.22500\times 10^{0}$ | $1.12170\times 10^{-1}$ | $5.24119\times 10^{-1}$ | 1.62896 | 1.58265
0.6 | $6.35111\times 10^{-2}$ | $6.58478\times 10^{-2}$ | $2.33667\times 10^{-3}$ | $3.54860\times 10^{-2}$ | $8.72001\times 10^{-1}$ | $7.28337\times 10^{-2}$ | $5.45434\times 10^{-1}$ | $4.28534\times 10^{-1}$ | $1.39327\times 10^{0}$ | $1.32452\times 10^{-1}$ | $5.34038\times 10^{-1}$ | 1.62759 | 1.59243
0.64 | $6.09453\times 10^{-2}$ | $6.31589\times 10^{-2}$ | $2.21356\times 10^{-3}$ | $3.50475\times 10^{-2}$ | $8.19001\times 10^{-1}$ | $7.44142\times 10^{-2}$ | $5.93835\times 10^{-1}$ | $4.70546\times 10^{-1}$ | $1.79139\times 10^{0}$ | $1.82965\times 10^{-1}$ | $5.53541\times 10^{-1}$ | 1.62486 | 1.61217
0.662445 | $5.95204\times 10^{-2}$ | $6.16558\times 10^{-2}$ | $2.13537\times 10^{-3}$ | $3.46337\times 10^{-2}$ | $7.93001\times 10^{-1}$ | $7.50572\times 10^{-2}$ | $6.21927\times 10^{-1}$ | $4.95258\times 10^{-1}$ | $2.05635\times 10^{0}$ | $2.18273\times 10^{-1}$ | $5.64303\times 10^{-1}$ | 1.62333 | 1.62333
0.72 | $5.59314\times 10^{-2}$ | $5.78406\times 10^{-2}$ | $1.90924\times 10^{-3}$ | $3.30086\times 10^{-2}$ | $7.33001\times 10^{-1}$ | $7.63046\times 10^{-2}$ | $6.97282\times 10^{-1}$ | $5.62469\times 10^{-1}$ | $2.90495\times 10^{0}$ | $3.38632\times 10^{-1}$ | $5.91365\times 10^{-1}$ | 1.61940 | 1.65214
0.76 | $5.35030\times 10^{-2}$ | $5.52369\times 10^{-2}$ | $1.73389\times 10^{-3}$ | $3.13901\times 10^{-2}$ | $6.98001\times 10^{-1}$ | $7.66518\times 10^{-2}$ | $7.52709\times 10^{-1}$ | $6.12854\times 10^{-1}$ | $3.67288\times 10^{0}$ | $4.55200\times 10^{-1}$ | $6.09759\times 10^{-1}$ | 1.61668 | 1.67225

Properties of marginally stable configurations that satisfy the adiabatic
index criterion $\langle\gamma\rangle=\gamma_{{\rm cr}}$, given different
values of the boundary temperature $b=k_{B}T(R)/mc^{2}$. The description for
each column can be found in the caption of Table 0(d) and $v(R)/c$ is the
velocity dispersion on the core boundary. The electronic version of the data
can be found at https://github.com/michaelwxfeng/truncated-Maxwell-Boltzmann.

$b$ | $w(0)$ | $\hat{E}=\hat{M}$ | $\hat{R}$ | $C=\hat{M}/\hat{R}$ | $Z(0)$ | $\epsilon_{c}(0)/mc^{2}$ | $\hat{\rho}(0)$ | $\hat{p}(0)$ | $v(0)/c$ | $v(R)/c$ | $\langle\gamma\rangle=\gamma_{{\rm cr}}$
---|---|---|---|---|---|---|---|---|---|---|---
$5.0$ | $6.48150\times 10^{-2}$ | $2.06015\times 10^{-1}$ | $2.59700\times 10^{0}$ | $7.93280\times 10^{-2}$ | $6.12906\times 10^{-1}$ | $4.79448\times 10^{-1}$ | $1.60344\times 10^{-1}$ | $1.71013\times 10^{-2}$ | $5.65650\times 10^{-1}$ | $3.36500\times 10^{-3}$ | 1.62358
$3.0$ | $1.08185\times 10^{-1}$ | $1.58573\times 10^{-1}$ | $2.00600\times 10^{0}$ | $7.90492\times 10^{-2}$ | $6.13491\times 10^{-1}$ | $4.80491\times 10^{-1}$ | $2.71859\times 10^{-1}$ | $2.89847\times 10^{-2}$ | $5.65552\times 10^{-1}$ | $5.01669\times 10^{-3}$ | 1.62356
$2.0$ | $1.62585\times 10^{-1}$ | $1.28444\times 10^{-1}$ | $1.63300\times 10^{0}$ | $7.86552\times 10^{-2}$ | $6.14252\times 10^{-1}$ | $4.81854\times 10^{-1}$ | $4.16754\times 10^{-1}$ | $4.44142\times 10^{-2}$ | $5.65434\times 10^{-1}$ | $8.39347\times 10^{-4}$ | 1.62354
$1.0$ | $3.27070\times 10^{-1}$ | $8.86277\times 10^{-2}$ | $1.14300\times 10^{0}$ | $7.75395\times 10^{-2}$ | $6.16612\times 10^{-1}$ | $4.86021\times 10^{-1}$ | $8.91191\times 10^{-1}$ | $9.48500\times 10^{-2}$ | $5.65059\times 10^{-1}$ | $5.65015\times 10^{-3}$ | 1.62347
$0.5$ | $6.62445\times 10^{-1}$ | $5.95204\times 10^{-2}$ | $7.93001\times 10^{-1}$ | $7.50572\times 10^{-2}$ | $6.21927\times 10^{-1}$ | $4.95258\times 10^{-1}$ | $2.05635\times 10^{0}$ | $2.18273\times 10^{-1}$ | $5.64303\times 10^{-1}$ | $3.66372\times 10^{-3}$ | 1.62333
$0.3$ | $1.12526\times 10^{0}$ | $4.27849\times 10^{-2}$ | $5.99001\times 10^{-1}$ | $7.14271\times 10^{-2}$ | $6.30406\times 10^{-1}$ | $5.09562\times 10^{-1}$ | $4.23693\times 10^{0}$ | $4.48081\times 10^{-1}$ | $5.63266\times 10^{-1}$ | $9.11920\times 10^{-3}$ | 1.62311
$0.2$ | $1.73635\times 10^{0}$ | $3.14517\times 10^{-2}$ | $4.77001\times 10^{-1}$ | $6.59364\times 10^{-2}$ | $6.44181\times 10^{-1}$ | $5.31996\times 10^{-1}$ | $8.65548\times 10^{0}$ | $9.11403\times 10^{-1}$ | $5.62044\times 10^{-1}$ | $7.05363\times 10^{-3}$ | 1.62278
$0.15$ | $2.39865\times 10^{0}$ | $2.41230\times 10^{-2}$ | $4.07001\times 10^{-1}$ | $5.92702\times 10^{-2}$ | $6.63548\times 10^{-1}$ | $5.61946\times 10^{-1}$ | $1.68117\times 10^{1}$ | $1.76500\times 10^{0}$ | $5.61211\times 10^{-1}$ | $9.52828\times 10^{-3}$ | 1.62239
$0.14$ | $2.60081\times 10^{0}$ | $2.24259\times 10^{-2}$ | $3.93001\times 10^{-1}$ | $5.70632\times 10^{-2}$ | $6.70637\times 10^{-1}$ | $5.72538\times 10^{-1}$ | $2.03949\times 10^{1}$ | $2.14055\times 10^{0}$ | $5.61128\times 10^{-1}$ | $1.01614\times 10^{-2}$ | 1.62225
$0.13$ | $2.84350\times 10^{0}$ | $2.06203\times 10^{-2}$ | $3.80001\times 10^{-1}$ | $5.42638\times 10^{-2}$ | $6.80073\times 10^{-1}$ | $5.86379\times 10^{-1}$ | $2.56589\times 10^{1}$ | $2.69362\times 10^{0}$ | $5.61190\times 10^{-1}$ | $8.83738\times 10^{-3}$ | 1.62209
$0.12$ | $3.14240\times 10^{0}$ | $1.86835\times 10^{-2}$ | $3.69001\times 10^{-1}$ | $5.06326\times 10^{-2}$ | $6.93337\times 10^{-1}$ | $6.05343\times 10^{-1}$ | $3.40099\times 10^{1}$ | $3.57521\times 10^{0}$ | $5.61577\times 10^{-1}$ | $5.28167\times 10^{-3}$ | 1.62188
$0.11$ | $3.52490\times 10^{0}$ | $1.65825\times 10^{-2}$ | $3.62001\times 10^{-1}$ | $4.58080\times 10^{-2}$ | $7.13594\times 10^{-1}$ | $6.33260\times 10^{-1}$ | $4.88977\times 10^{1}$ | $5.16148\times 10^{0}$ | $5.62734\times 10^{-1}$ | $6.38807\times 10^{-3}$ | 1.62159
$0.1$ | $4.05150\times 10^{0}$ | $1.42556\times 10^{-2}$ | $3.68001\times 10^{-1}$ | $3.87378\times 10^{-2}$ | $7.50107\times 10^{-1}$ | $6.81035\times 10^{-1}$ | $8.17923\times 10^{1}$ | $8.74006\times 10^{0}$ | $5.66189\times 10^{-1}$ | $8.76561\times 10^{-3}$ | 1.62117
$0.09$ | $5.08620\times 10^{0}$ | $1.13990\times 10^{-2}$ | $4.83001\times 10^{-1}$ | $2.36004\times 10^{-2}$ | $8.89250\times 10^{-1}$ | $8.44163\times 10^{-1}$ | $2.57654\times 10^{2}$ | $2.96813\times 10^{1}$ | $5.87872\times 10^{-1}$ | $5.72277\times 10^{-3}$ | 1.62023

## References

  * [1] S. Chandrasekhar, _The Dynamical Instability of Gaseous Masses Approaching the Schwarzschild Limit in General Relativity_ , _Astrophys. J._ 140 (1964) 417–433.
  * [2] Y. B. Zel’dovich and M. A. Podurets, _The Evolution of a System of Gravitationally Interacting Point Masses_ , _Sov. Astron._ 9 (1966) 742.
  * [3] J. R. Ipser and K. S. Thorne, _Relativistic, spherically symmetric star clusters. i. stability theory for radial perturbations_ , _Astrophys. J._ 154 (1968) 251.
  * [4] J. R. Ipser, _Relativistic, spherically symmetric star clusters. II. sufficient conditions for stability against radial perturbations_ , _Astrophys. J._ 156 (1969) 509.
  * [5] J. R. Ipser, _Relativistic, spherically symmetric star clusters. III. stability of compact isotropic models_ , _Astrophys. J._ 158 (1969) 17.
  * [6] S. L. Shapiro and S. A. Teukolsky, _Relativistic stellar dynamics on the computer. i - motivation and numerical method_ , _Astrophys. J._ 298 (1985) 34.
  * [7] S. L. Shapiro and S. A. Teukolsky, _Relativistic stellar dynamics on the computer - part two - physical applications_ , _Astrophys. J._ 298 (1985) 58.
  * [8] S. L. Shapiro and S. A. Teukolsky, _Relativistic stellar dynamics on the computer. IV - collapse of a star cluster to a black hole_ , _Astrophys. J._ 307 (1986) 575.
  * [9] W.-X. Feng, H.-B. Yu and Y.-M. Zhong, _Seeding Supermassive Black Holes with Self-interacting Dark Matter: A Unified Scenario with Baryons_ , _Astrophys. J. Lett._ 914 (2021) L26, [2010.15132].
  * [10] J. Pollack, D. N. Spergel and P. J. Steinhardt, _Supermassive Black Holes from Ultra-Strongly Self-Interacting Dark Matter_ , _Astrophys. J._ 804 (2015) 131, [1501.00017].
  * [11] J. Choquette, J. M. Cline and J. M. Cornell, _Early formation of supermassive black holes via dark matter self-interactions_ , _JCAP_ 07 (2019) 036, [1812.05088].
  * [12] H. Xiao, X. Shen, P. F. Hopkins and K. M. Zurek, _SMBH seeds from dissipative dark matter_ , _JCAP_ 07 (2021) 039, [2103.13407].
  * [13] D. N. Spergel and P. J. Steinhardt, _Observational evidence for selfinteracting cold dark matter_ , _Phys. Rev. Lett._ 84 (2000) 3760–3763, [astro-ph/9909386].
  * [14] R. Dave, D. N. Spergel, P. J. Steinhardt and B. D. Wandelt, _Halo properties in cosmological simulations of selfinteracting cold dark matter_ , _Astrophys. J._ 547 (2001) 574–589, [astro-ph/0006218].
  * [15] K.-J. Ahn and P. R. Shapiro, _Formation and evolution of the self-interacting dark matter halos_ , _Mon. Not. Roy. Astron. Soc._ 363 (2005) 1092–1124, [astro-ph/0412169].
  * [16] M. Rocha, A. H. G. Peter, J. S. Bullock, M. Kaplinghat, S. Garrison-Kimmel, J. Onorbe et al., _Cosmological Simulations with Self-Interacting Dark Matter I: Constant Density Cores and Substructure_ , _Mon. Not. Roy. Astron. Soc._ 430 (2013) 81–104, [1208.3025].
  * [17] M. Vogelsberger, J. Zavala and A. Loeb, _Subhaloes in self-interacting galactic dark matter haloes_ , _Mon. Not. Roy. Astron. Soc._ 423 (2012) 3740–3752, [1201.5892].
  * [18] M. Kaplinghat, S. Tulin and H.-B. Yu, _Dark Matter Halos as Particle Colliders: Unified Solution to Small-Scale Structure Puzzles from Dwarfs to Clusters_ , _Phys. Rev. Lett._ 116 (2016) 041302, [1508.03339].
  * [19] A. Robertson, R. Massey and V. Eke, _Cosmic particle colliders: simulations of self-interacting dark matter with anisotropic scattering_ , _Mon. Not. Roy. Astron. Soc._ 467 (2017) 4719–4730, [1612.03906].
  * [20] E. O. Nadler, A. Banerjee, S. Adhikari, Y.-Y. Mao and R. H. Wechsler, _Signatures of Velocity-Dependent Dark Matter Self-Interactions in Milky Way-mass Halos_ , _Astrophys. J._ 896 (2020) 112, [2001.08754].
  * [21] S. Tulin and H.-B. Yu, _Dark Matter Self-interactions and Small Scale Structure_ , _Phys. Rept._ 730 (2018) 1–57, [1705.02358].
  * [22] S. Balberg, S. L. Shapiro and S. Inagaki, _Selfinteracting dark matter halos and the gravothermal catastrophe_ , _Astrophys. J._ 568 (2002) 475–487, [astro-ph/0110561].
  * [23] J. Koda and P. R. Shapiro, _Gravothermal collapse of isolated self-interacting dark matter haloes: N-body simulation versus the fluid model_ , _Mon. Not. Roy. Astron. Soc._ 415 (2011) 1125, [1101.3097].
  * [24] R. Essig, S. D. Mcdermott, H.-B. Yu and Y.-M. Zhong, _Constraining Dissipative Dark Matter Self-Interactions_ , _Phys. Rev. Lett._ 123 (2019) 121102, [1809.01144].
  * [25] R. Huo, H.-B. Yu and Y.-M. Zhong, _The Structure of Dissipative Dark Matter Halos_ , _JCAP_ 06 (2020) 051, [1912.06757].
  * [26] S. Balberg and S. L. Shapiro, _Gravothermal collapse of selfinteracting dark matter halos and the origin of massive black holes_ , _Phys. Rev. Lett._ 88 (2002) 101301, [astro-ph/0111176].
  * [27] R. C. Tolman, _Static solutions of Einstein’s field equations for spheres of fluid_ , _Phys. Rev._ 55 (1939) 364–373.
  * [28] J. R. Oppenheimer and G. M. Volkoff, _On Massive neutron cores_ , _Phys. Rev._ 55 (1939) 374–381.
  * [29] B. K. Harrison, K. S. Thorne, M. Wakano and J. A. Wheeler, _Gravitation Theory and Gravitational Collapse_.  University of Chicago Press, 1965.
  * [30] Y. B. Zeldovich and I. D. Novikov, _Relativistic astrophysics. Vol.1: Stars and relativity_.  University of Chicago Press, 1971.
  * [31] R. Sorkin, _A criterion for the onset of instability at a turning point_ , _Astrophys. J._ 249 (1981) 254.
  * [32] R. W. Michie, _The dynamics of spherical stellar systems: III. the relative loss of stars WithDifferent mass_ , _Mon. Not. Roy. Astron. Soc._ 126 (1963) 331–341.
  * [33] I. King, _The structure of star clusters. i. an empirical density law_ , _Astron. J._ 67 (1962) 471.
  * [34] R. Ruffini and L. Stella, _On semi-degenerate equilibrium configurations of a collisionless self-gravitating Fermi gas_ , _Astron. Astrophys._ 119 (1983) 35–41.
  * [35] M. Merafina and R. Ruffini, _Systems of selfgravitating classical particles with a cutoff in their distribution function_ , _Astron. Astrophys._ 221 (1989) 4–19.
  * [36] M. Merafina and R. Ruffini, _The relativistic equilibrium configurations of systems of particles with a cutoff in their phase space distribution functions_ , _Astron. Astrophys._ 227 (1990) 415–421.
  * [37] R. C. Tolman and P. Ehrenfest, _Temperature equilibrium in a static gravitational field_ , _Phys. Rev._ 36 (1930) 1791–1798.
  * [38] O. Klein, _On the thermodynamical equilibrium of fluids in gravitational fields_ , _Rev. Mod. Phys._ 21 (1949) 531–533.
  * [39] S. L. Shapiro and S. A. Teukolsky, _Black holes, white dwarfs, and neutron stars : the physics of compact objects_.  John Wiley and Sons, 1983.
  * [40] C. W. Misner, K. S. Thorne and J. A. Wheeler, _Gravitation_.  W. H. Freeman, San Francisco, 1973.
  * [41] S. Weinberg, _Gravitation and Cosmology: Principles and Applications of the General Theory of Relativity_.  John Wiley and Sons, New York, 1972.
  * [42] J. Katz, _On the number of unstable modes of an equilibrium_ , _Mon. Not. Roy. Astron. Soc._ 183 (1978) 765–770.
  * [43] J. M. T. Thompson, _Stability predictions through a succession of folds_ , _Philos. Trans. Royal Soc. A_ 292 (1979) 1–23.
  * [44] J. L. Friedman, J. R. Ipser and R. D. Sorkin, _Turning-point method for axisymmetric stability of rotating relativistic stars_ , _Astrophys. J._ 325 (1988) 722.
  * [45] K. Takami, L. Rezzolla and S. Yoshida, _A quasi-radial stability criterion for rotating relativistic stars_ , _Mon. Not. Roy. Astron. Soc._ 416 (2011) L1–L5, [1105.3069].
  * [46] Z. Roupas, _Thermodynamical instabilities of perfect fluid spheres in General Relativity_ , _Class. Quant. Grav._ 30 (2013) 115018, [1411.0325].
  * [47] Z. Roupas and P.-H. Chavanis, _Relativistic Gravitational Phase Transitions and Instabilities of the Fermi Gas_ , _Class. Quant. Grav._ 36 (2019) 065001, [1809.07169].
  * [48] J. Butcher, _Numerical methods for ordinary differential equations in the 20th century_ , _J. Comput. Appl. Math._ 125 (2000) 1–29.
  * [49] M. Onoue, N. Kashikawa, Y. Matsuoka, N. Kato, T. Izumi, T. Nagao et al., _Subaru high-z exploration of low-luminosity quasars (SHELLQs). VI. black hole mass measurements of six quasars at 6.1 $\leq$ z $\leq$ 6.7_, _Astrophys. J._ 880 (2019) 77.
  * [50] X. Shen, P. F. Hopkins, L. Necib, F. Jiang, M. Boylan-Kolchin and A. Wetzel, _Dissipative dark matter on FIRE – I. Structural and kinematic properties of dwarf galaxies_ , _Mon. Not. Roy. Astron. Soc._ 506 (2021) 4421–4445, [2102.09580].
  * [51] A. H. G. Peter, M. Rocha, J. S. Bullock and M. Kaplinghat, _Cosmological Simulations with Self-Interacting Dark Matter II: Halo Shapes vs. Observations_ , _Mon. Not. Roy. Astron. Soc._ 430 (2013) 105, [1208.3026].
  * [52] M. Shibata and S. L. Shapiro, _Collapse of a rotating supermassive star to a supermassive black hole: Fully relativistic simulations_ , _Astrophys. J. Lett._ 572 (2002) L39–L44, [astro-ph/0205091].
  * [53] S. L. Shapiro and M. Shibata, _Collapse of a rotating supermassive star to a supermassive black hole: Analytic determination of the black hole mass and spin_ , _Astrophys. J._ 577 (2002) 904–908, [astro-ph/0209251].
  * [54] S. L. Shapiro, _Collapse of uniformly rotating stars to black holes and the formation of disks_ , _Astrophys. J._ 610 (2004) 913–919, [astro-ph/0404338].
  * [55] H. Nishikawa, K. K. Boddy and M. Kaplinghat, _Accelerated core collapse in tidally stripped self-interacting dark matter halos_ , _Phys. Rev. D_ 101 (2020) 063009, [1901.00499].
  * [56] M. Kaplinghat, M. Valli and H.-B. Yu, _Too Big To Fail in Light of Gaia_ , _Mon. Not. Roy. Astron. Soc._ 490 (2019) 231–242, [1904.04939].
  * [57] O. Sameie, H.-B. Yu, L. V. Sales, M. Vogelsberger and J. Zavala, _Self-Interacting Dark Matter Subhalos in the Milky Way’s Tides_ , _Phys. Rev. Lett._ 124 (2020) 141102, [1904.07872].
  * [58] F. Kahlhoefer, M. Kaplinghat, T. R. Slatyer and C.-L. Wu, _Diversity in density profiles of self-interacting dark matter satellite halos_ , _JCAP_ 12 (2019) 010, [1904.10539].
  * [59] C. A. Correa, _Constraining velocity-dependent self-interacting dark matter with the Milky Way’s dwarf spheroidal galaxies_ , _Mon. Not. Roy. Astron. Soc._ 503 (2021) 920–937, [2007.02958].
  * [60] H. C. Turner, M. R. Lovell, J. Zavala and M. Vogelsberger, _The onset of gravothermal core collapse in velocity-dependent self-interacting dark matter subhaloes_ , _Mon. Not. Roy. Astron. Soc._ 505 (2021) 5327–5339, [2010.02924].
  * [61] S. Y. Kim and A. H. G. Peter, _The Milky Way satellite velocity function is a sharp probe of small-scale structure problems_ , 2106.09050.
  * [62] F. Jiang, M. Kaplinghat, M. Lisanti and O. Slone, _Orbital Evolution of Satellite Galaxies in Self-Interacting Dark Matter Models_ , 2108.03243.
  * [63] Z. C. Zeng, A. H. G. Peter, X. Du, A. Benson, S. Kim, F. Jiang et al., _Core-collapse, evaporation and tidal effects: the life story of a self-interacting dark matter subhalo_ , 2110.00259.
  * [64] D. Yang and H.-B. Yu, _Self-Interacting Dark Matter and the Excess of Small-Scale Gravitational Lenses_ , 2102.02375.
  * [65] R. F. Tooper, _Adiabatic fluid spheres in general relativity._ , _Astrophys. J._ 142 (1965) 1541.
  * [66] G. Alberti, _Statistical Mechanics of Self-Gravitating Systems in General Relativity JURY_ , Ph.D. thesis, Physique Statistique des Systèmes Complexes (LPT), France, 2017.