# AT2018fyk: Candidate Tidal Disruption Event by a (Super)massive Black Hole
Binary

S. Wen<EMAIL_ADDRESS>National Astronomical Observatories, Chinese Academy of
Sciences, 20A Datun Road, Beijing 100101, China  Department of
Astrophysics/IMAPP, Radboud University, P.O. Box 9010, 6500 GL, Nijmegen, The
Netherlands  P.G. Jonker Department of Astrophysics/IMAPP, Radboud University,
P.O. Box 9010, 6500 GL, Nijmegen, The Netherlands  A.J. Levan Department of
Astrophysics/IMAPP, Radboud University, P.O. Box 9010, 6500 GL, Nijmegen, The
Netherlands  D. Li National Astronomical Observatories, Chinese Academy of
Sciences, 20A Datun Road, Beijing 100101, China  N. C. Stone Racah Institute
of Physics, The Hebrew University, Jerusalem, 91904, Israel  A. I. Zabludoff
University of Arizona and Steward Observatory, 933 N. Cherry Ave., Tucson, AZ
85721  Z. Cao Department of Astrophysics/IMAPP, Radboud University, P.O. Box
9010, 6500 GL, Nijmegen, The Netherlands  SRON, Netherlands Institute for
Space Research, Niels Bohrweg 4, 2333 CA, Leiden, The Netherlands  T. Wevers
Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218,
USA  European Southern Observatory, Alonso de Córdova 3107, Vitacura,
Santiago, Chile  D.R. Pasham Kavli Institute for Astrophysics and Space
Research, Massachusetts Institute of Technology, Cambridge, MA, USA  C. Lewin
Kavli Institute for Astrophysics and Space Research, Massachusetts Institute
of Technology, Cambridge, MA, USA  E. Kara  Kavli Institute for Astrophysics
and Space Research, Massachusetts Institute of Technology, Cambridge, MA, USA

###### Abstract

The tidal disruption event (TDE) AT2018fyk has unusual X-ray, UV, and optical
light curves that decay over the first $\sim$600d, rebrighten, and decay again
around 1200d. We explain this behavior as a one-off TDE associated with a
massive black hole (BH) _binary_. The sharp drop-offs from $t^{-5/3}$ power
laws at around 600d naturally arise when one BH interrupts the debris fallback
onto the other BH. The BH mass $M_{\bullet}$ derived from fitting X-ray
spectra with a slim disk accretion model and, independently, from fitting the
early UV/optical light curves, is smaller by two orders of magnitude than
predicted from the $M_{\bullet}$–$\sigma_{*}$ host galaxy relation, suggesting
that the debris is accreted onto the _secondary_ , with fallback cut off by
the primary. Furthermore, if the rebrightening were associated with the
primary, it should occur around 5000d, not the observed 1200d. The secondary’s
mass and dimensionless spin is $M_{\bullet,{\rm s}}=2.7^{+0.5}_{-1.5}\times
10^{5}M_{\odot}$ and $a_{\bullet,{\rm s}}>0.3$ (X-ray spectral fitting), while
the primary’s mass is $M_{\bullet,{\rm p}}=10^{7.7\pm 0.4}M_{\odot}$
($M_{\bullet}$-$\sigma_{*}$ relation). An intermediate mass BH secondary is
consistent with the observed UV/optical light curve decay, i.e., the
secondary’s outer accretion disk is too faint to produce a detectable emission
floor. The time of the first accretion cutoff constrains the binary separation
to be $(6.7\pm 1.2)\times 10^{-3}~{}{\rm pc}$. X-ray spectral fitting and
timing analysis indicate that the hard X-rays arise from a corona above the
secondary’s disk. The early UV/optical emission, suggesting a super-Eddington
phase for the secondary, possibly originates from shocks arising from debris
circularization.

Tidal disruption (1696), X-ray transient sources (1852), Accretion (14), Black
hole physics (159), galaxies: supermassive black hole binary

##  1 Introduction

When a star passes too close to a supermassive black hole (SMBH), the strong
gravitational field can cause the star to be torn apart in what is known as a
tidal disruption event (TDE; Hills, 1975; Rees, 1988). A significant fraction
of the resulting debris is then accreted onto the SMBH, causing a bright flare
of radiation that can last for months to years. The mass fallback of the
debris typically follows a power law, decaying at late times as $t^{-5/3}$
(Phinney, 1989). The energy dissipated by the fallback debris powers a
luminous flare that in principle enables the measurement of the SMBH mass and
even spin(Mockler et al., 2019; Ryu et al., 2020; Wen et al., 2020; Mummery et
al., 2023b), although existing efforts at parameter estimation all involve
significant assumptions. For this reason, TDEs have evolved into a powerful
tool for investigating the properties of SMBHs.

When two galaxies merge, the two SMBHs will sink to the centre of the newly
formed galaxy due to dynamical friction (Begelman et al., 1980). It is unclear
if the BH binary shrinks all the way down to $\approx$0.1–0.01 pc where
gravitational wave radiation determines the orbital evolution (Milosavljević &
Merritt, 2003). Proposed processes that allow the binary SMBHs to evolve down
to such small separations are: friction caused by the presence of a gaseous
disk, a triaxial stellar distribution, and the presence of a third SMBH of
similar size due to a subsequent merger (Blandford & Rees, 1974; Merritt &
Poon, 2004; Merritt & Milosavljević, 2005; Mayer et al., 2007; Lodato et al.,
2009; Vasiliev et al., 2015). The existence of merging binary SMBHs emitting
gravitational wave radiation is a major prediction underlying ESA-NASA’s LISA
mission (Amaro-Seoane et al., 2023) and Chinese TianQin/Taiji mission (Gong et
al., 2021). The underlying population of low-mass SMBH binaries relevant for
LISA remains unclear, although recent results from Pulsar Timing Arrays do
lend support for merging SMBHs at higher masses (Antoniadis et al., 2023;
Agazie et al., 2023).

When two SMBHs orbit each other at distances of a parsec or less, the TDE rate
is enhanced by a combination of the Kozai-Lidov mechanism (Ivanov et al.,
2005a) and chaotic three-body scatterings (Chen et al., 2011). The TDE rate
can be temporarily enhanced more than a hundred-fold, implying that multiple
TDEs can be seen from the same galaxy on timescales of decades or less (Chen
et al., 2011; Wegg & Nate Bode, 2011). Some theories predict that overall
about $10\%$ of the observed number of TDEs are associated with a supermassive
black hole binary (SMBHB) (Chen et al., 2011; Wegg & Nate Bode, 2011). This
could explain why TDEs prefer host galaxies that experienced a burst of star
formation that ended within the last $\sim$Gyr (Arcavi et al., 2014; French et
al., 2016; Graur et al., 2018); such “E+A” or “post-starburst” galaxies (PSBs)
often have disturbed morphologies indicative of galaxy-galaxy interactions and
mergers.

Aside from indirect evidence related to statistics of host galaxy properties,
the TDE light curves themselves can provide information on the question
whether SMBHs merge. If one SMBH in a binary pair with a separation $\sim{\rm
mpc}$ disrupts a star, the first stages of the TDE proceed as if the other
SMBH were not present. However, due to the time-dependent gravitational field
of the binary SMBH, some of the stellar debris is captured by the binary
companion, truncating the mass fallback and thus the radiation emitted by
accretion onto the SMBH responsible for the tidal disruption. The truncation
time scale is a measure of the orbital period of the binary SMBH (Liu et al.,
2009; Coughlin et al., 2019).

A stellar disruption by the SMBHB _secondary_ is also possible and would be
the exclusive channel for a TDE if the primary SMBH mass exceeds
$10^{8}M_{\odot}$, which would swallow instead of disrupt stars (Hills, 1975).
Although the probability of disruption by the secondary is low (proportional
to the mass ratio) compared to that of the primary (Coughlin et al., 2017),
the event rate associated with the secondary is shown to be significantly
increased by the eccentric Kozai-Lidov (EKL) mechanism (Kozai, 1962; Lidov,
1962; Naoz, 2016a), as reported by Mockler et al. (2023). Such events may be
more distinguishable observationally compared to cases of disruption by the
primary. They not only exhibit enhanced variability in the light curve due to
more pronounced perturbations from the primary, but may also produce huge
discrepancies in SMBH mass estimates, specifically between estimates (i) based
on emission from the TDE accretion disk, and (ii) those from host galaxy
scaling relationships.

The H$\alpha$ and H$\beta$ absorption indicates no ongoing star formation, and
the lack of significant H$\delta$ absorption suggests that the host galaxy of
AT2018fyk does not belong to the “E+A” galaxy class (Wevers et al., 2019).
Furthermore, the study of the emission lines shows that the host galaxy is not
an AGN (Wevers et al., 2023). AT2018fyk was initially classified as a TDE with
$M_{\bullet}=10^{7.7\pm 0.4}M_{\odot}$ from its host spectral properties
(Wevers et al., 2019, 2021) and later reclassified as a repeating partial TDE
(pTDE) due to the observed late-time rebrightening (Wevers et al., 2023). In
this paper, we propose an alternative explanation: that AT2018fyk is
consistent with stellar disruption by a SMBHB, with the secondary playing a
crucial role.

Throughout this paper, we assume a flat concordance cosmological model
($\Lambda$CDM) with parameters $H_{0}=67.4~{}\rm{km~{}s^{-1}~{}Mpc^{-1}}$ and
$\Omega_{m}=0.315$ (Planck Collaboration et al., 2020). The redshift of the
host of AT2018fyk is $z=0.059\pm 0.0005$, corresponding to a luminosity
distance of $D_{L}=274$ Mpc. The paper is organized as follows. We describe
the (archival) optical/UV emission data and X-ray emission data in Section 2.
We show the results from analysing the UV/optical and X-ray emission in
Section 3. In Section 4, we explain AT2018fyk as a SMBHB TDE. Finally, we
summarize our conclusions in Section 5.

##  2 Data and Analysis

###  2.1 X-ray Spectra

Table 1:  A log of the Swift and XMM-Newton observations of AT2018fyk used in this paper. The number of source counts in the photon 0.3–10 keV energy band is given. The spectra are background subtracted. For XMM epochs, the exposure time and counts are shown in order of cameras pn, MOS1, and MOS2, respectively.  Spectrum | date | Exposed time [ks] | Counts
---|---|---|---
Swift 1 | 2018/09/22-2018/10/17 | 16.9 | 446
Swift 2 | 2018/10/18-2018/11/21 | 19.5 | 742
Swift 3 | 2018/11/23-2019/01/08 | 13.6 | 375
Swift 4 | 2019/03/20-2019/05/20 | 25.7 | 2419
Swift 5 | 2019/05/23-2019/08/17 | 20.7 | 1930
Swift 6 | 2019/08/21-2019/10/23 | 9.6 | 432
Swift 7 | 2019/10/27-2020/01/09 | 18.3 | 1267
Swift 8 | 2022/03/25-2022/08/07 | 44.2 | 556
XMM 1 | 2018/12/09 | 23.1/26.8/25.3 | 10372/2219/1791
XMM 2 | 2019/10/27 | 30.2/41.8/42.7 | 22325/6443/6416
XMM 4 | 2022/05/20 | 6.2/12.6/12.5 | 1213/528/637

We reduce Swift X-ray data (Evans et al., 2009) obtained for AT2018fyk over
the period Sept. 22, 2018 to Aug. 7, 2022. We group and average the spectra
into 8 epochs. In addition, we use 3 out of 4 epochs of the XMM-Newton
(EPIC/PN) (Jansen et al., 2001) X-ray data (see Table 1 for a journal of the
data used).

We use the online Swift XRT tool to extract the time-averaged
spectra111https://www.swift.ac.uk/user-objects/. We centre on the source using
the default Single pass method. The Max attempts and Search radius are set at
default values of 10 and 1, respectively. We do not select the Super-soft
source option and use the default all valid Grade range. For the X-ray light
curve, we choose the binning method as “observation”, allowing for upper
limits and Bayesian bins.

For the XMM-Newton data, we run the default SAS v20 (20211130) tools under the
HEASOFT ftools software version 6.30 to extract source spectra and to filter
the data. The EPIC pn detector (Strüder et al., 2001; Turner et al., 2001)
provides the highest count rates but we also use the MOS1 and MOS2 data where
available. We filter the pn data for periods of enhanced background radiation
by deselecting episodes where the 10–12 keV detection rate of pattern 0 events
is higher than 0.4 counts s-1. We filter the MOS1 & MOS2 data for periods of
enhanced background radiation by deselecting episodes where the detection rate
of pattern 0 events at energies $>10$ keV is higher than 0.35 counts s-1. The
effective exposure time for each observation after filtering is given in Table
1. The observations are done with the pn and MOS1&2 detectors in Large Window
(XMM 1) or Full Frame (XMM 2 & 4) mode, respectively. We do not use the XMM 3
observation obtained on April 13, 2020 (MJD 58952.20), as the source was
undetected (Wevers et al., 2021).

###  2.2 NICER Photometry

NICER data was analyzed in the exact same manner as described in Wevers et al.
(2023) taking into account the possible contamination from the foreground
Oxygen contamination. Additional details about the data reduction and
screening thresholds can be found in Pasham et al. 2021.

###  2.3 UV/Optical Photometry

We reduce the UV/optical Telescope (UVOT; Roming et al. 2005) data using the
uvotsource task, extracting fluxes using the standard 7 arcsec radius circular
aperture. The emission in the UV bands has brightened by a factor of 10
compared to the host galaxy brightness, while the brightness of the B and V
filters remains consistent with the inferred host galaxy brightness (Wevers et
al., 2021). We therefore do not include these latter two filters in our
analysis.

###  2.4 Hubble Space Telescope Imaging

The location of AT2018fyk was imaged by the Hubble Space Telescope on August
19, 2021 as part of programme 16239 (PI Foley). At this epoch, observations
were obtained in F225W and F275W with exposure times of 780 and 710 seconds,
respectively. Since only two dither positions were obtained in each filter,
there remain numerous cosmic-ray hits which cannot be fully removed. However,
a source at the location of AT2018fyk is clearly detected in both drizzled
images obtained from the MAST archive (given the small number of dithers, a
re-processing of the data to drizzle to a smaller pixel scale is not
plausible). We measure the source magnitude in these epochs in a small (0.12
arcsecond) aperture centred on the position, obtaining magnitudes of
F225W=24.48 $\pm$ 0.04 and F275W = 22.70 $\pm$ 0.02, when corrected for
encircled energy. The use of a larger (0.32 arcsec) aperture gives magnitudes
approximately 0.4 and 0.2 magnitudes brighter (again after correction for
encircled energy), indicating there may be some extended underlying emission.

##  3 Results

###  3.1 X-ray spectral fitting

Throughout this paper, we fit the spectra using XSPEC version 12.12.0 (Arnaud,
1996), applying Poisson statistics (Cash 1979; C-stat in XSPEC). We quote
parameter errors at a $1\sigma$ $(68.3\%)$ confidence level (CL), using the
method $\rm{Statistic}=\rm{Statistic_{best-fit}}+\Delta\,Cstat$ (Arnaud, 1996)
and assuming $\Delta{\rm Cstat}=1.0$ and $\Delta{\rm Cstat}=2.3$ for single
and two parameter models, respectively (Wen et al., 2021). We use the slim
disk model, slimd (Wen et al., 2022), to model the disk spectrum and the
thermally comptonized continuum model, thcomp (Zdziarski et al., 2020), to
model the hard part of the spectrum. The priors of their parameters are listed
in Table 2. We set a prior on the inclination as $[2^{\circ},60^{\circ}]$,
where $0^{\circ}$ represents a face-on orientation. When the disk is nearly
edge-on ($>60^{\circ}$), a low accretion rate disk spectrum can be well fitted
by a high accretion rate disk model, due to self-shielding by the finite
aspect ratio of the disk (Wen et al., 2020). Edge-on TDEs typically exhibit
additional absorption due to the sight line crossing outflowing material or
intersecting the thick disk (Dai et al., 2018; Wen et al., 2020). However,
there is no evidence for the presence of additional absorption in addition to
the (Galactic) absorption parameterized by $N_{\rm H}=1.15\times
10^{20}\mathrm{cm^{-2}}$ (Wevers et al., 2021). Therefore, we enforce an upper
limit on the disk inclination of $60^{\circ}$. We provide more details about
the slimd model and its capability to constrain disk parameters in Appendix B.

