# Uniform Distribution Technique for Neutrino Beam Scan Simulation

D. A. Wickremasinghe  S. Ganguly  K. Yonehara  R. Zwaska  Fermi National Lab,
Batavia, IL, USA  Y. Yu  P. Snopok  Illinois Institute of Technology, Chicago,
IL, USA

###### Abstract

In Fermilab’s neutrino facilities such as the Neutrinos at the Main Injector
(NuMI) and the upcoming Long Baseline Neutrino Facility (LBNF), a proton beam
strikes high-power target, producing positively and negatively charged pions
and kaons. There is a need for detailed simulations in order to capture all
particle interactions and beam propagation from protons on target to short-
lived mesons decaying into muons and neutrinos. The generation of individual
beam simulations is a resource-intensive and time-consuming process. In this
paper, we describe a method through which many simulation samples with high
statistics can be generated to study the effects of beam scan across a target
for given beam configurations.

##  I Introduction

The neutrino mass hierarchy and CP asymmetry in the lepton sector can be
tested using long-baseline neutrino oscillation experiments that provide
precision measurements [1, 2, 3, 4]. Neutrino flux at the Near Detector
strongly correlates with proton beam misalignment. The alignment of the
primary proton beam, target, and focusing horn greatly impacts the neutrino
energy spectrum.

Neutrinos at the Main Injector (NuMI) [5] neutrino beam at Fermilab is
generated by aiming high energy protons into a carbon target. The NuMI
beamline as depicted in Fig. 1 produces an intense muon neutrino beam for NuMI
neutrino experiments. A 120-GeV proton beam from the main injector collides
with a fixed graphite target to provide neutrinos. Using two focusing horns
operating at 200 kA horn current, charged particles produced from proton
interactions with target nuclei are focused into a 675-m long, 2-m diameter
cylindrical decay pipe. As the mesons decay, they can produce neutrinos and
muons before being absorbed by the hadron absorber. High-energy muons produced
by meson decay could pass through muon monitors downstream of the hadron
absorber.

Figure 1: Beamline layout for NuMI.

A world-leading neutrino beam will be produced by the upcoming Long Baseline
Neutrino Facility (LBNF) [6] as shown in Fig. 2. The baseline beamline design
consists of a 1.2-MW, 120-GeV primary proton beam impinging on a cylindrical
graphite target that measures 1.8 m in length and 1.6 cm in diameter. Using
300 kA currents, three magnetic horns focus hadrons produced in the target,
which is supported inside horn 1. In addition to the target chase, there is a
194-m long helium-filled decay pipe and a hadron absorber.

In both NuMI and LBNF, neutrinos and anti-neutrinos can be created by
operating the focusing horns in forward or reverse current configurations.

Figure 2: An illustration of the upstream portion of the LBNF neutrino
beamline. A horn-protection baffle, three focusing horns, and the decay pipe
are shown inside the target chase respectively from left to right (the beam
direction).

We can estimate the effects of changes in the beam parameters by scanning the
primary proton beam across the target. Typical beam scan processes [7] involve
systematically changing proton beam position on target, beam spot size and
horn currents. In NuMI, these changes are visible in the response of each of
the three muon monitors (MM1–MM3) that are identical to each other. Due to
intervening shielding between two consecutive detectors, each detector sees a
different range of muon momenta. Note that these initial ranges are subsets
rather than distinct ranges.

Each muon monitor in the NuMI facility consists of a $9\times 9$ array of
ionization chambers as shown in Fig. 3. There are two parallel-plate
electrodes separated by a 0.3-cm gap in each ionization chamber filled with
helium gas. LBNF will monitor its beam quality through several detectors,
including the Muon Monitor System (MuMS). A NuMI-based approach will be used
in the conceptual design of MuMS.

Figure 3: Schematics of the NuMI muon monitors.

For the purpose of establishing a correlation between data and simulations,
beam scan data must be compared with simulation results. Combining data and
simulation also provides valuable beam diagnostics. High-statistics
simulations with “realistic” Gaussian distributions for various beam and horn
current settings generally require substantial computing resources and long
simulation times. A Gaussian distribution characterizes the intensity profile
of a Gaussian beam. In a Gaussian beam, the intensity is highest at the center
and gradually decreases while moving away from the center in a transverse
direction (perpendicular to the beam axis). Simulation time and computing
resources can be significantly reduced by simulating a single high-statistics
sample with a uniform distribution for each beam parameter where the uniform
beam’s intensity remains approximately the same from the center to the edges.
After that, a Gaussian weight can be calculated and applied to the sample in
post-processing. To understand neutrino beam variations from muon monitor data
and simulations, Monte Carlo samples with high statistics are essential. Using
a single high-statistics simulation sample with a uniform distribution, many
Gaussian distributions can be created that are similar to those found in real
beam scan studies. This proposed technique can generate a significant number
of Monte Carlo (MC) data samples by varying the incident beam parameters and
horn current settings. Section II provides a brief description of the
simulation tools used in the studies. Section III presents the details of the
simulation technique we have proposed. Section IV contains information about
the simulation data that have been generated. In Sec. V we compare uniform
beam simulation results with the nominal beam simulation results and validate
the uniform beam simulation technique. Section VI gives details on the
computer resources needed for these simulations. Section VII illustrates the
application of the technique with a few beam scan examples with NuMI
simulations. In Sec. VIII examples of uniform beam simulation technique are
illustrated for the LBNF simulation. Sec. IX highlights the usefulness of the
technique for machine learning applications. Sec. X briefly summarizes the
details of the current state and the future applications of the proposed
technique.

##  II Simulation Tools

For NuMI and LBNF neutrino studies with simulation, g4numi and g4lbnf use end-
to-end Monte Carlo simulation based on Geant4 [8, 9]. The simulations account
for beamline geometry, including the targets, focusing horns, and decay
volumes. Both secondary meson production and subsequent interactions with
beamline components are considered as part of the full sequence of particles
produced in the first collisions between protons and targets. When a neutrino
is produced in a particle decay, information about its momentum, position, and
ancestor is stored in the “dk2nu” format. A dk2nu ntuple is essentially a list
of neutrinos generated by a beam simulation. It is a flux ntuple format and
library with methods for analyzing flux ntuples (e.g., location weights).
DUNE, NOvA, and MINERvA are some of the Fermilab experiments that are using
this ntuple format. Neutrinos are associated with dk2nu objects that contain
detailed information about their kinematics and their hadron ancestors.

Based on the ancestry information provided in the g4numi/g4lbnf output, we
extract the momentum and position of each pion/kaon that decayed to a
neutrino. In the NuMI beamline, the intersection of Horn 1’s front face and
the primary proton beamline’s trajectory marks the origin of the coordinate
system used in this simulation which begins with a 120-GeV primary proton
beam. The Z-axis points towards the beamline and the Y-axis points vertically
up. Geometrical features such as the target hall, decay pipe, absorber and
muon monitors are included in the simulation that provides the location and
kinematics of each decay into a neutrino. All particles in g4numi have a
default energy threshold of 0.02 GeV. The simulation cuts off protons outside
the beamline. In the muon monitor simulation, muons with energies less than 1
GeV, as well as other particles, are omitted. Using a 1 GeV muon cut for the
second stage of the simulation (i.e., the muon monitor simulation) saves
computing time and space. The threshold of the hadron absorber in front of the
MM1 is about 5 GeV. Therefore, cutting muons below 1 GeV is safe in the decay
pipe. g4lbnf simulation uses the same tracking cuts as g4numi.

##  III Simulation Technique

