# Magnetic flux plays an important role during a BHXRB outburst in radiative
2T GRMHD simulations

M.T.P. Liska Center for Relativistic Astrophysics, Georgia Institute of
Technology, Howey Physics Bldg, 837 State St NW, Atlanta, GA 30332, USA
Institute for Theory and Computation, Harvard University, 60 Garden Street,
Cambridge, MA 02138, USA  N. Kaaz Center for Interdisciplinary Exploration &
Research in Astrophysics (CIERA), Physics & Astronomy, Northwestern
University, Evanston, IL 60202, USA  K. Chatterjee Black Hole Initiative at
Harvard University, 20 Garden Street, Cambridge, MA 02138, USA  Razieh Emami
Institute for Theory and Computation, Harvard University, 60 Garden Street,
Cambridge, MA 02138, USA  G. Musoke Canadian Institute for Theoretical
Astrophysics, 60 St George St, Toronto, ON M5S 3H8, Canada

###### Abstract

Black hole (BH) X-ray binaries cycle through different spectral states of
accretion over the course of months to years. Although fluctuations in the BH
mass accretion rate are generally recognized as the most important component
of state transitions, it is becoming increasingly evident that magnetic fields
play a similarly important role. In this article, we present the first
radiative two-temperature (2T) general relativistic magnetohydrodynamics
(GRMHD) simulations in which an accretion disk transitions from a quiescent
state at an accretion rate of $\dot{M}\sim 10^{-10}\dot{M}_{\rm Edd}$ to a
hard-intermediate state at an accretion rate of $\dot{M}\sim
10^{-2}\dot{M}_{\rm Edd}$. This huge parameter space in mass accretion rate is
bridged by artificially rescaling the gas density scale of the simulations. We
present two jetted BH models with varying degrees of magnetic flux saturation.
We demonstrate that in ‘Standard and Normal Evolution’ models, which are
unsaturated with magnetic flux, the hot torus collapses into a thin and cold
accretion disk when $\dot{M}\gtrsim 5\times 10^{-3}\dot{M}_{\rm Edd}$. On the
other hand, in ‘Magnetically Arrested Disk’ models, which are fully saturated
with vertical magnetic flux, the plasma remains mostly hot with substructures
that condense into cold clumps of gas when $\dot{M}\gtrsim 1\times
10^{-2}\dot{M}_{\rm Edd}$. This suggests that the spectral signatures observed
during state transitions are closely tied to the level of magnetic flux
saturation.

##  1 Introduction

Black hole X-Ray binaries (BHXRBs) are observed in various spectral state of
accretion (e.g Esin et al. 1997; Remillard & McClintock 2006; Fender et al.
2004; Belloni 2010; Belloni & Motta 2016). Most of the time, BHXRBs are
quiescent, with low luminosities and a hard X-Ray spectrum. However, according
to our current understanding, thermal-viscous instabilities (Lightman &
Eardley, 1974; Shakura & Sunyaev, 1976) in the outer accretion disk drive
periodic outbursts during which the gas supply falling into the black hole
increases by several orders of magnitude (e.g. Lasota 2001). When BHXRBs go
into outburst, their luminosity increases by several orders of magnitude
before and their spectrum changes shape. Similar cycles might exist in AGN
(e.g. Noda & Done 2018), albeit they occur on much longer timescales and can
thus not be observed. It is also well-known that the (quasi-periodic)
variability and jet production are tightly linked to the spectral state of a
black hole (e.g. Ingram & Motta 2019). However, due the complexity of the
physics involved, such as radiation-matter interaction and thermal decoupling
between ions and electrons, our ability to model outbursts is extremely
limited.

Spectral states are classified by two parameters (e.g. Kalemci et al. 2022).
The first parameter is the luminosity with respect to the Eddington limit
($L_{\rm Edd}=\frac{4\pi GM_{\rm Bh}c}{\kappa_{\rm es}}\sim 1.3\times
10^{38}\frac{M_{\rm Bh}}{M_{\odot}}$erg/s, where $\kappa_{\rm es}$ is the
opacity due to electron scattering), which represents the point at which
radiation pressure overcomes gravity. The second parameter is the hardness of
the spectrum, which is usually defined as the luminosity ratio between the
hard and soft X-ray bands, and quantifies how much the spectrum deviates from
a thermal blackbody. During an outburst, the BHXRB first transitions from the
quiescent state to the hard-intermediate state, during which the luminosity
rises rapidly and the spectrum remains dominated by hard Comptonized emission.
Then, the source transitions to the high-soft state characterised by a soft
spectrum with no jet. After spending (typically) weeks to months in the high-
soft state, the source transitions back through several (jetted) intermediate
states to the quiescent state. This transition from the high-soft state to the
quiescent state occurs at a much lower luminosity than the other way around,
for reasons that are not well understood (e.g. Begelman & Armitage 2014; Scepi
et al. 2021).

Most general relativistic magnetohydrodynamic (GRMHD) simulations to date
address accretion in the quiescent state. While BHXRBs indeed spend most of
their time in the quiescent state, they accrete most of their gas (and hence
grow most rapidly) in the hard-intermediate and high-soft states (e.g. Fabian
2012). However, simulating accretion disks in these luminous states is
numerically challenging due to the presence of dynamically important radiation
fields and thermal decoupling between ions and electrons. Presently, only a
handful of GRMHD codes are able to model radiation (e.g. Sadowski et al. 2013;
McKinney et al. 2013; Fragile et al. 2014; Ryan et al. 2017; White et al.
2023). In addition, since radiative cooling makes such accretion disks
thinner, one needs a much higher resolution to resolve them. For example, to
resolve on a spherical grid a disk that is two times thinner without static or
adaptive mesh refinement requires a factor $32$ more computational time. These
factors make such simulations extremely expensive and complex. Due to recent
algorithmic and computational advances, radiative GRMHD simulations of
accretion disks accreting above a few percent of the Eddington limit (i.e.,
very thin disks) came within the realm of possibility (e.g. Ohsuga & Mineshige
2011; Mishra et al. 2016; Morales Teixeira et al. 2018; Fragile et al. 2018;
Lančová et al. 2019; Mishra et al. 2020, 2022; Liska et al. 2022a, 2023).
These recent advances supplement earlier work that attempted to tackle the
physics driving accretion in the luminous states using an ad hoc cooling
function in place of first-principles radiation (e.g. Noble et al. 2009; Avara
et al. 2016; Hogg & Reynolds 2017, 2018; Scepi et al. 2023; Nemmen et al.
2023; Bollimpalli et al. 2023).

Recently, first-of-their-kind radiative GRMHD simulations of accretion disks
accreting at $L\sim 0.35L_{\rm Edd}$ demonstrated that in systems where no
vertical magnetic flux is present, a thin and cold accretion disk forms,
possibly explaining the high-soft state (Liska et al., 2022a). However, in the
presence of dynamically important large scale vertical magnetic flux (e.g.
‘MADs’, Narayan et al., 2003; Tchekhovskoy et al., 2011)), the accretion disk
is truncated and decouples within $r\sim 20r_{g}$ into a two-phase plasma of
cold and dense clumps surrounded by hot and dilute gas (Liska et al., 2022a).
The presence of cold plasma down to the innermost stable circular orbit (ISCO)
provides an interesting explanation for the observed relativistic broadened
iron-reflection lines in the hard-state (e.g. Reis et al. 2010), which is
thought to only feature hot gas unable to produce such lines. In between these
two extreme regimes, where vertical magnetic flux is present but does not
saturate the disk, Lančová et al. (2019) demonstrated that a hot plasma with
both inflowing and outflowing components sandwiches a thin accretion disk.
Such puffy disk models can potentially describe BHXRBs in the intermediate
spectral state, which launch relativistic jets but show no clear evidence of
significant disk truncation (e.g. Kara et al. 2019).

