# Reinforcement Learning in Ultracold Atom Experiments

Malte Reinschmidt  József Fortágh  Andreas Günther  a.guenther@uni-
tuebingen.de Center for Quantum Science, Physikalisches Institut, Eberhard
Karls Universität Tübingen, Auf der Morgenstelle 14, D-72076 Tübingen, Germany
Valentin Volchkov<EMAIL_ADDRESS>Max Planck Institute for
Intelligent Systems, Max-Planck-Ring 4, D-72072 Tübingen, Germany

###### Abstract

Cold atom traps are at the heart of many quantum applications in science and
technology. The preparation and control of atomic clouds involves complex
optimization processes, that could be supported and accelerated by machine
learning. In this work, we introduce reinforcement learning to cold atom
experiments and demonstrate a flexible and adaptive approach to control a
magneto-optical trap. Instead of following a set of predetermined rules to
accomplish a specific task, the objectives are defined by a reward function.
This approach not only optimizes the cooling of atoms just as an
experimentalist would do, but also enables new operational modes such as the
preparation of pre-defined numbers of atoms in a cloud. The machine control is
trained to be robust against external perturbations and able to react to
situations not seen during the training. Finally, we show that the time
consuming training can be performed in-silico using a generic simulation and
demonstrate successful transfer to the real world experiment.

Laser cooling and magneto-optical trapping of neutral atoms has been one of
the major breakthroughs in science in the last decades [1]. This enabled
access to ultracold atoms and quantum matter with many applications in quantum
simulation, computation and sensing [2]. Narrow-line cooling lies at the heart
of a new generation of magneto-optical traps (MOT) for non-alkali elements and
molecules [3], leading to novel quantum matter [4, 5] and vastly improving
optical clocks [6]. Nowadays, the MOT is the cornerstone for all applications
and devices based on ultracold atoms. Despite its simple design, operation can
be quite intricate, with the behavior of the atoms being determined by the
interplay of laser cooling, optical pumping, and light-induced interactions.
For example, precise temporal control of laser parameters and magnetic fields
are required for implementation of cooling schemes such as compressed MOT [7],
molasses [8], or Raman cooling [9]. With the advent of non-alkali and
molecular MOTs also the number of lasers used is rising. This increasing
complexity can create difficulties in controlling and optimizing the
production of ultracold samples; however, it can also provide major
advantages. The key challenge is to handle and make use of this complexity.
One promising approach is the application of machine learning: methods like
Bayesian optimization [10] have been used to improve control sequences of
different phases of an experiment [11], even creating a BEC [12, 13, 14, 15].
However, these approaches do not allow for direct feedback during a sequence
and cannot react to perturbations and drifts of external conditions.

Figure 1: Reinforcement Learning for MOT. a Schematic illustration of the
interaction between an RL agent and a MOT-system. The agent takes control
(action) via dedicated control parameters, which it freely adjusts during the
MOT operation. At each time step $n_{t}$, the action is based upon the
observables (e.g. fluorescence images) given to the agent and the knowledge
acquired during the training phase. The latter is based on a reward value $R$,
extracted at the end of the MOT cycle (episode) from an absorption image of
the atom cloud, which the agent aims to maximize. b Evolution of reward $R$
(normalized to reference sequences, as described in the Methods) during a
successful training. The insets show the control parameter sequences at four
points during the training phase. Animated fluorescence images corresponding
to the episodes in the insets are provided in Supplementary Movie 1. c
Sequence of control parameter ($\Delta_{\mathrm{C}}$) over the course of an
episode for an agent trained to maximize $N/T$. The horizontal line indicates
the value for maximal loading rate of the MOT. The inset shows the
fluorescence images observed by the agent at various time steps $n_{t}$ during
the episode.

In the last decade, reinforcement learning (RL) [16] has shown human-level
performances in a range of tasks that were previously thought to be impossible
for machines to complete. While impressive success was achieved in games and
simulated environments [17, 18, 19, 20, 21], application of RL to real world
systems remained challenging [22]. In recent years, control of various
physical systems was established through RL algorithms, leading to a paradigm
shift in how control is approached and implemented [23, 24, 25, 26, 27].
Generally speaking, in RL the control algorithm denoted as agent interacts
with an environment based on observations. As a result of its actions the
agent receives rewards, which are used in a training algorithm, allowing the
agent over time to learn how to maximize the cumulative reward. The task is
defined directly through the reward function, leaving the specific control
sequence solely to the agent. So-called deep RL [28] incorporates deep
artificial neural networks allowing for model-free control as a function of
many, potentially high dimensional inputs, such as images.

Here, we introduce deep RL to ultracold atom experiments by means of
controlling a MOT-system. We enhance an existing cold atom apparatus used to
cool rubidium (${}^{87}\mathrm{Rb}$) atoms by an RL agent, taking advantage of
the fluorescence light, which is a by-product of laser cooling. In most
similar apparatuses, fluorescence images are only a basic diagnostic tool,
whereas they can contain a complex and non-linear representation of the system
and are thus an ideal candidate for input to an RL agent. We use an off-the-
shelf implementation [29] of a popular RL algorithm with minimal modifications
to the code, running on a desktop computer. As proof-of-concept, we perform
training of the agent directly on the experiment, which rapidly converges to a
cooling scheme resembling well known molasses cooling. Furthermore, we train
the agent to be capable to react to sudden changes and observe generalization
to situations that are not included during the training. We demonstrate new
modes of operations of the MOT, such as loading predefined number of atoms by
engineering the corresponding reward function. Finally, we explore how an
agent trained on a simple, neural network-based MOT simulator can be
transferred to the real world experiment, opening the avenue to fast pre-
training, extensive parameter optimization and large number of control
parameters.

## Setup

### Experimental System

The demonstration of RL-based MOT control is implemented on a cold atom
apparatus, that employs many standard techniques (see Methods). Similar setups
can be found in many labs around the world [30]. Like all experiments with
cold atoms, these experiments are performed in recurring cycles, each starting
with the preparation of an atomic cloud in a MOT. We dedicate the present
manuscript to this preparation phase and conclude the experimental cycle
directly thereafter by detecting the collected atoms. For that, we use time-
of-flight absorption imaging [31], in which the expanding cold atom cloud is
illuminated with a resonant laser beam and the beam transmission measured on a
CCD camera. This image is used to determine the final state of the atomic
cloud, consisting of the number of atoms $N$ and their average kinetic energy,
which we denote by ’temperature’ $T$. In our specific system, after a loading
time of 1.5 s, the number of atoms typically amounts to $N\approx 1.5\times
10^{8}$ at a temperature of about 100 $\mu K$.

### Implementation of the control

While most of the experimental control is running on a separate system, the
agent algorithm is implemented on a desktop PC with access to fluorescence
images and control over specific experimental parameters (see Methods). We
limit the parameters controlled by the agent to the frequency detuning
$\Delta=\omega_{0}-\omega_{L}$ of the cooling laser at frequency $\omega_{L}$
with respect to the atomic resonance at $\omega_{0}$ and the magnetic field
gradient $B^{\prime}$ at the MOT center along the field’s symmetry axis (see
Methods). They are among the most influential parameters, both on the number
of atoms and the temperature in the MOT, and are typically varied over the
course of the MOT loading. In this study, the RL agent is allowed to
continuously adjust one or both of these parameters. The detuning $\Delta$ is
restricted to the range of 0 –8.25 Γ with $\Gamma$=2$\pi\times$6.063 MHz being
the cooling transition linewidth [32]. The magnetic field gradient
$B^{\prime}$ can be varied by the agent between 7.5 –22.5 G/cm.

Each MOT cycle of 1.5 s duration constitutes an episode and is divided into 25
sequential time steps $n_{t}$ of 60 ms length. Figure 1a illustrates the main
principle for controlling the MOT with an RL agent. The MOT is part of the
environment, which denotes the entire physical system, including the atoms,
the surrounding experimental set-up and the lab. Within each time step a
vector of observables consisting of fluorescence images and additional
parameters such as current time step index and values of control parameters is
passed as input to the agent. The agent processes the input using a trainable
neural network and determines the next action, i.e. computes the output values
for the controlled experimental parameters. Subsequently, the control
parameters are applied through the experimental control and the environment
evolves for the duration of the time step.

### Training of the agent

In order to train the agent for a specific task, the reward function needs to
be defined accordingly. In general, a reward can be generated after each time
step, in our implementation a reward is attributed only after the last time
step when the state of the atomic cloud has been obtained through absorption
imaging. The essence of RL is to train the agent’s policy to maximize the
reward. The RL algorithm used here is the deep deterministic policy gradient
(DDPG) [33], which has been developed for high dimensional observations (e.g.
images) and is capable of outputting continuous actions, i.e. control
parameters. The training takes place in the course of a number of episodes,
during which the policy is continuously updated according to incoming rewards.
During the training, exploration is introduced via noise added to the actions
of the agent. A single training consists of up to $10^{5}$ episodes, depending
on the number of control parameters and the specific definition of the reward
function. Figure 1b shows a successful training, in the course of which the
obtained reward increases and saturates, while the output values (control
parameters) converge to a stable sequence (see insets). Figure 1c shows a
typical control sequence within one episode, together with representative
fluorescence images.

A well known problem of learning algorithms is that some inputs may be ignored
and instead ”shortcuts” are memorized [34]. In our scenario, this approach is
particularly tempting as the experimental setups are typically very stable and
run reproducibly, thus operation with static control sequences is very common.
To prevent this we introduce a training perturbation, i.e. an offset, that is
randomly generated at the start of each episode and added to the agent’s
output value of the control parameter, effectively shifting the experimental
conditions. Such a perturbation is equivalent to a long-term drift of the
control parameter that may occur in experiments due to changes in the
environment (e.g. lab temperature). As the offset value is unknown to the
agent, it is forced to extract information from the fluorescence images in
order to compute the control parameter. This approach is similar to domain
randomization that induces variability during training on simulations, making
the agent more robust in real world environments [35].

Figure 2: Evaluation of agent trained to optimize $N/T$. a Control parameter
sequence over the course of episodes with mean value (colored solid line) and
standard deviation (transparent region). Data are obtained by an agent
evaluated at $\Delta_{\mathrm{Offset}}=0$. The horizontal black line indicates
the point of optimal loading $\Delta_{\text{opt}}$. b Same agent as in (a),
but evaluated at an offset $\Delta_{\mathrm{Offset}}=0.66\,\mathrm{\Gamma}$.
The optimal loading line has been shifted accordingly. c Same agent as in (b),
but evaluated with a sudden step in the $\Delta_{\mathrm{Offset}}$ at
$n_{t}=15$ with amplitude 0.82 $\Gamma$. d Control parameter sequence over the
course of episodes for an agent controlling two output parameters,
corresponding to the cooling laser detuning $\Delta_{\mathrm{C}}$ (blue) and
the magnetic field gradient $B^{\prime}$ (red). e Reward values (mean and
standard deviation) of three different agents as extracted during the
evaluation for different $\Delta_{\mathrm{Offset}}$ values:
${\color[rgb]{0,0,1}\circ}$/${\color[rgb]{0,0,1}\bullet{}}$ \- agent
controlling perturbed detuning ($\Delta_{\mathrm{Offset}}\neq 0$ during
training); $\blacklozenge{}$ \- agent controlling detuning
($\Delta_{\mathrm{Offset}}=0$ during training);
${\color[rgb]{1,0,0}\blacksquare{}}$ \- agent controlling $B^{\prime}$ and
perturbed detuning. The ${\color[rgb]{0,0,1}\bullet{}}$ measurements include a
sudden change in the $\Delta_{\mathrm{Offset}}$ value at $n_{t}=15$ during
evaluation, as illustrated in c. The timing and amplitude of this ”step” in
the offset was the same for all detuning offsets.

