# Nonlinear focusing of Terahertz laser beam using layered superconductor

H.V. Ovcharenko  V.N. Karazin Kharkov National University, 61077 Kharkov,
Ukraine  Z.A. Maizelis  V.N. Karazin Kharkov National University, 61077
Kharkov, Ukraine  A.Ya. Usikov Institute for Radiophysics and Electronics
NASU, 61085 Kharkov, Ukraine  S.S. Apostolov  V.N. Karazin Kharkov National
University, 61077 Kharkov, Ukraine  A.Ya. Usikov Institute for Radiophysics
and Electronics NASU, 61085 Kharkov, Ukraine  V.A. Yampol’skii  V.N. Karazin
Kharkov National University, 61077 Kharkov, Ukraine  A.Ya. Usikov Institute
for Radiophysics and Electronics NASU, 61085 Kharkov, Ukraine

###### Abstract

We theoretically study the propagation of a Terahertz (THz) Gaussian beam
through a thin sample of layered superconductor. We consider the beam axis and
the superconducting layers to be perpendicular to the sample interface, while
the electric field in the beam is perpendicular to the layers. We show that,
in such a geometry, the Josephson current between the superconducting layers
supports lensing of the beam instead of divergence on the Rayleigh range.
Moreover, due to the nonlinearity, the focal length and waist of the
transmitted beam depend on the incident beam intensity. These dependences
demonstrate nontrivial hysteresis behavior that can be observed in experiments
with THz lasers.

###### pacs:

05.45.Ac, 02.50.Ng, 02.60.Lj, 02.30.Hq

##  I Introduction

In recent years, the physical properties of layered superconductors have
attracted the attention of many research groups (see, e.g., review [1] and
references therein). The strongly anisotropic $\rm
Bi_{2}Sr_{2}CaCu_{2}O_{8+\delta}$, $\rm La_{2-\delta}Sr_{\delta}CuO_{4}$, and
$\rm La_{2-\delta}Ba_{\delta}CuO_{4}$ single crystals are the most prominent
examples of such structures [2, 3, 4, 5]. Numerous experiments on the $c$-axis
transport currents in layered high-$T_{c}$ superconductors justify the use of
a model in which the superconducting $\mathrm{CuO_{2}}$ layers are coupled
through the insulator layers by the intrinsic Josephson effect [1, 6, 7]. This
makes the layered superconductors to be anisotropic media not only
quantitatively, but also qualitatively. While the currents in the plane of
layers are of the same nature as in the bulk superconductors, the currents
across the layers are caused by the tunneling of the Cooper pairs and quasi-
particles.

The Josephson current flowing along the $c$-axis is coupled with the
electromagnetic field inside the insulating layers, thereby providing a
specific kind of elementary excitations called Josephson plasma waves (see,
e.g., [8, 7]). These waves are of considerable interest because of their THz
and sub-THz frequency ranges, which are still hardly reachable for both
electronic and optical devices. The frequencies of Terahertz waves are in the
region of resonance frequencies of molecules and are expected to have many
applications [9, 10].

Theoretical studies predict variety of interesting nonlinear phenomena in
layered superconductors even in the regime where the Josephson vortices are
not formed [6]. This becomes possible for the frequencies not far from the
Josephson plasma frequency due to a specific nonlinearity of equations
describing the electrodynamic properties of layered superconductors. In
particular, the nonlinearity results in the hysteresis response of the system
to the electromagnetic excitation [11, 12, 13, 14] and in sensitivity of the
system to the external DC-magnetic field [15, 16]. In all these theoretical
studies, the irradiation of the sample was considered either by the plane
waves or in the waveguide geometry.

Meanwhile, for experimental investigations of layered superconductor
properties, the $\mathrm{ZnTe}$ lasers are commonly used [1, 17] with the
spatially localized radiation. These lasers emit unit pulses in the near-
infrared range. Energy, produced by unit pulse, reaches nanojoules with
electric field of less than kilovolt per centimeter [18]. As was shown in [11,
19], such fields are well described by linear equations inside layered
superconductor. To investigate nonlinear effects, the stronger fields are
needed, and to achieve them the main three types of beam generation can be
used: the tilted pulse front method [20], the free electron lasers [21], and
the gaseous lasers [22].

Tilted pulse method [20, 23, 24] is based on the passage of external near-
infrared radiation through a non-linear $\mathrm{LiNbO_{3}}$ crystal. Due to
nonlinear optical effects, frequency of passed radiation is shifted to the THz
range, which is lower than frequency of incident light. The pulses from such
lasers are strongly localized in frequency (1-2$\%$), their duration varies
from one to hundreds picosecond, time between pulses is from $5\cdot 10^{-5}$
s to $10^{-3}$ s [20, 25], depending on external laser properties. Area of
beam cross-section also depends on external laser and varies from $4\,{\rm
mm}^{2}$ to $30\,{\rm cm}^{2}$. The tilted pulse method was used in Ref. [26]
to investigate response of layered superconductor
$\mathrm{La_{1.84}Sr_{0.16}CuO_{4}}$ to the terahertz external radiation, as
well as to measure frequency dependence of reflectivity and conductivity. To
achieve strong electric fields, laser was focused down to 1 ${\rm mm}^{2}$
beam cross-section, thus field increased up to 100 kV/cm. In Ref. [26], the
frequency of the radiation was 450 GHz, which is below the Josephson plasma
frequency of 2 THz.

The free electron laser sources use radiation by relativistic electrons moving
in non-homogeneous periodically changing magnetic field [21]. This radiation
is also localized in frequency (2 $\%$), however, duration of pulse is tens of
picoseconds, and time between pulses is of the order of $10^{-7}$ s [19, 27,
28]. 2 THz free electron laser with pulse duration of 25 ps was used in Ref.
[19] for excitation of the nonlinear Josephson plasma solitons, predicted in
Ref. [29] for frequencies close enough to the Josephson plasma frequency. The
laser beam was focused to 1 ${\rm mm}^{2}$ beam cross-section, producing 10
kV/cm field.

Recently, rapid development of $\mathrm{CO}$-based gaseous lasers with THz
frequencies shows perspective abilities to use them in investigation of
layered superconductors [22, 30, 31]. These lasers produce constant in time
radiation with the electric field up to fractions of kilovolts per centimeter.
They can be tuned to the pulse regime with the corresponding gain in
intensity.

Thus, the high amplitudes sufficient to observe nonlinear effects in layered
superconductors can be reached by the pulsed radiation. On the other hand, the
pulse duration should be long enough in order to establish the stationary
field distribution in the sample. As one can see, the pulse length reachable
in experiments is of the order of several centimeters for $\sim 100$ ps
pulses. In the present study, we consider samples of several millimeter
thicknesses for which the existing experimental setups can be used to observe
the strongly nonlinear effects even taking into account the multiple
reflections from the interfaces of the superconducting slab.

In this paper, we theoretically investigate the Gaussian beam which falls onto
a thin slab of layered superconductor. We consider the case of constant in
time beam amplitude, which corresponds to the long enough pulses in
experiments. Though the width of the beam does not change significantly within
the slab, the curvature of the wavefront becomes negative due to the
nonlinearity, resulting in the convergence of the beam. To characterize such
nonlinear focusing effect, we introduce focal length $F$ (length from right
sample interface to the point, at which width of beam becomes minimal). We
find the dependence of this parameter on the incident beam amplitude and
frequency. We show that, for the specially chosen frequency and amplitude of
the incident beam, the focusing effect can be strongly increased, and the
focal length can be significantly decreased down to the distance of several
centimeters. Thus, this focusing effect becomes available for the precise
experimental investigation.

The paper is organized as follows. In Section II, we describe the model of the
beam in vacuum and in layered superconductor and present the main equations
for the electromagnetic field in the system. In Section III, taking into
account boundary conditions at both interfaces of the slab, we find the
curvature of the wave front of the transmitted beam which determines the focal
length and the beam waist. In Section IV, we describe the numerical scheme
which we use to verify analytic results. Finally, in Section V, we analyze the
dependence of focusing characteristics on frequency and amplitude of the
incident beam.

##  II Model

###  II.1 Gaussian beam in the vacuum

We start our analysis from considering the behavior of the THz Gaussian beam
in the vacuum. For the beam propagating along the $x$-axis, the well-known
distribution of the electric field reads [32, 33],

$E(x,r,t)=E_{0}\exp\left(-\dfrac{r^{2}}{r_{b}^{2}}\right)\sin\Big{[}k_{v}\Big{(}x+\dfrac{\alpha
r^{2}}{2}\Big{)}-\omega t+\phi\Big{]}.$ (1)

Here $r=\sqrt{y^{2}+z^{2}}$ is the distance from the axis of beam, $k_{v}$ and
$\omega$ are the wave-number and frequency, related to each other by the
dispersion relation $k_{v}=\omega/c$, and $c$ is the speed of light. The
amplitude $E_{0}$, characteristic radius $r_{b}$, wave-front curvature
$\alpha$, and the Gouy phase $\phi$ are the main characteristics of the beam.
Being governed by the Maxwell equations, they vary along the beam path, thus
being the functions of $x$.

Presume the beam transmitted through the slab of layered superconductor, see
Fig. 1, has radius ${r_{b}(x=0)=r_{0}}$ and negative wavefront curvature
${\alpha(x=0)=-\alpha_{0}<0}$. The beam radius $r_{b}(x)$ obeys the following
relation [32]:

$r_{b}(x)=r_{\rm min}\left[1+\dfrac{4(x-F)^{2}}{k_{v}^{2}r_{\rm
min}^{4}}\right]^{1/2},$ (2)