However, none of this work addresses how and at which accretion rates these
high-luminosity accretion states form and what role magnetic fields play in
that process. In this work we present the first radiative two-temperature
GRMHD simulations spanning 8 orders of magnitude in mass accretion rate. These
simulations demonstrate a transition from a hot torus in the quiescent state
to either a magnetically truncated (e.g. Liska et al. 2022a) or puffy
accretion disk (e.g. Lančová et al. 2019) in the (hard-) intermediate state
depending on the amount of magnetic flux saturation. In Section 3 we describe
our radiative GRMHD code and numerical setup, before presenting our results in
Section 4 and concluding in Section 5.

##  2 Numerical Setup

To model the rise from quiescence to the hard-intermediate state we use the
GPU-accelerated GRMHD code H-AMR (Liska et al., 2018, 2022b). H-AMR evolves
the radiative two-temperature GRMHD equations (e.g. Sadowski et al. 2013;
Sa̧dowski et al. 2017) on a spherical grid. Similar to Liska et al. 2022a, we
model the radiation field as a second fluid using the M1 approximation and, in
addition, also evolve the photon number density to get a better estimate for
the radiation temperature (Sa̧dowski & Narayan, 2015). Radiative processes
such as Brehmstrahlung, Synchrotron, bound-free, iron line emission, and
scattering (including Comptonization) are included assuming a $M_{\rm
Bh}=10M_{\odot}$ black hole with solar abundances ($X=0.70$, $Y=0.28$,
$Z=0.02$). The associated energy-averaged grey opacities are provided in
McKinney et al. (2017) (equations C16, D7 and E1).

At each timestep, the dissipation rate is calculated by subtracting the
internal energy provided by the entropy equation from the internal energy
provided by the energy equation (e.g. Ressler et al. 2015). Subsequently, the
total energy dissipation is divided as a source term between the electron and
ions based on a reconnection heating model (Rowan et al., 2017). This deposits
a fraction $\delta_{e}\lesssim 0.5$ of the dissipation into the electrons,
which varies between $\delta_{e}\sim 0.2$ in less magnetized regions to
$\delta_{e}\sim 0.5$ in highly magnetized regions. Coulomb collisions
(Stepney, 1983) are taken into account through an implicit source term (e.g.
Sa̧dowski et al. 2017). To avoid the jet funnel becoming devoid of gas and to
keep the GRMHD scheme stable, we floor the density in the drift frame of the
jet (Ressler et al., 2015) such that the ratio of the density and magnetic
pressure $\frac{p_{b}}{\rho}\lesssim 12.5$.

We use a spherical grid in the Kerr-Schild foliation with coordinates
$x^{1}=\log(r)$, $x^{2}=\theta$, and $x^{3}=\varphi$ with a resolution of
$N_{r}\times N_{\theta}\times N_{\varphi}=420\times 192\times 192$ for our
SANE model and $N_{r}\times N_{\theta}\times N_{\varphi}=560\times 192\times
192$ for our MAD model. This adequately resolves the fastest growing MRI-
wavelength by $\gtrsim 16$ cells in all 3 dimensions. We place the outer
boundary at $R_{\rm out}=10^{3}r_{g}$ for our SANE model, and at $R_{\rm
out}=10^{4}r_{g}$ for our MAD model. We also maintain at least 5 cells within
the event horizon such that the inner boundary is causally disconnected from
the rest of the computational domain. We speed up the simulations
approximately $3$-fold by introducing $4$ levels of local adaptive
timestepping (Liska et al., 2022b). To prevent cell squeezing around the polar
axis from slowing down our simulations (e.g. Courant et al. 1928) we use $4$
levels of static mesh derefinement (Liska et al., 2018, 2022b) to reduce the
$\varphi-$resolution to $N_{\varphi}=[96,48,24,12]$ within
$\theta\lesssim[30^{\circ},15^{\circ},7.5^{\circ},3.75^{\circ}]$ from each
pole. This maintains a cell aspect ratio of roughly $|\Delta
r|:|\Delta\theta|:|\Delta\varphi|\sim 1:1:2$ throughout the grid, which is
sufficient to capture the 3-dimensional nature of the turbulence.

Figure 1: Panel A: The event horizon mass accretion rate $\dot{M}$ closely
follows the target mass accretion rate $\dot{M}_{\rm target}$ (red) for both
the SANE (black) and MAD (blue) models. Panel B:  The normalized magnetic flux
$\phi$ maintains saturation value in the MAD model and stays a factor $\gtrsim
2.0$ below saturation in the SANE model. Figure 2: The SANE (upper panels) and
MAD (lower panels) models at 3 different accretion rates. The left hemisphere
illustrates the electron temperature ($T_{e}$) while the right hemisphere
illustrates the density ($\rho$). The disk-jet boundary ($b^{2}/\rho=1$) is
demarcated by white line and the last scattering surface ($\tau_{es}=1$) is
demarcated by a magenta line. The inset in the left hemisphere gives the mass
accretion rate and luminosity in Eddington units. See our Youtube playlist for
a full animation of both figures that also includes the ion temperature. At
very low accretion rates (left panels), the electron temperature is determined
by the heating rate and adiabatic evolution of the electrons. At intermediate
accretion rates (middle panels), the electrons cool efficiently but there is
no noticeable change in the disk structure. At accretion rates of
$\dot{M}\gtrsim 10^{-2}\dot{M}_{\rm Edd}$ the torus collapses into a thin
accretion disk (SANE) sandwiched by a hot plasma or forms a magnetically
truncated accretion disk (MAD).

##  3 Physical Setup

To understand the effects of magnetic flux saturation on the transition from
the quiescent to the hard-intermediate state, we include two models in the
SANE (‘Standard and Normal Evolution’, Narayan & Yi, 1994) and MAD
(‘Magnetically Arrested Disk’, Narayan et al., 2003) regimes. We assume a
rapidly spinning black hole with spin parameter $a=0.9375$. Our SANE model
(XRB SANE) features a standard Fishbone & Moncrief torus (Fishbone & Moncrief,
1976) with inner radius $r_{\rm in}=6\,r_{g}$, radius of maximum pressure at
$r_{\rm max}=12\,r_{g}$, and outer radius $r_{\rm out}\sim 50r_{g}$. Our MAD
model (XRB MAD), on the other hand, features a much larger torus with $r_{\rm
in}=20\,r_{g}$ and $r_{\rm max}=41r_{g}$ whose outer edge lies at $r_{\rm
out}\sim 800\,r_{g}$. These torii are pretty standard choices in the GRMHD
community (e.g. Porth et al. 2019; Chatterjee et al. 2023). We thread the SANE
model with magnetic vector potential $A_{\Phi}\propto\rho-0.2$ and the MAD
model with magnetic vector potential $A_{\Phi}\propto\rho
r^{3}\sin^{3}\left(\theta\right)\exp\left(-\frac{r}{400\,r_{g}}\right)-0.2$.
Here $\rho$ is the gas density. In both cases this produces a field loop that
is approximately the size of the torus. Because the field loop in our MAD
torus is much larger than in our SANE torus, only the MAD torus contains
enough magnetic flux to get saturated and become MAD. In both cases, we
normalize the resulting magnetic field such that $\beta^{\rm max}=p^{\rm
max}_{\rm gas}/p^{\rm max}_{b}=100$ where $p^{\rm max}_{\rm gas}$ and $p^{\rm
max}_{b}$ are the maximum gas and magnetic pressure in the torus. For the
purpose of calculating the initial torus solution we set the adiabatic index
$\gamma=5/3$ for our SANE model and $\gamma=13/9$ for our MAD model. We
subsequently distribute, according to our heating prescription involving a
magnetic reconnection model (Rowan et al., 2017), the total pressure between
the ions and electrons before we self-consistently evolve their entropy and
adiabatic indices (e.g. Sa̧dowski et al. 2017).