### Evaluating the agents performance

In order to evaluate and compare the performance of trained agents, we define
a testing procedure, aiming for the agent’s ability to adapt to different
experimental conditions, as well as the reproducibility of the results when
presented with the same conditions. As such we test the agent by adding
different predetermined perturbation offset values (unknown to the agent) to
the control parameter at the start of an episode, repeating the test for each
of these values up to 200 times. We then compare the final results in terms of
the number of collected atoms, the achieved cooling, and the corresponding
reward. Additionally, we study to which policy the agent has converged by
analyzing the sequences and spread of the output values over the course of
multiple episodes for each of the different offset values [36].

## Results

### Standard MOT operation mode

In many cold atom experiments MOT operation aims for maximizing the number of
trapped atoms while reaching the lowest possible temperature, as it is one of
the keys towards quantum degeneracy. For a proof-of-principle demonstration we
define the reward as $R\propto N/T$ and focus first on a single control
parameter – the laser detuning $\Delta$. The physical value of the laser
detuning $\Delta$ is determined through the agent’s control value
$\Delta_{\mathrm{C}}$ and the perturbation offset $\Delta_{\mathrm{Offset}}$
(unknown to the agent) via
$\Delta=\Delta_{\mathrm{C}}-\Delta_{\mathrm{Offset}}$. The magnetic field
gradient is set to $B^{\prime}=15\,\mathrm{G/cm}$.

With the MOT principles being based on velocity dependent light forces, the
loading rate of atoms strongly depends on laser detuning and shows a maximum
at a specific value. In our experiment we measured this value to
$\Delta_{\mathrm{opt}}=$ 1.9 $\Gamma$. The temperature is known to
approximately scale inversely proportional to $\Delta$ [37]. The final
temperature is then mainly given by the value of $\Delta$ applied at the end
of the episode. With this in mind, the intuitive optimal policy consists of
first going to the point of highest loading rate, before rapidly ramping
$\Delta$ to the maximum value at the end of the episode to reduce the
temperature. This policy corresponds to what is commonly used in labs working
with MOTs and provides a reference for studying the application of RL to MOT
control (see Methods).

In Fig. 1b, we present the learning progress: the agent discovers early on
that going to large detunings at the end of an episode results in a massive
boost in the reward, as witnessed in the first insets. This is consistent with
the fact that reward is issued only in the last time step. During the second
half of the training, the agent fine-tunes the loading phase. The control
output sequences $\Delta_{\mathrm{C}}$ for a trained agent evaluated at two
different detuning offsets are shown in Fig. 2a and b. We observe, that the
sequence of the detuning $\Delta_{\mathrm{C}}$ over the course of an episode
for a trained agent matches our expected intuitive policy. While the first
fluorescence images of each episode are blank, as no atoms have yet been
loaded into the MOT, the agent develops an initial search pattern, jumping to
different values of $\Delta_{\mathrm{C}}$. Once the first atoms are trapped
and visible on the fluorescence image, the agents compensates for the
perturbation offset $\Delta_{\mathrm{Offset}}$ and sets the output value
corresponding to the highest loading rate, as shown in Fig. 2b. This allows
the agent to load as many atoms into the MOT as possible, before increasing
$\Delta_{\mathrm{C}}$ towards the end of the episode in order to reduce the
temperature and maximize the reward. This control policy is remarkably stable,
as seen in the small standard deviations in Fig. 2a and b. Fig. 2e summarizes
the performance in terms of the achieved rewards for several trained agents,
evaluated at different offsets. At large offsets, the drop in the obtained
reward is due to an increased number of time steps needed to locate the
optimal loading point.

For sake of consistency, we also trained an agent without perturbations on the
control values. This agent learns the optimal policy, however, it ignores the
changes related to $\Delta_{\mathrm{Offset}}$ in the observations during the
evaluation phase and plays back the same memorized action sequence. The reward
then simply mirrors the loading rate of the MOT at different offsets from the
optimal loading point, as can be seen in Fig. 2e.

### Generalization to unfamiliar situations

The capacity to adequately react to new situations, unseen during the
training, is one of the key interests in RL since it is a witness of the
robustness and usefulness of the trained agent. Thus, we tested our agent by
exposing it to sudden changes in the environment and observing the real-time
adaption of its behavior. In particular, the offset applied to the output
value $\Delta_{\mathrm{C}}$ during evaluation was suddenly increased at
$n_{t}=15$ by 0.82 $\Gamma$. Such a shift can be compared to a sudden change
in the environment, which may be caused by singular acoustic or electronic
events. As shown in the control parameter sequence of Fig. 2c, the agent was
able to compensate for this shift after just a single step, i.e. based on a
single fluorescence image. This demonstrates the ability of the trained agent
to generalize and adapt to conditions not present during the training. Note,
that the corresponding reward values (${\color[rgb]{0,0,1}\bullet{}}$ in Fig.
2e) cannot be directly compared to those without shift
(${\color[rgb]{0,0,1}\circ}$ in Fig. 2e), as for a short interval of time,
before the agent was given a chance to react, the detuning was being shifted
close enough to the resonance where loading not only stops, but atoms were
actively removed from the trap. As a result the total collected number of
atoms was reduced, leading to lower reward values.

Figure 3: Atom number on demand. a Evaluation for three different offsets in
terms of the targeted ($N_{\mathrm{T}}$) and achieved ($N$) number of atoms
for an agent trained to load an adjustable amount of atoms. The inset shows
the reward function, defined as gaussian profile around $N_{\mathrm{T}}$, with
standard deviation $\sigma_{N}=15\times 10^{6}$. The black line and grey area
correspond to the response $N=N_{\mathrm{T}}\pm\sigma_{N}$, with mean and
standard deviation, respectively. b Control parameter sequence for
$\Delta_{\text{offset}}=0$ and different vales of $N_{\mathrm{T}}$, ranging
from $5\times 10^{6}$ (red) to $110\times 10^{6}$ (blue). Qualitatively
similar sequences were observed for non-zero offsets.

### Increasing action space dimension

So far, the action space has been limited to a single control parameter. In
the following we expanded the action space to also include the current that
defined the gradient of the magnetic field, i.e. $B^{\prime}$, while the
definition of the reward as $R\propto N/T$ and other aspects during the
training remained unchanged. The perturbation offset applied to $\Delta$ was
found to be sufficient in preventing the learning of shortcuts, so no offset
was added to $B^{\prime}$.

With the controls extended to the gradient of the magnetic field, schemes
similar to compressed MOT [7] or polarization gradient cooling [37] could
expected to be discovered by the agent. However, these techniques require
control on time scales much shorter than the time step duration of 60 ms. As
consequence, no intuitive optimal policy is known for this scenario.

We evaluate the trained agent, similar as before, with the resulting control
sequences for zero offset shown in Fig. 2d. The control of the detuning
$\Delta$ is very similar to before. For $B^{\prime}$, the agent initially
maximizes this control parameter. Once the point of highest loading rate is
found, $B^{\prime}$ is lowered. In a manner reminiscent of a compressed MOT,
$B^{\prime}$ is increased to maximum at the end of the episode. In terms of
performance, one could expect an improvement due to the additional control
parameter. Indeed Fig. 2e shows generally better results for this agent.
Interestingly, for $\Delta_{\mathrm{Offset}}=-0.33\Gamma,-0.66\Gamma$ the
agent outperformed the reference-sequences in which $B^{\prime}$ was kept
constant. This suggests that MOT operation can be further improved, provided
the agents is granted control over further experimental parameters.

Figure 4: Sim-to-real transfer. a Schematic diagram of the training on the
simulation: the graphs illustrate the experimentally measured values for the
number of atoms $dN$ loaded during a time step and the temperature $T$ as
function of the detuning $\Delta$. At the bottom the generation of
fluorescence images is illustrated: a convolutional neural network takes the
number of atoms $N$ and $\Delta$ as input to produce realistically looking
images (see Methods). b Control parameter sequences of an agent trained on the
simulation and evaluated for $\Delta_{\text{offset}}=0.66\Gamma$ on the
simulation (left) and experiment (right). The red curves show the resulting
behavior if the offset is suddenly changed during the episode (cf. Fig. 2c).
The solid/dashed black line indicate the point of optimal loading with/without
the sudden change of the offset. c Reward values of agents trained on
simulation and evaluated on the experiment.
${\color[rgb]{0,0,1}\bullet{}}$/${\color[rgb]{1,0,0}\bullet{}}$ — reward
without/with sudden change of the offset.

### Reward function engineering

Now we go beyond the standard task of maximizing $N/T$ and showcase how the
engineering of the reward function enables the realization of new modes of
operation. As an example, one important task is the preparation of a sample
with an adjustable number of atoms.

Using RL, we were able to train an agent to achieve variable targeted atom
numbers in each episode, while acting only on the detuning of the cooling
laser and keeping all other parameters fixed ($B^{\prime}=15\,\mathrm{G/cm}$).
For this, we introduce a new parameter $N_{\mathrm{T}}$ for the target number
of atoms and pass it to the agent as additional observable. The reward
function is altered to depend on this parameter, using a Gaussian profile
centered around $N_{\mathrm{T}}$ with a width $\sigma_{N}$. To ensure that the
atoms remain cold, we include a factor $1/T$ in the reward function:

$R\propto\frac{1}{T}\exp{\left(-\frac{(N-N_{\mathrm{T}})^{2}}{2\sigma_{N}^{2}}\right)}$
(1)

For the training we limit the range of $N_{\mathrm{T}}$ to 0–120 million atoms
and set the acceptance range $\sigma_{N}$ equal to 15 million atoms.
$N_{\mathrm{T}}$ is randomly generated for each training episode, with a
uniform distribution across the entire range. As before, a random perturbation
offset was added to the output of the agent at each training episode to
prevent the learning of shortcuts. The agent must therefore learn to
manipulate the loading rate to achieve the requested atom number, while
simultaneously compensating the offset.

Intuitively, one possible control strategy could be to divide the episode in
two phases: a loading phase with the detuning remaining at the optimal loading
point but its length adjusted to meet the requested atom number, followed by a
cooling phase at larger detunings. Alternatively, the loading could take place
at variable loading rate but constant length and terminated with a jump to
maximal detuning for cooling at the episode’s end.

We evaluated the trained agent at given $N_{\mathrm{T}}$ values, ranging from
5–110 million atoms. When comparing the measured number of atoms to the target
values in Fig. 3a, we find, that they fall roughly within a region of
$\pm\sigma_{N}$ around $N_{\mathrm{T}}$. The agent is thus capable of adapting
to different number of atoms in consecutive episodes. More importantly, the
performance of the agent is robust with respect to the applied perturbation
during the evaluation, as demonstrated by the data obtained for non-zero
$\Delta_{\mathrm{Offset}}$. This robustness is a direct consequence of the
real-time feedback of the RL agent and would be impossible to obtain for
manual and static control sequences.