Table 2: Priors on the parameters of the slimd and Thcomp models  slimd |  $M_{\bullet}$ $[M_{\odot}]$ | $a_{\bullet}$ |  $\theta$ [deg] |  $\dot{m}$ [Edd] | Thcomp | $\tau$ |  $T_{e}$ [keV]
---|---|---|---|---|---|---|---
| [$10^{3},9\times 10^{7}$] | [-0.998, 0.998] | [2, 60] | [0.05, 100] |  | [-100, -0.01] | [0.1, 300]

####  3.1.1 Swift 1 Spectrum: Inner Accretion Disk

Wevers et al. (2021) showed that the source starts in a thermally dominated
X-ray spectral state, which later develops into a power-law dominated state
$\sim$100–200 days after the UV/optical peak. When analysing the time-averaged
Swift spectra, we found that the power-law component becomes increasingly
important after the date 2018/10/17. As a result, we designate the time-
averaged spectrum from 2018/09/22 to 2018/10/17 as Swift 1\. We fit the
spectrum with a slim disk model (slimd) (Wen et al., 2022), a thin disk model
(diskbb), and a simplified but general disk model valid in the Wien emission
limit (tdediscspec) (Mummery, 2021; Mummery et al., 2023b), each attenuated by
an absorption model TBAbs with a Galactic column density of NH = 1.15$\times
10^{20}$ cm-2 (the same as in Wevers et al. 2021). Table 3 lists the fit
results. See Figure 1 for the spectrum with the best-fit slim disk model
overlaid. For the slim disk model (slimd), the free parameters are the BH mass
($M_{\bullet}$), the dimensionless BH spin ($a_{\bullet}$), the disk
inclination ($\theta$), and the disk accretion rate $\dot{m}$ (in units of
Eddington accretion rate, $\dot{M}_{\rm Edd}=1.37\times 10^{21}M_{6}~{}{\rm
kg~{}s}^{-1}$, where $M_{6}=M_{\bullet}/(10^{6}M_{\odot})$). We find that the
spectrum can be well fitted with a $\rm{Cstat/d.o.f}\approx 1.0$ for all the
three models. These results show that the early X-ray spectrum originated from
the inner part of an accretion disk. The inferred values for the disk
bolometric luminosity and $M_{\bullet}$, as determined through the slim disk
model fits, are $1.4^{+4.0}_{-0.8}\times 10^{43}$ erg/s and
$0.5^{+1.5}_{-0.3}\times 10^{5}M_{\odot}$ at $1\sigma$ CL, respectively. The
disk bolometric luminosity is consistent with the Eddington luminosity of
$L_{\rm Edd}\approx 2.5\times 10^{43}$ erg/s for a BH mass of
$M_{\bullet}=2\times 10^{5}M_{\odot}$, indicating that the disk is in a super-
Eddington phase (also supported by the high accretion rate of Swift 1 in Table
3). The best-fit unabsorbed X-ray luminosity is $L_{\rm X}=7.5\times 10^{41}$
erg/s (for the band 0.3–10 keV), which implies that the ratio $L_{\rm
X}/L_{\rm d}$ is in the range of $\approx 1-13\%$.

Table 3: The best-fit parameter values and their $1\sigma$ errors from X-ray spectral fits to the Swift 1 time-averaged spectrum. We fit the spectrum with either the fit-function TBAbs $\times$ slimd , TBAbs $\times$ diskbb, or TBAbs $\times$ tdediscspec. For the TBAbs model, we fix the NH to the value 1.15$\times$ 1020 cm-2. Note that, the Eddington accretion rate is defined as $\dot{M}_{\rm Edd}=1.37\times 10^{21}M_{6}~{}{\rm kg~{}s}^{-1}$.  TBAbs $\times$ slimd |  $M_{\bullet}$ $[10^{5}M_{\odot}]$ | $a_{\bullet}$ | $\theta~{}[deg]$ |  $\dot{m}$ [Edd] | $Cstat/d.o.f$
---|---|---|---|---|---
| $0.5^{+1.5}_{-0.3}$ | $0.4^{+0.6}_{-0.8}$ | $54^{+6}_{-52}$ | $62^{+38}_{-61}$ | 87.86/90
TBAbs $\times$ diskbb |  | T [eV] |  | norm | $Cstat/d.o.f$
|  | $150\pm 5$ |  | $200\pm 40$ | 89.01/92
TBAbs $\times$ tdediscspec |  |  $R_{\rm peak}$ [$10^{12}$ cm] |  $T_{\rm peak}$ [$10^{5}$ K] | $\Gamma$ | $Cstat/d.o.f$
|  | $0.09\pm 0.02$ | $14\pm 3$ | $1.5\pm 1.2$ | 87.90/91

####  3.1.2 Evolution of Power-Law Component

We fit the Swift and XMM-Newton X-ray spectra simultaneously, with the fit-
function (TBAbs $\times$ Thcomp(slimd)) consisting of a slim disk model and a
thermally comptonized continuum (Zdziarski et al., 2020), absorbed by a
Galactic column density of NH = 1.15$\times$ 1020 cm-2. For the slim disk
component (slimd), we force $M_{\bullet}$, $a_{\bullet}$ and $\theta$ to have
time-independent values, while allowing for different accretion rates
$\dot{m}$ at each epoch. We do allow the (time-independent) $M_{\bullet}$,
$a_{\bullet}$ and $\theta$ parameters to float freely and hence the values can
be different from the values determined from the fit to the Swift 1 and XMM 2
spectra separately. We will discuss this in more detail in Section 3.1.3 and
Appendix A.

We first fit the spectra by allowing the different epochs to have a different
optical depth and a different electron temperature for the Thcomp model. We
get a good fit with the total ${\rm Cstat/d.o.f.}=6749.3/6713$. We then fit
the spectra by requiring all the epochs to have the same electron temperature.
We also get a good fit with the total ${\rm Cstat/d.o.f.}=6755.7/6722$. The
unfolded data and best-fit model spectrum for each epoch are shown in Figure
1. The best-fit parameters are listed in Table 4. The $\Delta{\rm AIC}=-11.6$
222$\Delta{\rm AIC}=\Delta Cstat+2\Delta K$, where $K$ is the number of free
parameters. Generally speaking, when comparing the best-fit statistics of two
models, $\Delta{\rm AIC}=5$ and $10$ are considered to be strong and very
strong evidence, respectively, against the weaker model. (Akaike, 1974),
obtained when comparing the fit where the electron temperatures are allowed to
be different with the fit where they are tied together, shows that using the
same corona electron temperature for all the epochs is strongly preferred by
the data. We can see from Figure 1 that the contribution of the hard component
to the overall spectrum starts to become larger about one month after the
discovery on 2018 September 8 (MJD = 58369.2), a trend which continues for the
next three months. After epoch 4, the parameters of the best-fit spectral
model stay roughly constant. This can, for instance, be seen from the fitted
optical depth of each epoch. The optical depth grows from 0.1 at Swift 2 to
about 0.3 at Swift 4\. We interpret these results as the formation of a corona
starting about one month after the peak of the UV emission. Its optical depth
increases until it reaches a constant value in about 100 days.

Wevers et al. (2021) have also found that a power-law X-ray component emerges,
with the source evolving from a soft state to a hard state, and they also
attribute this power-law component to originate in a corona. In contrast to
some X-ray binaries (XRBs) (Remillard & McClintock, 2006), where the power-law
strength increases simultaneously as the source luminosity decreases,
AT2018fyk exhibits a disk luminosity that either increases or remains constant
during the strengthening of the power-law component (refer to the lower panel
of Figure 4) 333Note that if we assume the decreasing UV/optical emission
associate directly with the source of the X-ray emission, then the source
luminosity does decrease at the same time the power-law X-ray component become
stronger.. It is noteworthy that a similar power-law X-ray component manifests
in the initial stages of certain TDEs, such as 3XMM J150052.0+015452 (Cao et
al., 2023), AT2019azh (Liu et al., 2022), AT2019ehz (van Velzen et al., 2021),
XMMSL1J0740 (Saxton et al., 2017), AT2020ocn (Pasham et al., 2020), and
AT2021ehb (Yao et al., 2022). Some of the hard components observed in TDEs
like 3XMM J150052.0+015452 and AT2020ocn have been explained by accounting for
the effects of inverse comptonization by a corona, as proposed by Cao et al.
(2023). AT2021ehb exhibits a closely analogous evolution of its hard
component, as illustrated in figure 15 of Yao et al. (2022). The behavior of
AT2021ehb’s hard X-ray spectral component is attributed to the gradual
formation of a magnetically dominated corona. The presence of a hard power-law
X-ray component in the early spectra of TDEs may be associated with a super-
Eddington accretion rate, diverging from the situation observed in certain
Active Galactic Nuclei (AGNs) where the accretion rate is sub-Eddington (Done
& Gierliński, 2005). If we assume that the power-law component comes from the
corona, then the observed variations in coronal strength between different
TDEs may reflect different magnitudes and topologies between the seed fields
provided by the stellar debris, as magnetic initial conditions can strongly
guide later evolution of accretion disc magnetic fields and coronal emission
(Liska et al., 2022).

Figure 1: Simultaneous slim disk fits to eight Swift and three XMM–Newton 0.3–10 keV spectra for AT2018fyk from early (Swift 1) to late (Swift 8) times. The data are obtained over a four-year period. All the spectra are background subtracted. For plotting purposes only, all the spectra are rebinned such that each bin has a $5\sigma$ significance. For the XMM-Newton spectra, we only plot the pn spectra. The fitted model and best-fit parameters are listed in Table 4.  Table 4: Best-fit values and the $1\sigma$ errors on the fit parameters. We fit the spectra simultaneously with a fit-function comprised of the model components ${\tt TBAbs\times Thcomp(slimd)}$. The Galactic column density, NH, is fixed to a value of 1.15$\times 10^{20}$ cm-2. The best-fit statistic is ${\rm Cstat/d.o.f.}=6755.7/6722$.  Spectrum | date |  $10^{5}M_{\bullet}$ $[M_{\odot}]$ | $a_{\bullet}$ |  $\theta$ (deg) |  $\dot{m}$ [Edd] | $\tau/(-1)$ |  $T_{e}$ [keV] | ${\rm Cstat}/\nu$
---|---|---|---|---|---|---|---|---
Swift 1 | 2018/09/22-2018/10/17 | $0.5^{+1.5}_{-0.3}$ | $0.4^{+0.6}_{-0.8}$ | $54^{+16}_{-52}$ | $62^{+38}_{-61.2}$ | - | - | 87.8/94
Swift 2 | 2018/10/18-2018/11/21 | $2.5^{+1.0}_{-1.0}$ | $0.9^{+0.1}_{-1.9}$ | $30_{-28}^{+30}$ | $1.0^{+0.7}_{-0.2}$ | $0.1^{+0.06}_{-0.06}$ | $300_{-10}$ | 144.9/142
Swift 3 | 2018/11/23-2019/01/08 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.69^{+0.07}_{-0.05}$ | $0.11^{+0.07}_{-0.07}$ | $=T_{e,2}$ | 124.3/118
Swift 4 | 2019/03/20-2019/05/20 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $1.5^{+5.1}_{-0.6}$ | $0.25^{+0.15}_{-0.15}$ | $=T_{e,2}$ | 334.1/337
Swift 5 | 2019/05/23-2019/08/17 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $1.5^{+5.1}_{-0.6}$ | $0.25^{+0.15}_{-0.15}$ | $=T_{e,2}$ | 294.6/313
Swift 6 | 2019/08/21-2019/10/23 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.67^{+0.14}_{-0.05}$ | $0.26^{+0.16}_{-0.16}$ | $=T_{e,2}$ | 164.4/177
Swift 7 | 2019/10/27-2020/01/09 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.95^{+0.43}_{-0.34}$ | $0.31^{+0.19}_{-0.19}$ | $=T_{e,2}$ | 254.3/302
Swift 8 | 2022/03/25-2022/08/07 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.23^{+0.05}_{-0.01}$ | $0.31^{+0.19}_{-0.19}$ | $=T_{e,2}$ | 183.5/225
XMM 1 | 2018/12/09 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.73^{+0.07}_{-0.03}$ | $0.06^{+0.04}_{-0.04}$ | $=T_{e,2}$ | 1084.9/1083
XMM 2 | 2019/10/27 | $4.5^{+1.0}_{-1.0}$ | $0.998_{-0.1}$ | $60_{-10}$ | $0.44^{+0.32}_{-0.03}$ | $0.23^{+0.14}_{-0.14}$ | $=T_{e,2}$ | 3157.2/2950
XMM 4 | 2022/05/20 | $=M_{\bullet,2}$ | $=a_{\bullet,2}$ | $=\theta_{2}$ | $0.22^{+0.01}_{-0.01}$ | $0.27^{+0.17}_{-0.17}$ | $=T_{e,2}$ | 925.7/1012

####  3.1.3 Fitted Black Hole Mass and Spin

We find that the fitted $M_{\bullet}$ differs moderately when fitting the
different epochs (see Table 4). We explore the variation in the best-fit
values of $M_{\bullet}$ and $a_{\bullet}$ between different epochs in Appendix
A. We found that the best-fit $M_{\bullet}$ of XMM 2 differs from the best-fit
value found for the Swift, XMM 1 and XMM 4 spectra. XMM 2 has the highest
number of counts in the spectrum. If we perform a joint fit that includes XMM
2, the fitted values of $M_{\bullet}$ and $a_{\bullet}$ are driven by the data
of XMM 2\. However, a strong power-law component is present in the XMM 2
spectrum, making the slimd parameters less robust. Therefore, we exclude the
XMM 2 data from the joint fit that has as a goal to constrain the
$M_{\bullet}$ and $a_{\bullet}$.

In order to get the strongest constraints on $M_{\bullet}$ and $a_{\bullet}$
possible with this data, we fit the spectra Swift 1-Swift 8, XMM 1 and XMM 4
simultaneously with the fit-function ${\tt TBAbs\times Thcomp(slimd)}$. Since
we showed that the Swift 1 spectrum is consistent with just a disk spectrum,
we disable the model component Thcomp from the fit-function for Swift 1 by
fixing the covering fraction factor to 0. For all the other epochs, we assume
the same $M_{\bullet}$, $a_{\bullet}$ and $\theta$, but allow for a different
$\dot{m}$ value. Except for Swift 1, we use the model Thcomp to describe the
hard component and assume the same electron temperature, but allow for a
different optical depth at each epoch. Figure 2 shows the constraints on
$M_{\bullet}$ and $a_{\bullet}$. We obtain the best-fit value of
$M_{\bullet}=(2.7^{+0.5}_{-1.5})\times 10^{5}M_{\odot}$ and $a_{\bullet}>0.3$
at a $1\sigma$ CL. As shown in the plot, the $M_{\bullet}$ inferred from the
X-ray spectra fitting also agrees with the $M_{\bullet}$ measured by TDEmass
(Ryu et al., 2020) and MOSFiT (Mockler et al., 2019), obtained through
modeling the UV/optical emission of this event (see Section 3.2.4 for details
on UV/optical emission modeling).

Figure 2: Constraints on $M_{\bullet}$ and $a_{\bullet}$ from simultaneously
fitting a fit-function ${\tt TBAbs\times Thcomp(slimd)}$ to the X-ray spectra
of all the Swift epochs, and the XMM 1 and XMM 4 spectra. The red cross shows
the best-fit value of the parameters, while the dark and light blue regions
denote the 1 and $2\sigma$ contours, respectively. We also plot the $1\sigma$
mass range inferred from modeling the UV/optical emission by the MOSFiT and
TDEmass packages.

We also estimate $M_{\bullet}$ with the Swift-1 fits using a thin disk and a
model for a spreading TDE disk (Mummery et al., 2023b). Assuming the disk
terminates at the innermost stable circular orbit (ISCO) radius $r_{\rm
ISCO}$, the thin disk model predicts the evolution of the maximum disk
temperature as $T_{\rm eff,max}\approx 78(\frac{r_{\rm
ISCO}}{r_{g}})^{-3/4}(\frac{\dot{m}}{M_{7}\eta_{-1}})^{1/4}~{}\rm{eV}$
(Shakura & Sunyaev, 1973), where $r_{g}=GM_{\bullet}/c^{2}$,
$M_{7}=M_{\bullet}/10^{7}M_{\odot}$, $\eta_{-1}=\eta/0.1$, $G$ is Newton’s
constant, $c$ is the speed of light, and $\eta$ is the radiative efficiency).
We get the mean effective temperature $T_{\rm eff,mean}\approx
2\sqrt{\frac{r_{\rm ISCO}}{r_{\rm out}}}T_{\rm eff,max}$ (for the case of
$r_{\rm out}\gg r_{\rm ISCO}$) by integrating the flux across radii. Due to
electron scattering in a hot atmosphere with finite temperature gradients, the
real X-ray flux will be higher than the corresponding black-body flux (Shimura
& Takahara, 1995), resulting in a higher effective temperature (by a factor of
$f_{c}$, the spectral hardening factor).

