# Post-Merger Jets from Supermassive Black Hole Coalescences as
Electromagnetic Counterparts of Gravitational Wave Emission

Chengchao Yuan11affiliation: Department of Physics, The Pennsylvania State
University, University Park, PA 16802, USA 22affiliation: Department of
Astronomy & Astrophysics, The Pennsylvania State University, University Park,
PA 16802, USA 33affiliation: Center for Multimessenger Astrophysics, Institute
for Gravitation and the Cosmos, The Pennsylvania State University, University
Park, PA 16802, USA<EMAIL_ADDRESS>Kohta Murase11affiliationmark:
22affiliationmark:  33affiliationmark:  44affiliation: Center for
Gravitational Physics, Yukawa Institute for Theoretical Physics, Kyoto
University, Kyoto, Kyoto 606-8502, Japan<EMAIL_ADDRESS>B. Theodore
Zhang11affiliationmark:  22affiliationmark:  33affiliationmark:  Shigeo S.
Kimura55affiliation: Frontier Research Institute for Interdisciplinary
Sciences; Astronomical Institute, Tohoku University, Sendai 980-8578, Japan
Peter Mészáros11affiliationmark:  22affiliationmark:  33affiliationmark:

###### Abstract

As a powerful source of gravitational waves (GW), a supermassive black hole
(SMBH) merger may be accompanied by a relativistic jet that leads to
detectable electromagnetic (EM) emission. We model the propagation of post-
merger jets inside a pre-merger wind bubble formed by disk-driven winds, and
calculate multi-wavelength EM spectra from the forward shock region. We show
that the non-thermal EM signals from SMBH mergers are detectable up to the
detection horizon of future GW facilities such as the Laser Interferometer
Space Antenna (LISA). Calculations based on our model predict slowly fading
transients with time delays from days to months after the coalescence, leading
to implications for EM follow-up observations after the GW detection.

jets, non-thermal, supermassive black holes

##  1 Introduction

Supermassive black hole (SMBH) mergers are ubiquitous in the history of the
Universe (Begelman et al., 1980; Kormendy & Ho, 2013; Kroupa et al., 2020) and
can produce powerful gravitational wave (GW) bursts when they coalesce (e.g.,
Thorne & Braginskii, 1976; Sesana et al., 2004), making them promising
candidates for GW detectors such as Laser Interferometer Space Antenna (LISA,
Amaro-Seoane et al., 2017; Baker et al., 2019) and pulsar timing arrays (PTAs,
e.g., Mingarelli et al., 2017; Taylor et al., 2019; Arzoumanian et al., 2020)
in single-source and/or stochastic GW background searches. The accretion
activity between the binary system and the surrounding disk can produce multi-
wavelength electromagnetic (EM) emission (e.g., Milosavljević & Phinney, 2005;
Moesta et al., 2012; Farris et al., 2015; Kelly et al., 2017; Haiman, 2017;
d’Ascoli et al., 2018), and the time-variable EM signatures from the
circumbinary disks could be detectable (e.g., Schnittman & Krolik, 2008;
Haiman et al., 2009; Tanaka & Menou, 2010). The spinning SMBH expected to form
after the SMBHs have coalesced may also lead to relativistic jets, in which
particle acceleration will take place. The resulting non-thermal emission from
the accelerated electrons may provide a promising post-merger EM counterpart
of the GW emission, and will not only provide complementary information on
SMBH mergers but also shed light on the physical processes in these systems
(e.g., Schnittman, 2011; Ravi, 2018; Mészáros et al., 2019). Yuan et al.
(2020) recently suggested that the SMBH mergers can also be high-energy
neutrino emitters, and demonstrated that they are also promising targets for
high-energy multi-messenger astrophysics (Murase & Bartos, 2019).

We study the EM emission produced in relativistic jets launched after the
coalescence of SMBHs. The physical picture is that the disk winds originating
from the circumbinary disk and mini-disks around each SMBH form a pre-merger
wind bubble, and jets powered by the Blandford-Znajek (BZ, Blandford & Znajek,
1977) mechanism are launched after the merger. The jets push ahead inside the
pre-merger disk wind material, resulting in the formation of forward and
reverse shocks. In the forward shock region, electrons are accelerated to high
energies with a power-law distribution as observed in afterglows of gamma-ray
bursts (GRBs) (e.g., Mészáros, 2006). These particles then produce broadband
non-thermal EM emission through synchrotron and synchrotron self-Compton (SSC)
processes.

This letter is organized as follows. In Sec. 2, we introduce the physical
conditions of the pre-merger wind bubble and model the propagation of jets.
The radiation processes and the resulting photon spectra, light curves and
detection horizons are presented in Sec. 3. In Sec. 4, we discuss implications
of our results. Throughout the letter, we use the conventional notation
$Q_{x}=Q/10^{x}$ and physical quantities are written in the centimeter-gram-
second units, unless otherwise specified.

##  2 Jet dynamics

We discuss here the physical conditions in a pre-merger circumbinary
environment and derive relevant quantities that describe the jet propagation.
We consider on-axis observers, which is sufficient for the purpose of this
work. The emission region is typically expected to be only mildly relativistic
on time scales of interest (the corresponding observation time after the jet
launch is $T\sim 10^{5}-10^{6}$ s).

