# Late Jets, Early Sparks: Illuminating the Pre-Maximum Bumps in Superluminous
Supernovae

Ore Gottlieb<EMAIL_ADDRESS>Center for Computational
Astrophysics, Flatiron Institute, 162 5th Avenue, New York, NY 10010, USA
Department of Physics and Columbia Astrophysics Laboratory, Columbia
University, Pupin Hall, New York, NY 10027, USA  Brian D. Metzger Department
of Physics and Columbia Astrophysics Laboratory, Columbia University, Pupin
Hall, New York, NY 10027, USA  Center for Computational Astrophysics, Flatiron
Institute, 162 5th Avenue, New York, NY 10010, USA

###### Abstract

Superluminous supernovae (SLSNe) radiate $\gtrsim 10-100$ times more energy
than ordinary stellar explosions, implicating a novel power source behind
these enigmatic events. One frequently discussed source, particularly for
hydrogen-poor (Type I) SLSNe, is a central engine such as a millisecond
magnetar or accreting black hole. Both black hole and magnetar engines are
expected to channel a fraction of their luminosity into a collimated
relativistic jet. Using 3D relativistic hydrodynamical simulations, we explore
the interaction of a relativistic jet, endowed with a luminosity $L_{\rm
j}\approx 10^{45}\,{\rm erg\,s^{-1}}$ and duration $t_{\rm eng}\approx$ days
compatible with those needed to power SLSNe, launched into the envelope of the
exploding star. The jet successfully breaks through the expanding ejecta and
its shocked cocoon powers ultraviolet/optical emission lasting several days
after the explosion and reaching a peak luminosity $\gtrsim 10^{44}$ erg s-1,
corresponding to a sizable fraction of $L_{\rm j}$. This high radiative
efficiency is the result of the modest adiabatic losses the cocoon experiences
due to the low optical depths of the enlarged ejecta at these late times,
e.g., compared to the more compact stars in gamma-ray bursts. The luminosity
and temperature of the cocoon emission match those of the “bumps” in SLSNe
light curves observed weeks prior to the optical maximum in many SLSNe.
Confirmation of jet breakout signatures by future observations (e.g., days to
weeks-long internal X-ray emission from the jet for on-axis observers,
spectroscopy confirming large photosphere velocities $v/c\gtrsim 0.1$, or
detection of a radio afterglow) would offer strong evidence for central
engines powering SLSNe.

##  1 Introduction

The collapse of a rotating massive star giving birth to a magnetized spinning
compact object $-$ either a Kerr black hole or a neutron star $-$ offers a
compelling scenario for generating some of the most energetic transients in
the Universe, including gamma-ray bursts (GRB) and superluminous supernovae
(SLSNe). While most GRBs with long durations are accompanied by stripped-
envelope supernovae with abnormally large kinetic energies (“broad-lined” Type
Ibc SNe; e.g., Stanek et al. 2003; Woosley & Bloom 2006; Liu & Modjaz 2016),
SLSNe are instead distinguished by their large radiated energies, which can
exceed those of ordinary stellar explosions by a factor of $\gtrsim 10-100$
(e.g., Gal-Yam 2018; Nicholl 2021 for reviews). Type I (i.e., hydrogen-poor)
SLSNe, in particular, may be associated with the birth of a millisecond (ms)
neutron star with a strong magnetic field (“ms magnetar”; e.g., Kasen &
Bildsten 2010; Woosley 2010; Nicholl et al. 2017; Vurm & Metzger 2021). The
magnetized wind which accompanies the first tens of seconds in the life of a
ms magnetar can feed a relativistic jet with the necessary luminosity and
baryon loading (Lorentz factor) to power a GRB (e.g., Wheeler et al. 2000;
Thompson et al. 2004; Bucciantini et al. 2008; Metzger et al. 2011), though
the accretion-powered jet of a Kerr black hole remains a more commonly invoked
scenario (e.g., Woosley 1993; MacFadyen & Woosley 1999; Gottlieb et al. 2023).

What distinguishes the central engines of GRBs versus those that power SLSNe?
Although the total energetics of both systems can be comparable, their main
distinction is the duration of the central engine activity (Metzger et al.,
2015). While GRB engines typically emit a sizable fraction of their energy
(Panaitescu & Kumar, 2002) over seconds to minutes (comparable to the duration
of the gamma-ray emission; e.g., Bromberg et al. 2012), the engines of SLSNe
must instead heat the supernova ejecta continuously from within over much
longer timescales of days to weeks. This heating being substantially delayed
after the time of the stellar explosion is crucial to explaining the high
radiative output of SLSNe. A delay ensures that less deposited energy is
“wasted” through adiabatic losses by occurring at later times after the
optical depth through the expanding ejecta has dropped. Fig. 1 shows estimates
of the engine duration and luminosity required to explain the rise times and
peak luminosities of a large sample of SLSNe-I (Gomez et al., 2024).

In magnetar scenarios, the engine duration is typically related to the
magnetic dipole spin-down timescale $t_{\rm sd}\propto B^{2}$ which varies
considerably with the surface magnetic field strength ($B\sim 10^{14}$ G to
power SLSNe vs. $B\gtrsim 3\times 10^{15}$ G to power a typical long GRB).
Even at fixed $B$, differences in the magnetar’s mass accretion rate (Metzger
et al., 2018) or the inclination angle of the magnetic dipole relative to the
spin axis (Margalit et al., 2018) can lead to large differences in the jet
power. So-called “ultra-long” GRBs (e.g., Levan et al. 2014) require engine
durations of hours, commensurate with their longer gamma-ray emission; these
events may represent an intermediate case in which the supernova’s luminosity
is boosted only modestly above the baseline level set by 56Ni decay (Greiner
et al., 2015; Metzger et al., 2015), the latter dominating in ordinary GRB
supernovae (e.g., Woosley et al. 2002).

While a variety of indirect evidence supports a compact object residing at the
center of at least some SLSNe-I (e.g., Lunnan et al. 2015; Nicholl et al.
2019; Eftekhari et al. 2019; Vurm & Metzger 2021; Omand & Jerkstrand 2023), a
central engine is not the only plausible power source for these events. A
commonly invoked alternative is shock interaction between the fast supernova
ejecta and a compact wind or other circumstellar medium (CSM) that surrounds
the progenitor star at the time of explosion (e.g., Chevalier & Irwin 2011;
Chatzopoulos et al. 2013; Lunnan et al. 2018; Fraser 2020; Aguilera-Dena et
al. 2020; Stevance & Eldridge 2021). Given the challenges in cleanly
distinguishing engine from CSM interaction models, alternative tests of the
engine hypothesis are sorely needed. Gamma-ray emission from a magnetar nebula
on timescales of months to years after the ejecta becomes optically-thin
offers one compelling test (e.g., Murase et al. 2016; Vurm & Metzger 2021);
however, the low predicted gamma-ray flux at these late times makes detection
challenging in all but the closest explosions (Acharyya et al., 2023). A more
promising test is one that probes the central engine at earlier times closer
to its maximal power.

If a rotating compact object powers both GRBs and SLSNe, it is natural to ask
whether SLSNe should also be accompanied by ultra-relativistic jets. While the
luminosities of SLSN engines are lower than GRB jets, the density of the
external medium through which a putative SLSN jet must propagate is also
lower. This is due to the larger ejecta radius $R_{\rm ej}\approx\,{v_{\rm
ej}}t\sim 10^{13}-10^{14}$ cm at such late times $t\sim$ days$-$weeks after
the explosion. Here, $\,{v_{\rm ej}}\approx 5000$ km s-1 is the speed of the
supernova ejecta, which for typical stripped-star progenitor radii
$R_{\star}\sim 10^{11}$ cm, achieves homologous expansion on a timescale
$\gtrsim R_{\star}/\,{v_{\rm ej}}\sim 10^{2}$ s. Analytic estimates show that
even jet luminosities $\,{L_{\rm j}}\simeq 10^{45}$ erg s-1, much lower than
for GRBs but comparable to those required of many SLSNe-I engines (Fig. 1),
are sufficient for jet breakout through the homologously expanding stellar
ejecta (e.g., Quataert & Kasen 2012; Margalit et al. 2018; Gottlieb & Nakar
2022). In particular, for a jet of half-opening angle $\theta_{j}$ and
luminosity $\,{L_{\rm j}}$, breakout can occur on a timescale (Eq. 45 in
Gottlieb & Nakar, 2022):

$t_{\rm bo}\approx 1.5\,\frac{E_{\rm KE}}{10^{51}\,\,{{\rm
erg}}}\frac{10^{45}\,\,{{\rm erg}}\,\,{{\rm s}}^{-1}}{\,{L_{\rm
j}}}\left(\frac{\theta_{j}}{10^{\circ}}\right)^{4}\,\,{{\rm d}}\,,$ (1)

where the prefactor assumes a non-magnetized jet propagation, as expected for
extended-duration weak jets, and other quantities including the ejecta kinetic
energy $E_{\rm KE}$ are normalized to characteristic values (Fig. 1). The jet
breakout times of days indicated by Eq. (1) are notably similar to the typical
engine durations (e.g., magnetar spin-down time) required to power the bulk of
the SLSNe emission.

Figure 1: Required engine peak luminosity $L_{\rm eng}$ and duration at this
luminosity, $t_{\rm eng}$ (e.g., magnetar spin-down luminosity and timescale,
respectively) to explain the observed peak luminosities, $L_{\rm pk}$, and
rise times, $t_{\rm pk}$, of a large sample of SLSNe-I light curves (Gomez et
al., 2024), for two different assumptions about the ejecta mass, $M_{\rm
ej}=10M_{\odot}$ (black stars) and $3M_{\odot}$ (blue squares), bracketing the
typical range inferred for Type I SLSNe (Nicholl et al., 2017). For
comparison, we show the approximate range of timescales and luminosities of
the early pre-maximum bumps observed in the light curves of SLSNe-I (shaded
blue ellipse; see e.g. Nicholl & Smartt 2016). A multiday-lasting jet from the
same central engine responsible for powering the bulk of the SLSNe is thus
capable of powering the early bumps, assuming the jet luminosity, $\,{L_{\rm
j}}$, is a fraction $\sim 1-10\%$ of $L_{\rm eng}$. The assumed jet properties
for the numerical simulation presented in this paper are illustrated with an
orange square, while left-right arrows indicate the degeneracy in engine
properties $L_{\rm eng}\propto t_{\rm eng}^{-1}$ that result in the same
hydrodynamic evolution for fixed ejecta properties. Details: For each event,
we estimate the ejecta kinetic energy, $E_{\rm KE}=M_{\rm ej}\,{v_{\rm
ej}}^{2}/2$, by equating the observed supernova rise-timescale to the
diffusion timescale for a centrally-concentrated energy source, $t_{\rm
diff}\simeq 0.4(\kappa M_{\rm ej}/\,{v_{\rm ej}}c)^{1/2}$ (Khatami & Kasen,
2019). Assuming the explosion energy $E_{\rm KE}$ comes exclusively from the
engine, i.e. $E_{\rm eng}=L_{\rm eng}t_{\rm eng}=E_{\rm KE}$, we estimate the
required engine timescale assuming the engine luminosity decays similar to
magnetic-dipole spin-down $L_{\rm pk}=(E_{\rm eng}/t_{\rm eng})(1+t_{\rm
pk}/t_{\rm eng})^{-2}$, for $t_{\rm eng}<t_{\rm pk}$. If the engine
contributes only a fraction of the explosion energy $E_{\rm KE}$, then a
greater engine duration is needed to explain the observed luminosities and
points will shift towards the lower-right-hand part of the diagram as
indicated by the arrow.

Do observations rule out jets in SLSNe? With one interesting exception (Levan
et al. 2013), no X-ray emission has been detected to faint limits from a
sample of 26 SLSNe-I on timescales of months to years after the explosion
(e.g., Margutti et al. 2018). However, this is not necessarily constraining
because the spectral energy distribution of any (internal) jetted emission is
uncertain; and even if the jet efficiently powers X-rays, they may only be
detectable at early times and for observers within the narrow opening angle of
the jet. In principle, radio afterglow emission from the jet as it interacts
with surrounding CSM may be detectable at late times even for slightly off-
axis observers; however, while again (with one or two interesting exceptions;
Eftekhari et al. 2019; Margutti et al. 2023111In PTF10hgi (Eftekhari et al.,
2019; Mondal et al., 2020), radio emission detected starting 7.5 years after
the explosion may originate from a magnetar nebula (e.g., Margalit & Metzger
2018; Omand et al. 2018), while the late radio emission in SN 2017ens
(Margutti et al., 2023) originates from the supernova ejecta colliding with
CSM.) no radio detections of SLSNe-I have been made (e.g., Coppejans et al.
2018; Hatsukade et al. 2021; Eftekhari et al. 2021). As we shall discuss, the
resulting constraints on the jet power are sensitive to the observer viewing
angle, CSM density, and uncertain microphysical parameters of the shock.

A cleaner signature of a successful jet in SLSNe-I may come at earlier times,
on timescales comparable to the engine duration of days to weeks (Margalit et
al., 2018). In ordinary GRB jets, the compact size of the shocked stellar
debris (the so-called “cocoon”; e.g., Gottlieb et al. 2018a) results in most
of the cocoon radiation being initially trapped and only subsequently released
after many expansion times. By contrast, for the more sustained jet-ejecta
interaction expected for longer-lived SLSNe jets interacting with
substantially larger ejecta, the lower optical depths (smaller trapping radius
in the ejecta frame) will alter the observed cocoon light curve, allowing it
to more faithfully track energy input from the engine.