The control parameter sequences presented in Fig. 3b were applied by the agent
for $\Delta_{\mathrm{Offset}}=0$. They show that the agent varies both: the
loading rate (i.e. distance to $\Delta_{\mathrm{opt}}$) as well as the
duration of the cooling phase (i.e. jump to larger detuning) depending on
$N_{\mathrm{T}}$. Furthermore, for increasing values of $N_{\mathrm{T}}$ the
agent limits the cooling in a trade-off for larger number of atoms. Indeed,
increasing the detuning reduces the confinement of the MOT and leads to an
increased loss of atoms. Due to the absolute nature of the acceptance range
$\sigma_{N}$ in the reward function (1), the loss of atoms affects the reward
only for larger $N_{\mathrm{T}}$. In principle, the agent could collect a
larger number of atoms in order to compensate for this loss and to better meet
the target by staying closer to $\Delta_{\mathrm{opt}}$ during the loading
phase. This shows that despite its overall success, the policy of the trained
agent is not yet optimal.

### Sim-to-real transfer

As with other implementations of RL in real world environments, a major
challenge is the required training time. Depending on the complexity of the
reward function and the number of parameters, the training performed for this
work took 12-48 h. Adding complexity to the reward function, or expanding the
system to include more parameters increases this time.

A typical approach to this problem in RL for robotic systems is sim-to-real
transfer. If given a simulation sufficiently close to the environment, general
behavior can be studied and hyperparameters of the algorithm being optimized.
Ideally, the algorithm trained on the simulation could then be transferred to
the experiment, either directly, or as a basis for an additional shorter
training to adapt to the real world.

To investigate sim-to-real transfer, we have implemented a simple
phenomenological model of the MOT, as depicted in Fig. 4a and detailed in the
Methods section. We modeled the internal state of the MOT assuming that both,
the number of atoms $dN$ loaded per time step as well as the temperature $T$
only depend on the current value of $\Delta$ and extracted their functional
relation, $dN(\Delta)$ and $T(\Delta)$, from the experiment. Furthermore, we
trained a generative neural network on experimental data to reproduce
fluorescence images corresponding to the simulated MOT state. As with the
training on the experiment, a perturbation offset was applied to the action in
order to force the agent to rely on the fluorescence images.

Figure 4b (left) shows the control sequence of an agent trained and evaluated
on the simulated environment, reproducing the behavior of experimentally
trained agents. We then evaluated the simulation-trained agent directly on the
real world experiment. The corresponding control presented in Fig. 4b (right)
qualitatively mimics the behavior of the agent evaluated on the simulation and
thereby indicates a successful sim-to-real transfer. More importantly, when
presented with a sudden shift of the offset value during the episode, the sim-
to-real agent adapted to the change in the same manner as the agent trained on
the experiment. The agent was therefore not only able to generalize from its
training on the simulations to the real world experiment, but was even able to
take this a step further and generalize to unknown situations that did not
occur during the training.

The absolute reward values achieved by the sim-to-real agent, shown in Fig.
4c, are inferior to those of the agents trained directly on the experiment.
The main reason for this is that the agent did not fully compensate for the
perturbation offset when evaluated on the experiment, as can be seen in Fig.
4b. This shows that there is still a considerable gap between the simulations
and the experiments. Possibly, one contributing factor is the smoothing of the
synthetic images, lacking high spatial frequency features of the fluorescence
images in the experiment (see Methods). Other factors may be the
simplifications in the simulation, which ignore light-induced interactions in
the MOT, multiple photon scattering and changing lab conditions. By addressing
these effects in the simulation, and improving the details in the synthetic
images, the sim-to-real gap may be further closed.

We publish the full code needed to run the simulation and encourage the
interested reader to try out the RL on a simulated MOT [38]. Furthermore, the
reward function can be easily modified to test new modes of operation of a MOT
on a computer. Even the phenomenological model can be retrained on custom
training data, allowing experimentalists to safely evaluate the possibilities
provided by RL without changes to their experimental setup.

## Discussion and Outlook

We have successfully applied RL control to ultracold atom experiments and
demonstrated its applicability by means of a MOT system. For this we have
implemented a well studied algorithm (DDPG), running on a desktop computer.
Our pioneering work will enable many labs around the world to adapt this
approach to their specific needs and extend their research. In detail we have
shown, that an RL agent for MOT control is able to find meaningful control
strategies during a training phase and is capable to react in real-time to
unexpected events, not seen during the training. The real-time feedback to the
agent through fluorescence images allows for operation under unstable
conditions and is a central distinction with respect to other machine
learning-optimized ultracold atom experiments. Providing additional
observations such as real-time values of the laser intensity and lab
temperature (that tend to fluctuate and drift over time) has the potential to
further improve the performance of such a system. In this work we have
facilitated the convergence of the training to comprehensible control
strategies by limiting the number of parameters and their range. In future
work, increasing the range, giving control to additional parameters (e.g.
laser intensities and polarizations) and reducing the duration of the time
steps could lead to discovery of new and possibly counter-intuitive operation
modes of a MOT [15]. At the same time, allowing the agent to control
additional parameters would enable the operation of MOTs for non-alkali
species and molecules that require multiple lasers and elaborate control
sequences. The policies discovered by RL agents would have the potential to
shed new light into the physics of the light-matter interaction, especially
when working with such complex systems.

In addition, we have demonstrated a paradigm change in ultracold atom physics
by using reward engineering to focus on the goal of a specific task, instead
of manually crafting control sequences. As an application of this principle,
we have trained an agent that is capable of producing samples with user
defined number of atoms. Conceivably, position dependent reward functions
could be used to realize complex spatio-temporal control of ultracold samples,
reminiscent to RL control of a plasma in a tokamak fusion generator [27]. The
promising results of transferring an agent trained on a simple simulation to
the real experiment opens multiple opportunities for further research: the
physical model used for the simulation can be improved by incorporating light-
induced losses, with the necessary parameters readily acquired on the
experiment. Domain randomization in the form of manipulation of synthetic
fluorescence images could further improve the sim-to-real transfer efficiency.

With the MOT being the workhorse in the field of ultracold atoms, it presents
a new, exciting platform for the development and testing of novel RL
algorithms. Indeed, RL algorithms such as soft actor-critic (SAC) [39] or
proximal policy optimization (PPO) [40] can be readily implemented and tested
in our simulation and brought to the experiment. To foster progress in this
new and exciting field, we made the code for training and evaluating an agent
based on a simulated MOT publicly available [38]. Finally, our work
constitutes the first step to the production of degenerate quantum gases using
RL control, provided that suitable real-time observables (e.g. ion signals
[41] or non destructive imaging methods [42, 43]) can be implemented.

## Methods

### MOT operation and implementation

MOT operation requires a quadrupole magnetic field and in the most common
configuration a pair of counter-propagating laser beams (cooling beams) in all
three spatial directions. Together, this makes for a position and velocity
dependant light force (radiation pressure), which cools and traps neutral
atoms from the gaseous phase at the center of the magnetic quadrupole field
where the beams are overlapping.

Our MOT system is based on ${}^{87}Rb$ atoms and described in [44]. The
quadrupole field is generated by a pair of coils in Anti-Helmholtz
configuration, providing a field gradient of $B^{\prime}=14\,\mathrm{G/cm}$ at
the quadrupole center, along the coil symmetry axis for about 1 A current in
the coils. The coil system is placed inside a UHV science chamber at a
pressure of $10^{-11}\,\mathrm{mbar}$ and operated via a voltage controlled
current source. The cooling beams with waists ($1/e^{2}$-beam radius) of 12 mm
and powers of about 15 mW are circularly polarized and overlap at the center
of the magnetic quadrupole between the coils. The MOT is loaded from a beam of
precooled Rb atoms emitted from a two-dimensional (2D)-MOT-source. The latter
consists of field coils and thermal rubidium emitter sources (dispensers),
providing a rubidium partial pressure of $>10^{-9}\,\mathrm{mbar}$ in the
corresponding 2D-MOT chamber. It is attached to the science chamber via a
differential pumping stage with a direct line of sight between the centers of
the two MOT systems. The 2D-MOT is loaded with rubidium atoms from the
background pressure and creates a transversally cooled beam of cold atoms
guided to the MOT in the science chamber.

As cooling and trapping relies on (position and velocity dependant) radiation
pressure due to light scattering, MOT operation requires a closed cooling
cycle transition for the atoms. That is a state pair, which the atoms cannot
leave during absorption or spontaneous emission. In ${}^{87}Rb$ such a
transition is found in the D2 manifold between the $5S_{1/2}$ ground state and
the $5P_{3/2}$ excited state, near 780 nm. Both states are split in two
respectively four hyperfine manifolds with total angular momentum $F=1,2$ and
$F^{\prime}=0,1,2,3$. We use the $5S_{1/2},F=2\rightarrow
5P_{3/2},F^{\prime}=3$ cooling transition, which is considered to be closed.
For MOT operation, the cooling light is red detuned to this transition. Atoms
can leave the otherwise closed cooling cycle only via off resonant excitation
to the $5P_{3/2},F^{\prime}=2$ state and adjacent decay to the F=1 ground
state. A repump laser resonant to the $5S_{1/2},F=1\rightarrow
5P_{3/2},F^{\prime}=2$ transition and overlapped with the cooling light is
used to bring back the atoms to the cooling cycle. The repump beam is
circularly polarised, has a power of 20 mW, and a beam waist equal to the
cooling beams. The 2D-MOT is operated with separate cooling and repump beams
at the same frequencies.

### Control Implementation

The control of the experimental apparatus is divided between the real-time
system Jäger Adwin Pro II and a desktop PC dedicated for the RL agent. While
the Adwin system is responsible for the time-of-flight imaging operations,
which require exact timing, the desktop PC controls the selected experimental
parameters ($\Delta$ and $B^{\prime}$) via a digital-to-analog converter (DAC)
card (AdLink DAQ-2502). The analog output line of the DAC card for the
detuning $\Delta$ controls directly the offset lock used for frequency
stabilization of the laser. $B^{\prime}$ is set via further analog output
lines that control the current in the MOT-coils.

An episode is initiated from the PC via enabling the magnetic field gradient,
effectively turning the MOT on. An additional output line from the DAC card to
the Adwin system is used to trigger the time-of-flight imaging sequence at the
end of an episode. At the conclusion of the imaging sequence, the Adwin system
returns into the idle state, ready to start loading atoms.

The fluorescence camera (UEye UI-3060CP) is connected to the PC and is
constantly recording images and transferring them into a reserved memory
register at a rate of 150 fps. During an episode, in each time step the latest
available image is read out from the memory and processed to be passed as part
of the observation to the agent.

### Fluorescence imaging

As atoms in the MOT scatter laser light during the cooling and trapping
process, the scattered light can be observed as fluorescence. This phenomenon
makes the atomic cloud in the MOT observable, with the spatial distribution of
atoms being reflected in the appearance of the fluorescence. Typically, the
extension of the fluorescing cloud is on the order of millimeters. The rate at
which fluorescence is emitted depends on the laser detuning, and is modulated
by the magnetic field gradient. This results in a highly non-linear mapping of
the 3-dimensional atomic density distribution onto the imaging plane. The
images then correspond to a projection along the imaging axis and show the
fluorescence intensity in the direction of the camera. Different processes
such as radiation trapping or imbalanced radiation pressure can lead to a
variety of non-trivial shapes of the atomic cloud, including ring structures
[45].
The fluorescence images have a dimension of 50x50 pixels with 8-bit
resolution. The exposure settings are chosen such that even weak signals can
be detected, which is important for the initial phase of the learning. This
leads to fast saturation when loading at maximum rate.

