\mathcal{E}^r_{\paramtheta^r, \boldsymbol{x}}(\rho) := U_k(\paramtheta^r_k)V(\boldsymbol{x}) \rho V^{\dagger}(\boldsymbol{x})U^\dagger_k(\paramtheta^r_k)
\end{align}

For two reuploads the density state is:

\begin{equation} \label{eqn:two_reupload_density_QNN}
\rho^2(\paramtheta, \paramalpha, \boldsymbol{x}) =
\sum_{k=1}^K \alpha_k U(\paramtheta^2_k)V_{\boldsymbol{x}}U(\paramtheta^1_k)\ketbra{\boldsymbol{x}}{\boldsymbol{x}}U^\dagger(\paramtheta^1_k)V_{\boldsymbol{x}}U^\dagger(\paramtheta^2_k)
\end{equation}

Where we drop the subscript $k$ on the sub-unitaries for compactness, and $V_{\boldsymbol{x}} := V({\boldsymbol{x}})$. Now, in the case of sub-unitaries which are decomposed into (commuting-)blocks, somewhat confusingly, $\paramtheta$ is a rank $4$ tensor, $\paramtheta := \{\theta^{r, b}_{k, j}\}$ where $r, b, k, j$ indexes the upload, block, sub-unitary and parameter respectively.