To make GRMHD simulations of BHXRB outbursts feasible, we artificially shorten
the relevant timescales by introducing a rescaling method that, as a function
of time, sets the accretion rate to a predetermined value. However, before we
apply this method, we first run the simulation for $t=10^{4}r_{g}/c$ in two-
temperature non-radiative GRMHD to get the accretion disk into a quasi-steady
state. We subsequently restart the simulation in radiative two-temperature
GRMHD and re-normalize the density ($\rho$) every full timestep with a factor
$\zeta$ such that the running average of the black hole mass accretion rate,
$\langle\dot{M}\rangle$,

$\langle\dot{M}\rangle=\langle\int-\sqrt{-g}\rho u^{r}d\theta
d\varphi\rangle|_{t-10^{3}r_{g}/c}^{t},$ (1)

is scaled to a time-dependent ‘target’ mass accretion rate,

$\dot{M}_{\rm Target}=10^{-10}\times 2^{\frac{t-10^{4}}{10^{4}}}\dot{M}_{\rm
Edd},$ (2)

via the rescaling factor,

$\zeta=\dot{M}_{\rm Target}/\langle\dot{M}\rangle|_{r=5r_{g}}$ (3)

Here $\dot{M}_{\rm Edd}=\frac{L_{\rm Edd}}{\eta_{\rm NT}c^{2}}$ is the
Eddington accretion rate and $\eta_{\rm NT}=0.178$ the Novikov & Thorne (1973)
radiative efficiency. We also rescale the internal energy density ($u_{g}$),
radiation energy density ($E_{\rm rad}$) and magnetic energy density ($b^{2}$)
with the same prefactor $\zeta$ as the density. This leads to a doubling of
the black hole mass accretion rate every $t=10^{4}r_{g}/c$. Note that this
approach automatically increases the total amount of event horizon magnetic
flux ($\Phi=\frac{\sqrt{4\pi}}{2}\int|B^{r}|d\theta d\varphi$) by a factor
$\sqrt{\zeta}$, such that the normalized magnetic flux
($\phi=\frac{\Phi|_{r=r_{\rm Bh}}}{\sqrt{\langle\dot{M}\rangle}|_{r=5r_{g}}}$)
remains constant. When $\phi\sim 40-50$ we expect that the disk turns MAD and
flux bundles get ejected from the black hole (e.g. Tchekhovskoy et al. 2011;
McKinney et al. 2012). We achieve inflow-outflow equilibrium over the mass
accretion rate doubling time ($\Delta t=10^{4}r_{g}/c$) up to approximately
$r\sim 15r_{g}$ for our SANE model and $r\sim 30r_{g}$ for our MAD model.

##  4 Results

We evolve both models for $t\sim 260,000-280,000\,r_{g}/c$, during which the
targeted mass accretion rate increases by $8$ orders of magnitude. As
illustrated in Figure 1a, the black hole mass accretion rate closely follows
the targeted mass accretion rate for both models. However, as illustrated in
Figure 1b, while the normalized magnetic flux threading the BH event horizon
stays constant in our MAD model, it increases by a factor of $\sim 3$ in our
SANE model. This is because a significant fraction of the initial gas
reservoir accretes or gets ejected in the form of winds, leading to a relative
larger increase in the dimensionless magnetic flux compared to the mass
accretion rate. The rapid variability of the magnetic flux observed in our MAD
model is a well-known characteristic of MADs (e.g. Tchekhovskoy et al. 2011;
McKinney et al. 2012) caused by flux bundles being ejected from the black hole
event horizon through magnetic reconnection (e.g. Ripperda et al. 2021).

The contour plots of density and electron temperature in Figure 2 illustrate 3
different stages of the artificially-induced state transition. An accompanying
animation also illustrating the ion temperature is included in the
supplementary materials and on our Youtube playlist. In the first stage
($\dot{M}\lesssim 10^{-6}\dot{M}_{\rm Edd}$), the radiative efficiency is low
and radiative cooling plays a negligible role. The ions in the plasma are
significantly hotter than the electrons because our heating prescription
typically injects only a fraction $\delta_{e}\sim 0.2-0.4$ of the dissipative
heating into the electrons (and a fraction $\delta_{i}\sim 0.6-0.8$ into the
ions). In the second stage ($\dot{M}\gtrsim 10^{-6}\dot{M}_{\rm Edd}$),
radiative cooling of the electrons becomes efficient leading to a drop in
electron temperature but no other structural change (see also Chatterjee et
al. 2023). In the third stage ($\dot{M}\gtrsim 10^{-2}\dot{M}_{\rm Edd}$),
Coulomb collisions become efficient. This allows the ions to cool by
transferring their energy to the radiation-emitting electrons and, eventually,
leads to a rapid collapse of the hot torus.

In our SANE model this collapse results in a geometrically thin accretion disk
surrounded by hot magnetic pressure supported gas outside of $r\gtrsim
3\,r_{g}$. Thus, the disk is only truncated very close to the black hole. The
production of hot non-thermal electrons within the ISCO was predicted by
Hankla et al. (2022). Interestingly, this hot coronal gas rather than the thin
accretion disk seems to be responsible for the majority of $\dot{M}$ (Fig. 3).
This is similar to the puffy accretion disks first presented in other
radiative GRMHD simulations (Lančová et al., 2019) and in pure MHD simulations
of weakly to moderately magnetized disks (Jacquemin-Ide et al., 2021). On the
other hand, in our MAD model the hot torus transitions into a two-phase medium
with cold optically thick patches of gas surrounded by hot, optically thin,
plasma. These cold patches of gas are visible for $20\,r_{g}\lesssim r\lesssim
100\,r_{g}$ and do not reach the event horizon. Since this work was performed
at a rather low resolution, we were forced to stop these simulations after the
cold plasma became under-resolved. Follow-up simulations featuring a much
higher resolution will be necessary to resolve this cold and slender plasma
evolves as we keep increasing the accretion rate.

Figure 3: A contourplot of density with velocity streamlines in black and the
last scattering surface in magenta for model XRB SANE. Similar to the puffy
accretion disk models presented in Lančová et al. (2019) the majority of gas
accretion seems to be driven by inflows outside of the disk’s midplane.

Nevertheless, these findings diverge from the magnetically truncated accretion
disk models detailed in Liska et al. (2022a), where a slender disk was
truncated at $r\sim 20r_{g}$ and cold patches of gas reached the event
horizon. The absence of this thin disk in our simulations may be attributed to
the considerably higher saturation of magnetic flux within our torus,
distinguishing it from the disk in Liska et al. 2022a. This discrepancy could
feasibly result in a significantly larger truncation radius. Consequently, if
the truncation radius in our MAD model lies much further out, it is plausible
that our simulation’s duration is insufficient to capture the formation of a
thin accretion disk. A larger truncation radius (assuming the cold clumps of
gas were sufficiently resolved, which might not be true) might consequently
also explain why no cold plasma reaches the event horizon. Namely, as proposed
in Liska et al. (2022a), magnetic reconnection can potentially evaporate the
cold clumps of gas before they reach the event horizon. This is less likely to
happen if the magnetic truncation radius moves further in and, hence, the cold
clumps have less time to evaporate.