Numerical simulations have demonstrated that binary SMBH mergers can produce
jet-like emissions driven by the Poynting outflow (e.g., Kelly et al., 2017).
We assume that a jet is launched after the coalescence and subsequently
propagates in the wind bubble formed by pre-merger disk winds. Fig. 1
schematically illustrates the configuration of the system. The disk wind
expands in the gaseous environment of the host galaxy. We focus on emissions
from the shock between the jet and the wind bubble. Initially, the
circumbinary disk can react promptly to the evolution of the binary system.
The ratio between the disk radius $R_{d}$ and the semi-major axis of the
binary system $a$ remains unchanged ($R_{d}/a\sim 2$), until the inspiral time
scale $t_{\rm GW}$ of the binary system (e.g., Shapiro & Teukolsky, 1983)
equals the viscosity time scale $t_{\rm vis}$ (e.g., Pringle, 1981), which is
known as the disk decoupling. After the disk becomes decoupled, the merger of
SMBHs in binary system occurs within the time interval $t_{m}\sim 3\times
10^{-2}~{}{\rm yr}~{}M_{\rm BH,6}\alpha_{-1}^{-8/5}h_{-1}^{-16/5}$, where
$M_{\rm BH}=10^{6}~{}M_{\rm BH,6}M_{\odot}$ is the mass of the binary system,
the dimensionless parameter $h$ is defined as $h=H/R_{d}$, $\alpha\sim 0.1$ is
the viscosity parameter, and $H$ is the disk scale height. The disk gas starts
to fill the cavity between the disk and the SMBHs in the viscosity timescale,

$t_{\rm vis}\sim 0.1~{}{\rm yr}~{}M_{\rm
BH,6}\alpha_{-1}^{-8/5}h_{-1}^{-16/5}$ (1)

after the coalescence (Farris et al., 2015). This leads to a time delay
($t_{\rm delay}\sim t_{\rm vis}$) of days to months between the GW burst and
the launch of post-merger jets, if $h\sim 0.1-0.3$ is assumed. However, for a
thick and highly magnetized disk with $h\sim\alpha\sim 1$, $t_{\rm delay}$
could be much shorter.

Figure 1: Schematic description of our model. Left panel: pre-merger disk
winds launched from the circumbinary disk. The green arrows illustrate the
disk-driven outflows that form a wind bubble. Mini-disks around each SMBH are
also shown. Right panel: post-merger jets launched by a merged SMBH. The
forward shock region is shown as the purple area. The cocoon is not depicted.

On the other hand, within the duration of these two short-term processes,
e.g., $t_{\rm m}$ and $t_{\rm vis}$, the disk wind radius may reach
$v_{d}(t_{\rm vis}-t_{m})\sim 10^{14}-10^{\rm 16}~{}\rm cm$ above the disk,
where $v_{d}$ is the disk wind velocity that is of the order of the escape
velocity, $v_{\rm esc}(R_{d,\rm dec})\approx\sqrt{2GM_{\rm BH}/R_{d,\rm dec}}$
for the circumbinary disk, and $R_{d,\rm dec}\approx 1.2\times 10^{13}~{}{\rm
cm}~{}M_{{\rm BH},6}\alpha_{-1}^{-2/5}h_{-1}^{-4/5}$ is the radius of the
circumbinary disk at the decoupling. In reality, not only the circumbinary
disk but also mini-disks around two SMBHs contribute, which would make the
wind bubble more complicated. For simplicity, we assume the density profile of
the winds at the decoupling to obtain the density distribution of the wind
bubble at larger distances,

$\varrho_{w}(r)=\frac{\eta_{w}(1+\chi)\dot{M}_{\rm BH}}{4\pi r^{2}v_{d}}\equiv
Dr^{-2},$ (2)

where $\dot{M}_{\rm BH}$ is the mass accretion rate onto the binary system,
$\chi\sim 1$ is introduced to take into account the contribution of mini-
disks, and $\eta_{w}$ represents the fraction of accreted mass converted to
the disk wind. According to the simulations, for SANE (Standard And Normal
Evolution) models, the parameter $\eta_{w}$ may vary from $10^{-4}$ to
$10^{-1}$ (Jiang et al., 2019b, a; Ohsuga et al., 2009) when the mass
accretion rate changes from sub-Eddington to super-Eddington. In MAD
(Magnetically Arrested Disk) models, $\eta_{w}$ can reach $10^{-2}$ to
$10^{-1}$ (Akiyama et al., 2019). With $v_{d}\sim v_{\rm esc}(R_{d,\rm dec})$,
we have $D\simeq 5.9\times 10^{11}~{}{\rm
g~{}cm^{-1}}~{}\tilde{\eta}_{w,-1.5}(\dot{m}/0.5)M_{\rm
BH,6}\beta_{d,-1}^{-1}$, where $\tilde{\eta}_{w}\equiv(1+\chi)\eta_{w}$,
$\beta_{d,-1}\equiv v_{d}/(0.1c)$, the parameter $\dot{m}$ is defined as the
ratio of $\dot{M}_{\rm BH}$, and the Eddington value $\dot{M}_{\rm Edd}\equiv
10L_{\rm Edd}/c^{2}$ (assuming a radiation efficiency of 0.1).

After the coalescence, a powerful jet driven by the spin energy of the newly
formed SMBH can appear, subsequently propagating in the pre-merger wind
bubble. Considering a sub-Eddington accretion rate with the MAD configuration,
we estimate the jet kinetic luminosity to be