A suggestive hint comes from the early maxima (“bumps”) of luminosity $\sim
10^{43}-10^{45}$ erg s-1 observed in some SLSNe-I optical light curves on
timescales of days to a week after the explosion, well before the main
supernova peak (Leloudas et al., 2012; Nicholl et al., 2015). Although such
pre-maximum bumps are not detected in all events due to observational
limitations, Nicholl & Smartt (2016) show that they are likely common features
in SLSNe. Subsequent studies have shown that early maxima are not present in
all SLSNe, and that their properties are heterogeneous, with a wide range of
rise and decline times (e.g., Angus et al. 2019). Qualitatively similar
“double peaked” light curves are seen in a fraction $\sim 3-9\%$ of ordinary-
luminosity stripped-envelope supernovae (e.g., Das et al. 2023 and references
therein), which however are frequently attributed to shock interaction between
the supernova ejecta and CSM from pre-explosion mass-loss.

It may not be a coincidence that the timescales of the observed bumps coincide
both with the first light from the jet breakout [Eq. (1)] and the engine
durations needed to power SLSNe emission. The luminosities of the bumps can
also be accommodated within the energy budget of the engine (light blue cloud
in Fig. 1). Kasen et al. (2016) proposed this early emission phase is powered
by ejecta re-heating due to a delayed spherical shock breakout through the
ejecta, distinct from the earlier supernova shock breakout and driven by the
pressure of a pulsar-like nebula inflated inside the ejecta by the magnetar
wind (see also Chen et al. 2016; Suzuki & Maeda 2016).

In this _Letter_ , we investigate whether this early emission phase can
instead be powered by a collimated jet$-$from an engine with similar
properties needed to power the full SLSNe light curve at later times$-$that
successfully breaks out of the supernova ejecta. In §2.1, we establish a
numerical simulation to track the jet propagation inside and outside the
ejecta over a few weeks. In §2.2, we outline the post-processing calculations
to estimate the resulting emission. In §3, we present our findings,
demonstrating their consistency with observations of the early bumps in
SLSNe-I. We discuss some implications of our results in §4.

##  2 Numerical calculation

###  2.1 Simulation setup

We carry out a 3D relativistic hydrodynamic simulation of a jet launched into
an expanding stellar envelope, using the code pluto (Mignone et al., 2007). We
follow the jet propagation in the ejecta, its activity after breakout, and
after it shuts off. We post-process the simulation output semi-analytically to
produce estimates of the post-breakout emission. As the simulations are
dimensionless, we can calibrate the output with dimensionless mass scale
$m_{s}$ and timescale (length-scale) $t_{s}$.

The jet operates for $\,{t_{\rm eng}}=4460\,t_{s}\,\,{{\rm s}}$ with a total
(two-sided) luminosity, $\,{L_{\rm j}}=2\times 10^{47}\,m_{s}/t_{s}\,{\rm
erg\,\,{{\rm s}}^{-1}}$. We consider as fiducial values
$\\{t_{s}=80,m_{s}=1\\}$, corresponding to jet duration $\,{t_{\rm
eng}}\approx 4.1\,\,{{\rm d}}$, and luminosity $L=2.5\times 10^{45}\,{\rm
erg}\,\,{{\rm s}}^{-1}$, typical of the engine properties needed to power
SLSNe-I (Fig. 1; see also Metzger et al. 2015, their Fig. 1). These jet
properties are illustrated as an orange square in Fig. 1, while left-right
arrows indicate the change in engine properties obtained by varying $t_{s}$
($L_{\rm eng}\propto t_{\rm eng}^{-1}$) that would result in the same
hydrodynamic evolution (for fixed ejecta properties $m_{s}$). This degeneracy
running approximately parallel to the SLSNe population suggests that the
hydrodynamic evolution predicted by our single simulation could apply to a
wide range of SLSNe engines.

Building on the results of Gottlieb et al. (2022a), we set the jets at the
inner boundary to be precessing about the spin-axis of the black hole, with an
inclination angle of $0.05$ rad and a period of $t_{s}\,\,{{\rm s}}$. The
initial jet Lorentz factor is $\Gamma_{0}=10$ with a maximum Lorentz
factor222We note that the precise ultra-relativistic Lorentz factor of the jet
does not affect the cocoon, which maintains $\Gamma\lesssim 3$ at all times
(Gottlieb et al., 2021). of 20 (equipartition between initial thermal and
kinetic energy), and its opening angle is $\theta_{0}=0.1$ rad. The jet runs
into a homologously expanding $M_{\rm ej}=10\,m_{s}\,{M_{\odot}}$ ejecta with
a front velocity of $\,{v_{\rm ej}}=10^{4}\,\,{{\rm km}}\,\,{{\rm s}}^{-1}$
and initial position $\,{r_{\rm ej,0}}=7\times 10^{6}\,t_{s}\,\,{{\rm km}}$.
The ejecta initial comoving mass density profile is:

$\rho(r)=\left(\frac{r}{9\times 10^{5}\,\,{{\rm
km}}}\right)^{-1.5}\left(1-\frac{r}{\,{r_{\rm
ej}}}\right)^{3}\frac{m_{s}}{t_{s}^{3}}\frac{\rm g}{\rm cm^{3}}\,.$ (2)

We use an ideal gas equation of state with a polytropic index of 4/3, as
appropriate for radiation-dominated gas. This applies at all relevant times
for our calculation as we focus on the outflow while radiation is coupled to
the gas. The 3D Cartesian grid is divided into 3 patches on the $\hat{x}$\-
and $\hat{y}$-axes: The inner patch is at the innermost $2\times
10^{5}t_{s}\,\,{{\rm km}}$ with 100 uniformly distributed cells. The outer
patches are stretched from $\lvert 2\times 10^{5}\,t_{s}\,\,{{\rm km}}\rvert$
to $\lvert 3\times 10^{8}\,t_{s}\,\,{{\rm km}}\rvert$ with 350 logarithmically
distributed cells. We divide the $\hat{z}$-axis (jet axis) into 4 patches:
starting from $z_{0}=2\times 10^{5}\,t_{s}\,\,{{\rm km}}$ to $3\times
10^{6}\,t_{s}\,\,{{\rm km}}$ with 300 uniformly distributed cells, followed by
3 logarithmically spaced patches with 600 cells to $3\times
10^{7}\,t_{s}\,\,{{\rm km}}$, 350 cells to $9.5\times 10^{7}\,t_{s}\,\,{{\rm
km}}$, and 480 cells to $6\times 10^{8}\,t_{s}\,\,{{\rm km}}$. The jet is
injected at the lower boundary, $z_{0}$, with a cylindrical injection radius
of $r_{0}=z_{0}\theta_{0}$. The integration is performed with a Harten-Lax-van
Leer (HLL) solver, Runge-Kutta time stepping, and piecewise parabolic
interpolation.

###  2.2 Emission calculation

We first consider the released thermal near-ultraviolet (NUV)/optical emission
from the gas that escaped the trapping radius in the last dynamical time. We
assume that each fluid element is propagating radially at the relevant times
of emission, such that we calculate the trapping radius along each line of
sight $\Omega$, and time $t$:

$\,{R_{\rm tr}}(t,\Omega)=R\left[\tau(t,\,{R_{\rm
tr}},\Omega)=\frac{1}{\beta(t,\,{R_{\rm tr}},\Omega)}\right]\,,$ (3)

where $\beta$ is the dimensionless velocity, and $\tau$ is the optical depth
along a radial line of sight:

$\tau(t,r,\Omega)=\int_{r}^{\infty}\kappa\rho(t,r^{\prime},\Omega)\left[1-\beta\left(t,r^{\prime},\Omega\right)\right]\gamma(t,r^{\prime},\Omega)\mathrm{d}r^{\prime}\,,$
(4)

where $\kappa=0.1\,\,{{\rm cm}}^{2}\,{\rm g}^{-1}$ is the opacity of the gas
(e.g., Kleiser et al. 2018), and $\gamma$ is its Lorentz factor.

For two given lab times, $t_{1},t_{2}$, we find the trapping radius in the
Lagrangian frame, $\,{R_{\rm tr}}^{L1},\,{R_{\rm tr}}^{L2}$. We calculate the
thermal energy between these radii as:

$\,{E_{\rm th}}\left(\bar{t},\Omega\right)=\int_{\,{R_{\rm
tr}}^{L2}}^{\,{R_{\rm
tr}}^{L1}}{3p(\bar{t},r,\Omega)\gamma(\bar{t},r,\Omega)r^{2}dr}\,,$ (5)

where $p$ is the comoving gas thermal pressure and
$\bar{t}=\frac{t_{1}+t_{2}}{2}$. As the emitting gas is sub-relativistic with
$\,{\beta_{\rm emit}}\equiv\beta(\,{R_{\rm tr}})\ll 1$, we neglect
relativistic effects for the luminosity and the relative arrival times of the
photons. The thermal energy is released during $t_{2}-t_{1}$, such that the
local luminosity in the lab frame is:

$L(\bar{t},\Omega)=\frac{\,{E_{\rm
th}}\left(\bar{t},\Omega\right)}{t_{2}-t_{1}}\,,$ (6)

and the bolometric luminosity is:

$L_{\rm bol}(\bar{t})=\int L(\bar{t},\Omega)d\Omega\,.$ (7)

Eq. (3) implies that the optical depth of the non-relativistic emitting shell
is always much greater than unity. Thus, the diffusing photons will interact
with free electrons before reaching the photosphere. These interactions will
govern the photon energy at the photosphere, defined as $\,{R_{\rm
ph}}(t,\Omega)=R\left[\tau(t,\,{R_{\rm ph}},\Omega)=1\right]$. At the relevant
times and temperatures, there is enough time for the gas to produce photons
through free-free emission to share the photon energy and reduce the
temperature from the trapping radius until reaching thermalization (see
discussion in Nakar & Sari, 2010). Therefore, the local temperature at the
photosphere is taken to be a blackbody:

$T(\bar{t},\Omega)=\left[\frac{L\left(\bar{t},\Omega\right)}{4\pi\sigma_{b}\,{R_{\rm
ph}}\left(\bar{t},\Omega\right)^{2}}\right]^{1/4}\,.$ (8)

Using the local temperature, we calculate the local spectral luminosity,
$L_{\nu}$, at various bands, and integrate it over all angles to obtain the
observed spectral luminosity. From the spectral luminosity, we infer the color
temperature as:

$T_{c}(\bar{t})=\frac{h\nu_{\rm max}(\bar{t})}{2.821\,k_{B}}\,,$ (9)

where $h$ and $k_{B}$ are the Planck and Boltzmann constants, respectively,
and $\nu_{\rm max}(\bar{t})={\rm max}\left[\nu L_{\nu}(\bar{t})\right]$.

Finally, as the jet interacts with the high-density CSM, it powers afterglow
synchrotron emission. We follow the semi-analytic solution in the thin-shell
approximation of Oren et al. (2004) for the propagation of a top-hat jet in a
wind medium, to calculate the resulting radio emission.

##  3 Emission

###  3.1 NUV/optical

Figure 2:  3D rendering of the logarithmic comoving mass density at the end of
the jet activity, $t=\,{t_{\rm eng}}$, when the ejecta front is at $\,{r_{\rm
ej}}\approx 3.5\times 10^{14}\,\,{{\rm km}}$. As the jet escapes from the
high-density expanding spherical ejecta (yellow), it generates a stratified
structure of the lower-density cocoon (red). The interaction of the jet and
the cocoon with the expanding ejecta close to the breakout radius generates
bright NUV/optical emission over a $\sim$ week timescale.

Figure 3: 2D $\hat{x}-\hat{z}$-plane cuts at various times of the azimuthally-
averaged logarithmic comoving thermal energy density $3p/c^{2}$ maps. The
white and black contours represent the photosphere and trapping radius,
respectively.

Figure 4:  (a): Bolometric luminosity peaks soon after breakout ($t_{b}$), and
then declines as $\sim t^{-2}$ (dashed black line). The dashed magenta lines
illustrate the rise of the main SLSN-I peak. (b): Color temperature of the
radiating gas. (c): Absolute (left vertical axis) and apparent (right vertical
axis) magnitude at various bands ($f_{\rm NUV}=1.2\times 10^{15}\,{\rm Hz}$,
$f_{g}=6.3\times 10^{14}\,{\rm Hz}$, $f_{r}=4.8\times 10^{14}\,{\rm Hz}$,
$f_{y}=3.0\times 10^{14}\,{\rm Hz}$, $f_{K_{s}}=1.4\times 10^{14}\,{\rm Hz}$)
demonstrate that the peak magnitude aligns with observations of SLSN early
bumps. (d): Luminosity-weighted average dimensionless velocity of the
radiating gas at the trapping radius demonstrates that the emitting gas is
always subrelativistic. The dimensionless front ejecta velocity is shown for
comparison in red.

The head of the weak jet propagates slowly inside the expanding ejecta with an
average velocity in the ejecta frame of $v_{h}-\,{v_{\rm ej}}\approx 0.01c$.
This implies that all the jet energy is practically deposited into the cocoon
for as long as the jet head remains within the ejecta. For our choice of
parameters: $t_{s}=80$ and $m_{s}=1$, the two-sided jet with a luminosity of
$\,{L_{\rm j}}=2.5\times 10^{45}\,\,{{\rm erg}}\,\,{{\rm s}}^{-1}$ breaks out
from the $M_{\rm ej}=10\,\,{M_{\odot}}$ ejecta after $t_{b}\approx
2.2\,\,{{\rm d}}$, indicating that the cocoon energy upon breakout is
$L_{c}\approx\,{L_{\rm j}}t_{b}\approx 5\times 10^{50}\,\,{{\rm erg}}$. Figure
2 illustrates a 3D rendering of the logarithmic comoving mass density,
showcasing the interaction of the jet-cocoon outflow (red) with the
homologously expanding ejecta (yellow). This interaction generates emission
near the optically thin breakout radius and shapes the large-scale stratified
structure observed beyond the confined expanding ejecta.