demonstrating convergence of the beam to the minimal reachable radius $r_{\rm
min}$ (the beam waist),

$r_{\rm
min}=r_{0}\left[1+\dfrac{\alpha^{2}_{0}k^{2}_{v}r_{0}^{4}}{4}\right]^{-1/2},$
(3)

at a distance $F$ from the right interface, that we call the focal length,

$F=-\dfrac{1}{\alpha_{0}}\left[1+\dfrac{4}{\alpha_{0}^{2}k^{2}_{v}r_{0}^{4}}\right]^{-1}.$
(4)

Figure 1: (Colour online) Schematics of the setup. The laser beam propagates
from left to right. The incident beam waist is $r_{0}$ at the left interface.
The transmitted beam waist is $r_{\min}$, it is located at the focal length
$F$ from the right interface. The $x$ axis is directed along the beam, the $z$
axis is perpendicular to the layers of the superconducting slab of the
thickness $d$. The incident, reflected, and transmitted beams propagate in the
vacuum regions, while forward and backward beams propagate inside the
superconducting slab (arrows show the directions of correspondent beams).

In the present paper, we investigate the beam transmitted through the slab of
layered superconductor. To be specific, we suppose that the waist point of
incident beam coincides with the slab interface, thus the beam curvature in
this point is zero. In common media, after such point, the beam should diverge
on the Rayleigh range. However, we show that even thin slab of layered
superconductor can make the curvature of the beam negative, thus producing the
focusing effect (although it does not change practically the radius of the
beam within the slab).

###  II.2 Gaussian beam in a slab of layered superconductor

####  II.2.1 Choice of the polarization

We focus our attention on the case when a laser beam falls perpendicularly
onto a surface of layered superconductor, which leads to the strong focusing
of the beam without its distortion. We are aiming to study the situation when
the Josephson tunneling current across the layer plays a decisive role in the
phenomenon under consideration. Therefore, we consider the simplest geometry
where the layers are perpendicular to the slab interface, see Fig. 1. In
another simple geometry, when the layers of the slab are parallel to the
interface, the electric field of the beam lays in the plane of superconducting
layers, and the nonlinear Josephson current is absent.

In the considered setup, the slab interface is anisotropic, which makes the
system highly sensitive to the polarization of the radiation. Actually, the
layered superconductor slab acts as a THz polarization filter in a wide
frequency range. Indeed, if the electric field in the incident wave is
directed along the layers of the slab, the superconducting currents flow only
along the layers, and the so-called ordinary waves are excited in the sample.
These ordinary waves, as was shown in [1], for characteristic frequencies not
far from Josephson plasma frequency $\omega_{J}$, evanesce inside the sample
on the characteristic depth of $\lambda_{ab}\sim 10^{-5}$ cm, which makes the
sample totally reflective for the realistic experimental situations. On the
contrary, the linearly polarized incident wave with the electric field
directed perpendicularly to the layers, generates Josephson current and
excites the so-called extraordinary waves in the slab, which are nonlinear.
Here we focus on the specific polarization with the electric field
perpendicular to the layers, which induces only extraordinary waves and
results in the nonlinear focusing of transmitted beam.

####  II.2.2 Multiple reflections and Gaussian profile

It is important to take into account the internal and external reflections on
both sample interfaces. Thus, we consider incident, reflected, and transmitted
beams in vacuum regions and the beams propagating in both directions in the
superconducting slab. It should be noted that, in general, the electromagnetic
field inside the slab could not be represented as the only two, forward and
backward, Gaussian beams. Firstly, the Gaussian profile could be lost due to
superposition of numerously reflected beams and, secondly, all these
reflections could pump energy into each other due to the nonlinearity of the
system. However, we study the thin slabs, and this assumption allows us
overcome the mentioned problem.

Namely, we assume here that the thickness of the slab is of the order of
effective wavelength inside the slab and that the radius of the beam is much
greater than the characteristic penetration length
$\lambda_{c}=c/(\omega_{J}\sqrt{\varepsilon})$ along the layers, where
$\varepsilon$ is the permittivity of the insulating layers in the slab. These
assumptions are fulfilled for the realistic experimental implications of the
THz beams (see, e.g., Ref. [19]), and allow one to consider the radius and,
thus, amplitude of each reflected beam inside the slab to be constants.
Indeed, using for the estimation Eq. (2) with the wavenumber $k_{s}$ of linear
waves in superconductor [1, 6],

$k_{s}^{2}=\varepsilon(\omega^{2}-\omega_{J}^{2})/c^{2},$ (5)

and characteristic parameters $\omega_{J}/2\pi=2$ THz,
$\omega-\omega_{J}=10^{-3}\omega_{J}$, we get $\Delta r/r\sim 10^{-4}$, where
$\Delta r$ is the variation of the thickness in the slab. This estimation is
also verified by numerical simulation, see Section IV for the details.
Therefore, superposition of all the forward (backward) beams of the same
radius can be regarded as a single forward (backward) beam.

####  II.2.3 Forward and backward nonlinear beams

Thus, taking into account the numerous reflections from the slab interfaces,
within the assumption of small sample thickness and nearly unchanged beam
width, we can seek the electric field in the slab in the form of two Gaussian
beams of the same radius $r_{0}$, propagating in forward (index “$+$”) and
backward (index “$-$”) directions,

$\displaystyle E_{s}(x,r,t)={\cal
E}_{0}\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\Big{(}E_{+}\sin\Phi_{+}+E_{-}\sin\Phi_{-}\Big{)},$
(6)

where the total phases $\Phi_{+}$ and $\Phi_{-}$ are

$\displaystyle\Phi_{\pm}=\pm
k_{\pm}(r)\Big{(}x+\dfrac{\alpha_{\pm}r^{2}}{2}\Big{)}-\omega t+\phi_{\pm}.$
(7)

Here ${\cal E}_{0}=\Phi_{0}/2\pi s\lambda_{c}$ is the characteristic scale of
electric and magnetic fields in layered superconductor, ${\Phi_{0}=\pi
c\hbar/e}$ is the magnetic flux quantum, $s$ is the period of the layered
superconductor structure. Parameters $\alpha_{\pm}$ and $\phi_{\pm}$ are
considered here to be constants, which is correct for the thin samples. This
means that, in the linear approximation, the curvature of wave-front would not
distinctly change. However, in the nonlinear regime, we should take into
account the radial dependency of the wave number $k_{\pm}(r)$, because
intensity of the electromagnetic field depends on the distance $r$ from the
beam axis. Therefore, the cross-terms depending on both $r$ and $x$ appear in
the phase that means the effective curvature of the wave-front changes along
the propagation of the beam, see Eq. (16) for details. In other words, due to
the nonlinearity, different regions of the wave-front have different speeds,
and this makes the beam converge after passing the slab. The latter is the
mechanism of the nonlinear focusing of the THz beam by the layered
superconductor.

###  II.3 Electrodynamic equations for layered superconductors

To find the nonlinear wave-number $k_{\pm}(r)$ we now briefly describe the
electrodynamic equations for a layered superconductor.

####  II.3.1 Equation for the gauge-invariant phase difference of the order
parameter

The field inside the superconducting slab is determined by the gauge-invariant
phase difference $\varphi$ of the order parameter in the neighboring
superconducting layers of the slab (see, e.g., Ref. [6]). This parameter is
actually discrete, but here we consider rather thick beams in comparison with
the period of the layered superconductor structure, $r_{0}\gg s$, when the
continual approximation is valid. In the considering geometry, the electric
field is oriented in perpendicular to the layers and induces the Josephson
tunneling current, $J_{z}=J_{c}\sin{\varphi}$ with $J_{c}$ being the maximal
value of non-dissipative Josephson current density. Meanwhile, the current
along the layers can be described by the London model,
$J_{x}=-c/(4\pi\lambda_{ab}^{2})\,A_{x}$, where $A_{x}$ is the $x$-component
of the vector potential in the slab. The calibration for the vector potential
can be chosen in such a way that the following relation for $A_{z}$ is valid,
$\varphi=-2\pi sA_{z}/\Phi_{0}$ (see, e.g., Ref. [6]).

Then, using the Maxwell equations and the relation of electromagnetic field to
the vector potential, one can express the electric and magnetic fields via the
phase difference $\varphi$,

$\displaystyle E_{s}$ $\displaystyle=$ $\displaystyle-{\cal
E}_{0}\dfrac{1}{\omega_{J}\sqrt{\varepsilon}}\dfrac{\partial\varphi}{\partial
t},$ (8) $\displaystyle\dfrac{\partial H_{s}}{\partial x}$ $\displaystyle=$
$\displaystyle-\dfrac{{\cal
E}_{0}}{\lambda_{c}}\Bigg{[}\dfrac{1}{\omega_{J}^{2}}\dfrac{\partial^{2}\varphi}{\partial
t^{2}}+\sin{\varphi}\Bigg{]},$ (9)

and derive the differential equation [34, 6], which is the continual version
of the well-known coupled sin-Gordon equations,

$\Big{(}1-\lambda_{ab}^{2}\dfrac{\partial^{2}}{\partial
z^{2}}\Big{)}\Bigg{[}\dfrac{1}{\omega_{J}^{2}}\dfrac{\partial^{2}\varphi}{\partial
t^{2}}+\sin{\varphi}\Bigg{]}-\lambda_{c}^{2}\Big{(}\dfrac{\partial^{2}\varphi}{\partial
x^{2}}+\dfrac{\partial^{2}\varphi}{\partial y^{2}}\Big{)}=0.$ (10)

Here the Josephson plasma frequency $\omega_{J}$ is related to the other
parameters of the layered superconductor, ${\omega_{J}=\sqrt{8\pi
esJ_{c}/(\hbar\varepsilon)}}$.