We have conducted tests in the simulated environment providing the agent only
with integrated fluorescence signals instead of the full image as observable.
The training in such a scenario did not converge to the optimal policy, in
particular, the agent was not able to adapt to any perturbation.

### Reference episodes

The reference episodes, in which the control sequence was manually set to the
intuitive optimal policy, i.e. loading at optimal detuning and ramping the
detuning to the maximal value in the last time step, allowed us to normalize
the evolution of the performance (i.e. reward $R$), correcting for possible
long term drifts of the experimental apparatus. In practice it means that in
each time step of a reference episode the actions of the agent were
overwritten with predefined values while keeping the timing constraints of the
agent, i.e. the 60 ms step size. The reference episodes were included during
training (see below) and when evaluating a trained agent.

Figure 5: Experimental fluorescence images from the test dataset (upper line)
and corresponding synthetic images (lower line) for different parameter sets
($N$, $\Delta$).

### RL Algorithm

For the realization of an RL-controlled MOT we have taken recourse to a
modular framework known as coach [29], containing implementations of a variety
of RL algorithms in Python. In this work, we have focused on the DDPG
algorithm, an actor-critic algorithm, that uses memory replay for training and
is able to learn directly from raw images [33]. We have found that almost all
parameters are close to the optimum at their default values (as defined in
coach). The parameters used for the simulation can be found in our preset
published as a part of the code [38].

Inspired by the performance of the Deep Q-Network [17] we have adapted the
stacking of the 4 most recent fluorescence images as input to the agent. This
allows the agent to compute the rate of change of the fluorescence signal,
which can be used to determine the loading rate of the MOT. All values passed
as observables are normalized to the range [0,1].

### Scheduling

The learning starts with a heat-up phase of 400 episodes during which the
control is fully determined by a Brownian motion-like exploration. During the
subsequent training phase, every 200 training episodes, a set of 10 reference
and 10 evaluation episodes with random perturbation offsets are recorded. The
performance of the agent is evaluated by disabling the exploration during the
evaluation episodes. During the training episodes, unlike in the original DDPG
implementation, we compute a new action in every time step (i.e. no action
repeat).

### MOT simulation

We have implemented a phenomenological, data-driven computer simulation for
the MOT in order to better understand the behavior of the agent in such a
novel environment, as well as to investigate sim-to-real transfer.

For our simulation we assume that the number of trapped atoms increases by a
certain number $dN$ in each time step, where $dN$ only depends on the applied
control parameter $\Delta_{n_{t}}$ at time step $n_{t}$.

$N(n_{t})=\sum_{i=1}^{n_{t}}dN\left(\Delta_{i}\right)\\\ $ (2)

This means that we neglect loss processes that are known for MOTs and
typically lead to saturation of the accumulated number of atoms. Furthermore,
we assume that at each time step the temperature of the trapped cloud is
dictated by the cooling laser light and thus only depends on the value of the
detuning.

$T(n_{t})=T\left(\Delta_{n_{t}}\right)$ (3)

This assumption is valid as long as the atomic cloud stays dilute and
transparent to the light. Thus, we neglect density dependent multiple
scattering. In order to obtain the most accurate simulation, we experimentally
determine the loading rate $dN$ as well as the temperature $T$ for a discrete
number of detuning values and create look-up-tables (LUT), which are used for
interpolations during the simulation (cf. Fig. 4a). The loading rate as
function of the detuning is measured by setting the desired value $\Delta$ for
a duration of 1 s and measuring the final number of atoms via absorption
imaging. The temperature LUT is recorded by loading atoms at the optimal
loading rate and subsequently rapidly ramping the detuning to a given value
before performing absorption imaging. The fluorescence images that are
captured right before the measurement of $dN$ and $T$ via absorption imaging
are stored in a database. Each measurement is repeated 5 times.

The full dataset contains on the order of 5000 samples. The environment based
on the simulation computes the intermediate values of $dN$ using a linear
interpolation. The temperature $T$ is obtained by interpolating between the
measured values using an exponential function. During the training on the
simulation, noise is introduced to the loading rate $dN$ in order to simulate
the physical fluctuation of the experimental conditions.

$dN\rightarrow dN\cdot\mathcal{N}\left(\mu,\sigma^{2}\right)$ (4)

with $\mathcal{N}$ being the normal distribution at mean value $\mu=1$ and
standard deviation $\sigma=0.1$.

We trained a convolutional neural network (CNN) to generate realistically
looking fluorescence images for the simulation on the same experimental
dataset. We have implemented the decoder-part of an autoencoder network [46,
47] as the generator and performed supervised training with $N$ and $\Delta$
as inputs and the fluorescence image as output. The network consists of a
fully connected layer, followed by 5 convolutional layers with rectified
linear unit as activation function and upsampling layers in between. Prior to
the training we set fluorescence images to 0 for samples with $N=0$, as well
as remove outliers from the dataset. A small, randomly chosen subset of the
dataset was held out during training and used for testing the generator. In
Fig. 5, we compare the recorded images to the output of the generator for
several parameter sets. The overall appearance of the generated images shows
high similarity to the real images, with some fine grained details not fully
captured by the CNN.

### Data availability

Data is fully available on request from the authors. Our code for training and
evaluating an RL agent based on a simulated MOT can be found under [38].

## References

  * [1] Phillips, W. D.  Nobel lecture: Laser cooling and trapping of neutral atoms.  _Rev. Mod. Phys._ 70, 721–741 (1998).
  * [2] Acín, A. _et al._ The quantum technologies roadmap: a european community view.  _New Journal of Physics_ 20, 080201 (2018).
  * [3] Schreck, F. & Druten, K. v.  Laser cooling for quantum gases.  _Nature Physics_ 17, 1296–1304 (2021).
  * [4] Kadau, H. _et al._ Observing the rosensweig instability of a quantum ferrofluid.  _Nature_ 530, 194–197 (2016).
  * [5] Norcia, M. A. _et al._ Two-dimensional supersolidity in a dipolar quantum gas.  _Nature_ 596, 357–361 (2021).
  * [6] Bothwell, T. _et al._ Resolving the gravitational redshift across a millimetre-scale atomic sample.  _Nature_ 602, 420–424 (2022).
  * [7] Petrich, W., Anderson, M. H., Ensher, J. R. & Cornell, E. A.  Behavior of atoms in a compressed magneto-optical trap.  _J. Opt. Soc. Am. B_ 11, 1332–1335 (1994).
  * [8] Lett, P. D. _et al._ Observation of atoms laser cooled below the doppler limit.  _Phys. Rev. Lett._ 61, 169–172 (1988).
  * [9] Kasevich, M. & Chu, S.  Laser cooling below a photon recoil with three-level atoms.  _Phys. Rev. Lett._ 69, 1741–1744 (1992).
  * [10] Shahriari, B., Swersky, K., Wang, Z., Adams, R. P. & de Freitas, N.  Taking the human out of the loop: A review of bayesian optimization.  _Proceedings of the IEEE_ 104, 148–175 (2016).
  * [11] Tranter, A. D. _et al._ Multiparameter optimisation of a magneto-optical trap using deep learning.  _Nature Communications_ 9, 4360 (2018).
  * [12] Wigley, P. B. _et al._ Fast machine-learning online optimization of ultra-cold-atom experiments.  _Scientific Reports_ 6, 25890 (2016).
  * [13] Barker, A. J. _et al._ Applying machine learning optimization methods to the production of a quantum gas.  _Machine Learning: Science and Technology_ 1, 015007 (2020).
  * [14] Davletov, E. T. _et al._ Machine learning for achieving bose-einstein condensation of thulium atoms.  _Physical Review A_ 102 (2020).
  * [15] Vendeiro, Z. _et al._ Machine-learning-accelerated bose-einstein condensation.  _Phys. Rev. Res._ 4, 043216 (2022).
  * [16] Sutton, R. S. & Barto, A. G.  _Reinforcement Learning: An Introduction_ (The MIT Press, 2018), second edn.
  * [17] Mnih, V. _et al._ Human-level control through deep reinforcement learning.  _Nature_ 518, 529–533 (2015).
  * [18] Silver, D. _et al._ A general reinforcement learning algorithm that masters chess, shogi, and go through self-play.  _Science_ 362, 1140–1144 (2018).
  * [19] Silver, D. _et al._ Mastering the game of go without human knowledge.  _Nature_ 550, 354–359 (2017).
  * [20] Vinyals, O. _et al._ Grandmaster level in starcraft ii using multi-agent reinforcement learning.  _Nature_ 575, 350–354 (2019).
  * [21] Bakhtin, A. _et al._ Human-level play in the game of Diplomacy by combining language models with strategic reasoning.  _Science_ 378, 1067 (2022).
  * [22] Dulac-Arnold, G. _et al._ Challenges of real-world reinforcement learning: definitions, benchmarks and analysis.  _Machine Learning_ 110, 2419–2468 (2021).
  * [23] Reuer, K. _et al._ Realizing a deep reinforcement learning agent discovering real-time feedback control strategies for a quantum system.  _arXiv:2210.16715_ (2022).
  * [24] Tünnermann, H. & Shirakawa, A.  Deep reinforcement learning for coherent beam combining applications.  _Opt. Express_ 27, 24223–24230 (2019).
  * [25] Konishi, M., Inubushi, M. & Goto, S.  Fluid mixing optimization with reinforcement learning.  _Scientific Reports_ 12, 14268 (2022).
  * [26] Bellemare, M. G. _et al._ Autonomous navigation of stratospheric balloons using reinforcement learning.  _Nature_ 588, 77–82 (2020).
  * [27] Degrave, J. _et al._ Magnetic control of tokamak plasmas through deep reinforcement learning.  _Nature_ 602, 414–419 (2022).
  * [28] Duan, Y., Chen, X., Houthooft, R., Schulman, J. & Abbeel, P.  Benchmarking deep reinforcement learning for continuous control.  In Balcan, M. F. & Weinberger, K. Q. (eds.) _Proceedings of The 33rd International Conference on Machine Learning_ , vol. 48 of _Proceedings of Machine Learning Research_ , 1329–1338 (PMLR, New York, New York, USA, 2016).
  * [29] Caspi, I., Leibovich, G., Novik, G. & Endrawis, S.  Reinforcement learning coach (2017).
  * [30] https://everycoldatom.com/.
  * [31] Ketterle, W., Durfee, D. & Stamper-Kurn, D.  Making, probing and understanding Bose-Einstein condensates.  In Inguscio, M., Stringari, S. & Wieman, C. E. (eds.) _Proceedings of the International School of Physics - Enrico Fermi_ , 67–176 (IOS Press, 1999).
  * [32] Gutterres, R. F. _et al._ Determination of the ${}^{87}\mathrm{Rb}$ $5p$ state dipole matrix element and radiative lifetime from the photoassociation spectroscopy of the ${\mathrm{rb}}_{2}$ ${0}_{g}^{-}{(P}_{3/2})$ long-range state.  _Phys. Rev. A_ 66, 024502 (2002).
  * [33] Lillicrap, T. P. _et al._ Continuous control with deep reinforcement learning.  _arXiv:1509.02971_ (2015).
  * [34] Geirhos, R. _et al._ Shortcut learning in deep neural networks.  _Nature Machine Intelligence_ 2, 665–673 (2020).
  * [35] Tobin, J. _et al._ Domain randomization for transferring deep neural networks from simulation to the real world.  In _2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 23–30 (2017).
  * [36] Note, that the agent does not differentiate between training or evaluation episodes, since its actions solely depend on the obtained observations.
  * [37] Dalibard, J. & Cohen-Tannoudji, C.  Laser cooling below the doppler limit by polarization gradients: simple theoretical models.  _J. Opt. Soc. Am. B_ 6, 2023–2045 (1989).
  * [38] https://github.com/mpi-is/rl-coach-for-mot.
  * [39] Haarnoja, T., Zhou, A., Abbeel, P. & Levine, S.  Soft actor-critic: Off-policy maximum entropy deep reinforcement learning with a stochastic actor.  In Dy, J. & Krause, A. (eds.) _Proceedings of the 35th International Conference on Machine Learning_ , vol. 80 of _Proceedings of Machine Learning Research_ , 1861–1870 (PMLR, 2018).
  * [40] Schulman, J., Wolski, F., Dhariwal, P., Radford, A. & Klimov, O.  Proximal policy optimization algorithms (2017).  eprint 1707.06347.
  * [41] Günther, A., Bender, H., Stibor, A., Fortágh, J. & Zimmermann, C.  Observing quantum gases in real time: Single-atom detection on a chip.  _Phys. Rev. A_ 80, 011604(R) (2009).
  * [42] Andrews, M. R. _et al._ Propagation of sound in a bose-einstein condensate.  _Phys. Rev. Lett._ 79, 553–556 (1997).
  * [43] Kaminski, F. _et al._ In-situ dual-port polarization contrast imaging of faraday rotation in a high optical depth ultracold 87rb atomic ensemble.  _Eur. Phys. J. D_ 66, 227 (2012).
  * [44] Fortagh, J., Ott, H., Grossmann, A. & Zimmermann, C.  A miniaturized magnetic guide for neutral atoms.  _App. Phys. B_ 70, 701 (2000).
  * [45] Bagnato, V. S. _et al._ Spatial distribution of atoms in a magneto-optical trap.  _Phys. Rev. A_ 48, 3771–3775 (1993).
  * [46] Hinton, G. E. & Salakhutdinov, R. R.  Reducing the dimensionality of data with neural networks.  _Science_ 313, 504–507 (2006).
  * [47] Masci, J., Meier, U., Cireşan, D. & Schmidhuber, J.  Stacked convolutional auto-encoders for hierarchical feature extraction.  In Honkela, T., Duch, W., Girolami, M. & Kaski, S. (eds.) _Artificial Neural Networks and Machine Learning – ICANN 2011_ , 52–59 (Springer Berlin Heidelberg, Berlin, Heidelberg, 2011).