In this section, we first describe a technique to reproduce NuMI beam scans
with simulated data. The pions generated from g4numi simulations are over-
sampled to decay into multiple muons in this study. Then each muon is
simulated through the absorber and the muon monitors. In order to predict the
muon flux at muon monitors for selected beam and horn current settings, the
output of the g4numi simulation is combined with the output of the muon
monitor simulation. Using multiple muons as decay products is intended to
decrease statistical uncertainty in the analysis of each muon monitor pixel.
The simulation samples used here have all been generated in forward horn
current (FHC) or neutrino mode, during which the horn current flows in the
inner conductor in the same direction as the beam. As the first step of the
proposed technique, a uniformly distributed single simulation data sample is
generated for the selected beam scan variable. With the real beam width and
beam positions, the uniformly distributed sample is then used to generate
Gaussian samples for the selected Gaussian functions.

Figure 4: Random throws along the horizontal beam positions to generate a
uniform proton beam sample.

Using this method, we can run the simulation scan through a set of beam
variable values comparable to the actual beam scans. In this process, the
weights ($w_{i}$) are calculated based on selected two-dimensional normalized
Gaussian functions as shown in Eq. (1):

$w_{i}=\frac{1}{2\pi\sigma_{x}\sigma_{y}}\cdot\exp\left\\{-\frac{(x_{i}-\mu_{x})^{2}}{2\sigma_{x}^{2}}-\frac{(y_{i}-\mu_{y})^{2}}{2\sigma_{y}^{2}}\right\\}$
(1)

The $w_{i}$ is the calculated weight at the $i^{th}$ beam position
($x_{i},y_{i}$). The beam spot size ($\sigma_{x}$ and $\sigma_{y}$) and the
beam centroid position at the target ($\mu_{x},\mu_{y}$) have been fixed for
selected beam settings. These weights are applied during the process of
filling the histograms on the observable variables. Steps involved in the
proposed simulation technique are as follows:

  * •

Prepare uniformly distributed simulation samples on a selected beam scan
variable for a selected variable range.

  * •

Calculate the event weights according to a selected Gaussian function.

  * •

Apply the calculated weights on the observable variables.

  * •

Repeat the weight calculation algorithm by varying the Gaussian mean through
the selected beam scan range.

Every event from the uniformly distributed sample is contributing to each
selected Gaussian sample based on their weight distribution.

To demonstrate the technique, a uniformly distributed sample of protons on the
target as shown in Fig. 4 is used to generate possible interactions with
matter.

Figure 5: Two-dimensional horizontal vs. vertical proton beam positions with
uniform distribution sample for recorded beam interactions on the NuMI target,
that has a neutrino candidate at the downstream neutrino detectors.

From this uniformly distributed sample, a two-dimensional distribution of
horizontal and vertical proton beam positions on the NuMI target is drawn as
shown in Fig. 5, for each recorded beam interaction that has a neutrino
candidate at the downstream near detector. The recorded events include beam
interactions with the carbon target, aluminum baffle and downstream matter.

Figure 6: Two-dimensional horizontal vs. vertical proton beam positions with
uniform distribution sample for recorded beam interactions on the LBNF target
that has a neutrino candidate at the downstream neutrino detector.

The same exercise has been performed with the g4lbnf simulation in the forward
horn current mode to draw two-dimensional proton beam horizontal vs. vertical
positions on the LBNF target as shown in Fig. 6. A detailed description of the
LBNF optimized beam design, performed using the g4lbnf simulation has been
documented in the conceptual design report for the optimized LBNF beamline
[6]. As shown in the Table. 1 below, g4numi and g4lbnf use different nominal
beam spot sizes. As per NuMI and LBNF specifications, these values have been
chosen in all our studies since they represent the recommended beam spot sizes
for high-power beam operations at the respective beamlines.

|  Beam spot size (cm) $\sigma$
---|---
g4numi |  0.15
g4lbnf |  0.27
Table 1: Nominal beam spot size values for the g4numi and g4lbnf simulations.

Next, the Gaussian weights for different proton beam configurations are
calculated using Eq. (1).

Figure 7: Example of calculated weights for different beam settings with
varying beam width (top) and mean (bottom).

The top plot in Fig. 7 shows two examples of the calculated weights for a
given beam setting with mean $\mu$ = 0.0 cm and widths $\sigma$ = 0.10 cm and
$\sigma$ = 0.15 cm. The bottom plot in Fig. 7 shows two examples of the
calculated weights with beam $\sigma$ fixed at 0.15 cm and varying $\mu$
between 0.0 cm and 0.2 cm. The third step is to weight each event in the
distribution to select corresponding proton beam interactions horizontally and
vertically as shown in Fig. 8. Five Gaussian distributions as shown in Fig. 8
for the proton beam in the horizontal direction are drawn by changing the
proton beam mean $\mu$ from $-0.2$ cm to $0.2$ cm while keeping the Gaussian
width $\sigma$ fixed at 0.15 cm. The vertical proton beam distribution is
asymmetric since the target geometry is asymmetrical in y. In each Gaussian
profile, the weights are distinct.

Figure 8: The horizontal (top) and vertical (bottom) proton beam positions for
recorded beam interactions that have a neutrino candidate at the downstream
neutrino detectors.

This demonstrates that starting with a single uniformly distributed sample of
proton beam in both horizontal and vertical directions, for which there is a
recorded beam interaction, as many Gaussian distributions in either horizontal
or vertical proton beam direction can be selected. In the simulation, pions
decay into multiple muons, which are propagated through the absorber and the
muon monitors. The corresponding weights have been applied to the observed
muon flux from the Gaussian proton beam slice with $\mu=0.0$ cm to obtain the
two dimensional event distributions on the muon monitors as shown in Fig. 9.

Figure 9: The observed two-dimensional muon distribution on NuMI muon monitor
1 for a given beam $\mu=0.0$ cm and $\sigma=0.15$ cm.

The uniform beam simulation allows us to generate as many Gaussian samples as
we require. It is possible to draw samples with means $\mu_{i}$ with i ranging
from $-\infty$ to $+\infty$. Figure 10 illustrates an example case with three
Gaussian samples with means of $\mu_{-1}$, $\mu_{0}$, and $\mu_{+1}$. The
$\mu_{-1}$ value in this example is $-0.1$ cm, the $\mu_{0}$ value is 0.0 cm,
and the $\mu_{+1}$ value is +0.1 cm. For all three samples, the width $\sigma$
has been fixed at 0.15 cm. Our reference Gaussian is the Gaussian sample with
mean $\mu_{0}$. By using the uniform beam method, we start with a static pool
of randomly thrown events interacting with matter to produce neutrinos.
Gaussian samples are generated based on weights by selecting events from the
same static pool. However, this is not a concern, because we give different
statistical weights to each event. These statistical weights determine with
what probability these events will contribute to the Gaussian samples. Within
one Gaussian sample, no two events will have the same statistical weight, as
the weight is calculated using Eq. (1). In the case of two Gaussian samples,
there might be a scenario in which the same event could be shared with the
same statistical weight between these two samples as shown in Fig. 10.

Figure 10: Two selected Gaussian profiles from the uniform simulation, showing
the volume of overlap between them.

In Fig. 10, we can see that by plotting two Gaussian samples with
$\mu_{-1}=-0.1$ cm and $\mu_{+1}=+0.1$ cm, an overlapping volume is created
between these two Gaussian samples, as shown by the shaded blue region. It is
possible for two events belonging to the above two Gaussian samples to have
the same weight in this overlap region, but they are really two different
events since they each have a unique beam position and will interact with
matter differently. Only the event at proton beam position $x_{0}$, shown by
the black dotted line will have the same interactions. So, it is only the
event shown by the black dotted line which will be shared by the Gaussian with
$\mu_{-1}=-0.1$ cm as well as the Gaussian with $\mu_{+1}=+0.1$ cm. Within one
generated Gaussian sample (for example the Gaussian with mean $\mu_{-1}$), we
are using the same event only once.

