# Light Curve Model for Luminous Red Novae and Inferences about the Ejecta of
Stellar Mergers

Tatsuya Matsumoto  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 Ave, New York, NY 10010, USA

###### Abstract

The process of unstable mass transfer in a stellar binary can result in either
a complete merger of the stars or successful removal of the donor envelope
leaving a surviving more compact binary. “Luminous red nova” (LRN) are the
class of optical transients believed to accompany such merger/common envelope
events. Past works typically model LRNe using analytic formulae for supernova
light curves which make assumptions (e.g., radiation dominated ejecta, neglect
of hydrogen recombination energy) not justified in stellar mergers due to the
lower velocities and specific thermal energy of the ejecta. We present a one-
dimensional model of LRN light curves, which accounts for these effects.
Consistent with observations, we find that LRNe typically possess two light
curve peaks, an early phase powered by initial thermal energy of the hot,
fastest ejecta layers and a later peak powered by hydrogen recombination from
the bulk of the ejecta. We apply our model to a sample of LRNe to infer their
ejecta properties (mass, velocity, and launching radius) and compare them to
the progenitor donor star properties from pre-transient imaging. We define a
maximum luminosity achievable for a given donor star in the limit that the
entire envelope is ejected, finding that several LRNe violate this limit.
Shock interaction between the ejecta and pre-dynamical mass-loss, may provide
an additional luminosity source to alleviate this tension. Our model can also
be applied to the merger of planets with stars or stars with compact objects.

XXX

##  1 Introduction

The direct interaction between the stars in a binary system is a common
occurrence (e.g., Sana et al. 2012), with implications for a number of topics
in stellar evolution (e.g., Podsiadlowski et al. 2004; de Mink et al. 2014; De
Marco & Izzard 2017) and energetic transients (e.g., Soker & Tylenda 2006;
Chevalier 2012; Smith 2014). Such interaction often begins after the
(typically more evolved) donor star overflows its Roche lobe and begins
transferring mass onto the more compact accretor star. Depending on the binary
mass ratio and eccentricity and the evolutionary state of the donor, this
mass-transfer can become unstable, resulting in a runaway growth of the mass
transfer rate (e.g., Soberman et al. 1997; Pavlovskii et al. 2017; Ge et al.
2020; Klencki et al. 2021).

Following an initial period of mass-loss, potentially lasting many orbits and
likely concentrated in the binary plane (e.g., from the $L_{2}$ Lagrange
point; Pejcha et al. 2016a, b, 2017; MacLeod & Loeb 2020) it is commonly
believed that the accretor will ultimately be dragged into the envelope of the
donor star (Paczynski, 1976; Taam et al., 1978; Iben & Livio, 1993; Ivanova et
al., 2013b). The final outcome of such a “common envelope” phase can be either
a complete merger of both stars (e.g., Kruckow et al. 2016), or the successful
ejection of the donor envelope and the survival of a tighter binary comprised
of the accretor and the core of the original donor (e.g., Ivanova & Nandez
2016; Sand et al. 2020; Law-Smith et al. 2020).

Stellar mergers and/or common envelope events have been invoked to give rise
to number of atypical stellar and binary populations, including cataclysmic
variables (Patterson, 1984), magnetic A-stars (Schneider et al., 2019),
magnetic white dwarfs (Tout et al., 2008; Nordhaus et al., 2011), sdB-stars
(e.g., Kramer et al. 2020), carbon-deficient giants (e.g., Bond 2019), and
blue stragglers (e.g., Davies et al. 2004; Ferreira et al. 2019; Wang et al.
2020), among other examples. The conditions required to remove the envelope,
and the final separations of the surviving binary, are also important for
predicting the population of merging binary compact objects by LIGO/Virgo
(e.g., Tauris et al. 2017; Belczynski et al. 2018; Vigna-Gómez et al. 2018).
However, despite increasingly sophisticated multi-dimensional simulations of
the common envelope phase over the past decade (e.g., Ricker & Taam 2012;
Ohlmann et al. 2016; Ivanova & Nandez 2016; MacLeod et al. 2018b; Prust &
Chang 2019; Chamandy et al. 2019; Law-Smith et al. 2020), the theoretical
community has yet to converge on a consistent physical picture which includes
all of the relevant physical processes (Clayton et al., 2017; Soker et al.,
2018; Wilson & Nordhaus, 2019; Reichardt et al., 2020). Observations remain
extremely important to unraveling this mystery.

The dynamical merger event may generate a short-lived transient powered by the
ejection of mass, primarily from the donor’s envelope (Soker & Tylenda, 2006).
“Luminous red novae” (LRN) are the class of optical transients with typical
durations of weeks to months and peak luminosities $\sim 10^{38}-10^{41}$ erg
s-1 between those of classical novae and supernovae (SNe) and
characteristically red colors (e.g., Martini et al. 1999; Munari et al. 2002;
Bond et al. 2003; Kulkarni et al. 2007; Tylenda et al. 2011; Kurtenkov et al.
2015; Smith et al. 2016; Blagorodnova et al. 2017; Cai et al. 2019;
Stritzinger et al. 2020; see Pastorello et al. 2019a for a recent review and
Fig. 1 for a sample of LRNe). The watershed event which confirmed the
association between LRNe and merging binary stars was the Galactic source
V1309 Sco (Mason et al., 2010), detected by OGLE as a contact binary for
almost a decade prior to the outburst (Tylenda et al., 2011). The observed
evolution of the binary light curve and orbital period leading up to the
merger, provided evidence for high rates of mass-loss, starting hundreds of
orbits prior to the final dynamical coalescence and concomitant abrupt rise to
optical maximum (Pejcha, 2014; Pejcha et al., 2017; MacLeod & Loeb, 2020).
Similar “precursor” emission, potentially powered by internal shocks in
outflowing streams from the $L_{2}$ point (e.g., Pejcha et al. 2016a, b), is
observed in the light curves of other LRNe, suggesting that pre-dynamical mass
loss is common if not ubiquitous in stellar mergers (Metzger & Pejcha 2017;
Blagorodnova et al. 2021).

The light curves of LRNe encode information on the merger ejecta, particularly
its mass, velocity, and launching radius (e.g., Ivanova et al. 2013a; MacLeod
et al. 2017; Blagorodnova et al. 2021). Accurate determinations of these
properties can, in combination with information on the progenitor binary from
pre-outburst imaging, be used to constrain the physical processes at work
leading up to and during the merger (e.g., MacLeod et al. 2017; Blagorodnova
et al. 2021). Ivanova et al. (2013a) pioneered this approach by applying an
analytic model for the luminosity and duration of Type II SNe (SNe II) light
curves from Popov (1993) to infer the ejecta properties for a sample of LRNe.
Using an updated version of the Popov (1993) formulae calibrated using SNe II
radiative transfer calculations (Sukhbold et al., 2016), Blagorodnova et al.
(2021) constrained the ejecta mass and launching radius for the AT2018bwo.
Combined with an application of binary stellar evolution models to the
observed yellow supergiant progenitor, the ejecta properties enabled them to
paint a consistent picture of the event (see also MacLeod et al. 2017 for a
similar analysis of M31-OT2015). Most LRN progenitor stars are moderately
evolved, with radii between $\sim 1-10$ times that of main sequence stars of
their mass (Fig. 2).

Despite the advances of these past works, several physical inconsistencies can
arise when applying SNe II models to LRN. SNe ejecta are dominated by
radiation pressure, an assumption implicit in SN II light curve models (Popov,
1993; Faran et al., 2019). By contrast, given the comparatively lower
velocities, higher densities, and lower thermal energies of LRN ejecta, gas
pressure can dominate over radiation pressure and sources of opacity other
than electron scattering may be important. Furthermore, while hydrogen
recombination could be an important source of luminosity for LRN (e.g.,
Ivanova et al. 2013a), recombination energy is irrelevant in SNe and is not
included in models like Popov (1993). SN light curves are powered by the
initial thermal energy of the recently-shocked ejecta from the explosion or
the radioactive decay of 56Ni. Although the early-time peaks seen in many LRN
light curves can be powered by the initial thermal energy of hot ejecta
(MacLeod et al., 2017; Metzger & Pejcha, 2017), the longer-lived plateau
emission (sometimes manifesting as a second luminosity peak) is powered by
hydrogen recombination energy or other forms of radially-distributed ejecta
heating (e.g., shock interaction between the fast merger ejecta and
circumbinary material from the pre-dynamical phase; Metzger & Pejcha 2017, or
from an accretion-powered jet; e.g., Soker 2020; Soker & Kaplan 2021).

In this paper we present a new model for LRN light curves, which self-
consistently accounts for the effects of arbitrary gas and radiation pressure,
thermal energy released by recombination, and a realistic distribution of
ejecta velocities. After presenting the light curve model in Section 2 we
describe our results in Section 3. After walking the reader through a few
example light curves (Sec. 3.1) and elucidating the physics giving rise to
multiple light curve peaks, we perform a parameter study covering the full
range of stellar merger properties (Sec. 3.2), apply our model to extract the
ejecta properties from a sample of individual LRNe (Sec. 3.3), and compare the
predictions of our light model to findings obtained using SN II models (Sec.
3.4). Section 4 discusses implications of our findings, including the need for
alternative energy sources in LRN from circumbinary interaction (Sec. 4.1);
the populations of binaries giving rise to LRNe (Sec. 4.2); and other
applications of our models to planet-star or star-compact object mergers (Sec.
4.3). We briefly summarize our conclusions in Section 5.

Figure 1: Bolometric light curves (bottom panel) and effective temperatures
(top panel) as a function of time since peak or discovery for a sample of
LRNe. See Table 1 for references to the data used to construct the figure.
Figure 2: (Left Panel) Progenitor star $R_{\star}$ versus progenitor star mass
$M_{\star}$ as inferred from pre-transient observations for a sample of LRNe
(Table 1). For comparison we show the terminal age main sequence (TAMS; black
line) and contours of surface escape velocity $v_{\rm
esc,\star}=(2GM_{\star}/R_{\star})^{1/2}$ (gray lines). (Right Panel) Ratio of
the ejecta velocity $v_{\rm obs}$ inferred from spectroscopy of the transient
emission to progenitor donor surface escape speed $v_{\rm esc,\star}$. As
expected, $v_{\rm obs}/v_{\rm esc,\star}$ is larger for more evolved
progenitors ($R_{\star}>R_{\rm TAMS})$ with well-developed core-envelope
structure, for which the ejecta launching radius $R_{0}$ can be $\ll
R_{\star}$.

##  2 Light Curve Model

Our light curve model follows that developed by Metzger et al. (2021) in the
different context of optical transients from the tidal disruption of a star by
a white dwarf, but with added physical effects such as radiation pressure and
Helium ionization (Kasen & Ramirez-Ruiz, 2010). The general approach is to
divide the radially-stratified ejecta into successive thin shells, calculate
the emission from each such mass/velocity shell separately by means of a time-
dependent one-zone model, and then obtain the total light curve by summing
over the contributions from each shell.

We begin by describing the one-zone model, as will ultimately be applied to
each shell of velocity $v$ and external mass $M_{\rm ej}(>v)$. Immediately
after the stellar merger or common envelope event, the evolution of the
unbound ejecta will be complex, potentially involving self-interaction between
different ejecta components (e.g., Metzger & Pejcha 2017; MacLeod et al.
2018b; Reichardt et al. 2019). However, we may reasonably assume that on a few
dynamical times the outflow geometry will become quasi-spherical and the
thermal and kinetic energies roughly equal. Hence we set the initial internal
energy of the shell as $E_{0}=dMv^{2}/2$ (here $dM$ is the shell’s mass) and
calculate the subsequent evolution starting from time $t_{0}=R_{0}/v$, where
$R_{0}$ is the initial (“ejection”) radius of the ejecta (typically equal to,
or smaller than, the radius of the donor star). After the initial thermal
energy has been converted into kinetic energy on a timescale $\sim t_{0}$, the
shell will thereafter expand freely with radius $R=vt$.

The internal energy of the shell can be written

$\displaystyle E$ $\displaystyle=\frac{3}{2}(1+\bar{x})Nk_{\rm
B}T+aT^{4}V+\sum_{i}NA_{i}x_{i}\varepsilon_{i}\ ,$ (1)

where $\bar{x}=\sum_{i}A_{i}x_{i}$, $A_{i}$, $x_{i}$, and $\varepsilon_{i}$
are the mean ionization degree, abundance fraction, degree of ionization, and
ionization energy of species $i$, respectively. We calculate $\bar{x}$ and
$x_{i}$ by solving Saha equation taking into account singly-ionized hydrogen
and helium, which we have checked capture recombination effects on the light
curve to high accuracy. The other quantities are the total number of nuclei
$N$, the Boltzmann constant $k_{\rm B}$, temperature $T$, radiation constant
$a$, and shell volume $V=4\pi R^{2}dR$ with the shell width $dR$.

The internal energy evolves according to the first law of thermodynamics:

$\displaystyle\frac{dE}{dt}=-\frac{3(\gamma_{3}-1)E}{t}-L\ ,$ (2)

where the effective adiabatic index is given by a function of density
$\rho=dM/V$ and temperature $T$ (Krishna Swamy, 1961; Kasen & Ramirez-Ruiz,
2010)

$\displaystyle\gamma_{3}(\rho,T)$ $\displaystyle\simeq
1+\frac{(1+\bar{x})(1+4\alpha)+\sum_{i}A_{i}\frac{\bar{x}(1-\bar{x})}{(2-\bar{x})}\bigl{(}\frac{3}{2}+\frac{\varepsilon_{i}}{k_{\rm
B}T}\bigl{)}}{(1+\bar{x})\big{(}\frac{3}{2}+12\alpha\big{)}+\sum_{i}A_{i}\frac{\bar{x}(1-\bar{x})}{(2-\bar{x})}\bigl{(}\frac{3}{2}+\frac{\varepsilon_{i}}{k_{\rm
B}T}\bigl{)}^{2}}\ ,$ (3)

and

$\displaystyle\alpha$ $\displaystyle\equiv\frac{P_{\rm rad}}{P_{\rm
gas}}=\frac{aT^{3}V}{3(1+\bar{x})Nk_{\rm B}}\ .$ (4)

is the ratio of radiation pressure $P_{\rm rad}=aT^{4}/3$ to gas pressure
$P_{\rm gas}=N(1+\bar{x})k_{\rm B}T/V$. The radiative luminosity of the shell
is given by

$\displaystyle L=\frac{E_{\rm rad}}{t_{\rm d}+t_{\rm lc}}\ ,$ (5)

where $E_{\rm rad}=aT^{4}V$ and the photon diffusion time through the column
associated with the external mass is

$\displaystyle t_{\rm d}$ $\displaystyle=\frac{\kappa M_{\rm ej}(>v)}{4\pi
cR}\ ,$ (6)

and the factor $t_{\rm lc}=R/c$ limits the energy escape time to the light
crossing time, where $c$ is the speed of light. It would be more accurate to
express Eq. (6) as an integral over the diffusion time through the external
layers $>v$ accounting for their respective opacity $\kappa$ and density
distribution. However, this approach would require following energy transport
between the shells via radiative diffusion; such a detailed treatment is
beyond the scope of this paper. Hence we adopt this simpler single-zone
expression for $t_{\rm d}$ because we find the calculated light curve turn
outs to be relatively insensitive to this assumption (a factor of 10 change in
the pre-factor of Eq. (6) changes the light curve duration by a factor
$\lesssim 2$).

We approximate the Rosseland mean opacity using the analytic formula similar
to that used by (Metzger & Pejcha, 2017)