Acknowledgments
The authors acknowledge fruitful discussions with G. Martius. We acknowledge
funding from the German Research Foundation (DFG) through the Research Unit
FOR 5413/1, Grant No. 465199066.

Author contributions
A.G., J.F., and V.V. planned and supervised the project. M.R. carried out the
experiments. V.V. carried out the simulations. M.R., A.G. and V.V. analyzed
the data. M.R., A.G. and V.V. wrote the manuscript. All authors discussed the
results and contributed to the manuscript.

Competing Interests
The authors declare no competing interests.


# Neutrino secret self-interactions: a booster shot for the cosmic neutrino
background

Anirban Das  0000-0002-7880-9454<EMAIL_ADDRESS>SLAC National
Accelerator Laboratory, Stanford University, 2575 Sand Hill Road, Menlo Park,
CA 94025, USA  Yuber F. Perez-Gonzalez  0000-0002-2020-7223 yuber.f.perez-
<EMAIL_ADDRESS>Institute for Particle Physics Phenomenology, Durham
University, South Road DH13EL, Durham, United Kingdom  Manibrata Sen
0000-0001-7948-4332<EMAIL_ADDRESS>Max-Planck-Institut für
Kernphysik, Saupfercheckweg 1, 69117 Heidelberg, Germany

###### Abstract

Neutrinos might interact among themselves through forces that have so far
remained hidden. Throughout the history of the Universe, such _secret_
interactions could lead to scatterings between the neutrinos from supernova
explosions and the non-relativistic relic neutrinos left over from the Big
Bang. Such scatterings can boost the cosmic neutrino background (C$\nu$B) to
energies of ${\cal O}$(MeV), making it, in principle, observable in
experiments searching for the diffuse supernova neutrino background. Assuming
a model-independent four-Fermi interaction, we determine the upscattered
cosmic neutrino flux, and derive constraints on such secret interactions from
the latest results from Super-Kamiokande. Furthermore, we also study prospects
for detection of the boosted flux in future lead-based coherent elastic
neutrino-nucleus scattering experiments.

††preprint: SLAC-PUB-17673, IPPP/22/24

##  I Introduction

The Standard Model of particle physics has come a long way since Pauli’s
proposal of the elusive neutrino in 1930. Not only have neutrinos been
detected, but also it has been proven conclusively that these elementary
particles have a tiny mass, thereby allowing different flavors of neutrinos to
oscillate among each other [1]. While three flavors of neutrinos have been
shown to exist, tremendous amount of theoretical and experimental efforts have
been underway to explore the existence of additional families of neutrinos.
Neutrinos are known to interact with other elementary particles as well as
themselves through weak interactions. However, what is not known conclusively
is whether there exists secret interactions among neutrinos.

The hypothesis for such secret self-interactions of neutrinos has a plethora
of consequences for the early Universe, astrophysical setups like compact
objects as well as terrestrial experiments (see Ref. [2] for a comprehensive
recent review). Laboratory bounds on such self-interactions can arise from
different sources - Higgs invisible decay bounds, constraints arising from Z
boson decay width, decays of $\pi,K,D$ mesons, as well as neutrinoless double
beta decay searches [3, 4, 5, 6, 7, 8, 9]. From a cosmological perspective, a
strong bound on these self-interactions exists from the successful predictions
of primordial abundances of elements during big-bang nucleosynthesis [10, 11],
as well as from its impact on the CMB [12, 13, 14, 15, 16, 17, 18] and
structure-formation [19]. Astrophysical probes of such strong self-
interactions can arise from stellar cooling arguments [20], core-collapse
supernovae (CCSN) [21, 22, 23, 24, 25, 26, 27], as well as the observation of
high-energy neutrinos in neutrino telescopes [28, 29, 30, 31, 32, 33, 34].
From a flavor point of view, it is important to appreciate that the most
stringent bounds exist in the electron and the muon sector, with substantially
weaker bounds for tau-neutrinos. Self-interactions of active neutrinos have
also been shown to assist in successful production of sterile neutrino dark
matter, alleviating astrophysical bounds arising from un-observation of X-rays
[6, 9, 35, 36, 37, 38].

One particular direction remains unexplored, nevertheless. The Universe is
filled with an abundance of thermalized relic neutrinos, which decoupled when
the Universe had a temperature of $\sim 1{\rm\ MeV}$, and have been free-
streaming since then. This cosmic neutrino background (C$\nu$B) consists non-
relativistic neutrinos***If the lightest neutrino has a mass $m_{0}\lesssim
0.6$ meV, it would still be relativistic today., with a temperature of around
$\mathcal{O}(0.1)$ meV. On the other hand, the Earth also receives a
considerable flux of relativistic neutrinos from different sources, for
example, an isotropic sea of MeV neutrinos from the diffuse supernova neutrino
background (DSNB), higher energy neutrinos from cosmogenic sources, and so on.
In the presence of large secret self-interactions, these relativistic
neutrinos can scatter with the C$\nu$B, and boost them to higher energies.
This can provide a healthy fraction of _relativistic_ C$\nu$B flux, with
energies comparable to those of the scattering neutrinos, and hence can be
detected in current as well as future neutrino detectors.

In this work, we demonstrate this idea with the help of the DSNB-C$\nu$B
scattering. As mentioned before, the DSNB consists of MeV neutrinos of all
flavors emanating from all CCSNe since the beginning of the star-formation
epoch [39]. After interacting, the neutrinos from the C$\nu$B get scattered to
MeV energies, and thus can be probed by experiments sensitive to the DSNB. The
boosted C$\nu$B spectra will add on to the DSNB spectra, and lead to an
increased number of events detected at lower energies. On the experimental
front, a lot of effort has been underway towards the detection of the DSNB.
Super-Kamiokande (SK) [40] in Japan, loaded with Gadolinium, leads the pack,
followed by upcoming experiments like the Jiangmen Underground Neutrino
Observatory (JUNO) in China [41], Deep Underground Neutrino Obervatory (DUNE)
in the USA [42], Hyper-Kamiokande (HK) — an upgrade of the SK detector [43],
among others. A number of other ideas like water-based liquid scintillator
detectors [44, 45] , paleo-detectors [46], as well as using archaeological
lead [47] for detection are being explored. The possibility of using coherent
elastic neutrino-nucleus scattering to detect the DSNB has also been explored
recently [48]. We utilize some of these experiments to probe the boosted
C$\nu$B spectra, and use the results to put constraints on secret neutrino
self-interactions.

The paper is organized as follows. In Section II, we introduce neutrino self-
interactions and establish the notation for the following sections. In section
III, we describe briefly the standard properties of the DSNB and the C$\nu$B.
We then describe the computation of the upscattered C$\nu$B flux in section
IV. Using the computed boosted flux, in section V, we determine the current
constraints from SK, and analyse the future prospects for detection in a lead-
based coherent elastic neutrino-nucleus scattering experiment. Finally, we
draw our conclusions in section VI. We use natural units in which
$\hbar=c=k_{\rm B}=1$ throughout this manuscript.

##  II Neutrino self-interactions

Neutrino self-interactions can arise due to exchange of a new scalar, and/or a
vector particle. Scalar mediated self-interactions can arise from the
following higher-dimensional operator involving the Higgs doublet $(H)$ and
the lepton doublet $(L)$ [4],

$\mathcal{L}\supset\frac{f^{\varphi}_{\alpha\beta}}{\Lambda^{2}}\left(\overline{L}_{\alpha}\,H\right)\left(H^{T}\,L^{C}_{\beta}\right)\varphi^{*},\,$
(II.1)

where $\varphi$ is a complex scalar with nonzero lepton-number to avoid
violating such a symmetry. After the electroweak symmetry breaking, one
arrives at the following neutrinophilic interaction,

$\mathcal{L}\supset
g^{\varphi}_{\alpha\beta}\,\overline{\nu}_{\alpha}\nu^{C}_{\beta}\varphi^{*},\,$
(II.2)

where $g^{\varphi}_{\alpha\beta}=f^{\varphi}_{\alpha\beta}v^{2}/\Lambda^{2}$.

Similarly, for vector-boson mediated self-interactions, we consider a gauge
boson $V$ which couples only to the active neutrinos at low energies through
the following operator [9],

$\mathcal{L}\supset\frac{f^{\rm
V}_{\alpha\beta}}{\Lambda^{2}}\left(\overline{L}_{\alpha}\,i\sigma_{2}H^{*}\right)\gamma^{\mu}\left(H^{T}\,i\sigma_{2}L_{\beta}\right)V_{\mu}\,,$
(II.3)

where the cut-off scale $\Lambda$ signifies the energy scale close to which
the effective theory breaks down. Once the Higgs attains a vacuum expectation
value, the above operator leads to active neutrino self-interactions of the
form,