To elucidate the different phases of the evolution and emission, Figure 3
presents azimuthally-averaged logarithmic comoving thermal energy density maps
at various times, where the black and white contours delineate the trapping
radius and photosphere, respectively. Figure 4 displays the resultant $L_{\rm
bol}$ (a), $T_{c}$ (b), absolute and apparent magnitudes (c), and
$\,{\beta_{\rm emit}}$ (d). using the numerical calculation described in §2.2.
We identify three phases of emission:

(1) The first emission phase begins upon breakout and continues for as long as
the jet engine is active at $t_{b}<t<t_{\rm eng}$. Fig. 3(a) depicts the
breakout of the collimated jet and its hot, pressurized cocoon from the
expanding ejecta. Upon breakout, photons above the trapping radius escape from
the top layers of the cocoon, producing the first light through shock breakout
emission. At this time, the emitting layers are mildly relativistic, as
indicated by the proximity of the trapping radius to the photosphere. Shortly
after the shock breakout signal, the gas from deeper layers starts to
contribute to the emission with $\,{\beta_{\rm emit}}\approx 0.2$ [Fig. 4(d)],
and the luminosity begins to fall off [Fig. 4(a)]. The high temperature after
breakout of $T_{c}\approx 2\times 10^{4}\,{\rm K}$ [Fig. 4(b)] results in a
peak in the NUV [Fig. 4(c)].

After the breakout, the jet head accelerates to ultra-relativistic velocities,
and the mildly relativistic cocoon expands in the angular direction to
$\theta\sim(\Gamma\beta)^{-1}$. In typical GRBs, the breakout takes place
within $\sim 100\,\,{{\rm s}}$, implying that the photons generated in the
jet-cocoon interaction are still trapped well within the outflow. Here, Figure
3(b), which depicts the jet-cocoon system at the time of jet engine shut-off
at $t_{\rm eng}=4.1\,\,{{\rm d}}$, suggests that the jet is optically thin,
and the late breakout results in a trapping radius located just outside of the
jet-cocoon interface (JCI; Gottlieb et al., 2021) and within the ejecta outer
radius. Consequently, the emission originates from the continuous heating in
the JCI, as indicated by the high thermal energy along the increasing trapping
radius in the angular direction away from the jet axis. Over time, the drop in
the temperature [Fig. 4(b)] results in a transition to softer emission [Fig.
4(c)].

(2) Figures 3(c,d) demonstrate that after the jet shuts off, the polar cavity
is loaded with baryon and turns optically thick. At the same time, the
emission from the jet-cocoon interaction ceases, resulting in a drop in the
bolometric luminosity at $t\approx t_{\rm eng}$. At this point, the thermal
energy along the angular direction of the trapping radius becomes roughly
homogeneous. The dominant heat source then becomes the newly shocked material
from the lateral spreading of the cocoon into the unshocked ejecta below the
trapping radius, with $\,{\beta_{\rm emit}}\lesssim\beta_{\rm ej}$, where
$\beta_{\rm ej}$ is the dimensionless front ejecta velocity [Fig. 4(d)]. As
the heat source is located below the trapping radius, the emission follows a
cooling-envelope model, with the light curve evolving as $L\sim t^{-2}$, as
expected from the expansion of the sub-relativistic gas in the cocoon
(Gottlieb et al., 2022b).

The high temperatures around the time of breakout cause the emission to first
peak in the NUV and then in the optical bands, as the latter overlaps with the
blackbody peak while the jet is active. As the luminosity declines more
rapidly than the temperature at $t\gtrsim t_{\rm eng}$, the NIR $K_{s}$ band
emission rises since it is located on the Rayleigh-Jeans portion of the
spectrum. Consequently, the decrease in the bolometric luminosity following
jet shutdown also reduces the optical emission.

(3) Comparing Figures 3(c) and 3(d), it is evident that the trapping radius
shrinks in the ejecta frame as the ejecta becomes optically thin. The final
phase of the emission will occur once the cocoon crosses the entire unshocked
ejecta, and the entire ejecta becomes optically thin. The timing of this
transition depends on the opacity, which may decrease with temperature,
depending on the composition of the supernova ejecta. In particular, if the
ejecta is primarily comprised of O-Ne-Mg nuclei, the electron scattering and
line-opacity drops substantially once $T\approx 6\times 10^{3}\,{\rm K}$, at
$t\approx 5\,\,{{\rm d}}$ due to recombination effects (e.g., Kleiser et al.
2018). However, if the line opacity of Fe-group nuclei (e.g., Fe-Co-Ni)
dominates, the opacity will remain relatively constant across this temperature
range, causing the optically thin transition to occur somewhat later. We note
that the relevant opacity for the emission is determined below the trapping
radius, where the temperature is $T_{t}=T_{c}\,{\beta_{\rm
emit}}^{-0.25}\approx 2T_{c}$, implying that $T_{t}\gtrsim 6\times
10^{3}\,{\rm K}$ at $t\leq 16\,{\rm d}$. This suggests that the numerical
calculation remains valid at the relevant times.

Ultimately, the rise of the main supernova emission, powered primarily by (a
more isotropic, non-jet-related) source of ejecta heating by the same central
engine (e.g., the absorption of X-rays and gamma-rays from the magnetar
nebula; Vurm & Metzger 2021), will take over the optical light curve. The
dashed magenta curves in the bolometric luminosity panel illustrate the rise
of the SLSN-I emission towards the main peak at $L_{\rm bol}\approx
10^{44}\,\,{{\rm erg}}\,\,{{\rm s}}^{-1}$ at $t\approx 50\,\,{{\rm d}}$ (e.g.,
Quimby et al., 2011; Gal-Yam, 2012). We present both a slow rise, $L_{\rm
bol}\sim t$, and a fast rise, $L_{\rm bol}\sim t^{2}$.

Our light curve calculation considers the integrated emission over all solid
angles. However, an observer will only detect the whole emission within
certain viewing angles relative to the jet axis, where most of the thermal
energy is generated. The first light – shock breakout emission [Figure 3(a)]
is quasi-isotropic with mild beaming due to the mildly relativistic velocities
of the front cocoon. After the initial signal, the emitting gas becomes
subrelativistic [Figure 4(d)], so it is no longer subject to relativistic
beaming effects. However, the emitted photons can only reach observers along
optically thin lines of sight, i.e. “geometric” beaming may still play a role.
Figure 3(b) shows that the phase in which the jet is active is only observable
within the JCI at $\theta_{\rm obs}\lesssim 2\theta_{j}$. Over time, as the
cocoon expands angularly, the accessible angle gradually increases. Therefore,
the light curve observed at angles farther from the jet axis may peak at shock
breakout, decay, and rise again at later times compared to the bolometric
luminosity shown in Figure 4(a). Full radiation transfer calculations are
needed to provide a self-consistent modeling of the light curves at different
viewing angles.

###  3.2 Radio

Ultimately, the jet will decelerate due to its interaction with the dense CSM,
powering an “orphan” radio afterglow synchrotron emission and an optical flash
from the reverse shock (Sari & Piran, 1999). However, the latter is likely to
be outshone by the thermal cocoon emission associated with the early bump.
While the lack of orphan afterglow detections in blind surveys (Ho et al.,
2018, see however Perley et al. 2024 for a recent potential detection), along
with the potential for high CSM densities around some SLSNe-I (e.g., Margutti
et al., 2023) that would amplify the signal, may argue against the generic
presence of jets in these systems, the large uncertainties prevent us from
reaching definitive constraints.

Figure 5:  Top-hat jet afterglow light curves at $\nu=6\,{\rm GHz}$ for
observers at $\theta_{\rm obs}=0.5\,{\rm rad}$ (solid lines) and on-axis
observers (dashed lines) at $z=0.25$, assuming jet Lorentz factor $\Gamma=100$
and opening angle $\theta_{j}=0.1\,{\rm rad}$, microphysics parameters:
electron equipartition parameter $\epsilon_{e}=0.1$, magnetic equipartition
parameter $\epsilon_{B}=10^{-4}$, and electron energy power-law index $p=2.2$,
and wind velocity $v_{w}=1000\,\,{{\rm km}}\,\,{{\rm s}}^{-1}$. The black dots
denote the most constraining upper limits measured for SLSN-I Gaia2016apd
(Coppejans et al., 2018), placing the source at $z=0.25$.

Figure 5 presents the semi-analytic thin-shell approximation solution for the
radio synchrotron emission of our model at various assumed progenitor wind
mass-loss rates, for both off-axis (solid lines) and on-axis (dashed lines)
viewing angles. We assume a characteristic redshift of $z=0.25$ for SLSNe-I
(e.g., Gomez et al., 2024), standard jet properties with $\Gamma=100$ and
$\theta_{j}=0.1\,{\rm rad}$, and microphysics parameters (detailed in the
caption) inferred from GW170817 afterglow and VLBI images (Mooley et al.,
2018). If the mass-loss rates are $M_{w}\gtrsim 10^{-5}\,{\rm yr}^{-1}$, the
emission is self-absorbed during the first $\sim 1$ week, leading to an
initial rise in the light curves for all observers.

The extended engine activity indicates that the jet’s initial radial extent is
$\sim(t_{\rm eng}-t_{b})c\sim 10^{16}\,\,{{\rm cm}}$. This suggests that not
all of the jet’s energy is transferred to the forward shock before the jet’s
width becomes comparable to the shock width, $\sim R/12\,\Gamma^{2}$, at an
observed time of $\sim 12\,(t_{\rm eng}-t_{b})\sim 1$ month. This duration
could be shorter if, for example, the jet’s radial extent shrinks beforehand
due to high CSM densities. As a result, the light curve, calculated under the
thin-shell approximation, may overestimate the radio emission during this
initial period, thus providing an upper limit.

Beyond the general uncertainties in the afterglow and jet parameters, the
strong dependence of the flux on the viewing angle (Nakar et al., 2002) and
the wind mass-loss rate complicates robust predictions that can support or
refute our model. The most stringent upper limits on SLSNe-I were obtained for
Gaia2016apd (Coppejans et al., 2018) at $z=0.1018$ (Yan et al., 2017). The
black dots represent these upper limits (adjusted for $z=0.25$), potentially
ruling out some scenarios with high mass-loss rates and small viewing angles
for this source, but not excluding off-axis observers and lower mass-loss
rates within the considered parameters.

##  4 Discussion

Cocoon emission following the birth of a compact object has been invoked in
several energetic transients, including collapsars (e.g., Nakar 2015; Nakar &
Piran 2017) and neutron star mergers (e.g., Lazzati et al. 2017; Gottlieb et
al. 2018a, b). However, in general, the cocoon emission in these systems is
radiatively inefficient, such that a significant fraction of the deposited jet
energy goes into adiabatic losses rather than observable radiation. By
contrast, the longer duration of the engines at work in SLSNe, which inject
the bulk of their energy over days to weeks instead of seconds to minutes, and
correspondingly much lower densities of the supernova ejecta at these late
times, reduce the impact of adiabatic losses. In particular, the trapping
radius moves inside the JCI on a timescale comparable to the engine duration.
This enables an order-unity fraction of the jet luminosity to emerge as
radiation.

Such a high radiative efficiency is important because the luminosities of the
early bumps seen in SLSN light curves $\sim 10^{43}-10^{45}$ erg s-1
correspond to a sizable fraction $\sim 1-100\%$ of the total energy budget of
the engine (Fig. 1). In the case of a magnetar engine, the fraction of the
magnetar spin-down energy placed into a jet depends on several uncertain
conditions, such as the inclination angle of the magnetic dipole (Margalit et
al. 2018, which determines the fraction of the striped pulsar wind that
undergoes magnetic reconnection at the wind termination shock; e.g., Lyubarsky
2003; Komissarov 2013) and the efficiency with which any remaining large-scale
field avoids instabilities and associated reconnection (Bucciantini et al.,
2009; Porth et al., 2013; Mösta et al., 2015). In the case of a black hole
engine, the power of the relativistic jet depends on the magnetic field and
spin of the black hole (Blandford & Znajek, 1977), while accretion disk
outflows can release their energy more isotropically and contribute to
powering the bulk of the SLSN emission (e.g., Dexter & Kasen 2013). If the jet
engine operates for $\,{t_{\rm eng}}\lesssim t_{b}$, then the jet will be
choked deep inside the ejecta, resulting in no powerful breakout and hence
absent or weaker associated UV/optical emission. A substantial failed jet
fraction is compatible with early maxima being detected in only 3/14 SLSNe-I
observed by DES (Angus et al., 2019).

For several days after the jet emerges from the star, the cocoon emission is
hot $T_{\rm c}\approx(1-2)\times 10^{4}$ K and peaks in the optical/NUV at an
absolute magnitude $M_{\rm AV}\approx-21.5$. Although multi-band measurements
of early SLSN bumps which constrain their temperatures are sparse, such data
on PS1-10pm indeed indicated high $\approx(2-3)\times 10^{4}$ K temperatures
during its early bump phase (Nicholl & Smartt, 2016). Similar temperature
evolution was observed from the early maximum in DES14x3taz, which cooled from
$\approx 2.5\times 10^{4}$ K at discovery to $\lesssim 10^{4}$ K over roughly
a week (Smith et al., 2016).