It is important to emphasize that, in spite of strong anisotropy of the
layered superconductor, the Gaussian beam propagating through the thin slab
nearly preserves its axial symmetry of the amplitude distribution in its
cross-section. Indeed, the terms with second derivatives over $y$ and $z$ in
Eq. (10) appear to be small in comparison to other terms, if the slab is thin,
$d\ll r_{0}$. We additionally verify this assumption comparing the cross-
sectional distribution in $y$ and $z$ calculated in numerical simulation with
correspondent analytical results, see Section IV and Fig. 2 for details.

####  II.3.2 Weak nonlinearity in the vicinity of Josephson plasma frequency

The nonlinearity in Eqs. (9) and (10) leads, in principal, to the generation
of higher harmonics both in space and time coordinates. As was reported in
Ref. [6], this generation can be neglected if the phase difference $\varphi$
is small and $\sin\varphi$ can be expanded into series up to the third order,
$\sin{\varphi}\approx\varphi-\varphi^{3}/6$. Usually, such approximation
provides only weak nonlinear effects. However, there is an important range of
frequencies, close enough to $\omega_{J}$, where the nonlinearity plays
crucial role [6]. The strong nonlinear effects can be observed when

$\varphi\sim\beta\equiv\sqrt{\omega^{2}/\omega_{J}^{2}-1}\ll 1,$ (11)

because the linear terms in Eq. (10) nearly cancel each other and the cubic
term becomes significant. In the present paper, we study this frequency range
and predict strong nonlinear focusing of the Gaussian beam.

In order to determine the nonlinear wave-number $k_{\pm}(r)$ in Eq. (7), we
expand all the radial dependencies in series over small $r/r_{0}$ up to the
second order, supposing that the Gaussian distribution for thin slab is
preserved. We find the phase difference $\varphi$ from Eq. (8) with the
electric field $E_{s}$ in the form of Eq. (6) and then substitute it into Eq.
(10). Neglecting higher spatial and temporal harmonics, the nonlinear wave-
numbers $k_{\pm}(r)$ can be related to the amplitudes $E_{\pm}$ of the forward
and backward beams via parameters $\kappa_{\pm}$,

$\displaystyle
k_{\pm}(r)=k_{s}+\dfrac{\kappa_{\pm}}{\lambda_{c}}-\dfrac{\gamma_{\pm}}{\lambda_{c}}\dfrac{r^{2}}{r_{0}^{2}},$
(12)
$\displaystyle\gamma_{\pm}=\kappa_{\pm}\dfrac{2\beta+\kappa_{\pm}}{\beta+\kappa_{\pm}},$
(13) $\displaystyle
8\kappa_{\pm}(2\beta+\kappa_{\pm})=\varepsilon(E_{\pm}^{2}+2E_{\mp}^{2}),$
(14)

where $k_{s}$ is the wavenumber of linear waves, Eq. (5). Note that the
wavenumbers for the forward and backward propagating beams are tangled with
each other via their amplitudes, which is natural for nonlinear problems.

The first term in Eq. (12) corresponds to the linear limit of small
amplitudes, while the second and third terms are provided by the nonlinearity.
The important effect of nonlinearity is the dependence of wave-numbers
$k_{\pm}(r)$ on the radial coordinate, which means that the curvature of the
beam changes along its path. Indeed, recombining terms with $r^{2}/r_{0}^{2}$
in Eq. (7) and neglecting terms $\propto(r/r_{0})^{4}$, we can write
$\Phi_{\pm}$ in the following form,

$\displaystyle\Phi_{\pm}=\pm k_{\pm}^{(0)}\Big{[}x+\dfrac{\alpha^{\rm
eff}_{\pm}(x)r^{2}}{2}\Big{]}-\omega t+\phi_{\pm},$ (15)

introducing the effective curvature $\alpha^{\rm eff}_{\pm}(x)$,

$\displaystyle\alpha^{\rm
eff}_{\pm}(x)=\alpha_{\pm}-\dfrac{2\gamma_{\pm}}{\lambda_{c}k^{(0)}_{\pm}r_{0}^{2}}x,\quad
k^{(0)}_{\pm}=k_{\pm}(r=0).\quad$ (16)

Note that, strictly speaking, Eqs. (12) and (15) are valid only for
$r/r_{0}\ll 1$. However, they are correct practically for all $r$, where
$\exp[-(r/r_{0})^{2}]\sim 1$, see Section IV and Fig. 2 for details.

For magnetic field (which is directed along the $y$ axis) we derive

$H_{s}(x,r,t)={\cal
E}_{0}\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\Big{(}H_{+}\sin\Phi_{+}+H_{-}\sin\Phi_{-}\Big{)},$
(17)

with dimensionless amplitudes $H_{\pm}$ related to $E_{\pm}$,

$H_{\pm}=\mp\sqrt{\varepsilon}E_{\pm}(\beta+\kappa_{\pm}).$ (18)

Now, having the expressions for the electromagnetic field in the
superconducting slab, we proceed with finding the parameters $E_{\pm}$,
$\alpha_{\pm}$, $\phi_{\pm}$ in superconductor and amplitudes, curvatures and
phases for the reflected and transmitted waves in the vacuum.

##  III Focusing of the Gaussian beam by the slab of layered superconductor

In this section, we derive an analytic expression for the curvature
$\alpha_{t}$ in the transmitted beam that defines the focal length $F$ and the
waist $r_{\min}$ of the beam. To that purpose we should relate the parameters
of the beams in vacuum and the slab of layered superconductor by matching the
tangential components of the electric and magnetic fields at the interfaces.
For the slab, we use Eqs. (6) and (17), while the corresponding expressions
for the vacuum regions can be written analogously to Eq. (1). Aiming to
determine the curvatures and phases of the incident, reflected, and
transmitted beams at the interfaces, we present the field only in the vicinity
of the slab. Near the slab, as well as within the slab, we can set the radius
of all beams to $r_{0}$. There exist the incident and reflected beams,

$E_{L}(x,r,t)=\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\big{(}E_{i}\sin\Phi_{i}+E_{r}\sin\Phi_{r}\big{)},$
(19)

in the left vacuum region (see Fig. 1), and there is the transmitted beam
only,

$E_{R}(x,r,t)=E_{t}\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\sin\Phi_{t},$
(20)

in the right vacuum region. Here the phases $\Phi_{t}$ of transmitted and
$\Phi_{r}$ of reflected beams are defined as follows,

$\Phi_{t,r}=\pm k_{v}\Big{(}x+\dfrac{\alpha_{t,r}r^{2}}{2}\Big{)}-\omega
t+\phi_{t,r},$ (21)

with signs $+$ and $-$ corresponding to indices $t$ and $r$, respectively.
Curvatures $\alpha_{r}$ and $\alpha_{t}$ as well as Gouy phases $\phi_{r}$ and
$\phi_{t}$ can be considered as constants only in the vicinity of the slab and
are to be determined from the boundary conditions. Once they are found, it can
be considered as initial values for the correspondent beams to determine the
focal length and the beam waist.

The total phase $\Phi_{i}$ of the incident beam we choose in the following
simple form,

$\Phi_{i}=k_{v}x-\omega t,$ (22)

which means that the incident beam waist position (where the curvature is
absent, $\alpha_{i}=0$) coincides with the left interface.

The corresponding magnetic field in the vacuum can be easily obtained from
Maxwell’s equations,

$\displaystyle H_{L}(x,r,t)$ $\displaystyle=$
$\displaystyle\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\big{(}-E_{i}\sin\Phi_{i}+E_{r}\sin\Phi_{r}\big{)},$
$\displaystyle H_{R}(x,r,t)$ $\displaystyle=$ $\displaystyle-
E_{t}\exp\Big{(}-\dfrac{r^{2}}{r_{0}^{2}}\Big{)}\sin\Phi_{t}.$ (23)

Now we can match the tangential components of the electric and magnetic fields
at the two interfaces between the vacuum and the slab. To obtain the closed
set of equation for the sought quantities, we expand the fields as the
functions of $r$ into series and keep the summands up to the second order in
$r/r_{0}$ only. Straightforward calculations yield the following equations
determining the amplitudes $E_{\pm}$ of the two beams in the slab,

$\displaystyle E_{+}+E_{-}+H_{+}+H_{-}=0,$
$\displaystyle\dfrac{E_{+}-E_{-}}{2}\sin\Big{[}k_{s}d+\dfrac{d}{2\lambda_{c}}(\kappa_{+}+\kappa_{-})\Big{]}=\dfrac{E_{i}}{{\cal
E}_{0}}.$ (24)

Recall that $H_{\pm}$ is related to $E_{\pm}$ by Eq. (17), and $\kappa_{\pm}$
is determined by Eq. (14). Thus, the values of $E_{\pm}$ appear to be
nonlinear functions of the amplitude $E_{i}$ of the incident beam. Moreover,
these functions can be multi-valued, which will be discussed in Section V.

Another result of matching the fields at the vacuum-slab interfaces is the
expression for the curvature of the wave-front of the transmitted beam, which
is valid for the amplitudes $E_{\pm}$ small enough as compared to $\beta$ (and
yet these amplitude can be high enough for the nonlinear effects to be
pronounced),

$\alpha_{t}(x=0)=\dfrac{2d}{k_{\upsilon}r_{0}^{2}}\dfrac{\gamma_{+}E_{+}-\gamma_{-}E_{-}}{E_{+}+E_{-}}.$
(25)

It can be shown, that this curvature is negative and thus the beam converges.
One can substitute this value into Eqs. (3) and (4),
$\alpha_{0}=\alpha_{t}(x=0)$, and thus find the waist $r_{\min}$ of the
transmitted beam and the focal length $F$, which determine the focusing
ability of the superconducting slab.