For the case of $a_{\bullet}=1.0$ (where $r_{\rm ISCO}/r_{g}=1.0$ and
$\eta=0.42$) and $\dot{m}\geq 1$ 444For highly super-Eddington accretion
rates, the bolometric luminosities increase over the Eddington limit by
factors $\sim\ln(\dot{m})$ (Abramowicz, 2005). The X-ray luminosities increase
slower (even constant), due to a decreasing $f_{c}$ (Wen et al., 2020).
Overall, we assume a constant $T_{\rm eff}$ for $\dot{m}>1$. , by adopting
$r_{\rm ISCO}/r_{\rm out}=1/20$ and $f_{c}=2.0$ (for the super-Eddington
cases; Davis & El-Abd 2019), we get $T_{\rm eff}\approx
50M_{7}^{-1/4}~{}\rm{eV}$. We note that this estimate is in good agreement
with observations from ASASSN-14li, which exhibited a temperature of $T\sim
55$eV during the first 200 days post-discovery (Miller et al., 2015), with
best-fit $M_{\bullet}=10^{7}M_{\odot}$ and $a_{\bullet}=0.998$ (inferred from
slim disk modeling; Wen et al. 2020), as well as CXO J1348, where a constant
black-body temperature of $\sim 90$ eV was detected early in the outburst,
with $M_{\bullet}=10^{5.5\pm 0.5}M_{\odot}$ (inferred from the
$M_{\bullet}$–$M_{\rm bulge}$ relation; Donato et al. 2014). Given the $T_{\rm
eff}\approx 50M_{7}^{-1/4}~{}\rm{eV}$, the 150 eV inferred from diskbb fitting
for epoch 1 (see Table 3) indicates that $M_{\bullet}$ is on the order of
$\sim 10^{5}M_{\odot}$. In addition, the TDE disk model-fits using tdediscspec
constrains $M_{\bullet}$ to be $(4.4\pm 1.0)\times 10^{5}M_{\odot}$, with a
systematic error of $3.5\times 10^{5}M_{\odot}$. This systematic error arises
from the uncertainty in the relationship between $R_{\rm Peak}$ and
$M_{\bullet}$. The disk bolometric luminosity is found from tdediscspec to be
$(2.2_{-1.7}^{+5.3})\times 10^{43}$ ergs/s. The fitted $M_{\bullet}$ and the
disk bolometric luminosity, obtained by applying Equations 13 and 19 of
Mummery et al. (2023b), are consistent with the outcomes of the slim disk
modeling.

In summary, multiple types of X-ray spectral fitting analysis reveal that the
thermal X-ray emission originates from a disk surrounding a BH with a mass of
approximately $2\times 10^{5}M_{\odot}$. Our fiducial multi-epoch X-ray
spectral modeling using the slim disk model finds that the BH mass and spin
are $(2.7^{+0.5}_{-1.5})\times 10^{5}M_{\odot}$ and $a_{\bullet}>0.3$ at
$1\sigma$ CL, respectively. The corona that becomes unambiguously detected
after the first epochs contributes to the X-ray spectra in the form of (hard)
power-law emission. Its contribution grows, eventually stabilizing after about
100 days. The X-ray emission from this event shares similarities with that of
IMBH TDEs, such as 3XMM J215022.4$-$055108 (Lin et al., 2018; Wen et al.,
2021) and 3XMM J150052.0+015452 (Lin et al., 2017; Cao et al., 2023), which
also exhibit high effective temperatures and prolonged super-Eddington decay
phases. Notably, in the case of 3XMM J150052.0+015452, its early hard X-ray
component can also be explained by incorporating a corona. In the following
section, we will present additional evidence supporting the association
between the power-law X-ray emission and a corona for AT2018fyk.

###  3.2 UV/optical emission analysis

####  3.2.1 Swift UVOT and HST observations

Figure 3 shows the evolution of the UV emission as a function of time. The
magnitudes are not corrected for any extinction from our own Galaxy or the
host galaxy. The UV observations at $\sim 600$ days and the HST observation on
day $1076$ are consistent with the host galaxy’s brightness, indicating that
the TDE flare became significantly dimmer during these periods. The w2 and m2
magnitudes decline more rapidly compared to w1 around 600 days, suggesting a
substantial reduction in the temperature of the source (see also the top panel
of Figure 4). The HST observation also supports a lower temperature on day
$1076$: the magnitude of short-wavelength band F225W is significantly fainter
than that of the long-wavelength band F275W, while m2 tends to be brighter
than w1 in the period $\hbox to0.0pt{$<$\hss}{}_{{}_{{}_{\textstyle\sim}}}500$
days and after $\sim 1200$ days. Additionally, X-ray observations by the XMM-
Newton in particular observations XMM 3, Swift, and eROSITA satellite between
the days 561 and 1163 show that the brightness of the source fell below the
detection limit (Wevers et al., 2021, 2023). This further supports the
conclusion that the TDE flare experienced a significant decrease in brightness
during those periods.

To determine the decay rate of the UV light curve, we use a power-law function
for the flux, ${\rm Flux}\propto t^{n}$, to fit the w1 observations during the
initial 500 days. We obtain a best-fit index of $n=-1.70$, in good agreement
with the classical mass fallback rate index of $-5/3$ (Phinney, 1989). Even
though the UV/optical luminosity fluctuates, it generally seems to follow the
mass fallback rate. We also fit the rebrightening apparent in the w1 light
curve with a power-law function for the flux, ${\rm
Flux}\propto(t-t_{0})^{n}$, where we fix the index at the best-fit value from
the first light curve fitting, $n=-1.70$. Here, $t_{0}$ represents the fitted
peak date of the rebrightening. We obtain a best-fit value of $t_{0}=1185$ d,
with $\chi^{2}/\nu=54.3/30$. Notably, this value of $t_{0}$ is in line with
the last X-ray non-detection from $\it eROSITA$ at $t=1163.7$ d, and the first
detection after the rebrightening in the UV was at $t=1216$ d (Wevers et al.,
2023). For comparison, we also show the theoretical partial disruption
fallback decay rate, which is given by ${\rm Flux}\propto t^{-9/4}$ (Coughlin
et al., 2019). The best-fit $\chi^{2}$ values for the early $\sim 500$ days
and the rebrightening $\hbox
to0.0pt{$>$\hss}{}_{{}_{{}_{\textstyle\sim}}}1200$ days are 3530.9 and 57.68,
respectively. The first one is notably larger than the corresponding value of
1962.7, while the second one close to that of 54.3 obtained from the $n=-1.70$
fit.

Figure 3: The Swift/UVOT and HST observations of AT2018fyk. The colored
circles represent the Swift/UVOT observations, while the two pink pluses
denote the HST observations. The dashed and dotted solid vertical lines denote
the times of the observations of XMM and eROSITA, respectively. Over the time
range between $t=583$ and 1164 d, the XMM 3 and the four eROSITA observations
yielded non-detections of the source X-ray flux (Wevers et al., 2023). The
black dashed lines depict the best power-law fit to the w1 light curve, while
the black dotted lines show the best fit of $m_{\rm w1}\propto(9/4)\log_{10}t$
to the w1 light curve. This figure shows that the UV luminosity decays as the
mass fallback rate for the first 500 days, decreases quickly and reaches the
host galaxy level by $\sim 561$ days, at least until the HST observation at
$t$=1076 d. This is followed by another power-law decay after the
rebrightening that started by $t\hbox
to0.0pt{$<$\hss}{}_{{}_{{}_{\textstyle\sim}}}1216$ d.

Past studies show that TDE UV/optical light curves (after peak) decay
initially as a power-law (typically $<350$ days), then reach a plateau. The
initial power-law decay is often consistent with theoretical fallback rates,
although the precise early-time emission mechanism remains unclear (Loeb &
Ulmer, 1997; Piran et al., 2015; Roth et al., 2020). Late-time emission in the
plateau phase is thought to originate from the outer edge of a circularized
accretion disk (Ulmer, 1999; van Velzen et al., 2019; Hammerstein et al.,
2023; Mummery et al., 2023a). However, in AT2018fyk, both the decay during the
first 500 days and the decay after the rebrightening are well-fit by power-
laws, with no evidence of a plateau. The outer disk UV emission is thus
probably fainter than the observed values. Using the slim disk model (Wen et
al., 2023), we can calculate the expected $w1$ magnitude from the outer disk
using the parameters derived from fitting the X-ray spectra. For a disk around
a BH with $M_{\bullet}=2\times 10^{5}M_{\odot}$, and an accretion rate of
$\dot{m}=2$, we get $m_{w1}\sim 24$, which is indeed much fainter than the
observed value of $16<m_{w1}<20$ and also far below the host brightness
($m_{w1}\sim 21$) (Wevers et al., 2019). As the brightness of the disk in $w1$
is fainter than that of the host, we predict that the observed $w1$ emission
will continue to decay following a power-law until it reaches the host level,
which is expected to occur approximately 5000 days after discovery. However,
if the star was disrupted by the secondary of a SMBHB, accretion to the
primary BH would eventually interrupt the power-law decay.

In general, both power-law decay phases are consistent with the idea that the
UV/optical luminosity is dominated by the same emission process. Observations
from Swift, HST, XMM, and eROSITA also support the fact that the TDE flare
subsides at the very least, or even shuts off, between 561 and 1163.7 days. We
attribute the rapid decrease in both UV and X-ray emission to a cutoff in the
gas fallback stream and subsequent cessation of accretion onto the BH.

####  3.2.2 Black-body fit to the UV/optical emission

Figure 4: Evolution of the black-body radius (top panel), temperature (middle
panel), and luminosity (bottom panel) over time. Black open circles denote the
early and rebrightening epochs, respectively. The red dashed lines denote the
best-fit power-law model for the early epoch, while the red solid lines
represent the model for the rebrightening epoch. In the bottom panel, the cyan
pluses show the disk bolometric luminosity as a function of time derived from
X-ray spectral fitting of the 8 Swift averaged spectra.

The SED of the early power-law decay UV/optical emission in TDEs can often be
well described by thermal black-body emission (van Velzen et al., 2011). We
therefore measure the black-body temperature and radius, as well as the
luminosity, by fitting a single-temperature black-body model to the UV/optical
data as a function of time. For all the UV/optical observations available, we
subtract the host galaxy contribution as determined from SED fitting in Wevers
et al. (2021). The host galaxy subtracted fluxes are also corrected for
Galactic extinction by employing the IDL code fm_unred.pro with $E(B-V)=0.01$
(Schlegel et al., 1998), assuming the Milky Way (MW) extinction curve (Gordon
et al., 2003). To carry out the black-body model fitting of the UV/optical
data, only observations with positive host-subtracted fluxes in more than 2
filters are used555The host contribution dominates the emission in some
observations, which might result in negative host-subtracted flux densities,
especially in optical bands.. We use the flux density at a specific
wavelength, $f_{\lambda}=4\pi R_{\rm bb}^{2}D_{L}^{-2}B_{\lambda}(T_{\rm
bb},\lambda)$, to do the fitting. Here, $R_{\rm bb}$ and $T_{\rm bb}$ are the
black-body radius and temperature, respectively, and $B_{\lambda}(T_{\rm
bb},\lambda)$ is the Planck law as a function of wavelength. We then calculate
the black-body luminosity with the Stefan-Boltzmann law, $L_{\rm
bb}=4\pi\sigma_{sb}R_{\rm bb}^{2}T_{\rm bb}^{4}$, where $\sigma_{sb}$ is the
Stefan-Boltzmann constant.

Figure 4 shows the black-body temperature, radius, and luminosity as a
function of time. The relatively large uncertainties are due to the fact that
the peak of the black-body spectra is not sampled by the observations. The
evolution of the temperature and radius is consistent with the results of
Wevers et al. (2021) for both the early and the rebrightening epochs, e.g.,
constant temperature and slowly decreasing radius.

We fit the luminosity with the equation $\log_{10}L_{\rm
bb}=-\alpha\log_{10}(t-t_{0})+C$, where $t_{0}$ is the peak day and $C$ is the
luminosity at $t-t_{0}=1$. Here, we set the zero point of $t$ at MJD 58369.2
days. Giving $t_{0}=0$ d for the early epochs, the best-fit parameters are
$\alpha=0.65\pm 0.09$ and $C=44.82\pm 0.20$. During the rebrightening epochs,
assuming $\alpha$ is the same as the early epochs, the best-fit $t_{0}$ and
$C$ are $1221.92\pm 41.19$ d and $44.14\pm 0.08$, respectively. The luminosity
decays as a power-law for both the early and rebrightening epochs, but with a
smaller index of 0.69 when compared with the power-law decay index of $w1$.
Integrating the best-fit over time, we obtain the total black-body energy
released for the early and rebrightening epochs; they are $1.1\times 10^{51}$
ergs (close to the estimation of $1.4\times 10^{51}$ ergs from Wevers et al.
2019) and $1.1\times 10^{50}$ ergs, respectively.

If we adopted a radiation efficiency of $\eta=0.01$ (see our analysis in the
next subsection), the radiated energy amounts to $0.06M_{\odot}$, and
$0.006M_{\odot}$ for the early and rebrightening epochs, respectively,
indicating $M_{\star}>0.13M_{\odot}$ assuming 50% of the star’s material is
unbound and 50% is accreted. This is in agreement within 1 $\sigma$ with the
star’s mass $M_{\star}=0.07^{+0.12}_{-0.04}M_{\odot}$ estimated from fitting
the UV light curves (see Table 5). For the early epochs, the peak black-body
luminosity and temperature are $10^{44.28\pm 0.75}$ erg/s (in agreement with
the estimate of $(3.0\pm 0.5)\times 10^{44}$ erg/s from Wevers et al. 2019)
and $10^{4.72\pm 0.17}$ K, respectively. Given the Eddington luminosity
$L_{\rm Edd}\approx 3.3\times 10^{43}$ erg/s (for the case of
$M_{\bullet}=2.7\times 10^{5}M_{\odot}$), the UV black-body luminosity is
about 10 times the Eddington luminosity. It takes about 200 days for the
UV/optical emission to decline to the Eddington luminosity. We also plot disk
bolometric luminosities derived from the slim disk modeling of the Swift X-ray
spectra. We fit each of the Swift spectra with the model ${\tt TBAbs\times
Thcomp\times slimd}$ separately, and then calculate the unabsorbed disk
bolometric luminosity with the disk parameters from the slimd component. As
shown in Figure 4, the BB luminosity of the UV/optical emission is brighter
than the slim disk bolometric luminosity for the first epoch.

It is possible to generate super-Eddington luminosities in the early UV
emission. Shock-powered emission (Piran et al., 2015) is not necessarily
Eddington-limited, and neither is reprocessing in an unbound disk wind
(Metzger & Stone, 2016) if the central engine is itself super-Eddington. The
cooling envelope model (Metzger, 2022) predicts that the bolometric luminosity
is about $L=L_{\rm Edd}+L_{\rm fb}$, where $L_{\rm Edd}$ is the disk Eddington
luminosity and $L_{\rm fb}$ is produced by the fall back debris. We present
observational evidence supporting the presence of a super-Eddington
luminosity. Firstly, Wevers et al. (2019) detected low-ionization Fe II
emission lines about a month after discovery. These emission lines require an
obscuring medium with significant particle density and optical depth, along
with an energy source to heat the gas. This obscuring medium can be attributed
to the super-Eddington outflow. The presence of these emission lines is
consistent with our assumed prior on the slim disk inclination angle to our
line of sight, $\theta<60^{\circ}$, since a sight line with
$\theta>60^{\circ}$ implying we would be looking through the outflow, might
have resulted in absorption lines. In addition, the X-ray flux exhibited a
gradual increase, eventually reaching a plateau around 200 days (Wevers et
al., 2023). The gradual increase in the X-ray flux (accretion rate) can be
plausibly explained by the gradual decrease in the super-Eddington outflow,
leaving more material in g/s to be accreted.