A few wide-field UV satellites with time-domain capabilities are planned over
the next decade, including the Ultraviolet Transient Astronomy Satellite
(ULTRASAT, launch date 2027; Sagiv et al. 2014), Ultraviolet Explorer (UVEX;
Kulkarni et al. 2021), and the Czech mission QUVIK (Werner et al., 2023). We
focus on ULTRASAT, which will reach a $5\sigma$ sensitivity of 22.4 AB
magnitude in the NUV band ($h\nu\approx$ 5 eV) across an instantaneous field
of view of $\approx 200$ deg2 for a 900 second (15 minutes) integration (Sagiv
et al., 2014).333 UVEX will possess sensitivity extending also into the FUV
and reaching an AB magnitude depth of 24.5 for a 900 s integration (Kulkarni
et al., 2021), roughly a factor of 6 deeper than ULTRASAT. However, UVEX’s
smaller instantaneous field of view $\approx 12$ deg2, results in a comparable
survey speed time to ULTRASAT. Most of the observing time will be spent on a
low cadence survey, cycling through $\Delta\Omega\sim 6800$ deg2
($f_{\Omega}=0.16$ of the whole sky) covering 10 fields of view per day (4-day
cadence per field).

For a source with $M_{\rm AB}\approx-21$ lasting a few days, we estimate a
ULTRASAT horizon distance of $D_{\rm lim}\approx 4$ Gpc. We assume the early
NUV emission to be detectable over $\theta_{\rm obs}\approx 0.5\,{\rm rad}$.
Given the volumetric rate of SLSNe-I of $\mathcal{R}\approx 40$ Gpc-3 yr-1
(Quimby et al., 2013; Frohmaier et al., 2021; Zhao et al., 2021) at $z\sim
0.2-1$, and assuming that all SLSNe produce early jet breakout emission
similar to that we have predicted, the number of detectable sources per year
can be estimated as $N_{\rm det}=(2\pi/3)\theta_{\rm obs}^{2}D_{\rm
lim}^{3}\mathcal{R}f_{\Omega}\sim 10^{2}\,{\rm yr}^{-1}$. The detection rates
will be lower if jets are less powerful than our baseline model, or only
successful in a fraction of explosions, as discussed above and supported by
SLSNe-I observations (Angus et al., 2019). More direct rate estimates could be
made based on a large uniform sample of pre-maximum SLSNe light curve
observations (Nicholl & Smartt, 2016). Given its greater depth and similar
cadence, many of the ULTRASAT sources will also have optical light curve data
from the Vera C. Rubin Observatory (Ivezić et al., 2019), enabling better
constraints on the temperature evolution.

An early photometric detection could also enable prompt spectra to be
obtained. We predict the earliest phases of the emission could reveal Doppler-
broadened velocities extending to $v/c\gtrsim 0.1$ (lower right-hand panel of
Fig. 4), at least for observers relatively close to the jet axis. These and
other jet breakout signatures, such as luminous X-ray emission for face-on
observers (possibly already observed in one SLSN; Levan et al. 2013) or radio
afterglow emission from the jet or its cocoon (Coppejans et al. 2018;
Eftekhari et al. 2021; Fig. 5), would provide the ultimate “smoking gun” of an
engine-powered origin for SLSNe. Our model for powering the early bumps in
SLSNe bears some resemblance to the proposal by Kasen et al. (2016), insofar
that we are also invoking shock heating of the supernova ejecta due to energy
input from a central engine which is substantially delayed after the
explosion. However, while Kasen et al. (2016) consider a spherical (or quasi-
spherical) shock breakout driven by the high pressure of a nebula (“bubble”)
inflated by a central engine (see also Suzuki & Maeda 2016), we have
considered the breakout of a tightly collimated jet, closer to the proposal of
Margalit et al. (2018). Furthermore, for the nebula-driven shock to reach the
ejecta surface, the energy of the central engine must exceed the original
kinetic energy of the explosion (Kasen et al., 2016), while collimated jet
breakout from homologous ejecta places a less stringent requirement on the
engine energy (e.g., Gottlieb & Nakar 2022). The successful breakout of a
relativistic jet also comes with additional consequences, such as a
synchrotron afterglow from relativistic ejecta, that in general will not be
present in the spherical breakout case (however, see Arons 2003; Bucciantini
et al. 2011). The outer layers of the ejecta accelerated by a spherical
breakout can produce their own synchrotron emission (e.g., Maeda et al. 2023),
but exhibit a different light curve and spectral evolution than a jet
afterglow, particularly for observers close to the jet axis.

The direct detection of internal high-energy radiation from the jet viewed on-
axis would provide an even more dramatic confirmation of the model. The
mechanism, radiative efficiency, and spectrum of any internal jet emission are
highly uncertain. Nevertheless, we focus on the soft X-ray band, motivated by
the X-ray detection of SCP 06F6 (Levan et al., 2013) and the recent launch of
the satellite mission Einstein Probe, which reaches a sensitivity in the soft
X-ray band across a large fraction of the sky of $F_{\rm lim}\approx 10^{-12}$
erg cm-2 s-1 for a $10^{5}$ s integration (Yuan et al., 2022), comparable to
the expected duration of SLSNe jets. Assuming that a fraction $f_{\rm X}$ of
the jet power $L_{\rm j}$ is placed into soft X-rays, the resulting detection
horizon is $D_{\rm lim}\approx(f_{\rm X}L_{\rm j}/4\pi F_{\rm lim})^{1/2}$.
Adopting a beaming fraction $f_{\rm b}\sim 10^{-2}$ similar to those of GRB
jets, the number of detectable on-axis SLSNe jets per year can be estimated
as:

$N_{\rm X}\approx\frac{4\pi}{3}D_{\rm lim}^{3}\mathcal{R}f_{\rm b}\sim
1.3\,{\rm yr^{-1}}\left(\frac{f_{\rm X}}{0.1}\right)^{3/2}\left(\frac{L_{\rm
X}}{10^{45}\,{\rm erg\,s^{-1}}}\right)^{3/2}$ (10)

It is thus possible that Einstein Probe or successor missions could detect on-
axis SLSNe-I jets over the coming years, though given the many uncertainties,
a lack of such detections is not necessarily constraining.

We thank Sebastian Gomez and Matt Nicholl for sharing data from their recent
compilation of SLSNe light curves, and for helpful comments on the manuscript.
OG is supported by the Flatiron Research Fellowship. BDM was supported in part
by the National Science Foundation (grant No. AST-2009255) and by the NASA
Fermi Guest Investigator Program (grant No. 80NSSC22K1574). The Flatiron
Institute is supported by the Simons Foundation.