##  IV Numerical simulation

Before analyzing the waist of the transmitted beam and the focal length as
functions of the system parameters, we describe our numerical simulation
scheme, worked out to verify the obtained analytic results. It should be
emphasized that, in our analytic approach, we make several important
assumptions, such as neglecting of higher spatial and temporal harmonics, and
invariance of the radius and the Gaussian profile of the beam inside the slab.
In order to check them to be correct, we perform a direct numerical simulation
of Eq. (10) for the phase difference $\varphi$ as a function of coordinates
$x,y,z$ inside the slab of layered superconductor and time $t$, while the
incident, reflected, and transmitted beams are accounted by the corresponding
boundary conditions for electric and magnetic fields at each point of the slab
interfaces. Additionally, we take into account the conditions of free
radiation for the reflected and transmitted beams, which is standard routine
in numerical simulations of spatially unbounded systems.

The size of the sample is chosen to be ${N=3}$ times greater than the diameter
of the beam, namely 2.1 cm $\times$ 2.1 cm for the $7$ mm wide laser. These
dimensions are high enough to neglect the edge effects due to the fact that
the electromagnetic field at the lateral edges is $\exp(N^{2})$ times weaker
than on the axis of the beam.

Technically, the numerical simulation of Eq. (10) is performed for several
values of the incident beam amplitude. For each value of $E_{i}$, the
following procedure is performed. The amplitude is gradually increased from
zero to $E_{i}$ slowly enough and then is kept constant, until the
electromagnetic field in the whole sample becomes well established. The
obtained distribution of the field is used to estimate deviations from the
expected analytic result and to calculate the focusing parameters of the
transmitted beam.

Figure 2: (Color online) Spatial distributions of the amplitude (panels $a$
and $c$) and the total phase (panels $b$ and $d$) of electric field in the
transmitted beam at the right interface. Panels $a$ and $b$ show the
distributions obtained by the numerical simulation as functions of spatial
coordinates $y$ and $z$ by the color gradient. Panels $c$ and $d$ present
these numerical distributions in two perpendicular cross-sections, vertical
(at $y=0$ as function of $z$) and horizontal (at $z=0$ as function of $y$),
plotted by the blue lines with solid circles and red lines with empty squares.
These curves are compared with the analytically obtained distributions (as
functions of $r$), plotted by the green lines and calculated as described in
Section III. The same blue and red colors are used for the dashed straight
lines in panels $a$ and $c$ to indicate respective cross-sections. The
distributions are calculated for the incident beam with amplitude $E_{i}=0.2$
kV/cm with pronounced nonlinear effects. Other parameters are:
$\omega_{J}/2\pi=2$ THz, $\omega/\omega_{J}-1=1.1\cdot 10^{-3}$,
$\lambda_{c}/\lambda_{ab}=15$, $r_{0}=3.5$ mm, $d=2.5$ mm, $s=2\cdot 10^{-7}$
m, $\varepsilon=15$.

Figure 2 shows the spatial distributions of the amplitude (panels $a$ and $c$)
and the total phase (panels $b$ and $d$) of electric field in the transmitted
beam at the right interface, $x=0$. The phase is counted from its value in the
central point of the cross-section,

$\Delta\Phi_{t}(y,z)=\Phi_{t}(y,z,t)-\Phi_{t}(y=0,z=0,t).$

Panels $a$ and $b$ present the distributions obtained by the numerical
simulation as functions of spatial coordinates $y$ and $z$ by color gradient.
Panels $c$ and $d$ are assigned with analytic results for accurate comparison
with numerically obtained distributions. Namely, the green curves present the
amplitude $E_{t}\exp(-r^{2}/r_{0}^{2})$ (panel $c$) and the total phase
$\Delta\Phi_{t}(r)=\Phi_{t}(r,t)-\Phi_{t}(r=0,t)$ (panel $d$) as functions of
the spatial coordinate $r$, calculated by analytic approach described in
Sections II and III. The blue and red lines with solid circles and empty
squares correspond to the numerically obtained distributions in two
perpendicular cross-sections, vertical (at $y=0$ as function of $z$) and
horizontal (at $z=0$ as function of $y$), respectively. The same colors are
used for the dashed straight lines in panels $a$ and $b$ to indicate
corresponding cross-sections.

As can be seen from panels $c$ and $d$ in Fig. 2, even for the intensive
enough incident beam with amplitude $E_{i}=0.2$ kV/cm, where nonlinear effects
are distinctly pronounced (see Section V and Fig. 4 for details), the field
distribution appears to be nearly isotropic and Gaussian. Note that the
seeming deviation of the amplitude from the Gaussian in the region of small
radii is due to a slightly narrower distribution of the beam in the
simulation, which results in the higher amplitude in the central part. The
total phase is well consistent with the quadratic dependence as assumed in the
model, while the deviation far from the beam axis are inessential there
because of the exponentially small amplitude, i.e. $\exp[-(r/r_{0})^{2}]\ll
1$.

Thus, the numerical simulation ensures that our assumptions made in the model
are reasonable and we can use analytic results to study the nonlinear focusing
of the Gaussian beam by a thin plate of layered superconductor.

##  V Analysis of the results

In this section, we apply both the analytic approach and numerical simulation
to study the dependence of focal length and the waist of the transmitted beam
on the frequency and amplitude of the incident radiation.

We start from the dependence on frequency, which is especially interesting in
the vicinity of the Josephson plasma frequency $\omega_{J}$. Figure 3 shows
the dependence of the focal length $F$ (upper blue curve) and the waist
$r_{\min}/r_{0}$ of the transmitted beam normalized to the initial radius
(lower red curve) as the functions of frequency detuning
$\omega/\omega_{J}-1$. One can see that, as the frequency $\omega$ of laser
beam approaches $\omega_{J}$, both the focal length and beam waist
pronouncedly oscillate. These oscillations emerge from the variation of the
wavelength in the sample, and are analogous to the Fabry-Pérot oscillations.
Indeed, the critical points of the curves in Fig. 3 correspond to the
frequencies, for which the thickness of the sample is equal to integer number
of wavelengths, $2\pi n/k_{\pm}^{(0)}$. However, the wavelength is not
proportional to $\omega^{-1}$ as in linear optics, but nonlinearly depends on
the frequency detuning $\omega/\omega_{J}-1$, significantly affecting not only
the amplitude of the transmitted beam, but also the focal length and the beam
waist.

It should be noted that the similar oscillations of the transmission
coefficients of plain waves and wave-guide modes when changing frequency
detuning $\omega/\omega_{J}-1$ were predicted in Refs. [12, 13]. So, the
nonlinear effects lead not only to rapid increase of the transmitted
amplitude, but also to strong focusing effect of the Gaussian beam.

Figure 3: (Colour online) Dependence of focal length $F$ (upper blue curve)
and normalized beam waist $r_{\min}/r_{0}$ (lower red curve) on frequency
detuning of the incident beam. Points $a$ and $b$ correspond to the local
maxima at $\omega/\omega_{J}-1\approx 0.75\cdot 10^{-4}$ and at
$\omega/\omega_{J}-1\approx 2\cdot 10^{-4}$ on the lower curve. The amplitude
of incident beam $E_{i}$ is 50 V/cm. The insets show the intensity
distribution in the transmitted beam for frequencies corresponding to points
$a$ and $b$, while dashed curves show $1/e^{2}$ width. Other parameters are
the same as in Fig. 2.

Now let us analyze the effect of the beam intensity, which also comes from the
nonlinearity of the problem. Figure 4 shows the dependence of the focal length
(upper panel) and transmitted beam waist (lower panel) on the amplitude of the
incident beam. The green points with vertical error bars in Fig. 4 show the
results of the numerical simulation described in Section IV. One can see that
the simulation points fit well the analytic curve, and the simulation errors
increase with the growth of the incident field amplitude.

It is important that, according to Eq. (16) with $\kappa_{\pm}$ from Eq. (14),
the dependence of curvature $\alpha_{t}$ on the incident beam amplitude
$E_{i}$ is implicit and strongly nonlinear. Therefore, as seen from Fig. 4,
the focal length appears to be non-monotonic as function of $E_{i}$. Moreover,
if one increases the amplitude even more, upto approximately 0.7 kV/cm for the
chosen parameters, the dependences become even multi-valued (see the right
panels). One can see that each dependence $F(E_{i})$ and $r_{\rm min}(E_{i})$
involves three branches. The first branch starts at $E_{i}=0$ and ends up at a
certain critical value $E_{i,1}$, while the third branch starts at a certain
critical value $E_{i,2}$ and goes to the greater values. The second
(intermediate) branch shown by dashed line between $E_{i,2}$ and $E_{i,1}$ is
unreachable as is usual for such multi-valued dependences. If one gradually
increases/decreases the amplitude of the incident beam, the parameters of the
transmitted beam can exhibit hysteresis behavior with jumps marked by arrows
in Fig. 4.

Figure 4: (Colour online) Dependence of focal length (upper panel) and
transmitted beam waist (lower panel) on the amplitude of incident beam. Green
points show the results of the numerical simulation. The left and right panels
show the dependences in different ranges of $E_{i}$. Solid and dotted curves
show reachable and unreachable branches. Arrows show possible hysteresis
variation of the focusing parameters, when the incident beam intensity
increases/decreases gradually. Frequency detuning is
$\omega/\omega_{J}-1=1.1\cdot 10^{-3}$, other parameters are the same as in
Fig. 2.