####  3.2.3 The origin of the UV/optical emission

The source of early power-law UV/optical emission of TDEs is still uncertain.
For instance, it may be produced by a reprocessing layer, which is powered by
the X-ray emission emanating from the inner accretion disk (Loeb & Ulmer,
1997; Guillochon et al., 2014; Metzger & Stone, 2016; Roth et al., 2016), from
the outer shock that forms at the intersection of the debris streams near
their orbital apocenters (Dai et al., 2015; Shiokawa et al., 2015; Piran et
al., 2015; Bonnerot & Lu, 2020; Lu & Bonnerot, 2019), or from an elliptical
accretion disk (Zanazzi & Ogilvie, 2020; Liu et al., 2021).

The bolometric disk luminosity of Swift 1 can be inferred either from slim
disk modeling ($1.4^{+4.0}_{-0.8}\times 10^{43}$ erg/s) or from the tdediscpec
(Mummery et al., 2023b) model $(2.2^{+5.3}_{-1.7}\times 10^{43}$ erg/s). These
models reach statistically consistent estimates, each of which is
significantly fainter than the peak UV luminosity (see also Figure 4). If the
UV/optical emission originates from a reprocessing layer, the disk would not
be initially bright enough to power the layer (though at later epochs this
conclusion no longer holds, see Figure 4). This indicates that the initial
UV/optical emission in this event comes from shock dissipation or the outer
regions of an elliptical accretion disk, and the lack of significant
temperature or radius evolution suggests that this may remain the UV/optical
power source at later epochs. In Section 3.3.2, we will present additional
evidence to support this conclusion.

Here, we show that the outer shock scenario can explain the brighter-than-the-
disk UV/optical emission. Due to general relativistic pericenter precession,
debris streams will undergo a self-intersection after each pericenter passage
(Rees, 1988; Dai et al., 2015; Shiokawa et al., 2015; Lu & Bonnerot, 2019;
Bonnerot & Lu, 2020; Andalman et al., 2022). Other shocks may also contribute
to the hydrodynamic evolution of the debris (Guillochon et al., 2014;
Steinberg & Stone, 2022), with the net result being dissipation of orbital
energy and the formation of a disk with outer radius equal at least to the
circularization radius $R_{\rm c}=2R_{\rm p}$, where $R_{\rm p}$ is the
pericenter of the disrupted star. The rate of energy release through the
circularization process is $\dot{E}_{\rm c}=\mu GM_{\bullet}\dot{M}/(2R_{\rm
c})$ (here, we add $\mu\leq 1$ to account for incomplete circularization,
which would result in an elliptical rather than a circular disk). The
circularization efficiency is then $\eta_{c}=\mu r_{\rm g}/(4R_{\rm p})$. For
super-Eddington accretion disks, the disk radiation efficiency, $\eta_{\rm
d}$, decreases significantly as the accretion rate increases, and typically
$\eta_{d}<0.1$ for $a_{\bullet}<0.9$, due to advective cooling (Abramowicz et
al., 1988; Sadowsk, 2009). For highly super-Eddington accretion disks, their
bolometric luminosities increase over the Eddington limit as
$\propto\ln\dot{m}$ (Abramowicz, 2005). As a result, the disk radiation
efficiency $\eta_{\rm d}$ drops by a factor $\sim(\ln\dot{m})/\dot{m}$. For
the case of $M_{\bullet}=2\times 10^{5}M_{\odot}$, we obtain
$\eta_{c}=0.002\mu\beta$ (where $\beta=R_{\rm p}/r_{\rm t}$ is the star’s
penetration parameter) and $\eta_{d}<0.0005$ (assuming $\dot{m}\sim 1500$).
For the case of $\mu=1$ and $\beta=1$, we find that $\eta_{c}>\eta_{d}$,
resulting in the UV/optical bolometric luminosity being larger than the disk
bolometric luminosity.

For the disruption of solar-like stars with $R_{p}=r_{t}$, the debris
circularization efficiency ($\eta_{c}\propto M_{\bullet}^{2/3}$) decreases
more slowly than the disk efficiency ($\eta_{d}\propto\sim M_{\bullet}^{3/2}$)
as $M_{\bullet}$ decreases. The shock paradigm predicts that TDEs with
$M_{\bullet}<\text{a few}\times 10^{5}~{}M_{\odot}$ will have a brighter peak
UV/optical bolometric luminosity than the disk bolometric luminosity. Future
observations of IMBH TDEs will help verify this prediction.

####  3.2.4 BH mass inferred from UV/optical emission

Currently, various models are employed to fit UV/optical light curves and to
estimate the black hole mass and other parameters of interest. Notable
examples include MOSFiT (Modular Open Source Fitter for Transients) (Mockler
et al., 2019), TDEmass (Ryu et al., 2020) and Redback (Sarin & Metzger, 2023).
TDEmass assumes that the UV emission originates from shocks near the apocenter
of the most bound debris (Piran et al., 2015), while Redback postulates that
the UV emission arises from the cooling envelope formed by the returning
debris (Metzger, 2022). MOSFit makes fewer physical assumptions but instead
posits a power-law relationship between various observables (Guillochon et
al., 2014). For our analysis, we utilize both TDEmass and MOSFiT to place
constraints on TDE parameters.

TDEmass posits that the dissipation of orbital energy and the production of
UV/optical emission occur due to shocks near the debris apocenter. Its
radiation efficiency is a function of the semimajor axis ($a_{\rm o}$) of the
debris’ orbit, $\eta_{c}\propto r_{\rm g}/a_{\rm o}$. Its temperature is
determined by $T^{4}=L/(2\pi a_{\rm o}^{2}\sigma_{\rm sb})$. As $a_{\rm o}$ is
a function of the specific energy of the debris, TDEmass does not aim to fit
the entire UV light curve but instead focuses on the peak emission, primarily
constraining two key parameters: $M_{\bullet}$ and $M_{\star}$. The estimation
of $M_{\bullet}$ relies heavily on the temperature observed at peak
luminosity, while $M_{\star}$ primarily depends on the peak luminosity. We use
the peak BB luminosity and temperature of the early epochs, $L=10^{44.28\pm
0.75}$erg/s and $T=10^{4.72\pm 0.17}$K, for the fit.

MOSFiT calculates the black-body UV/optical luminosity based on the mass
fallback rate, assuming a time-independent efficiency ($\eta_{\rm c}$). The
mass fallback rate $\dot{M}_{\rm f}(t)$, obtained from hydrodynamical
simulations of TDEs (Guillochon & Ramirez-Ruiz, 2013), depends on parameters
such as $M_{\bullet}$, $M_{*}$, and $\beta$. Instead of predicting the BB
radius and BB temperature from a specific theoretical model, MOSFiT fits these
parameters to the observed data, assuming that the BB radius decreases as a
power-law of the luminosity, $R_{\rm BB}\propto(L/L_{\rm
Edd})^{l}=\dot{m}^{l}\propto t^{l^{\prime}}$. This assumption aligns with
predictions from certain theoretical models; for instance, elliptical disk
models predict $R_{\rm BB}\propto t^{-1.1}$ (Liu et al., 2014), while cooling
envelope models predict $R_{\rm BB}\propto t^{-1}$ (Metzger, 2022), although
it is in disagreement with some radiation-hydrodynamics simulations (Steinberg
& Stone, 2022). An additional assumption added in MOSFiT is that the
luminosity does not follow the mass fallback rate directly, but only follows
convolution with a “viscous delay,” adding an additional free parameter to
account for potential delays in the circularization process.

MOSFiT sets a soft cap to prevent the UV black-body luminosity from exceeding
the Eddington luminosity. However, as we have shown, a super-Eddington UV
luminosity is not impossible. As a result, we remove the soft cap from the
code, allowing the luminosity to exceed the Eddington luminosity. As we have
shown, the decay of the $w1$ magnitude as a function of time follows the mass
fallback rate, a constant radiation efficiency may be correct for early UV
emission of this event. The free parameters are $M_{\bullet}$, $M_{\star}$,
$\beta$, viscous decay time $T_{\rm v}$, $\eta_{\rm c}$, BB radius
normalization factor $n_{\rm R}$, $l$, and the first fallback time $t_{0}$
(days since first detection). We restricted our fitting to the $w2$, $m2$, and
$w1$ light curves from the early epochs, and we applied priors for two
specific parameters: $t_{0}$ and $\eta_{c}$, with their ranges set to [-30,0]
days and [$10^{-4}$, 0.4], respectively. We use the Markov–Chain–Monte–Carlo
routine in MOSFiT to perform the fit and stop the calculation until it
converged, with the potential scale reduction factor (PSRF) $\leq 1.1$ (Gelman
& Rubin, 1992).

Table 5 shows our fit results. Both TDEmass and MOSFiT produce consistent
estimates for $M_{\bullet}$, which are also consistent with the result
obtained from our X-ray spectral modeling (as illustrated in Figure 2). In the
case of TDEmass, the determination of the stellar mass ($M_{\star}$) is
primarily driven by the peak luminosity, $L_{\rm p}$ (Ryu et al., 2020). This
reliance on $L_{\rm p}$ stems from the fact that TDEmass fixes the radiation
efficiency at $\eta_{\rm c}=r_{\rm g}/a_{\rm o}$, giving a relatively small
value of 0.0003 (see Table 5). This value is considerably smaller than the
circularization efficiency got by MOSFiT ($\eta_{c}=r_{g}/(4R_{p})=0.01$),
which accounts for the significantly larger fitted values of $M_{\star}$
obtained by TDEmass compared to those from MOSFiT. The $R_{\rm BB}$ of MOSFiT
can be recovered from the relationship, $R_{\rm BB}=n_{\rm R}a_{o}(L/L_{\rm
Edd})^{l}$. From the fit results, we see that $R_{\rm BB}$ evolves with time
very slowly, $R_{\rm BB}\propto t^{-0.1}$, and $R_{\rm BB}\sim 0.7a_{o}$,
consisting with the location of the emission assumed by TDEmass 666Note that
TDEmass assume a solid angle of $2\pi$, while MOSFiT assume a solid angle of
$4\pi$. There should be a factor of $\sqrt{2}$ between the two $R_{\rm BB}$.
The MOSFiT model fits the UV emission well, consistent with a stream self-
interaction shock scenario. We note that the small fitted $M_{\star}$ value is
consistent with the result from X-ray spectral fitting, which found between
$\approx 0.005-0.065M_{\odot}$ of debris accreted during the early epochs.
There are at least three possible reasons to account for the small fitted
$M_{\star}$ and total mass accretion: 1) the disruption is attributed to the
secondary of a binary system, leading to an orbital energy shift of the debris
caused by the primary (for further analysis, we refer to Section 4.3); 2) a
partial disruption; 3) the super-Eddington UV luminosity blows out a
significant amount of debris. In the first of these explanations, a shift in
the orbital energy results in less debris being bound to the secondary,
consequently leading to a smaller estimated value for $M_{\star}$.

Table 5: Best-fit parameters and their constraints derived from modelling the early UV emission with MOSFiT and TDEmass. We calculate the $\eta_{c}$ from the best-fit values of $M_{\bullet}$ and $M_{\star}$ for the case of TDEmass. See Section 3.2.4 for the meaning of the symbols. Parameter | TDEmass | MOSFiT
---|---|---
$M_{\bullet}$ / $10^{5}M_{\odot}$ | $5.4^{+10}_{-4.9}$ | $8.8^{+10}_{-7.7}$
$M_{*}$ / $M_{\odot}$ | $1.2_{-0.77}^{+88}$ | $0.07^{+0.12}_{-0.04}$
$\eta_{c}/10^{-2}$ | $0.03^{+0.01}_{-0.02}$ | $1.0^{+3.0}_{-1.0}$
$n_{\rm R}$ | - | $0.69^{+0.74}_{-0.21}$
$\beta$ | - | $0.9^{+0.01}_{-0.02}$
$t_{0}$ / days | - | $-29^{+2}_{-1}$
$T_{\rm v}$ / days | - | $0.03^{+0.32}_{-0.03}$
$l$ | - | $0.05^{+0.07}_{-0.04}$

###  3.3 X-ray versus UV/optical Timing

In this section, we investigate whether there are any time lags between the
X-ray and UV/optical emissions, to explore possible relationships or sources
of the UV/optical and hard X-ray emission.

####  3.3.1 X-ray Spectral Lags

Zhang (2022) reported tentative evidence of soft lags in one XMM-Newton
observation (obs id: 0853980201; XMM 2) of AT2018fyk. In a soft lag the
arrival time of lower energy photons lags behind that of higher energy
photons. In their frequency-resolved analysis, the soft lags were detected in
only the lowest frequency bin. To further investigate this result, we re-
computed the frequency-resolved lags from this observation, using the same
choice of energy bands and reference band (0.5–1 keV) as Zhang (2022). While
finding results in general agreement with their analysis, the number of
frequencies available in this lowest-frequency bin is small ($\sim$2), a
factor of 2-3 less than recommended for reliable time lag estimates (Uttley et
al. 2014). Given that the soft lags are detected in this bin only, and that
the count rate of the observation is very low ($<$1.5 counts/s in the 0.3–5
keV band), we do not find the lags to be statistically significant. As an
additional check, we computed the (non-frequency-resolved) lags using the
linear Interpolated Cross-Correlation Function (ICCF) technique of Peterson et
al. (1998), using the same energy bands/reference band as before. We again
find no significant time lag, with the maximum correlation coefficient (often
denoted $R_{\rm{max}}$) being less than 0.1 in all cases.

####  3.3.2 Cross-Correlations between UV, Soft X-ray, and Hard X-ray Light
Curves

Figure 5: The light curves and the cross correlation functions. Left panels:
The L1–L4 are as defined in the top right hand side of the left panels. From
top to bottom; the Swift $w1$ (L1), 1.5–5 keV (L2), 0.3-0.5 keV (L3), and the
NICER 0.3-0.5 keV (L4) light curves. Middle and right panels: The black solid
lines denote the correlation values of each time lag bin (with a 1-day width).
The red dotted and the blue dashed lines show the $2\sigma$ and $3\sigma$ CLs,
respectively. Positive time lag denotes the first light curve is leading. The
middle panels show the correlation values of the early light curves, while the
right panels show the correlation values of rebrightening light curves. These
plots show marginal evidence for a correlation between the hard X-ray emission
(L2) and the soft X-ray emission (L3 and L4). The centroids lags of L1 vs. L3,
L2 vs. L3 and L2 vs. L4 of the middle plot are $3.8_{-1.9}^{+1.6}$,
$-0.5_{-1.6}^{+0.9}$ and $1.4^{+19.3}_{-1.4}$ days, respectively.

We use the ICCF (Peterson et al., 1998) to explore whether there are
correlations among the UV (L1; $w1$), soft X-ray (0.3-0.5 keV; L3 Swift; L4
NICER), and hard X-ray (1.5-5 keV; L2 Swift) emission. We use the package
PYCCF to calculate the ICCF (Sun et al., 2018). We estimated the uncertainty
in the peak and the centroid of the ICCFs by using the Random Subset Selection
and the Flux Randomization procedures (Peterson et al., 1998, 2004).The study
(Welsh, 1999) shows that extended trends in light curves can introduce biases
in the ICCF analysis. Consequently, removing these trends by fitting the light
curves with a smooth function enhances the accuracy of recovering correlations
and lags between two considered time series. The Swift w1 light curve exhibits
a power law decay during the initial $\sim 500$ days and subsequent
rebrightening after $>1200$ days. We detrend the w1 light curve by subtracting
the best-fit power law values, as detailed in Section 3.2.1. However, for the
X-ray light curves, there are no clearly discernible long-term trends, as
depicted in the left panel of Figure 5. Consequently, we opt not to detrend
them during the ICCF analysis. We calculate the ICCF of the early light curves
and the rebrightening light curves separately. In order to determine the
confidence of the correlations, we generate $10^{5}$ random white noise X-ray
light curves using the method described in Timmer & Koenig (1995). The
randomly generated X-ray light curves require a power-law index, which relate
power and frequency. We adopt the power-law index as $1.49$ and $1.32\pm 0.3$
for the band 0.3-0.5 keV and 1.5-5 keV (Zhang, 2022), respectively 777Our
conclusions are unchanged for a different index of 0.5 or 2.5 for both bands..
We then cross-correlate the light curves with each of these simulated X-ray
light curves, to built a distribution of $10^{5}$ correlation values at each
time lag bin. Using these distributions, we extracted the $95\%$ ($2\sigma$)
and $99.7\%$ ($3\sigma$) confidence contours.