## References

  * Acharyya et al. (2023) Acharyya, A., Adams, C. B., Bangale, P., et al. 2023, ApJ, 945, 30, doi: 10.3847/1538-4357/acb7e6
  * Aguilera-Dena et al. (2020) Aguilera-Dena, D. R., Langer, N., Antoniadis, J., & Müller, B. 2020, ApJ, 901, 114, doi: 10.3847/1538-4357/abb138
  * Angus et al. (2019) Angus, C. R., Smith, M., Sullivan, M., et al. 2019, MNRAS, 487, 2215, doi: 10.1093/mnras/stz1321
  * Arons (2003) Arons, J. 2003, ApJ, 589, 871, doi: 10.1086/374776
  * Blandford & Znajek (1977) Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433
  * Bromberg et al. (2012) Bromberg, O., Nakar, E., Piran, T., & Sari, R. 2012, ApJ, 749, 110, doi: 10.1088/0004-637X/749/2/110
  * Bucciantini et al. (2011) Bucciantini, N., Metzger, B. D., Thompson, T. A., & Quataert, E. 2011, ArXiv e-prints.  https://arxiv.org/abs/1106.4668
  * Bucciantini et al. (2008) Bucciantini, N., Quataert, E., Arons, J., Metzger, B. D., & Thompson, T. A. 2008, MNRAS, 383, L25, doi: 10.1111/j.1745-3933.2007.00403.x
  * Bucciantini et al. (2009) Bucciantini, N., Quataert, E., Metzger, B. D., et al. 2009, MNRAS, 396, 2038, doi: 10.1111/j.1365-2966.2009.14940.x
  * Chatzopoulos et al. (2013) Chatzopoulos, E., Wheeler, J. C., Vinko, J., Horvath, Z. L., & Nagy, A. 2013, ApJ, 773, 76, doi: 10.1088/0004-637X/773/1/76
  * Chen et al. (2016) Chen, K.-J., Woosley, S. E., & Sukhbold, T. 2016, ApJ, 832, 73, doi: 10.3847/0004-637X/832/1/73
  * Chevalier & Irwin (2011) Chevalier, R. A., & Irwin, C. M. 2011, ApJ, 729, L6, doi: 10.1088/2041-8205/729/1/L6
  * Coppejans et al. (2018) Coppejans, D. L., Margutti, R., Guidorzi, C., et al. 2018, ApJ, 856, 56, doi: 10.3847/1538-4357/aab36e
  * Das et al. (2023) Das, K. K., Kasliwal, M. M., Sollerman, J., et al. 2023, arXiv e-prints, arXiv:2306.04698, doi: 10.48550/arXiv.2306.04698
  * Dexter & Kasen (2013) Dexter, J., & Kasen, D. 2013, ApJ, 772, 30, doi: 10.1088/0004-637X/772/1/30
  * Eftekhari et al. (2019) Eftekhari, T., Berger, E., Margalit, B., et al. 2019, ApJ, 876, L10, doi: 10.3847/2041-8213/ab18a5
  * Eftekhari et al. (2021) Eftekhari, T., Margalit, B., Omand, C. M. B., et al. 2021, ApJ, 912, 21, doi: 10.3847/1538-4357/abe9b8
  * Fraser (2020) Fraser, M. 2020, Royal Society Open Science, 7, 200467, doi: 10.1098/rsos.200467
  * Frohmaier et al. (2021) Frohmaier, C., Angus, C. R., Vincenzi, M., et al. 2021, MNRAS, 500, 5142, doi: 10.1093/mnras/staa3607
  * Gal-Yam (2012) Gal-Yam, A. 2012, Science, 337, 927, doi: 10.1126/science.1203601
  * Gal-Yam (2018) —. 2018, arXiv e-prints.  https://arxiv.org/abs/1812.01428
  * Gomez et al. (2024) Gomez, S., Nicholl, M., Berger, E., et al. 2024, arXiv e-prints, arXiv:2407.07946, doi: 10.48550/arXiv.2407.07946
  * Gottlieb et al. (2023) Gottlieb, O., Jacquemin-Ide, J., Lowell, B., Tchekhovskoy, A., & Ramirez-Ruiz, E. 2023, ApJ, 952, L32, doi: 10.3847/2041-8213/ace779
  * Gottlieb et al. (2022a) Gottlieb, O., Liska, M., Tchekhovskoy, A., et al. 2022a, ApJ, 933, L9, doi: 10.3847/2041-8213/ac7530
  * Gottlieb & Nakar (2022) Gottlieb, O., & Nakar, E. 2022, MNRAS, 517, 1640, doi: 10.1093/mnras/stac2699
  * Gottlieb et al. (2021) Gottlieb, O., Nakar, E., & Bromberg, O. 2021, MNRAS, 500, 3511, doi: 10.1093/mnras/staa3501
  * Gottlieb et al. (2018a) Gottlieb, O., Nakar, E., & Piran, T. 2018a, MNRAS, 473, 576, doi: 10.1093/mnras/stx2357
  * Gottlieb et al. (2018b) Gottlieb, O., Nakar, E., Piran, T., & Hotokezaka, K. 2018b, MNRAS, 479, 588, doi: 10.1093/mnras/sty1462
  * Gottlieb et al. (2022b) Gottlieb, O., Tchekhovskoy, A., & Margutti, R. 2022b, MNRAS, doi: 10.1093/mnras/stac910
  * Greiner et al. (2015) Greiner, J., et al. 2015, Nature, 523, 189, doi: 10.1038/nature14579
  * Hatsukade et al. (2021) Hatsukade, B., Tominaga, N., Morokuma, T., et al. 2021, ApJ, 922, 17, doi: 10.3847/1538-4357/ac20d5
  * Ho et al. (2018) Ho, A. Y. Q., Kulkarni, S. R., Nugent, P. E., et al. 2018, ApJ, 854, L13, doi: 10.3847/2041-8213/aaaa62
  * Ivezić et al. (2019) Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111, doi: 10.3847/1538-4357/ab042c
  * Kasen & Bildsten (2010) Kasen, D., & Bildsten, L. 2010, ApJ, 717, 245, doi: 10.1088/0004-637X/717/1/245
  * Kasen et al. (2016) Kasen, D., Metzger, B. D., & Bildsten, L. 2016, ApJ, 821, 36, doi: 10.3847/0004-637X/821/1/36
  * Khatami & Kasen (2019) Khatami, D. K., & Kasen, D. N. 2019, ApJ, 878, 56, doi: 10.3847/1538-4357/ab1f09
  * Kleiser et al. (2018) Kleiser, I., Fuller, J., & Kasen, D. 2018, MNRAS, 481, L141, doi: 10.1093/mnrasl/sly180
  * Komissarov (2013) Komissarov, S. S. 2013, MNRAS, 428, 2459, doi: 10.1093/mnras/sts214
  * Kulkarni et al. (2021) Kulkarni, S. R., Harrison, F. A., Grefenstette, B. W., et al. 2021, arXiv e-prints, arXiv:2111.15608, doi: 10.48550/arXiv.2111.15608
  * Lazzati et al. (2017) Lazzati, D., Deich, A., Morsony, B. J., & Workman, J. C. 2017, MNRAS, 471, 1652, doi: 10.1093/mnras/stx1683
  * Leloudas et al. (2012) Leloudas, G., Chatzopoulos, E., Dilday, B., et al. 2012, A&A, 541, A129, doi: 10.1051/0004-6361/201118498
  * Levan et al. (2013) Levan, A. J., Read, A. M., Metzger, B. D., Wheatley, P. J., & Tanvir, N. R. 2013, ApJ, 771, 136, doi: 10.1088/0004-637X/771/2/136
  * Levan et al. (2014) Levan, A. J., Tanvir, N. R., Starling, R. L. C., et al. 2014, ApJ, 781, 13, doi: 10.1088/0004-637X/781/1/13
  * Liu & Modjaz (2016) Liu, Y.-Q., & Modjaz, M. 2016, ArXiv e-prints.  https://arxiv.org/abs/1612.07321
  * Lunnan et al. (2015) Lunnan, R., Chornock, R., Berger, E., et al. 2015, ApJ, 804, 90, doi: 10.1088/0004-637X/804/2/90
  * Lunnan et al. (2018) Lunnan, R., Fransson, C., Vreeswijk, P. M., et al. 2018, Nature Astronomy, 2, 887, doi: 10.1038/s41550-018-0568-z
  * Lyubarsky (2003) Lyubarsky, Y. E. 2003, MNRAS, 345, 153, doi: 10.1046/j.1365-8711.2003.06927.x
  * MacFadyen & Woosley (1999) MacFadyen, A. I., & Woosley, S. E. 1999, ApJ, 524, 262, doi: 10.1086/307790
  * Maeda et al. (2023) Maeda, K., Suzuki, A., & Izzo, L. 2023, MNRAS, 522, 2267, doi: 10.1093/mnras/stad1075
  * Margalit & Metzger (2018) Margalit, B., & Metzger, B. D. 2018, ApJ, 868, L4, doi: 10.3847/2041-8213/aaedad
  * Margalit et al. (2018) Margalit, B., Metzger, B. D., Thompson, T. A., Nicholl, M., & Sukhbold, T. 2018, MNRAS, 475, 2659, doi: 10.1093/mnras/sty013
  * Margutti et al. (2018) Margutti, R., et al. 2018, ArXiv e-prints.  https://arxiv.org/abs/1801.03531
  * Margutti et al. (2023) Margutti, R., Bright, J. S., Matthews, D. J., et al. 2023, ApJ, 954, L45, doi: 10.3847/2041-8213/acf1fd
  * Metzger et al. (2018) Metzger, B. D., Beniamini, P., & Giannios, D. 2018, ArXiv e-prints.  https://arxiv.org/abs/1802.07750
  * Metzger et al. (2011) Metzger, B. D., Giannios, D., Thompson, T. A., Bucciantini, N., & Quataert, E. 2011, MNRAS, 413, 2031, doi: 10.1111/j.1365-2966.2011.18280.x
  * Metzger et al. (2015) Metzger, B. D., Margalit, B., Kasen, D., & Quataert, E. 2015, MNRAS, 454, 3311, doi: 10.1093/mnras/stv2224
  * Mignone et al. (2007) Mignone, A., Bodo, G., Massaglia, S., et al. 2007, ApJS, 170, 228, doi: 10.1086/513316
  * Mondal et al. (2020) Mondal, S., Bera, A., Chandra, P., & Das, B. 2020, MNRAS, 498, 3863, doi: 10.1093/mnras/staa2637
  * Mooley et al. (2018) Mooley, K. P., Deller, A. T., Gottlieb, O., et al. 2018, Nature, 561, 355, doi: 10.1038/s41586-018-0486-3
  * Mösta et al. (2015) Mösta, P., Ott, C. D., Radice, D., et al. 2015, Nature, 528, 376, doi: 10.1038/nature15755
  * Murase et al. (2016) Murase, K., Guetta, D., & Ahlers, M. 2016, Phys. Rev. Lett., 116, 071101, doi: 10.1103/PhysRevLett.116.071101
  * Nakar (2015) Nakar, E. 2015, ApJ, 807, 172, doi: 10.1088/0004-637X/807/2/172
  * Nakar & Piran (2017) Nakar, E., & Piran, T. 2017, ApJ, 834, 28, doi: 10.3847/1538-4357/834/1/28
  * Nakar et al. (2002) Nakar, E., Piran, T., & Granot, J. 2002, ApJ, 579, 699, doi: 10.1086/342791
  * Nakar & Sari (2010) Nakar, E., & Sari, R. 2010, ApJ, 725, 904, doi: 10.1088/0004-637X/725/1/904
  * Nicholl (2021) Nicholl, M. 2021, Astronomy and Geophysics, 62, 5.34, doi: 10.1093/astrogeo/atab092
  * Nicholl et al. (2019) Nicholl, M., Berger, E., Blanchard, P. K., Gomez, S., & Chornock, R. 2019, ApJ, 871, 102, doi: 10.3847/1538-4357/aaf470
  * Nicholl et al. (2017) Nicholl, M., Guillochon, J., & Berger, E. 2017, ApJ, 850, 55, doi: 10.3847/1538-4357/aa9334
  * Nicholl & Smartt (2016) Nicholl, M., & Smartt, S. J. 2016, MNRAS, 457, L79, doi: 10.1093/mnrasl/slv210
  * Nicholl et al. (2015) Nicholl, M., et al. 2015, ArXiv e-prints.  https://arxiv.org/abs/1503.03310
  * Omand & Jerkstrand (2023) Omand, C. M. B., & Jerkstrand, A. 2023, A&A, 673, A107, doi: 10.1051/0004-6361/202245406
  * Omand et al. (2018) Omand, C. M. B., Kashiyama, K., & Murase, K. 2018, MNRAS, 474, 573, doi: 10.1093/mnras/stx2743
  * Oren et al. (2004) Oren, Y., Nakar, E., & Piran, T. 2004, MNRAS, 353, L35, doi: 10.1111/j.1365-2966.2004.08247.x
  * Panaitescu & Kumar (2002) Panaitescu, A., & Kumar, P. 2002, ApJ, 571, 779, doi: 10.1086/340094
  * Perley et al. (2024) Perley, D. A., Ho, A. Y. Q., Fausnaugh, M., et al. 2024, arXiv e-prints, arXiv:2401.16470, doi: 10.48550/arXiv.2401.16470
  * Porth et al. (2013) Porth, O., Komissarov, S. S., & Keppens, R. 2013, MNRAS, 431, L48, doi: 10.1093/mnrasl/slt006
  * Quataert & Kasen (2012) Quataert, E., & Kasen, D. 2012, MNRAS, 419, L1, doi: 10.1111/j.1745-3933.2011.01151.x
  * Quimby et al. (2013) Quimby, R. M., Yuan, F., Akerlof, C., & Wheeler, J. C. 2013, MNRAS, 431, 912, doi: 10.1093/mnras/stt213
  * Quimby et al. (2011) Quimby, R. M., et al. 2011, Nature, 474, 487, doi: 10.1038/nature10095
  * Sagiv et al. (2014) Sagiv, I., Gal-Yam, A., Ofek, E. O., et al. 2014, AJ, 147, 79, doi: 10.1088/0004-6256/147/4/79
  * Sari & Piran (1999) Sari, R., & Piran, T. 1999, ApJ, 520, 641, doi: 10.1086/307508
  * Smith et al. (2016) Smith, M., Sullivan, M., D’Andrea, C. B., et al. 2016, ApJ, 818, L8, doi: 10.3847/2041-8205/818/1/L8
  * Stanek et al. (2003) Stanek, K. Z., et al. 2003, ApJ, 591, L17, doi: 10.1086/376976
  * Stevance & Eldridge (2021) Stevance, H. F., & Eldridge, J. J. 2021, MNRAS, 504, L51, doi: 10.1093/mnrasl/slab039
  * Suzuki & Maeda (2016) Suzuki, A., & Maeda, K. 2016, ArXiv e-prints.  https://arxiv.org/abs/1612.03911
  * Thompson et al. (2004) Thompson, T. A., Chang, P., & Quataert, E. 2004, ApJ, 611, 380, doi: 10.1086/421969
  * Vurm & Metzger (2021) Vurm, I., & Metzger, B. D. 2021, ApJ, 917, 77, doi: 10.3847/1538-4357/ac0826
  * Werner et al. (2023) Werner, N., Řípa, J., Thöne, C., et al. 2023, arXiv e-prints, arXiv:2306.15080, doi: 10.48550/arXiv.2306.15080
  * Wheeler et al. (2000) Wheeler, J. C., Yi, I., Höflich, P., & Wang, L. 2000, ApJ, 537, 810, doi: 10.1086/309055
  * Woosley (1993) Woosley, S. E. 1993, ApJ, 405, 273, doi: 10.1086/172359
  * Woosley (2010) —. 2010, ApJ, 719, L204, doi: 10.1088/2041-8205/719/2/L204
  * Woosley & Bloom (2006) Woosley, S. E., & Bloom, J. S. 2006, ARAA, 44, 507, doi: 10.1146/annurev.astro.43.072103.150558
  * Woosley et al. (2002) Woosley, S. E., Heger, A., & Weaver, T. A. 2002, Reviews of Modern Physics, 74, 1015, doi: 10.1103/RevModPhys.74.1015
  * Yan et al. (2017) Yan, L., Quimby, R., Gal-Yam, A., et al. 2017, ApJ, 840, 57, doi: 10.3847/1538-4357/aa6b02
  * Yuan et al. (2022) Yuan, W., Zhang, C., Chen, Y., & Ling, Z. 2022, in Handbook of X-ray and Gamma-ray Astrophysics, 86, doi: 10.1007/978-981-16-4544-0_151-1
  * Zhao et al. (2021) Zhao, W.-C., Xue, X.-X., & Cao, X.-F. 2021, New A, 83, 101506, doi: 10.1016/j.newast.2020.101506


where $w_{j,\,N_{j}+1}=0$. Expanding $x$ by using this Jordan basis in a way
that

$x=\sum_{j=1}^{N}\sum_{k=1}^{N_{j}}x_{j,\,k}w_{j,\,k}\;.$ (7.1)

Here, by the Putzer spectral method (cf. [2, Proof of Lemma B.1.]), we can
write

$e^{\mathbb{A}t}x=\sum_{j=1}^{N}\sum_{k=1}^{N_{j}}\sum_{i=1}^{k}\frac{t^{k-i}e^{-\lambda_{j}t}}{(k-i)!}x_{j,\,i}w_{j,\,k}\;.$

Then, define

$\displaystyle\eta$
$\displaystyle=\eta(x):=\min\left\\{\textup{Re}(\lambda_{j}):x_{j,\,k}\neq
0\text{ for some }k\in\llbracket 1,\,N_{j}\rrbracket\right\\}>0\;,$
$\displaystyle\nu$ $\displaystyle=\nu(x):=\max\left\\{N_{j}-k:x_{j,\,k}\neq
0\text{ for some }k\in\llbracket 1,\,N_{j}\rrbracket\right\\}\geq
0\;\text{and}$ $\displaystyle I$
$\displaystyle=I(x):=\\{j\in[1,\,N]:\textrm{Re}(\lambda_{j})=\eta\text{ and
}x_{j,\,N_{j}-\nu}\neq 0\\}\;,$

so that we can write

$\frac{\mathrm{e}^{\eta t}}{t^{\nu}}e^{\mathbb{A}t}x=\sum_{j\in
I}\frac{e^{-(\lambda_{j}-\eta)t}}{\nu!}x_{j,\,N_{j}-\nu}w_{j,\,N_{j}}+\varepsilon_{t}\;,$
(7.2)

where $\varepsilon_{t}\rightarrow 0$ as $t\rightarrow\infty$. Since
$\lambda_{j}-\eta$ is pure imaginary for $j\in I$, by writing

$\lambda_{j}-\eta=-i\theta_{j}\;\;\;\text{and}\;\;\;v_{j}=\frac{x_{j,\,N_{j}-\nu}}{\nu!}w_{j,\,N_{j}}\;,$

we can deduce from (7.2) that

$\lim_{t\rightarrow\infty}\left|\frac{\mathrm{e}^{\eta
t}}{t^{\nu}}e^{\mathbb{A}t}x-\sum_{k\in
I}\mathrm{e}^{i\theta_{k}t}v_{k}\right|=0\;.$ (7.3)

In [2, Lemma B.2.], using the Hartman-Grobman theorem mentioned above, it is
proven that we can replace $e^{\mathbb{A}t}x$ at (7.3) with $X_{t}(x)$.
Therefore, we can derive (2.15) with $\tau(x)=0$ in the case
$x\in\mathcal{U}$.

For general $x\in\mathbb{R}^{2d}$, define $\tau(x)>0$ as

$\tau(x)=\inf\\{t>0:X_{t}(x)\in\mathcal{U}\\}+1$

so that $X_{\tau(x)}(x)\in\mathcal{U}$. Note that $\tau(x)<\infty$ follows
from Theorem 2.3. Since $X_{t+\tau(x)}(x)=X_{t}(X_{\tau(x)}(x))$, we also have
(2.15) with the constants defined at $X_{\tau(x)}(x)\in\mathcal{U}$ in the
previous step. ∎

######  Remark 7.1.

In some special cases, the constants obtained in the previous lemma can be
simplified.

  1. (1)

For $x\in\mathcal{U}$ with $x_{j,\,k}\neq 0$ for all $j,\,k$ (which happens
almost everywhere in $\mathcal{U}$ with respect to the Lebesgue measure),
$\eta$ is the smallest real part of eigenvalue, and $\nu+1$ is the largest
size of Jordan block associated with an eigenvalue of smallest real part.

  2. (2)