$\displaystyle\kappa=\kappa_{\rm m}+\kappa_{\rm e}+\Big{(}\kappa_{\rm
H^{-}}^{-1}+\kappa_{\rm K}^{-1}\Big{)}^{-1}\ ,$ (7)

which accounts for the contribution from molecular opacity $\kappa_{\rm
m}\simeq 0.1\,Z\,\rm cm^{2}\,g^{-1}$, H- opacity $\kappa_{\rm H^{-}}\simeq
1.1\times 10^{-25}\,Z^{0.5}\rho^{0.5}T^{7.7}\,\rm cm^{2}\,g^{-1}$, Kramers
(bound-bound and bound-free) opacity $\kappa_{\rm K}\simeq 1.2\times
10^{26}\,Z(1+X)\rho T^{-7/2}\,\rm cm^{2}\,g^{-1}$, and electron scattering
$\kappa_{\rm es}\simeq 0.4\,\bar{x}\mu^{-1}\,\rm cm^{2}\,g^{-1}$. Here $Z$ is
the metallicity, $X$ is the mass fraction of hydrogen, and $\mu\equiv
dM/m_{\rm p}N$($\simeq 1.25$ for the solar metallicity) is the mean molecular
weight neglecting the number of electrons, where $m_{\rm p}$ is the proton
mass. We neglect dust opacity, which will become important once the dense
ejecta cools below the temperature for solid condensation; our model cannot
therefore be applied to the late phases of LRN emission when dust formation is
seen to occur (e.g., Wisniewski et al. 2008; Nicholls et al. 2013; Banerjee et
al. 2015; Kasliwal et al. 2017). Figure 3 compares our assumed opacity law
(Eq. 7) to that obtained using the OPAL opacity code for $Z=Z_{\odot}\simeq
0.02$.111https://opalopacity.llnl.gov/existing.html#type1GN Note that our
prescription of the electron scattering opacity in Eq. (7) is different from
that of Metzger & Pejcha (2017). Our formula is more accurate than theirs when
the recombination occurs at low densities $\rho\lesssim 10^{-11}{\,\rm
g\,cm^{-3}}$. We also adopted normalization of the Kramers opacity three times
larger than that used in Metzger & Pejcha (2017), to obtain better agreement
with the OPAL calculations.

Figure 3: Comparison of our assumed opacity law (dashed curve; Eq. 7) as a
function of gas temperature and density with results from the OPAL opacity
code (solid curves) for solar metallicity gas.

The total light curve is obtained by summing the luminosity contributions from
each mass/velocity shell (see Metzger et al., 2021). For a given mass profile
$M_{\rm ej}(>v)$ of the form discussed below, we divide the ejecta into thin
shells of velocity from $v$ to $v+dv$, which gives the shell width $dR=dvt$.
The mass contained in each shell is given by $dM=(dM_{\rm ej}/dv)dv$ and hence
its density is given by $\rho=dM/(4\pi R^{2}tdv)$. With these quantities one
can calculate the one-zone light curve for each shell $dL(t;v)$ and the total
light curve by summing over all such shells $L(t)=\int dL(t;v)$.

The ejecta velocity distribution depends on the physical mechanism responsible
for unbinding mass in stellar mergers and is uncertain theoretically. We
consider two general functional forms (exponential and power-law) for the
distribution of ejecta mass above a given velocity,

$\displaystyle M_{\rm ej}(>v)=M_{\rm ej}\begin{cases}e^{-(v/v_{\rm
ej})^{\beta_{\rm ex}}}&;\text{exponential}\ ,\\\ (v/v_{\rm ej})^{-\beta_{\rm
pl}}&;\text{power-law}\ ,\end{cases}$ (8)

where $\beta_{\rm ex}$, $\beta_{\rm pl}$, and $v_{\rm ej}$ are parameters and
the velocity range is $v_{\rm min}<v<v_{\rm max}$ ($v_{\rm min}=v_{\rm ej}$
for the power-law case). In what follows it will prove useful to introduce a
mean ejecta velocity $\bar{v}_{\rm E}$, which is defined by $M_{\rm
ej}\bar{v}_{\rm E}^{2}/2=\int v^{2}dM/2$ and related to $v_{\rm ej}$ entering
Eq. (8) according to

$\displaystyle\frac{\bar{v}_{\rm E}}{v_{\rm ej}}=\sqrt{\frac{\int_{v_{\rm
min}}^{v_{\rm max}}v^{2}dM}{v_{\rm ej}^{2}M_{\rm
ej}}}\simeq\begin{cases}\Gamma^{1/2}\left(\frac{\beta_{\rm ex}+2}{\beta_{\rm
ex}}\right)&;\text{exponential}\ ,\\\ \left(\frac{\beta_{\rm pl}}{\beta_{\rm
pl}-2}\right)^{1/2}&;\text{power-law}\ ,\end{cases}$ (9)

where $\Gamma(x)$ is the gamma function and we have assumed $v_{\rm min}\ll
v_{\rm ej}\ll v_{\rm max}$ in the exponential case, and $\beta_{\rm pl}>2$ and
$v_{\rm max}\gg v_{\rm min}$ in the power-law case.

To summarize, a given light curve model is defined by three main parameters:
total ejecta mass $M_{\rm ej}$, initial (ejection) radius $R_{0}$, and mean
ejecta velocity $\bar{v}_{\rm E}$ (equivalently, total kinetic energy $E_{\rm
K}=M_{\rm ej}\bar{v}_{\rm E}^{2}/2$), which is related to $v_{\rm ej}$ through
Eq. (9). In what follows, we will sometimes scale the ejecta mass to a fixed
fraction $\eta\equiv M_{\rm ej}/M_{\star}$ of the mass of the donor star
$M_{\star}$. We shall adopt a fiducial value of $\eta\simeq 0.25$. This is
somewhat larger than $\eta\simeq 0.1$ found in the common-envelope simulations
by MacLeod et al. (2018b) for a fixed accretor to donor mass ratio, 1/3. The
envelope ejection efficiency $\eta$ may also be even smaller $\lesssim 0.1$,
in the case of main sequence stellar mergers (e.g., Lombardi et al. 2002;
Glebbeek et al. 2013) or mergers of high mass-ratio binaries (e.g.,
Blagorodnova et al. 2021).

We likewise will frequently scale the mean ejecta velocity $\bar{v}_{\rm E}$
to the escape velocity of the donor at the outflow launching radius,

$v_{\rm esc}=\left(\frac{2GM_{\star}}{R_{0}}\right)^{1/2}\simeq 340\,{\rm
km\,s^{-1}\,}\left(\frac{M_{\star}}{3\,{\mathrm{M}_{\odot}}}\right)^{1/2}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)^{-1/2}\
,$ (10)

where $G$ is the gravitational constant. If $R_{0}$ equals the radius of the
donor star $R_{\star}$, then $v_{\rm esc}=v_{\rm esc,\star}$ is the surface
escape speed. However, for $R_{0}\lesssim R_{\star}$ then $v_{\rm esc}\gtrsim
v_{\rm esc,\star}.$ Ejecta speeds spanning $\bar{v}_{\rm E}\sim 0.3-3\,v_{\rm
esc,\star}$ are supported by LRN observations (Fig. 2).

The choice of velocity distribution (Eq. 8) introduces additional parameters,
$\\{\beta_{\rm ex}(\beta_{\rm pl}),v_{\rm min},v_{\rm max}\\}$. Typically, we
set $v_{\rm min}=0.1\,\bar{v}_{\rm E}$ and $v_{\rm max}=3\,\bar{v}_{\rm E}$
but the light curve is not sensitive to these choices as long as $v_{\rm
min}\ll\bar{v}_{\rm E}\ll v_{\rm max}$. The calculation also depends on the
assumed composition of the ejecta through the opacity; we assume solar
metallicity ($X=0.74$ and $Z=0.02$) throughout this work, motivated by the
moderately-evolved donor stars of observed LRNe progenitor binaries (Fig. 2).

##  3 Results

###  3.1 Example Light Curves

Figure 4 shows an example light curve, calculated for a total ejecta mass
$M_{\rm ej}=1\,{\mathrm{M}_{\odot}}$, initial radius
$R_{0}=10\,{\mathrm{R}_{\odot}}$, and velocity $\bar{v}_{\rm E}=300$ km s-1.
We employ the exponential velocity distribution, $M_{\rm ej}(>v)\propto
e^{-(v/v_{\rm ej})^{3}}$ (Eq. 9). A black line shows the total light curve,
while the colored lines show the separate contributions from distinct mass
shells (from $v=50$ km s-1 to $v=650$ km s-1 with an interval of $50\,\rm
km\,s^{-1}$ and shell width of $dv=5$ km s-1).

The light curve is characterized by an initial peak of high luminosity and
short duration $\lesssim 1$ week followed by longer duration plateau lasting
several months. Such double-peaked light curve behavior is common in LRNe
(MacLeod et al. 2017; Metzger & Pejcha 2017; Pastorello et al. 2021a; Fig. 1).
As we now describe, these two components have distinct physical origins.

Figure 4:  Example light curve for fiducial parameters: $M_{\rm
ej}={\mathrm{M}_{\odot}}$, $R_{0}=10\,{\mathrm{R}_{\odot}}$, $\bar{v}_{\rm
E}=300{\,\rm km\,s^{-1}}$ ($v_{\rm ej}\simeq 320\,\rm km\,s^{-1}$). Colored
curves denote one-zone light curves for shells with $v=50$ to $650{\,\rm
km\,s^{-1}}$ (violet to orange) with a shell width $dv=5\,\rm km\,s^{-1}$. The
total light curve is characterized by two peaks, early cooling peak ($\lesssim
10\,\rm d$) and late recombination plateau, produced by distinct physical
processes. The gray dashed curve shows an otherwise equivalent model but
neglecting energy release from recombination. The approximate effective
temperature of the emission is shown in the top panel (Eq. 19).  Figure 5:
Colored lines show trajectories in the density-temperature plane of mass
shells with $v=50$ to $650\,\rm km\,s^{-1}$ from our fiducial model (following
the same labeling scheme as in Fig. 4). Background contours indicate the
opacity $\kappa$ (left panel) and effective adiabatic index $\gamma_{3}$
(right panel). Dotted and dashed contours indicate different hydrogen
ionization fractions $\bar{x}$ and ratios of gas-to-radiation pressure. Prior
to hydrogen recombination, each shell decompresses and cools adiabatically
with $\gamma_{3}\in[4/3,5/3],$ moving from the upper right hand corner of the
plot to the lower left. However, once recombination starts at $T\simeq 10^{4}$
K, the energy released from recombination holds the temperature roughly
constant at $T\simeq 10^{4}$ K ($\gamma_{3}\simeq 1.1$), causing a convergence
of the mass shell trajectories. However, as the density continues to drop,
once $P_{\rm rad}>P_{\rm gas}$ ($\rho\lesssim 10^{-11}$ g cm-3 for $T\sim
T_{\rm i}\simeq 5000$ K), $\gamma_{3}$ again starts to rise, causing $T$ to
decrease faster and recombination to complete. The resulting decrease in the
opacity leads to the light curve peaking, thus motivating our analytic
estimate of the transient duration (Eq. 16).  Figure 6:  Time evolution of key
quantities as measured at the shell which dominates the luminosity at a given
epoch, for the fiducial model shown in Fig. 4. The top panel shows the
effective adiabatic index $\gamma_{3}$ (Eq. 3); for comparison dashed lines
show the usual values for radiation- (gas-) pressure-dominated gas,
$\gamma_{3}=4/3$ ($5/3$), absent recombination effects. A dash-dotted curve in
the top panel shows the value of $\gamma_{3}$ for the single fixed shell with
$v=\bar{v}_{\rm E}$.

####  3.1.1 Early Cooling Peak

At early stages of evolution, the fastest, lower-mass layers of the ejecta
dominate the transient luminosity. As we will confirm below, these layers are
typically dominated by radiation pressure and their energy source is “cooling”
emission produced by the initially hot ejecta (e.g., MacLeod et al. 2017),
similar to the mechanism responsible for powering the plateau of SNe II
(Arnett, 1980). The duration and luminosity of this early emission phase are
given by (e.g., Arnett, 1980; Padmanabhan, 2000)

$\displaystyle t_{\rm pk}\simeq\biggl{(}\frac{\kappa_{\rm es}M}{4\pi
cv}\biggl{)}^{1/2}\simeq 6.7{\,\rm
d\,}\biggl{(}\frac{M}{10^{-2}\,{\mathrm{M}_{\odot}}}\biggl{)}^{1/2}\biggl{(}\frac{v}{500\,\rm
km\,s^{-1}}\biggl{)}^{-1/2}\ ,$ (11) $\displaystyle L_{\rm pk}$
$\displaystyle\simeq$ $\displaystyle\frac{\frac{1}{2}Mv^{2}(R_{0}/R_{\rm
pk})}{t_{\rm pk}}\simeq 1.0\times 10^{39}{\,\rm erg\,s^{-1}}\times$ (12)
$\displaystyle\biggl{(}\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\biggl{)}\biggl{(}\frac{v}{500\,\rm
km\,s^{-1}}\biggl{)}^{2}\ ,$ $\displaystyle E_{\rm pk}$ $\displaystyle\simeq$
$\displaystyle L_{\rm pk}t_{\rm pk}\simeq 6\times 10^{44}{\,\rm erg\,}\times$
$\displaystyle\left(\frac{M}{10^{-2}\,{\mathrm{M}_{\odot}}}\right)^{1/2}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)\left(\frac{v}{500\rm\,km\,s^{-1}}\right)^{3/2},$

where $M$ is the ejecta mass of the fastest layer, $R_{\rm pk}=vt_{\rm pk}$,
and we have assumed that electron scattering of ionized hydrogen $\kappa_{\rm
es}\simeq 0.32$ cm2 g-1 dominates the opacity during the early peak.

The assumptions entering Eqs. (11)-(3.1.1) can now be checked. The ratio of
radiation to gas pressure, on the timescale $\sim t_{\rm pk}$, can be
estimated as

$\displaystyle\frac{P_{\rm rad}}{P_{\rm gas}}\biggl{|}_{t_{\rm pk}}$
$\displaystyle\simeq
5.1\left(\frac{M}{10^{-2}{\mathrm{M}_{\odot}}}\right)^{-1/4}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)^{3/4}\left(\frac{v}{500\,{\rm
km\,s^{-1}}}\right)^{3/2}$ $\displaystyle\simeq
0.98\,\eta_{0.25}^{-1/4}\left(\frac{v}{v_{\rm
esc}}\right)^{3/2}\left(\frac{M_{\star}}{3\,{\mathrm{M}_{\odot}}}\right)^{1/2},$
(14)

where (assuming $P_{\rm rad}\gg P_{\rm gas}$) we have used
$Mv^{2}/2=(4\pi/3)R_{0}^{3}aT_{0}^{4}$ to derive the initial temperature and
$T(t_{\rm pk})=T_{0}(R_{0}/R_{\rm pk})$. In the second line of Eq. (14) we
have set $M=M_{\rm ej}$ and scaled the ejecta mass to donor star mass
$M_{\star}$ via $\eta=M_{\rm ej}/M_{\star}$ ($\eta_{0.25}=\eta/0.25$) and the
ejecta velocity to the escape speed from the outflow launching radius (Eq.
10). Note that as long as the radiation pressure dominates initially, the
pressure ratio does not evolve significantly until the light curve
peaks.222Neglecting the ionization energy and radiation energy loss, we can
solve Eqs. (1), (2), and (3) to obtain the time evolution of the pressure
ratio ($\alpha=P_{\rm rad}/P_{\rm gas}$, see Eq. 14):
$\displaystyle\frac{t}{t_{0}}=\biggl{(}\frac{\alpha}{\alpha_{0}}\biggl{)}^{-1/3}\biggl{(}\frac{2\alpha+1}{2\alpha_{0}+1}\biggl{)}^{3}e^{\frac{32}{3}(\alpha_{0}-\alpha)}\
,$ where $\alpha_{0}$ is the initial ratio at $t_{0}=R_{0}/v$. We see that if
the radiation pressure dominates even a little over gas pressure,
$\alpha_{0}\gtrsim$ a few, the pressure ratio does not evolve until radiative
cooling becomes important. This is because the radiation entropy, which is
conserved under adiabatic expansion, is proportional to the pressure ratio.
This allows us to estimate the minimum velocity of the radiation-pressure
dominated shell, $v\simeq 500\,\rm km\,s^{-1}$ for our fiducial model with
$M\propto v^{3}e^{-(v/v_{\rm ej})^{3}}$, consistent with Fig. 4.