$\displaystyle L_{k,j}$ $\displaystyle=$ $\displaystyle\eta_{j}\dot{M}_{\rm
BH}c^{2}$ (3) $\displaystyle\simeq$ $\displaystyle 6.3\times 10^{44}{~{}\rm
erg\ s^{-1}}\eta_{j}(\dot{m}/0.5)M_{\rm BH,6},$

where $\eta_{j}\sim 0.3-1$ is the ratio of the accretion energy converted to
the jet energy (Tchekhovskoy et al., 2011).

Following the standard jet propagation theory (e.g., Bromberg et al., 2011;
Mizuta & Ioka, 2013), we write down the dimensionless parameter that
represents the ratio of the energy density of the jet and the rest-mass energy
density of the surrounding medium

$\tilde{L}\approx\frac{L_{k,\rm iso}}{{4\pi r^{2}\varrho_{w}c^{3}}}\simeq
63~{}\tilde{\eta}_{w,-1.5}^{-1}\eta_{j}\theta_{j,-0.5}^{-2}\beta_{d,-1},$ (4)

where $\theta_{j}$ is the jet opening angle, $L_{k,\rm iso}\approx
2L_{k,j}/\theta_{j}^{2}$ is the isotropic-equivalent luminosity. Since the
quantity $\tilde{L}$ lies in the regime $\theta_{j}^{-4/3}\ll\tilde{L}$, we
expect that the jet is “uncollimated” for our fiducial parameters111However,
jet collimation, which was assumed in Yuan et al. (2020), would be achievable
for the super-Eddington accretion accompanied by disk winds with $\eta_{w}\sim
0.1-0.3$.. This situation is similar to that in choked jet propagation in the
circumstellar material (Senno et al., 2016; Nakar, 2015), and $\beta_{h}$ is
evaluated from $\tilde{L}$. In the relativistic limit, the jet head Lorentz
factor is $\Gamma_{h}\approx\tilde{L}^{1/4}/\sqrt{2}$ (Senno et al., 2016),
and we have $\Gamma_{h}\sim 2$ in our fiducial case with
$\tilde{\eta}_{w}=10^{-1.5}$. Note that the jet head radius is
$R_{h}=c\beta_{h}\hat{T}\approx c\hat{T}$, and $\hat{T}$ is introduced to
represent time measured in the central engine frame, which can be converted to
the observation time $T$ via $T=(1+z)(1-\beta_{h})\hat{T}$ (that is
$T\approx(1+z)\hat{T}/[2\Gamma_{h}^{2}]$ in the relativistic limit) for on-
axis observers.

Furthermore, to ensure particle acceleration, we impose radiation constraints
requiring that the shock is collisionless, without being mediated by radiation
(Murase & Ioka, 2013; Senno et al., 2016). Here, ignoring effects of pair
production, we use the conservative condition,
$\tau_{T}\approx\varrho_{w}\sigma_{T}R_{h}/m_{p}<1$, where $\sigma_{T}$ is the
Thomson cross section. Numerically, this condition is satisfied at
$\hat{T}\gtrsim 10~{}\rm s$, which is much shorter than the duration of EM
emission.

Figure 2: Left panel: Non-thermal energy spectra expected for uncollimated
post-merger jets from a SMBH merger located at $z=1$. The solid and dashed
lines represent the synchrotron and SSC components. The dash-dotted lines show
the sensitivity curves for current and future detectors. Right panel: Multi-
wavelength light curves. The yellow and blue dashed vertical lines illustrate
respectively the characteristic times, e.g., $T_{\rm ssa}$, of 100 GHz and 5
GHz emissions. The used parameters are $\dot{m}=0.5$, $M_{\rm
BH}=10^{6}~{}M_{\odot}$, $\tilde{\eta}_{w}={10}^{-1.5}$, $\eta_{j}=1$,
$\theta_{j}=10^{-0.5}$, $s=2.0$, $\zeta_{e}=0.4$, $\epsilon_{e}=0.1$ and
$\epsilon_{B}=0.01$.

##  3 Electromagnetic emission from post-merger jets

With the jet dynamics presented in the previous section, we calculate the EM
spectra resulting from synchrotron and SSC emission. As in the standard theory
of GRB afterglows (e.g., Mészáros, 2006), we assume that electrons are
accelerated at the external forward shock with a power-law spectral index $s$.
The energy fractions of the downstream energy density converted to non-thermal
electron and magnetic field energy are defined as $\epsilon_{e}$ and
$\epsilon_{B}$, respectively. The upstream number density is given by
$n_{h,u}\approx\varrho_{w}(R_{h})/m_{p}\propto R_{h}^{-2}$, and
$B\approx[\epsilon_{B}32\pi\Gamma_{h}(\Gamma_{h}-1)n_{h,u}m_{p}c^{2}]^{1/2}$
is the downstream magnetic field strength.

In the relativistic limit ($\Gamma_{h}\gg 1$), the characteristic injection
frequency $\nu_{m}$ and the cooling frequency $\nu_{c}$ in the observer frame
are written respectively as,