When $\mathbb{A}$ is diagonalizable, then $\nu=0$.

  3. (3)

When all eigenvalues of $\mathbb{A}$ are positive real, we have $\theta_{k}=0$
for all $k\in I$ and all vectors $v_{k}$ are real vectors.

###  7.2. Strategy to estimate total-variation distance

Let $x\in\mathbb{R}^{2d}$ and recall the constants and vectors defined in
Lemma 2.10. The mixing time for the process $(X_{t}^{\epsilon}(x))_{t\geq 0}$
is defined as

$t_{\mathrm{mix}}^{\epsilon}(x):=\frac{1}{2\eta}\log\left(\frac{1}{2\epsilon}\right)+\frac{\nu}{\eta}\log\log\left(\frac{1}{2\epsilon}\right)+\tau(x)\;.$

For $t\geq\tau(x)$ define $v(t,\,x)\in\mathbb{R}^{2d}$ and
$D^{\epsilon}(t,\,x)\in[0,\,1]$ as

$\displaystyle v(t,\,x)$
$\displaystyle:=\frac{(t-\tau(x))^{\nu}}{\mathrm{e}^{\eta(t-\tau(x))}}\Sigma^{-1/2}\sum_{k=1}^{m}\mathrm{e}^{i\theta_{k}(t-\tau(x))}v_{k}\;,$
(7.4) $\displaystyle D^{\epsilon}(t,\,x)$
$\displaystyle:=\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{1}{\sqrt{2\epsilon}}v(t,\,x),\,\mathbb{I}_{2d}\right),\,\mathcal{N}(0,\,\mathbb{I}_{2d})\right)\;.$
(7.5)

Our main result regarding the estimate of the total variation distance to the
equilibrium is the following theorem.

######  Theorem 7.2.

For all $x\in\mathbb{R}^{2d}$ and $w\in\mathbb{R}$, it holds that

$\lim_{\epsilon\rightarrow
0}\left|\mathrm{d}_{\mathrm{TV}}\left(X_{t_{\mathrm{mix}}^{\epsilon}(x)+w}^{\epsilon}(x),\,\mu^{\epsilon}\right)-D^{\epsilon}(t_{\textup{mix}}^{\epsilon}(x)+w,\,x)\right|=0\;.$

We note here that Theorems 2.12 and 2.13 are a simple consequence of this
theorem. It remains to analyze the behavior of
$D^{\epsilon}(t_{\textup{mix}}^{\epsilon}(x)+w)$ for $w\in\mathbb{R}$ as
$\epsilon\rightarrow 0$, and this step is identical to [2] as for the
overdamped Langevin dynamics. Hence, we formally conclude the proof of these
theorems by assuming Theorem 7.2.

###### Proof of Theorems 2.12 and 2.13 .

The proof of Theorems 2.12 and 2.13 are identical to that of [2, Corollary
2.7] and [2, Corollary 2.9], respectively. ∎

The basic strategy to prove Theorem 7.2 is based on the coupling argument
developed in [2], but due to the degeneracy and lack of global contraction, we
need to considerably refine the argument. Our scheme of proof allows us to
tackle the degeneracy in (1.2) and could be extended to other processes with
similar degeneracy. We explain now the structure of this proof.

The first step is to control the total variation distance between the process
$(X_{t}^{\epsilon}(x))_{t\geq 0}$ and its Gaussian approximation
$(Z_{t}^{\epsilon}(x))_{t\geq 0}$ using the moment estimates obtained in
Proposition 5.3.

######  Proposition 7.3.

Let $(t^{\epsilon})_{\epsilon>0}$ be a sequence of positive real numbers such
that

$\alpha_{0}\log\frac{1}{\epsilon}\leq
t^{\epsilon}\leq\frac{1}{\epsilon^{\theta}}$ (7.6)

for some constant $\alpha_{0}>0$, where $\theta$ is the constant defined in
Notation 5.2. Then, for any $r>0$, we have

$\limsup_{\epsilon\rightarrow
0}\sup_{x\in\mathcal{D}_{r}}\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)=0\;.$
(7.7)

This proposition is proven in Section 7.3. The next step is to estimate the
distribution of the Gaussian process $(Z_{t}^{\epsilon}(x))_{t\geq 0}$.

######  Proposition 7.4.

Let $x\in\mathbb{R}^{2d}$ and let $(t^{\epsilon})_{\epsilon>0}$ be a sequence
of positive reals such that

$t^{\epsilon}\geq t_{\mathrm{mix}}^{\epsilon}(x)-c$ (7.8)

for some constant $c>0$ for all small enough $\epsilon$. Then, we have

$\lim_{\epsilon\rightarrow
0}\left|\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\,\mathcal{N}(0,2\epsilon\Sigma)\right)-D^{\epsilon}(t^{\epsilon},\,x)\right|=0\;.$

The remaining part of this article is devoted to the proof of Propositions
7.3, 7.4, and Theorems 2.8 and 7.2. Let us start with the proof of Theorem 7.2
assuming that Propositions 7.3, 7.4, and Theorem 2.8 hold.

###### Proof of Theorem 7.2.

By the triangle inequality, we have

$\displaystyle\bigg{|}\,\big{|}\mathrm{d}_{\mathrm{TV}}\left(X_{t_{\mathrm{mix}}^{\epsilon}(x)+w}^{\epsilon}(x),\,\mu^{\epsilon}\right)-D^{\epsilon}(t_{\textup{mix}}^{\epsilon}(x)+w,\,x)\big{|}$
$\displaystyle\;\;\;-\big{|}\mathrm{d}_{\mathrm{TV}}\left(Z_{t_{\mathrm{mix}}^{\epsilon}(x)+w}^{\epsilon}(x),\,\mathcal{N}(0,2\epsilon\Sigma)\right)-D^{\epsilon}(t_{\textup{mix}}^{\epsilon}(x)+w,\,x)\big{|}\bigg{|}$
$\displaystyle\leq\mathrm{d}_{\mathrm{TV}}\left(X_{t_{\mathrm{mix}}^{\epsilon}(x)+w}^{\epsilon}(x),\,Z_{t_{\mathrm{mix}}^{\epsilon}(x)+w}^{\epsilon}(x)\right)+\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(0,2\epsilon\Sigma),\,\mu^{\epsilon}\right)\;.$

Therefore, the proof is completed by Propositions 7.3 and 7.4 with
$t^{\epsilon}=t_{\mathrm{mix}}^{\epsilon}(x)+w$ (which satisfies (7.6) and
(7.8)) and Theorem 2.8. ∎

###  7.3. Proof of Proposition 7.3

We start by proving the following lemma which provides an upper-bound on the
total variation distance between the random variables $X_{t}^{\epsilon}(x)$
and $Z_{t}^{\epsilon}(x)$ for small times $t>0$

######  Lemma 7.5.

For all $x\in\mathbb{R}^{2d}$, we have that

$\mathrm{d}_{\mathrm{TV}}\left(X_{t}^{\epsilon}(x),Z_{t}^{\epsilon}(x)\right)^{2}\leq\frac{1}{2\epsilon}\int_{0}^{t}\mathbb{E}\left[\left|F(q_{s}^{\epsilon}(x))-F(q_{s}(x))-\textup{D}F(q_{s}(x))\left(q_{s}^{\epsilon}(x)-q_{s}(x)\right)\right|^{2}\right]\textup{d}s\;.$

###### Proof.

Write
$Z_{t}^{\epsilon}=(\widetilde{q}_{t}^{\epsilon},\,\widetilde{p}_{t}^{\epsilon})$
so that we can rewrite the SDE describing the process
$(Z_{t}^{\epsilon})_{t\geq 0}$ as

$\begin{cases}\textup{d}\widetilde{q}_{t}^{\epsilon}=\widetilde{p}_{t}^{\epsilon}\textup{d}t\\\
\textup{d}\widetilde{p}_{t}^{\epsilon}=-(F(q_{t})+\textup{D}F(q_{t})(\widetilde{q}_{t}^{\epsilon}-q_{t})+\gamma\widetilde{p}_{t}^{\epsilon})\textup{d}t+\sqrt{2\epsilon}\,\textup{d}B_{t}\;.\end{cases}$
(7.9)

Now, denote by $\mathbb{P}_{X}$ and $\mathbb{P}_{Z}$ the law of the processes
$(X_{s}^{\epsilon}(x))_{s\in[0,\,t]}$ and
$(Z_{s}^{\epsilon}(x))_{s\in[0,\,t]}$, respectively. Then, by Pinsker’s
inequality, we have

$\mathrm{d}_{\mathrm{TV}}\left(X_{t}^{\epsilon}(x),\,Z_{t}^{\epsilon}(x)\right)^{2}\leq-2\mathbb{E}_{\mathbb{P}_{X}}\left[\log\frac{\textup{d}\mathbb{P}_{Z}}{\textup{d}\mathbb{P}_{X}}\right]\;.$
(7.10)

By the Girsanov theorem, (1.2), and (7.9), we have

$\displaystyle\log\frac{\textup{d}\mathbb{P}_{Z}}{\textup{d}\mathbb{P}_{X}}=\;$
$\displaystyle\frac{1}{\sqrt{2\epsilon}}\int_{0}^{t}\left\langle
F(q_{s}^{\epsilon}(x))-F(q_{s}(x))-\textup{D}F(q_{s}(x))\left(q_{s}^{\epsilon}(x)-q_{s}(x)\right),\,\textup{d}B_{s}\right\rangle$
$\displaystyle-\frac{1}{4\epsilon}\int_{0}^{t}\left|F(q_{s}^{\epsilon}(x))-F(q_{s}(x))-\textup{D}F(q_{s}(x))\left(q_{s}^{\epsilon}(x)-q_{s}(x)\right)\right|^{2}\textup{d}s\;.$

Inserting this to (7.10), we can complete the proof. ∎

Our strategy is to first couple the processes $(X_{t}^{\epsilon}(x))_{t\geq
0}$ and $(Z_{t}^{\epsilon}(x))_{t\geq 0}$ for small times as explained in
Proposition 7.6, similarly to what was done in [2]. To that end, we need to
introduce the following small time scales. Recall the constants
$\alpha_{1},\,\alpha_{2}>0$ from Proposition 5.3 and the constant
$\alpha_{0}>0$ from the statement of Proposition 7.3. Let us take $\sigma$
small enough so that

$0<\sigma<\min\left(\frac{2\alpha_{1}\alpha_{0}}{3},\,\frac{2\alpha_{2}}{3}\right)\;.$
(7.11)

Additionally, let us take $\nu$ small enough from Proposition 5.3 such that

$0<\nu<\min\left(\frac{\sigma}{2},\alpha_{1}\alpha_{0}-\frac{3\sigma}{2}\right).$
(7.12)

Let us now define the time scale $\delta^{\epsilon}$ by
$\delta^{\epsilon}:=\epsilon^{\sigma}$. The following is the first step in the
approximation of $X_{t^{\epsilon}}^{\epsilon}(x)$ by
$Z_{t^{\epsilon}}^{\epsilon}(x)$ in the total variation distance.

######  Proposition 7.6.

Let $(t^{\epsilon})_{\epsilon>0}$ be a sequence of positive reals satisfying
(7.6). Then, for all $r>0$, we have

$\limsup_{\epsilon\rightarrow
0}\sup_{x\in\mathcal{D}_{r}}\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))\right)=0.$

###### Proof.

Let $r_{0}>0$ be small enough so that we have for all
$x\in\mathcal{D}_{r_{0}}$,

$T(x)=0\;\;\;\text{and\;\;\;}H(x)<\frac{1}{16(d+2)\kappa_{0}\rho}\;,$ (7.13)

where $T(x)$ is defined in Notation 5.2, and the constants $\kappa_{0}$ and
$\rho$ are the ones appearing in Lemma 4.2 and Assumption 2.11, respectively.

Since
$X_{t^{\epsilon}}^{\epsilon}(x)=X_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))$,
we have

$\displaystyle\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))\right)$
$\displaystyle\leq\int_{\mathbb{R}^{d}}\mathrm{d}_{\mathrm{TV}}\left(X_{\delta^{\epsilon}}^{\epsilon}(u),Z_{\delta^{\epsilon}}^{\epsilon}(u)\right)\mathbb{P}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u)$
$\displaystyle\leq\int_{\mathcal{D}_{r_{0}}}\mathrm{d}_{\mathrm{TV}}\left(X_{\delta^{\epsilon}}^{\epsilon}(u),Z_{\delta^{\epsilon}}^{\epsilon}(u)\right)\mathbb{P}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u)+\mathbb{P}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r_{0})\;.$
(7.14)

By Chebychev inequality and Proposition 4.9, we have

$\displaystyle\mathbb{P}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r_{0})$
$\displaystyle\leq\frac{\mathbb{E}\left[|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2}\right]}{r_{0}^{2}}\leq\frac{\kappa_{0}}{r_{0}^{2}}\left[H(x)\mathrm{e}^{-\lambda(t^{\epsilon}-\delta^{\epsilon})}+\frac{d\epsilon}{\lambda}\right]\;.$
(7.15)

Next we fix $u\in\mathcal{D}_{r_{0}}$. Then, by Lemma 7.5, we have

$\mathrm{d}_{\mathrm{TV}}\left(X_{\delta^{\epsilon}}^{\epsilon}(u),\,Z_{\delta^{\epsilon}}^{\epsilon}(u)\right)^{2}\leq\frac{1}{\epsilon}\left[A^{\epsilon}(1)+A^{\epsilon}(2)\right]\;,$
(7.16)