It should be noted that the accuracy of the simulation did not allow us to
check the hysteresis behavior though we tried several strategies to achieve
higher amplitudes with high enough accuracy in the simulations. So, it is an
open problem to find out whether correspondent hysteresis behavior can be
attained in the numerical simulation or/and in the real experiments. Yet, we
clearly see that our analytic approach well fits the results of numerical
simulation and, thus, can be used to predict the behavior of laser beam in
nonlinear layered superconducting slabs and describe its focusing parameters.

##  VI Conclusions

In this paper, we have studied theoretically propagation of the Gaussian laser
beam through the layered superconductor slab with layers perpendicular to the
slab interface. We have chosen polarization where the electric field is
perpendicular to the layers thus inducing the Josephson interlayer currents in
the sample. Solving differential equations for the electromagnetic field
distribution in the slab with appropriate boundary conditions, accounting for
the reflections from the interfaces, we have presented the dependences of the
focal length and the transmitted beam waist on the frequency and amplitude of
the incident beam in an implicit algebraic form. We also have performed the
numerical simulation to determine the field distribution in the slab and,
thus, verified the analytic results.

We have shown that, in the nonlinear regime, the laser beam acquires the
negative curvature of the wavefront after passing the layered superconductor
slab which leads to the beam focusing. The focusing effect strongly depends on
amplitude and frequency of the incident radiation, and becomes more pronounced
for the frequencies close to the Josephson plasma frequency and for high
enough amplitudes. Note that the results are shown for the small frequency
detunings and for the reachable amplitudes of order of 1 kV/cm. If the greater
amplitudes are used, the frequency may be detuned further from the Josephson
plasma frequency. Moreover, the analysis of the results admits that the focal
length and the transmitted beam waist can show hysteresis behavior, when the
amplitude of the incident laser beam is increased and then decreased
gradually.

##  VII Acknowledgment

We gratefully acknowledge support from the National Research Foundation of
Ukraine, Project No. 2020.02/0149 “Quantum phenomena in the interaction of
electromagnetic waves with solid-state nanostructures”.

## References

  * Laplace and Cavalleri [2016] Y. Laplace and A. Cavalleri, Advances in Physics: X 1, 387 (2016).
  * Curran  _et al._ [2018] P. J. Curran, H. A. Mohammed, S. J. Bending, A. E. Koshelev, Y. Tsuchiya, and T. Tamegai, Scientific Reports 8, 10914 (2018).
  * Nicoletti  _et al._ [2018] D. Nicoletti, D. Fu, O. Mehio, S. Moore, A. S. Disa, G. D. Gu, and A. Cavalleri, Phys. Rev. Lett. 121, 267003 (2018).
  * Hücker  _et al._ [2013] M. Hücker, M. v. Zimmermann, Z. J. Xu, J. S. Wen, G. D. Gu, and J. M. Tranquada, Phys. Rev. B 87, 014501 (2013).
  * Khanna  _et al._ [2016] V. Khanna, R. Mankowsky, M. Petrich, H. Bromberger, S. A. Cavill, E. Möhr-Vorobeva, D. Nicoletti, Y. Laplace, G. D. Gu, J. P. Hill, M. Först, A. Cavalleri, and S. S. Dhesi, Phys. Rev. B 93, 224522 (2016).
  * Savel’ev  _et al._ [2010] S. Savel’ev, V. A. Yampol’skii, A. L. Rakhmanov, and F. Nori, Rep. Prog. Phys. 73, 026501 (2010).
  * Hu and Lin [2010] X. Hu and S.-Z. Lin, Superconductor Science and Technology 23, 053001 (2010).
  * Averkov  _et al._ [2013] Y. O. Averkov, V. M. Yakovenko, V. A. Yampol’skii, and F. Nori, Phys. Rev. B 87, 054505 (2013).
  * Salén  _et al._ [2019] P. Salén, M. Basini, S. Bonetti, J. Hebling, M. Krasilnikov, A. Y. Nikitin, G. Shamuilov, Z. Tibai, V. Zhaunerchyk, and V. Goryashko, Physics Reports 836-837, 1 (2019).
  * Batov  _et al._ [2006] I. E. Batov, X. Y. Jin, S. V. Shitov, Y. Koval, P. Müller, and A. V. Ustinov, Applied Physics Letters 88, 262504 (2006).
  * Yampol’skii  _et al._ [2008] V. A. Yampol’skii, T. M. Slipchenko, Z. A. Mayzelis, D. V. Kadygrob, S. S. Apostolov, S. E. Savel’ev, and F. Nori, Phys. Rev. B 78, 184504 (2008).
  * Apostolov  _et al._ [2010] S. S. Apostolov, Z. A. Maizelis, M. A. Sorokina, V. A. Yampol’skii, and F. Nori, Phys. Rev. B 82, 144521 (2010).
  * Rokhmanova  _et al._ [2013] T. N. Rokhmanova, S. S. Apostolov, Z. A. Maizelis, V. A. Yampol’skii, and F. Nori, Phys. Rev. B 88, 014506 (2013).
  * Rokhmanova  _et al._ [2014] T. N. Rokhmanova, S. S. Apostolov, Z. A. Maizelis, V. A. Yampol’skii, and F. Nori, Phys. Rev. B 90, 184503 (2014).
  * Apostolov  _et al._ [2016] S. S. Apostolov, Z. A. Maizelis, N. M. Makarov, F. Pérez-Rodríguez, T. N. Rokhmanova, and V. A. Yampol’skii, Phys. Rev. B 94, 024513 (2016).
  * Kvitka  _et al._ [2021] N. Kvitka, S. S. Apostolov, N. M. Makarov, T. Rokhmanova, A. A. Shmat’ko, and V. A. Yampol’skii, Phys. Rev. B 103, 104512 (2021).
  * Lee  _et al._ [2021] H.-T. Lee, G.-S. Ji, J.-Y. Oh, C.-W. Seo, B.-W. Kang, K.-W. Kim, and H.-R. Park,  11, 651 (2021).
  * Vidal  _et al._ [2013] S. Vidal, J. Degert, M. Tondusson, E. Freysz, and J. Oberlé, Journal of the Optical Society of America B 31, 149 (2013).
  * Dienst  _et al._ [2013] A. Dienst, E. Casandruc, D. Fausti, L. Zhang, M. Eckstein, M. Hoffmann, V. Khanna, N. Dean, M. Gensch, S. Winnerl, W. Seidel, S. Pyon, T. Takayama, H. Takagi, and A. Cavalleri, Nature materials 12 (2013).
  * Wang  _et al._ [2020] L. Wang, G. Tóth, J. Hebling, and F. Kärtner, Laser & Photonics Reviews 14, 2000021 (2020).
  * Huang and Kim [2007] Z. Huang and K.-J. Kim, Phys. Rev. ST Accel. Beams 10, 034801 (2007).
  * Chevalier  _et al._ [2019] P. Chevalier, A. Amirzhan, F. Wang, M. Piccardo, S. G. Johnson, F. Capasso, and H. O. Everitt, Science 366, 856 (2019).
  * Nugraha  _et al._ [2019] P. S. Nugraha, G. Krizsán, C. Lombosi, L. Pálfalvi, G. Tóth, G. Almási, J. A. Fülöp, and J. Hebling, Optics Letters 44, 1023 (2019).
  * Hebling  _et al._ [2008] J. Hebling, K. L. Yeh, M. C. Hoffmann, B. Bartal, and K. A. Nelson, Journal of The Optical Society of America B-optical Physics 25 (2008).
  * Avetisyan  _et al._ [2017] Y. Avetisyan, A. Makaryan, V. Tadevosyan, and M. Tonouchi, Journal of Infrared, Millimeter, and Terahertz Waves 38, 1439 (2017).
  * Dienst  _et al._ [2011] A. Dienst, M. C. Hoffmann, D. Fausti, J. C. Petersen, S. Pyon, T. Takayama, H. Takagi, and A. Cavalleri, Nature Photonics 5, 485 (2011).
  * Piovella and Volpe [2021] N. Piovella and L. Volpe, Atoms 9, 28 (2021).
  * Huang  _et al._ [2018] R. Huang, W. Li, Z. Zhao, H. Li, J. Wang, T. Ma, Q. Huang, Z. He, Q. Jia, L. Wang, and Y. Lu, Particles 1, 267 (2018).
  * Rajaraman [1982] R. Rajaraman,  _Solitons and Instantons: An Introduction to Solitons and Instantons in Quantum Field Theory_ (Elsevier Science Ltd, 1982).
  * Bosco  _et al._ [2019] L. Bosco, M. Franckié, G. Scalari, M. Beck, A. Wacker, and J. Faist, Applied Physics Letters 115, 010601 (2019).
  * Wang  _et al._ [2018] F. Wang, J. Lee, D. J. Phillips, S. G. Holliday, S.-L. Chua, J. Bravo-Abad, J. D. Joannopoulos, M. Soljačić, S. G. Johnson, and H. O. Everitt, Proceedings of the National Academy of Sciences 115, 6614 (2018).
  * Simon [2016] D. S. Simon,  _A Guided Tour of Light Beams: From lasers to optical knots_ (IOP Publishing, 2016).
  * Guenther [1990] R. Guenther,  _Modern Optics_, Wiley international edition (Wiley, 1990).
  * Artemenko and Remizov [1997] S. N. Artemenko and S. V. Remizov, Journal of Experimental and Theoretical Physics Letters 66, 853 (1997).


# The compact circumstellar material of SN 2024ggi: Another supernova with a
pre-explosion effervescent zone and jet-driven explosion

Noam Soker  Department of Physics, Technion, Haifa, 3200003, Israel;
<EMAIL_ADDRESS>
###### Abstract