In Figure 4 we plot the time evolution of the bolometric luminosity (panels a
and b), density scale height (panels c and d), and outflow efficiencies (panel
c and d) as function of the mass accretion rate. While the luminosity
increases from $L=10^{-15}L_{\rm Edd}$ to $L=10^{-2}L_{\rm Edd}$, the
radiative efficiency increases by $3-5$ orders of magnitude. Similar to
results presented in the radiative GRMHD simulations of Ryan et al. (2017) and
Dexter et al. (2021), the MAD model is significantly more radiatively
efficient, especially at low accretion rates. This is caused by more efficient
Synchrotron cooling in the highly magnetized gas of a MAD. Around
$\dot{M}=5\times 10^{-3}\dot{M}_{\rm Edd}$ the SANE model collapses into a
thin accretion disk, and we observe a rapid order-of-magnitude rise in the
radiative efficiency to the NT73 (Novikov & Thorne, 1973) limit of
$\eta_{rad}\sim\eta_{NT}\sim 0.18$. Here
$\eta_{rad}=\frac{(\int\sqrt{-g}R^{r}_{t}d\theta
d\varphi)|_{r=200r_{g}}}{\langle\dot{M}\rangle|_{r=5r_{g}}}$ with
$R^{\mu}_{\nu}$ being the radiation stress-energy tensor. This collapse
manifests itself as a rapid decrease in the density scale height of the disk
($\frac{h}{r}=\langle\theta-\pi/2\rangle|_{\rho}$ with $|_{\rho}$ denoting a
density-weighted average). On the other hand, in our MAD model the radiative
efficiency asymptotes to $\eta_{rad}\sim 1.2\eta_{NT}$. This has been observed
in other radiative MADs (e.g. Curd & Narayan 2023) and could, pending future
analysis, potentially be explained by the presence of a dynamically important
magnetic field that injects energy into the accreting gas, which is not
accounted for in Novikov & Thorne (1973). In addition, there is only a
marginal factor $\sim 2$ decrease in the disk scale height after the formation
of cold plasma, because magnetic pressure is able to stabilize the accretion
disk against runaway thermal collapse (e.g Sa̧dowski 2016; Jiang et al. 2019).
The total wind ($\eta_{wind}=\frac{\int\sqrt{-g}T^{r}_{t}d\theta
d\varphi|_{r=5r_{g}}^{(b^{2}/\rho)<1}-\dot{M}|_{r=5r_{g}}}{\langle\dot{M}\rangle|_{r=5r_{g}}}$
with $T^{\mu}_{\nu}$ the stress energy tensor and $\frac{b^{2}}{\rho}=1$ the
wind-jet boundary) and jet ($\eta_{jet}=\frac{\int\sqrt{-g}T^{r}_{t}d\theta
d\varphi|_{r=5r_{g}}^{(b^{2}/\rho)>1}}{\langle\dot{M}\rangle|_{r=5r_{g}}}$)
driven outflow efficiency remain relatively constant throughout the evolution
in our MAD model. However, in our SANE model, the increase in the normalized
magnetic flux causes the jet to become significantly more efficient over time
(e.g. $\eta_{jet}\propto\phi^{2}$).

Figure 4: Panels (a,b): As the mass accretion rate rises the luminosity
increases from $L\sim 10^{-15}-10^{-13}L_{EDD}$ to $L\sim 10^{-2}L_{EDD}$.
This is driven by both an increase in the mass accretion rate and radiative
efficiency. Panels (c,d): The density scaleheight of the disk stays relatively
constant until the accretion rate exceeds $\dot{M}\sim 10^{-3}\dot{M}_{EDD}$
at which point the disk collapses. At this point Coulomb collisions are
efficient and both the ions can cool through the radiation emitting electrons.
Panels (e,f): The jet (blue), wind (black), radiative (purple) and NT73
(dashed green) efficiency as function of $\dot{M}$. Interestingly, the MAD
model maintains a significantly higher radiative efficiency, presumably due to
more efficient Synchrotron emission. Figure 5: Panels (a,b): In both models
the electron temperature $T_{e}$ starts cooling for $\dot{M}\gtrsim
10^{-6}\dot{M}_{EDD}$ due to more efficient radiative cooling. The ion
temperature $T_{i}$ is governed by the ability of the ions to transfer their
energy to the radiation-emitting electrons and doesn’t start cooling until
$\dot{M}\sim 10^{-3}\dot{M}_{EDD}$. Panels (c,d): We compare the radiative
emission ($t_{Em}$), Comptonization ($t_{Compt}$), and Coulomb collision
($t_{Coulomb}$) timescale against the accretion ($t_{acc}$) timescale to
illustrate the importance of these processes at different accretion rates.
Panels (e,f): The adiabatic index for the electrons $\gamma_{e}$ is
relativistic while the adiabatic index for the ions $\gamma_{i}$ is non-
relativistic. This leads to a semi-relativistic gas with adiabatic index
$\gamma_{g}\sim 1.55$ at $\dot{M}\lesssim 10^{-6}\dot{M}_{\rm Edd}$.

To better understand when, during an outburst, certain physical processes
become important, we plot in Figure 5(a,b) the radiative cooling timescale
($t_{rad}=\frac{\int\sqrt{-g}u_{e}d\theta
d\varphi|_{r=10r_{g}}}{\Lambda_{Em}|_{r=10r_{g}}}$ with $\Lambda_{Em}$ the
radiative emission rate and $u_{i,e}$ the ion/electron internal energy),
Comptonization timescale ($t_{Compt}=\frac{\int\sqrt{-g}u_{e}d\theta
d\varphi|_{r=10r_{g}}}{\Lambda_{Compt|_{r=10r_{g}}}}$ with $\Lambda_{Compt}$
the Compton scattering emission rate), Coulomb coupling timescale
($t_{Coulom}=\frac{\int\sqrt{g}\left(u_{e}+u_{i}\right)d\theta
d\varphi|_{r=10r_{g}}}{\Lambda_{Coulomb}|_{r=10r_{g}}}$ with
$\Lambda_{Coulomb}$ the Coulomb coupling rate), and accretion timescale
($t_{Acc}=\frac{\int\sqrt{-g}\rho rd\theta
d\varphi|_{r=10r_{g}}}{\int\sqrt{-g}\rho u^{r}d\theta
d\varphi|_{r=10r_{g}}}$). $\Lambda_{\rm Em}$, $\Lambda_{\rm Compt}$, and
$\Lambda_{\rm Coul}$ are derived from the opacities given in McKinney et al.
(2017) and the Coulomb coupling rate given in Sa̧dowski et al. (2017).
Evidently, the radiative and Compton cooling timescale become similar to the
accretion timescale around $\dot{M}\gtrsim 10^{-6.5}\dot{M}_{\rm Edd}$ and
$\dot{M}\gtrsim 10^{-4}\dot{M}_{\rm Edd}$ respectively. This manifests itself
in figure 5(c,d) as a decrease in the electron temperature
($T_{e}=\frac{\int\sqrt{-g}p_{e}d\theta
d\varphi|_{r=10r_{g}}}{\int\sqrt{-g}\rho d\theta d\varphi|_{r=10r_{g}}}$). The
ion temperature ($T_{i}=\frac{\int\sqrt{-g}p_{i}d\theta
d\varphi|_{r=10r_{g}}}{\int\sqrt{-g}\rho d\theta d\varphi|_{r=10r_{g}}}$) only
drops when the Coulomb coupling timescale becomes comparable to the accretion
timescale around $\dot{M}\sim 5\times 10^{-3}\dot{M}_{\rm Edd}$. Meanwhile the
plasma transitions in figure 5(e,f) from a quasi-relativistic adiabatic index
$\gamma\sim 1.5$ to a non-relativistic $\gamma\sim 5/3$. Even at accretion
rates that are typically associated with radiatively inefficient accretion
($\dot{M}\lesssim 10^{-7}\dot{M}_{\rm Edd}$), future work will need to test if
electron cooling (see also Dibi et al. 2012; Yoon et al. 2020) and/or a self-
consistent adiabatic index can change the spectral signatures compared to
equivalent non-radiative single-temperature GRMHD models (e.g. Mościbrodzka et
al. 2014).