where

$\displaystyle
A^{\epsilon}(1)=\int_{0}^{\delta^{\epsilon}}\mathbb{E}\left[\left|F(q_{s}^{\epsilon}(u))-F(q_{s}(u))\right|^{2}\right]ds\;,$
$\displaystyle
A^{\epsilon}(2)=\left\|\textup{D}F(q_{s}(u))\right\|\int_{0}^{\delta^{\epsilon}}\mathbb{E}\left[|X_{s}^{\epsilon}(u)-X_{s}(u)|^{2}\right]ds\;.$

Let us first look $A^{\epsilon}(1)$. By the Cauchy-Schwarz inequality and
Assumption 2.11, we have

$\displaystyle\big{|}F(q_{s}^{\epsilon}(u))-F(q_{s}(u))\big{|}^{2}$
$\displaystyle=\big{|}\int_{0}^{1}\textup{D}F(q_{s}(u)+\zeta(q_{s}^{\epsilon}(u)-q_{s}(u)))\,(q_{s}^{\epsilon}(u)-q_{s}(u))\mathrm{d}\zeta\bigg{|}^{2}$
$\displaystyle\leq\big{|}q_{s}^{\epsilon}(u)-q_{s}(u)\big{|}^{2}\int_{0}^{1}|\textup{D}F(q_{s}(u)+\zeta(q_{s}^{\epsilon}(u)-q_{s}(u)))|^{2}\mathrm{d}\zeta$
$\displaystyle\leq
C\,\big{|}q_{s}^{\epsilon}(u)-q_{s}(u)\big{|}^{2}\int_{0}^{1}\mathrm{e}^{2\rho|q_{s}(u)+\zeta(q_{s}^{\epsilon}(u)-q_{s}(u))|^{2}}\mathrm{d}\zeta$
$\displaystyle\leq
C\big{|}q_{s}^{\epsilon}(u)-q_{s}(u)\big{|}^{2}\mathrm{e}^{4\rho|q_{s}(u)|^{2}}\mathrm{e}^{4\rho|q_{s}^{\epsilon}(u)|^{2}}\;.$

By Theorem 2.3, we have for $u=(q,p)$
$|q_{s}(u)|^{2}\leq\kappa\left(|u|^{2}+U(q)\right)$ for all $s\geq 0$.
Inserting this bound in the exponential above and noticing that

$\big{|}q_{s}^{\epsilon}(u)-q_{s}(u)\big{|}^{2}\leq\big{|}X_{s}^{\epsilon}(u)-X_{s}(u)\big{|}^{2}\leq
2\left[\big{|}X_{s}^{\epsilon}(u)-Z_{s}^{\epsilon}(u)\big{|}^{2}+2\epsilon\big{|}Y_{s}(u)\big{|}^{2}\right]\;,$
(7.17)

and then applying Cauchy-Schwarz inequality, we get

$A^{\epsilon}(1)\leq
C\mathrm{e}^{4\rho\kappa\left(|u|^{2}+U(q)\right))}\int_{0}^{\delta^{\epsilon}}\mathbb{E}\left[\mathrm{e}^{8\rho|X_{s}^{\epsilon}(u)|^{2}}\right]^{1/2}\left(\mathbb{E}\left[\big{|}X_{s}^{\epsilon}(u)-Z_{s}^{\epsilon}(u)\big{|}^{4}\right]^{1/2}+2\epsilon\mathbb{E}\left[\big{|}Y_{s}\big{|}^{4}\right]^{1/2}\right)\mathrm{d}s\;.$

Since $u\in\mathcal{D}_{r_{0}}$, by (7.13), for all small enough $\epsilon>0$,

$8\rho<\frac{1}{2(d+2)\kappa_{0}\left(H(u)+\frac{d\epsilon}{\lambda}\right)}\leq\frac{1}{2(d+2)\kappa_{0}\left(H(u)\mathrm{e}^{-\lambda
s}+\frac{d\epsilon}{\lambda}\right)}\;,$

and therefore by Corollary 4.11, Proposition 5.3 and Lemma 5.7, we get

$A^{\epsilon}(1)\leq C(n,\,r_{0})\epsilon^{1-2\nu}\delta^{\epsilon}\;.$ (7.18)

On the other hand, by Theorem 2.3, the second inequality of (7.17),
Proposition 5.3 and Lemma 5.7, we also get

$A^{\epsilon}(2)\leq C(n,\,r_{0})\epsilon^{1-2\nu}\delta^{\epsilon}\;.$ (7.19)

Inserting (7.15), (7.16), (7.18) and (7.19) to (7.14) yields

$\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))\right)\leq
C(n,\,r_{0})(\delta^{\epsilon}\epsilon^{-2\nu})^{1/2}+\frac{\kappa_{0}}{r_{0}^{2}}\left[H(x)\mathrm{e}^{-\lambda(t^{\epsilon}-\delta^{\epsilon})}+\frac{d\epsilon}{\lambda}\right]\;.$

This completes the proof since
$\delta^{\epsilon}\epsilon^{-2\nu}=\epsilon^{\sigma-2\nu}$ and $\sigma>2\nu$
by (7.12). Besides, $\mathrm{e}^{-\lambda(t^{\epsilon}-\delta^{\epsilon})}$
converge to $0$ as $\epsilon\rightarrow 0$. ∎

Next we control the total variation distance between
$Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))$
and $Z_{t^{\epsilon}}^{\epsilon}(x)$.

######  Proposition 7.7.

Let $(t^{\epsilon})_{\epsilon>0}$ be a sequence of positive reals satisfying
(7.6). Then, for all $r>0$, it holds that

$\limsup_{\epsilon\rightarrow
0}\sup_{x\in\mathcal{D}_{r}}\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)=0\;.$

###### Proof.

Since
$Z_{t^{\epsilon}}^{\epsilon}(x)=Z_{\delta^{\epsilon}}^{\epsilon}(Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))$,
we have

$\displaystyle\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)$
$\displaystyle\leq\int_{\mathbb{R}^{d}}\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(u),\,Z_{\delta^{\epsilon}}^{\epsilon}(u^{\prime})\right)\mathbb{P}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u,\,Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u^{\prime})$
$\displaystyle\leq\int_{u,\,u^{\prime}\in\mathcal{D}_{r}}\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(u),\,Z_{\delta^{\epsilon}}^{\epsilon}(u^{\prime})\right)\mathbb{P}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u,\,Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u^{\prime})$
$\displaystyle\qquad+\mathbb{P}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)+\mathbb{P}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)\;.$
(7.20)

By the Markov inquality and the Cauchy-Schwarz inequality, we have

$\displaystyle\mathbb{P}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)$
$\displaystyle\leq\frac{\mathbb{E}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2})}{r^{2}}\leq\frac{2}{r^{2}}\left[\mathbb{E}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)-Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2})+\mathbb{E}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2})\right]\;,$
$\displaystyle\mathbb{P}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)$
$\displaystyle\leq\frac{1}{r^{2}}\mathbb{E}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2})\;.$

Since by the Cauchy-Schwarz inequality we have

$\mathbb{E}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|^{2})\leq
2|X_{t^{\epsilon}-\delta^{\epsilon}}(x)|^{2}+4\epsilon\mathbb{E}(|Y_{t^{\epsilon}-\delta^{\epsilon}}(x)|^{2})\;,$

we get from Theorem 2.3, Proposition 5.3 and Lemma 5.7 that

$\mathbb{P}(|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)+\mathbb{P}(|Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|>r)\leq\frac{C(r)}{r^{2}}\left[\epsilon+\textup{e}^{-\lambda(t^{\epsilon}-\delta^{\epsilon})}\right]\;.$
(7.21)

Next we compute the integral at the right-hand side of (7.20). By Lemma 6.8
and Proposition 5.3, this integral is bounded by

$\displaystyle\int_{u,\,u^{\prime}\in\mathcal{D}_{r}}\left[\frac{C(r)}{(\delta^{\epsilon})^{3/2}\epsilon^{1/2}}|u-u^{\prime}|+\delta^{\epsilon}\right]\mathbb{P}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u,\,Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)\in\mathrm{d}u^{\prime})$
$\displaystyle\leq\frac{C(r)}{(\delta^{\epsilon})^{3/2}\epsilon^{1/2}}\mathbb{E}\left[|X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)-Z_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)|\right]+\delta^{\epsilon}\leq\frac{C(r)}{(\delta^{\epsilon})^{3/2}}\left[\frac{\textup{e}^{-\alpha_{1}(t^{\epsilon}-\delta^{\epsilon})}}{\epsilon^{\nu}}+\epsilon^{\alpha_{2}}\right]+\delta^{\epsilon}\;.$

By this computation, (7.20), and (7.21), we can conclude that

$\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)\leq
C(r)\left[\epsilon+e^{-\lambda(t^{\epsilon}-\delta^{\epsilon})}+\frac{1}{(\delta^{\epsilon})^{3/2}}\left(\frac{\textup{e}^{-\alpha_{1}(t^{\epsilon}-\delta^{\epsilon})}}{\epsilon^{\nu}}+\epsilon^{\alpha_{2}}\right)\right]+\delta^{\epsilon}\;.$

The right-hand side tends to $0$ by (7.6), (7.11) and (7.12). ∎

Now we are ready to prove Proposition 7.3.

###### Proof of Proposition 7.3.

By the triangle inequality, we have

$\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon},\,Z_{t^{\epsilon}}^{\epsilon}\right)\leq\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon},\,Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x))\right)+\mathrm{d}_{\mathrm{TV}}\left(Z_{\delta^{\epsilon}}^{\epsilon}(X_{t^{\epsilon}-\delta^{\epsilon}}^{\epsilon}(x)),\,Z_{t^{\epsilon}}^{\epsilon}\right)$

and therefore the proof is completed by Propositions 7.6 and 7.7. ∎

###  7.4. Proof of Proposition 7.4

###### Proof of Proposition 7.4.

By (6.20) and Lemma 6.5 (1), we can write

$\displaystyle\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\,\mathcal{N}(0,2\epsilon\Sigma)\right)$
$\displaystyle=\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(X_{t^{\epsilon}}(x),\,2\epsilon\Sigma_{t^{\epsilon}}(x)),\,\mathcal{N}(0,2\epsilon\Sigma)\right)$
$\displaystyle=\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma_{t^{\epsilon}}(x)\right),\,\mathcal{N}(0,\Sigma)\right)\;.$

Therefore, by the triangle inequality and Lemma 6.5 (2),

$\displaystyle\left|\,\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\mathcal{N}(0,2\epsilon\Sigma)\right)-\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma\right),\,\mathcal{N}(0,\Sigma)\right)\right|$
$\displaystyle\leq\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma_{t^{\epsilon}}(x)\right),\,\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma\right)\right)=\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(0,\Sigma_{t^{\epsilon}}(x)),\,\mathcal{N}(0,\Sigma)\right)\;.$

By Lemma 6.2 and the fact that $t^{\epsilon}\rightarrow\infty$ (by (7.4)), the
right-hand side converges to $0$, and hence, it suffices to prove that

$\limsup_{\epsilon\rightarrow
0}\left|\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma\right),\,\mathcal{N}(0,\Sigma)\right)-D^{\epsilon}(t^{\epsilon},\,x)\right|=0\;.$
(7.22)

By Lemma 6.5 (3),

$\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\Sigma\right),\,\mathcal{N}(0,\Sigma)\right)=\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\Sigma^{-1/2}\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\mathbb{I}_{2d}\right),\,\mathcal{N}(0,\,\mathbb{I}_{2d})\right)\;.$

Hence, recalling the definitions (7.4) of $v(t,\,x)$ and (7.5) of
$D^{\epsilon}(t,\,x)$, and applying the triangle inequality along with Lemmas
6.5 (2) and 6.6, we can bound the absolute value at (7.22) from above by

$\displaystyle\left|\,\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\Sigma^{-1/2}\frac{X_{t^{\epsilon}}(x)}{\sqrt{2\epsilon}},\,\mathbb{I}_{2d}\right),\,\mathcal{N}(0,\,\mathbb{I}_{2d})\right)-\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{v(t^{\epsilon},\,x)}{\sqrt{2\epsilon}},\,\mathbb{I}_{2d}\right),\,\mathcal{N}(0,\,\mathbb{I}_{2d})\right)\right|$
$\displaystyle\leq\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}\left(\frac{1}{\sqrt{2\epsilon}}\left(\Sigma^{-1/2}X_{t^{\epsilon}}(x)-v(t^{\epsilon},\,x)\right),\,\mathbb{I}_{2d}\right),\,\mathcal{N}\left(0,\,\mathbb{I}_{2d}\right)\right)$
$\displaystyle\leq\frac{1}{\sqrt{4\pi\epsilon}}\left|\Sigma^{-1/2}X_{t^{\epsilon}}(x)-v(t^{\epsilon},\,x)\right|$
$\displaystyle=\frac{1}{\sqrt{4\pi\epsilon}}\left|\frac{(t^{\epsilon}-\tau(x))^{\nu}}{\mathrm{e}^{\eta(t^{\epsilon}-\tau(x))}}\Sigma^{-1/2}\left(\frac{\mathrm{e}^{\eta(t^{\epsilon}-\tau(x))}}{(t^{\epsilon}-\tau(x))^{\nu}}X_{t^{\epsilon}}-\sum_{k=1}^{m}\mathrm{e}^{i\theta_{k}(t^{\epsilon}-\tau(x))}v_{k}\right)\right|\;.$

The last line converges to $0$ as $\epsilon\rightarrow 0$ by Lemma 2.10 and
since by (7.8) we have