Equation (2) estimates the ratio of the overlapping volume $R_{i}$ between a
Gaussian sample with mean $\mu_{i}$ and a given beam width $\sigma$ as
compared to the reference Gaussian sample with mean $\mu_{0}$.

$R_{i}=\frac{2\int_{-\infty}^{x_{0}}\mathcal{N}(\mu_{i},\,\sigma^{2})dx}{\int_{-\infty}^{+\infty}\mathcal{N}(\mu_{0},\,\sigma^{2})dx}$
(2)

Figure 11 plots $R_{i}$ as a function of $\mu_{i}$ indicating that the
overlapping volume will be maximum when the mean of the $i$-th Gaussian aligns
with the mean of the central Gaussian slice (with $\mu_{0}$). As the mean of
the Gaussian $\mu_{i}$ moves away from $\mu_{0}$, the overlapping volume
between two consecutive Gaussian samples will decrease according to Eq. 2. The
larger $|\mu_{i}-\mu_{0}|$ is, the smaller the overlap will be, reducing the
weight of the shared event between two consecutive Gaussian samples.
Consequently, the statistical contribution of this shared event will decrease
according to Eq. 2.

Figure 11: Ratio of the overlapping volume of Gaussian beam profiles for
different $|\mu-\mu_{0}|$.

Using Figure 11, we can determine the overlap between two Gaussian samples at
a given separation between their means. By separating Gaussian means by
2$\sigma$, for example, the overlap volume is reduced from 100% to 30%,
compared to the case when the samples have the same mean.

##  IV Simulation Data

Validation and demonstration of the proposed technique has been performed by
generating several data samples. We have generated a uniformly simulated beam
sample of 1 billion protons on target (POT). Using FermiGrid, this sample was
the largest we were able to produce under the limit of available storage per
user. FermiGrid is a collection of computing resources that Fermilab Computing
Division makes available through grid protocols. It would be possible to
create multiple 1 billion uniform beam samples by utilizing supercomputers.

While real measurements utilize the three muon monitors to measure integrated
flux of muon beams that exceed 5 GeV and two other energies separated by at
least 7 GeV, simulated data samples can only be generated using muon monitors
1 and 2 as computing resources are limited. Simulating muon monitor 3
responses will require supercomputer resources in order to generate large
statistics.

In each Gaussian beam sample, 250 million POTs were generated with different
beam positions while maintaining a large beam spot size of 0.15 cm. Based on
the beam spot size ($\sigma$ = 0.15 cm) and the cross sectional area of the
target and the baffle, the vertical and horizontal limits of the uniform beam
simulation data are chosen to be 1.0 $\times$ 1.0 cm perpendicular to the beam
direction. The uniform beam window was chosen to cover approximately $7\sigma$
of the Gaussian beam, which allowed us to incorporate the interaction of the
protons with matter in the Gaussian beam tail. The target and a significant
portion of the baffle are covered by this limit. In order to prepare the
horizontal and vertical beam scan data samples, Gaussian beams with mean
values in the range of $[-0.2,0.2]$ cm are generated from uniform beam
simulation data. These beams cover a minimum of $5~{}\sigma$ beam profiles for
each beam configuration. The randomly thrown protons in the uniform beam
simulation interact with matter to produce secondary particles and are only
recorded if there is a neutrino candidate at the near detector from the hadron
decay. For each proton beam position along the horizontal direction, Gaussian
samples were created by applying weights with a beam width of 0.15 cm. As a
result, every Gaussian sample contains $\sim$250 million POTs. This is
estimated based on the generated 1 billion POTs in the uniform beam sample.
By applying the Gaussian weight to random throws of $N$ POTs, we calculated
statistically that we would select $\sim$ $\frac{1}{4}$ of the starting N POTs
in one Gaussian beam. In order to verify the calculation, we generated one
uniform beam profile with 1 billion POTs from which we drew Gaussian beam
profiles by applying weights. Our calculation shows that each of these
Gaussian beam profiles will contain 250 million POTs. Following this, we
compared these Gaussian beam profiles with Gaussian beams from nominal beam
simulations each containing 250 million POTs with the same beam width of 0.15
cm. Therefore, we would generate a uniform simulation sample just once with
four times as many POTs as a nominal simulation sample. On the other hand, we
would have to generate a new Gaussian with 250 million POTs every time if we
were to generate a new reference simulation sample with different beam
parameter settings.

Figure 12 shows a two-dimensional distribution of total hadron momentum vs.
hadron angle at decay for a Gaussian slice of 250 million POT with a beam mean
$\mu$ = 0.0 cm. By weighting the uniform simulation sample of 1 billion POT,
this Gaussian slice was created. The same two-dimensional distribution of
total hadron momentum vs. hadron angle at decay is plotted for a Gaussian
slice of 250 million particles with a beam mean $\mu$ = 0.0 cm using the
nominal simulation, as shown in Fig. 13. These two plots only show events for
which muon candidates with energies greater than 5 GeV are present in muon
monitor 1.

Figure 12: Distributions of hadron angle and hadron momentum in the
longitudinal direction at decay with uniform simulation for a Gaussian beam
with $\mu$ = 0.0 cm. Figure 13: Distributions of hadron angle and total hadron
momentum at decay with nominal simulation for a Gaussian beam of $\mu$ = 0.0
cm.

##  V Validation

Validation studies of the proposed beam simulation technique are an important
requirement to prove the capability of reproducing the nominal simulation
results. Therefore, two studies with the NuMI simulation have been carried out
to validate the uniform beam technique against the nominal simulation:

  * •

Testing statistical reproducibility.

  * •

Testing kinematics reproducibility.

In the validation studies, the uniform beam technique and multiple decays of
hadrons were combined to increase the statistics of the muon events in the
simulation. The uniqueness of the uniform beam simulation is in that we can
draw as many Gaussian beam profiles horizontally or vertically.

###  V.1 Testing Statistical Reproducibility

We generated five Gaussian beams using the nominal simulation, each containing
250 million POT, and one uniform simulation containing a total of 1 billion
POT to compare beam profiles. We drew five Gaussian weighted proton beam
samples from the uniform simulation sample with beam centers at -0.2 cm, -0.1
cm, 0 cm, 0.1 cm, and 0.2 cm along the horizontal direction and compared them
to the nominal Gaussian beams.

Figure 14: Comparison of the proton beam profiles drawn from uniform and
nominal simulations in the horizontal direction with the mean spanning from
$-0.2$ to $+0.2$ cm. All plots shown are area normalized.

Figure 14 shows a comparison of the proton beam profiles from the uniform beam
simulation (cyan) and Nominal beam simulation (red) at five different beam
positions along the horizontal direction. The beam profiles have been plotted
only for the incident protons that have corresponding muon candidates from the
decay of hadrons which were created by the proton-matter interactions. The
uniform and nominal simulations match up with all five beam profiles. The
bottom plot shows the ratio of two beam profiles for the mean of 0.0 cm based
on uniform and nominal simulations. A ratio of one confirms an excellent
agreement.

Figure 15 illustrates an example of NuMI muon monitor 1 response with a given
Gaussian beam generated from the uniform beam simulation. The same Gaussian
beam was also generated using the nominal beam simulation. A ratio between
those two simulations for each pixel has also been shown.

Figure 15: NuMI muon monitor 1 responses from two different simulations:
uniform (top) and nominal (middle). Ratio of the muon monitor 1 responses
between the uniform and the nominal simulations (bottom).

We plotted muon monitor responses only for incident protons decaying into muon
candidates. Statistical fluctuations can explain a maximum difference of
$\sim$4% between two simulations in some of the edge pixels of a muon monitor.
We demonstrated that the statistical fluctuations affect muon monitor pixel
responses even within the same nominal simulation. We generated five different
nominal simulation samples using five unique random seeds and then calculated
the ratio of the muon monitor 1 responses between two Nominal simulation
samples with two different seeds. It has been shown that even within the same
nominal simulation, some edge pixels fluctuated based on the random seed by
$\sim$ 4% at most.