Figure 5 show the light curves and the correlation values of each time lag
bin. The time lags of L1 vs. L2 and L1 vs. L3 are $17^{+23}_{-8}$ days and
$3.8^{+1.6}_{-1.9}$ days for the early emission, respectively. However, the
significance of the two time lags is less than $2\sigma$. During the
rebrightening phase, the positive lag between L1 and L3 weakens, and four high
correlation values appear for L1 and L2. This may be attributed to the bad
observed cadence and low detected X-ray count rates being suboptimal.
Recently, Guo et al. (2023) found a similar lag of $9.2^{+6.1}_{-6.0}$ between
the w1 and X-ray emission, limited to the initial 120 days. The positive time
lag between L1 and L3 here aligns with that from Guo et al. (2023) within
$1\sigma$ CL. As proposed by Guo et al. (2023), the positive lag disfavors the
assumption that the early UV/optical emission originates from a reprocessing
layer but supports the paradigm that the early UV/optical emission comes from
the outer shocks. Nevertheless, the lack of a significant correlation between
the UV and X-ray bands seems consistent with the X-ray source being centrally
compact and dominated by light bending close to the SMBH. These findings are
in line with the early analysis suggesting that the UV/optical emission
originates from the outer shocks.

The correlation analysis between hard and soft X-ray emissions reveals a
significance above the $3\sigma$ confidence level for the L2 vs. L3 case
during the initial 600 days. The calculated time lag is $-0.5^{+0.9}_{-1.6}$
days, consistent with a lag of 0 days. This 0-time lag is also observed in the
L2 vs. L3 case during the rebrightening phase and the L2 vs. L4 case in both
the early and rebrightening phases. However, it is noteworthy that the
significance of the correlation diminishes in the rebrightening phase and
between L2 and L4, possibly due to poor cadence during that period and
variations in the NICER observation cadence. Despite this, the persistent
0-time lag suggests a correlation between the hard and soft components. This
finding aligns with our analysis in Section 3.1.2, supporting the hypothesis
that the hard component originates from the corona above the disk, which
generates the observed soft X-ray emission. A correlation between the hard and
soft X-ray emission is also shown in some AGN system, e.g. Narrow-line Seyfert
1 galaxy 1H 0707-495 (Robertson et al., 2015), IRAS 13224-3809 (Chiang et al.,
2015). We also conducted the ICCF by de-trending the X-ray light curves
through subtracting the best-fit values of the power-law model. The results
also align with the findings here, e.g. no significant correlation between
$w1$ and soft- or hard X-ray emissions, and a $3\sigma$ confidence level
correlation between soft and hard X-ray emissions for case of L2 vs. L3 in the
initial 600 days.

##  4 Discussion: A TDE in a Binary SMBH system?

In this section, we interpret our previous data analysis as evidence for the
production of a TDE in a binary SMBH system and explore its implications.

###  4.1 Discrepancy with Black Hole Mass from Galaxy Scaling Relation

The lack of a high equivalent width H-$\delta$ absorption line implies that
the galaxy does not belong to the E + A galaxy class (Wevers et al., 2019).
Wevers (2020) measured the velocity dispersion of the host galaxy of the
AT2018fyk to be $\sigma_{*}=159\pm 1\rm{km~{}s^{-1}}$. Using the
M–$\sigma_{*}$ relation of McConnell & Ma (2013), Wevers (2020) estimated the
host BH mass as $M_{\bullet}=10^{7.7\pm 0.4}M_{\odot}$. Adopting the
M–$\sigma_{*}$ relation of Ferrarese & Ford (2005), the host BH mass is
$M_{\bullet}=10^{7.7\pm 0.4}M_{\odot}$, whereas it is $M_{\bullet}=10^{8.0\pm
0.4}M_{\odot}$ for the M–$\sigma_{*}$ relation of Kormendy & Ho (2013). A
measurement from the $M_{\bullet}$–$M_{\rm bulge}$ relation also supports the
inference of a high-mass SMBH, for example, Wevers et al. (2019) found that
$M_{\bullet}=10^{7.3\pm 0.4}M_{\odot}$ using the $M_{\bullet}$–$M_{\rm bulge}$
relation of Häring & Rix (2004).

In contrast, the $M_{\bullet}$ inferred from our fiducial X-ray spectral fits
is orders of magnitude lower, $(2.7^{+0.5}_{-1.5})\times 10^{5}M_{\odot}$.
Even considering a (non-fiducial) joint fit with all the epochs, the X-ray
parameter estimation still rules out $M_{\bullet}>7\times 10^{5}M_{\odot}$ at
$>3\sigma$ CL. Likewise, the $M_{\bullet}$ inferred from the UV emission is
$10^{4.7}-10^{6.3}M_{\odot}$, further supporting the small SMBH mass inferred
from the X-ray spectral fits. The discrepancy between the inferred
$M_{\bullet}$ from X-ray spectral fits and that derived from the M-$\sigma$
relation is consistent with the idea that the thermal X-ray flare after the
disruption is not generated from a disk around the host galaxy’s primary SMBH,
but rather from a disk around a BH with much smaller $M_{\bullet}$, orbiting
close to the central SMBH.

Arguments based on galaxy scaling relations are of course approximate. The
M–$\sigma_{*}$ measurement may over estimate the inferred $M_{\bullet}$, e.g.
Koss et al. (2022) found that the M–$\sigma_{*}$ relation of Kormendy & Ho
(2013) overestimates $M_{\bullet}$ by 0.7 dex (this would be 0.4 dex if using
the relationship of Ferrarese & Ford 2005 and McConnell & Ma 2013) by
comparing to the direct $M_{\bullet}$ measurement from megamasers
($M_{\bullet}<10^{8.5}M_{\odot}$). Such an overestimate is also reported by
Greene et al. (2016). In general, however, even these possible systematic
biases are unlikely to resolve the $\sim 2.5$ order of magnitude discrepancy
between X-ray continuum fitting measurements of $M_{\bullet}$ and estimates
from galaxy scaling relations.

###  4.2 Rapid UV Light Curve Decline Suggesting Accretion Cutoff

A star with mass $M_{*}$ and radius $R_{*}$ will be tidally disrupted as it
approaches a SMBH within the tidal radius,
$r_{t}=R_{*}(M_{\bullet}/M_{*})^{1/3}$ (Rees, 1988; Phinney, 1989). After the
disruption at the tidal radius, the specific energy E across the stellar
radius varies from $E_{o}-\Delta E$ to $E_{o}+\Delta E$, where $E_{o}$ is the
original orbital binding energy of the star, and $\Delta
E=GM_{\bullet}R_{*}/r_{t}^{2}$ is the spread in specific energy across the
stellar radius. Typically, $E_{o}$ is negligible compared to $\Delta E$, as
most stars tidally disrupted in realistic galactic nuclei approach the SMBH on
nearly zero energy orbits whose apocenters lie at parsec scales (Magorrian &
Tremaine, 1999; Wang & Merritt, 2004). The bound debris with $E<0$ moves in
high eccentricity Keplerian orbits, and falls back to the tidal radius after
one Keplerian period, $T=2\pi GM_{\bullet}(-2E)^{-3/2}$. The debris returning
to the tidal radius generates a mass fallback rate,

$\dot{M}_{\rm f}=133\dot{M}_{\rm
Edd}M_{6}^{-3/2}m_{*}^{2}r_{*}^{-3/2}\left(\frac{t}{t_{\rm
f}}+1\right)^{-5/3},$ (1)

where $m_{\star}=M_{\star}/M_{\odot}$ and $r_{\star}=R_{\star}/R_{\odot}$. The
return time of the most bound debris is

$t_{\rm f}=40{~{}\rm days}\times M_{6}^{1/2}m_{*}^{-1}r_{*}^{3/2}.$ (2)

As the (initial) viscous time scale is much shorter than this orbital period,
debris will be accreted quickly once it returns to pericenter, so long as
circularization is effective (e.g. Steinberg & Stone 2022, though see also Ryu
et al. 2023). If this is the case, an accretion disk will be formed and its
accretion rate will initially follow the mass fallback rate, deviating at
later times once the viscous time grows (Cannizzo et al., 1990; van Velzen et
al., 2019).

In a SMBH binary system, the star can in principle be disrupted by both the
primary, of mass $M_{\bullet,\rm p}$, and the secondary, of mass
$M_{\bullet,\rm s}$ (Liu et al., 2009; Coughlin et al., 2017). After
disruption, the most tightly bound debris falls back to the pericenter of the
SMBH that caused the TDE as in a single SMBH system; gravitational torques
from the companion SMBH cause only small perturbations on the slow secular
timescale (Kozai, 1962; Lidov, 1962; Liu et al., 2014). For less tightly bound
debris, with an apocentre larger than a critical radius, fluid elements will
move on chaotic orbits (Liu et al., 2009; Coughlin et al., 2017) due to the
lack of a hierarchy in the gravitational fields of the two SMBHs. As a result,
this gas will fail to fall back to the pericenter after disruption, producing
an accretion cutoff at a time $T_{\rm c}$ (Liu et al., 2009). Following
Coughlin et al. (2017), we take the Hill sphere as this critical radius. As we
have shown that the TDE flare is due to the smaller BH, the corresponding
Keplerian time $T_{\rm c}$ can be estimated as (Coughlin & Armitage, 2018)

$T_{\rm c}=22{~{}\rm day}\times M_{\rm p,8}^{-1/2}\left(\frac{a}{10^{-3}{\rm
pc}}\right)^{3/2},$ (3)

here, $M_{\rm p,8}=M_{\bullet,\rm p}/(10^{8}M_{\odot})$ and $a$ is the
separation of the binary.

The cutoff time $T_{\rm c}$ can be directly measured by observations, and used
to constrain a combination of $a$ and $M_{\rm p}$. AT2018fyk was discovered on
2018 September 8 (MJD 58369.23), with a non-detection reported in the
$g^{\prime}$ band ($g^{\prime}>$ 17.4 mag) on 2018 August 29 (MJD 58359). As
shown in Figures 3 & 5, there is a rapid decline for both UV/optical and X-ray
light curves on 2020 March 22 (MJD 58930.15), indicating an accretion cutoff
starts approximately on this day. The last observation before the sharp cutoff
is on 2020 Jan 09 (MJD 58857.6). The upper and lower limit of $T_{\rm c}$ are
therefore 571 and 488 days, respectively. When the disk is cut off from an
external resupply of debris, the disk accretion rate will drop quickly, over a
timescale much smaller888More specifically, the accretion rate in an isolated
disk evolves as $\dot{M}\propto(t_{\rm v}/t)^{4/3}$ (Cannizzo et al., 1990;
Liu et al., 2014). Here, the viscous time $t_{\rm v}\approx 0.034~{}{\rm
days}\times\alpha^{-1}\beta^{-3/2}(H/R_{c})^{-2}$. The initial ($H/R\sim 0.5$)
viscous time at $R_{\rm c}$ is thus $t_{\rm v}\sim 1~{}\rm{days}$. As can be
seen from Figure 4, the BB luminosity drop by a factor of $\sim 10$ during the
long declining phase, but if the disk were isolated, the time for the
accretion rate to drop by a factor of 10 would be a far shorter $5.6t_{v}$
$\sim 8$ days, showing that the initial evolution is fallback-regulated rather
than viscously regulated. than the uncertainty on $T_{\rm c}$, and we
therefore neglect this complication. We thus take $T_{\rm c}=530\pm 40$ days.
Together with $M_{\bullet,\rm p}=10^{7.7\pm 0.4}M_{\odot}$, we obtain
$a=(6.7\pm 1.2)\times 10^{-3}{~{}\rm pc}$.

Assuming a circular orbit, the period of this binary is $T_{b}=2\pi
a^{3/2}/\sqrt{GM_{\bullet,\rm p}}=7.4\pm 2.2{\rm~{}yr}$. When we use a
different accretion cutoff assumption, as suggested by Liu et al. (2009), we
obtain a similar binary period, $T_{b}\approx T_{c}/\epsilon=4.8\pm 2.7$ yr,
with semi-empirical constant $\epsilon=0.32\pm 0.18$. This result is
consistent with the $7.4\pm 2.2{\rm~{}yr}$ above. The merger time of a binary
can be estimated as (Peters, 1964)

$\displaystyle T_{\rm m}=$ $\displaystyle 580{\rm~{}yr}\times
M_{p,8}^{-3}\left(\frac{a}{10^{-3}\rm pc}\right)^{4}\frac{M_{\bullet,\rm
p}^{2}}{M_{\bullet,\rm s}(M_{\bullet,\rm p}+M_{\bullet,\rm s})}$
$\displaystyle\approx$ $\displaystyle 1.8\times 10^{9}~{}{\rm
yr}~{}\left(\frac{T_{\rm c}}{530~{}\rm
d}\right)^{8/3}\left(\frac{M_{\bullet,\rm p}}{10^{7.7}M_{\sun}}\right)^{-2/3}$
(4) $\displaystyle\times\left(\frac{M_{\bullet,\rm
s}}{10^{5.4}M_{\sun}}\right)^{-1}.$

In the second equation, we have approximated $M_{\bullet,\rm p}\gg
M_{\bullet,\rm s}$. With the values of $M_{\bullet,\rm p}$, $M_{\bullet,\rm
s}$ and $T_{\rm c}$ as given above, we obtain $T_{\rm
m}=1.8^{+1.8}_{-0.9}\times 10^{9}{\rm yr}$. This result passes a non-trivial
consistency test for our model, which would have suffered from fine-tuning if
we had determined $T_{\rm m}\ll$ the typical time between SMBHB formation
(e.g. the galaxy merger timescale).

###  4.3 Debris Bound to Secondary versus Primary

Here, we explore the consequences of the disruption being done by the
secondary SMBH in the binary. We define an effective energy, $E_{\rm eff}$,
instead of the energy $E$, to take the impact from the presence of the primary
SMBH into account, $E_{\rm eff}=V_{\rm E}^{2}/2-GM_{\bullet,\rm s}/r_{s,t}$,
where $V_{\rm E}$ is the debris velocity relative to the secondary at the
tidal radius, $r_{s,t}$ is the tidal radius of the secondary. $V_{\rm E}$ can
be determined by, $V_{\rm E}^{2}/2=V^{2}/2+V_{o}^{2}/2+VV_{o}\cos\phi$. Here,
$V\approx\sqrt{2E+2GM_{\bullet,\rm s}/r_{s,t}+2GM_{\bullet,\rm p}/a}$ is the
total velocity of the debris, $V_{o}=\sqrt{GM_{\bullet,\rm p}/a}$ is the
orbital velocity of the secondary and $\phi$ is the difference in the
direction between $V$ and $V_{o}$. Setting $f=2V\cos\phi/V_{o}+3$ (typically,
$V>V_{o}$), we get $V_{\rm E}^{2}/2=E+GM_{\bullet,\rm
s}/r_{s,t}+fGM_{\bullet,\rm p}/2a$. Following Coughlin et al. (2017), we
define an energy ratio,

$K_{E}=\frac{GM_{s}R_{*}}{r_{s,t}^{2}}\frac{2a}{GM_{\bullet,\rm p}}.$ (5)

The effective energy of the most bound debris to the secondary can be written
as,

$E_{\rm eff}\approx-\Delta E+\frac{fGM_{\bullet,\rm p}}{2a}=(f/K_{E}-1)\Delta
E.$ (6)

If $f>K_{E}$, no debris will return to the pericenter of secondary (this is
only roughly correct, as debris that can go through the Roche lobe of the
secondary can also approach the secondary again). After the disruption, the
bound debris with effective energy $E_{\rm eff}<0$ (and an apocenter less than
the Hill sphere) moves in approximately Keplerian orbits of the secondary. Due
to the impact from the primary, the fallback time scale and peak accretion
rate are modified by a factor of $(1-f/K_{E})^{-3/2}$ and $(1-f/K_{E})^{5/2}$,
respectively, compared with the disruption by single BH. We note that, this
analysis is consistent with the simulation of Coughlin & Armitage (2018), who
constructed 30 disruptions subjected to the secondary, where 23 yield no
accretion to the secondary, and 2 cases show fall back later or earlier than
the single BH disruption case.