I examine recent theoretical studies and observations of the recent core-
collapse supernova (CCSN) SN 20224ggi and find that the likely explanation for
its dense, compact circumstellar material is an effervescent model, where
parcels or streams of gas are uplifted by stellar convection and pulsation and
fall back. The effervescent zone exists alongside the regular wind from the
red supergiant (RSG) progenitor of SN 2024ggi. I find that an extended wind-
acceleration zone encounters some difficulties in accounting for the required
CSM mass. Recent modelling finds the explosion energy of SN 2024ggi to be
$E_{\rm exp}>10^{51}{~{}\rm erg}$, and up to $E_{\rm exp}\simeq 2\times
10^{51}{~{}\rm erg}$. I examine this explosion energy against a recent study
of the delayed neutrino explosion mechanism and find that this mechanism has
difficulties in accounting for the required energy. This suggests that the
explosion was caused by the jittering jets explosion mechanism (JJEM), adding
to other recent pieces of evidence supporting the JJEM, particularly point-
symmetric CCSN remnants.

###### keywords:

stars: massive – stars: mass-loss – supernovae: general; supernova:
individual: SN 2024ggi

##  1 Introduction

Studies fit the type IIP core-collapse supernova (CCSN) SN 2024ggi lightcurve
and spectroscopic properties with a compact pre-explosion circumstellar
material (CSM; e.g., Chen et al. 2024a, b; Jacobson-Galán et al. 2024; Pessi
et al. 2024; Shrestha et al. 2024; Xiang et al. 2024; Zhang et al. 2024). By
compact CSM, I refer to CSM with which the ejecta collides within days after
the explosion, unlike SN 1987A, where a collision occurs years after the
explosion. In principle, compact dense CSM of CCSN progenitors might result
from an enhanced mass loss rate that starts years to weeks before explosion,
possibly accompanied by a pre-explosion outburst (e.g., Foley et al. 2007;
Pastorello et al. 2007; Smith et al. 2010; Margutti et al. 2014; Ofek et al.
2014; Svirski & Nakar 2014; Tartaglia et al. 2016; Yaron et al. 2017; Wang,
Wang, & Dai 2019; Bruch et al. 2020; Prentice et al. 2020; Strotjohann et al.
2021; Jacobson-Galán et al. 2022), an extended long-lived dense zone above the
stellar photosphere (e.g., Dessart et al. 2017), an extended accelerated zone
of the wind (e.g., Moriya et al. 2017, 2018) or by a long-lived extended dense
zone of uprising and falling gas parcels or streams above the stellar
photosphere (e.g., Soker 2021, 2023; Fuller & Tsuna 2024; see further
discussion by Fuller & Tsuna 2024). SN 2024ggi did not experience an outburst
within years before the explosion (e.g., Shrestha et al. 2024).

The properties of the compact CSM of SN 2024ggi are not much different from
those of SN 2023ixf (for SN 2023ixf, e.g., Berger et al. 2023; Bostroem et al.
2023; Grefenstette et al. 2023; Jacobson-Galan et al. 2023; Kilpatrick et al.
2023; Singh Teja et al. 2023; Smith et al. 2023), which also did not have any
pre-explosion outburst that could have formed the compact CSM (e.g., Jencson
et al. 2023; Neustadt, Kochanek, & Rizzo Smith 2023; Soraisam et al. 2023). In
Soker (2023), I argued that the best explanation for the compact pre-explosion
CSM of SN 2023ixf is an effervescent Zone. In this study, I make a similar
claim for the compact CSM of SN 2024ggi.

Most studies of SN 2024ggi mentioned above ignore the effervescent zone model
and the similar model by Fuller & Tsuna (2024), which includes many more
details; in both models, convective motion and pulsation uplift the bound
material, but Fuller & Tsuna (2024) calculate the uplifting of material by
shocks that the convection excites above the photosphere. In this study, I use
the term effervescent zone model to include both the model studied by Soker
(2023) and Fuller & Tsuna (2024). In some aspects, but not all, the
hydrostatic model of Dessart et al. 2017 can also be grouped with the
effervescent zone model.

In section 2, I examine the proposed wind acceleration zone that Chen et al.
(2024a) construct for the dense compact CSM of SN 2024ggi. I reiterate my
claim that the effervescent zone should not be ignored. In section 3 I comment
on the interesting finding by Chen et al. (2024a) that their best fit gives an
explosion energy of $E_{\rm exp}\simeq 2\times 10^{51}{~{}\rm erg}$. This has
implications for the likely explosion mechanism of SN 2024ggi that I
emphasize. In my summary in section 4, I comment on the relation between the
effervescent zone model and the explosion mechanism.

##  2 The compact CSM of SN 2024ggi

A most recent study and analysis of the compact CSM of SN 2024ggi is that by
Chen et al. (2024a). They consider the wind acceleration model of Moriya et
al. (2018) for a compact CSM that extends to $R_{\rm CSM}\simeq 6\times
10^{14}{~{}\rm cm}$, and that has a mass loss rate of $\dot{M}_{\rm CSM}\simeq
0.001M_{\odot}{~{}\rm yr}^{-1}$ for a wind terminal velocity of $v_{\rm
w}=10{~{}\rm km}s^{-1}$. The wind velocity is lower near the star as it is
accelerated, such that the total wind mass within the radius of $R_{\rm CSM}$
is $M_{\rm CSM}\simeq 0.4M_{\odot}$. There are two challenges to this CSM
model, as I discuss next.

The first challenge is to explain what accelerates the wind. The reason is
that the momentum in the radiation is much smaller than that of the wind.
Momentum balance limits the mass loss of the wind that radiation can
accelerate. It is given by the equation $\dot{M}_{\rm wc}v_{\rm
w}\simeq\eta_{\rm w}L/c$, where $v_{\rm w}$ is the terminal wind speed, $L$
the stellar luminosity, and $\eta_{\rm w}$ is the average number of times that
a photon transfers its momentum to the wind in the outward radial direction,

$\displaystyle\begin{aligned} \dot{M}_{\rm wc}&\simeq 4\times 10^{-4}\eta_{\rm
w}\left(\frac{L}{2\times 10^{5}L_{\odot}}\right)\\\ &\times\left(\frac{v_{\rm
w}}{10{~{}\rm km}{~{}\rm s}^{-1}}\right)^{-1}M_{\odot}{~{}\rm
yr}^{-1}.\end{aligned}$ (1)

The mass loss rate of the wind-accelerated zone of Chen et al. (2024a)
requires $\eta_{\rm w}\simeq 2.5$. This requires several scattering of the
photons, which is unclear how to achieve to distances of $>10{~{}\rm AU}$ from
the star. With the same model Moriya & Singh (2024) estimate the mass loss
rate to build the CSM of SN 2023ixf to be $0.001-0.01M_{\odot}{~{}\rm
yr}^{-1}$, which requires the highly-unlikely value of up to $\eta_{\rm
w}\simeq 25$.

The demand for $\eta_{\rm w}\simeq 2.5$ at distances of tens of AU from the
star is challenging for the wind-acceleration zone model. Note that the number
$\eta_{\rm w}$ does not refer to the total scattering of photons but to the
scattering in the radial direction, such that the photon deposits a radial
outward momentum. A simple optically-thick medium does not increase $\eta_{\rm
w}$ above 1. One needs scattering from different sides of the star, namely,
that photons bounce from one side to another.

The terminal wind velocity that Chen et al. (2024a) and Moriya & Singh (2024)
take in their model, $v_{\rm w}=10{~{}\rm km}s^{-1}$, is much lower than the
escape velocity from red supergiant (RSGs) progenitors of SNe II. It is even
more problematic as this is lower than the escape velocity at the outer
boundary of the compact dense CSM

$\displaystyle\begin{aligned} v_{\rm w}&=10{~{}\rm km}{~{}\rm s}^{-1}<v_{\rm
esc}(R_{\rm CSM})=13.3\\\ &\times\left(\frac{M_{\rm
RGB}}{4M_{\odot}}\right)^{1/2}\left(\frac{R_{\rm CSM}}{6\times 10^{14}{~{}\rm
cm}}\right)^{-1/2}{~{}\rm km}{~{}\rm s}^{-1}.\end{aligned}$ (2)

The RSG mass is likely higher than the above scaling of $M_{\rm
RGB}=4M_{\odot}$. This implies that the acceleration zone must continue beyond
$R_{\rm CSM}=6\times 10^{14}{~{}\rm cm}$.

The most important and robust parameter of the CSM that Chen et al. (2024a)
and Moriya & Singh (2024) obtain in their modeling is the CSM mass. If the
terminal wind velocity is taken to be higher, as should be for such RSG stars,
by a factor of $k_{\rm w}$, then the mass loss rate should also be $k_{\rm w}$
higher. This implies a wind radial momentum larger by a factor of $k^{2}_{\rm
w}$, and so is the value of $\eta_{\rm w}$. This would make the value of
$\eta_{\rm w}\gg 10$ unrealistically high.

This section concludes that the wind-accelerated zone model for a compact and
dense CSM of SN 2023ixf and SN 2024ggi encounters severe challenges. I
consider this model unlikely to explain the CSM of SN 2024ggi and SN 2023ixf.
However, the total CSM mass and the explosion energy that Chen et al. (2024a)
and Moriya & Singh (2024) obtain in their modeling hold.

##  3 The implication of $E_{\rm exp}\ga 2\times 10^{51}{~{}\rm erg}$