$\mathcal{L}\supset g^{\rm
V}_{\alpha\beta}\,\overline{\nu}_{\alpha}\gamma^{\mu}\nu_{\beta}V_{\mu}\,,$
(II.4)

where $g^{\rm V}_{\alpha\beta}=f^{\rm V}_{\alpha\beta}v^{2}/\Lambda^{2}$.

In this paper, we shall mostly be concerned with scattering energy ($\lesssim
10{\rm\ MeV}$) below the mediator mass. In that case, the neutrino self-
interaction can be described by four-Fermi interaction and the scattering
cross sections can always be approximated as

$\displaystyle\sigma\approx\frac{(g_{\alpha\beta}^{\varphi,V})^{4}s}{16\pi
M_{\varphi,V}^{4}}\,\equiv G_{X}^{2}m_{\nu}E_{\nu}\,.$ (II.5)

Here $s=2m_{\nu}E_{\nu}$ is the center-of-mass energy when the neutrino from
the DSNB has energy $E_{\nu}$ and the C$\nu$B is almost at rest, and the
effective coupling
$G_{X}=\sqrt{2/\pi}(g_{\alpha\beta}^{\varphi,V})^{2}/(4M_{\varphi,V}^{2})$
being determined by the nature of the interaction. This $E_{\nu}$-scaling in
Eq.(II.5) is inspired by the exact expression of a four-Fermi interaction
cross section when one particle is at rest.

We also focus on flavor universal couplings such that
$g_{ee}=g_{\mu\mu}=g_{\tau\tau}$, while the flavor off-diagonal couplings are
set to zero. Note that in such a case, the couplings are identical in the
flavor and the mass bases†††The distinction between flavor and mass basis is
crucial since at the time of scattering, the neutrinos from the DSNB and well
as the C$\nu$B have decohered and behave as pure mass-eigenstates. In Fig. 1,
we show the contours of a few values of $G_{X}^{2}m_{\nu}$, that we shall be
using in this paper, in the plane of $g^{\varphi}-M_{\varphi}$ for the scalar
mediator. The results for the vector mediators are similar, see [9] for
further details. Other experimental constraints are also shown for reference.
Heretofore, we stay agnostic of the nature of the mediator when computing the
boosted C$\nu$B spectra. The formalism presented below is generic, and can be
used for either case. We will specify the nature of the mediator only when
translating experimental probes to constraints on the mediator mass-coupling
parameter space.

Figure 1: Isocontours of $G_{X}^{2}m_{\nu}$ for the scalar case in the plane
$g_{\alpha\beta}^{\varphi}$ vs $M_{\varphi}$. We present constraints from
Invisible $Z$ decay (dark orange region), meson decays
$\mathfrak{m}^{\pm}\rightarrow\ell\nu\varphi$ (purple), and double-beta decay
(green), taken from [2].

##  III DSNB and C$\nu$B

Each CCSN emits a near-thermal flux of neutrinos that can be approximated by a
Fermi-Dirac distribution with a flavor-dependent temperature $T_{\nu}$ [39],

$F_{\nu}(E)=\frac{E_{\nu}^{\mathrm{tot}}}{6}\,\frac{120}{7\pi^{4}}\,\frac{E_{\nu}^{2}}{T_{\nu}^{4}}\,\frac{1}{e^{E_{\nu}/T_{\nu}}+1}\,,$
(III.6)

Table 1: The star formation rate parameters used in Eq.(III.7). Parameter |  Value
---|---
$\rho_{*}$ |  $0.02\,M_{\odot}\mathrm{~{}yr^{-1}Mpc^{-3}}$
$a$ |  $3.4$
$b$ |  $-0.3$
$c$ |  $-2.5$
$z_{1}$ |  $1$
$z_{2}$ |  $4$
$\zeta$ |  $-10$

where $E_{\nu}^{\mathrm{tot}}$ is the total energy in all neutrinos. We assume
$T_{\nu_{e}}=6.6{\rm\ MeV},T_{\bar{\nu}_{e}}=7{\rm\ MeV},$ and
$T_{\nu_{x}}=10{\rm\ MeV}$ unless otherwise stated. We will explore the effect
of having different temperatures when determining current constrains from SK.
The net DSNB flux received at earth can be obtained by integrating the fluxes
from each individual SN that are distributed over a range of redshifts. The
redshift distribution of the SNe is proportional to the star formation rate
which can be fitted with a broken power law [49],

$\begin{split}R_{\mathrm{CCSN}}(z)&=\frac{\rho_{*}}{143\,M_{\odot}}\\\
&\times\left((1+z)^{a\zeta}+\left(\frac{1+z}{B}\right)^{b\zeta}+\left(\frac{1+z}{C}\right)^{c\zeta}\right)^{1/\zeta}\,,\end{split}$
(III.7)

where $B=(1+z_{1})^{1-a/b},C=(1+z_{1})^{(b-a)/c}(1-z_{2})^{1-b/c}$. The
numerical values of the various parameters are given in Table 1. The DSNB flux
at a redshift $z$ can be obtained by convolving the CCSN rate with the
neutrino flux from an individual SN including the expansion of the Universe,

$\Phi(E_{\nu},z)=(1+z)^{2}\int_{z}^{z_{\mathrm{max}}}\frac{dz^{\prime}}{H(z^{\prime})}F_{\nu}[E_{\nu}(1+z^{\prime})]R_{\mathrm{CCSN}}(z^{\prime})\,.$
(III.8)

We multiply $E_{\nu}$ by an appropriate factor of $(1+z^{\prime})$ as the the
neutrinos were emitted at that redshift. The Hubble parameter is given by
$H(z)=H_{0}\sqrt{\Omega_{\Lambda}+\Omega_{M}(1+z)^{3}}$ with
$H_{0}=67\,\mathrm{km~{}s^{-1}Mpc^{-1}}$ as per the Planck 2018 data [50]. We
integrate the SN formation rate up to redshift ${z_{\mathrm{max}}}=6$. As an
illustration, in Fig. 2, we show the DSNB spectrum for $\overline{\nu}_{e}$ at
$z=0$.

Figure 2: $\overline{\nu}_{e}$ flux as function of neutrino energy for the
standard DSNB (dark orange band) and the upscattered C$\nu$Bflux for four
different values of $G^{2}_{X}m_{\nu}=10^{-9}~{}{\rm MeV^{-3}}$ (purple),
$10^{-10}~{}{\rm MeV^{-3}}$ (emerald), $10^{-11}~{}{\rm MeV^{-3}}$ (dark
yellow), $10^{-12}~{}{\rm MeV^{-3}}$ (lilac). We present current constrains on
$\overline{\nu}_{e}$ searches from Borexino [51] (blue) and SK-IV [40] (red).

The C$\nu$B is the population of relic neutrinos that were frozen out after
decoupling from the thermal plasma. The weak interactions of the SM neutrinos
became slower than the Hubble expansion rate around a photon temperature
$T_{\gamma}\simeq 1{\rm\ MeV}$. Since then, neutrinos have been free-streaming
and cooling down as a result of the expansion of the Universe. They are
believed to have a Fermi-Dirac distribution today with a temperature
$T_{\nu}\simeq 1.95\,\mathrm{K}$. Hence the number density of C$\nu$B is given
by

$n_{\nu}=\frac{3}{4}\frac{\zeta(3)}{\pi^{2}}gT_{\nu}^{3}\,,$ (III.9)

where $g$ is the number of degrees of freedom. Looking at their temperature
today, it is clear that these neutrinos have very low energy ($T_{\nu}\simeq
1.95\,\mathrm{K}\approx 0.17\,\mathrm{meV}$) compared to the DSNB. Their
temperature $T_{\nu}$ today is slightly lower than the photons as the
electrons and positrons annihilated after neutrino decoupling and deposited
their energy into the photon bath heating it up approximately by a factor of
$(11/4)^{1/3}$. Although we do not know the absolute masses of the SM
neutrinos, the atmospheric mass difference $\Delta
m_{\mathrm{atm}}^{2}=2.5\times 10^{-3}\,\mathrm{eV^{2}}$ tells us that at
least one of the neutrino mass states is nonrelativistic today [52].

##  IV Boosted flux

As discussed before, in the presence of neutrino self-interactions, the
supernova relic neutrinos scatter off the C$\nu$B, and transfer their energy
to the latter. This allows a fraction of the non-relativistic C$\nu$B to get
boosted to higher energies, comparable to those of the DSNB. This can have
interesting consequences for experiments searching for the DSNB.

Figure 3: Total $\overline{\nu}_{e}$ flux for $E_{\nu}>17.8$ MeV as function
of $G^{2}_{X}m_{\nu}$ for two different sets of SN temperatures,
$T_{\bar{\nu}_{e}}=7{\rm\ MeV},$ and $T_{\nu_{x}}=10{\rm\ MeV}$ (green) and
$T_{\bar{\nu}_{e}}=6{\rm\ MeV},$ and $T_{\nu_{x}}=7{\rm\ MeV}$ (light blue).
The full lines correspond to the total DSNB+boosted C$\nu$B fluxes, while the
dashed curves depict the contribution from the upscattered cosmic neutrinos.
The dark lilac band corresponds to the SK-IV limit at the 90% C.L.

Neutrinos originating from a SN at a redshift
$z^{\prime}\in(0-{z_{\mathrm{max}}})$ will scatter off the C$\nu$B at some
later redshift $z^{\prime\prime}\leq z^{\prime}$. The scattered DSNB flux as a
function of energy $E_{D}$ at redshift $z$ can then be estimated as

$\frac{d\Phi_{D}(z)}{dE_{D}}=(1+z)^{2}\int_{z}^{z_{\mathrm{max}}}\frac{dz^{\prime}}{H(z^{\prime})}R_{\mathrm{CCSN}}(z^{\prime})F_{\nu}[E_{D}(1+z^{\prime})]\times\exp\left[{-\int_{z}^{z^{\prime}}\frac{dz^{\prime\prime}}{(1+z^{\prime\prime})H(z^{\prime\prime})}n_{\nu}(z^{\prime\prime})\sigma}\right]\,,$
(IV.10)

where $n_{\nu}(z^{\prime\prime})$ is the number density of the C$\nu$B at the
redshift $z^{\prime\prime}$. At the same time, the non-relativisitic neutrinos
from the C$\nu$B gets upscattered to an energy $E_{C}$. Because the neutrinos
from the DSNB carry $\mathcal{O}(10{\rm\ MeV})$ energy, the maximum energy
transferred to C$\nu$B, ignoring the mass of neutrino, is simply
$E_{C}^{\mathrm{max}}\approx E_{D}$. The upscattered C$\nu$B flux can be
computed as

$\frac{d\psi_{C}}{dE_{C}}=\int_{E_{C}}^{\infty}dE_{D}\int_{0}^{z_{\mathrm{max}}}\frac{dz}{H(z)}n_{\nu}(z)\,\frac{d\sigma}{dE_{C}}\,\frac{d\Phi_{D}(z)}{dE_{D}},$
(IV.11)

In the heavy mediator limit, i.e., when the mediator mass is heavier than the
energy transferred, the differential cross section can be simplified to
$d\sigma/dE_{C}=\sigma/E_{C}^{\mathrm{max}}$. Eq. IV.11 will be used to
estimate the upscattered C$\nu$B flux. The details of modelling the self-
interactions enter through the cross-section $\sigma$.