If the debris that goes beyond the Hill sphere of secondary can be accreted by
primary, the most bound debris to the primary should have an energy of
$E\approx-(f+(9q)^{2/3})\Delta E/K_{E}=-\kappa\Delta E$ (note that $\Delta
E=GM_{\bullet,\rm s}R_{*}/r_{s,t}^{2}$, $\kappa=f/K_{E}+(9q)^{2/3}/K_{E}$ and
$q=M_{\bullet,\rm s}/M_{\bullet,\rm p}$). We assume the debris has been
accreted after it orbits the primary by one Keplerian orbit. The fallback rate
to the primary can be approximated as,

$\dot{M}_{p}=\frac{M_{*}\kappa}{3t_{\rm f}}\left(\frac{t+t_{\rm f}}{t_{\rm
f}}\right)^{-5/3}=\kappa^{5/2}q^{1/2}\dot{M}_{\rm p,f}.$ (7)

where $t_{\rm f}$ is the returning time of the most bound debris to the
primary, and $t_{\rm f}=t_{\rm s,f}\kappa^{-3/2}q^{-1}$ ($t_{\rm s,f}$ is the
fallback of the most bound debris to the secondary, defined as in Eq. 2).
Here, we estimate the peak accretion rate to the primary for case of
$\kappa=1$ (indicating no accretion to the secondary) and
$M_{\star}=M_{\odot}$. We get $\dot{m}_{p}\approx 0.02$ in Eddington mass
accretion rate units and $t_{\rm f}\approx 5000$ days for case of
$M_{\bullet}=10^{7.7}M_{\odot}$. This accretion rate results in a bolometric
luminosity of $L_{\rm d}\approx 1.4\times 10^{44}$ erg/s, which is brighter
than the rebrightening UV black body luminosity $L_{\rm BB}\approx
10^{43}$erg/s, assuming a radiative efficiency of $\eta=0.1$. Consequently, it
is possible that the rebrightening of the UV/optical emission originates from
accretion onto the primary. However, we have also demonstrated that the
rebrightened UV/optical emission closely follows the mass fallback law of
$t^{-5/3}$. If this UV/optical emission were primarily attributed to accretion
onto the primary, the fallback timescale should align with the rebrightening
date, at approximately 1200 days. Nevertheless, the predict fallback timescale
is significantly longer than this. Based on the aforementioned analysis, we
favor an explanation in which the rebrightening epochs are a consequence of
accretion onto the secondary. The UV/optical brightness during the early
stages of the rebrightening epochs surpasses that observed near the cutoff at
approximately $\sim 500$ days. This can be explain as, the debris’s angular
momentum, relative to the secondary, starts to differ from the material
initially comprising the disk, resulting in an higher accretion rate (Coughlin
et al., 2017). For TDE subjected to the secondary in a binary system with a
small $q$, simulations (Coughlin & Armitage, 2018) indicate that, in general,
there is negligible accretion onto the primary. If any accretion does occur,
it typically happens at very late times, and the rate is very low. The
associated timescale for the debris traveling directly to the primary falls in
the range of $T_{\rm acc}=a/v_{\rm escape}\approx 6$ years for this event.

By using Jacobi’s constant (Roth, 1952), which describes the total energy of
the debris, we can constrain $\kappa$ by estimating the mass accreted by the
secondary before the first light curves cutoff. As the star is disrupted by
the secondary at $r_{s,t}$, debris with Jacobi’s constant bigger than the
Jacobi’s constant at the Lagrange points L1 would not be able to travel
through the Roche lobe of the secondary. For small $q$, the Roche lobe can be
approximated as the Hill sphere radius, $R_{H}=(q/3)^{1/3}a$, and Jacobi’s
constant at L1 can be written in the form,

$C_{\rm L1}\approx\frac{GM_{\bullet,\rm
p}}{a}(1-R_{H}/a)^{2}+\frac{2GM_{\bullet,\rm
p}}{a-R_{H}}+\frac{2GM_{\bullet,\rm s}}{R_{H}}.$ (8)

The Jacobi’s constant of debris with a different energy E at the disrupted
position can be estimated as,

$C_{\rm E}\approx\frac{GM_{\bullet,\rm
p}}{a}(1-r_{s,t}/a)^{2}+\frac{2GM_{\bullet,\rm
p}}{a-r_{s,t}}+\frac{2GM_{\bullet,\rm s}}{r_{s,t}}-V_{\rm E}^{2}.$ (9)

We assume that debris with $C_{E}>C_{\rm L1}$ is all accreted before the first
light curves cutoff. As a result, the mass ratio between the debris accreted
by secondary and the mass of the star is,

$f_{m}=\frac{C_{E}-C_{\rm L1}}{4\Delta E}\approx\frac{1-\kappa}{2}.$ (10)

$f_{m}$ and $\kappa$ can be constrained from observation. By fitting the X-ray
spectra, we get that the mass accreted by secondary is $\approx
0.005$—$0.065M_{\odot}$ at $1\sigma$ CL, during the first 561 days. Note that
this mass is in good agreement with the
$M_{\star}=0.07^{+0.12}_{-0.04}M_{\odot}$ inferred from the MOSFiT light curve
fits. Assuming a solar mass star, we get $f_{m}$ is about 0.01–0.13 and
$\kappa$ is about 0.984–0.74. Taking $\kappa=0.74$, the values of
$\dot{m}_{\rm p}$ and $t_{\rm f}$ above change from 0.02 to 0.01 and from 5000
days to 7900 days, respectively.

###  4.4 Binary TDE Event Rate

An isolated BH in a galactic nucleus will tidally disrupt stars due to
diffusion in orbital angular momentum space (Frank & Rees, 1976). In this
picture, stars undergo a random walk in angular momentum from uncorrelated
two-body scatterings, eventually reaching pericenters less than the tidal
disruption radius. The outcomes of this random walk can be quantified in terms
of an angular momentum “loss cone,” which is emptied on a short dynamical
timescale and refilled on a slow relaxational timescale. The TDE rate in such
a system is thus equal to the loss-cone feeding rate produced by the two-body
relaxation process (Lightman & Shapiro, 1977; Cohn & Kulsrud, 1978).
Theoretical loss cone calculations generally find average TDE rates $\sim{\rm
few}\times 10^{-5}-{\rm few}\times 10^{-4}~{}{\rm yr}^{-1}$ (Magorrian &
Tremaine, 1999; Wang & Merritt, 2004; Stone & Metzger, 2016; Pfister et al.,
2020; Stone et al., 2020; Teboul et al., 2024), in some tension with (lower)
recent TDE rate estimates (Yao et al., 2023).

The formation and evolution of a bound SMBHB in a galaxy’s centers would
significantly change TDE event rates. Shortly after the formation of a hard
SMBHB, interactions between the binary and a bound stellar cusp will enhance
the TDE rate to as high as 1 yr-1, through a combination of the EKL mechanism
(Kozai, 1962; Lidov, 1962; Ivanov et al., 2005b) and strong three-body
scatterings (Chen et al., 2011; Wegg & Nate Bode, 2011). Mockler et al. (2023)
also recently showed that around the smaller black hole, a modest binary
eccentricity will allow the EKL mechanism (Ford et al., 2000; Naoz, 2016b) to
produce an even larger number of TDEs from the smaller secondary BH.

The combined effect of secular and strong interactions between the binary and
the surrounding dense stellar cusp will produce a number $N_{\rm b}\approx
7\times 10^{4}q^{(2-\gamma)/(6-2\gamma)}M_{7}^{11/12}$ TDEs (Chen et al.,
2011), where $1.5\leq\gamma\leq 2$ is density profile index. The $N_{\rm b}$
estimation of Chen et al. (2011) is in general agreement with the simulations
of Coughlin et al. (2017), who found $N_{\rm b}\approx 800-8000$ for case of
$q=0.2$ and $M_{\bullet,\rm p}=10^{6}M_{\odot}$ on circular SMBHBs. The
probability of stellar disruption by the secondary is about a factor of $q/2$
to that of primary in circular SMBHBs (Coughlin et al., 2017). Supposing that
a galaxy on average experiences N SMBHB episodes during a Hubble time
($10^{10}$ years), and assuming an observationally motivated TDE rate of
$7.4\times 10^{-5}~{}{\rm yr}^{-1}$ (Yao et al., 2023), then during one duty
cycle of galaxy merger, the relative numbers of TDEs by single BHs ($N_{\rm
s}$), SMBHB primaries ($N_{\rm b}$), and SMBHB secondaries ($N_{\rm bs}$) are
$N_{\rm s}:N_{\rm b}:N_{\rm bs}\approx
21:2q^{0.1}M_{7}^{11/12}N:q^{1.1}M_{7}^{11/12}N$ (assuming $\gamma=1.75$). For
the case of AT2018fyk, $q\approx 0.004$ and $M_{\bullet,\rm p}\approx
10^{7.7}M_{\odot}$, gives $N_{\rm s}:N_{\rm b}:N_{\rm bs}\approx 2100:500N:N$.
It is worth noting that this estimation introduces some uncertainty due to the
following factors: 1) Isolated SMBHs having higher TDE rates, leading to a
lower binary TDE rate; and 2) the EKL mechanism, resulting in more TDEs around
the secondary (Mockler et al., 2023), thereby contributing to a higher binary
TDE rate due to the secondary.

The TDE Catalog currently has 116 TDEs candidates (Goldtooth et al., 2023). If
we conservatively take $N=1$, there would be probabilities of $\approx 100\%$
and $\approx 4\%$ to detect at least 1 TDE from the primary and the secondary,
respectively, of an SMBHB out of a sample of this size. EKL mechanism can
produce TDEs around the secondary that are comparable to those around the
primary in the case of q=0.1 (Mockler et al., 2023). Assuming that the EKL
mechanism can increase the TDE rate by an order of magnitude for other q
values, the probability of detecting 1 TDE caused by the secondary out of 116
samples would then increase from $4\%$ to $31\%$. Currently, there are two
TDEs, SDSS J120136.02+300305.5 (Liu et al., 2014) and OGLE16aaa (Shu et al.,
2020), explained as SMBHB TDE originating from the primary. The low detected
rate of binary TDEs may be due to the fact that, it is not easy to distinguish
the binary TDE from the single SMBH TDE, as the perturbation from the
secondary is not strong enough to produce a significant difference in the
light curves. However, some factors may make the observed rated of binary TDE
due to secondary higher than predicted, e.g. the direct capture of the star by
the primary SMBH with $M_{\bullet}>10^{7}M_{\odot}$ due to relativistic
effects, would result in a substantial reduction in the rate of TDEs by the
primaries in SMBHBs (Coughlin & Nixon, 2022).

###  4.5 Distinguishing between SMBHB TDE and repeating pTDE scenarios

The rebrightening of this event has also been explained as a repeating pTDE
with $M_{\bullet}=10^{7.7}M_{\odot}$ (Wevers et al., 2023). In this scenario,
a star, on a $\sim$1200-day orbit around the SMBH, is periodically stripped of
mass during each pericenter passage, powering its unique light curves. As
shown above, the binary TDE hypothesis can also explain the observations well
and presents certain advantages over the repeating pTDE scenario. These
include: (i) the detection of a BH mass, through multiple distinct methods,
that is two orders of magnitude below predictions from host galaxy scaling
relations; (ii) the large-amplitude fluctuations and sharp cutoff in the light
curves; and (iii) the alignment of late-time mass fallback behavior with the
$t^{-5/3}$ law rather than the $t^{-9/4}$ law expected for pTDEs (Coughlin &
Nixon, 2019). However, none of these arguments decisively rules out the
possibility of a repeating pTDE. Galaxy scaling relations possess intrinsic
scatter (Kormendy & Ho, 2013) and major outliers do exist (van den Bosch et
al., 2012), a repeating pTDE model can still explain the light curve cutoff
(Wevers et al., 2023), and the light curve fluctuations make that the light
curve decay could be consistent with the $t^{-9/4}$ law (see Fig. 3). Here, we
list possible (future) observations that can distinguish between these two
scenarios, either for this TDE or for others:

  1. 1.

The time interval of the repeating rebrightening: Binary TDEs predict chaotic
rebrightening times, while the pTDE predicts the same rebrightening time
interval.

  2. 2.

The blue/redshift of the emission lines (Gaskell, 1983; De Rosa et al., 2019):
the binary TDE model predicts periodic blue/red-shift evolution of the
emission lines due to the Doppler effect when the secondary moves with respect
to our line of sight, while the pTDE does not.

  3. 3.

The ultimate fading of UV/optical emission: the binary TDE model predicts that
the emission will drop to the host level in about 5000 days, whereas pTDE
predicts that the flare would die out once the core stops being disrupted,
which can be on either much shorter or much longer timescales than 5000 days.

  4. 4.

Future gravitational wave detections with pulsar timing arrays (PTAs; e.g.
Agazie et al. 2023; Antoniadis et al. 2023) may help resolve SMBHB hosts of
TDEs with unusual light curve features, but this will only be possible for
nearby sources, high-mass primaries, and secondaries that are as high-mass as
possible. While the chirp mass that we infer for AT2018fyk is too low to ever
be detected as a resolvable PTA source (Sesana et al., 2009), a binary with
$M_{\bullet,\rm p}\sim 10^{9}M_{\odot}$ and $M_{\bullet,\rm s}\sim
10^{7}M_{\odot}$ would be more promising.

##  5 Conclusions

We conducted a comprehensive analysis of AT2018fyk, scrutinizing both its
X-ray and UV/optical emission. For the X-ray spectra, we employed a model
comprising a thermal disk and a thermally comptonized continuum. This allowed
us to analyse the source of the soft and hard X-ray spectral component and
place constraints on parameters such as $M_{\bullet}$ and $a_{\bullet}$
through X-ray spectral fitting. To understand the temporal behavior of the
source, we fitted the $w1$ light curve to determine its decay rate. We
examined the UV/optical light curves to investigate the evolution of the
black-body temperature, radius, and luminosity. We investigate the origins of
the early-time UV/optical emission and derived $M_{\bullet}$ from these data.
Next, we studied the cross-correlation between the X-ray, UV/optical, and soft
and hard X-ray emission to uncover potential connections and emission
mechanisms. Finally, our interpretation of AT2018fyk centers on it being a
supermassive binary TDE. We placed constraints on the binary parameters and
also list possible observations that can distinguish between the binary TDE
scenario and the repeating pTDE scenario. Our findings are:

  1. 1.

The fit to the X-ray spectra shows that the soft X-ray spectrum can be
explained as a disk spectrum, while the hard X-ray spectrum can be explained
as emission from a corona above the disk.

  2. 2.

The $M_{\bullet,\rm s}$ and $a_{\bullet}$ are constrained to be
$M_{\bullet,\rm s}=(2.7^{+0.5}_{-1.5})\times 10^{5}M_{\odot}$ and
$a_{\bullet}>0.3$ at the $1\sigma$ CL from slim disk modeling.

  3. 3.

Both the early and rebrightened $w1$ light curves decays as a power-law,
following the theoretical mass fallback rate. The $M_{\bullet,\rm s}$ is
estimated to be $(5.4^{+10}_{-4.9})\times 10^{5}M_{\odot}$ and
$(8.8^{+10}_{-7.7})\times 10^{5}M_{\odot}$, as determined from the early time
UV emission using the TDEmass and MOSFiT light curve modeling tools,
respectively. These values are consistent with that inferred from our X-ray
spectral modeling.

  4. 4.

The constant fitted BB temperature, slowly decreasing BB radius, and the
higher BB luminosity compared to the expected disk luminosity suggest that the
UV/optical emission originates from the circulation of debris streams. Our
analysis shows that the early BB luminosity has super-Eddington values.

  5. 5.

The evidence for a correlation between the soft and hard X-ray emission
supports the conclusion that the hard X-ray emission comes from the corona,
and the lack of significant correlation between the X-ray and UV/optical
emission support the idea that the UV/optical emission originates from the
circulation of the debris, whereas the X-ray emission comes from the inner
accretion disk.

  6. 6.