With their accelerated-wind zone model for the CSM, Chen et al. (2024a)
estimated the explosion energy of SN 2024ggi to be $E_{\rm exp}\simeq 2\times
10^{51}{~{}\rm erg}$, and Moriya & Singh (2024) estimate the explosion energy
of SN 2023ixf to be $E_{\rm exp}\simeq 2-3\times 10^{51}{~{}\rm erg}$. Their
fitting might apply to the effervescent zone model if it has a similar density
profile. Jacobson-Galán et al. (2024) adopted an explosion energy of $E_{\rm
exp}\simeq 1.2\times 10^{51}{~{}\rm erg}$ for SN 2024ggi. Such explosion
energies challenge the delayed neutrino explosion mechanism of CCSNe, which
cannot account for explosion energies of $E_{\rm exp}\ga 2\times
10^{51}{~{}\rm erg}$, and even for $E_{\rm exp}\ga 10^{51}{~{}\rm erg}$ it has
severe difficulties to account for (e.g., Fryer et al. 2012; Ertl et al. 2016;
Sukhbold et al. 2016; Gogilashvili, Murphy, & Mabanta 2021). To emphasize this
limitation of the delayed neutrino explosion mechanism, I critically examine a
most recent study by Nakamura et al. (2024), who write, “We find that
neutrino-driven explosions occur for all models within 0.3 s after bounce.”

Nakamura et al. (2024) conduct three-dimensional magnetohydrodynamical
simulations of the delayed neutrino explosion mechanism for stellar models
with different zero-age main sequence masses from $9M_{\odot}$ to
$24M_{\odot}$. Their stellar models are from Sukhbold et al. (2016). The
diagnostic explosion energy of their results is the sum of all numerical cells
where the total energy is positive, and the materials have a positive radial
velocity. The integrand is the sum of kinetic, gravitational, internal, and
magnetic energy. In Figure 1, I present their figure of the evolution of the
diagnostic explosion energy, $E_{\rm Diag}$ with time.

Figure 1: A figure adapted from Nakamura et al. (2024). The different lines
show the evolution of the diagnostic energy in their simulations of the
delayed neutrino explosion mechanism of stellar models differ in their zero-
age main sequence mass, as indicated on the inset on the left (foe stands for
fifty-one erg, i.e., $10^{51}{~{}\rm erg}$). The diagnostic energy includes
only grid cells with positive energy that moves radially outwards. With three
horizontal black lines, I emphasize cases where the diagnostic energy does not
increase at the end of their simulation, or, for $M_{\rm ZAMS}=14M_{\odot}$,
even decreases at the end. In the inset attached to the three cases with the
largest diagnostic energies, I write the estimated binding energy of the
stellar parts with negative energy, namely the binding energy $\Delta E_{\rm
OV}$ (overburden; based on Mezzacappa et al. 2015), as well as the expected
final explosion energy $E_{\rm exp}=E_{\rm Diag}-\Delta E_{\rm OV}$: $E_{\rm
exp}(14M_{\odot})\simeq 0.1\times 10^{51}{~{}\rm erg}$ and $E_{\rm
exp}(24M_{\odot})\simeq 0.1\times 10^{51}{~{}\rm erg}$. For the cases of
$M_{\rm ZAMS}=16,18,19,20$, the overburden is larger than the diagnostic
energy, and I expect no explosion (or one with very low energy of $\ll 0.1$
foe) according to the simulations by Nakamura et al. (2024). I mark this only
for $E_{\rm exp}(20M_{\odot})<0$. Using the overburden

The final explosion energy is the diagnostic energy minus the binding energy
of the rest of the envelope and core (termed overburden), namely, $R_{\rm
exp}=E_{\rm Diag}-\Delta E_{\rm OV}$. Because Nakamura et al. (2024) do not
give the binding energy of the rest of their stellar models (overburden), I
take overburden energies from the simulations of Mezzacappa et al. (2015) who
give this for 4 models at $t\simeq 1.2{~{}\rm s}$ when the shock radius is
larger than those at the end of the simulations by Nakamura et al. (2024) each
at $t=0.3{~{}\rm s}$. This implies that the binding energy (overburden) values
should be larger even. Also, Mezzacappa et al. (2015) and Nakamura et al.
(2024) do not use the same models. Therefore, the binding energy I give in the
insets in Figure 1 are approximate values. These values, nonetheless, clearly
point to the very low final explosion energies (to distinguish from the
diagnostic explosion energy) that Nakamura et al. (2024) obtain. Those models
that explode have energies of $E_{\rm exp}\lesssim 0.1\times 10^{51}{~{}\rm
erg}$ while some others seem to have negative energy, namely, they do not
explode, e.g., the simulations for $M_{\rm ZAMS}=16,18,19,20$ lead to no
explosion.

For the purpose of the present study, relevant is the limit of $E_{\rm
exp}<1.2\times 10^{51}{~{}\rm erg}$ of most simulation results of the delayed
neutrino explosion mechanism. In their very long simulation to 7 seconds after
bounce, Bollig et al. (2021) find the explosion energy of a model with $M_{\rm
ZAMS}=19M_{\odot}$ to be $\simeq 10^{51}{~{}\rm erg}$. Even the diagnostic
explosion energies of Nakamura et al. (2024) do not reach these values.
Mezzacappa et al. (2015) obtain higher explosion energies than Nakamura et al.
(2024), but still of $E_{\rm exp}<10^{51}{~{}\rm erg}$. Burrows, Wang, &
Vartanyan (2024) obtain higher explosion energies, but still most models have
$E_{\rm exp}<1.2\times 10^{51}{~{}\rm erg}$. Their results do not agree with
those of Nakamura et al. (2024), e.g., Burrows, Wang, & Vartanyan (2024)
obtain clear explosion with $E_{\rm exp}\simeq 10^{51}{~{}\rm erg}$ in models
with negative value of $E_{\rm exp}$ in the simulations of Nakamura et al.
(2024), e.g., $M_{\rm ZAMS}=16,18,19,20$. Such qualitative disagreements
between groups simulating the delayed neutrino explosion mechanism are a
severe problem for this model.

The conclusion from the short discussion above is that the modelings of
SN2024ggi and SN2023ixf suggest that the explosion mechanism of their
progenitors was not the delayed neutrino explosion mechanism. Adding to this
energetic consideration the challenge of the delayed neutrino mechanism to
explain the point-symmetric morphologies of some CCSN remnants identified in
the last two years (e.g., Bear & Soker 2024; Soker 2024a, b, c), I conclude
that the most likely explosion mechanism of the progenitors of SN 2024ggi and
SN 2023ixf was the jittering jets explosion mechanism (JJEM). The JJEM can
also account for the polarization of CCSNe and its correlation with the
explosion energy (for the correlation, see Nagao et al. 2024b), as well as for
bipolar CCSNe (e.g., Nagao et al. 2024a).

##  4 Summary

The two recent relatively near CCSNe SN 2023ixf and SN 2024ggi allow good
fitting to their light curves and spectra, suggesting a dense compact CSM in
both cases. There are no indications of pre-explosion outbursts in either of
these CCSNe. This limits the possibility of the compact CSM being modeled with
a relatively long-lived CSM. Moriya & Singh (2024) and Chen et al. (2024a)
suggest this CSM be the wind acceleration zone of SN 2023ixf and SN 2024ggi,
respectively. In section 2 I concentrated on the recent model by Chen et al.
(2024a) and argued that the acceleration-wind zone model of the CSM requires
an unlikely large number of photon scattering events in the radial direction,
i.e., $\eta_{\rm w}\gg 1$. Instead, I proposed the effervescent zone model for
a pre-explosion dense compact CSM of SN 2024ggi, as I suggested for SN 2023ixf
(Soker, 2023), or the similar model of Fuller & Tsuna (2024). I encourage
future studies of ejecta-CSM interaction of CCSNe, including SN 2023ixf and SN
2024ggi, not to ignore the effervescent zone model and similar ones (e.g.,
Soker 2023; Fuller & Tsuna 2024).

Moriya & Singh (2024) and Chen et al. (2024a) also find their best fitting to
give explosion energy of $E_{\rm exp}\simeq 2\times 10^{51}{~{}\rm erg}$ for
both these CCSNe, and Jacobson-Galán et al. (2024) adopted an explosion energy
of $E_{\rm exp}\simeq 1.2\times 10^{51}{~{}\rm erg}$ for SN 2024ggi. In
section 3, I argued that these energies are a tough challenge to the delayed
neutrino explosion mechanism. The energetic challenge adds to the extremely
tough challenge of the delayed neutrino mechanism to explain point symmetric
CCSN remnants, as was pointed out in 2023-2024 (e.g., Soker 2024b). I
concluded in section 3 that the JJEM best explains the explosion energy of SN
2024ggi.