We see that the radiation pressure will dominate in the high-velocity tail of
the ejecta $v\gg v_{\rm esc}$ and that Eqs. (11)-(3.1.1) indeed determine the
early-time peak emission from these layers (and hence the early-peak of the
transient as a whole). However, for the bulk of the ejecta for which
$v\lesssim v_{\rm esc}$, gas pressure will dominate and the early cooling
emission from these layers will be reduced compared to $L_{\rm pk}$ in Eq.
(12) because of the greater adiabatic losses experienced by the ejecta for
$\gamma=5/3$ ($T\propto R^{-2}$) than $\gamma=4/3$. Previous models for LRNe
employ analytic formula for the luminosity and duration (e.g., Popov 1993)
which implicitly assume radiation pressure dominates, as is justified in SNe
due to their higher ejecta velocities and specific internal energy.

The ratio of Kramers to electron scattering opacity during the early peak can
also be checked:

$\displaystyle\frac{\kappa_{\rm es}}{\kappa_{\rm K}}\biggl{|}_{t_{\rm pk}}$
$\displaystyle\simeq$ $\displaystyle
5.3\left(\frac{M}{10^{-2}\,{\mathrm{M}_{\odot}}}\right)^{-3/8}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)^{7/8}\left(\frac{v}{500\,{\rm
km\,s^{-1}}}\right)^{3/2}\ ,$ (15)

where again we take $\kappa_{\rm es}\simeq 0.32$ cm2 g-1. This confirms that
$\kappa_{\rm es}\gtrsim\kappa_{\rm K}$, as assumed in Eqs. (11)-(3.1.1).

####  3.1.2 Late Recombination Plateau

The longer plateau-like emission phase in Fig. 4 is not powered by the
ejecta’s initial thermal energy, but instead by hydrogen recombination. After
forming its initial luminosity peak, the ejecta cool down and hydrogen starts
to recombine at temperatures $T\sim 10^{4}\,\rm K$. Once recombination sets
in, the additional energy supply causes the adiabatic index to decrease to
$\gamma_{3}-1\simeq 0.1$, thus slowing the temperature evolution greatly
compared to that of ordinary adiabatic expansion absent recombination
$(\gamma_{3}\lesssim 5/3)$. By the time a given shell completes recombination
at temperatures $T_{\rm i}\simeq 4000-6000\,\rm K$ most of the radiation
energy has been emitted, generating a second peak in its light curve. The sum
of the second peak from all mass shells powers the $\sim 100$ day plateau in
Fig. 4.

The duration of the plateau is roughly given by the timescale over which the
recombination completes at $T\gtrsim T_{\rm i}$. We find that recombination
completes at roughly a constant characteristic density $\rho_{\rm i}\sim
10^{-11}\,\rm g\,cm^{-3}$. As explained in Fig. 5, this density corresponds to
roughly where $P_{\rm rad}\sim P_{\rm gas}$ for $T\sim T_{\rm i}$. The
increasing adiabatic index once $P_{\rm rad}\gtrsim P_{\rm gas}$ triggers the
completion of ionization of each mass shell, causing the opacity to decrease
and leading to the light curve peak.

The duration of the recombination plateau is thus set by the timescale for the
density to reach $\rho_{\rm i}$,

$\displaystyle t_{\rm pl}\sim\biggl{(}\frac{3M_{\rm ej}}{4\pi\rho_{\rm
i}\bar{v}_{\rm E}^{3}}\biggl{)}^{1/3}\simeq 140{\,\rm d\,}\rho_{\rm
i,-11}^{-1/3}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{1/3}\biggl{(}\frac{\bar{v}_{\rm
E}}{300\,{\rm km\,s^{-1}}}\biggl{)}^{-1}\ ,$ (16)

where $\rho_{\rm i,-11}=\rho_{\rm i}/10^{-11}\,\rm g\,cm^{-3}$. Since the
plateau is powered by the hydrogen recombination ($\varepsilon_{\rm H}=13.6$
eV per hydrogen atom), the radiated energy and luminosity can be written

$\displaystyle E_{\rm pl}$ $\displaystyle\sim$ $\displaystyle f_{\rm
ad}\frac{\varepsilon_{\rm H}XM_{\rm ej}}{m_{\rm p}}\simeq 5.8\times
10^{45}{\,\rm erg\,}f_{\rm ad,0.3}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}\ ,$ (17) $\displaystyle L_{\rm pl}$
$\displaystyle\sim$ $\displaystyle\frac{E_{\rm pl}}{t_{\rm pl}}\simeq
4.8\times 10^{38}{\,\rm erg\,s^{-1}\,}f_{\rm ad,0.3}\rho_{\rm
i,-11}^{1/3}\times$ (18) $\displaystyle\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{2/3}\biggl{(}\frac{\bar{v}_{\rm
E}}{300\,{\rm km\,s^{-1}}}\biggl{)}\ ,$

where the dimensionless parameter $f_{\rm ad}=0.3f_{\rm ad,0.3}<1$ accounts
for inefficiency in radiating the recombination energy due to adiabatic
losses. Helium recombination also contributes, but since it occurs at higher
temperature (at hence much higher optical depth), its effect on the light
curve is small compared to that of hydrogen recombination.

In Appendix A we present a one-zone toy light curve model for a hydrogen-rich
explosion, in the spirit of Popov (1993) but assuming gas-pressure instead of
radiation dominates and accounting for recombination energy through an
(assumed constant) effective value of the adiabatic index $\gamma_{3}$; for
values of $\gamma_{3}\sim 1$ expected during recombination the analytic
expressions for the plateau luminosity (Eq. A14), duration (Eq. A13), and
radiated energy (Eq. A15), are quite similar to those Eqs. (16)-(18) above.
This is remarkable because the conceptual motivation of the two derivations
are different: the plateau duration derived in Appendix A (Eq. A13) is
dictated by the photon diffusion timescale through the recombining ejecta
shell, while Eq. (16) was derived assuming recombination occurs at roughly
fixed density.

The top panel of Fig. 4 shows the evolution of the effective temperature,
which we have defined by

$\displaystyle T_{\rm eff}\equiv\biggl{(}\frac{L}{4\pi\sigma_{\rm SB}R_{\rm
lum}^{2}}\biggl{)}^{1/4}\ ,$ (19)

where $\sigma_{\rm SB}$ is the Stefan-Boltzmann constant and $R_{\rm lum}$ is
the radius of the shell contributing the majority of the luminosity at any
epoch. During the initial peak, $T_{\rm eff}\simeq 10^{4}\,\rm K$ before
cooling during the recombination plateau to an almost fixed valued $T_{\rm
eff}\simeq 3000-4000\,\rm K$, consistent with the “red” colors of LRNe. This
temperature stability occurs for the same reason as in SNe II, namely that the
recombining shell recedes monotonically in the velocity coordinate but the
expansion effect compensates and leads to the almost spatially fixed location
of the shell (see Fig. 6). We emphasize however, that$-$unlike in SNe$-$the
plateau is powered not by radioactivity or the thermal energy released from
the initial explosion, but rather from the recombination energy itself. We
confirm this by showing an otherwise identical light curve model but
calculated with recombination energy artificially excluded, as a gray line in
Fig. 4, demonstrating the absence in this case of an extended plateau phase.
This further implies that standard relationships applied to model SN light
curves (e.g., Popov 1993) should not be applied to LRNe.

Figure 6 shows the time evolution of key quantities for the fiducial model, as
measured at the shell dominating the luminosity at each epoch ($R_{\rm lum}$).
As expected, $P_{\rm rad}\gg P_{\rm gas}$ during the early peak when the
highest velocity layers dominate the emission, while $P_{\rm rad}\lesssim
P_{\rm gas}$ during the later recombination peak. Also note the mean
ionization fraction $\bar{x}$ is low $\lesssim 10^{-2}$ during the prolonged
recombination plateau, and that electron scattering opacity $\kappa_{\rm
es}\simeq 0.32\,\bar{x}$ cm2 g-1 becomes comparable to H- opacity. During the
recombination phase, the effective adiabatic index $\gamma_{3}\lesssim 4/3$
due to the energy released by recombination offsetting losses due to adiabatic
expansion. Note that although the index $\gamma_{3}$ of $R_{\rm lum}$ is close
to $4/3$, this is not related to radiation pressure. Rather, during the brief
epoch when a shell dominates the total luminosity, its adiabatic index quickly
increases from $\gamma_{3}\simeq 1.1$ to $5/3$ passing through an value
coincidentally close to $4/3$ (as illustrated by a dash-dotted curve in Fig.
6).

####  3.1.3 Variations About the Fiducial Model

Figure 7: Impact on the fiducial model light curve (Fig. 4) of changing the
shape of the ejecta velocity distribution (Eq. 8), for a fixed mean velocity
$\bar{v}_{\rm E}$. The left panel shows the effect of varying the parameter
$\beta_{\rm ex}$ which enters the exponential velocity distribution $M_{\rm
ej}(>v)\propto e^{-(v/v_{\rm ej})^{\beta_{\rm ex}}}$, while the right panel
shows the impact of varying the power-law index $\beta_{\rm pl}$ for $M_{\rm
ej}(>v)\propto v^{-\beta_{\rm pl}}$. Although the light curve shape changes
for different velocity distributions, the qualitative features of an early
peak followed by a second peak or plateau are robust to the details.

Figure 7 shows light curve models calculated for ejecta properties otherwise
identical to the fiducial model (exponential model with $\beta_{\rm ex}=3$)
but varying the shape of the velocity distribution at fixed $\bar{v}_{\rm E}$
(Eq. 8). For the exponential velocity profile, the duration of the plateau
shortens with increasing values of $\beta_{\rm ex}$, but the qualitative shape
does not change significantly unless $\beta_{\rm ex}\lesssim 2$. Likewise, in
the case of a power-law profile, the light curve shape is not qualitatively
dependent on the power-law index $\beta_{\rm pl}$. Motivated thus, we
hereafter fix the velocity profile to the fiducial exponential $\beta_{\rm
ex}=3$.

Figure 8: Comparison of the light curve for the fiducial model ($M_{\rm
ej}=1\,{\mathrm{M}_{\odot}}$; $R_{0}=10\,{\mathrm{R}_{\odot}}$; $\bar{v}_{\rm
E}=300$ km s-1; Fig. 4) to otherwise similar models in which: the ejecta mass
is increased to $3\,{\mathrm{M}_{\odot}}$ (blue); the initial radius is
increased to $R_{0}=30\,{\mathrm{R}_{\odot}}$ (red); the mean velocity is
increased to $\bar{v}_{\rm E}=900$ km s-1 (orange). The trends in the
luminosity and duration of the first and second peaks follow Eqs. (11), (12)
and Eqs. (16), (18), respectively.

Figure 8 compares the light curve for the fiducial model ($M_{\rm
ej}=1\,{\mathrm{M}_{\odot}}$; $R_{0}=10\,{\mathrm{R}_{\odot}}$; $\bar{v}_{\rm
E}=300$ km s-1) to otherwise similar models but for which $\\{M_{\rm
ej},R_{0},\bar{v}_{\rm E}\\}$ are tripled from their fiducial values. As
expected from Eqs. (11)-(3.1.1), larger values of $R_{0}$ and $\bar{v}_{\rm
E}$ increase the luminosity of the early peak. As expected from Eqs.
(16)-(18), the properties of the recombination plateau are roughly independent
of $R_{0}$. Likewise, the plateau gets shorter but more luminous with
increasing velocity, while both the duration and luminosity increase for
larger ejecta mass.

As can be seen in Fig. 1, many LRNe light curves exhibit multi-peak light
curves in which the first peak is substantially brighter than the longer
plateau to follow (e.g., V838 Mon: Munari et al. 2002; M31-OT2011: Kurtenkov
et al. 2015; MacLeod et al. 2017; AT2019zhd: Pastorello et al. 2021a;
AT2020kog: Pastorello et al. 2021b) and others in which a distinct first peak
is more subtle or non-existent (e.g., V1309 Sco: Mason et al. 2010; Tylenda et
al. 2011; AT2020hat: Pastorello et al. 2021b). This difference may indicate
the presence or absence of a high-velocity ejecta component or different
progenitor’s properties. To understand the latter effect, using Eqs. (12) and
(18), we see that the ratio of the cooling envelope peak luminosity to that of
the recombination plateau can be estimated as

$\displaystyle\frac{L_{\rm pk}}{L_{\rm pl}}$ $\displaystyle\sim$
$\displaystyle 3\,f_{\rm ad,0.3}^{-1}\left(\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\right)^{-2/3}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)\left(\frac{\bar{v}_{\rm
E}}{300\,{\rm km\,s^{-1}}}\right)$ (20) $\displaystyle\sim$ $\displaystyle
4\,f_{\rm
ad,0.3}^{-1}\eta_{0.25}^{-2/3}\left(\frac{M_{\star}}{3\,{\mathrm{M}_{\odot}}}\right)^{-1/6}\left(\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\right)^{1/2}\left(\frac{\bar{v}_{\rm
E}}{v_{\rm esc}}\right)\ ,$

where we have taken $\rho_{\rm i,-11}=1$ and assumed that the outer layers of
the ejecta responsible for powering the first peak have a velocity twice as
large as that of the bulk, $v\simeq 2\,\bar{v}_{\rm E}$. Thus, a range of
values $L_{\rm pk}/L_{\rm pl}\sim 1-10$ can be obtained.

###  3.2 Parameter Study of Stellar Mergers Across the HR Diagram

We now explore the light curve behavior across a larger parameter space of
ejecta mass $M_{\rm ej}\sim 10^{-3}-10\,{\mathrm{M}_{\odot}}$ and initial
radius $R_{0}\sim 10^{-2}-100\,{\mathrm{R}_{\odot}}$ covering those expected
from a wide range of stellar and even sub-stellar merger events. For each
model we set the mean ejecta velocity to the escape speed from $R_{0}$,
$\bar{v}_{\rm E}=v_{\rm esc}$ (Eq. 10), where the stellar mass is related to
the ejecta mass by $M_{\rm ej}=\eta M_{\star}$ with $\eta=0.25$. For each
light curve model, we calculate the total radiated energy $E_{\rm rad}=\int
Ldt$. We also define a characteristic transient duration, $t_{90}$, as the
interval over which 90% of the total energy is radiated, $0.9\,E_{\rm
rad}=\int^{t_{90}}Ldt$. Finally, from the radiated energy and duration we
define a mean luminosity $L_{90}\equiv E_{\rm rad}/t_{90}$.