Figure 16: NuMI muon monitor 1 responses from nominal simulations with two
different seeds (top and middle, respectively). Ratio of the muon monitor 1
responses between nominal simulations with two different seeds (bottom).

It is consistent with the difference in edge pixel response between uniform
and nominal simulations on muon monitor 1. Hence, the $\sim$ 4% difference is
not an artifact of the uniform simulation technique, rather stems from
statistical fluctuations. Increasing statistics will reduce this fluctuation.
Figure 16 top distribution shows muon monitor 1 pixel responses for two
nominal simulation samples created with different random seeds, while the
bottom distribution shows the ratio between the two simulations. There is a
maximum difference of $\sim$ 4% in one of the corner left edge pixels.

Finally, Fig. 17 illustrates the NOvA Near Detector neutrino energy spectrum
plotted from nominal and uniform simulations. As can be seen from the bottom
plot, in the region of interest, the ratio of the two spectra stays close to 1
for the energy range of interest.

Figure 17: Neutrino energy spectra at the NOvA Near Detector from nominal and
uniform beam simulations (top). Ratio of the neutrino energy spectra at the
Near Detector between nominal and uniform simulations (bottom).

###  V.2 Testing Kinematics Reproducibility

Figure 18: A comparison of the muon momenta along the X (top),Y (middle) and Z
(bottom) axes produced at the hadron decay with uniform (red) and Nominal
(blue) simulations for a Gaussian beam of $\mu$ = 0.0 cm.

We have tested the reproducibility of the kinematics of the secondary
particles after applying the uniform beam technique. Pion decay kinematics and
decay pipe size determine the size of the muon beam at the Muon Monitors. In
this study, we have used 0.15 cm width Gaussian beam samples, with mean set at
$\mu$ = 0.0 cm from uniform and nominal simulations to compare the muon
momentum distributions at the hadron decay along $x$, $y$, and $z$ directions
at the muon monitor 1. Figure 18 shows the distributions of the momenta of the
muons produced at the hadron decay. In the red distributions, the momenta are
plotted using the uniform beam simulation, whereas in the blue distributions,
they are plotted using the nominal simulation. Muon production momenta are
shown to be consistent between uniform and nominal simulations, which is
validated by the ratio 1 between the momenta distributions from these
simulations.

##  VI Computing Resources

As a result of using uniform beam simulations to generate large simulation
samples, the amount of computing resources needed to run the simulations can
be reduced considerably. Here we show a comparison between the amount of
computing resources used to generate nominal Gaussian and uniform beam
simulation samples, which is shown in Table. 2. The table shows that running a
nominal simulation to generate one random sample takes about an hour, whereas
running a uniform simulation to generate four times as many statistics takes
almost the same time. The fractions of protons interacting with matter are
smaller in uniform beam simulation than in nominal simulation. For uniform
simulations, we create Gaussian beam post-processing. Using nominal
simulation, however, it would take approximately 1000 hours of CPU time to
generate 1000 random samples. Nevertheless, if the proposed weighting
technique is applied to a high POT single uniform beam simulation sample, the
CPU time required to generate 1000 random Gaussian samples can be
significantly reduced.

|  POT |  $<$CPU Time$>$ (min) |  $<$Memory$>$ (MiB)
---|---|---|---
Nominal |  250M |  59 |  51.2
Uniform |  1000 M |  67 |  52.7
Table 2: A comparison of the average CPU time and memory usage required to
generate one sample of uniform beam simulation versus one sample of nominal
beam simulation.

##  VII Applications of the Technique

Simulation samples were generated for selected proton beam position
configurations with different horn current settings to test the uniform beam
simulation technique. This section presents the simulation technique for
scanning targets by generating uniform beams along horizontal and vertical
positions at the target for 200 kA and 180 kA horn currents, respectively.
Simulated samples are generated using the combined g4numi and muon monitor
packages.

###  VII.1 Total Muon Scan With and Without the Target

To check the alignment of the target-baffle system, horizontal and vertical
target scans are generally performed with the primary proton beam at the
beginning of a run period. To align the horn, the target is removed and the
beam is scanned along the horn cross-hair to observe the responses on the
hadron monitor. A similar scan is presented here using our proposed technique
to generate 100 Gaussian beam profiles with $\sigma=0.15$ cm along the
horizontal axis from -0.8 cm to +0.8 cm.

Figure 19: Total number of muon on muon monitor 1 from the beam scan with
including the target (top) and removing the target (bottom).

The demonstration has been performed with and without the target in the
beamline. The horn current has been set to +200kA for this study. Muon monitor
1 shows a variation in the total number of muons as a function of beam
position on the target. Fig. 19 shows the total number of muons detected at
muon monitor 1 during the beam scan with (left) and without (left) the target.
By simulating the interactions of the beam with the baffle, target, and other
beamline materials, we are able to observe muons at muon monitor 1. A similar
study of the target scan with the NuMI beam has been described in chapter 10:
Beam Monitoring Measurements of the Accelerator Systems and Instrumentation
for the NuMI Neutrino Beam of reference 10.[10].

###  VII.2 Horizontal and Vertical Beam Scan

To demonstrate horizontal and vertical beam scans, we have generated Monte
Carlo samples each with 1 billion protons on target (POT) for +200 kA horn
current setting. We generate horizontal beam scan samples by throwing protons
uniformly along the horizontal axis from -1.0 cm to +1.0 cm while keeping the
vertical beam as a random Gaussian with a 0.15 cm beam width and zero
centroid. In order to generate a vertical beam scan sample, the axis
selections for the uniform beam and Gaussian beam are changed. Each iteration
uses a centroid separation of 0.02 cm and a Gaussian width of 0.02 cm. The
Gaussian width of the proton beam has been set to 0.15 cm for the current
scan. To verify the uniform beam simulation technique, we compare the results
of the uniform simulation beam scan with the results of the nominal simulation
beam scan. Fig. 20 shows the responses of the muon flux centroid along the
horizontal (top) and the vertical (bottom) axes at the muon monitor 1 (red)
and 2 (blue) as a function of the proton beam positions for +200 kA horn
current plotted with the uniform beam simulation.

Figure 20: The muon flux centroid estimated with uniform beam simulation and
nominal simulation respectively, on muon monitor 1 (uniform: red, nominal:
green) and muon monitor 2 (uniform: blue, nominal: purple) in a horizontal
beam scan (top) and vertical beam scan (bottom) for a +200 kA horn current as
a function of the proton beam position at the target.

The equivalent plots with the nominal simulation are shown at muon monitor 1
(green) and 2 (purple). The muon flux centroid on each muon monitor is
calculated by taking the average of the projected distribution on the
horizontal axis from 81 muon monitor pixels of ionization chamber signals. The
distributions have been fitted to a linear function. The fitted slopes from
the uniform beam simulation are compared with the slopes calculated using
nominal simulation. Another example of the muon monitor response to the
primary beam position variation can be demonstrated as the ratio of the muon
flux at muon monitor 1 for proton beam at $\mu$ = -0.2 cm to the beam at $\mu$
= +0.2 cm in the horizontal direction as shown in fig. 21. By observing this
ratio we can see the muon flux behavior at the muon monitor 1 due to the horn
focusing mechanism at work [11]. The same flux behavior is observed in both
uniform and nominal simulations. It can be seen that the muon flux moves from
right to left on muon monitor 1 as the primary beam moves from left to right.
This is due to the fact that hadrons are over-focused at muon monitor 1 hence
produces the negative slope as shown in fig. 20 (top).