Fig. 2 depicts the upscattered $\bar{\nu}_{e}$ component of the C$\nu$B for
different values of $G^{2}_{X}m_{\nu}$. The corresponding $\bar{\nu}_{e}$ from
the DSNB is also shown for comparison. We find that for
$G^{2}_{X}m_{\nu}>10^{-9}\,{\rm MeV}^{-3}$, the boosted cosmic neutrino flux
dominates over the SN relic neutrino flux for energies $E_{\nu}>1{\rm\ MeV}$.
Since $E_{D}\gg E_{C}$ and the scattering particles have the same mass, almost
all of the SN relic neutrino energy gets transferred to the cosmic neutrino
flux; hence the upscattered flux is independent of energy for low energies,
and follows the DSNB spectra for $E_{\nu}>$ few MeVs. The corresponding limits
from searches of $\bar{\nu}_{e}$ in these energy ranges from SK [40] and
Borexino [51] are also shown. While constraints from Borexino are weaker, and
can be probed with larger values of $G^{2}_{X}m_{\nu}$, the current bounds
from SK can already probe such a boosted $\bar{\nu}_{e}$ flux.

##  V Detection and limits

Fig. 3 shows the current limit of 90% C.L. from SK-IV searches of the DSNB as
function of $G^{2}_{X}m_{\nu}$. SK will be sensitive to the net neutrino flux
(solid lines) consisting of the original DSNB flux, as well as the boosted
C$\nu$B flux (dashed lines). The green and the blue lines correspond two
different choices of $T_{\bar{\nu}_{e}}=7{\rm\
MeV},\,T_{\nu_{\mu,\tau}}=10{\rm\ MeV}$, and $T_{\bar{\nu}_{e}}=6{\rm\
MeV},\,T_{\nu_{\mu,\tau}}=7{\rm\ MeV}$. We find that the former case,
characterized by slightly larger temperatures than what is expected from SN
simulations, is already ruled out by SK-IV for $G^{2}_{X}m_{\nu}\gtrsim
10^{-11}{\rm MeV}^{-3}$. This can be directly translated to bounds on the
mediator mass and couplings for a specific model of neutrino self-
interactions. The case corresponding to lower values of
$(T_{\bar{\nu}_{e}},T_{\nu_{\mu,\tau}})$ stays out of reach of the current SK
bounds, but can be probed by future surveys. Indeed, SK intends to constrain
$\Phi_{\overline{\nu}_{e}}\lesssim 1~{}{\rm cm^{-2}s^{-1}}$ for
$E_{\nu}\gtrsim 17.3$ MeV in the next ten-years, after doping their detector
with Gadolinium [40]. Furthermore, we find that there is little sensitivity to
values of $G^{2}_{X}m_{\nu}<10^{-11}{\rm MeV}^{-3}$. In this case, the boosted
component of the C$\nu$B flux is too tiny to make any appreciable difference
to the DSNB flux.

Figure 4: Expected CEvNS differential recoil rate for a Pb target as function
of the nuclear recoil energy. We present the recoil rates produced by the
upscattered C$\nu$B flux for four different values of
$G^{2}_{X}m_{\nu}=10^{-9}~{}{\rm MeV^{-3}}$ (purple), $10^{-10}~{}{\rm
MeV^{-3}}$ (emerald), $10^{-11}~{}{\rm MeV^{-3}}$ (dark yellow),
$10^{-12}~{}{\rm MeV^{-3}}$ (lilac). The black line represents the rate for
solar + atmospheric neutrinos, while the dark orange band corresponds to the
prediction for the DSNB. Figure 5: Expected CEvNS number of events for a Pb
target as function of the experimental threshold energy, produced by the
upscattered C$\nu$B flux for four different values of
$G^{2}_{X}m_{\nu}=10^{-9}~{}{\rm MeV^{-3}}$ (purple), $10^{-10}~{}{\rm
MeV^{-3}}$ (emerald), $10^{-11}~{}{\rm MeV^{-3}}$ (dark yellow),
$10^{-12}~{}{\rm MeV^{-3}}$ (lilac). The black line represents the events for
solar + atmospheric neutrinos, while the dark orange band corresponds to the
prediction for the DSNB. Figure 6: Total neutrino flux times energy as
function of neutrino momentum for the upscattered C$\nu$B flux for four
different values of $G^{2}_{X}m_{\nu}=10^{-9}~{}{\rm MeV^{-3}}$ (purple),
$10^{-10}~{}{\rm MeV^{-3}}$ (emerald), $10^{-11}~{}{\rm MeV^{-3}}$ (dark
yellow), $10^{-12}~{}{\rm MeV^{-3}}$ (lilac). For comparison, we present the
standard C$\nu$B (dashed lilac), relic neutrinos from BBN (orange dashed and
dot-dashed), solar neutrinos (blue dashed), and low-energy atmospheric
neutrinos (red dashed).

Recent works have demonstrated the feasibility of observing the DSNB via the
Coherent Elastic Neutrino-Nucleus Scattering (CEvNS) in future experiments
[48, 46]. Specifically, the RES-NOVA proposal will exploit the enhancement of
the CEvNS cross section for a heavy nucleus like lead, together with the high
purity of archaeological Pb [47], in such a way that it will be possible to
search for the DSNB. Therefore, we can expect that such facilities could also
search for the boosted C$\nu$B flux from our scenario. We estimate the recoil
rates and number of events in a generic experiment using Lead as follows. The
recoil rate is computed via

$\displaystyle\frac{dR}{dE_{R}}={\cal E}\int dE_{C}\frac{d\sigma_{\rm
CEvNS}}{dE_{R}}\frac{d\psi_{C}}{dE_{C}}\,,$ (V.12)

where ${\cal E}$ is the total experimental exposure (observation time $\times$
number of targets), $E_{R}$ is the nuclear recoil energy, and $d\sigma_{\rm
CEvNS}/dE_{R}$ is the CEvNS cross section, which for completeness we quote
next,

$\displaystyle\frac{d\sigma_{\rm
CEvNS}}{dE_{R}}=\frac{G_{F}^{2}m_{X}}{4\pi}{\cal
Q}_{W}^{2}\left(1-\frac{m_{X}E_{R}}{2E_{C}^{2}}\right){\cal F}(Q^{2}),$ (V.13)

being $G_{F}$ the Fermi constant, $m_{X}$ the nuclear target mass, ${\cal
Q}_{W}=N-Z(1-4\sin^{2}\theta_{W})$ the weak vector nuclear charge, with $N,Z$
the number of neutrons, nucleons, respectively, and ${\cal F}(Q^{2})$ the
nuclear form factor, taken here to be the Helm form factor [53]. The number of
CEvNS events from the upscattered C$\nu$B flux is therefore obtained by
integrating the recoil rate in the recoil energy,

$\displaystyle{\cal N}_{C}=\int_{E_{\rm th}}dE_{R}\frac{dR}{dE_{R}}\,,$ (V.14)

where the integration is performed from the experimental energy threshold
$E_{\rm th}$.

We present in Fig. 4 the recoil rate and number of events for the boosted
C$\nu$B flux for four different values of $G^{2}_{X}m_{\nu}$ as function of
the recoil and threshold energies, respectively. We also show the total
solar+atmospheric neutrinos recoil rate (black) and the DSNB contribution
(purple band). For the recoil rate, we observe that the upscattered C$\nu$B
contribution is above the standard CEvNS from the DSNB at very low recoil
energies ($E_{R}\lesssim 1$ keV) for $G^{2}_{X}m_{\nu}\gtrsim 10^{-10}~{}{\rm
MeV^{-3}}$. Nevertheless, the enormous contribution coming from solar
neutrinos makes the search in such ranges of recoil energies rather difficult.
For larger recoil energies, we find that the recoil events from the DSNB
become larger than the boosted contribution. However, in the region
$E_{R}\in[2-3]$ keV the DSNB rate is of the same order as the atmospheric
neutrino rate. Meanwhile, for $G^{2}_{X}m_{\nu}\gtrsim 10^{-9}~{}{\rm
MeV^{-3}}$, the C$\nu$Brate is only a factor of $\sim 70\%$ the DSNB value.
Thus, in principle, detailed searches in this energy range could aid in
exploring the parameter space of secret self-interactions.

Such a behavior is reflected in the number of events dependent on the energy
threshold, Fig. 5. For relatively large energy thresholds, above the solar
neutrinos, we observe that the DSNB contribution is about half the atmospheric
one. Indeed, for a threshold of $E_{\rm th}\sim 3$ keV, larger than the
expected 1 keV threshold in the RES-NOVA proposal [47], the events for the
fiducial value of the DSNB would be $\sim 0.015$ per ton$\cdot$y while the
atmospheric background corresponds to $\sim 0.03$ per ton$\cdot$y. Meanwhile,
the events from the boosted C$\nu$B neutrinos would be $\sim 0.007$ per
ton$\cdot$y. We can estimate that at least an exposure of ${\cal E}\sim 1$
kton$\cdot$y would be necessary to have $S/\sqrt{B}\sim 1$ for
$G^{2}_{X}m_{\nu}\gtrsim 10^{-9}~{}{\rm MeV^{-3}}$, including the standard
DSNB contribution as background. For smaller $G^{2}_{X}m_{\nu}$, we would
require much larger exposures, ${\cal E}\sim 20$ kton$\cdot$y, for
$G^{2}_{X}m_{\nu}\gtrsim 10^{-11}~{}{\rm MeV^{-3}}$. Surely, these are naíve
estimates; for much more realistic estimates, we would require introducing
additional backgrounds and other experimental characteristics of RES-NOVA. We
leave such an analysis for future work.

##  VI Conclusions

Currently, only neutrinos hold the guaranteed key to the gateway of physics
beyond the Standard Model. Search for secret neutrino self-interaction is a
promising avenue that might lead us to discovery of such BSM physics. Current
laboratory bounds on these interactions allow for the possibility of large
couplings, which can even be stronger than the weak interaction by several
orders of magnitude. This presents a perfect opportunity to study the effect
that strongly self-interacting neutrinos can have in different areas.

In this paper, we study a simple yet important outcome of strongly self-
interacting neutrinos. The Universe is filled with a thermal bath of non-
relativistic neutrinos (C$\nu$B), which have decoupled when the Universe was
quite hot, at a temperature of around $1$ MeV. Another ubiquitous source of
MeV-energy (relativistic) neutrinos is all the core-collapse supernovae that
have happened in the past. In the presence of large non-standard self-
interactions, the DSNB will scatter with the C$\nu$B, and upscatter the latter
to relativistic energies. We compute this upscattered flux. The results, shown
in a larger picture in Fig. 6 for four different cases, demonstrate how the
upscattered C$\nu$B flux compare to the original one, as well as the fluxes
from other sources. The upscattered flux is comparable to, and in some cases
even larger than, the DSNB flux, and could be probed by current and upcoming
experiments. Using the most recent search of the DSNB by the Super-Kamiokande
collaboration $(22.5\times 5823)$ kton$\cdot$days of data, we put model-
independent constraints on the combination of the coupling, and the mass of
the mediator which mediates the self-interactions. We emphasize that this
constitutes a direct test of the impact of the secret neutrino self-
interactions.