If the star is disrupted by the secondary of a binary SMBH, the rapid decline
in the light curve occurring approximately 561 days after its discovery can be
explained as the outcome of accretion onto the secondary experiencing a
cutoff. Utilizing this specific point in time, in conjunction with
$M_{\bullet,\rm s}$ and $M_{\bullet,\rm p}$, we constrain the separation of
this binary to be $(6.7\pm 1.2)\times 10^{-3}$ pc, and the corresponding
merger time to be $1.8^{+1.8}_{-0.9}\times 10^{9}$ years.

The behavior of X-ray and UV/optical emission is consistent with that of IMBH
TDEs, such as a high disk temperature, long super-Eddington time (Lin et al.,
2017, 2018; Wen et al., 2021; Cao et al., 2023), long UV power-law decay time,
and higher early UV/optical bolometric luminosity than the disk bolometric
luminosity. The $M_{\bullet}$ measured both from modeling the X-ray emission
and independently from modeling the UV emission, strongly supports that this
is an IMBH TDE. Notably, the $M_{\bullet}$ derived from TDE characteristics
appears to be incommensurate with that determined through the
$M_{\bullet}-\sigma_{\star}$ relation. We posit that this incongruity stems
from the disruption of a star by the secondary of a SMBHB. Considering the
protracted mass fallback time required for accretion onto the primary, we
attribute the rebrightening phase to accretion onto the secondary.
Consequently, we propose that both the early and rebrightening phases of
radiation emission can be ascribed to the accretion process onto the
secondary. The strong perturbation from the primary induces substantial
fluctuations and gaps in the mass fallback rate, thereby manifesting as
irregularities in the X-ray and UV/optical light curves.

While the properties of AT2018fyk appear to be consistent with a SMBHB TDE,
the inferred chirp mass is too low to be detected as a single PTA source.
Conversely, we cannot dismiss the possibility of a repeating pTDE based on
current observations alone. Additional observations are necessary to
invalidate the SMBHB, the pTDE hypothesis, or both. The SMBHB model predicts a
velocity change in the secondary BH of $\sim 5.5\times 10^{3}\rm{km~{}s^{-1}}$
over the $\approx$7.4 year binary orbital period. We attempted to detect
regular red-shifts or blue-shifts in the emission lines using optical spectra.
Unfortunately, only optical spectra taken during the peak of the initial flare
are available, and the object has since become too dim to for spectroscopic
ground based observations, making it challenging to verify this prediction
currently. However, in principle, if there is a third flare or outburst,
obtaining optical spectra might allow us to resolve any time-resolved offsets.

We thank F. Liu and E. R. Coughlin for the useful discussion about binary
TDEs. We thank B. Mockler for guidance in utilizing the MOSFiT package. This
work made use of data supplied by the UK Swift Science Data Centre at the
University of Leicester. This research has made use of data obtained through
the High Energy Astrophysics Science Archive Research Center Online Service,
provided by the NASA/Goddard Space Flight Center. Our work here is partly
based on observations obtained with XMM-Newton, an ESA science mission with
instruments and contributions directly funded by ESA Member States and NASA.
SW thanks the RU Department of Astrophysics/IMAPP for postdoctoral support
during the early part of this work. SW is also supported by the NSFC (grant
12333004), the Strategic Priority Research Program of the Chinese Academy of
Sciences (grant XDB0550200), and the Strategic Pioneer Program on Space
Science, Chinese Academy of Sciences (grant XDA15310300). AIZ acknowledges
support from NASA ADAP grant 80NSSC21K0988. PGJ has received funding from the
European Research Council (ERC) under the European Union’s Horizon 2020
research and innovation programme (Grant agreement No. 101095973).

##  Appendix A Constraints of $M_{\bullet}$ and $a_{\bullet}$ using different
X-ray spectra

We fit the spectra Swift 1, XMM 1, XMM 2 and XMM 4, separately. For Swift 1
data, we fit the spectrum with a fit-function ${\tt TBAbs\times slimd}$, while
for the three XMM spectra, we fit them with a fit-function ${\tt TBAbs\times
Thcomp(slimd)}$. For the model Thcomp, we fix the electron temperature to a
value of 300 keV. We derive the $M_{\bullet}$ and $a_{\bullet}$ contours for
each spectrum with the command steppar in XSPEC (Arnaud, 1996).

Figure 6: Constraints on $M_{\bullet}$ and $a_{\bullet}$ obtained by fitting
different spectra individually. The dark and light blue regions denote the 1
and $2\sigma$ contours, respectively. XMM 2 yields a larger value for
$M_{\bullet}$ than the three other epochs at the $>2\sigma$ level.

Figure 6 shows the fitting constraints on $M_{\bullet}$ and $a_{\bullet}$. The
fits to the individual XMM spectra do not constrain $a_{\bullet}$ as much as
the fit to the Swift 1 spectrum does, even though the individual XMM spectra
have much more photons. This is because in the XMM epochs the corona modifies
the spectrum, which is cloaking some information from the inner part of the
disk. The constraints on $M_{\bullet}$ and $a_{\bullet}$ from XMM 2 are at
least $2\sigma$ away from that of Swift 1\. This tension may be explained by
the possible soft lags in the XMM 2 as report by Zhang (2022). This tension
becomes smaller for XMM 4\. This may be because the outflowing material
causing the soft lag disappears over time, or just because the there are not
enough counts to give a robust constraints on $M_{\bullet}$ and $a_{\bullet}$
in the XMM 4 spectrum.

##  Appendix B Introduction of the slimd model.

The slimd model (Wen et al., 2022) is based on the work of Wen et al. (2020)
and Wen et al. (2021), who used a general relativistic stationary slim disk
model (Abramowicz et al., 1996; Sądowski, 2009) employing the ray-tracing code
from Psaltis & Johannsen (2012) to calculate the synthetic spectrum as seen by
an observer. The model assumptions include: 1. A circular and equatorial inner
disk. 2. A zero-torque inner boundary condition, large optical depth, and
absence of self-irradiation. 3. Generation of X-rays occur at the finite-
height photospheric surface of the disk, presumed to match its scale height.
The X-ray spectrum is determined by BH mass and spin, the mass accretion rate,
the system inclination, and the spectral hardening factor $f_{c}$. The value
of the latter parameter is determined through radiative transfer calculations
by Davis & El-Abd (2019). Studies have indicated that different choices of the
hardening factor have a limited impact on the estimation of mass and spin (see
figures 14 and 15 in Wen et al. 2020). To study the constraining power on disk
parameters of the thermal spectrum, we fit the Swift 1 spectrum with the slimd
model.

Figure 7: Constraints on different parameters obtained from fitting the Swift
1. We fit the spectrum with the function ${\tt TBAbs\times slimd}$ by setting
$N_{H}=1.15\times 10^{20}\rm{cm^{-2}}$ while allowing the other parameters
(including $M_{\bullet}$, $a_{\bullet}$, $\dot{m}$, and $\theta$) to float
freely. The dark and light blue regions denote the 1 and $2\sigma$ contours,
respectively.

Figure 7 shows the degeneracy between disk parameters and the BH mass and
spin. The BH mass and spin strongly (anti)correlate with the mass accretion
rate and system inclination, respectively. There is no constraint on the
inclination, and a poor constraint on the accretion rate. Although the
uncertainty in the best-fit value of $\theta$ and $\dot{m}$ is large, the
$M_{\bullet}$ and $a_{\bullet}$ values are relatively well constrained, as
shown in the upper left panel of figure 6. Overall, the BH mass is better
constrained than the BH spin. Due to the degeneracy between mass and spin, the
constraint on spin is weak. However, the degeneracy between mass and spin can
be partially alleviated by simultaneously fitting multi-epoch spectra with
varying accretion rates (Wen et al., 2022). These results indicate that the
X-ray spectra are mainly determined by the BH mass and accretion rate and to a
lesser degree BH spin. Spectra with super-Eddington accretion rates have a
stronger power in constraining $M_{\bullet}$ and $a_{\bullet}$, because then
the flux is insensitive to the actual value of the mass accretion rate (Wen et
al., 2022).