$\begin{split}\nu_{m}&\approx\frac{3\Gamma_{h}\gamma_{m}^{2}eB}{4\pi(1+z)m_{e}c}\\\
&\simeq 3.4\times 10^{3}{\ \rm GHz}\
\epsilon_{e,-1}^{2}{\zeta_{e,-0.4}^{2}}\epsilon_{B,-2}^{1/2}\eta_{j}^{1/2}T_{4}^{-1}\\\
&~{}~{}~{}\times(\dot{m}/0.5)^{1/2}M_{\rm
BH,6}^{1/2}\theta_{j,-0.5}^{-1}\end{split}$ (5)

and

$\begin{split}\nu_{c}&\approx\frac{3\Gamma_{h}\gamma_{c}^{2}eB}{4\pi(1+z)m_{e}c}\\\
&\simeq 4.6\times 10^{2}{\ \rm GHz\
}(1+z)^{-2}(1+Y)^{-2}\epsilon_{B,-2}^{-3/2}\\\
&~{}~{}~{}\times\tilde{\eta}_{w,-1.5}^{-2}\theta_{j,-0.5}^{-1}(\dot{m}/0.5)^{-3/2}T_{4}\beta_{d,-1}^{2}M_{\rm
BH,6}^{-3/2},\end{split}$ (6)

where $\gamma_{m}=\epsilon_{e}\zeta_{e}(\Gamma_{h}-1)m_{p}/m_{e}$ is the
electron minimum Lorentz factor, and $\gamma_{c}=6\pi
m_{e}c/[(1+Y)T^{\prime}\sigma_{T}B^{2}]$ is the cooling Lorentz factor. Here,
$\zeta_{e}=g_{s}/f_{e}=1/[f_{e}\ln(\gamma_{M}/\gamma_{m})]\sim 0.3-0.4$ is
constrained by the particle-in-cell simulations (Park et al., 2015) (where
$f_{e}$ is the fraction of accelerated electrons and the maximum Lorentz
factor of electrons is $\gamma_{M}={(6\pi
e)^{1/2}/[\sigma_{T}B(1+Y)]^{1/2}}$), $Y$ is the Compton parameter, and
$T^{\prime}=\hat{T}/\Gamma_{h}\approx 2\Gamma_{h}T/(1+z)$ is the comoving
time. For example, at $T=10^{4}~{}\rm s$, we have $Y\simeq 2.4$, corresponding
to the fast cooling regime. It changes to the slow cooling regime on a time
scale from days to weeks. We obtain the peak synchrotron flux (e.g., Wijers &
Galama, 1999)

$\begin{split}F_{\nu,\rm syn}^{\rm
max}&\approx\frac{(1+z){(0.6f_{e}n_{h,u}R_{h}^{3})}\Gamma_{h}e^{3}B}{\sqrt{3}m_{e}c^{2}d_{L}^{2}}\\\
&\simeq 0.24~{}{\rm
mJy}~{}(1+z)g_{s,-1.2}\zeta_{e,-0.4}^{-1}(\dot{m}/0.5)^{3/2}\\\
&~{}~{}~{}\times\eta_{j}^{1/2}\epsilon_{B,-2}^{1/2}\tilde{\eta}_{w,-1.5}\beta_{d,-1}^{-1}\theta_{j,-0.5}^{-1}M_{\rm
BH,6}^{3/2}d_{L,28}^{-2}.\end{split}$ (7)

The low-frequency synchrotron emission is subject to synchrotron self-
absorption (SSA). The SSA optical depth is written as $\tau_{\rm
ssa}(\nu)={\xi_{s}}{en_{h,u}R_{h}}(\nu/\nu_{n})^{-p}/[{B\gamma_{n}^{5}}]$,
where $\nu$ is the observed frequency, $\xi_{s}\sim 5-10$ depends on the
electron spectral index, $\gamma_{n}=\min\left[\gamma_{m},\
\gamma_{c}\right]$, $\nu_{n}=\gamma_{n}^{2}eB/[(1+z)m_{e}c]$, $p=5/3$ for
$\nu<\nu_{n}$ and $p=(4+s)/2$ or $p=3$ for $\nu>\nu_{n}$ depending on the slow
or fast cooling regime (e.g., Panaitescu & Kumar 2000; Murase et al. 2014).
The critical time scales set by $\tau_{\rm ssa}=1$ for $\nu<\nu_{n}$ and for
$\nu>\nu_{n}$ are $T_{\rm ssa}\simeq 5.4\times 10^{5}~{}{\rm
s}~{}\xi_{s,1}^{3/10}(1+z)^{1/2}(1+Y)^{1/2}\epsilon_{B,-2}^{3/5}\left(\frac{\nu}{1\
\rm
GHz}\right)^{-1/2}\tilde{\eta}_{w,-1.5}^{11/10}\beta_{d,-1}^{-11/10}(\dot{m}/0.5)^{9/10}M_{\rm
BH,6}^{9/10}$ and $T_{\rm ssa}\simeq 3.5\times 10^{5}~{}{\rm
s}~{}\xi_{s,1}^{1/2}(1+z)^{-1/2}(1+Y)^{-1/2}M_{\rm
BH,6}^{1/2}\left(\frac{\nu}{100\ \rm
GHz}\right)^{-3/2}(\dot{m}/0.5)^{1/2}M_{\rm
BH,6}^{1/2}\tilde{\eta}_{w,-1.5}^{1/2}\beta_{d,-1}^{-1/2}$, respectively.
Thus, we expect that EM emission at 5 GHz and 100 GHz reaches a peak about a
few days after the jet launch ($T_{\rm ssa}\simeq 7.1\times 10^{5}~{}\rm s$
and $T_{\rm ssa}\simeq 3.1\times 10^{5}~{}\rm s$, respectively, in our
fiducial case with $\xi_{s}=8.7$).