Now unfortunately it is not clear in general how to efficiently train the above models in def:density_qnns_reuploading_def even if each unitary is of a commuting-block form of Ref. [15]. This is because the arbitrary nature of the data encoding unitary, $V(\boldsymbol{x})$ prohibits a fixed commutation relation between unitaries in subsequent uploads. However, there is one potential solution to this, which is to train in a layerwise fashion. For example, one first trains the first uploaded unitaries, $U(\paramtheta_k^1)$. Then, defining $V'_{\boldsymbol{x}} := V_{\boldsymbol{x}}U(\paramtheta_k^1)V_{\boldsymbol{x}}$ as the new (fixed) `encoding' unitary, we train the next upload $U(\paramtheta_k^2)$, and so in. So a scheme would allow gradient scaling with an additional factor of $\mathcal{O}(R)$ overhead to that of the single reuploading model,  prop:gradient_scaling_dqnn. In the above, we assume that the unitaries applied are the same in successive reuploads - i.e. if the random variable we employ chooses sub-unitary, $U_k$, to apply to the initial state for the first upload, the subsequent sub-unitaries applied are also $U_k$, although with potentially different parameters, $\paramtheta_k^1 \neq \paramtheta_k^2$. In this case the distribution $\{\alpha_k\}$ represents the distribution of reuploading `sequences'. One could also account for different distributions over each reupload. For example, we apply one upload to the state and get $U(\paramtheta^1_k)\ketbra{\boldsymbol{x}}{\boldsymbol{x}}U^\dagger(\paramtheta^1_k)$ with probability $\alpha_k$. If we choose a different (but independent) distribution, $\{\beta_{k'}\}$, to select $U(\paramtheta^2_{k'})$, the resulting state would be $U(\paramtheta^2_{k'})V_{\boldsymbol{x}}U(\paramtheta^1_k)\ketbra{\boldsymbol{x}}{\boldsymbol{x}}U^\dagger(\paramtheta^1_k)V_{\boldsymbol{x}}U^\dagger(\paramtheta^2_{k'})$ with probability $\alpha_k\beta_{k'}$, and so on.

§ EXPERIMENT DETAILS

§.§ Hyperparameter optimisation

For hyperparameter optimisation, we use the hyperparameter tuning package $\texttt{optuna}$ [69].

For the experiments in the main text relating to $\orthoqnn$s and their density counterparts, we use the following choice of hyperparameters:

* $\texttt{batch\_size} \in \{32, 48, 64, 80, \dots, 256\}$

* $\texttt{learning\_rate} \in [1\times 10^{-4}, 1\times 10^{-2}]$

* $\texttt{optimiser} \in \{\texttt{Adam}, \texttt{SGD}, \texttt{RMSProp}\}$

* $\texttt{regularisation} \in \{\ell_2, \texttt{None}\}$

The (initial) learning rate is sampled uniformly on a logarithmic scale, and the regularisation (either $\ell_2$ norm regularisation or no regularisation) is added to the loss function, which in all cases is the cross entropy loss between model outputs and the true labels. The specific plot (<ref>) shows the average over the best $\frac{x}{32}$ hyperparameter runs, those which achieve at least $70\%$ accuracy during training. For each model we have <ref> i a) ($\orthoqnn$): $\frac{32}{32}$, <ref> i b) (Density $\orthoqnn$): $\frac{29}{32}$, <ref> i c) (Density $\orthoqnn$ (compressed): $\frac{31}{32}$. Increasing the criterion to $90\%$ accuracy, the relative success drops to $\orthoqnn$: $\frac{26}{32}$,  Density $\orthoqnn$: $\frac{20}{32}$, Density $\orthoqnn$ (compressed): $\frac{22}{32}$ respectively.


# Unstable Mass Transfer from a Main-Sequence Star to a Supermassive Black
Hole
and
Quasi-Periodic Eruptions

Itai Linial Racah Institute, Hebrew University of Jerusalem, Jerusalem, 91904,
Israel  Institute for Advanced Study, 1 Einstein Drive, Princeton, New Jersey,
08540, USA  Re’em Sari Racah Institute, Hebrew University of Jerusalem,
Jerusalem, 91904, Israel Itai Linial
itai.linial<EMAIL_ADDRESS>
###### Abstract

We discuss the formation and evolution of systems comprised of a low-mass
($M_{\star}\lesssim 4\,\rm M_{\odot}$) main sequence star, orbiting a
$10^{5}-10^{7}\,\rm M_{\odot}$ supermassive black hole with an orbital period
of order $\sim$hours, and a mild eccentricity ($e\approx 0.1-0.2$),
episodically shedding mass at each pericenter passage. We argue that the
resulting mass transfer is likely unstable, with Roche lobe overflow initially
driven by gravitational wave emission, but then being accelerated by the
star’s expansion in response to its mass loss, undergoing a runaway process.
We show that such systems are naturally produced by two-body gravitational
encounters within the inner parsec of a galaxy, followed by gravitational wave
circularization and inspiral from initially highly eccentric orbits. We argue
that such systems can produce recurring flares similar to the recently
identified class of X-ray transients known as Quasi-Periodic Eruptions,
observed at the centers of a few distant galaxies.

Stellar dynamics (1596) , Tidal disruption (1696) , X-ray transient sources
(1852) , Roche lobe overflow (2155) , Supermassive black holes (1663) ,
Gravitational waves (678)

##  1 Introduction

Almost all galaxies harbor a supermassive black hole (SMBH) located at their
center, surrounded by a dense nuclear star cluster, containing millions to
billions of stars packed into a volume of a few cubic parsec (e.g., Tremaine
et al., 2002; Volonteri et al., 2021).

These dense stellar environments are fertile grounds for a variety of
dynamical processes, including the tidal disruption of stars as they pass
nearby the central black hole (Rees, 1988), the disintegration of stellar
binaries and the subsequent expulsion of hypervelocity stars (Hills, 1988), as
well the formation of extreme-mass ratio inspirals (EMRIs) of compact objects
towards the SMBH through the emission of gravitational waves (GWs) (Amaro-
Seoane, 2018).

Over the past decades, observations have revealed a plethora of transient
phenomena associated with otherwise dormant galactic nuclei across the entire
electromagnetic spectrum, and more recently even in neutrino. An important
class of nuclear transients are tidal disrutption events (TDEs), which are now
being routinely observed thanks to wide-field transient surveys such as ZTF,
ASASSN, PAN-STARRS and in the near future, the Vera Rubin Observatory (for a
recent review, see Gezari, 2021). These events occur when a star passes close
to the SMBH and is torn apart by tidal forces, leading to high rates of
infalling material, generating bright panchromatic flares.

A related, but not as well-studied phenomena involves the inspiral of a star
towards the SMBH through GW emission, until mass is gradually stripped from
the star towards the black hole (Dai & Blandford, 2013; Linial & Sari, 2017).
We refer to these systems as “stellar-inspirals” or “stellar-EMRIs” (to
distinguish them from a compact object EMRI often discussed in the literature,
e.g., Amaro-Seoane 2018). In this paper, we consider dynamical channels that
bring stars to tight orbits around the SMBH, with hour-day periods. On these
tight orbits, tides overcome the star’s self gravity around pericenter,
stripping a small fraction of its mass once an orbit, powering a bright,
highly variable, electromagnetic emission. We suggest that these flares may
have already been observed as Quasi-Periodic Eruptions (QPEs) - a recently
identified class of bright soft X-ray flares associated with galactic nuclei
(Miniutti et al., 2019; Giustini et al., 2020; Arcodia et al., 2021;
Chakraborty et al., 2021).

The structure of the paper is as follows. We introduce the timescales
governing the relevant dynamical processes within the nucleus of a galaxy in
section 2. The formation channels that place stars on orbits that tidally
interact with the SMBH are then described in section 3. In section 4, the fate
of stellar inspirals and the corresponding mass transfer phase is discussed.
Section 5 is dedicated to the connection between QPEs and stellar inspirals,
tying together the observables with our model. Our summary, conclusions and
further discussion appear in section 6.

Shortly before the completion of this paper, we became aware of a recent work
by Lu & Quataert (2022), who have similarly considered a main-sequeunce star
around an SMBH, undergoing unstable mass transfer as the origin of QPEs.

##  2 Governing timescales

Consider a cluster of stars engulfing a supermassive black hole of mass
$M_{\bullet}$ in the center of a galaxy. Assuming that all stars are of equal
mass $M_{\star}$, a total of $N_{\star,\rm tot}=M_{\bullet}/M_{\star}$ are
enclosed within the sphere of influence, where the SMBH dominates the
potential. The radius of influence is related to the velocity dispersion of
the galactic bulge, $\sigma_{h}$ as

$r_{h}=\frac{GM_{\bullet}}{\sigma_{h}^{2}}\,,$ (1)

and by assuming the M-Sigma relation (e.g., Tremaine et al., 2002; Marconi &
Hunt, 2003; McConnell & Ma, 2013), one finds

$r_{h}\approx 1.0\,M_{\bullet,6}^{0.5}\,\rm pc\,,$ (2)

where $M_{\bullet,6}=M_{\bullet}/(10^{6}\,\rm M_{\odot})$.

###  2.1 Two-body relaxation

The stellar cluster settles onto an equilibrium configuration on the two-body
relaxation timescale. At every radius $r$

$\tau_{\rm
2B}(r)\approx\frac{1}{\Omega_{k}(r)N_{\star}(r)}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{2}\ln{\Lambda}^{-1}\,,$
(3)

where $\Omega_{k}(r)=\sqrt{GM_{\bullet}/r^{3}}$ is the orbital frequency of a
circular orbit of radius $r$, $N(r)$ is the number of stars enclosed within a
radius $r$ and $\ln{\Lambda}$ is the Coloumb logarithm,
$\ln{\Lambda}\approx\ln{M_{\bullet}/M_{\star}}$. The relaxation time at
$r_{h}$ is approximately

$\tau_{\rm
2B}(r_{h})\approx\frac{GM_{\bullet}^{2}}{\sigma_{h}^{3}M_{\star}}\ln{\Lambda}^{-1}\approx
10^{9}M_{\bullet,6}^{1.25}\;\rm yr\,,$ (4)

where we used again the M-Sigma relation and $M_{\star}=\rm M_{\odot}$ was
assumed.

For galaxies with $M_{\bullet}\lesssim 10^{7}\,\rm M_{\odot}$, $\tau_{\rm
2B}(r_{h})$ is shorter than the Hubble time, and their nuclei are expected to
have settled onto an equilibrium configuration, given by the single-mass, zero
particle-flux Bahcall-Wolf (BW) cusp, with $n_{\star}(r)\propto r^{-7/4}$
(Bahcall & Wolf, 1976, 1977).

A small fraction of the stars whose semi-major axis is $r$ have a pericenter
distance, $r_{p}$ much smaller than $r$. The orbital angular momentum of such
stars is smaller by a factor $\sqrt{2r_{p}/r}$ relative to a star with a
circular orbit of the same semi-major axis. Assuming a BW profile, the number
of stars scales as $N_{\star}(r)\propto r^{5/4}$, increasing sufficiently
rapidly with $r$ such that scatterings at apocenter dominate the torque
applied on stars with eccentric orbits. The angular momentum relaxation time
for highly eccentric orbits with $r_{p}\ll r$ is therefore given by

$\tau_{\rm 2B}^{J}(r,r_{p})\approx\tau_{\rm
2B}(r)\left(\frac{r_{p}}{r}\right)\,,$ (5)

where the above expression is proportional to $(r_{p}/r)$ rather than
$\sqrt{r_{p}/r}$ owing to the diffusive nature of the uncorrelated two-body
relaxation process. The angular momentum relaxation timescale for eccentric
orbits is shorter than the energy relaxation times by a factor $\tau_{\rm
2B}^{J}(r,r_{p})/\tau_{\rm 2B}^{E}(r,r_{p})\approx(r/r_{p})^{1/4}$, assuming a
BW profile (e.g., Sari & Fragione, 2019).

###  2.2 Gravitational-wave emission

A star on a circular orbit loses both angular momentum and energy to
gravitational waves on a similar timescale, given by

$\tau_{\rm
GW}=\frac{E}{|\dot{E}|}=2\frac{J}{|\dot{J}|}\approx\frac{R_{g}}{c}\left(\frac{r}{R_{g}}\right)^{4}\left(\frac{M_{\bullet}}{M\star}\right)\,,$
(6)

where $R_{g}=GM_{\bullet}/c^{2}$ is the gravitational radius of the SMBH.

The orbital dissipation rate of stars on highly eccentric orbits due to
gravitational wave emission can be approximated as

$\tau_{\rm GW}^{E}(r,r_{p})\approx\tau_{\rm
GW}(r)\left(\frac{r_{p}}{r}\right)^{7/2}\,,$ (7)

shorter than the GW angular-momentum dissipation timescale by approximately
$\tau_{\rm GW}^{E}/\tau_{\rm GW}^{J}\approx(r_{p}/r)$ (Peters, 1964).

###  2.3 Two-body vs. Gravitational wave regime

The gravitational wave evolution timescales is comparable to the two-body
timescale on orbits that satisfy

$\tau_{\rm GW}^{E}(r,r_{p})\approx\tau_{\rm 2B}^{J}(r,r_{p})\,,$ (8)

namely,

$r_{p}=R_{g}\left(\frac{N_{\star}(r)\ln{\Lambda}}{N_{\star}}\right)^{-2/5}\,,$
(9)

and for a BW profile

$r_{p}=R_{g}\left(\frac{r}{r_{h}}\right)^{-1/2}\,,$ (10)

where order unity factors were omitted. Note that for roughly circular orbits,
where $r_{p}\lesssim r$, the two timescales are comparable for $r\approx
r_{\rm GW}$, defined as

$r_{\rm GW}\equiv(R_{g}^{2}r_{h})^{1/3}\,,$ (11)

with the subscript indicating the transition to a GW dominated regime.

Put differently, a star with a semi-major axis $r\gg r_{\rm GW}$, evolves
primarily due to two-body scatterings as long as its orbit is not too
eccentric, $r_{p}\gg R_{g}(r/r_{h})^{-1/2}$. If $r_{p}\ll
R_{g}(r/r_{h})^{-1/2}$, the orbit circularizes due to GW emission on a
timescale shorter than the angular momentum relaxation time, $\tau_{\rm
2B}^{J}(r,r_{p})$.

##  3 Stellar-EMRI Production

Stars of radius $R_{\star}$ and mass $M_{\star}$ experience tidal forces of a
magnitude comparable to their self-gravity as they approach distances of order

$r_{\rm tidal}\approx
R_{\star}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{1/3}\,,$ (12)

from the SMBH. In this section we review two dynamical channels that place
stars on orbits with pericenter distance $r_{p}\approx r_{\rm tidal}$, where
mass is peeled from the star due to strong tides. Note that the pericenter
passage time on such orbits is approximately $\sqrt{r_{\rm
tidal}^{3}/GM_{\bullet}}\approx\sqrt{R_{\star}^{3}/GM_{\star}}$, i.e., similar
to the star’s dynamical timescale.

We emphasize the distinction between tidal disruption events (TDEs) and
stellar inspirals (or stellar EMRIs) in this context. TDEs typically involve
stars on nearly parabolic orbits, originating at the outskirts of the sphere
of influence, which undergo complete tidal disruption as they approach a
distance of order $r_{\rm tidal}$ from the SMBH. “Stellar EMRIs” occur when
stars approach the black hole on a tight inspiraling orbit, such that at the
onset of mass transfer, the orbit is roughly circular ($e\lesssim 0.5$), with
hours-day periods, where a small fraction of the stellar mass is stripped at
each pericenter passage. We note that “repeating TDEs”, producing flares
separated by much longer timescales (e.g., Payne et al. 2021 with $\sim
110\,\rm d$, Liu et al. 2022 with $\sim 200\rm d$ and Wevers et al. 2022b with
$\sim 10^{3}\,\rm d$) may represent an intermediate class of events, that like
stellar-EMRIs, are long-lived and are not fully destroyed at the onset of mass
transfer, but like TDEs, have very eccentric orbits. We briefly address these
systems in the discussion section.

###  3.1 Single-single scattering and GW capture

The combination of two-body scatterings and gravitational wave emission
naturally leads to the formation of stellar inspirals, that end up shedding
mass to the black hole on mildly eccentric orbits. This process is illustrated
in figure 1.

Figure 1: Schematic phase space describing dynamical processes occurring
around SMBHs in galactic nuclei - the formation of tidal disruption events
(TDEs) and stellar-inspirals. The plot represents the sphere of influence of
an SMBH with $M_{\bullet}=10^{6}\,\rm M_{\odot}$. The vertical axis is the
orbital pericenter distance (closely related to the orbital angular momentum),
and the horizontal axis is the orbital semi-major axis (related to the orbital
energy). The diagonal black solid line corresponds to circular orbits, and the
dashed horizontal red line is the tidal radius, where tidal forces exerted by
the central supermassive black hole exceed the self gravity of a sun-like
star. TDEs are generated as a result of diffusion in angular momentum of stars
orbiting the black hole at around the radius of influence, due to two-body
encounters (the right blue funnel). Stellar-inspirals are generated primarily
by stars on smaller initial radii, with $r_{h}(R_{g}/r_{\rm tidal})^{2}$.
These stars are stochastically driven to low angular momentum, until
circularization through gravitational wave emission becomes dominant (where
the blue funnel crosses the yellow line). The star’s orbit shrinks until it
begins to shed mass onto the black hole, on a tight, mildly eccentric orbit
(represented by the stars in the inset). Additionally, Hills’ mechanism
produces stars on the dashed black line, after which they continue to evolve
due to two-body scatterings and GW emission. The magenta track represents the
trajectory of the remaining star from a binary split at a radius at which GW
emission was marginally more rapid than two-body evolution.

At every radius $r$, stars on semi-circular orbits are scattered to highly
eccentric orbits through two-body scatterings at a rate

$\mathcal{R}_{\rm 2B}(r)\approx\frac{N_{\star}(r)}{\tau_{\rm
2B}(r)}\approx\frac{N_{\star}^{2}(r)}{P(r)}\left(\frac{M_{\star}}{M_{\bullet}}\right)^{2}\approx\\\
\approx\frac{\sigma_{h}^{3}}{GM_{\bullet}}\left(\frac{r}{r_{h}}\right)\approx
10^{-4}\,\left(\frac{r}{r_{h}}\right)\;M_{\bullet,6}^{-0.25}\;\rm yr^{-1}\,,$
(13)

where we assumed a BW profile and the M-Sigma relation. Stars with initial
semi-major axis, $r_{i}$ in the range $r_{c}<r_{i}<r_{h}$, where

$r_{c}\approx r_{h}(R_{g}/r_{\rm tidal})^{2}\,,$ (14)

reach $r_{p}=r_{\rm tidal}$ through two-body interactions (with $r_{i}$
essentially unchanged), before GW emission becomes important. These stars
therefore comprise the source of TDEs, with their integrated production rate
dominated at around $r_{h}$, namely

$\mathcal{R}_{\rm TDE}=\mathcal{R}_{\rm 2B}(r_{h})\approx
10^{-4}\;M_{\bullet,6}^{-0.25}\;\rm yr^{-1}\,.$ (15)

Formation of TDEs through two-body relaxation is depicted in figure 1 as the
vertical funnel illustrating the phase-space trajectory of such stars, up to
the point where $r_{p}\approx r_{\rm tidal}$.

On the other hand, stars with $r_{\rm GW}<r_{i}<r_{c}$ undergo a different
evolutionary track. These stars similarly evolve due to two-body angular
momentum relaxation, however, before they are disrupted by tides, their orbit
begins to circularize due to loss of orbital energy to gravitational wave
emission. This process is shown in figure 1 as the blue tracks that break to
the left at the intersection with the yellow separatrix. The transition to a
GW dominated evolution occurs once their pericenter first approaches
$r_{p,i}(r_{i})=R_{g}(r_{i}/r_{h})^{-1/2}$ ($\gtrsim r_{\rm tidal}$). The
orbital eccentricity at this phase is very close to 1, with
$1-e_{i}\approx(r_{i}/r_{\rm GW})^{-3/2}\ll 1$. As the orbit continues to
dissipate due to GW emission, the pericenter eventually approaches $r_{\rm
tidal}$, marking the onset of mass transfer. The residual eccentricity $e$
when mass transfer begins satisfies (e.g., Peters, 1964)

$g(e)=\frac{r_{\rm tidal}}{r_{p,i}}=\left(\frac{r_{i}}{r_{c}}\right)^{1/2}\,,$
(16)

where

$g(e)=\frac{2e^{12/19}}{1+e}\left[\frac{304+121e^{2}}{425}\right]^{870/2299}\,,$
(17)

where we approximated $e_{i}\approx 1$. Since the production rate of stars
with initial semi-major axis of order $r_{i}$ scales as $\mathcal{R}_{\rm
2B}(r_{i})\propto r_{i}$ (equation 13), the eccentricity’s probability density
is given by $p(e)=2g(e)g^{\prime}(e)$, plotted in figure 2. $p(e)$ peaks at
around $e_{\star}=0.081$ and vanishes as $e^{5/9}$ as $e\to 0$. System with
very small eccentricity at the onset of mass transfer, i.e., $e\ll 1$, are
unlikely because these would originate from stars with initially small semi-
major axis, whose scattering rate is reduced. Very eccentric systems with
$1-e\ll 1$ are also rare, because they originate from stars with $r_{p,i}$
very close to $r_{\rm tidal}$ when gravitational waves take over two-body
scatterings, requiring fine tuning.

Figure 2: The probability density of the orbital eccentricity at the onset of
mass transfer, for systems that circularize due to gravitational wave
emission, formed through single-single scatterings. The dashed black line
gives the low eccentricity asymptotic, $p(e)\propto e^{5/19}$. The
distribution peaks around $e=0.081$.

The two-body+GW formation channel produces stellar-EMRIs at a rate of

$\mathcal{R}_{\rm MS-EMRI,single}\approx\mathcal{R}_{\rm
2B}(r_{c})\approx\mathcal{R}_{\rm TDE}\left(\frac{R_{g}}{r_{\rm
tidal}}\right)^{2}\approx\\\ \approx
10^{-7}\;M_{\bullet,6}^{1.1}\,M_{\star,1}^{-0.9}\;\rm yr^{-1}\,,$ (18)

where the M-Sigma and main-sequence relations were assumed.

###  3.2 Hills’ mechanism

A fraction $f_{b}$ of the stars at the sphere of influence are paired in
binaries. A binary with semi-major axis $a_{b}$ is tidally ionized when it
approaches a distance

$r_{\rm tidal,b}\approx r_{\rm
tidal}\left(\frac{a_{b}}{R_{\star}}\right)=a_{b}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{1/3}\,,$
(19)

where we assumed that both binary components are of similar mass $M_{\star}$.
One of the stars escapes as hypervelocity star, whereas the remaining bound
star will have a post-disruption semi-major axis of

$r_{\rm Hills}\approx r_{\rm
tidal,b}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{1/3}=a_{b}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{2/3}\,,$
(20)

independently of the binary’s initial semi-major axis around the SMBH (as long
as it is much greater than $r_{\rm Hills}$). The remaining star then evolves
according to the evolution sketched out in the previous section - it will
produce a TDE if $r_{\rm Hills}>r_{c}$, or as a stellar-EMRI otherwise. This
formation channel is demonstrated in figure 1 as the magenta trajectory,
showing the aftermath of a tidally disintegrated binary. The remaining bound
star initially has an orbit located on the black dashed line, and then
continues to evolve due to GWs/two-body scattering, depending on its position
relative to the yellow $\tau_{\rm 2B}=\tau_{\rm GW}$ separatrix. The magenta
track corresponds to the case where GWs are (marginally) dominant, and the
orbit circularizes and inspirals towards the SMBH, until it intersects the
tidal radius (red dashed line) with some small residual eccentricity.

Note that the Hills’ mechanism contributes to the formation of stellar EMRIs,
as long as

$r_{\rm tidal}\left(\frac{M_{\bullet}}{M_{\star}}\right)^{1/3}<r_{c}\,,$ (21)

or

$s_{\rm Hills,GW}\equiv\left(\frac{R_{\star}}{r_{\rm
GW}}\right)\left(\frac{M_{\bullet}}{M_{\star}}\right)^{4/9}<1\,,$ (22)

namely, the remaining star following the tidal splitting of a near-contact
binary ($a_{b}\approx R_{\star}$), has $r_{\rm Hills}$ smaller than the
critical semi-major axis for which GWs dominates over two-body relaxation
(geometrically, the intersection of the dashed black line and the yellow
separatrix lies above $r_{\rm tidal}$ in figure 1). Scaling the left hand-side
of the inequality by typical values, we get

$s_{\rm Hills,GW}=0.8M_{\bullet,6}^{-0.39}M_{\star,1}^{0.35}<1\,,$ (23)

hence, the Hills’ mechanism is not expected to contribute to GW driven stellar
EMRIs for galaxies with $M_{\bullet}\lesssim 10^{5}\rm\;M_{\odot}$. The lowest
possible eccentricity at the onset of mass transfer satisfies $g(e_{\rm
Hills,min})=s_{\rm Hills,GW}$, where the post-binary-disruption eccentricity
was approximated to be 1. Since $s_{\rm Hills,GW}$ is never too small ($\ll
1$), and unless finely tuned, never too close to $1$, the minimal residual
eccentricity is mild, i.e., $e\sim\mathcal{O}(s_{\rm Hills,GW})$. In what
follows, we use $e=0.1$ as our fiducial eccentricity value. We note that our
results are quite insensitive to the specific value of $e$ we use, as long as
it is not very small or very close to 1. Residual eccentricity of $e=0.1$
corresponds to stars whose pericenter was $r_{p,i}=r_{\rm tidal}/g(e)\gtrsim
3\times r_{\rm tidal}$ when gravitational waves took over two-body
scatterings.

For such systems, tidal heating is not expected to qualitatively affect the
evolution we describe (see Press & Teukolsky 1977; Generozov et al. 2018 for
related works). We leave a more detailed analysis of the role of the
dissipation of tidal energy, during the phase of GW inspiral to a future work.

Assuming that the distribution of initial binary separations $a_{b}$ is log-
uniform, the Hills mechanism produces stellar-EMRIs at a rate

$\mathcal{R}_{\rm MS-EMRI,Hills}\approx f_{b}\times\mathcal{R}_{\rm
TDE}\approx\\\ \approx
10^{-5}\;\left(\frac{f_{b}}{0.1}\right)M_{\bullet}^{-0.25}\;\rm yr^{-1}\,,$
(24)

up to a small logarithmic correction accounting for the number of decades of
$a_{b}$ that contribute to EMRI production.

Comparing equations 18 and 24 we see that the Hill’s mechanism is typically
the dominant source of stellar-EMRI production. This result is due to the fact
that at around $r\approx r_{h}$, which is where most stars are located,
single-single scatterings usually result in a TDE before the orbit can
circularize through GWs, while sufficiently tight binaries that are scattered
at $r_{h}$ are ionized and leave a star on a tighter orbit that can eventually
produce a low-eccentricity stellar-EMRI.

We discussed the dynamical processes which form both TDEs (with $a\approx
r_{h}$, and $1-e\approx r_{\rm tidal}/r_{h}\approx 10^{-5}$) as well as
stellar-EMRIs (with $r_{p}\approx r_{\rm tidal}\lesssim a$ and $e\approx
0.1$). Stars on orbit that may produce repeating TDEs, with say $P\approx
100\,\rm d$ and $1-e\approx 10^{-2}$ (reminiscent of Payne et al. 2021) may
also form, albeit at lower rates. Referring to figure 2, the eccentricity
distribution of stellar-EMRIs at the onset of mass transfer approaches a
constant ($\approx 0.7$) as $e\to 1$, and therefore, systems with $e>0.99$
comprise roughly $0.7\%$ of the total stellar-EMRI population formed through
single-single scatterings. Unlike the mildly eccentric inspirals previously
discussed, highly eccentric systems have $r_{p}/r_{\rm tidal}$ that is very
close to 1 when GWs become important, possibly causing the star to undergo
substantial tidal heating, that may change the subsequent evolution.

##  4 The fate of stellar inspirals

As the orbit continues to shrink due to GWs, tides gradually increase relative
to the star’s self gravity, until mass transfer ensues. Since the remaining
eccentricity is non-negligible, the resulting mass transfer is highly variable
- mass is easily stripped from the star near pericenter, while near apocenter
the star remains mostly intact.

This behavior can also be described by generalizing the common Roche potential
formalism, usually defined for circular orbits. In an eccentric orbit, one can
consider the ‘instantaneous’ Roche lobe as equivalent to that of an object
orbiting the black hole on a circular orbit of the same radius. Since
corotation cannot be maintained for an eccentric orbit, the star cannot be
stationary in the local corotating frame. Additionally, the time dependent
tidal field acting on a (mildly) eccentric orbit excites stellar modes and
tidally heats the star, possibly altering the nature of the resulting mass
transfer and the stellar structure. Despite the difficulty in incorporating
these effects to the standard Roche lobe picture, they likely amount to order-
unity corrections to the following basic picture - mass transfer begins when
the star first overfills its instantaneous Roche lobe at pericenter, while the
star is well contained within its Roche lobe during the rest of the orbit, in
between subsequent pericenter passages. The Roche lobe size at pericenter is
defined as

$R_{\rm RL}=r_{p}\left(\frac{M_{\star}}{M_{\bullet}}\right)^{1/3}\,,$ (25)

up to some order-unity corrections.

###  4.1 Mass transfer stability

Perhaps the most fundamental aspect determining the star’s evolution after the
onset of mass transfer is the stability of the resulting mass transfer.

Stability depends on several factors, and in particular on the response of the
stellar radius to loss of mass $\varepsilon=d\ln{R_{\star}}/d\ln{M_{\star}}$,
and on the conservation of orbital angular momentum of mass lost from the
star. In the limit $M_{\star}/M_{\bullet}\ll 1$, The stability criterion can
be stated as (e.g., Soberman et al., 1997; Linial & Sari, 2017; King, 2022)

$\frac{5}{6}+\frac{\varepsilon}{2}-\alpha>0\,,$ (26)

where $\varepsilon=d\ln{R_{\star}}/d\ln{M_{\star}}$ and $\alpha$ is the
fraction of orbital angular momentum lost from the system due to mass loss,
i.e., angular momentum carried by the transferred mass that is ultimately not
returned to the orbit. In the above expression, we neglected corrections of
order $(M_{\star}/M_{\bullet})^{1/3}\approx 10^{-2}$.

To evaluate $\varepsilon$, we briefly repeat the arguments and derivations
used in section 5.3 of Linial & Sari (2017). The ratio between a star’s
cooling (Kelvin-Helmholtz) time and the GW evolution timescale at the onset of
mass transfer assuming a semi-circular orbits ($r\gtrsim r_{p}\approx r_{\rm
tidal}$) is given by

$\frac{\tau_{\rm GW}}{\tau_{\rm
KH}}\approx\frac{c^{5}R_{\star}^{5}L_{\star}}{G^{4}M_{\bullet}^{2/3}M_{\star}^{13/3}}\,,$
(27)

where $L\star$ is the star’s bolometric luminosity. Assuming power-law main-
sequence relations, $R_{\star}\propto M_{\star}^{\beta_{\rm R,MS}}$ and
$L_{\star}\propto M_{\star}^{\beta_{\rm L,MS}}$ with $\beta_{R}=0.8$ and
$\beta_{L}=3.5$ we have

$\frac{\tau_{\rm GW}}{\tau_{\rm KH}}\approx
0.02\;M_{\bullet,6}^{-2/3}\,\left(\frac{M_{\star}}{M_{\odot}}\right)^{3.1}\,.$
(28)

Since mass loss evolves on a timescale
$M_{\star}/\dot{M}_{\star}\approx\tau_{\rm GW}$ (if mass transfer is stable)
or shorter (for unstable mass transfer), we find that for main-sequence stars
with $M_{\star}\lesssim 4\,M_{\bullet,6}^{0.2}$, such that the condition
$\tau_{\rm GW}<\tau_{\rm KH}$ is satisfied, the star responds adiabatically to
mass loss. Taking the star’s effective adiabatic index to be $\gamma_{\rm
ad}$, the corresponding mass radius relation is given by
$\varepsilon=(2-\gamma_{\rm ad})/(4-3\gamma_{\rm ad})$ as the star loses mass
and evolves away from the main-sequence. For low mass stars, gas pressure
dominates such that $\gamma_{\rm ad}\approx 5/3$ and $\varepsilon=-1/3$,
whereas in more massive stars, radiation pressure becomes increasingly
important, and $\gamma_{\rm ad}$ approaches $4/3$, such that $\varepsilon$
becomes a large negative number. In summary, for sufficiently low mass stars,
we expect $\varepsilon\lesssim-1/3$, corresponding to the adiabatic response
of the star obtained for $\tau_{\rm GW}<\tau_{\rm KH}$.

The remaining piece of the puzzle is the value of $\alpha$ \- the fraction of
angular momentum carried by the lost mass that is ultimately transported back
to the orbit. Precise determination of $\alpha$ and the stability of the mass
transfer may require detailed numerical simulations involving the stellar
structure, relativistic orbital motion etc., but within the scope of this
paper, we address two important effects which help constrain the value of
$\alpha$.

The extreme ratio between $M_{\star}$ and $M_{\bullet}$ implies a small
potential difference between the inner and outer Lagrange points and their
corresponding equipotential surfaces (e.g., Linial & Sari, 2017). Since the
star is already overfilling its Roche lobe to some extent, stellar material
can easily escape through the L2 Lagrange point. The proximity between the L1
and L2 surfaces implies that once the star begins to shed mass towards the
black hole, a comparable amount of mass, $f_{L2}\approx 0.5$ leaks outwards
from the star, in the opposite direction. Matter lost through L2 follows an
orbit which is external to the star’s orbit, and the star will likely torque
the resulting circumbinary ring to larger radii. This mass loss component
certainly carries away orbital angular momentum from the system, contributing
roughly $\sim 0.5$ to the value of $\alpha$, and possibly even more if the
ring is torqued to large radii.

The proximity of the orbit to the black hole’s horizon plays an important role
in the determination of $\alpha$. The specific orbital angular momentum is
given by

$J=\sqrt{GM_{\bullet}(1+e)r_{\rm tidal}}\,,$ (29)

where $r_{p}\approx r_{\rm tidal}$. Taking $e\approx 0.1$, we find

$\left.\frac{J_{\rm ISCO}}{J}\right|_{\rm MS}\approx
0.5\;M_{\bullet}^{1/3}\,M_{\star,1}^{-0.4}\,,$ (30)

where $J_{\rm ISCO}=2\sqrt{3}GM_{\bullet}/c$ is the specific angular momentum
of the innermost stable circular (ISCO) orbit of a Schwartzschild SMBH, and
the main-sequeunce relation $\beta_{\rm MS,R}=0.8$ was assumed. Thus, as
material is accreted towards the black hole and approaches the ISCO, it still
retains some $\sim 50\%$ of its initial angular momentum, which is added to
the SMBH’s spin, implying that a significant fraction of the angular momentum
is not conserved in the orbit. Combining these results, we approximate
$\alpha$ as

$\alpha_{\rm L2+spin-up}\approx f_{\rm L2}+(1-f_{\rm L2})\frac{J_{\rm
ISCO}}{J}\approx 0.75\,.$ (31)

In summary, the combination of the adiabatic radial response of main-sequence
stars to mass loss, in combination with the lack of conservation of orbital
angular momentum due to black hole spin-up and L2 mass loss imply an unstable
mass transfer, i.e., the stability criterion of equation 26 is violated.

The above discussion suggests that stability is certainly not robustly
achieved in these systems. The rest of the paper deals mostly with
implications of unstable mass transfer. For completeness, we refer the reader
to Linial & Sari (2017) for treatment of the evolution of an MS-SMBH system
undergoing stable mass transfer on a circular orbit.

###  4.2 Mass transfer stability for low-mass white dwarfs

Mass transfer from white dwarfs on highly eccentric orbits towards SMBHs has
been invoked by several authors, who have considered the consequences of
stable mass-transfer driven by GW emission (King, 2020; Zhao et al., 2021;
King, 2022).

We emphasize that the arguments applied to the case of MS-SMBH, become even
more severe for WD-SMBH systems. While for main-sequence star $r_{\rm
tidal}/R_{g}\approx 10-100$, the tidal radius for WDs is just marginally
outside the ISCO for the SMBH masses considered in these studies (e.g., King,
2020; Zhao et al., 2021). This implies that almost all of the material’s
initial angular momentum is lost to spinning up the black hole, namely - the
mass transfer is highly non-conservative. Furthermore, even if these systems
are assumed to evolve stably, the implied values of $\dot{M}_{\rm WD}$ are
only achieved if the WD overfills its Roche lobe at pericenter to an extent
that would lead to substantial loss of mass (and with it, angular momentum)
through L2. These two effects imply $\alpha\lesssim 1$, and the stability
criterion of equation 26 cannot be satisfied (here again, $\varepsilon<-1/3$
as for main-sequence stars).

###  4.3 Evolution under unstable mass transfer

Unstable mass transfer implies that the system enters a runaway evolution,
where the amount of mass lost from the star increases at every pericenter
passage.

The small amount of mass initially stripped from the star has negligible
effect on the stellar structure and on the orbital evolution. The orbit
continues to shrink due to GW emission, such that the star overfills its Roche
lobe to a greater extent, gradually increasing the orbit-averaged mass loss
rate $\dot{M}_{\star}$. We define the relative Roche lobe overfilling extent
as (see equations 12 and 25)

$\xi=\frac{R_{\star}-R_{\rm RL}}{R_{\star}}=\frac{r_{\rm tidal}-r_{p}}{r_{\rm
tidal}}\,.$ (32)

Its time derivative is given by

$\dot{\xi}=\left(\frac{\dot{R}_{\star}}{R_{\star}}-\frac{\dot{R}_{\rm
RL}}{R_{\rm RL}}\right)(1+\xi)=\\\
=\left((\varepsilon+5/3-2\alpha)\frac{\dot{M}_{\star}}{M_{\star}}-\frac{\dot{r}_{p}}{r_{p}}\right)(1+\xi)\,.$
(33)

In the above expression, $\dot{M}_{\star}<0$, $(\varepsilon+5/3-2\alpha)<0$
(unstable mass transfer), as well as $\dot{r}_{p}<0$ (GWs cause the pericenter
to shrink), and therefore, $\dot{\xi}>0$. We further note that the (orbit-
averaged) mass loss rate is a sensitive function of $\xi$, with

$\dot{M}_{\star}\approx-\frac{M_{\star}}{P}\xi^{k}\,,$ (34)

where $k\approx 3.5-5$ (see derivation in appendix B). This strong dependence
is a result of the density stratification near the stellar surface. Using
equation 34, the equation for $\dot{\xi}$ is rewritten as

$\dot{\xi}\approx\\\ \left(P^{-1}\xi^{k}+\tau_{\rm GW}^{-1}\right)(1+\xi)\,,$
(35)

where we used $\dot{r}_{p}/r_{p}\equiv-\tau_{\rm GW}^{-1}$, and omitted order
unity prefactors.

Thus, when $\xi$ is very small, it increases roughly linearly with time,
$\xi\approx t/\tau_{\rm GW}$, where $t=0$ marks the onset of mass transfer.
However, as it approaches a value of $\xi_{\rm
linear/runaway}\approx(P/\tau_{\rm GW})^{1/k}$, changes in the star’s radius
and the corresponding shrinkage of the Roche lobe due to mass loss begin to
dominate over the GW term, causing $\xi$ to evolve in a runaway fashion.

At this phase, the star sheds ever increasing amounts of mass at every
pericenter passage, to the extent that GW emission, which was the initial
trigger of this process, stops playing a role in the subsequent evolution. The
star would continue to evolve in this runaway manner even if the emission of
gravitational waves was somehow artificially stopped. The system continues to
evolve on ever decreasing timescales, with

$\tau_{\rm runaway}(\xi)\approx\frac{\xi}{\dot{\xi}}\approx P\xi^{1-k}\,,$
(36)

representing the timescale spent at any value of $\xi$. As $\xi$ approaches 1,
the evolution time becomes comparable to the orbital period (and by
definition, the star’s dynamical time), at which point nearly half the stellar
mass is spilled over just a few orbits. The elapsed time from the onset of
mass transfer to the final disruption of the star is of order

$\tau_{\rm total}\approx\tau_{\rm GW}\left(\frac{P}{\tau_{\rm
GW}}\right)^{1/k}\,,$ (37)

or roughly

$\frac{\tau_{\rm total}}{P}\approx\left(\frac{\tau_{\rm
GW}}{P}\right)^{1-1/k}\approx
10^{6}\,M_{\bullet,6}^{-0.5}\,M_{\star}^{-2.3}\,,$ (38)

orbits, where we assumed $k=4$. We note that while the star does inflate and
evolve away from its initial main-sequence structure, throughout most of the
runaway evolution, when $\xi\ll 1$, the star has only lost a small fraction of
its initial mass, $\delta M_{\star}\approx\xi M_{\star}$, and therefore its
overall structure does not evolve significantly, at least until
$\xi\sim\mathcal{O}(1)$.

###  4.4 Emission resulting from transferred mass

The mass transfer from the star to the black hole can fuel bright
electromagnetic emission. For a given $\dot{M}_{\star}$, the system will
produce a time-averaged bolometric luminosity

$\left<L\right>\approx\left|\dot{M}_{\star}\right|c^{2}\left(\frac{R_{g}}{r_{\rm
diss}}\right)\,,$ (39)

where $r_{\rm ISCO}\leq r_{\rm diss}\leq r_{\rm tidal}$ is the radius at which
the orbital energy of the infalling material is dissipated and converted to
radiation (assuming a sufficiently efficient radiative process). As
$\dot{M}_{\star}$ increases with time, the averaged emission gradually
brightens. Normalized by typical values

$\left<L\right>(\xi)\approx 10^{41}\;\left(\frac{10R_{g}}{r_{\rm
diss}}\right)\left(\frac{\xi}{10^{-2}}\right)^{4}\,M_{\star,1}^{0.3}\;\rm
erg\,s^{-1}\,,$ (40)

where the main-sequence mass-radius relation $\beta_{\rm MS,R}=0.8$ has been
assumed.

##  5 Quasi-Periodic Eruptions from Unstable Mass Transfer from a Main-
Sequence Star

Quasi-periodic eruptions (QPEs) are enigmatic high amplitude bursts of X-ray
radiation with a recurrence period of a few to several hours, recently
discovered near the central supermassive black holes (SMBHs) of a few distant
galaxies (Miniutti et al., 2019; Giustini et al., 2020; Arcodia et al., 2021;
Chakraborty et al., 2021). Each QPE flare has a duration that is roughly
$5-20\%$ of the recurrence period, and emits roughly $10^{45}-10^{46}\;\rm
ergs$ in soft X-rays. Different models have been proposed in the literature to
explain this phenomena. Several authors have invoked stable mass-transfer from
a compact star or white dwarf onto a low-mass SMBH (e.g., King, 2020, 2022;
Zhao et al., 2021). Mass transfer from a main-sequence star has been recently
proposed in Krolik & Linial (2022), where shocks self-intersecting streams of
material have been considered as the source of radiation. Main-sequence stars
have also been considered by Metzger et al. (2021), who suggested that QPE
flares arise from interaction between two counter-rotating main-sequence stars
on circular orbits near the tidal radius, shedding an excess of mass at each
conjunction of the two stars. Xian et al. (2021) considered collisions between
a pre-existing AGN disk and a star on an inclined orbit to explain the timing
of the flares seen in GSN-069 (Miniutti et al., 2019).

In this section, we will argue that stellar inspirals may naturally explain
some of the fundamental observational properties of the known QPE systems.

###  5.1 Orbital period and QPE recurrence time

Stellar EMRIs orbit the SMBH on a mildly eccentric orbits with semi-major axis
$\approx r_{\rm tidal}$, and have an orbital period $P\approx\tau_{\rm dyn}$,
that depends only on the mean stellar density. Assuming a standard main-
sequence stellar structure and $e=0.1$, and by setting $r_{p}=r_{\rm tidal}$,
the stellar mass and orbital periods are related by
$M_{\star}=0.9\,M_{\odot}\,P_{10}^{1.4}$ where $P_{10}=P/10\,\rm h$ (see
appendix A).

Mass transfer from stellar inspirals varies considerably along an orbit. The
resulting emission would be strongly modulated over the orbital period $P$, as
long as the dissipation mechanism is capable of rapidly dissipating the
orbital energy of the mass transfer stream, and as long as the heated material
can readily cool to allow for large variability in the emitted flux (e.g.,
Krolik & Linial, 2022). One natural dissipation mechanism capable of producing
the high-amplitude, short duty cycle and observed temperature of the QPE
flares is due to shocks occurring when streams of mass transfer self-intersect
near the ISCO of the SMBH, as discussed in Krolik & Linial (2022). On the
other hand, the long viscous timescale of a standard accretion disk
($\tau_{\rm visc}\gg P(r_{p})$) implies that the high-amplitude, short duty
QPE-like flares, cannot originate directly from accretion onto the SMBH,
challenging the picture portrayed in Metzger et al. (2021). However, in their
scenario, mass transfer streams originating from the two counter-rotating
stars might mix and have a reduced angular momentum, possibly leading to the
formation of an accretion disk with $r<r_{\rm tidal}$, that may have a
sufficiently short viscous time. Interpreting the recurrence time of QPEs as
the orbital period of stellar inspirals, we find that the range of observed
periods $P=2.4-18.5\rm\,h$, corresponds to stellar masses in the range
$0.1-2.2\,\rm M_{\odot}$ \- which comprises the majority of stellar objects in
the galactic center. Our estimates of $M_{\star}$ for the known QPE systems
are summarized in table 1.

Table 1: Properties of the known QPE systems. $P$ is the average QPE period, and $\left<L\right>$ is the QPE luminosity averaged over a full period. The Eddington ratio is based on the SMBH mass estimates from Wevers et al. (2022a). Given these quantities, we estimate $M_{\star}$ \- the mass of the donor main-sequence star, the apparent timescale $\tau_{\rm app}$ and the luminosity evolution time $\tau_{L}$, as well as the Roche lobe overflow extent, $\xi$. Name | $P\;(\rm hr)$ | $\left<L\right>\;(\rm erg\,s^{-1})$ | $\frac{\left<L\right>}{L_{\rm Edd}}$ | $M_{\star}\;(\rm M_{\odot})$ | $\log_{10}\xi$ | $\tau_{L}\;\rm(yr)$ | $\tau_{\rm app}\;\rm(yr)$
---|---|---|---|---|---|---|---
GSN-069 (Miniutti et al., 2019) | $9$ | $5\times 10^{41}$ | $\sim 10^{-2}$ | $0.8$ | $-2.4$ | $14$ | $\sim 8\times 10^{3}$
RX J1301.9+2747 (Giustini et al., 2020) | $4.5$ | $\sim 10^{41}$ | $10^{-3}-10^{-4}$ | $0.3$ | $-2.4$ | $16$ | $\sim 10^{4}$
eRO-QPE1 (Arcodia et al., 2021) | $18.5$ | $4\times 10^{42}$ | $10^{-2}-10^{-1}$ | $2.3$ | $-1.4$ | $38$ | $\sim 3\times 10^{3}$
eRO-QPE2 (Arcodia et al., 2021) | $2.4$ | $2\times 10^{41}$ | $10^{-2}-10^{-1}$ | $0.12$ | $-2.4$ | $5$ | $\sim 3\times 10^{3}$
2MASXJ0249 (Chakraborty et al., 2021) | $2.5$ | $\sim 10^{41}$ | $10^{-2}$ | $0.13$ | $-2.6$ | $12$ | $\sim 10^{4}$

###  5.2 QPE energetics

The observed period-averaged QPE luminosity of roughly $\left<L\right>\approx
10^{41}-10^{42}\,\rm erg\,s^{-1}$ can be explained by an average mass transfer
rate $\dot{M}_{\star}=\left<L\right>(r_{\rm diss}/R_{g})/c^{2}\approx
10^{-4}\,\rm M_{\odot}\,yr^{-1}$ onto the SMBH, where $r_{\rm diss}$ was
defined as the radius at which the matter’s orbital energy is dissipated, and
we used a fiducial value of $10R_{g}$. Here we assumed that the observed X-ray
luminosity is similar to the total bolometric luminosity of the system. This
corresponds to an apparent mass-loss timescale, $\tau_{\rm
app}=M_{\star}/\dot{M}_{\star}\approx 5000\,\rm yr$ for $M_{\star}\approx
1\;\rm M_{\odot}$ (see table 1 for specific values for each system). We
interpret this implied mass loss rate as originating from a stellar EMRI. From
the observerd luminosity, we use equation 40 to conclude $\xi=(P/\tau_{\rm
app})^{1/k}\approx 10^{-2}$ (table 1).

We note that the apparent evolution timescale $\tau_{\rm app}$ is several
orders of magnitude shorter than the gravitational wave inspiral time of a
solar mass main-sequence star orbiting an SMBH at $r_{\rm tidal}$, $\tau_{\rm
GW}\approx 10^{5}\,M_{\star,1}^{0.9}M_{\bullet,6}^{-2/3}\,\rm yr$ (e.g.,
Linial & Sari, 2017). This suggests that the system cannot be undergoing a
stable mass transfer driven by GW emission, as that would imply $\tau_{\rm
app}\approx\tau_{\rm GW}$. We conclude that the system is currently undergoing
an unstable mass transfer, where the stellar response to mass loss causes the
star to expand, therefore $\xi$ to increase, further accelerating the mass
loss, in a runaway process. This determination is in agreement with the
analysis given in section 4.1, where we have demonstrated that the mass
transfer is not conservative enough to be stable.

Our inference of unstable mass transfer has strong consequences on the
estimate of the system’s remaining lifetime. Significant increase in the
accretion rate would occur once the mass of the star changes by only a small
fraction $\xi$, rather than by order unity as is in the case of stable mass
transfer. Thus, the evolution time is approximately
$\tau_{L}\approx\xi\tau_{\rm app}\approx 10\,P_{10}^{1.3}\,L_{42}^{-0.7}\,\rm
yr$. We therefore predict that over the next few years, the known QPEs will
continue to brighten, doubling their luminosity within about a decade. These
systems will eventually cease to episodically flare and brighten, once they
will approach mass transfer rates comparable to the Eddington rate, with
$L_{\rm Edd}\approx 10^{44}\,M_{\bullet,6}\;\rm erg\,s^{-1}$, at which point
radiation pressure regulates the black hole’s accretion rate, and the variable
mass transfer will no longer modulate the emitted radiation. The emitted
radiation could strongly affect the star itself, potentially enhancing the
mass transfer rate, coupling the star’s evolution to the radiation that the
transfer of mass produces. This effect, discussed in Krolik & Linial (2022)
could significantly affect the long term evolution of the system, possibly
regulating the runaway mass-loss. Recently, Lu & Quataert (2022) also
addressed a possible self-regulating effect, where interactions between the
star and the accretion flow it produces, ultimately determine the long term
steady-sate evolution of the system.

We emphasize the following hierarchy between the relevant timescales - the
star’s dynamical time is comparable to the orbital period $P$. The system
brightens over a timescale $\tau_{L}$ corresponding to $\xi^{1-k}\approx
10^{6}$ orbits. The apparent timescale, $\tau_{\rm app}$ overestimates the
actual evolution time by a factor $\xi^{-1}\approx 100$. Finally, the apparent
timescale is shorter than $\tau_{\rm GW}$ by several orders of magnitude.

###  5.3 Rates

In section 3.2 we estimated the stellar-EMRI formation rate through a
combination of Hills’ mechanism, two-body angular momentum relaxation and GW
emission (equation 24). The time spent at each logarithmic luminosity range is
approximately

$\tau_{\rm runaway}(\left<L\right>)\approx
P\left(\frac{\left<L\right>P}{M_{\star}(R_{g}/r_{\rm
diss})c^{2}}\right)^{1/k-1}\,,$ (41)

where the expression in the parentheses is the per-flare radiated energy
relative to the total available energy budget. Normalized by representative
values

$\tau_{\rm runaway}(\left<L\right>)\approx 10^{2}\,P_{10}^{1.3}L_{\rm
41}^{-0.75}\left(\frac{r_{\rm diss}}{10R_{g}}\right)^{0.75}\;\rm yr\,,$ (42)

Combining the formation rate with the expected lifetime implies that the
average number of QPE systems per galaxy is of order

$N_{\rm QPE}(\left<L\right>)\approx\mathcal{R}_{\rm MS-
EMRI,Hills}\times\tau_{\rm runaway}\approx\\\ 10^{-4}\;P_{10}^{1.3}L_{\rm
41}^{-0.75}M_{\bullet,6}^{-0.25}\left(\frac{r_{\rm
diss}}{10R_{g}}\right)^{0.75}\left(\frac{f_{b}}{0.1}\right)\,.$ (43)

The eROSITA survey detected eRO-QPE1 at $z=0.0505$, covering a volume
containing roughly $10^{5}$ galaxies hosting an SMBH with $M_{\bullet}\approx
10^{6}\,M_{\odot}$ (Arcodia et al., 2021; Metzger et al., 2021). The discovery
of 2 QPE systems within this volume can therefore be explained by our
theoretical abundance estimate.

Even though the lifetime of brighter systems is shorter, they could be seen to
further distances. Hence, in a flux limited survey, bright systems would
dominate the number of detections. The exact luminosity upper cutoff is
unclear - and is possibly related to the Eddington luminosity of the SMBH. The
observed luminosities of the current systems are about 1-3 orders of magnitude
smaller than the SMBH Eddington limit (Arcodia et al., 2021; Metzger et al.,
2021), as summarized in table 1. This discrepancy might be explained if either
QPEs cannot brighten past their current observed levels and the luminosity
stops tracing the mass transfer rate, or alternatively, if mass-loss is self-
regulated through the interaction between the star and the accretion flow or
the emitted flux (e.g., Lu & Quataert, 2022; Krolik & Linial, 2022).

##  6 Summary and discussion

We studied systems the formation and evolution of systems containing a main-
sequence star orbiting an SMBH on mildly eccentric orbits, episodically
shedding mass at every pericenter passage. We highlight the following key
properties of our model.

  * •

Combination of two-body scatterings, the Hills’ mechanism and GW emission
produces stars on mildly eccentric orbits around SMBHs, at rates of $\sim
10^{-7}-10^{-5}\;\rm yr^{-1}$.

  * •

These stars inspiral towards the SMBH due to GW emission until they begin to
shed mass, once $r_{p}\approx r_{\rm tidal}$, corresponding to orbital periods
of hour-day timescale (essentially dictated by the star’s dynamical
timescale).

  * •

The mass transfer rate is initially governed by the continued GW inspiral and
gradually increases with time, but after $\sim 10^{6}$ orbits, $\sim 10^{-3}$
of the stellar mass has been shed and the stellar response to mass loss begins
to dominate the enhancement of mass transfer, over GWs. Mass transfer then
increases in a runaway fashion, on ever decreasing timescales.

  * •

Mass transfer rates increase from zero to as high as $M_{\star}/\tau_{\rm
dyn}\approx 10^{3}\,\rm M_{\odot}\,yr^{-1}$, with all intermediate rates
realized along the process. As the system continues to evolve, the time spent
around any mass transfer rate becomes shorter, the higher $\dot{M}_{\star}$
is.

We have shown that such systems can naturally produce some of the key features
of the known QPE systems, and in particular, their recurrence timescale,
emission energy budget, and their rates. We address a few additional aspect
and caveats related to the association between stellar-EMRIs and QPEs.

###  6.1 Long term evolution of known QPE systems

A few of the known QPE systems already have years-long observational
baselines, which can be compared against the predictions of our model. We
emphasize that if mass transfer is self-regulated through interaction between
the resulting radiated flux or the accretion flow, the system may not
necessarily gradually brighten, and we might underestimate the remaining
system lifetime (see end of section 5.3). Nonetheless, we address the
properties of the known QPEs in this context. The most recently identified QPE
candidate, 2MASXJ0249, was discovered in archival XMM-Newton observations,
showing two QPE flares in 2006 (Chakraborty et al., 2021). Based on the 2006
flare data, we estimate the system’s evolution timescale as $\tau_{L}\approx
12\,\rm yr$, implying that the system is expected to brighten and then stop
flaring on this timescale. Indeed, followup observations of the system in
2021, 15 years later, revealed no QPE activity, in agreement with our model’s
prediction. The period and black hole mass of this system are remarkably
similar to those of eRO-QPE2. The current luminosity of eRO-QPE2 is slightly
higher than the 2006 luminosity of 2MASXJ0249, implying a shorter evolution
time of only 5 years. We predict a similar fate and evolution track for the
two systems - we expect eRO-QPE2 to disappear within the next decade, as
2MASXJ0249 did between the 2006 and 2021 observation epochs.

The QPE system with the longest observational baseline is RX J1301.9+2747,
which was identified in 2019 using XMM-Newton (Giustini et al., 2020).
Archival observations show that this object was flaring already in the year
2000, with a similar recurrence time of $\sim 4.5\,\rm hr$. We estimate an
evolution time of roughly $\tau_{L}\approx 16\,\rm yr$, implying that the
system was expected to roughly double in brightness between 2000 and 2019. The
observed fluxes are consistent with being constant, given the uncertainties in
the spectral fitting performed in Giustini et al. (2020). As additional
observations are accumulated over the next few years, we anticipate that this
system will brighten, until it will eventually approach its Eddington
luminosity, at roughly $\sim 10^{44}\,\rm erg\,s^{-1}$, and then disappear.

The first identified QPE system, GSN-069 was detected by XMM-Newton in 2018.
This system has displayed a long term gradual decay in X-ray emission since
2010, possibly indicating a long-lived tidal disruption event (Miniutti et
al., 2019). This gradual long term decay is not explained by our model. The
properties of the QPE flares detected in 2018-2019 imply an evolution time of
about a decade. However, no QPE activity was detected in a sufficiently long
observation performed in 2014. Our model suggests that the QPE was dimmer by a
factor of order of unity in 2014, while the gradually evolving X-ray
background was brighter by a factor of $\sim 2$. The overall amplitude of the
QPE flares with respect to quiescence was therefore a few times lower than in
2018-2019.

###  6.2 Quasi-periodicity.

QPEs are only quasi-periodic, with a typical scatter in the flare timing of
order $\delta t/P\approx 0.1$. One possible origin of this quasi-periodicity
may be stellar oscillations. The strong tidal forces exerted on the star can
naturally excite different modes of stellar oscillations. Low order radial
modes result in periodic changes in the mean stellar density, and as a result
- in $r_{\rm tidal}$. If the stellar oscillations are not in phase with the
orbital motion, the peak of mass transfer may occur slightly before or after
pericenter. The resulting flares will therefore occur quasi-periodically, with
a scatter that depends on the stellar density fluctuations, the orbital
eccentricity and on the (average) tidal radius penetration, $\xi$.

###  6.3 Formation of narrow-line region.

The accretion of a total mass $\delta M_{\star}$ onto the black hole is
capable of producing a total of $N_{\rm ph}\approx(\delta
M_{\star}/m_{e})\alpha^{-2}\approx 10^{63}\,(\delta M_{\star}/\rm M_{\odot})$
hydrogen-ionizing photons, where $m_{e}$ and $\alpha$ are the electron mass
and the fine structure constant. The known QPE systems have shed about $\delta
M_{\star}=\xi M_{\star}\approx 10^{-2}\,\rm M_{\odot}$, which could have
ionized a total of $\sim 10^{5}\,M_{\odot}$ in neutral hydrogen. The elapsed
time since the onset of mass transfer is roughly $\Delta t\approx\tau_{\rm
GW}(P/\tau_{\rm GW})^{0.25}\approx 3500\,\rm yr$, implying a maximal extent of
$\sim 1\,\rm kpc$ of the ionized region (although most of the radiation is
emitted on the much shorter timescale, $\tau_{L}$, reducing the size of the
ionized region). Whether or not the QPE emission can sustain a narrow-line
region depends on the recombination rate of the ionized gas. If this rate is
low enough, our model could potentially explain the recently identified
narrow-line regions in the QPE host galaxies (Wevers et al., 2022a), however,
more quantitative analysis of this problem is beyond the scope of this paper.

###  6.4 Additional dynamical effects

In our calculation of the QPE production rate we made several simplifying
assumptions regarding galactic nuclei dynamics - we neglected resonant
relaxation, stellar collisions and deviations from spherical symmetry.
Additionally, we assumed that all stars engulfing the SMBH are of equal mass.
More massive objects (i.e., stellar black holes) tend to sink towards the SMBH
due to mass segregation (e.g., Alexander & Hopman, 2009; Linial & Sari, 2022).
The concentration of stellar black holes in small radii will tend to scatter
low-mass main-sequence stars away from the SMBH, altering the stellar-inspiral
rates we report here.

We would like to thank Nicholas Stone, Riccardo Arcodia, Tsvi Piran, Julian
Krolik, Eliot Quataert and Brian Metzger for interesting discussions that
contributed to different stages of this work. IL thanks support from the Adams
Fellowship, the Rothschild Fellowship and the Gruber Foundation.

##  Appendix A Orbital period of stellar EMRIs

At the onset of mass transfer, a stellar inspiral has some residual
eccentricity $e$ and a pericenter distance

$r_{p}\approx r_{\rm tidal}=\alpha_{\rm
RL}R_{\star}\left(M_{\bullet}/M_{\star}\right)^{1/3}\,,$ (A1)

where $\alpha_{\rm RL}\approx 2$ is an order unity prefactor that may weakly
depend on $e$. The semi-major axis $a=r_{p}/(1-e)$, corresponds to an orbital
period

$P=2\pi\sqrt{\frac{a^{3}}{GM_{\bullet}}}=2\pi\alpha_{\rm
RL}^{3/2}(1-e)^{-3/2}\tau_{\rm dyn}\,,$ (A2)

showcasing the fact that at the Roche limit, the pericenter passage time is
similar to the star’s dynamical timescale, $\tau_{\rm
dyn}=\sqrt{R_{\star}^{3}/GM_{\star}}$. Taking the main-sequence star relation
$R_{\star}\propto M_{\star}^{0.8}$, the stellar mass is related to the period
as

$M_{\star}(P)\approx 0.9\,M_{\odot}\,\left(\frac{P}{10\,\rm
hr}\right)^{1.4}\,,$ (A3)

where we used the fiducial $e=0.1$.

##  Appendix B Mass transfer rate from eccentric inspirals

Consider a star of mass $M_{\star}$ and radius $R_{\star}$ on a marginally
mass-shedding orbit $r_{p}\lesssim r_{\rm tidal}$. We follow Linial & Sari
(2017) to obtain the mass transfer rate as a function of the Roche-Lobe
overfilling extent, $\xi$ (see equation 32).

Since the orbit is slightly eccentric, the transfer of mass from the star is
varying over the course of an orbital period. The star only overfills its
Roche lobe (and transfers mass) during a fraction of the orbit

$\frac{\Delta t}{P}=\xi^{1/2}\sqrt{\frac{(1-e)^{3}}{2\pi^{2}e}}\,,$ (B1)

valid for $\xi\ll e$. In the opposite regime (not considered here), where
$\xi>e$, part of the star is always exceeding the Roche lobe. During this
interval, the instantaneous Roche-lobe overfilling varies as

$\xi_{\rm i}(t)\approx\xi\left(1-\left(\frac{t}{\Delta
t/2}\right)^{2}\right)\,,$ (B2)

where $t=0$ was set to be the pericenter passage time, and the above
expression is valid for $-\Delta t/2<t<\Delta t/2$.

The instantaneous mass-transfer rate through the L1 nozzle is given by

$\dot{M}_{\star,\rm i}(\xi_{\rm i})\approx\rho vS\,,$ (B3)

where $\rho$ is the typical flow density, $v$ is the velocity of the flow
through the nozzle, and $S$ is the cross section of the flow. Assuming a
polytropic stellar structure with index $n$, we estimate

$\displaystyle\rho\approx K_{\rho}(n)\frac{M_{\star}}{R_{\star}^{3}}\xi_{\rm
i}^{n}\,,$ (B4) $\displaystyle v\approx
c_{s}=\frac{1}{\sqrt{n}}\sqrt{\frac{GM_{\star}}{R_{\star}}}\xi_{\rm
i}^{1/2}\,,$ (B5) $\displaystyle S\approx 2\pi R_{\star}^{2}\xi_{\rm i}\,,$
(B6)

where $K_{\rho}(n)$ is a dimensionless prefactor (see appendix B of Linial et
al., 2021) and $c_{s}$ is the material’s initial sound speed. All the
quantities were evaluated on the L1 Roche equipotential surface, located
beneath the stellar surface at pericenter, at depth $\xi R_{\star}$. Put
together

$\dot{M}_{\star,\rm i}(\xi_{\rm i})\approx
K_{\rho}(n)\frac{2\pi}{\sqrt{n}}\frac{M_{\star}}{\tau_{\rm
dyn}}\xi_{i}^{n+3/2}\,,$ (B7)

where $\tau_{\rm dyn}=\sqrt{R_{\star}^{3}/GM_{\star}}$.

The orbit averaged mass transfer rate is given by

$\dot{M}_{\star}=\frac{1}{P}\int_{-P/2}^{P/2}\dot{M}_{\star,\rm i}(\xi_{\rm
i}(t))\,dt=\alpha_{\rm MT}\frac{M_{\star}}{P}e^{-1/2}\xi^{n+2}\,,$ (B8)

where $\alpha_{\rm MT}$ is a dimensionless prefactor, which we approximate as
$\alpha_{\rm MT}\approx 9$ for $n=3/2$ and $\alpha_{\rm M}\approx 2$ for
$n=3$.

## References

  * Alexander & Hopman (2009) Alexander, T., & Hopman, C. 2009, ApJ, 697, 1861, doi: 10.1088/0004-637X/697/2/1861
  * Amaro-Seoane (2018) Amaro-Seoane, P. 2018, Living Reviews in Relativity, 21, 4, doi: 10.1007/s41114-018-0013-8
  * Arcodia et al. (2021) Arcodia, R., Merloni, A., Nandra, K., et al. 2021, Nature, 592, 704, doi: 10.1038/s41586-021-03394-6
  * Bahcall & Wolf (1976) Bahcall, J. N., & Wolf, R. A. 1976, ApJ, 209, 214, doi: 10.1086/154711
  * Bahcall & Wolf (1977) —. 1977, ApJ, 216, 883, doi: 10.1086/155534
  * Chakraborty et al. (2021) Chakraborty, J., Kara, E., Masterson, M., et al. 2021, ApJ, 921, L40, doi: 10.3847/2041-8213/ac313b
  * Dai & Blandford (2013) Dai, L., & Blandford, R. 2013, MNRAS, 434, 2948, doi: 10.1093/mnras/stt1209
  * Generozov et al. (2018) Generozov, A., Stone, N. C., Metzger, B. D., & Ostriker, J. P. 2018, MNRAS, 478, 4030, doi: 10.1093/mnras/sty1262
  * Gezari (2021) Gezari, S. 2021, ARA&A, 59, 21, doi: 10.1146/annurev-astro-111720-030029
  * Giustini et al. (2020) Giustini, M., Miniutti, G., & Saxton, R. D. 2020, A&A, 636, L2, doi: 10.1051/0004-6361/202037610
  * Hills (1988) Hills, J. G. 1988, Nature, 331, 687, doi: 10.1038/331687a0
  * King (2020) King, A. 2020, MNRAS, 493, L120, doi: 10.1093/mnrasl/slaa020
  * King (2022) —. 2022, MNRAS, 515, 4344, doi: 10.1093/mnras/stac1641
  * Krolik & Linial (2022) Krolik, J. H., & Linial, I. 2022, arXiv e-prints, arXiv:2209.02786.  https://arxiv.org/abs/2209.02786
  * Linial et al. (2021) Linial, I., Fuller, J., & Sari, R. 2021, MNRAS, 501, 4266, doi: 10.1093/mnras/staa3969
  * Linial & Sari (2017) Linial, I., & Sari, R. 2017, MNRAS, 469, 2441, doi: 10.1093/mnras/stx1041
  * Linial & Sari (2022) —. 2022, arXiv e-prints, arXiv:2206.14817.  https://arxiv.org/abs/2206.14817
  * Liu et al. (2022) Liu, Z., Malyali, A., Krumpe, M., et al. 2022, arXiv e-prints, arXiv:2208.12452.  https://arxiv.org/abs/2208.12452
  * Lu & Quataert (2022) Lu, W., & Quataert, E. 2022, arXiv e-prints, arXiv:2210.08023.  https://arxiv.org/abs/2210.08023
  * Marconi & Hunt (2003) Marconi, A., & Hunt, L. K. 2003, ApJ, 589, L21, doi: 10.1086/375804
  * McConnell & Ma (2013) McConnell, N. J., & Ma, C.-P. 2013, ApJ, 764, 184, doi: 10.1088/0004-637X/764/2/184
  * Metzger et al. (2021) Metzger, B. D., Stone, N. C., & Gilbaum, S. 2021, arXiv e-prints, arXiv:2107.13015.  https://arxiv.org/abs/2107.13015
  * Miniutti et al. (2019) Miniutti, G., Saxton, R. D., Giustini, M., et al. 2019, Nature, 573, 381, doi: 10.1038/s41586-019-1556-x
  * Payne et al. (2021) Payne, A. V., Shappee, B. J., Hinkle, J. T., et al. 2021, ApJ, 910, 125, doi: 10.3847/1538-4357/abe38d
  * Peters (1964) Peters, P. C. 1964, Physical Review, 136, 1224, doi: 10.1103/PhysRev.136.B1224
  * Press & Teukolsky (1977) Press, W. H., & Teukolsky, S. A. 1977, ApJ, 213, 183, doi: 10.1086/155143
  * Rees (1988) Rees, M. J. 1988, Nature, 333, 523, doi: 10.1038/333523a0
  * Sari & Fragione (2019) Sari, R., & Fragione, G. 2019, ApJ, 885, 24, doi: 10.3847/1538-4357/ab43df
  * Soberman et al. (1997) Soberman, G. E., Phinney, E. S., & van den Heuvel, E. P. J. 1997, A&A, 327, 620.  https://arxiv.org/abs/astro-ph/9703016
  * Tremaine et al. (2002) Tremaine, S., Gebhardt, K., Bender, R., et al. 2002, ApJ, 574, 740, doi: 10.1086/341002
  * Volonteri et al. (2021) Volonteri, M., Habouzit, M., & Colpi, M. 2021, Nature Reviews Physics, 3, 732, doi: 10.1038/s42254-021-00364-9
  * Wevers et al. (2022a) Wevers, T., Pasham, D. R., Jalan, P., Rakshit, S., & Arcodia, R. 2022a, arXiv e-prints, arXiv:2201.11751.  https://arxiv.org/abs/2201.11751
  * Wevers et al. (2022b) Wevers, T., Coughlin, E. R., Pasham, D. R., et al. 2022b, arXiv e-prints, arXiv:2209.07538.  https://arxiv.org/abs/2209.07538
  * Xian et al. (2021) Xian, J., Zhang, F., Dou, L., He, J., & Shu, X. 2021, ApJ, 921, L32, doi: 10.3847/2041-8213/ac31aa
  * Zhao et al. (2021) Zhao, Z. Y., Wang, Y. Y., Zou, Y. C., Wang, F. Y., & Dai, Z. G. 2021, arXiv e-prints, arXiv:2109.03471.  https://arxiv.org/abs/2109.03471


# Optimizing Real-Time Performances for Timed-Loop Racing under F1TENTH
††thanks: We thank Zhaojian Li from Michigan State University and Marko
Bertogna from University of Modena for fruitful discussions on the topic. We
also thank Nvidia for providing the Jetson Xavier AGX used in this paper.

Nitish A Gupta  North Carolina State University
Department of Computer Science
<EMAIL_ADDRESS>Kurt Wilson  University of Central Florida
Department of Computer Science
<EMAIL_ADDRESS>Zhishan Guo  North Carolina State University
Department of Computer Science
<EMAIL_ADDRESS>
###### Abstract

Motion planning and control in autonomous car racing are one of the most
challenging and safety-critical tasks due to high speed and dynamism. The
lower-level control nodes are expected to be highly optimized due to resource
constraints of onboard embedded processing units, although there are strict
latency requirements. Some of these guarantees can be provided at the
application level, such as using ROS2’s Real-Time executors. However, the
performance can be far from satisfactory as many modern control algorithms
(such as Model Predictive Control) rely on solving complicated online
optimization problems at each iteration. In this paper, we present a simple
yet effective multi-threading technique to optimize the throughput of online-
control algorithms for resource-constrained autonomous racing platforms. We
achieve this by maintaining a systematic pool of worker threads solving the
optimization problem in parallel which can improve the system performance by
reducing latency between control input commands. We further demonstrate the
effectiveness of our method using the Model Predictive Contouring Control
(MPCC) algorithm running on Nvidia’s Xavier AGX platform.

##  I Introduction

High-speed autonomous racing is a challenging and important application within
autonomous driving, wherein safety is a crucial consideration. In recent
years, this application has gained increased attention from the real-time
systems community. Autonomous navigation at high-speed demands tighter bounds
on response times of the localization and control algorithms. Several
competitions such as F1TENTH [1] and Indy Autonomous Challenge (IAC) [2] are
held once to twice every year to validate and improve the functional limits of
the existing software stack. The F1TENTH competition focuses on racing 1/10th-
scale vehicles in a dynamic head-to-head racing environment. At this scale,
the problem of autonomous navigation narrows down to designing an efficient
and real-time safe trajectory following an algorithm in a dynamic environment.
To this end, different control algorithms [3] are employed by the community
such as Reactive, Model Predictive Control (MPC) [4], raceline-optimization
[5] or Learning-based algorithms [6]. In this paper, we focus our attention on
an online optimization-based MPC algorithm, particularly, Model Predictive
Contouring Control (MPCC). MPCC algorithm provides high-performance autonomous
navigation accounting for non-linearities in dynamics efficiently at high
speed. This is mainly because the MPCC solver accounts for model deviations at
each step and applies relevant corrections dynamically. Therefore, the closed-
loop online optimization as in MPCC particularly tends to be more reliable
than other algorithms mentioned earlier. In the context of autonomous racing,
this algorithm was originally demonstrated in [7] on 1:43 scale RC cars.

Figure 1: F1TENTH Platform

As compared to Autonomous Vehicle (AV) applications for urban driving,
autonomous vehicle racing platform pose a different set of challenges such as
frequent edge cases, resource constraints, and critical reaction/response
time. Thus it becomes highly difficult to deploy algorithms that rely on
online optimization under these constraints. However, control algorithms such
as MPCC rely on solving a constrained convex optimization problem over a
finite horizon to provide a dynamic and efficient navigation strategy.

To apply such a technique for high-speed autonomous navigation on embedded
platforms such as Nvidia Xavier NX/AGX, the latency between actuation command
signals needs to be optimized from both the application and system
perspective. We discuss this problem in the next sections and present a simple
yet elegant approach for reducing this latency of online optimization control
algorithms for timed-loop autonomous racing platforms. In particular, we
propose a multi-threading based algorithm on top of the ROS2 implementation of
MPCC, such that the online optimization can be performed by multiple threads
in a synchronized fashion. This enables us to improve the publish latency of
the control algorithm while optimizing the CPU utilization in a resource-
constrained environment.

##  II Background

This section briefly explains ROS2 architecture and then provides the details
of Model Predictive Contouring Control, as well as the F1TENTH hardware and
software stack.

###  II-A Architecture of ROS2

ROS2 is a popular robotics framework that leverages the publisher-subscriber
paradigm to provide real-time inter-process communication via the data
distribution service layer in a distributed computing environment. The
simplest form of computing unit in ROS2 is called nodes which may in turn
contain several callbacks. Each callback in itself is an executable function
that could be either time-triggered or data-triggered. The message transfer in
ROS is handled through a topic. A node can publish and/or subscribe to one or
several topics. A subscriber callback of a node is triggered whenever the
topic it subscribes to receives a message from a publisher. A node publishes a
topic using a timer which is used to set the rate at which the timer callback
will be triggered and thus the message will be published. Beneath the
abstraction, all callbacks within a single node are scheduled for execution by
an executor. The default ROS2 executor behavior has been thoroughly studied
[8, 9] and several customized schedulers also have been proposed such as
dynamic-priority [10] and fixed-priority [11] schedulers for executors.

###  II-B Model Predictive Contouring Control

Model Predictive Control is a popular optimal control algorithm used for
process control and trajectory following in linear and non-linear systems. MPC
achieves the prediction of multiple control inputs over a finite receding
horizon by solving a constrained convex optimization problem. For optimal
performance, the cost function formulation must capture the process goals and
the constraints from the model and environment. In the case of autonomous
racing, the cost function can be represented by the progress along the race
line while maintaining a safe distance from the track borders and obstacles.
At each time step, the first predicted input is applied and the optimization
is repeated over the new horizon.

In this paper, we consider a reconfigured and optimized implementation of
Model Predictive Contouring Control by Liniger et al. [7]. However, instead of
a 1/43rd scale, we model it for a 1/10th scale by updating the vehicle
dynamics for F1TENTH using the methodology described in [12]. The state
representation is based on the vehicle bicycle model and control input
consists of duty cycle $d$ and steering angle $\delta$.

$\dot{X}=v_{x}\cos{\varphi}-v_{y}\sin{\varphi}$ (1)
$\dot{Y}=v_{x}\sin{\varphi}+v_{y}\cos{\varphi}$ (2) $\dot{\varphi}=\omega$ (3)
$\dot{v_{x}}=(F_{r,x}-F_{f,y}\sin{\delta}+mv_{y}\omega)/m$ (4)
$\dot{v_{y}}=(F_{r,y}+F_{f,y}\cos{\delta}-mv_{x}\omega)/m$ (5)
$\dot{\omega}=(F_{f,y}l_{f}\cos{\delta}-F_{r,y}l_{r})/I_{z}$ (6)

In the above differential equations, $X$, $Y$, and $\varphi$ is vehicle’s
position and heading angle relative to the inertial frame. Similarly, the
$v_{x}$, $v_{y}$, and $\omega$ represent linear and angular velocities. The
$F_{f}$ and $F_{r}$ represent tire forces along the lateral and longitudinal
direction which are modeled using a simplified Pacejka Tire Model. Finally,
$m$, $l_{r}$, and $l_{f}$ is vehicle mass and front/rear distance from
vehicle’s $COG$, respectively. The dynamics modeling part is beyond the scope
of this paper and can be found in [7, 12]. Considering the goal of maximizing
the progress measured along the parameterized centerline and satisfying
several equality and inequality constraints, a convex optimization problem is
formulated. The online optimization is solved using a high-performance
Quadratic Programming (QP) solver over a finite horizon $N$, such that all
constraints are satisfied. Only the first input from the entire horizon
solution is applied as input commands to the vehicle. This process is then
repeated in the next cycle over the entire horizon. This ensures the stability
of the control algorithm and keeps up with slight inconsistencies in vehicle
dynamics modeling. The horizon length parameter, N, controls the tracking
accuracy as it provides farsightedness to the vehicle and allows for better
planning.

###  II-C F1TENTH Platform

The F1TENTH platform is built upon a 1/10th scale chassis which hosts a
differential 4-wheel drive and Ackermann steering mechanism. The steering and
throttle commands are sent to the motors through the VESC 6 MkV controller
which also allows for motor feedback. An embedded platform such as Nvidia
Jetson AGX Xavier is used as an onboard processing unit. The sensor stack
consists of a Hokuyo UST-10LX Lidar, a 9-DoF IMU, and an optional ZED Depth
camera. The Lidar has a detection range of 10 meters, a scan range of 270∘ at
0.25∘ angular resolution, and a scan speed of 40 Hz. The Lidar along with IMU
sensor and odometry data enables mapping and localization of the ego vehicle
using algorithms such as Google’s cartographer [13] or slam-toolbox [14]. Fig.
2 shows the essential modules and data flow in a typical F1TENTH platform. The
communication between the sensors, embedded board, and speed controllers is
facilitated using ROS2.

Figure 2: F1TENTH Software Stack

The sensor stack and localization module enable state estimation of the ego
and opponent vehicle, which is then communicated to the MPCC algorithm for
calculating a control input command. The MPCC controller solves and publishes
the control input command consisting of $(d,\delta)$, to the lower level VESC
and steering controller nodes for actuation of the vehicle. Upon actuation,
the state estimation repeats and transitions the system to the next state.
However, the applied target input command is not achieved instantly due to the
physical constraints of the vehicle. Thus the updated state is not necessarily
a linear interpolation of the applied control input, but instead represents an
integration of model dynamics over small finite intervals. As we will see in
the next section, we try to capture this behavior using the iterative Runga-
Kutta (RK4) numerical method.

##  III Problem Statement

As the main focus of this paper is optimizing the control performance, we
consider a simplified representation of the framework explained in Section II
as two main nodes implemented in ROS2. As shown in Fig. 3, the first node is
the mpcController node that initializes and solves an online convex
optimization problem. The mpcController node subscribes to vehicle state on
the nextState topic and publishes the vehicle control input commands to
inputCmd topic. The second node on the control loop is mpcIntegrator node
which represents the model state integrator. The mpcIntegrator node subscribes
to control inputs commands on mpcInput topic and transitions the vehicle to
the next state which is then published to nextState topic. On a physical
F1TENTH platform, the role of the state integrator is distributed between
actuating the vehicle and then estimating the next state through the
localization module and other sensors such as IMU and Odometry. As the
combined sensor stack on the F1TENTH platform publishes data at a high
frequency ranging anywhere from 40 Hz to 100 Hz, the state updates are
processed faster than the control optimization. This behavior is captured in
our implementation by executing a state simulator thread within mpcIntegrator
node at 1000 Hz or 1 ms interval. The state integrator updates the states of
the vehicle using the most recently available control input command, thus it
is important that the mpcController node publishes new control inputs at a
high frequency to maintain data freshness.

Figure 3: Simplified MPCC Control Loop

##  IV Implementation

The mpcController node operates on the current vehicle state to solve a convex
optimization problem over a finite horizon $N$. The controller then publishes
the first control input from the solved horizon. As the autonomous racing
control problem is highly dynamic, the time interval between each new command
message is crucial to optimize to keep the vehicle as reactive as possible.
From the application perspective, the optimization solver is only dependent on
the vehicle’s current state. On the other hand, the mpcIntegrator node updates
the system state at a much higher frequency and publishes on the nextState
topic. This gives us an opportunity to parallelize the controller node using
multi-threading, wherein each worker thread is solving optimization based on
the most recently available vehicle state. However, the control input commands
that are based on two very similar states do not provide much benefit to
overall system throughput. Similarly, a new command message based on an old
vehicle state needs to be discarded when a command message based on a newer
state has already been published. Such cases need to be handled in a
systematic and synchronized manner to ensure reliability and data freshness.

Algorithm 1 Multi-threaded Controller Node

1:$last\\_queued\\_timestamp\leftarrow 0$

2:$last\\_output\\_time\leftarrow 0$

3:declare $atomic\\_msg$

4:declare $result\\_msg$

5:declare $worker\\_thread\\_cv$

6:declare $publish\\_thread\\_cv$

7:function Process_State_Message(State, Time)

8: /* if enough time has passed */

9: if $State.timestamp-last\\_queued\\_timestamp>MIN\\_GAP$ then

10: atomic_msg.exchange(State)

11: /* notify the worker threads */

12: worker_thread_cv.notify_one()

13: end if

14:end function

15:function workerthread

16: while should_run do

17: worker_thread_cv.wait()

18: State = atomic_msg.exchange(nullptr)

19: /* update the shared timer */

20: last_queued_timestamp = State.timestamp

21: result = mpc(State)

22: /* make sure another worker did not complete with a newer state */

23: if State.timestamp $<$ last_output_time then

24: break

25: end if

26: last_output_time = State.timestamp;

27: result_msg = result;

28: publish_thread_cv.notify_one()

29: end while

30:end function

31:function publishthread

32: while true do

33: publish_thread_cv.wait()

34: publish(result_msg)

35: end while

36:end function

We create a thread pool to run the mpcController on different nextState
messages in parallel. The state subscriber callback of mpcController node
checks the timestamp of each incoming message and compares it to the timestamp
of the last state message accepted by the thread pool. If the timestamp of an
incoming message is too close to the most recently accepted state message by
any thread, the received message is discarded. We call this interval
$MIN\\_GAP$, and the above check occurs on line-9 of Algorithm 1. This
condition prevents the thread pool from being saturated with state messages
from consecutive timestamps and ensures a well-spaced interval between
published control input messages. The unassigned worker threads are required
to wait until a state message arrives with a qualified timestamp. The next
received state message is then passed to one of the available worker threads
in the pool. The assigned thread updates the shared timestamp and starts
working on the MPCC optimization call. Once the thread completes the call, it
checks whether another worker has already published a result based on newer
state data. If not, it wakes the publisher thread with the new control input
message to publish on mpcInput topic. A global state variable
$last\\_output\\_time$, maintains the timestamp of the last used state data to
ensure that no control message is sent based on older state messages. This
check occurs on line-23 of the algorithm. This extra check must occur since a
thread working on state data could be delayed due to operating system
preemptions. During this time another thread could have finished with newer
state data and thus invalidating the control message obtained from the delayed
worker thread.

To ensure that the worker threads do not interfere with each other, we pin
each thread to a CPU core and assign each thread a high priority.

##  V Preliminary Evaluation

In this section, we present an initial evaluation to show the efficacy of our
proposed optimization method. We model the C++ implementation of MPCC
algorithm using the latest version of ROS2 Humble, based on the representation
described in Section III. We performed our experiments on Nvidia Jetson AGX
Xavier running Ubuntu 20.04 with Linux Real-Time kernel. We set $MIN\\_GAP$
parameter set to 10 ms, that is no worker thread can be started within 10 ms
of another worker. Therefore, the target publish interval is 10 ms. The
evaluation results are shown in Fig. 4(a) and Fig. 4(b).

The experiments show that adding additional workers significantly reduces the
maximum and average control publish interval length. This enables vehicle
actuation on more recent and frequent input command signal. Adding a second
thread brings the absolute maximum rate from 50 ms down to 30 ms. The three-
thread scenario keeps the maximum publishing rate much closer to the target
rate than with just one thread, demonstrating the benefits of synchronized
parallel processing. While the maximum observed publishing rates are close for
two and three threads, the variance of the publishing interval is
significantly decreased when using three threads.

(a) Distribution of the MPC control output interval as more threads are added
to the worker pool. Adding threads both reduces the mean and variance of the
interval

(b) Ranges of the MPC control output intervals. Using three threads allows the
control interval to stay much closer to the enforced minimum of 10ms. Note
that there are cases where the interval was less than 10ms; this can happen
when two threads started 10ms apart but completed around the same time and in
order

Figure 4: Results from testing the threaded MPC algorithm on a Nvidia Jetson
AGX Xavier

##  VI Conclusion

In this paper, we presented a challenge that is commonly faced while
developing optimization algorithms for autonomous vehicle racing. We then
presented a multi-threading-based solution to optimize the control throughput
performance of online optimization algorithms on resource-constrained
platforms such as F1TENTH. In the future, we plan to investigate other ways of
optimizing the performance by utilizing the shared GPU memory on embedded
platforms for performing expensive solver matrix calculations. Finally, we
plan to demonstrate our proposed algorithm during an upcoming F1TENTH
competition.

## References

  * [1] M. O’Kelly et al., “F1TENTH: An Open-source Evaluation Environment for Continuous Control and Reinforcement Learning,” Proceedings of Machine Learning Research, vol. 123, 2020.
  * [2] A. Wischnewski et al., “Indy autonomous challenge - autonomous race cars at the handling limits,” in 12th International Munich Chassis Symposium 2021 (P. Pfeffer, ed.), (Berlin, Heidelberg), pp. 163–182, Springer Berlin Heidelberg, 2022.
  * [3] J. Betz, H. Zheng, A. Liniger, U. Rosolia, P. Karle, M. Behl, V. Krovi, and R. Mangharam, “Autonomous vehicles on the edge: A survey on autonomous vehicle racing,” IEEE Open Journal of Intelligent Transportation Systems, vol. 3, pp. 458–488, 2022.
  * [4] E. Thilén, “Robust model predictive control for autonomous driving,” 2017\.
  * [5] M. O’Kelly, H. Zheng, A. Jain, J. Auckley, K. Luong, and R. Mangharam, “Tunercar: A superoptimization toolchain for autonomous racing,” in 2020 IEEE International Conference on Robotics and Automation (ICRA), pp. 5356–5362, IEEE, 2020.
  * [6] P. R. Wurman, S. Barrett, K. Kawamoto, J. MacGlashan, K. Subramanian, T. J. Walsh, R. Capobianco, A. Devlic, F. Eckert, F. Fuchs, et al., “Outracing champion gran turismo drivers with deep reinforcement learning,” Nature, vol. 602, no. 7896, pp. 223–228, 2022.
  * [7] A. Liniger et al., “Optimization-based autonomous racing of 1: 43 scale rc cars,” Optimal Control Applications and Methods, vol. 36, no. 5, pp. 628–647, 2015.
  * [8] D. Casini et al., “Response-time analysis of ros 2 processing chains under reservation-based scheduling,” in ECRTS, 2019.
  * [9] Y. Tang et al., “Response time analysis and priority assignment of processing chains on ros2 executors,” in RTSS, IEEE, 2020.
  * [10] A. A. Arafat et al., “Response time analysis for dynamic priority scheduling in ros2,” in DAC, ACM/IEEE, 2022.
  * [11] H. Choi et al., “Picas: New design of priority-driven chain-aware scheduling for ros2,” in RTAS, IEEE, 2021.
  * [12] D. Zahrádka, “Optimization-based control of the f1/10 autonomous racing car,” 2020.
  * [13] W. Hess et al., “Real-time loop closure in 2d lidar slam,” in ICRA, IEEE, 2016.
  * [14] S. Macenski and I. Jambrecic, “Slam toolbox: Slam for the dynamic world,” Journal of Open Source Software, vol. 6, no. 61, p. 2783, 2021.