Figure 21: Ratio of muon flux in each pixel on muon monitor 1 for proton beam
$\mu$ = -0.2 cm and $\mu$ = +0.2 cm respectively in the horizontal direction
for a +200 kA horn current. Uniform beam simulation (top). Nominal beam
simulation (bottom).

In order to repeat the beam scan studies for a different horn current setting
of +180 kA, we repeat the same analysis procedure described earlier. Fig. 22
shows the muon flux centroid along the horizontal (left) and vertical (right)
axes at muon monitor 1 as a function of the proton beam positions for a horn
current of +180 kA, using both the uniform beam and nominal simulations.

Figure 22: The uniform simulation estimates of the muon flux centroid along
the horizontal (top) and vertical (bottom) projections on muon monitor 1 for
180 kA horn current setting as a function of the proton beam position at the
target.

Two different horn current settings produced consistent slopes in both
horizontal and vertical scans of both simulations.

###  VII.3 Beam Spot Size Scans

In order to establish a correlation between neutrino and muon beam profiles,
it is essential to understand how beam spot size affects muon monitor
observations.

|
---|---
|
Figure 23: The beam spot size selections: $\sigma=0.08$ cm (top left) and
$\sigma=0.15$ cm (top right) and corresponding muon monitor responses for
$\sigma=0.08$ cm (bottom right) and $\sigma=0.15$ cm (bottom left).

With the uniform beam simulation technique, we have produced Gaussian beams by
gradually changing the Gaussian width $\sigma$. We have generated a simulation
sample with 1 billion POT of uniformly distributed proton throws along
horizontal and vertical directions ranging from -1.0 cm to +1.0 cm.

Figure 24: The muon flux centroid estimated along the horizontal projection on
Muon monitor 1 and 2 for 200 kA for two different beam spot sizes: 0.15 cm and
0.08 cm.

We have studied the muon flux distributions at monitor 1 for beam spot size
change of 0.08 cm, and 0.15 cm. Fig. 23 shows the beam spot sizes of 0.08 cm
(top left) and 0.15 cm (top right) and the corresponding muon monitor
responses in the contour plots (bottom left and right) respectively. Fig. 24
shows the muon flux centroid estimation for different proton beam positions
with changing beam spot sizes. According to these observations, muon flux
centroids and standard deviations do not change as beam spot size changes.
Simulation studies indicate that changing the beam spot size will not affect
muon monitor measurements.

##  VIII Application in LBNF simulation

Neutrino beam instrumentation (NBI) for LBNF provides measurements of the
secondary beam for commissioning, alignment, monitoring, and hardware
protection. This instrumentation complements the primary beam instrumentation
as well as the neutrino detectors. It includes a hadron alignment detector
system (HADeS) for measuring the remaining secondary particles in the Decay
Pipe and a muon monitor detector system (MuMS) for monitoring muons downstream
of the Hadron Absorber. HaDeS is based on the NuMI hadron monitor [7], which
will be used to align the beamline components by analyzing the residual beam
after scanning the primary beam across the target, baffle, and horn. The NuMI
hadron monitor consists of a linear array of parallel-plate ionization
chambers of the size of 1 m2, $7\times 7$ array. HaDeS will implement LBNF-
specific modifications in terms of higher channel counts to account for
smaller pixel sizes compared to NuMI. During beam commissioning, beam
alignment, and diagnosing failures, the HaDeS will be lowered into the beam
for low-intensity alignment scans. A number of purposes are accomplished by
the MuMS, including monitoring of long-term target degradation under beam
irradiation by comparing muon fluxes at different energies. In addition to
monitoring the condition of a target in detail, this monitor bypasses the need
to wait for neutrinos or other tertiary beam monitors. Similar detector
technology will be used for MuMS as for HaDeS, following NuMI’s lead.

The LBNF beam simulation (g4lbnf) is a dedicated Geant4 beamline simulation
based on the optimized design of the beamline. In the nominal g4lbnf
simulation, a Gaussian beam is used with beam $\sigma$ set at 0.27 cm. We have
implemented the uniform beam simulation technique in g4lbnf and drawn five
Gaussian distributions for the proton beam in the horizontal direction from
the uniform beam by changing the proton beam mean $\mu$ from -0.2 cm to +0.2
cm while keeping the Gaussian width $\sigma$ fixed at 0.27 cm as shown in Fig.
25 top. The plots shown below as a demonstration of the uniform method have
been generated using a simulation sample that contains 1000 jobs, each
containing 10000 POT. This is just an example. It should be noted that this is
only an example. POT can be simulated in different amounts, depending on the
need of the user. The bottom two plots here show the distributions of the
charged particles in a virtual detector at the location of the HaDeS for the
Gaussian beams with $\mu$ = -0.2 cm, and, 0.0 cm respectively. Only those
events for which a beam interaction has been recorded are included in the
Gaussian distribution.

Figure 25: The horizontal proton beam positions (top) for recorded beam
interactions that have a neutrino candidate at the downstream neutrino
detectors. Beam $\sigma=0.27$ cm. Distribution of charges particles at the
LBNF HaDeS location (middle) for a Gaussian proton beam with $\mu$ = -0.2 cm.
Equivalent distribution (bottom) for a Gaussian proton beam with $\mu$ = 0.0
cm.

During beam-based alignment, low intensity beam will be scanned across the
face of the HaDeS to determine beam direction. In order to determine if the
pixel map of the hadron monitor is accurate, we can compare the HaDeS
distributions from simulation with data.

---
Figure 26: Top: the Gaussian slices of the LBNF proton beam with mean at -0.2
cm, 0.0 cm and +0.2 cm from left to right. Bottom: muon beam profiles at the
location of alcove 1 for the corresponding proton beams.

In Fig. 26, the top plot illustrates how Gaussian distributions are drawn
horizontally after applying weights, only for the events for which beam
interactions are recorded. On a virtual detector placed downstream of the
absorber at the location of MuMS in muon alocve 1, the bottom plot shows the
muon beam profiles after applying the same corresponding weights for each of
the Gaussian beams.

Each MuMS system detects muons that are focused differently depending on the
energy of the parent pion. The two-dimensional event distribution on the
alcove 1 tracking plane from a Gaussian proton beam slice with mean $\mu$ =
0.0 cm is shown in Fig. 27, left. The right plot shows the muon flux centroid
estimation in a horizontal proton beam scan in the alcove 1 location. The horn
current is 300 kA. The ability to perform beam scan studies using the uniform
beam simulation technique will be valuable to understand impacts of the
changing LBNF beam parameters on the beam response measurements.

Figure 27: The observed two dimensional muon distribution on a tracking plane
at the LBNF alcove 1 location (top) for a given proton beam with $\mu=0.0$ cm
and $\sigma=0.27$ cm. Muon flux centroid estimated at the alcove 1 location as
a function of the horizontal proton beam position (bottom) at the target.

##  IX Application of machine Learning: predict neutrino flux

The application of Machine learning (ML) algorithms is becoming increasingly
critical in particle accelerators and beamlines. Due to their efficiency and
flexibility, artificial neural networks are ideally suited for complex
simulations. A multi-layered artificial neural network (ANN) consists of an
input layer, an output layer, and multiple hidden layers, each layer
containing a number of nodes that are connected by randomly assigned weights.
Every node in the layers produces an output according to its “activation
function”. The neural networks are then trained through repeated iterations to
reduce estimated errors in the output by adjusting the weights of the
connections between the nodes. Besides modeling nonlinear behavior, neural
networks can also adapt to changes in the system over time.