We numerically calculate the electron distribution and the resulting
synchrotron and SSC spectra of the forward shock, following the method used in
Murase et al. (2011) and Zhang et al. (2020). We solve the continuity equation
that describes the evolution of the electron spectra and calculate the
synchrotron/SSC components, in which the trans-relativistic regime can be
consistently treated as in Zhang et al. (2020). Combining the obtained radio,
millimeter, optical and X-ray light curves with the sensitivities of
corresponding detectors, we discuss the possibility of follow-up observations
of the EM counterpart.

The left panel of Fig. 2 shows the snapshots of synchrotron and SSC spectra at
$T=10^{4}-10^{6}~{}\rm s$ for an on-axis source located at $z=1$. We assume
$s=2.0$, $\epsilon_{e}=0.1$, and $\epsilon_{B}=0.01$. The solid and dashed
lines correspond to the synchrotron and SSC components. Very high-energy
gamma-ray emission at $\gtrsim 1$ TeV energies is suppressed due to the Klein-
Nishina effect (e.g., Murase et al., 2011; Zhang et al., 2020), and the
$\gamma\gamma$ annihilation with the extragalactic background light (EBL). For
the EBL correction, $\gamma\gamma$ optical depth data from Model C in Finke et
al. (2010) is used. To show how the EM signal evolves with time, we illustrate
the gamma-ray (1 GeV), X-ray (1 keV), UV (1 eV) and radio (5 GHz and 100 GHz)
light curves in the right panel. In particular, before the characteristic time
$T_{\rm ssa}$ (shown as the vertical yellow and blue lines). The radio
emission is suppressed by the SSA process, which is implemented by multiplying
$(1-e^{-\tau_{\rm ssa}})/\tau_{\rm ssa}$.

Figure 3: Detection horizons for multi-wavelength detectors, e.g., SKA, VLA,
EVLA, ALMA, HST, JWST, LSST and Chandra. The horizontal dotted line shows the
100 GHz detection window for ALMA assuming a source located at $z=3$. Similar
to Fig. 2, the dotted vertical lines are the characteristic times of 5 GHz and
100 GHz signals.