## References

  * Abramowicz (2005) Abramowicz, M. 2005, in Growing Black Holes: Accretion in a Cosmological Context, ed. A. Merloni, S. Nayakshin, & R. A. Sunyaev (Berlin, Heidelberg: Springer Berlin Heidelberg), 257–273
  * Abramowicz et al. (1996) Abramowicz, M. A., Chen, X. M., Granath, M., & Lasota, J. P. 1996, ApJ, 471, 762, doi: 10.1086/178004
  * Abramowicz et al. (1988) Abramowicz, M. A., Czerny, B., Lasota, J. P., & Szuszkiewicz, E. 1988, ApJ, 332, 646, doi: 10.1086/166683
  * Agazie et al. (2023) Agazie, G., Anumarlapudi, A., Archibald, A. M., et al. 2023, ApJ, 951, L8, doi: 10.3847/2041-8213/acdac6
  * Akaike (1974) Akaike, H. 1974, IEEE Transactions on Automatic Control, 19, 716
  * Amaro-Seoane et al. (2023) Amaro-Seoane, P., Andrews, J., Arca Sedda, M., et al. 2023, Living Reviews in Relativity, 26, 2, doi: 10.1007/s41114-022-00041-y
  * Andalman et al. (2022) Andalman, Z. L., Liska, M. T. P., Tchekhovskoy, A., Coughlin, E. R., & Stone, N. 2022, MNRAS, 510, 1627, doi: 10.1093/mnras/stab3444
  * Antoniadis et al. (2023) Antoniadis, J., Arumugam, P., Arumugam, S., et al. 2023, arXiv e-prints, arXiv:2306.16214, doi: 10.48550/arXiv.2306.16214
  * Arcavi et al. (2014) Arcavi, I., Gal-Yam, A., Sullivan, M., et al. 2014, ApJ, 793, 38, doi: 10.1088/0004-637X/793/1/38
  * Arnaud (1996) Arnaud, K. A. 1996, in Astronomical Society of the Pacific Conference Series, Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes, 17
  * Begelman et al. (1980) Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980, Nature, 287, 307, doi: 10.1038/287307a0
  * Blandford & Rees (1974) Blandford, R. D., & Rees, M. J. 1974, MNRAS, 169, 395, doi: 10.1093/mnras/169.3.395
  * Bonnerot & Lu (2020) Bonnerot, C., & Lu, W. 2020, MNRAS, 495, 1374, doi: 10.1093/mnras/staa1246
  * Cannizzo et al. (1990) Cannizzo, J. K., Lee, H. M., & Goodman, J. 1990, ApJ, 351, 38, doi: 10.1086/168442
  * Cao et al. (2023) Cao, Z., Jonker, P. G., Wen, S., Stone, N. C., & Zabludoff, A. I. 2023, MNRAS, 519, 2375, doi: 10.1093/mnras/stac3539
  * Cash (1979) Cash, W. 1979, ApJ, 228, 939, doi: 10.1086/156922
  * Chen et al. (2011) Chen, X., Sesana, A., Madau, P., & Liu, F. K. 2011, ApJ, 729, 13, doi: 10.1088/0004-637X/729/1/13
  * Chiang et al. (2015) Chiang, C.-Y., Walton, D. J., Fabian, A. C., Wilkins, D. R., & Gallo, L. C. 2015, MNRAS, 446, 759, doi: 10.1093/mnras/stu2087
  * Cohn & Kulsrud (1978) Cohn, H., & Kulsrud, R. M. 1978, ApJ, 226, 1087, doi: 10.1086/156685
  * Coughlin & Armitage (2018) Coughlin, E. R., & Armitage, P. J. 2018, MNRAS, 474, 3857, doi: 10.1093/mnras/stx3039
  * Coughlin et al. (2019) Coughlin, E. R., Armitage, P. J., Lodato, G., & Nixon, C. J. 2019, Space Sci. Rev., 215, 45, doi: 10.1007/s11214-019-0612-z
  * Coughlin et al. (2017) Coughlin, E. R., Armitage, P. J., Nixon, C., & Begelman, M. C. 2017, MNRAS, 465, 3840, doi: 10.1093/mnras/stw2913
  * Coughlin & Nixon (2019) Coughlin, E. R., & Nixon, C. J. 2019, ApJ, 883, L17, doi: 10.3847/2041-8213/ab412d
  * Coughlin & Nixon (2022) —. 2022, ApJ, 936, 70, doi: 10.3847/1538-4357/ac85b3
  * Dai et al. (2015) Dai, L., McKinney, J. C., & Miller, M. C. 2015, ApJ, 812, L39, doi: 10.1088/2041-8205/812/2/L39
  * Dai et al. (2018) Dai, L., McKinney, J. C., Roth, N., Ramirez-Ruiz, E., & Miller, M. C. 2018, ApJ, 859, L20, doi: 10.3847/2041-8213/aab429
  * Davis & El-Abd (2019) Davis, S. W., & El-Abd, S. 2019, ApJ, 874, 23, doi: 10.3847/1538-4357/ab05c5
  * De Rosa et al. (2019) De Rosa, A., Vignali, C., Bogdanović, T., et al. 2019, New A Rev., 86, 101525, doi: 10.1016/j.newar.2020.101525
  * Donato et al. (2014) Donato, D., Cenko, S. B., Covino, S., et al. 2014, ApJ, 781, 59, doi: 10.1088/0004-637X/781/2/59
  * Done & Gierliński (2005) Done, C., & Gierliński, M. 2005, MNRAS, 364, 208, doi: 10.1111/j.1365-2966.2005.09555.x
  * Evans et al. (2009) Evans, P. A., Beardmore, A. P., Page, K. L., et al. 2009, Monthly Notices of the Royal Astronomical Society, 397, 1177, doi: 10.1111/j.1365-2966.2009.14913.x
  * Ferrarese & Ford (2005) Ferrarese, L., & Ford, H. 2005, Space Sci. Rev., 116, 523, doi: 10.1007/s11214-005-3947-6
  * Ford et al. (2000) Ford, E. B., Kozinsky, B., & Rasio, F. A. 2000, ApJ, 535, 385, doi: 10.1086/308815
  * Frank & Rees (1976) Frank, J., & Rees, M. J. 1976, MNRAS, 176, 633, doi: 10.1093/mnras/176.3.633
  * French et al. (2016) French, K. D., Arcavi, I., & Zabludoff, A. 2016, ApJ, 818, L21, doi: 10.3847/2041-8205/818/1/L21
  * Gaskell (1983) Gaskell, C. M. 1983, in Liege International Astrophysical Colloquia, Vol. 24, Liege International Astrophysical Colloquia, ed. J.-P. Swings, 473–477
  * Gelman & Rubin (1992) Gelman, A., & Rubin, D. B. 1992, Statistical Science, 7, 457, doi: 10.1214/ss/1177011136
  * Goldtooth et al. (2023) Goldtooth, A., Zabludoff, A. I., Wen, S., et al. 2023, PASP, 135, 034101, doi: 10.1088/1538-3873/acb9bc
  * Gong et al. (2021) Gong, Y., Luo, J., & Wang, B. 2021, Nature Astronomy, 5, 881, doi: 10.1038/s41550-021-01480-3
  * Gordon et al. (2003) Gordon, K. D., Clayton, G. C., Misselt, K. A., Landolt, A. U., & Wolff, M. J. 2003, ApJ, 594, 279, doi: 10.1086/376774
  * Graur et al. (2018) Graur, O., French, K. D., Zahid, H. J., et al. 2018, ApJ, 853, 39, doi: 10.3847/1538-4357/aaa3fd
  * Greene et al. (2016) Greene, J. E., Seth, A., Kim, M., et al. 2016, ApJ, 826, L32, doi: 10.3847/2041-8205/826/2/L32
  * Guillochon et al. (2014) Guillochon, J., Manukian, H., & Ramirez-Ruiz, E. 2014, ApJ, 783, 23, doi: 10.1088/0004-637x/783/1/23
  * Guillochon & Ramirez-Ruiz (2013) Guillochon, J., & Ramirez-Ruiz, E. 2013, ApJ, 767, 25, doi: 10.1088/0004-637X/767/1/25
  * Guo et al. (2023) Guo, H., Sun, J., Li, S.-L., et al. 2023, arXiv e-prints, arXiv:2312.06771, doi: 10.48550/arXiv.2312.06771
  * Hammerstein et al. (2023) Hammerstein, E., van Velzen, S., Gezari, S., et al. 2023, ApJ, 942, 9, doi: 10.3847/1538-4357/aca283
  * Häring & Rix (2004) Häring, N., & Rix, H.-W. 2004, ApJ, 604, L89, doi: 10.1086/383567
  * Hills (1975) Hills, J. G. 1975, Nature, 254, 295, doi: 10.1038/254295a0
  * Ivanov et al. (2005a) Ivanov, P. B., Polnarev, A. G., & Saha, P. 2005a, MNRAS, 358, 1361, doi: 10.1111/j.1365-2966.2005.08843.x
  * Ivanov et al. (2005b) —. 2005b, MNRAS, 358, 1361, doi: 10.1111/j.1365-2966.2005.08843.x
  * Jansen et al. (2001) Jansen, F., Lumb, D., Altieri, B., et al. 2001, A&A, 365, L1, doi: 10.1051/0004-6361:20000036
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Koss et al. (2022) Koss, M. J., Trakhtenbrot, B., Ricci, C., et al. 2022, ApJS, 261, 6, doi: 10.3847/1538-4365/ac650b
  * Kozai (1962) Kozai, Y. 1962, AJ, 67, 591, doi: 10.1086/108790
  * Lidov (1962) Lidov, M. L. 1962, Planet. Space Sci., 9, 719, doi: 10.1016/0032-0633(62)90129-0
  * Lightman & Shapiro (1977) Lightman, A. P., & Shapiro, S. L. 1977, ApJ, 211, 244, doi: 10.1086/154925
  * Lin et al. (2017) Lin, D., Guillochon, J., Komossa, S., et al. 2017, Nature Astronomy, 1, 0033, doi: 10.1038/s41550-016-0033
  * Lin et al. (2018) Lin, D., Strader, J., Carrasco, E. R., et al. 2018, Nature Astronomy, 2, 656, doi: 10.1038/s41550-018-0493-1
  * Liska et al. (2022) Liska, M. T. P., Musoke, G., Tchekhovskoy, A., Porth, O., & Beloborodov, A. M. 2022, ApJ, 935, L1, doi: 10.3847/2041-8213/ac84db
  * Liu et al. (2021) Liu, F. K., Cao, C. Y., Abramowicz, M. A., et al. 2021, ApJ, 908, 179, doi: 10.3847/1538-4357/abd2b6
  * Liu et al. (2009) Liu, F. K., Li, S., & Chen, X. 2009, ApJ, 706, L133, doi: 10.1088/0004-637X/706/1/L133
  * Liu et al. (2014) Liu, F. K., Li, S., & Komossa, S. 2014, ApJ, 786, 103, doi: 10.1088/0004-637X/786/2/103
  * Liu et al. (2022) Liu, X.-L., Dou, L.-M., Chen, J.-H., & Shen, R.-F. 2022, ApJ, 925, 67, doi: 10.3847/1538-4357/ac33a9
  * Lodato et al. (2009) Lodato, G., Nayakshin, S., King, A. R., & Pringle, J. E. 2009, MNRAS, 398, 1392, doi: 10.1111/j.1365-2966.2009.15179.x
  * Loeb & Ulmer (1997) Loeb, A., & Ulmer, A. 1997, ApJ, 489, 573, doi: 10.1086/304814
  * Lu & Bonnerot (2019) Lu, W., & Bonnerot, C. 2019, MNRAS, 492, 686, doi: 10.1093/mnras/stz3405
  * Magorrian & Tremaine (1999) Magorrian, J., & Tremaine, S. 1999, MNRAS, 309, 447, doi: 10.1046/j.1365-8711.1999.02853.x
  * Mayer et al. (2007) Mayer, L., Kazantzidis, S., Madau, P., et al. 2007, Science, 316, 1874, doi: 10.1126/science.1141858
  * McConnell & Ma (2013) McConnell, N. J., & Ma, C.-P. 2013, ApJ, 764, 184, doi: 10.1088/0004-637X/764/2/184
  * Merritt & Milosavljević (2005) Merritt, D., & Milosavljević, M. 2005, Living Reviews in Relativity, 8, 8, doi: 10.12942/lrr-2005-8
  * Merritt & Poon (2004) Merritt, D., & Poon, M. Y. 2004, ApJ, 606, 788, doi: 10.1086/382497
  * Metzger (2022) Metzger, B. D. 2022, ApJ, 937, L12, doi: 10.3847/2041-8213/ac90ba
  * Metzger & Stone (2016) Metzger, B. D., & Stone, N. C. 2016, MNRAS, 461, 948, doi: 10.1093/mnras/stw1394
  * Miller et al. (2015) Miller, J. M., Kaastra, J. S., Miller, M. C., et al. 2015, Nature, 526, 542, doi: 10.1038/nature15708
  * Milosavljević & Merritt (2003) Milosavljević, M., & Merritt, D. 2003, ApJ, 596, 860, doi: 10.1086/378086
  * Mockler et al. (2019) Mockler, B., Guillochon, J., & Ramirez-Ruiz, E. 2019, ApJ, 872, 151, doi: 10.3847/1538-4357/ab010f
  * Mockler et al. (2023) Mockler, B., Melchor, D., Naoz, S., & Ramirez-Ruiz, E. 2023, arXiv e-prints, arXiv:2306.05510, doi: 10.48550/arXiv.2306.05510
  * Mummery (2021) Mummery, A. 2021, Monthly Notices of the Royal Astronomical Society: Letters, 507, L24, doi: 10.1093/mnrasl/slab088
  * Mummery et al. (2023a) Mummery, A., van Velzen, S., Nathan, E., et al. 2023a, MNRAS, doi: 10.1093/mnras/stad3001
  * Mummery et al. (2023b) Mummery, A., Wevers, T., Saxton, R., & Pasham, D. 2023b, MNRAS, 519, 5828, doi: 10.1093/mnras/stac3798
  * Naoz (2016a) Naoz, S. 2016a, ARA&A, 54, 441, doi: 10.1146/annurev-astro-081915-023315
  * Naoz (2016b) —. 2016b, ARA&A, 54, 441, doi: 10.1146/annurev-astro-081915-023315
  * Pasham et al. (2020) Pasham, D., Gendreau, K., Remillard, R., et al. 2020, The Astronomer’s Telegram, 13864, 1
  * Pasham et al. (2021) Pasham, D. R., Ho, W. C. G., Alston, W., et al. 2021, Nature Astronomy, 6, 249, doi: 10.1038/s41550-021-01524-8
  * Peters (1964) Peters, P. C. 1964, Physical Review, 136, 1224, doi: 10.1103/PhysRev.136.B1224
  * Peterson et al. (1998) Peterson, B. M., Wanders, I., Horne, K., et al. 1998, PASP, 110, 660, doi: 10.1086/316177
  * Peterson et al. (2004) Peterson, B. M., Ferrarese, L., Gilbert, K. M., et al. 2004, ApJ, 613, 682, doi: 10.1086/423269
  * Pfister et al. (2020) Pfister, H., Volonteri, M., Dai, J. L., & Colpi, M. 2020, MNRAS, 497, 2276, doi: 10.1093/mnras/staa1962
  * Phinney (1989) Phinney, E. S. 1989, in IAU Symposium, Vol. 136, The Center of the Galaxy, ed. M. Morris, 543, doi: 10.1017/S0074180900187054
  * Piran et al. (2015) Piran, T., Svirski, G., Krolik, J., Cheng, R. M., & Shiokawa, H. 2015, ApJ, 806, 164, doi: 10.1088/0004-637x/806/2/164
  * Planck Collaboration et al. (2020) Planck Collaboration, Aghanim, N., Akrami, Y., et al. 2020, A&A, 641, A6, doi: 10.1051/0004-6361/201833910
  * Psaltis & Johannsen (2012) Psaltis, D., & Johannsen, T. 2012, ApJ, 745, 1, doi: 10.1088/0004-637X/745/1/1
  * Rees (1988) Rees, M. J. 1988, Nature, 333, 523, doi: 10.1038/333523a0
  * Remillard & McClintock (2006) Remillard, R. A., & McClintock, J. E. 2006, ARA&A, 44, 49, doi: 10.1146/annurev.astro.44.051905.092532
  * Robertson et al. (2015) Robertson, D. R. S., Gallo, L. C., Zoghbi, A., & Fabian, A. C. 2015, MNRAS, 453, 3455, doi: 10.1093/mnras/stv1575
  * Roming et al. (2005) Roming, P. W. A., Kennedy, T. E., Mason, K. O., et al. 2005, Space Sci. Rev., 120, 95, doi: 10.1007/s11214-005-5095-4
  * Roth (1952) Roth, H. 1952, AJ, 57, 42, doi: 10.1086/106702
  * Roth et al. (2016) Roth, N., Kasen, D., Guillochon, J., & Ramirez-Ruiz, E. 2016, ApJ, 827, 3, doi: 10.3847/0004-637x/827/1/3
  * Roth et al. (2020) Roth, N., Rossi, E. M., Krolik, J., et al. 2020, Space Sci. Rev., 216, 114, doi: 10.1007/s11214-020-00735-1
  * Ryu et al. (2020) Ryu, T., Krolik, J., & Piran, T. 2020, ApJ, 904, 73, doi: 10.3847/1538-4357/abbf4d
  * Ryu et al. (2023) Ryu, T., Krolik, J., Piran, T., Noble, S. C., & Avara, M. 2023, ApJ, 957, 12, doi: 10.3847/1538-4357/acf5de
  * Sadowsk (2009) Sadowsk, A. 2009, ApJS, 183, 171, doi: 10.1088/0067-0049/183/2/171
  * Sarin & Metzger (2023) Sarin, N., & Metzger, B. D. 2023, arXiv e-prints, arXiv:2307.15121, doi: 10.48550/arXiv.2307.15121
  * Saxton et al. (2017) Saxton, R. D., Read, A. M., Komossa, S., et al. 2017, A&A, 598, A29, doi: 10.1051/0004-6361/201629015
  * Schlegel et al. (1998) Schlegel, D. J., Finkbeiner, D. P., & Davis, M. 1998, ApJ, 500, 525, doi: 10.1086/305772
  * Sesana et al. (2009) Sesana, A., Vecchio, A., & Volonteri, M. 2009, MNRAS, 394, 2255, doi: 10.1111/j.1365-2966.2009.14499.x
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337
  * Shimura & Takahara (1995) Shimura, T., & Takahara, F. 1995, ApJ, 445, 780, doi: 10.1086/175740
  * Shiokawa et al. (2015) Shiokawa, H., Krolik, J. H., Cheng, R. M., Piran, T., & Noble, S. C. 2015, ApJ, 804, 85, doi: 10.1088/0004-637x/804/2/85
  * Shu et al. (2020) Shu, X., Zhang, W., Li, S., et al. 2020, Nature Communications, 11, 5876, doi: 10.1038/s41467-020-19675-z
  * Sądowski (2009) Sądowski, A. 2009, ApJS, 183, 171, doi: 10.1088/0067-0049/183/2/171
  * Steinberg & Stone (2022) Steinberg, E., & Stone, N. C. 2022, arXiv e-prints, arXiv:2206.10641, doi: 10.48550/arXiv.2206.10641
  * Stone & Metzger (2016) Stone, N. C., & Metzger, B. D. 2016, MNRAS, 455, 859, doi: 10.1093/mnras/stv2281
  * Stone et al. (2020) Stone, N. C., Vasiliev, E., Kesden, M., et al. 2020, Space Sci. Rev., 216, 35, doi: 10.1007/s11214-020-00651-4
  * Strüder et al. (2001) Strüder, L., Briel, U., Dennerl, K., et al. 2001, A&A, 365, L18, doi: 10.1051/0004-6361:20000066
  * Sun et al. (2018) Sun, M., Grier, C. J., & Peterson, B. M. 2018, PyCCF: Python Cross Correlation Function for reverberation mapping studies, Astrophysics Source Code Library, record ascl:1805.032.  http://ascl.net/1805.032
  * Teboul et al. (2024) Teboul, O., Stone, N. C., & Ostriker, J. P. 2024, MNRAS, 527, 3094, doi: 10.1093/mnras/stad3301
  * Timmer & Koenig (1995) Timmer, J., & Koenig, M. 1995, A&A, 300, 707
  * Turner et al. (2001) Turner, M. J. L., Abbey, A., Arnaud, M., et al. 2001, A&A, 365, L27, doi: 10.1051/0004-6361:20000087
  * Ulmer (1999) Ulmer, A. 1999, ApJ, 514, 180, doi: 10.1086/306909
  * Uttley et al. (2014) Uttley, P., Cackett, E. M., Fabian, A. C., Kara, E., & Wilkins, D. R. 2014, A&A Rev., 22, 72, doi: 10.1007/s00159-014-0072-0
  * van den Bosch et al. (2012) van den Bosch, R. C. E., Gebhardt, K., Gültekin, K., et al. 2012, Nature, 491, 729, doi: 10.1038/nature11592
  * van Velzen et al. (2011) van Velzen, S., Körding, E., & Falcke, H. 2011, Monthly Notices of the Royal Astronomical Society: Letters, 417, L51, doi: 10.1111/j.1745-3933.2011.01118.x
  * van Velzen et al. (2019) van Velzen, S., Stone, N. C., Metzger, B. D., et al. 2019, ApJ, 878, 82, doi: 10.3847/1538-4357/ab1844
  * van Velzen et al. (2021) van Velzen, S., Gezari, S., Hammerstein, E., et al. 2021, ApJ, 908, 4, doi: 10.3847/1538-4357/abc258
  * Vasiliev et al. (2015) Vasiliev, E., Antonini, F., & Merritt, D. 2015, ApJ, 810, 49, doi: 10.1088/0004-637X/810/1/49
  * Wang & Merritt (2004) Wang, J., & Merritt, D. 2004, ApJ, 600, 149, doi: 10.1086/379767
  * Wegg & Nate Bode (2011) Wegg, C., & Nate Bode, J. 2011, ApJ, 738, L8, doi: 10.1088/2041-8205/738/1/L8
  * Welsh (1999) Welsh, W. F. 1999, PASP, 111, 1347, doi: 10.1086/316457
  * Wen et al. (2023) Wen, S., Jonker, P. G., Stone, N. C., Van Velzen, S., & Zabludoff, A. I. 2023, MNRAS, 522, 1155, doi: 10.1093/mnras/stad991
  * Wen et al. (2021) Wen, S., Jonker, P. G., Stone, N. C., & Zabludoff, A. I. 2021, ApJ, 918, 46, doi: 10.3847/1538-4357/ac00b5
  * Wen et al. (2022) Wen, S., Jonker, P. G., Stone, N. C., Zabludoff, A. I., & Cao, Z. 2022, ApJ, 933, 31, doi: 10.3847/1538-4357/ac70c5
  * Wen et al. (2020) Wen, S., Jonker, P. G., Stone, N. C., Zabludoff, A. I., & Psaltis, D. 2020, ApJ, 897, 80, doi: 10.3847/1538-4357/ab9817
  * Wevers (2020) Wevers, T. 2020, MNRAS, 497, L1, doi: 10.1093/mnrasl/slaa097
  * Wevers et al. (2019) Wevers, T., Pasham, D. R., van Velzen, S., et al. 2019, MNRAS, 488, 4816, doi: 10.1093/mnras/stz1976
  * Wevers et al. (2021) Wevers, T., Pasham, D. R., van Velzen, S., et al. 2021, ApJ, 912, 151, doi: 10.3847/1538-4357/abf5e2
  * Wevers et al. (2023) Wevers, T., Coughlin, E. R., Pasham, D. R., et al. 2023, The Astrophysical Journal Letters, 942, L33, doi: 10.3847/2041-8213/ac9f36
  * Yao et al. (2022) Yao, Y., Lu, W., Guolo, M., et al. 2022, ApJ, 937, 8, doi: 10.3847/1538-4357/ac898a
  * Yao et al. (2023) Yao, Y., Ravi, V., Gezari, S., et al. 2023, ApJ, 955, L6, doi: 10.3847/2041-8213/acf216
  * Zanazzi & Ogilvie (2020) Zanazzi, J. J., & Ogilvie, G. I. 2020, MNRAS, 499, 5562, doi: 10.1093/mnras/staa3127
  * Zdziarski et al. (2020) Zdziarski, A. A., Szanecki, M., Poutanen, J., Gierliński, M., & Biernacki, P. 2020, MNRAS, 492, 5234, doi: 10.1093/mnras/staa159
  * Zhang (2022) Zhang, W. 2022, MNRAS, 511, 19, doi: 10.1093/mnrasl/slab133