Furthermore, we also compute the expected event rate for the net boosted
C$\nu$B using coherent scattering on a generic lead-based detector. We find
that for a reasonable choice of parameters, the upscattered C$\nu$B flux can
be comparable to the DSNB flux at lower recoil energies, while for larger
energies the DSNB flux dominates. A naïve estimate indicates that an exposure
of $\sim 400$ ton$\cdot$y, such as the expected in RES-NOVA-3, would be needed
to start exploring the parameter space in our scenario. For smaller values of
the mediator mass (which would be in tension with BBN in standard
cosmologies), one can get stronger self-interaction cross-sections. In such a
scenario, the boosted C$\nu$B flux can dominate over the DSNB, and hence
already by probed by experiments. However, theoretically it is difficult to
motivate such large cross-sections in standard $\Lambda$CDM cosmology, and one
needs to invoke further exotic scenarios.

A couple of additional comments are in order. One might imagine that the
boosted C$\nu$B component may show up in a detector like PTOLEMY. However, the
event rate in PTOLEMY from such a boosted component is tiny, and hence not
observable. Secondly, the upscattering of the C$\nu$B can also take place from
high energy cosmogenic neutrinos, coming from far-away sources such as active
galactic nuclei. The resultant boosted C$\nu$B component can be an important
source of flux at higher energies, and have consequences for neutrino
telescopes. One can also imagine similar boosting mechanisms to take place
from solar and atmospheric neutrinos, however these sources are rather nearby,
and hence boosting may not be efficient. Nonetheless, we believe that
searching for boosted relic neutrinos with present and next-generation
neutrino telescopes will help shed valuable light on the existence of secret
self-interactions of neutrinos.

## Acknowledgements

AD was supported by the U.S. Department of Energy under contract number DE-
AC02-76SF00515. This work has made use of the Hamilton HPC Service of Durham
University.

## References

  * Tanabashi  _et al._ [2018] M. Tanabashi _et al._ (Particle Data Group), Phys. Rev. D 98, 030001 (2018).
  * Berryman  _et al._ [2022] J. M. Berryman _et al._ , in  _2022 Snowmass Summer Study_ (2022) arXiv:2203.01955 [hep-ph] .
  * Brune and Päs [2019] T. Brune and H. Päs, Phys. Rev. D 99, 096005 (2019), arXiv:1808.08158 [hep-ph] .
  * Berryman  _et al._ [2018] J. M. Berryman, A. De Gouvêa, K. J. Kelly, and Y. Zhang, Phys. Rev. D 97, 075030 (2018), arXiv:1802.00009 [hep-ph] .
  * Blinov  _et al._ [2019] N. Blinov, K. J. Kelly, G. Z. Krnjaic, and S. D. McDermott, Phys. Rev. Lett. 123, 191102 (2019), arXiv:1905.02727 [astro-ph.CO] .
  * De Gouvêa  _et al._ [2020] A. De Gouvêa, M. Sen, W. Tangarife, and Y. Zhang, Phys. Rev. Lett. 124, 081802 (2020), arXiv:1910.04901 [hep-ph] .
  * Deppisch  _et al._ [2020] F. F. Deppisch, L. Graf, W. Rodejohann, and X.-J. Xu, Phys. Rev. D 102, 051701 (2020), arXiv:2004.11919 [hep-ph] .
  * Brdar  _et al._ [2020] V. Brdar, M. Lindner, S. Vogl, and X.-J. Xu, Phys. Rev. D 101, 115001 (2020), arXiv:2003.05339 [hep-ph] .
  * Kelly  _et al._ [2020] K. J. Kelly, M. Sen, W. Tangarife, and Y. Zhang, Phys. Rev. D 101, 115031 (2020), arXiv:2005.03681 [hep-ph] .
  * Huang  _et al._ [2018] G.-y. Huang, T. Ohlsson, and S. Zhou, Phys. Rev. D 97, 075009 (2018), arXiv:1712.04792 [hep-ph] .
  * Grohs  _et al._ [2020] E. Grohs, G. M. Fuller, and M. Sen, JCAP 07, 001, arXiv:2002.08557 [astro-ph.CO] .
  * Lancaster  _et al._ [2017] L. Lancaster, F.-Y. Cyr-Racine, L. Knox, and Z. Pan, JCAP 07, 033, arXiv:1704.06657 [astro-ph.CO] .
  * Forastieri  _et al._ [2019] F. Forastieri, M. Lattanzi, and P. Natoli, Phys. Rev. D 100, 103526 (2019), arXiv:1904.07810 [astro-ph.CO] .
  * Das and Ghosh [2021] A. Das and S. Ghosh, JCAP 07, 038, arXiv:2011.12315 [astro-ph.CO] .
  * Brinckmann  _et al._ [2021] T. Brinckmann, J. H. Chang, and M. LoVerde, Phys. Rev. D 104, 063523 (2021), arXiv:2012.11830 [astro-ph.CO] .
  * Roy Choudhury  _et al._ [2021] S. Roy Choudhury, S. Hannestad, and T. Tram, JCAP 03, 084, arXiv:2012.07519 [astro-ph.CO] .
  * Das and Ghosh [2022] A. Das and S. Ghosh, PoS EPS-HEP2021, 124 (2022), arXiv:2109.03263 [hep-ph] .
  * Venzor  _et al._ [2022] J. Venzor, G. Garcia-Arroyo, A. Pérez-Lorenzana, and J. De-Santiago, (2022), arXiv:2202.09310 [astro-ph.CO] .
  * Hannestad [2005] S. Hannestad, JCAP 02, 011, arXiv:astro-ph/0411475 .
  * Escudero  _et al._ [2019] M. Escudero, D. Hooper, G. Krnjaic, and M. Pierre, JHEP 03, 071, arXiv:1901.02010 [hep-ph] .
  * Kolb and Turner [1987] E. W. Kolb and M. S. Turner, Phys. Rev. D 36, 2895 (1987).
  * Kachelriess  _et al._ [2000] M. Kachelriess, R. Tomas, and J. W. F. Valle, Phys. Rev. D 62, 023004 (2000), arXiv:hep-ph/0001039 .
  * Farzan [2003] Y. Farzan, Phys. Rev. D 67, 073015 (2003), arXiv:hep-ph/0211375 .
  * Farzan and Palomares-Ruiz [2014] Y. Farzan and S. Palomares-Ruiz, JCAP 06, 014, arXiv:1401.7019 [hep-ph] .
  * Das  _et al._ [2017] A. Das, A. Dighe, and M. Sen, JCAP 05, 051, arXiv:1705.00468 [hep-ph] .
  * Dighe and Sen [2018] A. Dighe and M. Sen, Phys. Rev. D 97, 043011 (2018), arXiv:1709.06858 [hep-ph] .
  * Shalgar  _et al._ [2021] S. Shalgar, I. Tamborra, and M. Bustamante, Phys. Rev. D 103, 123008 (2021), arXiv:1912.09115 [astro-ph.HE] .
  * Hooper [2007] D. Hooper, Phys. Rev. D 75, 123001 (2007), arXiv:hep-ph/0701194 .
  * Ioka and Murase [2014] K. Ioka and K. Murase, PTEP 2014, 061E01 (2014), arXiv:1404.2279 [astro-ph.HE] .
  * Kelly and Machado [2018] K. J. Kelly and P. A. N. Machado, JCAP 10, 048, arXiv:1808.02889 [hep-ph] .
  * Murase and Shoemaker [2019] K. Murase and I. M. Shoemaker, Phys. Rev. Lett. 123, 241102 (2019), arXiv:1903.08607 [hep-ph] .
  * Barenboim  _et al._ [2019] G. Barenboim, P. B. Denton, and I. M. Oldengott, Phys. Rev. D 99, 083515 (2019), arXiv:1903.02036 [astro-ph.CO] .
  * Bustamante  _et al._ [2020] M. Bustamante, C. Rosenstrøm, S. Shalgar, and I. Tamborra, Phys. Rev. D 101, 123024 (2020), arXiv:2001.04994 [astro-ph.HE] .
  * Esteban  _et al._ [2021] I. Esteban, S. Pandey, V. Brdar, and J. F. Beacom, Phys. Rev. D 104, 123014 (2021), arXiv:2107.13568 [hep-ph] .
  * Kelly  _et al._ [2021] K. J. Kelly, M. Sen, and Y. Zhang, Phys. Rev. Lett. 127, 041101 (2021), arXiv:2011.02487 [hep-ph] .
  * Alonso-Álvarez and Cline [2021] G. Alonso-Álvarez and J. M. Cline, JCAP 10, 041, arXiv:2107.07524 [hep-ph] .
  * Sen [2021] M. Sen, J. Phys. Conf. Ser. 2156, 012018 (2021), arXiv:2110.09907 [hep-ph] .
  * Benso  _et al._ [2022] C. Benso, W. Rodejohann, M. Sen, and A. U. Ramachandran, Phys. Rev. D 105, 055016 (2022), arXiv:2112.00758 [hep-ph] .
  * Beacom [2010] J. F. Beacom, Ann. Rev. Nucl. Part. Sci. 60, 439 (2010), arXiv:1004.3311 [astro-ph.HE] .
  * Abe  _et al._ [2021] K. Abe _et al._ (Super-Kamiokande), Phys. Rev. D 104, 122002 (2021), arXiv:2109.11174 [astro-ph.HE] .
  * An  _et al._ [2016] F. An _et al._ (JUNO), J. Phys. G 43, 030401 (2016), arXiv:1507.05613 [physics.ins-det] .
  * Abi  _et al._ [2020] B. Abi _et al._ (DUNE), (2020), arXiv:2002.03005 [hep-ex] .
  * Abe  _et al._ [2018] K. Abe _et al._ (Hyper-Kamiokande), (2018), arXiv:1805.04163 [physics.ins-det] .
  * Askins  _et al._ [2020] M. Askins _et al._ (Theia), Eur. Phys. J. C 80, 416 (2020), arXiv:1911.03501 [physics.ins-det] .
  * Orebi Gann [2015] G. D. Orebi Gann (THEIA Interest Group), (2015), arXiv:1504.08284 [physics.ins-det] .
  * Baum  _et al._ [2022] S. Baum, F. Capozzi, and S. Horiuchi, (2022), arXiv:2203.12696 [hep-ph] .
  * Pattavina  _et al._ [2020] L. Pattavina, N. Ferreiro Iachellini, and I. Tamborra, Phys. Rev. D 102, 063001 (2020), arXiv:2004.06936 [astro-ph.HE] .
  * Suliga  _et al._ [2022] A. M. Suliga, J. F. Beacom, and I. Tamborra, Phys. Rev. D 105, 043008 (2022), arXiv:2112.09168 [astro-ph.HE] .
  * Kistler  _et al._ [2013] M. D. Kistler, H. Yuksel, and A. M. Hopkins, (2013), arXiv:1305.1630 [astro-ph.CO] .
  * Aghanim  _et al._ [2018] N. Aghanim _et al._ (Planck), (2018), arXiv:1807.06209 [astro-ph.CO] .
  * Agostini  _et al._ [2021] M. Agostini _et al._ (Borexino), Astropart. Phys. 125, 102509 (2021), arXiv:1909.02422 [hep-ex] .
  * Esteban  _et al._ [2020] I. Esteban, M. C. Gonzalez-Garcia, M. Maltoni, T. Schwetz, and A. Zhou, JHEP 09, 178, arXiv:2007.14792 [hep-ph] .
  * Helm [1956] R. H. Helm, Phys. Rev. 104, 1466 (1956).