$\limsup_{\epsilon\rightarrow\infty}\frac{(t^{\epsilon})^{\nu}}{\epsilon^{1/2}\mathrm{e}^{\eta
t^{\epsilon}}}<\infty\;.$

∎

##  8\. Proof of Theorem 2.8

It only remains to prove Theorem 2.8 to conclude this work. The proof of this
theorem is also based on Propositions 7.3 and 7.4. In addition to these
propositions, we shall need the following a priori variance bound on the
stationary density $\mu_{\epsilon}$.

######  Lemma 8.1.

There exists a constant $c>0$ independent of $\epsilon$ such that

$\int_{\mathbb{R}^{2d}}|x|^{2}\mu_{\epsilon}(\mathrm{d}x)\leq c\epsilon\;.$

###### Proof.

This proof is similar to the proof provided in [2] and is inspired from [29,
p.122, Section 5]. Let $L>0$. Then, for all $t\geq 0$, by the stationarity of
$\mu_{\epsilon}$ and Jensen’s inequality we get

$\int_{\mathbb{R}^{2d}}(|x|^{2}\land
L)\mu_{\epsilon}(\mathrm{d}x)=\int_{\mathbb{R}^{2d}}\mathbb{E}\left[|X_{t}^{\epsilon}(x)|^{2}\land
L\right]\mu_{\epsilon}(\mathrm{d}x)\leq\int_{\mathbb{R}^{2d}}\left(\mathbb{E}\left[|X_{t}^{\epsilon}(x)|^{2}\right]\land
L\right)\mu_{\epsilon}(\mathrm{d}x)\;.$ (8.1)

We bound the last expectation using Proposition 4.9. Then letting
$t\rightarrow\infty$ and applying the dominated convergence theorem, we get

$\int_{\mathbb{R}^{2d}}(|x|^{2}\land
L)\mu_{\epsilon}(\mathrm{d}x)\leq\int_{\mathbb{R}^{2d}}\left(c_{3}\epsilon\land
L\right)\mu_{\epsilon}(\mathrm{d}x)\;.$

Finally, letting $L\rightarrow\infty$ and applying the monotone convergence
theorem completes the proof. ∎

The main ingredients in the proof of Theorem 2.8 are Propositions 7.3, 7.4,
and the following lemma providing a rough a priori bound
$\frac{1}{\epsilon^{\theta}}$ on the mixing time.

######  Lemma 8.2.

Let $t_{\epsilon}=\frac{1}{\epsilon^{\theta}}$ where $\theta$ is the constant
defined in Notation 5.2. Then, for all $x\in\mathbb{R}^{d}$, we have that

$\lim_{\epsilon\rightarrow
0}\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,\mu^{\epsilon}\right)=0\;.$

###### Proof.

We fix $x\in\mathbb{R}^{2d}$. Denote by $p_{t}^{\epsilon}(\cdot,\,\cdot)$ the
transition kernel of the process $(X_{t}^{\epsilon})_{t\geq 0}$. Then, we have

$\displaystyle\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,\mu^{\epsilon}\right)$
$\displaystyle=\frac{1}{2}\int_{\mathbb{R}^{2d}}\left|p_{t^{\epsilon}}^{\epsilon}(x,\,z)-\mu^{\epsilon}(z)\right|\textup{d}z$
$\displaystyle=\frac{1}{2}\int_{\mathbb{R}^{2d}}\left|\int_{\mathbb{R}^{2d}}\left(p_{t^{\epsilon}}^{\epsilon}(x,\,z)-p_{t^{\epsilon}}^{\epsilon}(y,\,z)\right)\mu^{\epsilon}(y)\,\textup{d}y\right|\textup{d}z$
$\displaystyle\leq\frac{1}{2}\int_{\mathbb{R}^{2d}}\int_{\mathbb{R}^{2d}}\left|p_{t^{\epsilon}}^{\epsilon}(x,\,z)-p_{t^{\epsilon}}^{\epsilon}(y,\,z)\right|\mu^{\epsilon}(y)\,\textup{d}y\textup{d}z$
$\displaystyle=\int_{\mathbb{R}^{2d}}\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,X_{t^{\epsilon}}^{\epsilon}(y)\right)\mu^{\epsilon}(\mathrm{d}y)\;.$

Take $R$ large enough so that $x\in\mathcal{D}_{R}$. By Lemma 8.1 and the
Markov inequality, there exists $c>0$ such that

$\mu^{\epsilon}(\mathcal{D}_{R}^{c})\leq\frac{\int_{\mathbb{R}^{2d}}|x|^{2}\mu_{\epsilon}(\mathrm{d}x)}{R^{2}}\leq\frac{c\epsilon}{R^{2}}$

and therefore we get

$\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,\mu^{\epsilon}\right)\leq\int_{\mathcal{D}_{R}}\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,X_{t^{\epsilon}}^{\epsilon}(y)\right)\mu^{\epsilon}(\mathrm{d}y)+\frac{c\epsilon}{R^{2}}\;.$

By Proposition 7.3,

$\limsup_{\epsilon\rightarrow
0}\int_{\mathcal{D}_{R}}\left|\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,X_{t^{\epsilon}}^{\epsilon}(y)\right)-\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\,Z_{t^{\epsilon}}^{\epsilon}(y)\right)\right|\mu^{\epsilon}(\mathrm{d}y)=0$

since both $x$ and $y$ at the integrand belong to the set $\mathcal{D}_{R}$.
On the other hand, by Lemma 6.7, we have that (note that
$x\in\mathcal{D}_{R}$)

$\lim_{\epsilon\rightarrow
0}\int_{\mathcal{D}_{R}}\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\,Z_{t^{\epsilon}}^{\epsilon}(y)\right)\mu^{\epsilon}(\mathrm{d}y)=0$

and the proof is completed. ∎

Now we prove Theorem 2.8.

###### Proof of Theorem 2.8.

Let $t^{\epsilon}=\frac{1}{\epsilon^{\theta}}$ where $\theta$ is the constant
defined in Notation 5.2. Then, by the triangle inequality, for any
$x\in\mathbb{R}^{2d}$

$\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(0,2\epsilon\Sigma),\,\mu^{\epsilon}\right)\leq\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(0,2\epsilon\Sigma),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)+\mathrm{d}_{\mathrm{TV}}\left(Z_{t^{\epsilon}}^{\epsilon}(x),\,X_{t^{\epsilon}}^{\epsilon}(x)\right)+\mathrm{d}_{\mathrm{TV}}\left(X_{t^{\epsilon}}^{\epsilon}(x),\,\mu^{\epsilon}\right)\;.$

By Proposition 7.3 and Lemma 8.2, the second and the third term at the right-
hand side converge to $0$ as $\epsilon\rightarrow 0$. On the other hand, by
the triangle inequality, the first term at the right-hand side is bounded from
above by

$\left|\mathrm{d}_{\mathrm{TV}}\left(\mathcal{N}(0,2\epsilon\Sigma),\,Z_{t^{\epsilon}}^{\epsilon}(x)\right)-D^{\epsilon}(t^{\epsilon})\right|+D^{\epsilon}(t^{\epsilon})\;.$

By Proposition 7.3, the first term tends to $0$ as $\epsilon\rightarrow 0$,
while it is easy to check that $D^{\epsilon}(t^{\epsilon})\rightarrow 0$ for
$t^{\epsilon}=\frac{1}{\epsilon^{\theta}}$ as $\epsilon\rightarrow 0$. This
completes the proof. ∎

Acknowledgement. M. R. and S. L. were supported by the Samsung Science and
Technology Foundation (Project Number SSTF-BA1901-03). I. S. was supported by
the Samsung Science and Technology Foundation (Project Number SSTF-BA1901-03)
and National Research Foundation(NRF) of Korea grant funded by the Korea
government(MSIT) (No. 2022R1F1A106366811, 2022R1A5A600084012, and
2023R1A2C100517311). The authors are also supported by BK21 SNU Mathematical
Sciences Division.

## References

  * [1] G. Barrera.  Limit behavior of the invariant measure for Langevin dynamics.  Probability and Mathematical Statistics, 42(1):143–162, 2022.
  * [2] G. Barrera and M. Jara.  Thermalisation for small random perturbations of dynamical systems.  The Annals of Applied Probability, 30(3):1164 – 1208, 2020.
  * [3] Paul Cuff, Jian Ding, Oren Louidor, Eyal Lubetzky, Yuval Peres, and Allan Sly.  Glauber dynamics for the mean-field Potts model.  Journal of Statistical Physics, 149(3):432–477, 2012.
  * [4] L. Devroye, A. Mehrabian, and T. Reddad.  The total variation distance between high-dimensional Gaussians with the same mean.  arXiv:1810.08693, 2018.
  * [5] Dor Elboim and Dominik Schmid.  Mixing times and cutoff for the tasep in the high and low density phase.  arXiv preprint arXiv:2208.08306, 2022.
  * [6] M. I. Freidlin and A. D. Wentzel.  Random perturbations of dynamical systems.  Grundlehren der mathematischen Wissenschaften. Springer, 1984.
  * [7] Shirshendu Ganguly, Eyal Lubetzky, and Fabio Martinelli.  Cutoff for the east process.  Communications in mathematical physics, 335:1287–1322, 2015.
  * [8] Shirshendu Ganguly and Insuk Seo.  Information percolation and cutoff for the random-cluster model.  Random Structures & Algorithms, 57(3):770–822, 2020.
  * [9] Hubert Lacoin.  The cutoff profile for the simple exclusion process on the circle.  2016\.
  * [10] P. Lancaster and M. Tismenetsky.  The theory of matrices second edition, with applications.’1 academic press.  Inc.. San Diego, 1985.
  * [11] J. Lee and I. Seo.  Non-reversible metastable diffusions with Gibbs invariant measure i: Eyring-Kramers formula.  Probability Theory and Related Fields, 182:849–903, 2022.
  * [12] J. Lee and I. Seo.  Non-reversible metastable diffusions with Gibbs invariant measure ii: Markov chain convergence.  Journal of Statistical Physics, 189(25):1–34, 2023.
  * [13] S. Lee, M. Ramil, and I. Seo.  Eyring-Kramers formula for underdamped Langevin dynamics.  In preparation, 2023.
  * [14] T. Lelièvre, M. Ramil, and J. Reygner.  A probabilistic study of the kinetic Fokker–Planck equation in cylindrical domains.  Journal of Evolution Equations, 22(2):38, 2022.
  * [15] David A Levin, Malwina J Luczak, and Yuval Peres.  Glauber dynamics for the mean-field Ising model: cut-off, critical power law, and metastability.  Probability Theory and Related Fields, 146:223–265, 2010.
  * [16] David Asher Levin, Yuval Peres, and Elizabeth Lee Wilmer.  Markov chains and mixing times.  American Mathematical Soc., 2009.
  * [17] Eyal Lubetzky and Yuval Peres.  Cutoff on all Ramanujan graphs.  Geometric and Functional Analysis, 26(4):1190–1216, 2016.
  * [18] Eyal Lubetzky and Allan Sly.  Cutoff phenomena for random walks on random regular graphs.  2010\.
  * [19] Eyal Lubetzky and Allan Sly.  Explicit expanders with cutoff phenomena.  2011\.
  * [20] Eyal Lubetzky and Allan Sly.  Cutoff for the Ising model on the lattice.  Inventiones mathematicae, 191(3):719–755, 2013.
  * [21] Eyal Lubetzky and Allan Sly.  Cutoff for general spin systems with arbitrary boundary conditions.  Communications on Pure and Applied Mathematics, 67(6):982–1027, 2014.
  * [22] Eyal Lubetzky and Allan Sly.  Information percolation and cutoff for the stochastic Ising model.  Journal of the American Mathematical Society, 29(3):729–774, 2016\.
  * [23] Eyal Lubetzky and Allan Sly.  Universality of cutoff for the Ising model.  The Annals of Probability, 45(6A):3664–3696, 2017.
  * [24] J. C Mattingly, A. M Stuart, and D. J Higham.  Ergodicity for SDEs and approximations: locally Lipschitz vector fields and degenerate noise.  Stochastic processes and their applications, 101(2):185–232, 2002\.
  * [25] Mathieu Merle and Justin Salez.  Cutoff for the mean-field zero-range process.  2019\.
  * [26] P. Monmarche and M. Ramil.  Overdamped limit at stationarity for non-equilibrium Langevin diffusions.  Electron. Commun. Probab., 27:1–8, 2022.
  * [27] Danny Nam and Allan Sly.  Cutoff for the Swendsen-Wang dynamics on the lattice.  arXiv preprint arXiv:1805.04227, 2018.
  * [28] D. Le Peutrec and L. Michel.  Sharp spectral asymptotics for nonreversible metastable diffusion processes.  Probability and Mathematical Physics, 1(1):3–53, 2020.
  * [29] E. Priola, A. Shirikyan, L. Xu, and J. Zabczyk.  Exponential ergodicity and regularity for equations with Lévy noise.  Stochastic Processes and their Applications, 122(1):106–133, 2012\.
  * [30] Justin Salez.  Universality of cutoff for exclusion with reservoirs.  The Annals of Probability, 51(2):478–494, 2023.
  * [31] D. W Stroock and S. Karmakar.  Lectures on topics in stochastic differential equations, volume 68.  Springer Berlin, 1982.
  * [32] Seoyeon Yang.  Cutoff and dynamical phase transition for the general multi-component ising model.  Journal of Statistical Physics, 190(9):151, 2023.
  * [33] Z. Zhang and D. Chen.  A new criterion on existence and uniqueness of stationary distribution for diffusion processes.  Advances in Difference Equations, 2013(1):1–6, 2013.