In addition to being useful for Fermilab’s current neutrino beamlines, the
development of machine learning algorithms for beam quality monitoring,
anomaly detection, neutrino beam systematic studies, and neutrino beam quality
monitoring will also be useful for future beamlines such as LBNF. In order to
understand some of the rare “anomaly” scenarios in the NuMI beamline, we need
to generate simulated data. For example, g4numi simulation is used to generate
simulated data for ML applications in the NuMI beamline. We can predict some
“hard-to-measure” beam parameters ML algorithms based on simulations such as
Horn tilting angle, target offset. The muon monitors are sensitive to changes
in the primary beam and variations in horn current. Machine learning
applications can be built to monitor the quality of the NuMI neutrino beam
based on the unique responses of muon monitors. The fig. 28. below illustrates
how different ML algorithms can be used to predict different beam parameters
and neutrino fluxes based on muon monitor pixel information.

Figure 28: An example of application of Machine Learning Algorithms with Muon
Monitor pixel inputs.

A uniform beam simulation sample was used for the training data samples, while
three independent uniform beam simulation samples were used as the test data
samples. There are 6000 points in the training data sample, each containing
beam parameters, muon events at pixels, and neutrino flux at the near
detector.

Muon monitor observations depend on a number of parameters. Taking pixel
information from the muon monitor and disentangling it from other correlations
is important. Simulated data should be divided into training and test sets. To
understand correlations between muon monitor observations and beamline
incidents, we need MC samples with large statistics, $\sim$300 million POTs in
each sample. In the past, the bottleneck has been obtaining large simulated
samples running the nominal g4numi simulation, due to the lengthy sample
preparation process. Creating a nominal g4numi simulation sample with more
than $\sim$ 300 million events takes $\sim$ 1 day.

Since we need to scan over parameter space, we need many large samples. Each
time we change a beam parameter, we need to run the simulation again with
large statistics for each change. Solutions such as GPU processing have been
considered, but it becomes non-trivial as Geant4 does not support GPU
processing. CPU processing on the DOE supercomputer National Energy Research
Scientific Computing Center (NERSC) with 16000 cores would be an easier
option. However, by using uniform beam simulation, we can create large
simulation samples more quickly and efficiently.

A Monte Carlo data sample with uniform beam simulation has been generated for
ML applications to demonstrate the power of uniform beam simulation. The
uniform beam simulation and multiple decay simulation have been combined to
increase simulation efficiency.

Figure 29: Top: Horizontal beam position prediction, bottom: horn current
prediction using a linear regression model.

In this technique, we generate a uniformly distributed single simulation data
sample on the grid with large statistics for the selected beam variable range.
The uniformly distributed sample is then used to generate Gaussian beam
profiles for different selected beam parameters. Beam positions, widths etc.
can be varied post-processing once we have generated the uniform sample. In
this way, computational overhead is greatly reduced. Based on the uniform beam
simulation data, a linear regression model was tested to predict the
horizontal proton beam position and the horn current. To train the linear
model, we used muon monitor pixels as inputs. The model shows high accuracy in
predicting beam position (left) and horn current (right) after training as
shown in fig. 29.

The ML algorithm and uniform beam simulations have also been shown to be
successful in predicting neutrino flux. We have predicted neutrino flux at the
NOvA Near detector using simulated muon events from 81 pixels from each muon
monitor. A neural network has been used to make prediction of the NuMI
neutrino flux. Fig. 30 shows the predicted and true neutrino energy spectra at
the NOvA Near detector where the prediction matches the training sample well.

Figure 30: A comparison of predicted (blue) and true (orange) neutrino energy
spectra at the NOvA Near detector.

The ML predictions provide an additional layer of monitoring of neutrino beam
behavior and horn current behavior. The results show potential for monitoring
beam performance, developing trends, and identifying issues during regular
beam operations. The uniform beam simulation has enabled ML models to be used
to predict rare incidents and anomalies due to its efficiency and cost-
effectiveness in generating large simulation samples.

##  X Summary and Outlook

In this paper, we have presented a simulation technique that can be applied to
generate as many random beam simulation options desired for beam scan studies.
In contrast to generating multiple Gaussian samples with very large
statistics, this technique produces large simulation samples with a uniform
distribution of beam parameters. From a test sample of 1 billion POT with
uniformly distributed beam X and Y with g4numi simulation, each generated
Gaussian beam with 0.15 cm beam width will contain 540 million weighted
interactions. For the same number of interactions, a sample of 250 million POT
must be generated every time from a nominal sample. Therefore, instead of
generating multiple Gaussian beams with 250 million POT every time, one can
generate a uniform simulation with 1 billion POT once and generate as many
Gaussian beams from it as desired. To generate a single Gaussian beam with 250
million POT with nominal simulation takes more than 2 hrs of wall time on the
Fermigrid parallel processing resource at Fermilab computing facility. Grid-
submitted jobs are subject to varying and lengthy waiting times. By using
uniform beam simulation with high statistics, one can generate weights for
many beam simulation combinations without requiring a large computing resource
and a long simulation time. In order to demonstrate the application of the
uniform simulation technique to NuMI studies, a number of examples have been
presented.

A similar exercise has been performed with g4lbnf simulation, and the studies
have been extended to show that the technique can be used to study LBNF beam
scans through simulation. The results of the uniform simulations are shown to
be consistent with the nominal simulation. A combination of uniform beam
simulations and multiple decays has been shown to be applicable to multiple
neutrino beamline simulations. The high statistics of Gaussian throws from the
uniform simulation technique can be leveraged not only for neutrino
experiments but also for various simulation requirements related to the
upcoming Mu2e experiment at Fermilab or its upgrade, known as Mu2e-II, which
is the next generation muon conversion experiment. Due to its efficiency and
cost-effectiveness in generating large simulation samples, the utilization of
the uniform beam simulation technique has facilitated the use of machine
learning models for predicting rare incidents and anomalies. Consequently,
this technique will play a crucial role in the application of machine learning
algorithms across diverse research contexts.

###### Acknowledgements.

This work is supported by the Fermi Research Alliance, LLC manages and
operates the Fermi National Accelerator Laboratory pursuant to Contract number
DE-AC02-07CH11359 with the United States Department of Energy. This work is
partially supported by the U.S. Department of Energy grant DE-SC0019264.

## References

  * [1] P. Adamson, et al., First measurement of electron neutrino appearance in NOvA, Phys. Rev. Lett. 116 (15) (2016) 151806.  arXiv:1601.05022, doi:10.1103/PhysRevLett.116.151806.
  * [2] K. Abe, et al., Letter of Intent: The Hyper-Kamiokande Experiment — Detector Design and Physics Potential — (9 2011).  arXiv:arXiv:1109.3262.
  * [3] K. Abe, et al., Measurements of neutrino oscillation in appearance and disappearance channels by the T2K experiment with 6.6×1020 protons on target, Phys. Rev. D 91 (7) (2015) 072010.  arXiv:1502.01550, doi:10.1103/PhysRevD.91.072010.
  * [4] R. Acciarri, et al., Long-Baseline Neutrino Facility (LBNF) and Deep Underground Neutrino Experiment (DUNE): Conceptual Design Report, Volume 2: The Physics Program for DUNE at LBNF (12 2015).  arXiv:1512.06148.
  * [5] G. Barenboim, et al., Detector R and D for Future Neutrino Experiments with the NuMI Beamline (11 2002).  arXiv:hep-ex/0304017.
  * [6] P. Adamson, et al., The conceptual design report for the optimized neutrino beamline (2017).