Figure 9 shows the distribution of the characteristic luminosity, duration,
and radiated energy. We also plot the ratio of the radiation energy to the
total hydrogen recombination energy $f_{\rm ad}\equiv E_{\rm
rad}/(\varepsilon_{\rm H}XM_{\rm ej}/m_{\rm p}$). For comparison we show
tracks of stellar radius versus terminal-age main sequence (TAMS), from the
MIST database of MESA stellar evolution models (Dotter, 2016; Choi et al.,
2016; Paxton et al., 2011, 2013, 2015, 2019), as well as a hypothetical state
of post-main sequence evolution in which the stellar radius is 10 times larger
than its TAMS value (“10TAMS”). These bracket the state of LRN progenitor
stars from pre-transient imaging studies (Fig. 2). For matter ejected from the
stellar surface, we expect $R_{0}\sim R_{\star}$, while for matter ejected
from deeper within the donor star $R_{0}<R_{\star}$.

Across most of the parameter space, we find that the total radiated energy is
dominated by the later recombination plateau, consistent with a “by eye”
integration of LRN light curves (Fig. 1). The parameter dependencies of the
luminosity, duration, and radiated energy are roughly consistent with those
expected from Eqs. (16)-(18). In particular, for $\bar{v}_{\rm E}=v_{\rm
esc}\propto M_{\star}^{1/2}R_{0}^{-1/2}$ and $M_{\rm ej}\propto M_{\star}$,
these formulae predict $t_{90}\sim t_{\rm pl}\propto
M_{\star}^{1/6}R_{0}^{1/2}$, $L_{90}\sim L_{\rm pl}\propto
M_{\star}^{7/6}R_{0}^{-1/2}$, and $E_{\rm rad}\sim E_{\rm pl}\propto
M_{\star}$, consistent with the contours in Fig. 9. These relations no longer
hold when the ejecta is already partially ionized at the initial time (upper
left and bottom regions in Fig. 9). Figure 10 shows the ratios of the
transient luminosity and duration as predicted by our analytic formulae (Eqs.
11,12) to those obtained from the full numerical calculations (Fig. 9),
demonstrating good agreement (to within a factor $\lesssim 3$) across a band
encompassing the entire TAMS sequence.

For mergers that occur when the donor is near the end of the main sequence
($R\sim R_{\rm TAMS}$) or after evolving moderately up the giant branch
($R\sim 10R_{\rm TAMS}$) the predicted luminosities are in the range $\sim
10^{39}-10^{40}$ erg s-1, broadly consistent with the bulk of the LRN
population. Luminosities up to $\sim 10^{41}$ erg s-1 may be possible in
extreme cases for efficient ejection of the stellar envelope $\eta\sim 1$ from
deep within the potential of the stellar core ($R_{0}\ll R_{\rm TAMS}$). The
minimum launching radius of the donor star envelope for an evolved star is
determined by the size of the final He core, which we show as a white dash-
dotted curve in Fig. 9 (we adopt the analytic $M_{\star}-M_{\rm He}$
relationship from Woosley 2019 and employ the $M_{\rm He}-R_{\star}$
relationship from Eq. 4 of Schaerer & Maeder 1992; note that $\eta\lesssim
0.5-0.9$ in this case). Even under the most optimistic assumptions, the model
is strained to reach luminosities $L>10^{41}$ erg s-1, which can therefore be
taken as a theoretical upper limit (absent other energy sources not included
in our fiducial model; Section 4.1).

Typical values of the radiative efficiency of the recombination energy are
$f_{\rm ad}=E_{\rm rad}/(\varepsilon_{\rm H}XM_{\rm ej}/m_{\rm p})\simeq
0.1-0.3$. Across most of the parameter space, the luminosity of the initial
thermal peak $L_{\rm pk}$ will stick out above that of the recombination
plateau $L_{\rm pl}$, even though the radiated energy is dominated by the
plateau. The parameter space for which $L_{\rm pk}>L_{\rm pl}$ is delineated
by a solid black line in Fig. 9, roughly consistent with our analytic estimate
in Eq. (20). Note however that having $L_{\rm pk}>L_{\rm pl}$ does not
necessarily imply the early cooling peak will be detected, because its
duration is significantly shorter than the plateau, $t_{\rm pk}/t_{\rm
pl}\lesssim 0.01-0.1$, and could be missed depending on the discovery epoch
and observational cadence of the survey (Fig. 1).

Figure 9: Characteristic luminosity, timescale, radiated energy, and ratio of
radiated energy to recombination energy calculated for a range of stellar
(donor) mass and radius for our fiducial velocity profile (Eq. 8). Ejecta mass
is linked to the stellar mass via $M_{\rm ej}=\eta M_{\star}$ ($\eta=0.25$)
and velocity is set to the escape velocity. White curves denote the mass and
radius relation for TAMS, ten times larger radius of TAMS (“10TAMS”), and
brown dwarf and giant planet. A white dash-dotted curve shows the radius of He
core for the corresponding ZMAS mass ($M_{\star}$), which gives a rough lower
limit of the mass-ejection radius $R_{0}$ from an evolved star. A black solid
line denotes the contour along which the peak luminosities of the initial
cooling phase and later recombination plateau are equal (see Fig. 4), as
calculated directly from the light curve models. A gray dotted curve shows the
contour on which the mean ionization degree of the typical shell
($v=\bar{v}_{\rm E}$) is $\bar{x}=0.9$ (partially ionized), at the initial
time.  Figure 10: The same parameter space shown in Fig. 9 but showing the
ratio of the plateau luminosity and duration predicted by our analytical
formulae (Eqs. 16,18; for fixed $f_{\rm ad,0.3}=1$ and $\rho_{\rm i,-11}=1$)
to those obtained from the full numerical model (Fig. 9). The analytic
formulae accurately predict the luminosity and duration within a factor of
$\lesssim 3$ along the parameter space surrounding the TAMS sequence.

###  3.3 Application to Individual LRNe

Figure 11 shows the mean luminosity $L_{90}$ and duration $t_{90}$ for a
sample of LRNe (Table 1; Fig. 1) as a function of the donor mass $M_{\star}$
inferred from pre-transient imaging observations when available.333We note
that the measured duration $t_{90}$ and accordingly the radiated energy
$E_{\rm rad}$ may be underestimated in some events because of incomplete light
curve coverage (some events were discovered only at the the end of the plateau
phase and others were not observed until the end of the plateau; see Fig. 1).
Observed quantities are shown by circles or dot-dashed lines (when no
information on $M_{\star}$ is available). For comparison in Fig. 11 we show
our light curve predictions (Fig. 9) assuming $M_{\rm ej}=\eta M_{\star}$ and
$\bar{v}_{\rm E}=v_{\rm esc,\star}=(2GM_{\star}/R_{\star})^{1/2}$ for stars on
the TAMS (red curves) and with radii 10 times their TAMS radius (10TAMS; blue
curves), for a range of values $\eta=10^{-2}$ (thinnest lines) to $\eta=1$
(thickest lines). For fixed $\eta$, the theoretical relations roughly follow
the trends $L_{90}\propto M_{\star}$ and $t_{90}\propto M_{\star}^{1/2}$
predicted by Eqs. (16,18).

Figure 11 shows that, in broad brush terms, our calculations can reproduce the
observed luminosities $\sim 10^{38}-10^{40}$ erg s-1 of most LRNe, and also
the general trend of higher luminosity and longer duration for more massive
stellar progenitors (e.g., Kochanek et al. 2014; Blagorodnova et al. 2021).
However, there are several exceptions, including OT2011 in NGC4490
($L_{90}\sim 2\times 10^{41}\,\rm erg\,s^{-1}$; Pastorello et al. 2019a),
AT2014ej ($L_{90}\sim 2\times 10^{41}\,\rm erg\,s^{-1}$; Stritzinger et al.
2020), and AT2017jfs ($L_{90}\sim 3\times 10^{41}$ erg s-1; Pastorello et al.
2019b), with luminosities that approach or exceed even the $\eta=1$ TAMS
prediction for very massive stars. We return to this tension regarding the
implied ejecta masses of some LRNe below.

We can determine the ejecta properties for individual LRNe by inverting our
analytic formulae for the transient luminosity and duration (Eqs. 16, 18),

$\displaystyle M_{\rm ej}$ $\displaystyle\simeq
1.6\,{\mathrm{M}_{\odot}}\,f_{\rm ad,0.3}^{-1}\biggl{(}\frac{t_{\rm
pl}}{100{\,\rm d}}\biggl{)}\biggl{(}\frac{L_{\rm pl}}{10^{39}\,\rm
erg\,s^{-1}}\biggl{)}\ ,$ (21) $\displaystyle\bar{v}_{\rm E}$
$\displaystyle\simeq 480{\,\rm km\,s^{-1}}\rho_{\rm i,-11}^{-1/3}f_{\rm
ad,0.3}^{-1/3}\biggl{(}\frac{t_{\rm pl}}{100\,\rm
d}\biggl{)}^{-2/3}\biggl{(}\frac{L_{\rm pl}}{10^{39}{\,\rm
erg\,s^{-1}}}\biggl{)}^{1/3}\ ,$ (22)

where we note that our calculations reveal that $f_{\rm ad,0.3}$ increases
moderately $\sim 0.5-2$ with increasing stellar mass and radius (Fig. 9). Eq.
(22) provides an independent velocity estimate from that obtained from
spectroscopy of the transient, $v_{\rm obs}$, the consistency between which
provides a check on the model (see below).

The ejecta launching radius $R_{0}$ can likewise be inferred from the
observables in two ways:
(1) assuming $\bar{v}_{\rm E}=v_{\rm esc}=(2GM_{\star}/R_{0})^{1/2}$ using the
measured $M_{\star}$ when available,

$\displaystyle R_{0}$ $\displaystyle\simeq$ $\displaystyle
5.4\,{\mathrm{R}_{\odot}}\,f_{\rm ad,0.3}^{2/3}\rho_{\rm i,-11}^{2/3}\times$
(23) $\displaystyle\biggl{(}\frac{t_{\rm pl}}{100{\,\rm
d}}\biggl{)}^{4/3}\biggl{(}\frac{L_{\rm pl}}{10^{39}\,\rm
erg\,s^{-1}}\biggl{)}^{-2/3}\biggl{(}\frac{M_{\star}}{3\,{\mathrm{M}_{\odot}}}\biggl{)}\
;$

(2) assuming $v_{\rm obs}=v_{\rm esc}$,

$\displaystyle R_{0}$ $\displaystyle\simeq
13\,{\mathrm{R}_{\odot}}\,\biggl{(}\frac{v_{\rm obs}}{300{\,\rm
km\,s^{-1}}}\biggl{)}^{-2}\biggl{(}\frac{M_{\star}}{3\,{\mathrm{M}_{\odot}}}\biggl{)}\
.$ (24)

Figure 12 compare the inferred ejecta properties $\\{M_{\rm ej},\bar{v}_{\rm
E},R_{0}\\}$ to the stellar properties or transient properties
$\\{M_{\star},R_{\star},v_{\rm obs}\\}$ for events for which the latter
information is available (see also Table 1). At a basic level most of our
results for the baseline model (shown as circles or triangles in Fig. 12) make
physical sense. Inferred ejecta masses range from $\lesssim
10^{-2}\,{\mathrm{M}_{\odot}}$ in V4332 Sgr to $\sim 10\,{\mathrm{M}_{\odot}}$
in AT2018bwo, for example. Likewise, the implied ejecta launching radii are
all safely within the donor’s interior, $R_{0}<R_{\star}$. We return to what
inferences can be drawn about the LRN population based on these data in
Section 4.2.

Upon closer inspection, several inconsistencies or tensions are also
noticeable in Fig. 12. Firstly, the implied ejecta masses from the most
massive binary systems (e.g., NGC4490-OT2011 and AT2017jfs) are enormous
$\gtrsim 100\,{\mathrm{M}_{\odot}}$ and greatly violate the theoretical limit
$M_{\rm ej}>M_{\star}(\eta>1)$ in which the entire mass of the donor is
ejected, while a comparable number of events hover around $\eta\sim 1$ (see
also Soker 2020). A second point of tension is that the inferred ejecta
velocities (estimated from Eq. 22) systematically exceed those observed,
$\bar{v}_{\rm E}\sim(2-3)\,v_{\rm obs}$. Equivalently, the ejection radius
$R_{0}$ one infers by equating $\bar{v}_{\rm E}$ to the local escape velocity
(Eq. 23) is smaller than that implied if the observed spectroscopic expansion
velocity equals the escape speed at the ejection radius (squares in bottom
panel of Fig. 12).

Keeping in mind that our formulae (Eqs. 21 and 22) exhibit errors of a factor
of a few, the latter mentioned velocity/radius discrepancy could in principle
have several origins: (1) the observed spectroscopic velocities underestimate
the true energy-weighted ejecta speed, for instance due how excitation and
optical-depth affect what regions of the ejecta given lines probe; (2)
inaccuracies in our light curve model cause an over-estimate of the plateau
duration, which results in a large velocity for a given observed duration (see
Eq. 16);this could result from one or more of our simplifying assumptions,
e.g. spherical symmetry; (3) our model is missing physics, such as additional
sources of the ejecta heating beyond that from hydrogen recombination. The
latter possibility would also have implications for the inferred ejecta mass.

To illustrate this with an example, consider a hypothetical scenario in which
the ejecta is heated during the plateau phase at a specific rate $\sim 10$
times larger than from hydrogen recombination (e.g., due to shock interaction
with circumbinary material; Section 4.1). In such a case the true ejecta mass
would be $\sim 10$ times smaller than inferred from Eq. (21) from the observed
transient duration and luminosity; hence, for a fixed observed plateau
duration $t_{\rm pl}$ (Eq. 16) the implied ejecta speed $\bar{v}_{\rm
E}\propto M_{\rm ej}^{1/3}$ would be $\sim 10^{1/3}\sim 2$ times smaller than
previously estimated neglecting the additional heating source, bringing it
into closer agreement with $v_{\rm obs}$ (and likewise, the implied $R_{0}$
obtained assuming $\bar{v}_{\rm E}=v_{\rm esc}$ would be $\sim 4$ times
smaller).

To explore this possibility further, the top left panel of Fig. 12 includes a
second estimate of the ejecta mass (shown as squares and diamonds), which
makes use of our expression for the plateau duration (Eq. 16) but supplants
Eq. (21) (which assumes recombination as the sole energy source) with the
assumption that $\bar{v}_{\rm E}=v_{\rm obs}$, i.e. that the ejecta speed
equals that observed from the transient spectroscopy. Comparing the circles
and squares (or triangles and diamonds), we observe that the required value of
$M_{\rm ej}$ decreases significantly, bringing even the most luminous LRNe
near or into the physically allowed region ($\eta\lesssim 1$). Under the same
assumption ($\bar{v}_{\rm E}=v_{\rm obs}$), the implied mass ejection radius
increases but remains below the progenitor star surface for most events
(squares in the bottom panel of Fig. 1).

To summarize, the ejecta masses one infers from the optical light curve
assuming hydrogen recombination as the only energy source (as in our model)
represent only upper limits. If other sources of heating not accounted in our
model are present, then (1) $M_{\rm ej}$ will be lower than this upper limit
and; (2) the ejecta velocities $\bar{v}_{\rm E}$ implied by the observed LRN
duration are brought into better agreement with those observed
spectroscopically.