##  5 Discussion and Conclusion

In this article we addressed for the first time the transition from the
quiescent to the hard intermediate state using radiative two-temperature GRMHD
simulations. By rescaling the black hole mass accretion rate across 8 orders
of magnitude, these simulations demonstrated that radiative cooling and
Coulomb coupling become increasingly important and eventually lead to a
transition to a two-phase medium. While the hot torus in SANE models
transitions to a thin accretion disk surrounded by a sandwich-like corona,
reminiscent of a puffy (Lančová et al., 2019) or magnetically elevated
(Begelman & Silk, 2017) disk, the MAD torus transitions to a magnetically
truncated disk and forms cold clumps of gas embedded in a hot corona (see also
Bambic et al. 2023). This is similar to previous work (Liska et al., 2019,
2022a), but in this case an outer thin accretion disk is absent and the cold
plasma does not reach the event horizon. While a thin accretion disk forms
around $\dot{M}\sim 5\times 10^{-3}\dot{M}_{\rm Edd}$ for our SANE model,
which agrees well with Dexter et al. (2021), no cold plasma is observed in our
MAD model before $\dot{M}\gtrsim 1\times 10^{-2}\dot{M}_{\rm Edd}$. As
described in the appendix, we find remarkably similar conclusions for two
analogous models applicable to a $M=6.5\times 10^{9}M_{\odot}$ AGN. Pending
future ray-tracing analysis, we expect the MAD and SANE model have vastly
different spectral and time-variability signatures. We also plan to address
the structure and dynamics of the thin and truncated disks as we keep
increasing the density scale with a dedicated simulation campaign performed at
a much higher resolution.

The goal of this article was to study the transition from the quiescent state
to the hard intermediate state, both of which feature radio jets. Thus, we
have not considered models that do not produce any jets such as models with a
purely toroidal field (e.g. Liska et al. 2022a) and instead only considered a
jetted SANE model with $\phi\sim 5-15$ and a jetted MAD model with $\phi\sim
40-55$. By rescaling both the gas density and magnetic energy density
proportionally to the target mass accretion rate, this work implicitly assumes
that the accretion-rate normalized magnetic flux ($\phi$) remains constant
within a factor $\sim 2$. Conventional thinking might suggest that since the
hard-intermediate state is associated with the most powerful jets, and recent
polarimetric Event Horizon Telescope observations (Event Horizon Telescope
Collaboration et al., 2021) strongly imply that AGN accretion disks in the
quiescent state are MAD, the hard-intermediate state would be MAD as well.
However, the jet power is set by the total amount of magnetic flux threading
the black hole ($P_{\rm jet}\propto\Phi^{2}$), and thus a SANE jet at a much
higher accretion rate can easily outperform a MAD jet at a lower accretion
rate. Thus, an interesting possibility to be explored in future work would
include a model where the magnetic flux does not increase proportional to
$\Phi\propto\sqrt{\dot{M}}$ but is truncated at a maximum value $\Phi\propto
min(\sqrt{\zeta}\Phi_{0},\Phi_{max})$. This would cause the disk to transition
from a MAD disk in the quiescent state to a SANE disk in the hard-intermediate
state where, at least initially, the magnetic pressure is still dynamically
important (e.g. Begelman & Silk 2017; Dexter & Begelman 2019; Lančová et al.
2019). In upcoming work, we will employ ray-tracing calculations to compare
both our existing models and future models featuring a truncated magnetic flux
against multi-wavelength observations, which offer constraints on the
truncation radius and the size/geometry of coronal structures in actual
astrophysical systems (e.g. Ingram & Done 2011; Plant et al. 2014; Fabian et
al. 2014; García et al. 2015; Kara et al. 2019).

There are several theoretical and observational arguments that support this
‘truncated flux’ scenario. First, for systems to remain MAD during a 2-4
orders of magnitude increase in $\dot{M}$, they would need to advect $1-2$
orders of magnitude additional magnetic flux onto the BH (e.g
$\Phi\propto\sqrt{\dot{M}}$ in a MAD). Especially when the outer disk becomes
geometrically thin it is unclear if this is physically possible since
theoretical arguments suggest thin disks might not be able to advect magnetic
flux loops (e.g. Lubow et al. 1994). Second, observations suggests that the
disk truncation radius in the hard intermediate state appears (e.g. Reis et
al. 2010; Kara et al. 2019) to be rather small ($r_{t}\lesssim 5r_{g}$). This
is inconsistent with recent GRMHD simulations which demonstrated that even
when the disk only contained a factor $\sim 1.5$ of excess magnetic flux
(above the MAD limit), this led to a truncation radius $r_{t}\sim 20r_{g}$
(e.g. Liska et al. 2019, 2022a). Third, low-frequency quasi-periodic
oscillations which are ubiquitous in the hard-intermediate state (e.g. Ingram
& Motta 2019), are most likely seeded by a precessing disk which tears of from
a larger non-precessing disk (e.g. Stella & Vietri 1998; Ingram et al. 2009,
2016; Musoke et al. 2022). This has been observed in radiative and non-
radiative GRMHD simulations where a tilted thin accretion disk is threaded by
a purely toroidal magnetic field (e.g. Liska et al. 2022b; Musoke et al. 2022;
Liska et al. 2023) and in similar GRMHD simulations where the accretion disk
is threaded by a below saturation level vertical magnetic field (e.g. Liska et
al. 2021). However, there are no numerical simulations that have shown any
disk tearing or precession where the disk is saturated by vertical magnetic
flux (e.g. Fragile et al. 2023). The main problem is that for a disk to tear
(and precess), the warping of space-time needs to substantially exceed the
viscous torques holding the disk together (e.g. Nixon & King 2012; Nealon et
al. 2015; Doǧan et al. 2018; Doğan & Nixon 2020; Raj et al. 2021). However,
the viscous torques stemming from equipartition strength magnetic fields
within the truncation radius might be too strong for a disk to tear.

While our simulations incorporate the effects of radiation and thermal
decoupling between ions and electrons, they still rely on a rather simplistic
heating prescription for electrons extracted from particle-in-cell models
(Rowan et al., 2017). Since, absent any Coulomb collisions, the cooling rate
in a given magnetic field will be determined by the temperature and density of
the radiation emitting electrons, the radiative efficiency at lower accretion
rates can become sensitive to the used heating prescription (e.g. Chael et al.
2018). For example, in our models roughly a fraction $\delta_{e}\sim 20-40\%$
of the dissipation ends up in the electrons. If this electron heating fraction
would be smaller/bigger, we expect that the radiative efficiency to drop/rise
and the collapse to a two-phase medium to occur later/earlier. Similarly,
other microphysical effects, typically not captured by the ideal MHD
approximation, such as thermal conduction between the corona and disk (e.g.
Meyer & Meyer-Hofmeister 1994; Liu et al. 1999; Meyer-Hofmeister & Meyer 2011;
Cho & Narayan 2022), a non-unity magnetic Prandtl number (e.g. Balbus & Henri
2008), could alter the transition rate to a two-phase medium.