URL https://docs.dunescience.org/cgi-bin/ShowDocument?docid=4559

  * [7] R. Zwaska, et al., Beam-Based Alignment of the NuMI Target Station Components at FNAL, Nucl. Instrum. Meth. A 568 (2006) 548–560.  arXiv:physics/0609106, doi:10.1016/j.nima.2006.08.031.
  * [8] S. Agostinelli, et al., GEANT4–a simulation toolkit, Nucl. Instrum. Meth. A 506 (2003) 250–303.  doi:10.1016/S0168-9002(03)01368-8.
  * [9] J. Allison, et al., Geant4 developments and applications, IEEE Trans. Nucl. Sci. 53 (2006) 270.  doi:10.1109/TNS.2006.869826.
  * [10] R. M. Zwaska, Accelerator Systems and Instrumentation for the NuMI Neutrino Beam, Ph.D. thesis, Texas U. (2005).  doi:10.2172/879065.
  * [11] K. Yonehara, S. Ganguly, D. A. Wickremasinghe, P. Snopok, Y. Yu, Exploring the focusing mechanism of the numi horn magnets (2023).  arXiv:2305.08695.


Journal of Web Engineering, Vol. 0, No. 0 (2003) 000–000
©Rinton Press

1

RUN-TIME APPLICATION MIGRATION

USING CHECKPOINT/RESTORE IN USERSPACE

Aleksandar Tošić

University of Primorska Faculty of Mathematics, Natural Sciences and
Information Technologies, Koper

InnoRenew CoE, Livade 6a, Izola

Koper, Slovenia

<EMAIL_ADDRESS>
Received (20.7.2023)
Revised (revised date)

This paper presents an empirical study on the feasibility of using
Checkpoint/Restore In Userspace (CRIU) for run-time application migration
between hosts, with a particular focus on edge computing and cloud
infrastructures. The paper provides experimental support for CRIU in Docker
and offers insights into the impact of application memory usage on checkpoint
size, time, and resources. Through a series of tests, we find that the time to
checkpoint is linearly proportional to the size of the memory allocation of
the container, while the restore is less so. Our findings contribute to the
understanding of CRIU’s performance and its potential use in edge computing
scenarios. To obtain accurate and meaningful findings, we monitored system
telemetry while using CRIU to observe its impact on the host machine’s CPU and
RAM. Although our results may not be groundbreaking, they offer a good
overview and a technical report on the feasibility of using CRIU on edge
devices. This study’s findings and experimental support for CRIU in Docker
could serve as a useful reference for future research on performance
optimization and application migration using CRIU.

Keywords: Checkpoint/Restore, Edge computing, Run-time container migration

Communicated by: to be filled by the Editorial

##  1 Introduction

In recent years, there has been a significant shift in the way computing is
done. Instead of relying solely on centralized cloud-based systems, there has
been a move towards decentralized edge computing systems. Edge computing
systems bring computing resources closer to the end-users, reducing latency,
improving performance, and enabling new use cases that were previously
unfeasible. However, with edge computing comes a new set of challenges, such
as resource constraints, high variability, and unreliable connectivity.

Checkpoint/Restore In Userspace (CRIU) is a powerful tool that allows the
migration of running applications from one host to another while preserving
their states. This tool has significant potential for use in edge computing
systems, where migration of running applications may be necessary due to
resource constraints or other factors. The ability to migrate running
applications while preserving their state can significantly improve the
reliability, availability, and fault-tolerance of edge computing systems. A
few potential use case pertaining to the aforementioned use case are:

  1. 1.

Fault Tolerance and Disaster Recovery: CRIU can be used to provide fault
tolerance and disaster recovery capabilities in edge computing systems. By
migrating running applications to a different host in the event of a failure
or disaster, CRIU can help ensure that critical services remain available and
reduce the impact of downtime.

  2. 2.

Mobility: CRIU can be used to provide mobility capabilities in edge computing
systems. By allowing running applications to be migrated between different
hosts, CRIU can enable new use cases such as mobile edge computing, where
services follow users as they move around.

  3. 3.

Dynamic Scaling: CRIU can be used to enable dynamic scaling of services in
edge computing systems. By migrating running applications between hosts based
on changes in demand, CRIU can help ensure that the system remains responsive
and that resources are used efficiently.

  4. 4.

Energy Efficiency: CRIU can be used to enable energy-efficient computing in
edge computing systems. By migrating running applications to hosts that are
currently idle, CRIU can help reduce overall energy consumption and improve
the sustainability of the system.

One of the main points of interests is the advent of decentralized protocols
for autonomous container management and orchestration. Protocols such as
Caravela [5], and Nion Network [6] proposed a decentralized network of nodes
that perform container migrations in an effort to balance resource utilization
across the entire network. Nion network secures the decision making process
through an efficient and scalable blockchain consensus mechanism. Unlike
typical blockchain protocols, the blocks are viewed as states descriptions of
all the containers under management. Adding new blocks to the chain can be
considered as the state transition, while the state transition function is
nodes performing planed migrations of containerized applications.

Clearly, such networks are subject to Byzantine behaviour and while consensus
mechanisms blockchains typically use tolerate such faults, the applications do
not. These faults are traditionally addressed by introducing redundancy where
possible and/or backup. In this use-case, redundancy would require more nodes
to run the containerized application and serve as ready backups. Evidently,
this does raise concern over the efficiency of the entire network due to
additional compute resource utilization. On the other hand, backups only
impose additional storage requirements but have to be frequent to avoid loss
loosing states.

CRIU can be used to regularly snapshot applications and commit these snapshots
to other nodes. In an event of a fault on the host machine(node), the protocol
can detect it, and restore the application on another node. However, given the
experimental nature of CRIU, the exact resource requirements for creating
checkpoints and restoring them needs to be well understood. Moreover, studying
the telemetry while performing such operations can give a better understanding
into the limitations of run-time migrations.

In this paper, we explore the potential benefits of using CRIU in edge
computing systems using docker experimental features that support CRIU [3, 7].
We present an empirical study of the performance impact of CRIU on Docker
containers and investigate the feasibility of using CRIU in edge devices.
Specifically, we examine the impact of the size of the memory allocation of
the host system on the time to checkpoint and restore, as well as the
resources consumed during the process.

Our results demonstrate the potential benefits of using CRIU in edge computing
systems, providing insights into the performance impact of using this tool in
these settings. We believe that this research will contribute to the growing
body of work exploring the use of checkpointing and migration in edge
computing systems and will pave the way for further research into the
practical applications of CRIU in this domain.

##  2 Related work

CRIU, which stands for ”Checkpoint/Restore In Userspace,” is an open-source
software tool for performing live migrations of Linux processes. It was first
introduced in 2011 by Pavel Emelyanov as a way to checkpoint and restore
individual processes in Linux user space, with the goal of enabling faster
live migrations for cloud computing workloads. Since then, CRIU has gained
popularity in both cloud and edge computing settings, where it allows for
efficient and flexible management of distributed applications.

Since it’s inception the project has matured considerably. One of the earliest
studies analyzing the performance of checkpoint and restore was published by
Chen Yang [1]. Yang analyzed the time it took to checkpoint and restore a
process and unidentified process memory allocation to be the most impactful
factor. However, Yang’s study is very limiting as only small processes were
used, and no other telemetry such as CPU utilization were analyzed while
performing migrations.

Oh, SeungYong, and JongWon Kim experimented the use of CRIU for live
migrations of containerized applications [4]. In their study, they analyze the
performance of CRIU with a focus on latency analysis in transferring
checkpoints between hosts and were less concerned about the resources and time
needed to perform these tasks.

A similar study in which CRIU was used to improve the availability of docker
swarm, a cluster orchestration solution provided by Docke [2]. Their findings
report that higher availability can be achieved by periodically checkpointing
the states of containers within the swarm and restoring them to provide higher
availability. Similarly to the study of Oh, SeungYong, and JongWon Kim, their
analysis is focused on the memory use and storage on the network file system
(NFS), and did not study the impact on the computing resource utilization
while performing CRIU. Moreover, their study is very limiting with very few
scenarios and container configurations analyzed.