Table 1: LRNe Sample: Light curve observables, stellar progenitors, and inferred ejecta properties. Event | $L_{90}$ | $t_{90}$ | $E_{\rm rad}$ |  $v_{\rm obs}$‡ | $M_{\star}$ | $R_{\star}$ | $M_{\rm ej}^{\P}$ | $\bar{v}_{\rm E}^{\P}$ | Ref.
---|---|---|---|---|---|---|---|---|---
| [erg/s] | [d] | [erg] | [km/s] | [${\mathrm{M}_{\odot}}$] | [${\mathrm{R}_{\odot}}$] | [${\mathrm{M}_{\odot}}$] | [km/s] |
V4332 Sgr | 1.3e+37 | 15 | 1.7e+43 | $200\pm 20$ | 1 | - | 3.0e-03 | 410 | 1,2†
V838 Mon | 3.6e+39 | 79 | 2.4e+46 | $500\pm 50$ | 5-10 | 4-6 | 4.2e+00 | 860 | 3,4,5†
V1309 Sco | 8.9e+37 | 24 | 1.9e+44 | $150\pm 15$ | 1-3 | 3.5 | 3.2e-02 | 550 | 6,7†
NGC 4490-OT2011 | 2.2e+41 | 180 | 3.3e+48 | $800\pm 100$ | 30 | - | 5.8e+02 | 2000 | 8,9†
M31-OT2015 | 1.5e+39 | 50 | 6.4e+45 | $300\pm 20$ | 3-5.5 | 30-40 | 1.1e+00 | 880 | 10,11,12,13†
M101-OT2015 | 1.4e+40 | 270 | 3.1e+47 | $500\pm 50$ | 18 | 220 | 5.5e+01 | 590 | 14
SN Hunt248 | 6.2e+40 | 130 | 6.9e+47 | $1200\pm 100$ | 30 | 500 | 1.2e+02 | 1600 | 9†,15,16
AT 2014ej | 2.2e+41 | 77 | 1.4e+48 | $900\pm 80$ | - | - | 2.5e+02 | 3400 | 17
AT 2017jfs | 2.8e+41 | 170 | 4.0e+48 | $700\pm 70$ | - | - | 6.9e+02 | 2200 | 18
AT 2018bwo | 1.3e+40 | 60 | 6.5e+46 | $500\pm 65$ | 12-16 | 100 | 1.1e+01 | 1600 | 19
AT 2018hso | 4.5e+40 | 180 | 7.2e+47 | $500\pm 50$ | - | - | 1.2e+02 | 1100 | 20
AT 2019zhd | 7.1e+38 | 33 | 2.0e+45 | $280\pm 30$ | - | - | 3.5e-01 | 900 | 21
AT 2020kog | 6.6e+39 | 110 | 6.4e+46 | $250\pm 30$ | - | - | 1.1e+01 | 830 | 22
AT 2020hat | 3.0e+40 | 110 | 2.9e+47 | $470\pm 50$ | - | - | 5.1e+01 | 1400 | 22
‡ Observed velocities are taken from Table 5 of Blagorodnova et al. (2021) for
all events.
¶Obtained by Eqs. (21) and (22) for $\rho_{\rm i,-11}=f_{\rm ad,0.3}=1$.
† References for the light curve data.
Ref. 1: Martini et al. (1999), 2: Tylenda et al. (2005), 3: Munari et al.
(2002), 4: Bond et al. (2003),
5: Tylenda (2005), 6: Mason et al. (2010), 7: Tylenda et al. (2011), 8: Smith
et al. (2016), 9: Pastorello et al. (2019b)
10: Kurtenkov et al. (2015), 11: Lipunov et al. (2017), 12: Williams et al.
(2015), 13: MacLeod et al. (2017)
14: Blagorodnova et al. (2017), 15: Kankare et al. (2015), 16: Mauerhan et al.
(2015), 17 : Stritzinger et al. (2020)
18 : Pastorello et al. (2019a), 19: Blagorodnova et al. (2021), 20: Cai et al.
(2019), 21: Pastorello et al. (2021a)
22: Pastorello et al. (2021b)
Figure 11: Luminosity $L_{90}$ (left) and duration $t_{90}$ (right) of stellar
merger transients compared to LRNe with progenitor donor masses $M_{\star}$
inferred from pre-transient imaging (Table 1). Red and blue shaded regions
show the predictions of our light curve model as a function of donor mass
$M_{\star}$ for a range of $\eta=M_{\rm ej}/M_{\star}=0.01$ (thinnest lines)
to $\eta=1$ (thickest lines), and where we have assumed $\bar{v}_{\rm
E}=v_{\rm esc,\star}=(2GM_{\star}/R_{\star})^{1/2}$ for stars on the TAMS or
10 times their TAMS radii. In the left panel, a dashed curve denotes the
stellar luminosity at TAMS (Hurley et al., 2000), while a dotted line shows
the luminosity at the tip of the red giant branch (Paczyński 1970; Cummings et
al. 2018). Figure 12: For the same sample of LRNe from Fig. 11, here we show:
(top left panel) Ejecta mass $M_{\rm ej}$ estimated from light curve
observables $\\{L_{90}$, $t_{90}\\}$ versus donor mass from pre-transient
imaging, $M_{\star}$. Circles show $M_{\rm ej}$ calculated from our light
curve model assuming hydrogen recombination as the only energy source (Eq.
21), while squares show $M_{\rm ej}$ calculated from Eq. (16) assuming the
mean ejecta speed equals that observed spectroscopically, $\bar{v}_{\rm
E}=v_{\rm obs}$; triangles and diamonds show the same as circles and squares
but for events without donor-mass constraint (horizontal lines denote the
corresponding mass estimates for an event); for comparison, lines of constant
$\eta=M_{\rm ej}/M_{\star}=0.01,0.1,1$ are shown as solid black lines; (top
right) Ejecta velocity estimated from light curve observables (Eq. 22) versus
spectroscopic velocity, $v_{\rm obs}$ (Circles). Squares show the case when
the estimated velocity is the same exactly as $v_{\rm obs}$.; (bottom left)
Various radii versus donor mass. Circles show the mass ejection radius $R_{0}$
estimated by equating the mean velocity $\bar{v}_{\rm E}$ inferred light curve
modeling to the local escape speed $v_{\rm esc}=(2GM_{\star}/R_{0})^{1/2}$
(Eq. 23), while squares show $R_{0}$ calculated assuming $v_{\rm esc}=v_{\rm
obs}$ (Eq. 24). Stars show the radius of the donor progenitor star,
$R_{\star}$, inferred from pre-transient imaging.

###  3.4 Comparison to SNe II Analytic Formulae

The facts that recombination energy dominates the luminosity and gas pressure
can exceed radiation pressure during the plateau phase (Fig. 6), both violate
implicit assumptions of analytic SNe light curve models (Popov, 1993; Sukhbold
et al., 2016). Because several past works use these formulae to infer the
ejecta properties from observed LRNe (e.g., Ivanova et al. 2013a; MacLeod et
al. 2017; Blagorodnova et al. 2021), or predict observables in a population
analysis (e.g., Howitt et al. 2020), it is important to quantify the implied
inaccuracies.

Figure 13 shows the ratio of the transient luminosity and duration as
extracted from our light curve calculation (Fig. 9) to those predicted by
Popov (1993). For the latter we employ results calibrated using SN II
radiative transfer calculations (Sukhbold et al., 2016; Blagorodnova et al.,
2021),

$\displaystyle L_{\rm Popov}$ $\displaystyle=1.3\times 10^{39}{\,\rm
erg\,s^{-1}\,}\times$ $\displaystyle\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{1/3}\biggl{(}\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\biggl{)}^{2/3}\biggl{(}\frac{\bar{v}_{\rm
E}}{300\,\rm km\,s^{-1}}\biggl{)}^{5/3}\ ,$ (25) $\displaystyle t_{\rm Popov}$
$\displaystyle=51{\,\rm d\,}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{1/3}\biggl{(}\frac{R_{0}}{10\,{\mathrm{R}_{\odot}}}\biggl{)}^{1/6}\biggl{(}\frac{\bar{v}_{\rm
E}}{300\,\rm km\,s^{-1}}\biggl{)}^{-1/3}\ .$ (26)

Remarkably, Fig. 13 shows that the calibrated Popov formulae (Eqs. 25,26)
agree with our model predictions across most of the relevant parameter space
to within a factor $\lesssim 3$. However, this agreement is largely a
coincidence because the models are derived under completely different
assumptions (Popov’s formulae “know” nothing about hydrogen recombination
energy, for instance). The Popov formulae also exhibit different parameter
dependencies from our model (Eqs. 16,18), however they (again, coincidentally)
become similar for $\bar{v}_{\rm E}\sim v_{\rm
esc}\propto(M_{\star}/R_{0})^{1/2}$: $L_{\rm Popov}\propto M_{\rm
ej}^{7/6}R_{0}^{-1/6}$ and $t_{\rm Popov}\propto M_{\rm ej}^{1/6}R_{0}^{1/3}$.

As pointed out by MacLeod et al. (2017); Blagorodnova et al. (2021), using the
Popov (1993) expressions to infer ejecta properties from LRN observables
typically result in large errors on the former. To see why, consider inverting
Eqs. (25), (26) to obtain the ejecta mass and launching radius as functions of
plateau luminosity, duration, and ejecta velocity:

$\displaystyle M_{\rm ej}$ $\displaystyle\simeq
18\,{\mathrm{M}_{\odot}}\,\biggl{(}\frac{t_{\rm pl}}{100\,\rm
d}\biggl{)}^{4}\biggl{(}\frac{\bar{v}_{\rm E}}{300\,\rm
km\,s^{-1}}\biggl{)}^{3}\biggl{(}\frac{L_{\rm pl}}{10^{39}\,\rm
erg\,s^{-1}}\biggl{)}^{-1}\ ,$ (27) $\displaystyle R_{0}$ $\displaystyle\simeq
1.7\,{\mathrm{R}_{\odot}}\,\biggl{(}\frac{t_{\rm pl}}{100\,\rm
d}\biggl{)}^{-2}\biggl{(}\frac{\bar{v}_{\rm E}}{300\,\rm
km\,s^{-1}}\biggl{)}^{-4}\biggl{(}\frac{L_{\rm
pl}}{10^{39}\rm\,erg\,s^{-1}}\biggl{)}^{2}\ .$ (28)

The inferred values $M_{\rm ej},R_{0}$ thus depend sensitively on the
observables, particularly the ejecta velocity. This is problematic because the
measured velocity (e.g., as typically inferred from spectral line widths at
particular epochs) may not coincide precisely with the bulk of the ejecta. By
contrast, the sensitivity of inferred ejecta quantities to LRN observables as
predicted by our model (Eqs. 21,22,23), particularly if the donor star
properties are independently constrained, allow for a more robust extraction
of $\\{M_{\rm ej},\bar{v}_{\rm E},R_{0}\\}$, compared to Eqs. (27,28).

Figure 14 compares the inferences on $\\{M_{\rm ej},R_{0}\\}$ made employing
our model (Table 1) to those obtained by instead applying the modified Popov
(1993) formulae (Eqs. 27,28). We see that the ejecta masses implied by the
Popov formulae can over- or under-estimate those of our model (circles) by
more than 1-2 orders of magnitude. In particular, the Popov formulae
underestimate the implied $M_{\rm ej}$ for the most luminous LRNe (largest
$M_{\rm ej}$ events), masking the $M_{\rm ej}>M_{\star}$ tension implied by
our modeling (Sec. 3.3). In applying the Popov formulae here, we have
estimated the ejecta mass and radius using the observed spectroscopic velocity
$\bar{v}_{\rm E}=v_{\rm obs}$, which typically differs from the surface escape
speed (Fig. 2). The sensitive dependence of the Popov formulae velocity thus
accounts for the larger deviations of the ejecta properties from our estimates
in Fig. 14 than naively suggested by Fig. 13 which instead assumed
$\bar{v}_{\rm E}=v_{\rm esc}$.

Figure 13: The same as Fig. 9 but showing the ratios of the transient
luminosity and duration from our numerical calculations to the analytic
formulae for Type II supernovae (Eqs. 25, 26; Popov 1993; Sukhbold et al.
2016). Along the gray dashed line contour along the radiation pressure is
equal to the gas pressure at the time of ejection for the shell with
$v=\bar{v}_{\rm E}$. The assumption of radiation-dominated ejecta (e.g., Popov
1993) is violated to the left of this line. Figure 14: For the same sample of
LRNe in Fig. 12, here we compare the ejecta masses and launching radii
inferred from our light curve model, $M_{\rm ej}/R_{0}$ on the vertical axis
(using the two methods described in Fig. 12 with the same symbols notation) to
those obtained using the Popov (1993) supernova light curve formulae on the
horizontal axis (Eqs. 27, 28).

##  4 Discussion

###  4.1 Shock Interaction as a Luminosity Source

Hydrogen recombination typically dominates the radiated energy in our baseline
model because the ejecta’s initial thermal energy is mostly lost to adiabatic
expansion before it can be radiated. However, some of this thermal energy can
be regenerated if the merger ejecta were to collide with slower material at
larger radii $\gg R_{\star}$ surrounding the binary (Metzger & Pejcha, 2017).
This pre-existing slower material could arise from an earlier episode of mass-
transfer leading up to the dynamical merger event, for instance a brief
thermal-timescale mass-transfer phase (e.g., Blagorodnova et al. 2021) or
other physical processes generating mass-loss from the binary’s $L_{2}$ point
(Pejcha et al., 2016a, b, 2017). The presence of pre-dynamical ejecta is
supported observationally by the slow early rise in the pre-maximum light
curves of LRNe such as V1309 Sco (Pejcha, 2014; Pejcha et al., 2017) and the
redder than expected early-time colors of V838 Mon and AT2018bwo (Blagorodnova
et al., 2021).

The characteristic mass $M_{\rm pre}$ of the slowly-expanding (velocity
$v_{\rm pre}$) material surrounding the binary at the time of merger is
uncertain, but in many cases is expected to be smaller than the ejecta mass
released at later times, $M_{\rm ej}$ (e.g., at most $M_{\rm pre}\sim
0.1M_{\star}$ is predicted via $L_{2}$ mass-loss; MacLeod et al. 2017; Pejcha
et al. 2017; MacLeod & Loeb 2020; however, see Blagorodnova et al. 2021).
Assuming $M_{\rm pre}\ll M_{\rm ej}$ and $v_{\rm pre}\ll\bar{v}_{\rm E}$, and
that the shock-interaction is momentum conserving, the thermal energy released
by the collision is given by

$E_{\rm sh}\simeq\frac{1}{2}M_{\rm pre}\bar{v}_{\rm E}^{2}\ .$ (29)

This exceeds the energy released by hydrogen recombination, $\varepsilon_{\rm
H}XM_{\rm ej}/m_{\rm p}$, for

$\frac{M_{\rm pre}}{M_{\rm ej}}\gtrsim\frac{2X\varepsilon_{\rm H}}{m_{\rm
p}\bar{v}_{\rm E}^{2}}\sim 0.02\left(\frac{\bar{v}_{\rm E}}{300\,{\rm
km\,s^{-1}}}\right)^{-2}$ (30)

Thus, if even a few percent of the ejecta mass is present on radial scales
$R_{\rm pre}\sim\bar{v}_{\rm E}t_{\rm pl}\sim 10^{14}-10^{15}$ cm prior to the
merger, the shock will compete with recombination in powering the light curve.
A related, alternative heating mechanism for the LRN ejecta is a fast outflow
or jet powered by accretion onto the accretor star or donor star core (e.g.,
Soker & Kashi 2016; Soker 2020; Soker & Kaplan 2021).

If the duration of the transient is still set by hydrogen recombination in
this scenario, then the transient duration is still given by $t_{90}\sim
t_{\rm pl}\propto M_{\rm ej}^{1/3}\bar{v}_{\rm E}^{-1}$ (Eq. 16). However,
instead of Eq. (17) for the plateau luminosity, the shock scenario predicts
$E_{\rm pl}\propto E_{\rm sh}\propto\bar{v}_{\rm E}^{2}$ and hence