In addition, it was recently demonstrated that the physics driving accretion
in luminous black holes (e.g. with $L\gtrsim 0.01L_{\rm Edd}$), which are
misaligned with the black hole spin axis, is fundamentally different. Namely,
dissipation of orbital energy is driven by nozzle shocks induced by strong
warping (Kaaz et al., 2022; Liska et al., 2023) instead of magneto-rotational
instability (MRI) driven turbulence (e.g. Balbus & Hawley 1991, 1998). These
nozzle shocks form perpendicular to the line of nodes, where the disk’s
midplane intersects the equatorial plane of the black hole, and increase the
radial speed of the gas by $2-3$ orders of magnitude in luminous systems that
are substantially misaligned. This could, at a given accretion rate, lead to a
decrease in the disk’s density, potentially delaying the formation of a thin
accretion disk. We expect to address outbursts of warped accretion disks in
the coming years.

Numerically, this article has also introduced a method to study outbursts by
artificially rescaling the density as a function of time. This solves the
issue that the physical processes in the outer disk that drive such drastic
fluctuations in the mass accretion rate occur over timescales that are too
long to simulate (real outbursts typically take weeks to months, while our
simulations last for $t\sim 10-15s$). Future applications of this method might
include (i) ultra-luminous accretion disks, which decay from super-Eddington
to sub-Eddington accretion rates; (ii) the transition from the hard-
intermediate state to the high-soft state, where the magnetic flux threading
the black hole drops while the accretion rate remains constant; and (iii) the
transition from the high-soft state to the quiescent state, characterised by a
gradual drop in the mass accretion rate.

##  6 Acknowledgements

We thank Sera Markoff, Sasha Tchekhovskoy, and Ramesh Narayan for insightful
discussions. An award of computer time was provided by the Innovative and
Novel Computational Impact on Theory and Experiment (INCITE) and ASCR
Leadership Computing Challenge (ALCC) programs under awards PHY129 and AST178.
This research used resources of the Oak Ridge Leadership Computing Facility,
which is a DOE Office of Science User Facility supported under Contract DE-
AC05-00OR22725. ML was supported by the John Harvard, ITC and NASA Hubble
Fellowship Program fellowships. NK was supported by an NSF Graduate Research
Fellowship. RE was supported by the NASA ATP grant numbers 21-ATP21-0077, NSF
AST-1816420, and HST-GO-16173.001-A. KC was supported by the Black Hole
Initiative (BHI) fellowship program. GM was supported by a Netherlands
Research School for Astronomy (NOVA) Virtual Institute of Accretion (VIA) and
the Canadian Institute for Theoretical Astrophysics (CITA) postdoctoral
fellowships.

We present two radiative two-temperature GRMHD models (M87 SANE and M87 MAD)
where we change the black hole mass from a typical BHXRB of $M_{\rm
Bh}=10M_{\odot}$ to a large AGN such as M87 with $M_{\rm Bh}=6.5\times
10^{9}M_{\odot}$. Figures 6, 7, 8,9 and 10 in the Appendix correspond to
figures 1, 2, 3, 4 and 5 in the main article. Interestingly, the evolution of
our AGN models looks very similar to our BHXRB models. The most striking
difference between BHXRBs and AGN is a slightly lower radiative efficiency at
lower accretion rates, which can be explained by a weaker Synchrotron emission
opacity coefficient (e.g. McKinney et al. 2017) and is reflected in a longer
emission time $t_{Em}$. In addition, after the plasma condenses into a two-
phase medium, the temperature of the cold phase gas in AGN ($T_{e}\sim
10^{5}K$) is much lower than in BHXRBs ($T_{e}\sim 10^{7}K$). This is a well
known fact in the analytic theory of radiatively efficient AGN accretion disks
(e.g. Shakura & Sunyaev 1973; Novikov & Thorne 1973), which are less dense and
hence more radiation pressure dominated than their BHXRB analogues.

Figure 6: Same as figure 1, but for a $M=6.5\times 10^{9}M_{\odot}$ AGN.
Figure 7: Same as figure 2, but for a $M=6.5\times 10^{9}M_{\odot}$ AGN.
Figure 8: Same as figure 3, but for a $M=6.5\times 10^{9}M_{\odot}$ AGN.
Figure 9: Same as figure 4, but for a $M=6.5\times 10^{9}M_{\odot}$ AGN.
Figure 10: Same as figure 5, but for a $M=6.5\times 10^{9}M_{\odot}$ AGN.