It is useful to discuss the detection horizon $d_{\rm lim}$ for some detectors
such as the Square Kilometre Array (SKA), Very Large Array (VLA), Expanded VLA
(EVLA), Atacama Large Millimeter Array (ALMA), Hubble Space Telescope (HST),
James Webb Space Telescope (JWST), Large Synoptic Survey Telescope (LSST) and
the high-resolution camera on the Chandra X-ray Observatory (Chandra) 222For
information on these facilities see, e.g., : VLA (http://www.vla.nrao.edu),
EVLA (http://www.aoc.nrao.edu/evla/), SKA (https://www.skatelescope.org), ALMA
(https://public.nrao.edu/telescopes/alma/), HST
(https://www.nasa.gov/mission_pages/hubble/main/index.html), JWST
(https://stsci.edu/jwst), LSST (https://www.lsst.org/scientists/scibook) and
Chandra (https://cxc.cfa.harvard.edu/cdo/about_chandra/) as functions of the
observation time $T$. Given the observed flux $F_{\nu}(\nu_{\gamma},T,z)$ at
the observer time $T$ from an on-axis source located at redshift $z$, the
horizon can be calculated iteratively via

$d_{\rm lim}(\nu_{\gamma},T)=d_{L}\left(\frac{\frac{1}{\Delta T_{\rm
exp}}\int_{T}^{T+\Delta T_{\rm exp}}F_{\nu}(\nu_{\gamma},t,z)dt}{F_{\rm
lim}(\nu_{\gamma},\Delta T_{\rm exp})}\right)^{1/2},$ (8)

where $F_{\rm lim}(\nu_{\gamma},\Delta T_{\rm exp})$ is the detector
sensitivity normalized to the exposure time $\Delta T_{\rm exp}$. For example,
specifying the detection frequency $\nu=100~{}\rm GHz$, the sensitivity of
ALMA is approximately $34~{}\mu\rm Jy$ for one-hour integration, e.g., $\Delta
T_{\rm exp}=1$ hour. Fig. 3 indicates the detection horizons for SKA (5 GHz,
$\Delta T_{\rm exp}=10$ hr), SKA (1 GHz333At 1 GHz, the SKA field-of-view can
reach $\gtrsim$1 deg2, $\Delta T_{\rm exp}=10$ hr), VLA (5 GHz, $\Delta T_{\rm
exp}=1$ hr), ALMA (100 GHz, $\Delta T_{\rm exp}=1$ hr), JWST (1 eV, $\Delta
T_{\rm exp}=10$ ks), HST (1 eV, $\Delta T_{\rm exp}=10$ ks), LSST (r-band,
point source exposure time $\Delta T_{\rm exp}=30$ s in the 3-day revisit
time), and Chandra (1 keV, $\Delta T_{\rm exp}=100$ ks). The vertical black
and blue dotted lines respectively illustrate the times $T_{\rm ssa}$ at which
photons at 100 GHz and 5 GHz bands start to survive from the synchrotron self-
absorption.

From Fig. 3, we expect that ALMA, SKA and EVLA can detect SMBH mergers in the
radio bands respectively out to redshifts of $z\sim 4-6$. Remarkably, the
optical and X-ray signals from the mergers in the range $1\lesssim z\lesssim
2$ can also be identified through targeted searches by Chandra, HST and JWST
in a long duration after the merger. In addition, we can estimate the
observation time for each detector if the luminosity distance of the merger is
specified. For example, a source located at $z=3$ would remain detectable by
ALMA for roughly 20-30 days (see the black dotted horizontal line). One caveat
is that this calculation is carried out in the ideal case where the detectors
can point to the position of the source and start the observation immediately
after the EM signal reaches the Earth. We discuss the sky coverage and a
detection strategy in the following Sec. 4.

##  4 Summary and Discussion

We investigated broadband non-thermal EM emission from electrons accelerated
at the external forward shock expected in post-merger jets from the
coalescence of SMBHs. In our model, the jets can be launched at $t_{\rm
delay}\sim t_{\rm vis}\sim(0.003-0.1)M_{\rm BH,6}$ yr after the coalescence.
The time lag is primarily determined by the scale height of the circumbinary
disk and the viscosity parameter. We found that, for a moderate accretion rate
($\dot{m}\sim 0.5$), the multi-wavelength emission from such a system may
persist at detectable levels for months after the jet launch, depending on the
facilities and the luminosity distance. Moreover, according to our model, the
sources with moderate $\dot{m}=0.5$ can be detected up to $z\sim 5-6$,
covering the range that LISA-like GW detectors have the best detection chance,
e.g., $z\sim 1-2$, in which $(1-10)f_{b}$ mergers per year are expected (Menou
et al., 2001; Enoki et al., 2004; Arun et al., 2009; Amaro-Seoane et al.,
2012; Dal Canton et al., 2019; Kara et al., 2019). Here $f_{b}\sim
1/(2\Gamma_{h}^{2})$ is the beaming factor in our model. Because the jet head
Lorentz factor is as low as $\Gamma_{h}\lesssim 2$, the EM emission from the
forward shock region is not highly beamed and we expect $f_{b}\sim 0.1-1$.
This makes the binary SMBH mergers interesting targets for future multi-
messenger studies. If super-Eddington accretion (e.g., $\dot{m}\sim 10$)
occurs, as was optimistically assumed in Yuan et al. (2020), even LSST and
Chandra could detect EM signals from the sources in the redshift range
$4\lesssim z\lesssim 6$. We showed the case of $s=2.0$ for the demonstration.
If a larger spectral index, e.g., $s\sim 2.2-2.4$, is used, as expected from
observations of GRB afterglows, the radio detection would be more promising
whereas a higher accretion rate would be required for successful optical and
X-ray observations.

The density of the premerger bubble, which was assumed to be a wind profile,
is subject to large uncertainties. The extrapolation in the density
distribution would be applicable up to an outer wind radius of $\sim
10^{14}-10^{16}$ cm. The density predicted by equation 2 would drop below that
of the central molecular zone (indicated as the circumnuclear environment in
Fig.1), which may lead to the increase of radio emission. In addition, a
cocoon formed along with the jet, depending on uncertain details of the
medium, could produce thermal photons which may not only lead to detectable
signals but also serve as seed photons for inverse-Compton emission. We
focused on the more secure EM emission from the forward shock region as the
jet propagates in the wind. In this sense, our prediction for the fluxes are
conservative.

EM emission from the external reverse shock and internal shocks can also be
expected (e.g., Mészáros & Rees, 1999; Kobayashi & Zhang, 2003, for the
reverse shock emission in GRBs). Qualitatively, the ratio between the peak
fluxes of the reverse and forward shock emission depends on the value of
$\Gamma_{j}$, and the reverse shock contribution might be important for
$\Gamma_{j}\gg\Gamma_{h}$.

Previous studies based on general relativistic three dimensional
magnetohydrodynamics simulations have shown that the circumbinary disk and the
corona can emit light in UV/EUV bands (e.g., d’Ascoli et al., 2018), while
X-ray and infrared emission from the post-merger circumbinary disk are
expected to last for years (Milosavljević & Phinney, 2005; Schnittman &
Krolik, 2008). In the pre-merger phase, the orbits of dual SMBH cores may be
identified by radio facilities such as VLBI (e.g., Rodriguez et al., 2006).
Blind searches could identify radio or UV/EUV sources from the binary SMBH
systems, which would provide complementary constraints on the source location,
the accretion rate and the ambient gaseous environment.

Our model can provide a guidance, including the onset times and the detection
windows, in developing detection strategies for future EM follow-up
observations, once GW signals are detected. Considering the large
uncertainties in the localization with GW detectors, an initial follow-up
using large field-of-view (FOV) telescopes, like SKA and LSST, would be
necessary to more precisely localize the position of the source. After that,
we can use the putative positional information from the initial follow-up
imaging to guide the observation of narrower FOV telescopes. In particular,
for high-redshift mergers in the range $z\sim 2-5$, EM follow-up observations
rely more on radio detectors, and the detection is possible a few weeks after
the merger. SKA needs the source localization before follow-up observations by
VLA and ALMA. On the other hand, if the merger is close enough (e.g., $z\sim
1$), LISA observations staring from a few weeks before the merger can localize
the merger with a median precision of $\sim$1 deg2 (Mangiagli et al., 2020).
In this case, LISA and LSST can jointly guide other X-ray and optical
facilities in the very early stage. Amid these two regimes, e.g., $z\sim 1-2$,
detections in the optical and X-ray bands using HST, JWST and Chandra would be
promising if the source is localized by SKA.

C.C.Y. and P.M. acknowledge support from the Eberly Foundation. The work of
K.M. is supported by NSF Grant No. AST-1908689, and KAKENHI No. 20H01901 and
No. 20H05852. B.T.Z. acknowledges the IGC fellowship. S.S.K. acknowledges the
JSPS Research Fellowship, JSPS KAKENHI Grant No. 19J00198.

## References

  * Akiyama et al. (2019) Akiyama, K., Alberdi, A., Alef, W., et al. 2019, ApJ, 875, L5, doi: 10.3847/2041-8213/ab0f43
  * Amaro-Seoane et al. (2012) Amaro-Seoane, P., Aoudia, S., Babak, S., et al. 2012, Classical and Quantum Gravity, 29, 124016, doi: 10.1088/0264-9381/29/12/124016
  * Amaro-Seoane et al. (2017) Amaro-Seoane, P., Audley, H., Babak, S., et al. 2017, arXiv:1702.00786.  https://arxiv.org/abs/1702.00786
  * Arun et al. (2009) Arun, K., Babak, S., Berti, E., et al. 2009, Classical and Quantum Gravity, 26, 094027, doi: 10.1088/0264-9381/26/9/094027
  * Arzoumanian et al. (2020) Arzoumanian, Z., Baker, P. T., Blumer, H., et al. 2020, ApJ, 905, L34, doi: 10.3847/2041-8213/abd401
  * Baker et al. (2019) Baker, J., Barke, S. F., Bender, P. L., et al. 2019, 51, 243.  https://arxiv.org/abs/1907.11305
  * Begelman et al. (1980) Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980, Nature, 287, 307, doi: 10.1038/287307a0
  * Blandford & Znajek (1977) Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433, doi: 10.1093/mnras/179.3.433
  * Bromberg et al. (2011) Bromberg, O., Nakar, E., Piran, T., et al. 2011, ApJ, 740, 100, doi: 10.1088/0004-637X/740/2/100
  * Dal Canton et al. (2019) Dal Canton, T., Mangiagli, A., Noble, S. C., et al. 2019, ApJ, 886, 146, doi: 10.3847/1538-4357/ab505a
  * d’Ascoli et al. (2018) d’Ascoli, S., Noble, S. C., Bowen, D. B., et al. 2018, ApJ, 865, 140, doi: 10.3847/1538-4357/aad8b4
  * Enoki et al. (2004) Enoki, M., Inoue, K. T., Nagashima, M., & Sugiyama, N. 2004, ApJ, 615, 19, doi: 10.1086/424475
  * Farris et al. (2015) Farris, B. D., Duffell, P., MacFadyen, A. I., & Haiman, Z. 2015, MNRAS, 447, L80, doi: 10.1093/mnrasl/slu184
  * Finke et al. (2010) Finke, J. D., Razzaque, S., & Dermer, C. D. 2010, ApJ, 712, 238, doi: 10.1088/0004-637X/712/1/238
  * Haiman (2017) Haiman, Z. 2017, Phys. Rev. D, 96, 023004, doi: 10.1103/PhysRevD.96.023004
  * Haiman et al. (2009) Haiman, Z., Kocsis, B., Menou, K., Lippai, Z., & Frei, Z. 2009, Classical and Quantum Gravity, 26, 094032, doi: 10.1088/0264-9381/26/9/094032
  * Jiang et al. (2019a) Jiang, Y.-F., Blaes, O., Stone, J. M., & Davis, S. W. 2019a, ApJ, 885, 144, doi: 10.3847/1538-4357/ab4a00
  * Jiang et al. (2019b) Jiang, Y.-F., Stone, J. M., & Davis, S. W. 2019b, ApJ, 880, 67, doi: 10.3847/1538-4357/ab29ff
  * Kara et al. (2019) Kara, E., Margutti, R., Keivani, A., et al. 2019, BAAS, 51, 112.  https://arxiv.org/abs/1903.05287
  * Kelly et al. (2017) Kelly, B. J., Baker, J. G., Etienne, Z. B., Giacomazzo, B., & Schnittman, J. 2017, Phys. Rev. D, 96, 123003, doi: 10.1103/PhysRevD.96.123003
  * Kobayashi & Zhang (2003) Kobayashi, S., & Zhang, B. 2003, ApJ, 582, L75, doi: 10.1086/367691
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Kroupa et al. (2020) Kroupa, P., Subr, L., Jerabkova, T., & Wang, L. 2020, MNRAS, 498, 5652, doi: 10.1093/mnras/staa2276
  * Mangiagli et al. (2020) Mangiagli, A., Klein, A., Bonetti, M., et al. 2020, Phys. Rev. D, 102, 084056, doi: 10.1103/PhysRevD.102.084056
  * Menou et al. (2001) Menou, K., Haiman, Z., & Narayanan, V. K. 2001, ApJ, 558, 535, doi: 10.1086/322310
  * Mészáros (2006) Mészáros, P. 2006, Rept. Prog. Phys., 69, 2259, doi: 10.1088/0034-4885/69/8/R01
  * Mészáros et al. (2019) Mészáros, P., Fox, D. B., Hanna, C., & Murase, K. 2019, Nature Reviews Physics, 1, 585, doi: 10.1038/s42254-019-0101-z
  * Mészáros & Rees (1999) Mészáros, P., & Rees, M. J. 1999, MNRAS, 306, L39, doi: 10.1046/j.1365-8711.1999.02800.x
  * Milosavljević & Phinney (2005) Milosavljević, M., & Phinney, E. S. 2005, ApJ, 622, L93, doi: 10.1086/429618
  * Mingarelli et al. (2017) Mingarelli, C. M. F., Lazio, T. J. W., Sesana, A., et al. 2017, Nature Astronomy, 1, 886, doi: 10.1038/s41550-017-0299-6
  * Mizuta & Ioka (2013) Mizuta, A., & Ioka, K. 2013, ApJ, 777, 162, doi: 10.1088/0004-637X/777/2/162
  * Moesta et al. (2012) Moesta, P., Alic, D., Rezzolla, L., Zanotti, O., & Palenzuela, C. 2012, ApJ, 749, L32, doi: 10.1088/2041-8205/749/2/L32
  * Murase & Bartos (2019) Murase, K., & Bartos, I. 2019, Ann. Rev. Nucl. Part. Sci., 69, 477, doi: 10.1146/annurev-nucl-101918-023510
  * Murase & Ioka (2013) Murase, K., & Ioka, K. 2013, Phys. Rev. Lett., 111, 121102, doi: 10.1103/PhysRevLett.111.121102
  * Murase et al. (2014) Murase, K., Thompson, T. A., & Ofek, E. O. 2014, MNRAS, 440, 2528, doi: 10.1093/mnras/stu384
  * Murase et al. (2011) Murase, K., Toma, K., Yamazaki, R., & Mészáros, P. 2011, ApJ, 732, 77, doi: 10.1088/0004-637X/732/2/77
  * Nakar (2015) Nakar, E. 2015, Astrophys. J., 807, 172, doi: 10.1088/0004-637X/807/2/172
  * Ohsuga et al. (2009) Ohsuga, K., Mineshige, S., Mori, M., & Kato, Y. 2009, PASJ, 61, L7, doi: 10.1093/pasj/61.3.L7
  * Panaitescu & Kumar (2000) Panaitescu, A., & Kumar, P. 2000, ApJ, 543, 66, doi: 10.1086/317090
  * Park et al. (2015) Park, J., Caprioli, D., & Spitkovsky, A. 2015, Phys. Rev. Lett., 114, 085003, doi: 10.1103/PhysRevLett.114.085003
  * Pringle (1981) Pringle, J. 1981, ARA&A, 19, 137, doi: 10.1146/annurev.aa.19.090181.001033
  * Ravi (2018) Ravi, V. 2018, in Astronomical Society of the Pacific Conference Series, Vol. 517, Science with a Next Generation Very Large Array, ed. E. Murphy, 781
  * Rodriguez et al. (2006) Rodriguez, C., Taylor, G. B., Zavala, R. T., et al. 2006, ApJ, 646, 49, doi: 10.1086/504825
  * Schnittman (2011) Schnittman, J. D. 2011, Classical and Quantum Gravity, 28, 094021, doi: 10.1088/0264-9381/28/9/094021
  * Schnittman & Krolik (2008) Schnittman, J. D., & Krolik, J. H. 2008, ApJ, 684, 835, doi: 10.1086/590363
  * Senno et al. (2016) Senno, N., Murase, K., & Mészáros, P. 2016, Phys. Rev. D, 93, 083003, doi: 10.1103/PhysRevD.93.083003
  * Sesana et al. (2004) Sesana, A., Haardt, F., Madau, P., & Volonteri, M. 2004, ApJ, 611, 623, doi: 10.1086/422185
  * Shapiro & Teukolsky (1983) Shapiro, S. L., & Teukolsky, S. A. 1983, Black holes, white dwarfs, and neutron stars: The physics of compact objects
  * Tanaka & Menou (2010) Tanaka, T., & Menou, K. 2010, ApJ, 714, 404, doi: 10.1088/0004-637X/714/1/404
  * Taylor et al. (2019) Taylor, S. R., Burke-Spolaor, S., Baker, P. T., et al. 2019, BAAS, 51, 336.  https://arxiv.org/abs/1903.08183
  * Tchekhovskoy et al. (2011) Tchekhovskoy, A., Narayan, R., & McKinney, J. C. 2011, MNRAS, 418, L79, doi: 10.1111/j.1745-3933.2011.01147.x
  * Thorne & Braginskii (1976) Thorne, K. S., & Braginskii, V. 1976, ApJ, 204, L1, doi: 10.1086/182042
  * Wijers & Galama (1999) Wijers, R., & Galama, T. 1999, ApJ, 523, 177, doi: 10.1086/307705
  * Yuan et al. (2020) Yuan, C., Murase, K., Kimura, S. S., & Mészáros, P. 2020, Phys. Rev. D, 102, 083013, doi: 10.1103/PhysRevD.102.083013
  * Zhang et al. (2020) Zhang, B. T., Murase, K., Veres, P., & Mészáros, P. 2020, arXiv: 2012.07796.  https://arxiv.org/abs/2012.07796