$\displaystyle L_{\rm pl}^{\rm sh}$ $\displaystyle\sim\frac{E_{\rm sh}}{t_{\rm
pl}}\simeq 7\times 10^{39}{\rm erg\,s^{-1}}\rho_{\rm i,-11}^{1/3}\times$
$\displaystyle\left(\frac{M_{\rm pre}}{0.1M_{\rm
ej}}\right)\left(\frac{\bar{v}_{\rm E}}{300\,{\rm
km\,s^{-1}}}\right)^{3}\left(\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\right)^{2/3}.$ (31)

The value of $L_{\rm pl}^{\rm sh}$ can now reach $\gtrsim 10^{41}$ erg s-1 for
events like AT2017jfs which possess a relatively high velocity $v_{\rm
obs}\sim 700$ km s${}^{-1},$ and which otherwise required unphysically large
ejecta masses when considering recombination energy alone (Fig. 12). Depending
on its opacity $\kappa$, the optical depth through the circumbinary material,

$\displaystyle\tau_{\rm pre}$ $\displaystyle\sim$ $\displaystyle\frac{\kappa
M_{\rm pre}}{4\pi R_{\rm pre}^{2}}$ (32) $\displaystyle\sim$ $\displaystyle
1.6\left(\frac{\kappa}{0.1\,{\rm cm^{2}\,g^{-1}}}\right)\left(\frac{M_{\rm
pre}}{0.1\,{\mathrm{M}_{\odot}}}\right)\left(\frac{R_{\rm pre}}{10^{15}\,{\rm
cm}}\right)^{-2}\ ,$

is sufficiently large to redden the early transient emission (e.g.,
Blagorodnova et al. 2021). Future work is required to incorporate the time-
dependent shock luminosity, given an assumed density profile of the pre-merger
cirumbinary medium (e.g., Metzger & Pejcha 2017), self-consistently into our
light curve model framework (Section 2).

One might also expect shock interaction to manifest through observational
signatures other than just the light curve, such as optical emission lines and
thermal X-rays from hot shocked gas and non-thermal radio synchrotron emission
from relativistic electrons accelerated at the shock (as seen from internal
shocks in classical novae, for example; Chomiuk et al. 2021). However, the
equatorial geometry of the slow circumbinary medium (e.g., Pejcha et al.
2016a) implies that the polar regions of the merger ejecta and hence the
optical photosphere will quickly “wrap around” the equatorial shock, blocking
signatures of the latter from most viewing angles (see Andrews & Smith 2018
for an example of this in the context of shock-powered Type II supernovae). At
late times after recombination completes, the opacity will drop, the
photosphere recedes and residual shock signatures may become visible; however,
the high optical-wavelength opacity once dust formation begins around the same
time may keep these inner regions opaque and continue to obscure shock-related
signatures.

###  4.2 Inferences about the LRN Population

A goal of measuring the ejecta properties for a sample of LRNe is to deduce
information about the physics of the merger process as well as the
demographics of merging binary systems.

Putting aside the ambiguities that arise in measuring individual ejecta masses
in the face of additional sources of ejecta heating, the top panel of Fig. 12
reveals an interesting trend that is likely to remain robust in spite of these
uncertainties: the fraction of the stellar envelope being ejected $\eta=M_{\rm
ej}/M_{\star}$ systematically increases with more massive donor stars. This
arises because the average luminosity of LRNe increases sharply with stellar
mass, $L\propto M_{\star}^{\alpha}$ with $\alpha\sim 2-3$ (e.g., Kochanek et
al. 2014), while recombination powered models predict a sub-linear dependence
of $L_{\rm pl}\propto M_{\rm ej}^{2/3}$ (Eq. 18). What could give rise to the
increasing $\eta(M_{\star})$ trend?

The fraction of the donor’s star envelope unbound during the merger process is
likely to increase with secondary mass and hence the mass ratio of the binary
(e.g., Nordhaus et al. 2010). The fact that the average mass ratio of
interacting binaries is expected to increase with stellar mass (Moe & Di
Stefano, 2017) could therefore contribute to $\eta$ increasing with
$M_{\star}$. Another contributing factor is that more massive LRN progenitors
appear to be further evolved off the main sequence (Fig. 2): more evolved
stars possess a more pronounced core-envelope structure, facilitating the
ejection of a higher fraction of the envelope as the accretor releases more
gravitational energy spiraling into the compact core. A smaller ejection
radius also increases the ejecta velocity, which acts to increase the
transient luminosity (Eq. 18). Indeed, the most luminous events tend to have
larger observed velocities (see Table 1).

The current LRN sample could also be biased by a selection effect: since the
intrinsic nuclear luminosity of the star at all phases of evolution increases
steeply with stellar mass (dotted and dashed lines in Fig. 11), low-$M_{\rm
ej}$ (low-$\eta$) events could be selected against in their discovery because
of the lower contrast between the transient and star luminosity. Moreover,
since massive progenitors are rare, their outbursts are more typically
detected as extragalactic transients; hence, it is natural to miss low-
luminosity events from massive stars. These would imply the existence of a
hidden population of low-$M_{\rm ej}$ mergers (“minor mergers”) from massive
stars. This hidden population would need to be taken into account when
comparing the $M_{\star}$-dependent rate of LRN to the initial mass function
of stars or to rate predictions from binary population synthesis modeling
(e.g., Kochanek et al. 2014; Howitt et al. 2020).

###  4.3 Stellar Mergers with Planets and Compact Objects

Our model can be applied to predict the optical light curves from other types
of mergers involving mass ejection from hydrogen-rich bodies.

One example is the merger of a giant planet or brown dwarf with a more massive
host companion star (e.g., Metzger et al. 2012). When the average density of
the planet is larger than that of the star, the energy released as the planet
plunges slowly below the stellar surface can give rise to mass ejection (and a
longer-term evolution of the star’s luminosity due to energy released deeper
in the stellar envelope; Metzger et al. 2017; MacLeod et al. 2018a).
Alternatively, if the planet is less dense than the star (e.g., “inflated” Hot
Jupiters; Baraffe et al. 2010), then it can be tidally disrupted above the
surface, forming an accretion disk around the star. Outflows from the
resulting super-Eddington accretion flow are likely to unbind a significant
fraction of its mass from the planet (e.g., Metzger et al. 2012). In both
cases then, a fraction of the planet/brown dwarf mass $\sim
10^{-3}-10^{-2}\,{\mathrm{M}_{\odot}}$ is ejected from radii $R_{0}\sim
R_{\star}$ at close to the surface escape speed.

Figure 11 shows that for $M_{\star}\sim{\mathrm{M}_{\odot}}$ and $\eta\sim
10^{-3}-10^{-2}$, the predicted transient lasts days to weeks with a
luminosity $L\sim 10^{37}-10^{38}$ erg s-1. These light curve properties are
comparable to those of classical novae (e.g., Strope et al. 2010), with which
planet-star mergers could easily be confused. Based on the decay rate of
observed gas giant planets orbits due to tidal dissipation, Metzger et al.
(2012) estimate a planet-star merger rate of $\sim 0.1-1$ per year in the
Milky Way, consistent (within an order of magnitude) of the stellar merger
rate (Kochanek et al., 2014) and $\sim 100$ times lower than the rate of
classical novae (e.g., Kawash et al. 2021). Given their low luminosities, and
low relative frequency to classical novae, it is unsurprising that no
confirmed planet-star mergers have yet been discovered.

Our model can also be applied to mergers of hydrogen-rich stars with compact
objects, such as white dwarfs (WD; e.g., Shara & Shaviv 1977; Michaely & Shara
2021; Metzger et al. 2021), neutron stars (NS), or stellar- or intermediate-
mass black holes (BH; e.g., Perets et al. 2016; Fragione et al. 2020; Kremer
et al. 2021, 2022), which give rise to up to several solar masses of ejecta
expanding at velocities of hundreds to thousands of km s-1. Due to the deeper
potential well of the WD/NS/BH, the kinetic energies of the ejecta from these
events can exceed those of ordinary star mergers; however, the small ejection
radii$-$less than the tidal disruption radius, itself typically comparable to
the stellar radius$-$imply that most of the initial thermal energy will be
lost to adiabatic expansion, as in the stellar mergers we have studied. As a
result, the total radiated optical energy from the unbound ejecta will also be
limited to that from its hydrogen recombination, i.e. $\lesssim f_{\rm
ad}\varepsilon_{\rm H}X{\mathrm{M}_{\odot}}/m_{\rm p}\sim 10^{46}$ erg per
solar mass of ejecta. The predicted luminosity from star-WD/NS/BH mergers are
thus unlikely to exceed $\sim 10^{39}-10^{40}$ erg s-1 (Metzger et al., 2021),
less than estimates that do not account for adiabatic losses prior to
recombination (e.g., Michaely & Shara 2021). Of course, if as may be the case
also in stellar mergers (Sec. 4.1), star-WD/NS/BH mergers are accompanied by
an extra central luminosity source which re-energizes the unbound ejecta
(e.g., internal shocks due to ongoing outflows from the central engine; e.g.,
Dexter & Kasen 2013; Margalit & Metzger 2016; Gottlieb et al. 2022) then the
transient luminosities could be considerably higher than from recombination
power alone.

##  5 Conclusions

We have presented a new model for the light curves of LRNe from stellar
mergers, which self-consistently accounts for the range of relevant opacities,
an arbitrary degree of radiation and gas pressure, and the energy released
from hydrogen recombination. We apply our model to a sample of 14 LRNe (Fig.
1, Table 1) to infer the ejecta properties (mass, velocity, and launching
radius) and compare to the mass/radios of the progenitor donor star when
available (Fig. 2). Our conclusions can be summarized as follows:

  * •

Our model naturally reproduces the observation (Fig. 1) that LRN light curves
typically possess two luminosity peaks: (1) a short initial peak typically
lasting days powered by the initial thermal energy (“cooling emission”) from
the outermost fastest layers of the ejecta; and (2) a longer peak or plateau
lasting weeks or months, powered by energy released from hydrogen
recombination from the bulk of the mass. While radiation pressure can dominate
in the ejecta during the first peak, gas pressure typically dominates the
latter plateau, at odds with the common assumption of supernova light curve
models (e.g., Popov 1993).

The potential for two light curve peaks was recognized previously as a result
of distinct physical stages during the merger (MacLeod et al., 2017; Metzger &
Pejcha, 2017). Here, we emphasize (again) that such a double-peaked feature
results due to the stratification of even a single spherically symmetric
ejecta shell into layers with distinct velocities and initial thermal energy
content. Multiple distinct mass-loss episodes (e.g., Ivanova et al. 2013a;
Pastorello et al. 2019a) are not required to generate multiple light curve
peaks.

  * •

The luminosity and duration of the early and late light curve peaks follow
from simple analytic estimates (Eqs. 11,12 and Eqs. 16,18, respectively).
These formulae reproduce within a factor of a few the average luminosity
$L_{90}$ and duration $t_{90}$ (over which 90% of the total optical energy is
received) across a wide parameter space, corresponding to mergers spanning the
main sequence and post-main sequence phases of stars of mass $\sim
0.1\,{\mathrm{M}_{\odot}}$ to $\sim 100\,{\mathrm{M}_{\odot}}$ (Figs. 10,11).
A one-zone calculation, similar in methodology to Popov (1993) but assuming
gas pressure-dominated conditions and accounting for recombination energy,
reproduces similar scalings (Appendix A).

The total radiated energy $E_{\rm rad}$ is dominated by the recombination-
powered plateau; hence the measured value of $E_{\rm rad}$ is proportional to
the ejecta mass (up to an adiabatic loss factor $f_{\rm ad}\sim 0.1$; Fig. 9).

  * •

We determine the ejecta properties $\\{M_{\rm ej},\bar{v}_{\rm E},R_{0}\\}$ of
individual LRN events, first under the assumption that hydrogen recombination
is the sole source of ejecta heating after launch, and compare them to the
masses and radii of observed stellar progenitors of LRNe (Figs. 11, 12, Table
1). The derived ejecta properties differ from those derived by applying SN II
light curve models such as Popov (1993) (e.g., Figs. 13, 14). The latter also
come with larger uncertainties due to the greater sensitivity of the formulae
to the observables (Eqs. 26,25). Remarkably, our results roughly agree with
those by Popov (1993) under the assumption that the ejecta velocity equals the
stellar escape velocity (Figs. 13). However, such information is not typically
available in the Popov framework, and instead using the observed spectroscopic
velocity results in drastically different ejecta property estimates compared
to those obtained applying our model (Fig. 14).

  * •

In broad terms, for the range of expected ejecta masses and velocities of
stars at different stages of stellar evolution, we predict merger light curve
durations $\sim$ weeks to months and luminosities $\sim 10^{38}-10^{40}$ erg
s-1, consistent with observed LRNe, as well as a general trend of increasing
luminosity with increasing stellar mass (Fig. 11).

  * •

However, several LRNe arising from mergers involving the most massive stars
$\gtrsim 10\,{\mathrm{M}_{\odot}}$, exceed the luminosity upper limit
corresponding to ejection of the entire donor star envelope ($M_{\rm
ej}/M_{\star}=\eta>1$; circles in Figs. 11 and 12). The inferred mean ejecta
speed $\bar{v}_{\rm E}$ of many events also exceed those inferred by transient
spectroscopy by a factor of $\sim 3$.

Both the luminosity and velocity tension can be alleviated by postulating the
existence of an additional energy source for the ejecta exceeding that from
hydrogen recombination (squares in Fig. 12). This reduces the required ejecta
mass to explain the radiated energy as well as the ejecta velocity needed to
explain the light curve duration $t_{\rm pl}\propto M_{\rm
ej}^{1/3}\bar{v}_{\rm E}^{-1}.$

  * •

A promising source for the required “extra” heating is shock-interaction with
circumbinary material generated prior to the dynamical phase of the merger
(e.g., Metzger & Pejcha 2017). The strong velocity dependence of shock-
generated heating $\propto v^{2}$ and its associated plateau luminosity (Eq.
31), coupled with the observation that the ejecta velocities of LRNe increase
with stellar mass (Fig. 2), would together act to steepen the $L-M_{\star}$
relationship relative to the fiducial recombination-only prediction (Fig. 11),
bringing it closer into alignment with LRN observations which find $L\propto
M_{\star}^{2-3}$ (Kochanek et al. 2014; Fig. 11).

LRN donor star progenitors with larger masses appear systematically more
evolved off the main sequence and generate merger ejecta with velocities
exceeding the surface escape speed (Fig. 2); this could arise due to the
better developed core-envelope structure of evolved stars, which enables the
secondary to inspiral deeper into the gravitational potential well.

  * •

Although the robust conclusions we are able to draw are limited by the modest
statistics afforded by the current sample of stellar merger transients, this
population is expected to grow considerably in the years ahead. Howitt et al.
(2020) predict that the Vera Rubin Observatory will discover $\sim 500$ LRNe
per year in nearby galaxies, for which pre-transient imaging will enable
progenitor star constraints, allowing the methodology presented here to be
greatly expanded.

  * •

Our model can also be generalized to other merger events involving the
ejection of hydrogen-rich material, such as planet-star mergers and star-
compact object mergers. However, the relatively modest luminosities achieved
by hydrogen recombination power alone, coupled with the comparatively low
rates of these events, will likely render these transient signals challenging
to detect.

We are in debt to Nadia Blagorodnova for providing useful guidance on
assembling LRN light curve data. We thank Ondrej Pejcha for helpful
conversations and comments on an early draft of the text. This work is
supported in part by JSPS Overseas Research Fellowships (T.M.). B.D.M.
acknowledges support from the National Science Foundation (grant number
AST-2009255).

##  Appendix A Analytical light curve model