## References

  * Avara et al. (2016) Avara, M. J., McKinney, J. C., & Reynolds, C. S. 2016, MNRAS, 462, 636, doi: 10.1093/mnras/stw1643
  * Balbus & Hawley (1991) Balbus, S. A., & Hawley, J. F. 1991, ApJ, 376, 214, doi: 10.1086/170270
  * Balbus & Hawley (1998) —. 1998, Reviews of Modern Physics, 70, 1, doi: 10.1103/RevModPhys.70.1
  * Balbus & Henri (2008) Balbus, S. A., & Henri, P. 2008, ApJ, 674, 408, doi: 10.1086/524838
  * Bambic et al. (2023) Bambic, C. J., Quataert, E., & Kunz, M. W. 2023, arXiv e-prints, arXiv:2304.06067, doi: 10.48550/arXiv.2304.06067
  * Begelman & Armitage (2014) Begelman, M. C., & Armitage, P. J. 2014, ApJ, 782, L18, doi: 10.1088/2041-8205/782/2/L18
  * Begelman & Silk (2017) Begelman, M. C., & Silk, J. 2017, MNRAS, 464, 2311, doi: 10.1093/mnras/stw2533
  * Belloni (2010) Belloni, T. M. 2010, States and Transitions in Black Hole Binaries, ed. T. Belloni, Vol. 794, 53, doi: 10.1007/978-3-540-76937-8_3
  * Belloni & Motta (2016) Belloni, T. M., & Motta, S. E. 2016, in Astrophysics and Space Science Library, Vol. 440, Astrophysics of Black Holes: From Fundamental Aspects to Latest Developments, ed. C. Bambi, 61, doi: 10.1007/978-3-662-52859-4_2
  * Bollimpalli et al. (2023) Bollimpalli, D. A., Fragile, P. C., & Kluźniak, W. 2023, MNRAS, 520, L79, doi: 10.1093/mnrasl/slac155
  * Chael et al. (2018) Chael, A., Rowan, M., Narayan, R., Johnson, M., & Sironi, L. 2018, MNRAS, 478, 5209, doi: 10.1093/mnras/sty1261
  * Chatterjee et al. (2023) Chatterjee, K., Chael, A., Tiede, P., et al. 2023, Galaxies, 11, 38, doi: 10.3390/galaxies11020038
  * Cho & Narayan (2022) Cho, H., & Narayan, R. 2022, ApJ, 932, 97, doi: 10.3847/1538-4357/ac6d5c
  * Courant et al. (1928) Courant, R., Friedrichs, K., & Lewy, H. 1928, Mathematische Annalen, 100, 32, doi: 10.1007/BF01448839
  * Curd & Narayan (2023) Curd, B., & Narayan, R. 2023, MNRAS, 518, 3441, doi: 10.1093/mnras/stac3330
  * Dexter & Begelman (2019) Dexter, J., & Begelman, M. C. 2019, MNRAS, 483, L17, doi: 10.1093/mnrasl/sly213
  * Dexter et al. (2021) Dexter, J., Scepi, N., & Begelman, M. C. 2021, ApJ, 919, L20, doi: 10.3847/2041-8213/ac2608
  * Dibi et al. (2012) Dibi, S., Drappeau, S., Fragile, P. C., Markoff, S., & Dexter, J. 2012, Monthly Notices of the Royal Astronomical Society, 426, 1928, doi: 10.1111/j.1365-2966.2012.21857.x
  * Doğan & Nixon (2020) Doğan, S., & Nixon, C. J. 2020, MNRAS, 495, 1148, doi: 10.1093/mnras/staa1239
  * Doǧan et al. (2018) Doǧan, S., Nixon, C. J., King, A. R., & Pringle, J. E. 2018, MNRAS, 476, 1519, doi: 10.1093/mnras/sty155
  * Esin et al. (1997) Esin, A. A., McClintock, J. E., & Narayan, R. 1997, ApJ, 489, 865, doi: 10.1086/304829
  * Event Horizon Telescope Collaboration et al. (2021) Event Horizon Telescope Collaboration, Akiyama, K., Algaba, J. C., et al. 2021, ApJ, 910, L12, doi: 10.3847/2041-8213/abe71d
  * Fabian (2012) Fabian, A. C. 2012, Annual Review of Astronomy and Astrophysics, 50, 455, doi: 10.1146/annurev-astro-081811-125521
  * Fabian et al. (2014) Fabian, A. C., Parker, M. L., Wilkins, D. R., et al. 2014, MNRAS, 439, 2307, doi: 10.1093/mnras/stu045
  * Fender et al. (2004) Fender, R. P., Belloni, T. M., & Gallo, E. 2004, MNRAS, 355, 1105, doi: 10.1111/j.1365-2966.2004.08384.x
  * Fishbone & Moncrief (1976) Fishbone, L. G., & Moncrief, V. 1976, ApJ, 207, 962
  * Fragile et al. (2023) Fragile, P. C., Chatterjee, K., Ingram, A., & Middleton, M. 2023, MNRAS, 525, L82, doi: 10.1093/mnrasl/slad099
  * Fragile et al. (2018) Fragile, P. C., Etheridge, S. M., Anninos, P., Mishra, B., & Kluźniak, W. 2018, ApJ, 857, 1, doi: 10.3847/1538-4357/aab788
  * Fragile et al. (2014) Fragile, P. C., Olejar, A., & Anninos, P. 2014, ApJ, 796, 22, doi: 10.1088/0004-637X/796/1/22
  * García et al. (2015) García, J. A., Steiner, J. F., McClintock, J. E., et al. 2015, ApJ, 813, 84, doi: 10.1088/0004-637X/813/2/84
  * Hankla et al. (2022) Hankla, A. M., Scepi, N., & Dexter, J. 2022, MNRAS, 515, 775, doi: 10.1093/mnras/stac1785
  * Hogg & Reynolds (2017) Hogg, J. D., & Reynolds, C. S. 2017, ApJ, 843, 80, doi: 10.3847/1538-4357/aa774b
  * Hogg & Reynolds (2018) —. 2018, ApJ, 854, 6, doi: 10.3847/1538-4357/aaa6c6
  * Ingram & Done (2011) Ingram, A., & Done, C. 2011, MNRAS, 415, 2323, doi: 10.1111/j.1365-2966.2011.18860.x
  * Ingram et al. (2009) Ingram, A., Done, C., & Fragile, P. C. 2009, MNRAS, 397, L101, doi: 10.1111/j.1745-3933.2009.00693.x
  * Ingram et al. (2016) Ingram, A., van der Klis, M., Middleton, M., et al. 2016, MNRAS, 461, 1967, doi: 10.1093/mnras/stw1245
  * Ingram & Motta (2019) Ingram, A. R., & Motta, S. E. 2019, New A Rev., 85, 101524, doi: 10.1016/j.newar.2020.101524
  * Jacquemin-Ide et al. (2021) Jacquemin-Ide, J., Lesur, G., & Ferreira, J. 2021, A&A, 647, A192, doi: 10.1051/0004-6361/202039322
  * Jiang et al. (2019) Jiang, Y.-F., Blaes, O., Stone, J. M., & Davis, S. W. 2019, ApJ, 885, 144, doi: 10.3847/1538-4357/ab4a00
  * Kaaz et al. (2022) Kaaz, N., Liska, M. T. P., Jacquemin-Ide, J., et al. 2022, arXiv e-prints, arXiv:2210.10053.  https://arxiv.org/abs/2210.10053
  * Kalemci et al. (2022) Kalemci, E., Kara, E., & Tomsick, J. A. 2022, in Handbook of X-ray and Gamma-ray Astrophysics, 9, doi: 10.1007/978-981-16-4544-0_100-1
  * Kara et al. (2019) Kara, E., Steiner, J. F., Fabian, A. C., et al. 2019, Nature, 565, 198, doi: 10.1038/s41586-018-0803-x
  * Lančová et al. (2019) Lančová, D., Abarca, D., Kluźniak, W., et al. 2019, ApJ, 884, L37, doi: 10.3847/2041-8213/ab48f5
  * Lasota (2001) Lasota, J.-P. 2001, New A Rev., 45, 449, doi: 10.1016/S1387-6473(01)00112-9
  * Lightman & Eardley (1974) Lightman, A. P., & Eardley, D. M. 1974, The Astrophysical Journal, 187, L1, doi: 10.1086/181377
  * Liska et al. (2018) Liska, M., Hesp, C., Tchekhovskoy, A., et al. 2018, MNRAS, 474, L81, doi: 10.1093/mnrasl/slx174
  * Liska et al. (2021) —. 2021, MNRAS, 507, 983, doi: 10.1093/mnras/staa099
  * Liska et al. (2019) Liska, M., Tchekhovskoy, A., Ingram, A., & van der Klis, M. 2019, MNRAS, 487, 550, doi: 10.1093/mnras/stz834
  * Liska et al. (2023) Liska, M. T. P., Kaaz, N., Musoke, G., Tchekhovskoy, A., & Porth, O. 2023, ApJ, 944, L48, doi: 10.3847/2041-8213/acb6f4
  * Liska et al. (2022a) Liska, M. T. P., Musoke, G., Tchekhovskoy, A., Porth, O., & Beloborodov, A. M. 2022a, ApJ, 935, L1, doi: 10.3847/2041-8213/ac84db
  * Liska et al. (2022b) Liska, M. T. P., Chatterjee, K., Issa, D., et al. 2022b, ApJS, 263, 26, doi: 10.3847/1538-4365/ac9966
  * Liu et al. (1999) Liu, B. F., Yuan, W., Meyer, F., Meyer-Hofmeister, E., & Xie, G. Z. 1999, ApJ, 527, L17, doi: 10.1086/312383
  * Lubow et al. (1994) Lubow, S. H., Papaloizou, J. C. B., & Pringle, J. E. 1994, Monthly Notices of the Royal Astronomical Society, 267, 235, doi: 10.1093/mnras/267.2.235
  * McKinney et al. (2017) McKinney, J. C., Chluba, J., Wielgus, M., Narayan, R., & Sadowski, A. 2017, MNRAS, 467, 2241, doi: 10.1093/mnras/stx227
  * McKinney et al. (2012) McKinney, J. C., Tchekhovskoy, A., & Blandford, R. D. 2012, MNRAS, 423, 3083, doi: 10.1111/j.1365-2966.2012.21074.x
  * McKinney et al. (2013) —. 2013, Science, 339, 49, doi: 10.1126/science.1230811
  * Meyer & Meyer-Hofmeister (1994) Meyer, F., & Meyer-Hofmeister, E. 1994, A&A, 288, 175
  * Meyer-Hofmeister & Meyer (2011) Meyer-Hofmeister, E., & Meyer, F. 2011, A&A, 527, A127, doi: 10.1051/0004-6361/201015478
  * Mishra et al. (2022) Mishra, B., Fragile, P. C., Anderson, J., et al. 2022, arXiv e-prints, arXiv:2209.03317.  https://arxiv.org/abs/2209.03317
  * Mishra et al. (2016) Mishra, B., Fragile, P. C., Johnson, L. C., & Kluźniak, W. 2016, MNRAS, 463, 3437, doi: 10.1093/mnras/stw2245
  * Mishra et al. (2020) Mishra, B., Kluźniak, W., & Fragile, P. C. 2020, MNRAS, 497, 1066, doi: 10.1093/mnras/staa1848
  * Morales Teixeira et al. (2018) Morales Teixeira, D., Avara, M. J., & McKinney, J. C. 2018, MNRAS, 480, 3547, doi: 10.1093/mnras/sty2044
  * Mościbrodzka et al. (2014) Mościbrodzka, M., et al. 2014, A&A, 570, A7, doi: 10.1051/0004-6361/201424358
  * Musoke et al. (2022) Musoke, G., Liska, M., Porth, O., van der Klis, M., & Ingram, A. 2022, arXiv e-prints, arXiv:2201.03085.  https://arxiv.org/abs/2201.03085
  * Narayan et al. (2003) Narayan, R., Igumenshchev, I. V., & Abramowicz, M. A. 2003, PASJ, 55, L69, doi: 10.1093/pasj/55.6.L69
  * Narayan & Yi (1994) Narayan, R., & Yi, I. 1994, ApJ, 428, L13, doi: 10.1086/187381
  * Nealon et al. (2015) Nealon, R., Price, D. J., & Nixon, C. J. 2015, MNRAS, 448, 1526, doi: 10.1093/mnras/stv014
  * Nemmen et al. (2023) Nemmen, R., Vemado, A., Almeida, I., Garcia, J., & Motta, P. 2023, arXiv e-prints, arXiv:2305.11429, doi: 10.48550/arXiv.2305.11429
  * Nixon & King (2012) Nixon, C. J., & King, A. R. 2012, MNRAS, 421, 1201, doi: 10.1111/j.1365-2966.2011.20377.x
  * Noble et al. (2009) Noble, S. C., Krolik, J. H., & Hawley, J. F. 2009, ApJ, 692, 411, doi: 10.1088/0004-637X/692/1/411
  * Noda & Done (2018) Noda, H., & Done, C. 2018, Monthly Notices of the Royal Astronomical Society, 480, 3898, doi: 10.1093/mnras/sty2032
  * Novikov & Thorne (1973) Novikov, I. D., & Thorne, K. S. 1973, in Black Holes (Les Astres Occlus), ed. C. Dewitt & B. S. Dewitt, 343–450
  * Ohsuga & Mineshige (2011) Ohsuga, K., & Mineshige, S. 2011, ApJ, 736, 2, doi: 10.1088/0004-637X/736/1/2
  * Plant et al. (2014) Plant, D. S., Fender, R. P., Ponti, G., Muñoz-Darias, T., & Coriat, M. 2014, MNRAS, 442, 1767, doi: 10.1093/mnras/stu867
  * Porth et al. (2019) Porth, O., Chatterjee, K., Narayan, R., et al. 2019, ApJS, 243, 26, doi: 10.3847/1538-4365/ab29fd
  * Raj et al. (2021) Raj, A., Nixon, C. J., & Doğan, S. 2021, ApJ, 909, 81, doi: 10.3847/1538-4357/abdc24
  * Reis et al. (2010) Reis, R. C., Fabian, A. C., & Miller, J. M. 2010, MNRAS, 402, 836, doi: 10.1111/j.1365-2966.2009.15976.x
  * Remillard & McClintock (2006) Remillard, R. A., & McClintock, J. E. 2006, ARA&A, 44, 49, doi: 10.1146/annurev.astro.44.051905.092532PDF: http://arjournals.annualreviews.org/doi/pdf/10.1146/annurev.astro.44.051905.092532
  * Ressler et al. (2015) Ressler, S. M., Tchekhovskoy, A., Quataert, E., Chandra, M., & Gammie, C. F. 2015, MNRAS, 454, 1848, doi: 10.1093/mnras/stv2084
  * Ripperda et al. (2021) Ripperda, B., Liska, M., Chatterjee, K., et al. 2021, arXiv:2109.15115 [astro-ph, physics:gr-qc, physics:physics].  https://arxiv.org/abs/2109.15115
  * Rowan et al. (2017) Rowan, M. E., Sironi, L., & Narayan, R. 2017, ApJ, 850, 29, doi: 10.3847/1538-4357/aa9380
  * Ryan et al. (2017) Ryan, B. R., Gammie, C. F., Fromang, S., & Kestener, P. 2017, ApJ, 840, 6, doi: 10.3847/1538-4357/aa6a52
  * Sa̧dowski (2016) Sa̧dowski, A. 2016, MNRAS, 459, 4397, doi: 10.1093/mnras/stw913
  * Sa̧dowski & Narayan (2015) Sa̧dowski, A., & Narayan, R. 2015, MNRAS, 454, 2372, doi: 10.1093/mnras/stv2022
  * Sa̧dowski et al. (2017) Sa̧dowski, A., Wielgus, M., Narayan, R., et al. 2017, MNRAS, 466, 705, doi: 10.1093/mnras/stw3116
  * Scepi et al. (2021) Scepi, N., Begelman, M. C., & Dexter, J. 2021, MNRAS, 502, L50, doi: 10.1093/mnrasl/slab002
  * Scepi et al. (2023) —. 2023, arXiv e-prints, arXiv:2302.10226, doi: 10.48550/arXiv.2302.10226
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337
  * Shakura & Sunyaev (1976) —. 1976, MNRAS, 175, 613, doi: 10.1093/mnras/175.3.613
  * Sadowski et al. (2013) Sadowski, A., Narayan, R., Tchekhovskoy, A., & Zhu, Y. 2013, MNRAS, 429, 3533, doi: 10.1093/mnras/sts632
  * Stella & Vietri (1998) Stella, L., & Vietri, M. 1998, ApJ, 492, L59, doi: 10.1086/311075
  * Stepney (1983) Stepney, S. 1983, MNRAS, 202, 467, doi: 10.1093/mnras/202.2.467
  * Tchekhovskoy et al. (2011) Tchekhovskoy, A., Narayan, R., & McKinney, J. C. 2011, MNRAS, 418, L79, doi: 10.1111/j.1745-3933.2011.01147.x
  * White et al. (2023) White, C. J., Mullen, P. D., Jiang, Y.-F., et al. 2023, ApJ, 949, 103, doi: 10.3847/1538-4357/acc8cf
  * Yoon et al. (2020) Yoon, D., Chatterjee, K., Markoff, S. B., et al. 2020, MNRAS, 499, 3178, doi: 10.1093/mnras/staa3031