The existence of a massive effervescent zone requires strong pulsation and
envelope convection (e.g., Soker (2023), as was shown in detail by Fuller &
Tsuna (2024) for a similar CSM model (I use effervescent to include both). The
JJEM requires vigorous pre-explosion core convection (e.g., Papish & Soker
2011; Gilkis & Soker 2014, 2016; Shishkin & Soker 2021, 2022; Wang, Shishkin,
& Soker 2024). The delayed neutrino explosion mechanism also requires pre-
collapse convection to induce perturbations in the collapsing core that aid
explosion (e.g., Couch et al. 2015; Müller 2020; Bollig et al. 2021; Yamada et
al. 2024). I conclude that the correct modeling of convection, crucially in
the inner core but also important in the envelope, is mandatory to explore the
explosion and evolution of CCSNe.

## Acknowledgments

I thank Luc Dessart for his helpful comments. This research was supported by a
grant from the Pazy Foundation.

## References

  * Bear & Soker (2024) Bear E., Soker N., 2024, arXiv, arXiv:2403.07625.
  * Berger et al. (2023) Berger E., Keating G. K., Margutti R., Maeda K., Alexander K. D., Cendes Y., Eftekhari T., et al., 2023, arXiv, arXiv:2306.09311.
  * Bollig et al. (2021) Bollig R., Yadav N., Kresse D., Janka H.-T., Müller B., Heger A., 2021, ApJ, 915, 28.
  * Bostroem et al. (2023) Bostroem K. A., Pearson J., Shrestha M., Sand D. J., Valenti S., Jha S. W., Andrews J. E., et al., 2023, arXiv, arXiv:2306.10119.
  * Bruch et al. (2020) Bruch, R. J., Gal-Yam, A., Schulze, S., et al. 2020, arXiv:2008.09986
  * Burrows, Wang, & Vartanyan (2024) Burrows A., Wang T., Vartanyan D., 2024, ApJL, 964, L16.
  * Chen et al. (2024a) Chen T.-W., Yang S., Srivastav S., Moriya T. J., Smartt S. J., Rest S., Rest A., et al., 2024a, arXiv, arXiv:2406.09270.
  * Chen et al. (2024b) Chen X., Kumar B., Er X., Guo H., Yang Y.-P., Lin W., Fang Y., et al., 2024b, arXiv, arXiv:2405.07964.
  * Couch et al. (2015) Couch S. M., Chatzopoulos E., Arnett W. D., Timmes F. X., 2015, ApJL, 808, L21.
  * Dessart et al. (2017) Dessart, L., John Hillier, D., & Audit, E. 2017, A&A, 605, A83
  * Ertl et al. (2016) Ertl T., Janka H.-T., Woosley S. E., Sukhbold T., Ugliano M., 2016, ApJ, 818, 124.
  * Foley et al. (2007) Foley, R. J., Smith, N., Ganeshalingam, M., Li, W., Chornock, R., & Filippenko, A. V. 2007, ApJ, 657, L105
  * Fryer et al. (2012) Fryer C. L., Belczynski K., Wiktorowicz G., Dominik M., Kalogera V., Holz D. E., 2012, ApJ, 749, 91.
  * Fuller & Tsuna (2024) Fuller J., Tsuna D., 2024, arXiv, arXiv:2405.21049.
  * Gilkis & Soker (2014) Gilkis A., Soker N., 2014, MNRAS, 439, 4011.
  * Gilkis & Soker (2016) Gilkis A., Soker N., 2016, ApJ, 827, 40.
  * Gogilashvili, Murphy, & Mabanta (2021) Gogilashvili M., Murphy J. W., Mabanta Q., 2021, MNRAS, 500, 5393.
  * Grefenstette et al. (2023) Grefenstette B. W., Brightman M., Earnshaw H. P., Harrison F. A., Margutti R., 2023, arXiv, arXiv:2306.04827.
  * Jacobson-Galán et al. (2024) Jacobson-Galán W. V., Davis K. W., Kilpatrick C. D., Dessart L., Margutti R., Chornock R., Foley R. J., et al., 2024, arXiv, arXiv:2404.19006. doi:10.48550/arXiv.2404.19006
  * Jacobson-Galán et al. (2022) Jacobson-Galán W. V., Dessart L., Jones D. O., Margutti R., Coppejans D. L., Dimitriadis G., Foley R. J., et al., 2022, ApJ, 924, 15.
  * Jacobson-Galan et al. (2023) Jacobson-Galán W. V., Dessart L., Margutti R., Chornock R., Foley R. J., Kilpatrick C. D., Jones D. O., et al., 2023, ApJL, 954, L42.
  * Jencson et al. (2023) Jencson J. E., Pearson J., Beasor E. R., Lau R. M., Andrews J. E., Bostroem K. A., Dong Y., et al., 2023, arXiv, arXiv:2306.08678.
  * Kilpatrick et al. (2023) Kilpatrick C. D., Foley R. J., Jacobson-Galán W. V., Piro A. L., Smartt S. J., Drout M. R., Gagliano A., et al., 2023, arXiv, arXiv:2306.04722.
  * Margutti et al. (2014) Margutti, R., Milisavljevic, D., Soderberg, A. M., et al. 2014, ApJ, 780, 21
  * Mezzacappa et al. (2015) Mezzacappa A., Bruenn S. W., Lentz E. J., Hix W. R., Harris J. A., Bronson Messer O. E., Endeve E., et al., 2015, arXiv, arXiv:1501.01688.
  * Moriya et al. (2018) Moriya T. J., Förster F., Yoon S.-C., Gräfener G., Blinnikov S. I., 2018, MNRAS, 476, 2840.
  * Moriya & Singh (2024) Moriya T. J., Singh A., 2024, arXiv, arXiv:2406.00928.
  * Moriya et al. (2017) Moriya, T. J., Yoon, S.-C., Gräfener, G., & Blinnikov, S. I. 2017, MNRAS, 469, L108
  * Müller (2020) Müller B., 2020, LRCA, 6, 3.
  * Nagao et al. (2024a) Nagao T., Maeda K., Mattila S., Kuncarayakti H., Kawabata M., Taguchi K., Nakaoka T., et al., 2024a, arXiv, arXiv:2406.14057.
  * Nagao et al. (2024b) Nagao T., Patat F., Cikota A., Baade D., Mattila S., Kotak R., Kuncarayakti H., et al., 2024b, A&A, 681, A11.
  * Nakamura et al. (2024) Nakamura K., Takiwaki T., Matsumoto J., Kotake K., 2024, arXiv, arXiv:2405.08367.
  * Neustadt, Kochanek, & Rizzo Smith (2023) Neustadt J. M. M., Kochanek C. S., Rizzo Smith M., 2023, arXiv, arXiv:2306.06162.
  * Ofek et al. (2014) Ofek, E. O., Sullivan, M., Shaviv, N. J., et al. 2014, ApJ, 789, 104
  * Papish & Soker (2011) Papish O., Soker N., 2011, MNRAS, 416, 1697.
  * Pastorello et al. (2007) Pastorello, A., Smartt, S. J., Mattila, S., et al. 2007, Nature, 447, 829
  * Pessi et al. (2024) Pessi T., Cartier R., Hueichapan E., de Brito Silva D., Prieto J. L., Muñoz R. R., Medina G. E., et al., 2024, arXiv, arXiv:2405.02274. doi:10.48550/arXiv.2405.02274
  * Prentice et al. (2020) Prentice, S. J., Maguire, K., Boian, I., et al. 2020, arXiv:2009.10509
  * Shishkin & Soker (2021) Shishkin D., Soker N., 2021, MNRAS, 508, L43.
  * Shishkin & Soker (2022) Shishkin D., Soker N., 2022, MNRAS, 513, 4224.
  * Shrestha et al. (2024) Shrestha M., Bostroem K. A., Sand D. J., Hosseinzadeh G., Andrews J. E., Dong Y., Hoang E., et al., 2024, arXiv, arXiv:2405.18490. doi:10.48550/arXiv.2405.18490
  * Singh Teja et al. (2023) Singh Teja R., Singh A., Dutta A., Basu J., Anupama G. C., Sahu D. K., Swain V., et al., 2023, arXiv, arXiv:2306.10284
  * Smith et al. (2010) Smith, N., Miller, A., Li, W., et al. 2010, AJ, 139, 1451
  * Smith et al. (2023) Smith N., Pearson J., Sand D. J., Ilyin I., Bostroem K. A., Hosseinzadeh G., Shrestha M., 2023, arXiv, arXiv:2306.07964.
  * Soker (2021) Soker N., 2021, ApJ, 906, 1.
  * Soker (2023) Soker N., 2023, RAA, 23, 081002.
  * Soker (2024a) Soker N., 2024a, NewA, 107, 102154.
  * Soker (2024b) Soker N., 2024b, OJAp, 7, 31.
  * Soker (2024c) Soker N., 2024c, arXiv, arXiv:2403.08544.
  * Soraisam et al. (2023) Soraisam M. D., Szalai T., Van Dyk S. D., Andrews J. E., Srinivasan S., Chun S.-H., Matheson T., et al., 2023, arXiv, arXiv:2306.10783
  * Strotjohann et al. (2021) Strotjohann, N. L., Ofek, E. O., Gal-Yam, A., et al. 2021, arXiv:2010.11196
  * Sukhbold et al. (2016) Sukhbold T., Ertl T., Woosley S. E., Brown J. M., Janka H.-T., 2016, ApJ, 821, 38.
  * Svirski & Nakar (2014) Svirski, G., & Nakar, E. 2014, ApJ, 788, L14
  * Tartaglia et al. (2016) Tartaglia, L., Pastorello, A., Sullivan, M., et al. 2016, MNRAS, 459, 1039
  * Wang, Shishkin, & Soker (2024) Wang N. Y. N., Shishkin D., Soker N., 2024, arXiv:2401.06652.
  * Wang, Wang, & Dai (2019) Wang S.-Q., Wang L.-J., Dai Z.-G., 2019, RAA, 19, 063.
  * Xiang et al. (2024) Xiang D., Mo J., Wang X., Wang L., Zhang J., Lin H., Chen L., et al., 2024, arXiv, arXiv:2405.07699. doi:10.48550/arXiv.2405.07699
  * Yamada et al. (2024) Yamada S., Nagakura H., Akaho R., Harada A., Furusawa S., Iwakami W., Okawa H., et al., 2024, PJAB, 100, 190.
  * Yaron et al. (2017) Yaron, O., Perley, D. A., Gal-Yam, A., et al. 2017, Nature Physics, 13, 510
  * Zhang et al. (2024) Zhang J., Dessart L., Wang X., Zhai Q., Yang Y., Li L., Lin H., et al., 2024, arXiv, arXiv:2406.07806. doi:10.48550/arXiv.2406.07806