We present a one-zone model (Dexter & Kasen, 2013; Matsumoto et al., 2016) for
recombination powered transients analogous to Popov (1993) but for which the
ejecta pressure is dominated by gas instead of radiation and accounting for
the energy released by recombination. We again consider a spherical and
homologously expanding ejecta. The initial ejecta temperature is roughly given
by equating the kinetic energy $M_{\rm ej}v^{2}/2$ with internal energy
$3M_{\rm ej}k_{\rm B}T/2m_{\rm p}$:

$\displaystyle T_{0}\simeq\frac{m_{\rm p}v^{2}}{3k_{\rm B}}\simeq 3.6\times
10^{6}{\,\rm K\,}\biggl{(}\frac{v}{300\,\rm km\,s^{-1}}\biggl{)}^{2}\ .$ (A1)

The ejecta is dominated by gas pressure for initial densities larger than a
critical value,

$\displaystyle\rho_{0}\simeq 1.4{\,\rm g\,cm^{-3}\,}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}\biggl{(}\frac{R_{0}}{{\mathrm{R}_{\odot}}}\biggl{)}^{-3}\gtrsim
1.5\times 10^{-3}{\,\rm g\,cm^{-3}\,}\biggl{(}\frac{v}{300{\,\rm
km\,s^{-1}}}\biggl{)}^{6}\ .$ (A2)

The time evolution of temperature is determined by the first law of
thermodynamics:

$\displaystyle\frac{dE}{dt}=-P\frac{dV}{dt}-L\ ,$ (A3)

where $E$, $P$, $V$, and $L$ are the internal energy, pressure, volume, and
luminosity, respectively. Since the gas-pressure dominates the radiation
pressure, the internal energy and pressure are given by $E\sim 3M_{\rm
ej}k_{\rm B}T/2m_{\rm p}$ and $P=(\gamma_{3}-1)E/V$ with the adiabatic index
$\gamma_{3}=5/3$ as long as the recombination is negligible. The radiative
luminosity is given by the diffusion approximation:

$\displaystyle L\simeq 4\pi R^{2}\frac{c}{3\kappa\rho}\frac{d(aT^{4})}{dR}\sim
4\pi R^{2}\frac{acT^{4}}{3\kappa\rho R}\ ,$ (A4)

where $\rho=M_{\rm ej}/V$ is the density. Combining Eq. (A3) with (A4), we
have

$\displaystyle\frac{dT}{dR}$ $\displaystyle=-\frac{2T}{R}-{\cal A}R^{4}T^{4}\
,$ (A5)

where ${\cal A}={32\pi^{2}acm_{\rm p}}/({27k_{\rm B}\kappa M_{\rm ej}^{2}}v)$.
Assuming the opacity is constant, we integrate this equation to obtain

$\displaystyle\frac{T}{T_{0}}$
$\displaystyle=\biggl{[}1+\xi_{0}\biggl{(}1-\frac{R_{0}}{R}\biggl{)}\biggl{]}^{-1/3}\biggl{(}\frac{R}{R_{0}}\biggl{)}^{-2}\
,$ (A6) $\displaystyle\xi_{0}$ $\displaystyle=3{\cal
A}R_{0}^{5}T_{0}^{3}\simeq 2.0\times
10^{-10}\,\kappa_{0.32}^{-1}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{-2}\biggl{(}\frac{R_{0}}{{\mathrm{R}_{\odot}}}\biggl{)}^{5}\biggl{(}\frac{v}{300{\,\rm
km\,s^{-1}}}\biggl{)}^{5}\ ,$ (A7)

where $\kappa=0.32\,\kappa_{0.32}\,\rm cm^{2}\,g^{-1}$. For a typical gas-
pressure-dominated ejecta, $\xi_{0}\ll 1$ and the temperature and luminosity
decrease as $T\propto t^{-2}$ and $L\propto t^{-4}$, respectively.

When the temperature drops to a value $T_{\rm i}\sim 10^{4}$ K, recombination
begins in the outermost layers of the ejecta. As time proceeds, the
recombination layer recedes back into the ejecta shell. We denote the location
of the recombination front by a dimensionless radial coordinate $x\equiv
R_{\rm i}/R$ and assume that the gas temperature is fixed to $T_{\rm i}$ at
the layer. Considering now only the volume inside the recombination front
$4\pi(xR)^{3}/3$, the diffusion luminosity (Eq. A4) becomes

$\displaystyle L\simeq
4\pi(xR)^{2}\frac{c}{3\kappa\rho}\frac{d(aT^{4})}{dR}\sim
4\pi(xR)^{2}\frac{acT_{\rm i}^{4}}{3\kappa\rho xR}\ .$ (A8)

Then Eq. (A3) gives the time evolution of the recombination front $x$:

$\displaystyle\frac{dx}{dR}=-\frac{\lambda x}{R}-\frac{(1-\lambda){\cal
A}T_{\rm i}^{3}R^{4}}{3x}\ ,$ (A9)

where we define $\lambda\equiv(\gamma_{3}-1)/\gamma_{3}$. Assuming $\lambda$
and $\kappa$ are constant during the recombination phase, we can integrate
this equation analytically to obtain

$\displaystyle x$ $\displaystyle=\biggl{[}1-\frac{2(1-\lambda)\xi_{\rm
i}}{9(2\lambda+5)}\biggl{\\{}\biggl{(}\frac{R}{R_{\rm
i}}\biggl{)}^{2\lambda+5}-1\biggl{\\}}\biggl{]}^{1/2}\biggl{(}\frac{R}{R_{\rm
i}}\biggl{)}^{-\lambda}\ ,$ (A10) $\displaystyle\xi_{\rm i}$
$\displaystyle=3{\cal A}R_{\rm i}^{5}T_{\rm i}^{3}\simeq 1.0\times
10^{-11}\,\kappa_{0.32}^{-1}T_{\rm i,4}^{1/2}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{-2}\biggl{(}\frac{R_{0}}{{\mathrm{R}_{\odot}}}\biggl{)}^{5}\biggl{(}\frac{v}{300{\,\rm
km\,s^{-1}}}\biggl{)}^{4}\ ,$ (A11)

where $T_{\rm i}=10^{4}\,T_{\rm i,4}\,\rm K$ and the ejecta radius at the
beginning of the recombination is given by Eq. (A6),

$\displaystyle R_{\rm i}\simeq R_{0}\biggl{(}\frac{T_{0}}{T_{\rm
i}}\biggl{)}^{1/2}\simeq 19\,R_{0}\,T_{\rm
i,4}^{-1/2}\biggl{(}\frac{v}{300{\,\rm km\,s^{-1}}}\biggl{)}\ .$ (A12)

During the recombination phase, the light curve shows a plateau similar to
that predicted by Popov (1993), with its duration given by setting $x=0$:

$\displaystyle t_{\rm pl}=\frac{R_{\rm
i}}{v}\biggl{(}\frac{9(2\lambda+5)}{2(1-\lambda)\xi_{\rm
i}}\biggl{)}^{\frac{1}{2\lambda+5}}\overset{\lambda=0}{\simeq}150{\,\rm
d\,}\kappa_{0.32}^{1/5}T_{\rm i,4}^{-3/5}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{2/5}\biggl{(}\frac{v}{300{\,\rm
km\,s^{-1}}}\biggl{)}^{-4/5}\ .$ (A13)

In the last equality, we have set $\lambda=0$ because the temperature does not
evolve significantly during the plateau. A characteristic luminosity may be
given by the peak of the plateau, which is obtained by setting $dL/dt=0$,

$\displaystyle L_{\rm pl}=\frac{16\pi^{2}acT_{\rm i}^{4}R_{\rm i}^{4}}{9\kappa
M_{\rm
ej}}\biggl{(}\frac{2\lambda+5}{13}\biggl{)}^{1/2}\biggl{(}\frac{9(4-\lambda)(2\lambda+5)}{13(1-\lambda)\xi_{\rm
i}}\biggl{)}^{\frac{4-\lambda}{2\lambda+5}}\overset{\lambda=0}{\simeq}6.0\times
10^{38}{\,\rm erg\,s^{-1}\,}\kappa_{0.32}^{-1/5}T_{\rm
i,4}^{8/5}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}^{3/5}\biggl{(}\frac{v}{300\,\rm
km\,s^{-1}}\biggl{)}^{4/5}\ .$ (A14)

The total radiated energy is then given by

$\displaystyle E_{\rm pl}\sim L_{\rm pl}t_{\rm
pl}\overset{\lambda=0}{\simeq}4.4\,\frac{k_{\rm B}T_{\rm i}M_{\rm ej}}{m_{\rm
p}}\simeq 7.8\times 10^{45}{\,\rm erg\,}T_{\rm i,4}\biggl{(}\frac{M_{\rm
ej}}{{\mathrm{M}_{\odot}}}\biggl{)}\ .$ (A15)

For $\lambda=0$ ($\gamma_{3}=1$) and $T_{\rm i,4}\sim 1$ we see that the
plateau energy matches that from hydrogen recombination (Eq. 17) and the
parameter dependencies for $L_{\rm pl}\propto M_{\rm ej}^{3/5}v^{4/5}$ (Eq.
A14), $t_{\rm pl}\propto M_{\rm ej}^{2/5}v^{-4/5}$ (Eq. A13) and their
normalizations, are similar to those we derived in Eqs. (16)-(18) under
different assumptions in Sec. 3.1, namely $L_{\rm pl}\propto M_{\rm
ej}^{2/3}v$ and $t_{\rm pl}\propto M_{\rm ej}^{1/3}v^{-1}$.