A more recent study by Adityas et. al. analyzed the resource utilization of
the host machine while performing checkpoint and restore [8]. In their study,
they report the expected memory allocation of the container to have the most
impact on the time it takes to perform checkpoint and restore. Interestingly,
their study also reports on the CPU consumption while performing these tasks.
However, instead of extensively testing various containers, the authors chose
a scenario based approach. They identified four three main, namely _One way
migration between two hosts one direction at a time_ , _two way migrations
using one service_ , and _two way migrations using three services_. While the
preliminary results reported do offer some insight, the study is very limited
to the aforementioned scenarios. The study would benefit a more detailed
analysis of multiple services with various memory allocations and a much
larger sample size.

Previous research suggests there is the use of CRIU is ubiquitous in edge and
micro-service architectures. However, the lack of an in-depth performance
analysis that highlights the limitations is evident. Our study aims to address
the gap and pave the way for future use of CRIU for run-time application
migrations on the edge.

##  3 Experimental Results

When CRIU is used to checkpoint a running application, it performs the
following steps:

  1. 1.

Image creation: CRIU creates an image of the process that needs to be
checkpointed. This image contains information about the process’s memory, open
files, network connections, and other relevant information. The image is
created in a checkpoint file.

  2. 2.

Freezing: Once the image is created, CRIU freezes the process. Freezing is the
process of pausing the execution of the process so that its state does not
change while the checkpoint is being created.

  3. 3.

Image dumping: CRIU dumps the image of the process into the checkpoint file.
This process involves saving the state of the process’s memory, CPU registers,
file descriptors, network connections, and other relevant information.

  4. 4.

Thawing: After the image has been dumped, CRIU unfreezes the process and
allows it to resume execution.

The checkpoint file created by CRIU contains all the information necessary to
restore the process’s state at a later time. When the process needs to be
restored, CRIU reads the checkpoint file and recreates the process from the
saved state. This process involves creating a new process, restoring the
process’s memory, CPU registers, file descriptors, network connections, and
other relevant information, and then resuming the process’s execution.

The ability to store and restore the state of a running application is
possible due to the design of the Linux kernel. Linux uses a virtual memory
system to manage memory, which allows processes to have their own virtual
address space. This means that each process has its own view of memory, and
the kernel can manage the physical memory separately for each process. When a
process is checkpointed, its virtual memory state is saved, including the
contents of the process’s memory and its CPU registers. This allows the
process to be restored to the same state later, even if it is running on a
different host or at a different time.

Memory is obviously the main factor that impacts the performance of CRIU. To
get a better understanding relationship between the size of virtual memory
applications require and CRIU performance, a program written in Rust, was
deployed as a docker container, the program allocates a given amount of memory
and runs using the Ubuntu 20.04 as the base image. After the initial
allocation, an external process checkpoints the container and restores it
after while measuring the timings for individual tasks. Figure 3 illustrates
the timings (checkpoint and restore) observed for various memory allocation
sizes.

Fig. 1. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

Fig. 1. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

We can observe that the time it takes to perform a checkpoint is proportional
to the amount of memory the container had allocated. While this is not
surprising, it does set some boundaries on how checkpointing can be used.

Moreover, the analysis shows that CPU utilization does not seem to scale with
the size of virtual memory allocated during the checkpoint and restore phases.
This suggests that CPU utilization may not be a limiting factor in
checkpointing applications, and the trade-off between checkpointing and
performance may be acceptable for stateless applications with limited
resources available on edge devices.

Note that the container has to be paused in order to dump the pages of virtual
memory and maintain consistency. Checkpoints of applications, which require
high availability and constant up-time(i.e. servers) do not seem feasible.
However, micro-services and stateless applications seem reasonable. The
absence of global state in applications greatly reduces the memory used making
CRIU a viable approach for scalability, and fault tolerance.

However, edge computing solutions have to consider limited resources available
on edge devices. While the time it takes to perform a checkpoint and restore
may be a trade-off worth taking for stateless applications, the CPU
utilization on the host device may be a limiting factor.

Fig. 2. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

Fig. 2. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

To gain some insight into the CPU utilization of the host dockerd service and
CRIU we monitor system telemetry when performing checkpoint and restore tasks.
Figure 3 illustrates the CPU utilization in % both of the user running the
checkpoint as well as the overall system. We can observe that while the memory
allocation of the application directly impact the time to perform these tasks,
the CPU utilization does not seem to scale with the size of the virtual memory
allocated. Moreover, the same can be concluded for the restore phase observed
in Figure 3.

Fig. 3. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

Fig. 3. CRIU checkpoint and restore timings relative to the amount of memory
used by the docker container.

Overall, the analysis highlights the potential benefits and limitations of
using CRIU for checkpointing applications and provides insights into the
trade-offs between performance, memory usage, and fault tolerance.

##  4 Conclusion

In conclusion, we have investigated the capabilities and limitations of the
CRIU tool for checkpointing and restoring containers in the context of edge
computing. Our findings demonstrate that CRIU provides a robust solution for
stateless applications and microservices where fault tolerance and scalability
are essential. However, the size of the virtual memory allocated by the
application greatly impacts the time required to perform checkpoint and
restore operations. Thus, edge computing solutions should consider the trade-
off between checkpointing frequency and application performance. Additionally,
while the CPU utilization does not seem to scale with the size of the virtual
memory allocated, it remains a limiting factor in the deployment of
checkpointing on edge devices with limited resources.

Future work can explore techniques to reduce the memory footprint of the
container, enabling more frequent checkpoints without affecting the
application’s performance. Furthermore, it would be interesting to investigate
the use of CRIU for more complex stateful applications and the integration of
CRIU with other container orchestration tools. Overall, our study contributes
to the understanding of the capabilities and limitations of CRIU in the
context of edge computing and provides insights for the deployment of
checkpointing solutions in resource-constrained environments.

## References

  * [1] Yang Chen.  Checkpoint and restore of micro-service in docker containers.  In 2015 3rd International Conference on Mechatronics and Industrial Informatics (ICMII 2015), pages 915–918. Atlantis Press, 2015.
  * [2] Cheng-Hao Huang and Che-Rung Lee.  Enhancing the availability of docker swarm using checkpoint-and-restore.  In 2017 14th International Symposium on Pervasive Systems, Algorithms and Networks & 2017 11th International Conference on Frontier of Computer Science and Technology & 2017 Third International Symposium of Creative Computing (ISPAN-FCST-ISCC), pages 357–362. IEEE, 2017.
  * [3] Dirk Merkel.  Docker: lightweight linux containers for consistent development and deployment.  Linux journal, 2014(239):2, 2014.
  * [4] SeungYong Oh and JongWon Kim.  Stateful container migration employing checkpoint-based restoration for orchestrated container clusters.  In 2018 International Conference on Information and Communication Technology Convergence (ICTC), pages 25–30. IEEE, 2018.
  * [5] André Pires, José Simão, and Luís Veiga.  Distributed and decentralized orchestration of containers on edge clouds.  Journal of Grid Computing, 19(3):1–20, 2021.
  * [6] Aleksandar Tošić, Jernej Vičič, Michael David Burnard, and Michael Mrissa.  A blockchain-based edge computing architecture for the internet of things.  2022\.
  * [7] Ranjan Sarpangala Venkatesh, Till Smejkal, Dejan S. Milojicic, and Ada Gavrilovska.  Fast in-memory criu for docker containers.  In Proceedings of the International Symposium on Memory Systems, MEMSYS ’19, page 53–65, New York, NY, USA, 2019. Association for Computing Machinery.
  * [8] Adityas Widjajarto, Deden Witarsyah Jacob, and Muharman Lubis.  Live migration using checkpoint and restore in userspace (criu): Usage analysis of network, memory and cpu.  Bulletin of Electrical Engineering and Informatics, 10(2):837–847, 2021.