## References

  * Andrews & Smith (2018) Andrews, J. E., & Smith, N. 2018, MNRAS, 477, 74, doi: 10.1093/mnras/sty584
  * Arnett (1980) Arnett, W. D. 1980, ApJ, 237, 541, doi: 10.1086/157898
  * Banerjee et al. (2015) Banerjee, D. P. K., Nuth, Joseph A., I., Misselt, K. A., et al. 2015, ApJ, 814, 109, doi: 10.1088/0004-637X/814/2/109
  * Baraffe et al. (2010) Baraffe, I., Chabrier, G., & Barman, T. 2010, Reports on Progress in Physics, 73, 016901, doi: 10.1088/0034-4885/73/1/016901
  * Belczynski et al. (2018) Belczynski, K., Askar, A., Arca-Sedda, M., et al. 2018, A&A, 615, A91, doi: 10.1051/0004-6361/201732428
  * Blagorodnova et al. (2017) Blagorodnova, N., Kotak, R., Polshaw, J., et al. 2017, ApJ, 834, 107, doi: 10.3847/1538-4357/834/2/107
  * Blagorodnova et al. (2021) Blagorodnova, N., Klencki, J., Pejcha, O., et al. 2021, A&A, 653, A134, doi: 10.1051/0004-6361/202140525
  * Bond (2019) Bond, H. E. 2019, ApJ, 887, 12, doi: 10.3847/1538-4357/ab4e13
  * Bond et al. (2003) Bond, H. E., Henden, A., Levay, Z. G., et al. 2003, Nature, 422, 405, doi: 10.1038/nature01508
  * Cai et al. (2019) Cai, Y. Z., Pastorello, A., Fraser, M., et al. 2019, A&A, 632, L6, doi: 10.1051/0004-6361/201936749
  * Chamandy et al. (2019) Chamandy, L., Tu, Y., Blackman, E. G., et al. 2019, MNRAS, 486, 1070, doi: 10.1093/mnras/stz887
  * Chevalier (2012) Chevalier, R. A. 2012, ApJ, 752, L2, doi: 10.1088/2041-8205/752/1/L2
  * Choi et al. (2016) Choi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102, doi: 10.3847/0004-637X/823/2/102
  * Chomiuk et al. (2021) Chomiuk, L., Metzger, B. D., & Shen, K. J. 2021, ARA&A, 59, doi: 10.1146/annurev-astro-112420-114502
  * Clayton et al. (2017) Clayton, M., Podsiadlowski, P., Ivanova, N., & Justham, S. 2017, MNRAS, 470, 1788, doi: 10.1093/mnras/stx1290
  * Cummings et al. (2018) Cummings, J. D., Kalirai, J. S., Tremblay, P. E., Ramirez-Ruiz, E., & Choi, J. 2018, ApJ, 866, 21, doi: 10.3847/1538-4357/aadfd6
  * Davies et al. (2004) Davies, M. B., Piotto, G., & de Angeli, F. 2004, MNRAS, 349, 129, doi: 10.1111/j.1365-2966.2004.07474.x
  * De Marco & Izzard (2017) De Marco, O., & Izzard, R. G. 2017, PASA, 34, e001, doi: 10.1017/pasa.2016.52
  * de Mink et al. (2014) de Mink, S. E., Sana, H., Langer, N., Izzard, R. G., & Schneider, F. R. N. 2014, ApJ, 782, 7, doi: 10.1088/0004-637X/782/1/7
  * Dexter & Kasen (2013) Dexter, J., & Kasen, D. 2013, ApJ, 772, 30, doi: 10.1088/0004-637X/772/1/30
  * Dotter (2016) Dotter, A. 2016, ApJS, 222, 8, doi: 10.3847/0067-0049/222/1/8
  * Faran et al. (2019) Faran, T., Goldfriend, T., Nakar, E., & Sari, R. 2019, ApJ, 879, 20, doi: 10.3847/1538-4357/ab218a
  * Ferreira et al. (2019) Ferreira, T., Saito, R. K., Minniti, D., et al. 2019, MNRAS, 486, 1220, doi: 10.1093/mnras/stz878
  * Fragione et al. (2020) Fragione, G., Metzger, B. D., Perna, R., Leigh, N. W. C., & Kocsis, B. 2020, MNRAS, 495, 1061, doi: 10.1093/mnras/staa1192
  * Ge et al. (2020) Ge, H., Webbink, R. F., & Han, Z. 2020, ApJS, 249, 9, doi: 10.3847/1538-4365/ab98f6
  * Glebbeek et al. (2013) Glebbeek, E., Gaburov, E., Portegies Zwart, S., & Pols, O. R. 2013, MNRAS, 434, 3497, doi: 10.1093/mnras/stt1268
  * Gottlieb et al. (2022) Gottlieb, O., Tchekhovskoy, A., & Margutti, R. 2022, arXiv e-prints, arXiv:2201.04636.  https://arxiv.org/abs/2201.04636
  * Howitt et al. (2020) Howitt, G., Stevenson, S., Vigna-Gómez, A., et al. 2020, MNRAS, 492, 3229, doi: 10.1093/mnras/stz3542
  * Hurley et al. (2000) Hurley, J. R., Pols, O. R., & Tout, C. A. 2000, MNRAS, 315, 543, doi: 10.1046/j.1365-8711.2000.03426.x
  * Iben & Livio (1993) Iben, Icko, J., & Livio, M. 1993, PASP, 105, 1373, doi: 10.1086/133321
  * Ivanova et al. (2013a) Ivanova, N., Justham, S., Avendano Nandez, J. L., & Lombardi, J. C. 2013a, Science, 339, 433, doi: 10.1126/science.1225540
  * Ivanova & Nandez (2016) Ivanova, N., & Nandez, J. L. A. 2016, MNRAS, 462, 362, doi: 10.1093/mnras/stw1676
  * Ivanova et al. (2013b) Ivanova, N., Justham, S., Chen, X., et al. 2013b, A&A Rev., 21, 59, doi: 10.1007/s00159-013-0059-2
  * Kankare et al. (2015) Kankare, E., Kotak, R., Pastorello, A., et al. 2015, A&A, 581, L4, doi: 10.1051/0004-6361/201526631
  * Kasen & Ramirez-Ruiz (2010) Kasen, D., & Ramirez-Ruiz, E. 2010, ApJ, 714, 155, doi: 10.1088/0004-637X/714/1/155
  * Kasliwal et al. (2017) Kasliwal, M. M., Bally, J., Masci, F., et al. 2017, ApJ, 839, 88, doi: 10.3847/1538-4357/aa6978
  * Kawash et al. (2021) Kawash, A., Chomiuk, L., Rodriguez, J. A., et al. 2021, ApJ, 922, 25, doi: 10.3847/1538-4357/ac1f1a
  * Klencki et al. (2021) Klencki, J., Nelemans, G., Istrate, A. G., & Chruslinska, M. 2021, A&A, 645, A54, doi: 10.1051/0004-6361/202038707
  * Kochanek et al. (2014) Kochanek, C. S., Adams, S. M., & Belczynski, K. 2014, MNRAS, 443, 1319, doi: 10.1093/mnras/stu1226
  * Kramer et al. (2020) Kramer, M., Schneider, F. R. N., Ohlmann, S. T., et al. 2020, A&A, 642, A97, doi: 10.1051/0004-6361/202038702
  * Kremer et al. (2022) Kremer, K., Lombardi, James C., J., Lu, W., Piro, A. L., & Rasio, F. A. 2022, arXiv e-prints, arXiv:2201.12368.  https://arxiv.org/abs/2201.12368
  * Kremer et al. (2021) Kremer, K., Lu, W., Piro, A. L., et al. 2021, ApJ, 911, 104, doi: 10.3847/1538-4357/abeb14
  * Krishna Swamy (1961) Krishna Swamy, K. S. 1961, ApJ, 134, 1017, doi: 10.1086/147235
  * Kruckow et al. (2016) Kruckow, M. U., Tauris, T. M., Langer, N., et al. 2016, A&A, 596, A58, doi: 10.1051/0004-6361/201629420
  * Kulkarni et al. (2007) Kulkarni, S. R., Ofek, E. O., Rau, A., et al. 2007, Nature, 447, 458, doi: 10.1038/nature05822
  * Kurtenkov et al. (2015) Kurtenkov, A. A., Pessev, P., Tomov, T., Barsukova, E. A., et al. 2015, A&A, 578, L10, doi: 10.1051/0004-6361/201526564
  * Law-Smith et al. (2020) Law-Smith, J. A. P., Everson, R. W., Ramirez-Ruiz, E., et al. 2020, arXiv e-prints, arXiv:2011.06630.  https://arxiv.org/abs/2011.06630
  * Lipunov et al. (2017) Lipunov, V. M., Blinnikov, S., et al. 2017, MNRAS, 470, 2339, doi: 10.1093/mnras/stx1107
  * Lombardi et al. (2002) Lombardi, James C., J., Warren, J. S., Rasio, F. A., Sills, A., & Warren, A. R. 2002, ApJ, 568, 939, doi: 10.1086/339060
  * MacLeod et al. (2018a) MacLeod, M., Cantiello, M., & Soares-Furtado, M. 2018a, ApJ, 853, L1, doi: 10.3847/2041-8213/aaa5fa
  * MacLeod & Loeb (2020) MacLeod, M., & Loeb, A. 2020, ApJ, 895, 29, doi: 10.3847/1538-4357/ab89b6
  * MacLeod et al. (2017) MacLeod, M., Macias, P., Ramirez-Ruiz, E., et al. 2017, ApJ, 835, 282, doi: 10.3847/1538-4357/835/2/282
  * MacLeod et al. (2018b) MacLeod, M., Ostriker, E. C., & Stone, J. M. 2018b, ApJ, 868, 136, doi: 10.3847/1538-4357/aae9eb
  * Margalit & Metzger (2016) Margalit, B., & Metzger, B. D. 2016, MNRAS, 461, 1154, doi: 10.1093/mnras/stw1410
  * Martini et al. (1999) Martini, P., Wagner, R. M., Tomaney, A., et al. 1999, AJ, 118, 1034, doi: 10.1086/300951
  * Mason et al. (2010) Mason, E., Diaz, M., Williams, R. E., Preston, G., & Bensby, T. 2010, A&A, 516, A108, doi: 10.1051/0004-6361/200913610
  * Matsumoto et al. (2016) Matsumoto, T., Nakauchi, D., Ioka, K., & Nakamura, T. 2016, ApJ, 823, 83, doi: 10.3847/0004-637X/823/2/83
  * Mauerhan et al. (2015) Mauerhan, J. C., Van Dyk, S. D., Graham, M. L., et al. 2015, MNRAS, 447, 1922, doi: 10.1093/mnras/stu2541
  * Metzger et al. (2012) Metzger, B. D., Giannios, D., & Spiegel, D. S. 2012, MNRAS, 425, 2778, doi: 10.1111/j.1365-2966.2012.21444.x
  * Metzger & Pejcha (2017) Metzger, B. D., & Pejcha, O. 2017, MNRAS, 471, 3200, doi: 10.1093/mnras/stx1768
  * Metzger et al. (2017) Metzger, B. D., Shen, K. J., & Stone, N. 2017, MNRAS, 468, 4399, doi: 10.1093/mnras/stx823
  * Metzger et al. (2021) Metzger, B. D., Zenati, Y., Chomiuk, L., Shen, K. J., & Strader, J. 2021, arXiv e-prints, arXiv:2108.04305.  https://arxiv.org/abs/2108.04305
  * Michaely & Shara (2021) Michaely, E., & Shara, M. M. 2021, MNRAS, 502, 4540, doi: 10.1093/mnras/stab339
  * Moe & Di Stefano (2017) Moe, M., & Di Stefano, R. 2017, ApJS, 230, 15, doi: 10.3847/1538-4365/aa6fb6
  * Munari et al. (2002) Munari, U., Henden, A., Kiyota, S., et al. 2002, A&A, 389, L51, doi: 10.1051/0004-6361:20020715
  * Nicholls et al. (2013) Nicholls, C. P., Melis, C., Soszynski, I., et al. 2013, MNRAS, 431, L33, doi: 10.1093/mnrasl/slt003
  * Nordhaus et al. (2010) Nordhaus, J., Spiegel, D. S., Ibgui, L., Goodman, J., & Burrows, A. 2010, MNRAS, 408, 631, doi: 10.1111/j.1365-2966.2010.17155.x
  * Nordhaus et al. (2011) Nordhaus, J., Wellons, S., Spiegel, D. S., Metzger, B. D., & Blackman, E. G. 2011, Proceedings of the National Academy of Science, 108, 3135, doi: 10.1073/pnas.1015005108
  * Ohlmann et al. (2016) Ohlmann, S. T., Röpke, F. K., Pakmor, R., & Springel, V. 2016, ApJ, 816, L9, doi: 10.3847/2041-8205/816/1/L9
  * Paczyński (1970) Paczyński, B. 1970, Acta Astron., 20, 47
  * Paczynski (1976) Paczynski, B. 1976, in Structure and Evolution of Close Binary Systems, ed. P. Eggleton, S. Mitton, & J. Whelan, Vol. 73, 75
  * Padmanabhan (2000) Padmanabhan, P. 2000, Theoretical astrophysics. Vol.1: Astrophysical processes
  * Pastorello et al. (2019a) Pastorello, A., Mason, E., Taubenberger, S., et al. 2019a, A&A, 630, A75, doi: 10.1051/0004-6361/201935999
  * Pastorello et al. (2019b) Pastorello, A., Chen, T. W., Cai, Y. Z., et al. 2019b, A&A, 625, L8, doi: 10.1051/0004-6361/201935511
  * Pastorello et al. (2021a) Pastorello, A., Fraser, M., Valerin, G., et al. 2021a, A&A, 646, A119, doi: 10.1051/0004-6361/202039952
  * Pastorello et al. (2021b) Pastorello, A., Valerin, G., Fraser, M., et al. 2021b, A&A, 647, A93, doi: 10.1051/0004-6361/202039953
  * Patterson (1984) Patterson, J. 1984, ApJS, 54, 443, doi: 10.1086/190940
  * Pavlovskii et al. (2017) Pavlovskii, K., Ivanova, N., Belczynski, K., & Van, K. X. 2017, MNRAS, 465, 2092, doi: 10.1093/mnras/stw2786
  * Paxton et al. (2011) Paxton, B., Bildsten, L., Dotter, A., et al. 2011, ApJS, 192, 3, doi: 10.1088/0067-0049/192/1/3
  * Paxton et al. (2013) Paxton, B., Cantiello, M., Arras, P., et al. 2013, ApJS, 208, 4, doi: 10.1088/0067-0049/208/1/4
  * Paxton et al. (2015) Paxton, B., Marchant, P., Schwab, J., et al. 2015, ApJS, 220, 15, doi: 10.1088/0067-0049/220/1/15
  * Paxton et al. (2019) Paxton, B., Smolec, R., Schwab, J., et al. 2019, ApJS, 243, 10, doi: 10.3847/1538-4365/ab2241
  * Pejcha (2014) Pejcha, O. 2014, ApJ, 788, 22, doi: 10.1088/0004-637X/788/1/22
  * Pejcha et al. (2016a) Pejcha, O., Metzger, B. D., & Tomida, K. 2016a, MNRAS, 455, 4351, doi: 10.1093/mnras/stv2592
  * Pejcha et al. (2016b) —. 2016b, MNRAS, 461, 2527, doi: 10.1093/mnras/stw1481
  * Pejcha et al. (2017) Pejcha, O., Metzger, B. D., Tyles, J. G., & Tomida, K. 2017, ApJ, 850, 59, doi: 10.3847/1538-4357/aa95b9
  * Perets et al. (2016) Perets, H. B., Li, Z., Lombardi, James C., J., & Milcarek, Stephen R., J. 2016, ApJ, 823, 113, doi: 10.3847/0004-637X/823/2/113
  * Podsiadlowski et al. (2004) Podsiadlowski, P., Langer, N., Poelarends, A. J. T., et al. 2004, ApJ, 612, 1044, doi: 10.1086/421713
  * Popov (1993) Popov, D. V. 1993, ApJ, 414, 712, doi: 10.1086/173117
  * Prust & Chang (2019) Prust, L. J., & Chang, P. 2019, MNRAS, 486, 5809, doi: 10.1093/mnras/stz1219
  * Reichardt et al. (2020) Reichardt, T. A., De Marco, O., Iaconi, R., Chamandy, L., & Price, D. J. 2020, MNRAS, 494, 5333, doi: 10.1093/mnras/staa937
  * Reichardt et al. (2019) Reichardt, T. A., De Marco, O., Iaconi, R., Tout, C. A., & Price, D. J. 2019, MNRAS, 484, 631, doi: 10.1093/mnras/sty3485
  * Ricker & Taam (2012) Ricker, P. M., & Taam, R. E. 2012, ApJ, 746, 74, doi: 10.1088/0004-637X/746/1/74
  * Sana et al. (2012) Sana, H., de Mink, S. E., de Koter, A., et al. 2012, Science, 337, 444, doi: 10.1126/science.1223344
  * Sand et al. (2020) Sand, C., Ohlmann, S. T., Schneider, F. R. N., Pakmor, R., & Röpke, F. K. 2020, A&A, 644, A60, doi: 10.1051/0004-6361/202038992
  * Schaerer & Maeder (1992) Schaerer, D., & Maeder, A. 1992, A&A, 263, 129
  * Schneider et al. (2019) Schneider, F. R. N., Ohlmann, S. T., Podsiadlowski, P., et al. 2019, Nature, 574, 211, doi: 10.1038/s41586-019-1621-5
  * Shara & Shaviv (1977) Shara, M. M., & Shaviv, G. 1977, MNRAS, 179, 705, doi: 10.1093/mnras/179.4.705
  * Smith (2014) Smith, N. 2014, ARA&A, 52, 487, doi: 10.1146/annurev-astro-081913-040025
  * Smith et al. (2016) Smith, N., Andrews, J. E., Van Dyk, S. D., et al. 2016, MNRAS, 458, 950, doi: 10.1093/mnras/stw219
  * Soberman et al. (1997) Soberman, G. E., Phinney, E. S., & van den Heuvel, E. P. J. 1997, A&A, 327, 620.  https://arxiv.org/abs/astro-ph/9703016
  * Soker (2020) Soker, N. 2020, ApJ, 893, 20, doi: 10.3847/1538-4357/ab7dbb
  * Soker et al. (2018) Soker, N., Grichener, A., & Sabach, E. 2018, ApJ, 863, L14, doi: 10.3847/2041-8213/aad736
  * Soker & Kaplan (2021) Soker, N., & Kaplan, N. 2021, Research in Astronomy and Astrophysics, 21, 090, doi: 10.1088/1674-4527/21/4/90
  * Soker & Kashi (2016) Soker, N., & Kashi, A. 2016, MNRAS, 462, 217, doi: 10.1093/mnras/stw1686
  * Soker & Tylenda (2006) Soker, N., & Tylenda, R. 2006, MNRAS, 373, 733, doi: 10.1111/j.1365-2966.2006.11056.x
  * Stritzinger et al. (2020) Stritzinger, M. D., Taddia, F., Fraser, M., et al. 2020, A&A, 639, A104, doi: 10.1051/0004-6361/202038019
  * Strope et al. (2010) Strope, R. J., Schaefer, B. E., & Henden, A. A. 2010, AJ, 140, 34, doi: 10.1088/0004-6256/140/1/34
  * Sukhbold et al. (2016) Sukhbold, T., Ertl, T., Woosley, S. E., Brown, J. M., & Janka, H. T. 2016, ApJ, 821, 38, doi: 10.3847/0004-637X/821/1/38
  * Taam et al. (1978) Taam, R. E., Bodenheimer, P., & Ostriker, J. P. 1978, ApJ, 222, 269, doi: 10.1086/156142
  * Tauris et al. (2017) Tauris, T. M., Kramer, M., Freire, P. C. C., et al. 2017, ApJ, 846, 170, doi: 10.3847/1538-4357/aa7e89
  * Tout et al. (2008) Tout, C. A., Wickramasinghe, D. T., Liebert, J., Ferrario, L., & Pringle, J. E. 2008, MNRAS, 387, 897, doi: 10.1111/j.1365-2966.2008.13291.x
  * Tylenda (2005) Tylenda, R. 2005, A&A, 436, 1009, doi: 10.1051/0004-6361:20052800
  * Tylenda et al. (2005) Tylenda, R., Soker, N., & Szczerba, R. 2005, A&A, 441, 1099, doi: 10.1051/0004-6361:20042485
  * Tylenda et al. (2011) Tylenda, R., Hajduk, M., Kamiński, T., et al. 2011, A&A, 528, A114, doi: 10.1051/0004-6361/201016221
  * Vigna-Gómez et al. (2018) Vigna-Gómez, A., Neijssel, C. J., Stevenson, S., et al. 2018, MNRAS, 481, 4009, doi: 10.1093/mnras/sty2463
  * Wang et al. (2020) Wang, L., Kroupa, P., Takahashi, K., & Jerabkova, T. 2020, MNRAS, 491, 440, doi: 10.1093/mnras/stz3033
  * Williams et al. (2015) Williams, S. C., Darnley, M. J., Bode, M. F., & Steele, I. A. 2015, ApJ, 805, L18, doi: 10.1088/2041-8205/805/2/L18
  * Wilson & Nordhaus (2019) Wilson, E. C., & Nordhaus, J. 2019, MNRAS, 485, 4492, doi: 10.1093/mnras/stz601
  * Wisniewski et al. (2008) Wisniewski, J. P., Clampin, M., Bjorkman, K. S., & Barry, R. K. 2008, ApJ, 683, L171, doi: 10.1086/591735
  * Woosley (2019) Woosley, S. E. 2019, ApJ, 878, 49, doi: 10.3847/1538-4357/ab1b41