# Development of large-volume 130TeO2 bolometers for the CROSS $2\beta$ decay
search experiment

F.T. Avignone III  A.S. Barabash  V. Berest  L. Bergé  J.M. Calvo-Mozota  P.
Carniti  M. Chapellier  I. Dafinei  F.A. Danevich  L. Dumoulin  F. Ferella  F.
Ferri  A. Gallas  A. Giuliani  C. Gotti  P. Gras  A. Ianni  L. Imbert  H.
Khalife  V.V. Kobychev  S.I. Konovalov  P. Loaiza  P. de Marcillac  S.
Marnieros  C.A. Marrache-Kikuchi  M. Martinez  S. Nisi  C. Nones  E. Olivieri
A. Ortiz de Solórzano  Y. Peinaud  G. Pessina  D.V. Poda  Ph. Rosier  J.A.
Scarpaci  V.I. Tretyak  V.I. Umatov  M.M. Zarytskyy  and A. Zolotarova

###### Abstract

We report on the development of thermal detectors based on large-size
tellurium dioxide crystals (45$\times$45$\times$45 mm), containing tellurium
enriched in 130Te to about 91%, for the CROSS double-beta decay experiment. A
powder used for the crystals growth was additionally purified by the
directional solidification method, resulting in the reduction of the
concentration of impurities by a factor 10, to a few ppm of the total
concentration of residual elements (the main impurity is Fe). The purest part
of the ingot (the first $\sim$200 mm, about 80% of the total length of the
cylindrical part of the ingot) was determined by scanning segregation profiles
of impurities (both on the surface and in center of the ingot) and used for
the 130TeO2 powder production with no evidence of re-contamination. The
crystal growth was verified with precursors produced from powder with natural
Te isotopic composition, and two small-size (20$\times$20$\times$10 mm)
samples were tested at a sea-level laboratory showing high bolometric and
spectrometric performance together with acceptable 210Po content (below 10
mBq/kg). This growth method was then applied for the production of six large
cubic 130TeO2 crystals for the CROSS experiment and 4 of them were taken
randomly to be characterized at the Canfranc underground laboratory, in the
CROSS-dedicated low-background cryogenic facility. Two 130TeO2 samples were
coated with a thin metal film in form of Al layer (on 4 sides) or AlPd grid
(on a single side) to investigate the possibility to tag surface events by
pulse-shape discrimination. Similarly to the small natural precursors, large-
volume 130TeO2 bolometers show high performance and even better internal
purity (210Po activity $\sim$ 1 mBq/kg, while activities of 228Th and 226Ra
are below 0.01 mBq/kg), satisfying requirements for the CROSS and,
potentially, next-generation experiments.

##  1 Introduction

The search and investigation of spontaneous nuclear transitions like double-
beta ($2\beta$) processes, which are expected to occur rarely or even to be
forbidden according to the Standard Model (SM) of particle physics, have
represented a great interest for nuclear and particle physics until recently
[1, 2, 3, 4]. The mostly searched $2\beta$ processes (driven by high
transition energy and a relatively large presence in nature of the isotopes of
interest [5]) are decays in which two electrons appear in the final state
accompanied or not with the emission of two anti-neutrinos, corresponding to
two-neutrino ($2\nu$) and neutrinoless ($0\nu$) $2\beta$ decays respectively.
The former process has been detected for a dozen nuclei with half-lives of
1018–1024 yr [6], while the latter one remains unobserved and the most
stringent half-life limits reach 1024–1026 yr level [2]. The observation of
$0\nu 2\beta$ decay would help to scrutinize many open questions of particle
physics such as the lepton-number conservation, the absolute value of neutrino
mass and its origin (Majorana or Dirac particle), the admixture of right-
handed currents in weak interactions, the existence of majorons and bosonic
neutrinos, the imbalance between matter and anti-matter (via leptogenesis) and
other effects beyond the SM [2].

Among the candidates studied, 130Te stands out as one of the few extensively
researched $2\beta$-active isotopes, thanks to: a) a relatively high
transition energy ($Q_{2\beta}$ = 2527 keV [7]); b) a considerable presence in
natural tellurium ($\delta$ = 34% [8]); c) the possibility to grow large-
volume, high-quality, and radiopure tellurium dioxide (TeO2) crystals, i.e. a
dielectric material with a high tellurium content suitable for the use in low-
temperature detectors of rare events [9, 10, 11, 12]. Thanks to the above-
listed features, TeO2 is among a few materials with the longest history of
applications to $2\beta$ searches. Indeed, initial attempts were made in the
early 1990s using small samples. Subsequently, these efforts were expanded to
larger arrays in experiments such as MiBETA, Cuoricino, and CUORE-0.
Currently, this research has reached its pinnacle with the CUORE ton-scale
experiment [13], which is operating at the Gran Sasso underground laboratory
(LNGS) in Italy. (The evolution of TeO2-based bolometric $2\beta$ search
experiments is overviewed in [14].)

All past and present TeO2-based $2\beta$ experiments exploited and currently
used crystals with natural Te content, with the exception of two samples
enriched in 130Te used in addition to natural ones in the MiBETA [15] and
Cuoricino [16] experiments. Due to poor performance and radiopurity of
enriched samples compared to natural crystals, as well as taking into account
the still high natural isotopic abundance of 130Te, 130Te-enriched crystals
were abandoned from the CUORE program. However, the interest in 130TeO2 has
raised again in view of the CUORE upgrade [17], particularly in the CUPID
program [18, 19]. Two new enriched samples have been produced and tested at
LNGS showing high performance and low radioactivity [20]. Due to challenging
particle identification with TeO2 bolometers [21], mainly exploiting the
detection of $\gamma(\beta)$-induced Cherenkov radiation proposed in [22], and
difficulties in further reduction of $\gamma$ background in the CUORE
cryostat, CUPID has selected lithium molybdate (Li2MoO4) scintillating
crystals [23]. This compound allows for a more efficient scintillation-
assisted particle identification and a lower $\gamma$ background, by moving
the region of interest ($Q_{2\beta}$ of 100Mo is 3034 keV [7]) above the end-
point of natural radioactivity (2615 keV $\gamma$ of 208Tl from the 232Th
family).

Nevertheless, 130Te remains one of a few $2\beta$-active isotopes highly
suitable for multi-isotope bolometric study [24]. Furthermore, TeO2 material
has a few advantages compared with Li2MoO4: a) higher detection efficiency for
the same volume (higher density and effective atomic number); b) no
hygroscopic properties; c) two orders of magnitude longer half-life of the
ordinary $2\nu 2\beta$ decay, which makes the background contribution from
random coincidences of the associated events negligible for 130Te-enriched
thermal detectors, whereas this is the major component of background expected
in CUPID with 100Mo-enriched scintillating bolometers [25]. Taking into
account the possibility to improve the present, already extremely low,
$\gamma$ background in CUORE by even more careful selection of external
screens with lower contamination and/or by applying innovative techniques for
active suppression of backgrounds induced by $\gamma$ (e.g. proposed by BINGO
[26]) and $\alpha$ (e.g. by detection of $\gamma$($\beta$)-induced Cherenkov
radiation from TeO2 [22]), the choice of 130TeO2 for future high-sensitivity
$2\beta$ searches is still attractive and viable.

With this in mind, TeO2 together with Li2MoO4 detector materials have been
chosen for the implementation of the CROSS program on the development of
metal-coated bolometers with the capability of pulse-shape discrimination
between bulk and near surface interactions [27, 28]. The validation of the
CROSS technology for active rejection of surface-originated $\alpha$ and
$\beta$ backgrounds is planned to be realized in a small-scale demonstrator
with tens of Li2100MoO4 and a few 130TeO2 bolometers. For the sake of
completeness, we refer to early attempts in the development of active
background rejection for bolometric $2\beta$ detectors, including TeO2, which
are summarized in [21] and we underline a couple of currently ongoing projects
on this topic. As in CROSS, TeO2 and Li2MoO4 are both included in the BINGO
project, which is developing an innovative detector assembly system with high-
performance bolometric photodetectors and a cryogenic active veto based on BGO
crystals with bolometric photodetection [26]. Recently, the BINGO detector
system with the nylon wire serving to hold a similar size Li2MoO4 crystal
(twice lighter than TeO2) have been tested underground in the CROSS cryostat
and good detector performance has been reported [29]. Also, the idea of
scintillating polymer-based holder allowing active background rejection of
radioactivity induced by close materials of the detector structure has been
proposed. First prototypes with small and large TeO2 bolometers were tested in
such a 3D-printed detector holder [30, 31], as a part of R&D for next-
generation $2\beta$ experiments with this detector material.

Therefore, interest in TeO2 low-temperature detectors is still high for
$2\beta$ search applications. With this aim, we report on the development of
large-volume 130TeO2 crystals for the CROSS experiment. We present then the
first results of bolometric measurements both with newly developed small-size
natural crystals and large-size 130Te-enriched crystals, operated in dilution
refrigerators of above-ground and underground laboratories respectively.

##  2 Development of TeO2 crystals from 130Te-enriched tellurium

###  2.1 Initial 130Te-enriched material

The initial tellurium powder enriched in 130Te isotope to around 93% (provided
by the University of South Carolina) was assayed as received from the customer
for impurities and oxygen content. Using the inductively coupled plasma mass-
spectrometry (ICP-MS), the total concentration of impurities excluding oxygen
were measured at 25 ppm; with the duplicate sample, a consistent result of 27
ppm was obtained. The detected impurities are listed in table 1. The oxygen
content is found on the level of 2%.

Table 1: Impurities detected in samples which were produced during an R&D on purification of the 130Te-enriched material. We quote results for initial, as-purchased, raw material, after purification (average for pure section of the ingot produced by solidification, see text), and after the conversion to a 130TeO2 powder. No results are given for cases when the signal is below the detection limit. Uncertainties are at 20% level. Element | Concentration (ppm) of element in 130Te-enriched samples
---|---
| Raw material | After solidification |  130TeO2 powder
Al | 0.6 |  | 0.4
Ca | 1.0 | 0.3 | 0.3
Cr | 0.6 | 0.2 | 0.5
Cu | 2.8 |  | 0.3
Fe | 8.6 | 1.8 | 1.7
K | 0.4 |  |
Li | 0.1 |  |
Mg | 2.5 | 0.2 |
Mn | 0.1 |  |
Na | 0.4 |  | 0.2
Ni | 1.4 | 0.1 | 0.3
Si | 2.7 |  | 1.6
Sn |  |  | 0.2
Ti | 0.1 |  |
Zn | 3.1 | 0.5 | 0.2
Total | 24.6 | 3.1 | 5.7

###  2.2 Purification by the solidification method

The purification of the enriched material has been performed in Canada (5N
company, Garand, Montreal, Quebec). At the initial step — hydrogenation — a
deoxidized ingot was produced from the tellurium oxide powder. The material
was remelted under a stream of hydrogen to obtain a solid ingot, as
illustrated in figure 1.

Figure 1: View of a 130Te-enriched sample before (top) and after (bottom) the
hydrogenation step. A 1 kg load of the material is used for this demonstration
test. Figure 2: View of the 130Te-enriched tellurium ingot produced by the
directional solidification method. Red lines show the positions of the probes
(in the center and on the surface) used in the study of segregation profiles
of the detected impurities, exploited to define the poorest part of the
present ingot (below 180 mm, see details in text).

Then the directional solidification method was applied. The tellurium ingot
from the hydrogenation step was remelted under a partial atmosphere of
hydrogen and solidified directionally to get impurity segregation and to
obtain a material with a 5N purity. After the process, the ingot was cut at 6
positions over the length, as shown in figure 2, and each position was sampled
in the center and on the surface of the ingot. Each sample from the different
positions in the ingot was assayed to get segregation profiles, showing
signals of impurities above the detection limit at positions of 180–210 mm.
This is illustrated in figure 3 for the main impurity (Fe) of the 130Te-
enriched material. Therefore, the impure section above 180 mm was removed,
while the rest of the ingot was ground and the powder was homogenized for a
final sampling. The final assay of the purified tellurium, done using ICP-MS
measurements, is presented in table 1. The total impurity in the 130Te-
enriched powder is detected at 3 ppm, showing no hint of a re-contamination of
the material to be used for the crystal growth.

Figure 3: Segregation profile of iron in samples taken from the center (blue)
and on the surface (orange) of the 130Te-enriched tellurium ingot produced by
the directional solidification.

For the first run of purification a yield of 68% was achieved. This figure has
been improved to $\sim$80% in the following batch (of 7 kg) with further
improvement of the deoxidation step and a better definition of the impure
section after the directional solidification. A further purification of the
rejected material is considered later in the program.

###  2.3 Conversion to 130TeO2 powder

Tellurium oxidation assay was also performed by the 5N company using the
material purified via the directional solidification and following a standard
process. A 200 g load of tellurium was ground into a powder to carry out the
oxidation test. The obtained TeO2 powder was found visually white (see in
figure 4) and differential scanning calorimetry assay revealed no unoxidized
tellurium (detection limit was 100 ppm). A particles size distribution
analysis was performed on the final powder using laser scattering particle
size analyzer (Horiba, LA-950); the mean particle size is measured as 43
$\mu$m. Impurity levels investigated with ICP-MS are presented in table 1. The
total level of impurities is low and meets the objective of 10 ppm limit. The
main impurities are iron at 1.7 ppm and silicon at 1.6 ppm.

Figure 4: TeO2 powder made of highly purified tellurium enriched in 130Te to
93%.

Thus, all the purified material has been oxidized using the standard process.
The yield is close to 80% and the total purification factor is about 10, with
a residual average impurity concentration of 1.4 ppm in the final purified
130TeO2 powder.

###  2.4 Crystal growth

We identified a company in the US, Gooch & Housego (Clevelend, Ohio), for the
Czochralski growth of large-volume (about 90 cm3) 130Te-enriched TeO2 crystals
to be included in the CROSS demonstrator. In order to validate the crystal
growth process, two natural TeO2 samples (24 g mass each) were produced by
this company from non-purified materials and tested as Cherenkov-radiation-
emitting bolometers in Orsay to demonstrate acceptable bolometric performance
and the level of internal 210Po contamination (less than 10 mBq/kg), as
detailed in the next section. After the validation of the TeO2 crystal
producer, we proceeded with the growth of 130Te-enriched crystals and six
130TeO2 samples 45$\times$45$\times$45 mm each were fabricated using purified
starting materials.

In order to investigate the isotopic composition of Te in the grown 130TeO2
crystals, we took fragments (about 80 mg in total) of crystals coming from
different locations in the ingot, dissolved them in a mixture of HCl and HNO3
acids, and carried out the ICP-MS measurement in isotopic mode. The results on
the Te isotopic concentration in the 130Te-enriched crystals are given in
table 2. The enrichment in the isotope tellurium-130 is found at $\sim$91%,
which is still very high compared to the enrichment of the initial powder
($\sim$93%). Thus, the adoption of the Czochralski technique (i.e. a small
seed) avoids an important dilution of the enriched material in the growth
seed. This result would be impossible using the Bridgman method, e.g. employed
for the production of the natural TeO2 crystals for the CUORE experiment, as
in this case the seed is very large and dominate the final crystal
composition. This result is innovative and very important for future
bolometric $2\beta$ decay search experiments based on enriched tellurium.

Table 2: ICP-MS-measured isotopic composition of tellurium in 130Te-enriched crystals developed in the present work. Isotopic abundance of natural tellurium is given for comparison. Te isotopes | Isotopic abundance (%) in tellurium samples
---|---
| Enriched in 130Te | Natural Te [8]
120Te | 0.00020(4) | 0.0009(1)
122Te | 0.0038(8) | 0.0255(12)
123Te | 0.013(3) | 0.0089(3)
124Te | 0.072(15) | 0.0474(14)
125Te | 0.11(3) | 0.0707(15)
126Te | 0.30(5) | 0.1884(25)
128Te | 8.1(2) | 0.3174(8)
130Te | 91.4(5) | 0.3408(62)

##  3 Operation of TeO2-based low-temperature detectors

###  3.1 Above-ground test of small-size natural precursors

In order to validate the selected producer of 130TeO2 crystals for CROSS, we
carried out a low-temperature characterization of two TeO2 small-size samples
(20$\times$20$\times$10 mm, 24 g each), grown from high-purity powder with
natural Te isotopic composition (without the additional purification foreseen
in CROSS). The detector construction was carried out in the same way as CROSS
prototypes used in studies of metal-coated bolometers [27]. Specifically, each
crystal was equipped with a temperature sensor (NTD-Ge, neutron-transmutation-
doped Ge [32]) and a heater (heavily doped Si chip [33]); both elements were
glued with a bi-component epoxy (Araldite®) on a 20$\times$20-mm side of a
TeO2 crystal. The crystals were mounted on a Cu plate with the help of two Cu
rods and PTFE (polytetrafluoroethylene) spacers. The assembly was secured by
nuts. To provide electrical connections, we wire-bonded sensors and heaters
with Au wires to the gold-plated-on-Kapton contacts. The assembled TeO2
bolometric module is shown in figure 5 (Left). The lateral part of the
assembly was covered by a Cu cylinder; the internal side of this screen and
the Cu plate facing TeO2 crystals were covered by a reflective film
(Vikuiti™).

Figure 5: Bolometers based on two 20$\times$20$\times$10-mm TeO2 crystals
(Left) and a $\oslash$44-mm Ge slab with a 0.2 mm thickness (Right, from
[34]), used in the present work for the detection of Cherenkov radiation
induced by particle interactions in the TeO2 thermal detectors.

On the top of the TeO2 assembly we placed a low-threshold optical bolometer
[35] to attempt to detect the tiny Cherenkov radiation induced by $\beta$ /
$\gamma$ / $\mu$ interactions in the TeO2 bolometers [36], thus allowing for
selection and study of $\alpha$ events. This photodetector (NTLLD2 in [35])
was made of an electronic grade pure Ge slab ($\oslash 44\times 0.18$ mm) and
instrumented with a smaller size NTD-Ge sensor epoxy-glued on the wafer
surface. The Ge wafer was coated by a 70-nm SiO layer to reduce the
reflectivity of the surface and thus to improve photon collection. A set of
aluminum concentric electrodes were deposited on the Ge surface and
interconnected with Al wires, allowing to apply a voltage bias and to amplify
thermal signals exploiting the Neganov-Trofimov-Luke effect in semiconductors.
The Ge slab was kept in a Cu housing with the help of sapphire balls
($\varnothing 1.5$ mm) and polypropylene supporting elements; the NTD-Ge
sensor was polarized through Au bonding wires. The mounted optical bolometer
is depicted in figure 5 (Right). The photodetector was facing the TeO2
crystals by the side carrying the Al electrodes, while the opposite side was
facing another bolometer (ZnO). Its recently published characterization [34]
contains more details on the experiment and results of the light detector
performance, which are only briefly described below.

The fully assembled device was tested with the help of a pulse-tube-based
dilution refrigerator [37] placed at IJCLab (Orsay, France) inside a 10-cm-
thick lead shielding. The detector module was spring-suspended inside the
cryostat, having a flexible thermal link to the coldest part of the unit to
mitigate vibrational noise induced by the pulse tube and pumps [38]. The
detectors were cooled down to 15 mK. The NTD-Ge sensors were polarized with
$O$(1 nA) currents to get working resistances below 10 M$\Omega$, in
compliance with the electronics used (e.g. the set of the load resistances)
and to avoid a large contribution of the Johnson noise. The continuous streams
of data were acquired by 16-bit National Instruments cards (NI USB-6218 BNC)
with a sampling rate of 5 kS/s; an anti-aliasing Bessel filter with a cut-off
frequency of 675 Hz was used in the readout. The data processing was realized
using an optimum filter technique [39] implemented in a MATLAB-based
application [40] with data-driven templates of signal and noise of the
detectors used to construct the transfer function of the filter.

Figure 6: A partial energy spectrum measured by a 24-g TeO2 bolometer in a
low-temperature test (15 mK; 96 h of data taking) in an above-ground
laboratory. The most intense $\gamma$ peaks, corresponding to 214Pb and 214Bi,
are labeled. A 609-keV peak of 214Bi together with a fit to a Gaussian plus a
flat background component are shown in the Inset. The energy resolution of the
peak is measured as 1.6(4) keV FWHM.

Measurements were carried out with a 232Th calibration source for 24 h and
without the use of additional radioactive sources, the nominally background
data were taken for a total of 166 hours. Due to a comparatively high
environmental $\gamma$-radioactivity (mainly from decays of the 226Ra sub-
chain [41]), and thanks to high spectroscopic properties of the TeO2
bolometers, background data can be easily calibrated because of the presence
of many sharp peaks in a wide energy range as seen in figure 6. We observed
that both TeO2 thermal detectors show excellent energy resolution, comparable
to results achieved with bare (i.e. non-coated) samples of the same-size used
in the CROSS program [27, 42]. The detectors are characterized by high
sensitivity defined as the voltage signal per unit of deposited energy
(0.2–0.5 $\mu$V/keV) and a low noise represented by the baseline width at the
optimum filter output (0.25–0.60 eV RMS). As for the time properties of the
thermal response, both TeO2 bolometers show typically slow signals with the
ascending part of about 12 ms and with the descending part about 60 ms. The
performances of the TeO2 thermal detectors are summarized in table 3.

Table 3: Performance of two 24-g TeO2 bolometers and a Ge light detector exploiting the signal amplification based on the Neganov-Trofimov-Luke effect (60 V bias on Al electrode). We report the working point of the NTD thermistor (resistance $R_{NTD}$ at a given current $I_{NTD}$), pulse shape time constants (rise- and decay-time parameters, $\tau_{r}$ and $\tau_{d}$), detector sensitivity represented by the signal voltage amplitude per unit of the deposited energy ($A_{s}$), RMS baseline width (RMSnoise). Detector | $R_{NTD}$ | $I_{NTD}$ | $\tau_{r}$ | $\tau_{d}$ | $A_{s}$ | RMSnoise
---|---|---|---|---|---|---
ID | (M$\Omega$) | (nA) | (ms) | (ms) | ($\mu$V/keV) | (keV)
TeO-63 | 8.8 | 1.3 | 12.7 | 67 | 0.47 | 0.59
TeO-81 | 5.0 | 1.5 | 12.4 | 61 | 0.22 | 0.25
LD-NTL [34] | 4.2 | 1.3 | 1.3 | 11 | 38 | 0.017
Figure 7: Light signals detected by the Ge optical bolometer versus heat
signals measured by the 24-g TeO2 bolometer (TeO-81) over 96 h of exposure to
natural radioactivity at a sea-level laboratory. The major part of events
below $\sim$2.6 MeV corresponds to $\gamma$($\beta$) interactions in the
crystal medium, while muon-induced events represent dominant contribution
above this energy. A structure of high-energy TeO2-detected events which
exhibit no light signal correspond to $\alpha$ decays of 210Po inside the TeO2
detector medium, indicating an internal contamination by this radionuclide.

During the TeO2 detector calibrations we kept the Al electrodes of the Ge
bolometer unpolarized, thus allowing the detection of Te X-rays induced by
environmental radioactivity in the crystals, which was used to calibrate the
maximum distribution of the muon-induced events [34]. Knowing the energy scale
of the Ge bolometer, we found a factor $\sim$10 improvement of the signal-to-
noise ratio by applying a 60 V bias on the Al electrode, resulting to a high
performance reported in table 3 and in [34]. Thanks to the low noise achieved
with the Ge bolometer, we were able to detect the TeO2-emitted Cherenkov
light, providing highly efficient particle identification as illustrated in
figure 7. As one can see in figure 7, most of the events populate a band
including $\beta$-, $\gamma$\- and $\mu$-induced interactions, which are
associated with a registered Cherenkov light of the order of 40 eV per MeV
energy deposited in the TeO2 bolometer. The population of $\alpha$ events,
which do not produce the Cherenkov radiation (due to energy threshold of 400
MeV), is also observed in the data (see figure 7). Moreover, these
$\alpha$-events have even negative values of light signals due to a small
cross-talk, also seen in the data for heater-induced events (not shown in
figure 7). The observed $\alpha$ events are mainly ascribed to a 210Po
contamination in the crystal bulk corresponding to an activity of 2.8(4)
mBq/kg and 9.1(8) mBq/kg in the samples TeO-63 and TeO-81 respectively. This
contamination is typical for TeO2, as well as for different crystal
scintillators [43], but does not originate by the 210Pb content as for many
other materials. It is worth noting that we detected only tens of 210Po
events, while the presence of other $\alpha$-active radionuclides from U/Th
chains is not evident, and thus their activities are expected to be at least
an order of magnitude lower. Therefore, one can conclude that the radiopurity
of both samples is very high, especially keeping in mind that the starting
material did not pass an extra purification.

Accordingly, the results of the low-temperature test presented above prove
high bolometric and spectrometric properties together with a low level of
radioactive contamination of the TeO2 crystals, and validate the crystal
producer for the development of the 130Te-enriched crystals for the CROSS
experiment.

###  3.2 Underground test of large-size 130Te-enriched crystals

To test the newly developed 130Te-enriched crystals as bolometers, we took
four of the six samples with a size of 45$\times$45$\times$45 mm each,
produced for the CROSS project. Two chosen samples (TeO-20 and TeO-24) were
used uncoated, while the other two crystals were metal-coated as follows: i)
four-side coated with Al layer (TeO-Al); ii) one-side coated with an AlPd bi-
layer in the form of grid (TeO-AlPd). Each crystal was equipped with an NTD-Ge
thermistor and a Si:P heater, coupled with UV-curable (PERMABOND 620) and
epoxy (Araldite® Rapid) glue respectively. The samples were assembled in the
Cu structure with 3D-printed PLA (polylactic acid) spacers — the holder _Slim_
version designed for the CROSS experiment [44]— fasting two crystals per
flour. In addition, a thin Ge bolometer (45$\times$45$\times$0.3 mm) has been
mounted on each TeO2 crystal to detect Cherenkov radiation. All Ge slabs were
SiO-coated to reduce the reflectivity of the surface. A set of Al concentric
electrodes, not optimized for the adopted squared geometry of the light
detectors, was deposited following the procedure described in [35]. The
purpose of the electrodes is to amplify the thermal signal via the Neganov-
Trofimov-Luke effect. In total, the tower was composed of 5 floors; the rest
of the array was made of six lithium molybdate based bolometers and the same
type of Ge light detectors attached to them [45]. The part of the tower
containing the 130TeO2 bolometers and construction elements of the detector
structure are depicted in figure 8.

Figure 8: Fragment of a 5-floor array of dual-readout bolometers showing the
last two floors with large-volume 130TeO2 bolometers (Left). Two 130TeO2
crystals are coated with Al film on 4 sides and AlPd grid on a single side
respectively (Right top). Thin square-shaped Ge bolometers with concentric Al
electrodes (light detectors with Neganov-Trofimov-Luke effect signal
amplification) are placed above the 130TeO2 thermal detectors (Right bottom).

The detector array was installed inside a low-background facility with a
pulse-tube-based dilution refrigerator, built in the Canfranc underground
laboratory (Spain) for the CROSS experiment [46, 47]. In order to mitigate the
vibrational noise induced by the pulse tube cryocooler, and to gain space of
the available experimental volume, the set-up was recently upgraded with a
pendulum-like detector suspension with three magnetic dampers installed at the
1 K stage [48]. Thanks to the underground location of the facility and to
massive external and internal shields, the background level from natural
radioactivity is well reduced providing suitable conditions for
characterisation of large-volume low-temperature detectors [48, 49]. After
reaching a base temperature of around 12 mK, we carried out the
characterization of the bolometers at warmer temperatures (17, 22, and 27 mK)
regulated on the detector plate. The detector was controlled and read-out
using room-temperature electronics [50, 51]. The amplified voltage output of
the NTD-Ge thermistors was digitized in continuous streams by a 24-bit ADC
with a sampling rate of 2 or 10 kHz; the cut-off frequency of the active
Bessel filter was set as 0.3 or 2.5 kHz respectively. The optimum-filter-based
data processing was done with the same tool, and in a similar way, as
described above (section 3.1).

At low temperature, we lost no contacts but we faced an issue with the readout
of a single 130TeO2 channel (TeO-24), therefore we grounded it during the
experiment. A similar issue happened for another channel (TeO-Al) after the
electronics upgrade 111Currently, 11 over 14 front-end boxes and twelve
12-channel DAQ boards have been installed., thus we were not able to study
this detector at the warmest temperature (27 mK).

After the characterization of the NTD-Ge sensors at the coldest detector plate
temperature chosen (17 mK), we found that, for the optimal working point in
terms of highest detector sensitivity, the measured resistances of the
thermistors are too high — $O$(100 M$\Omega$) — for the optimal operation of
the existing electronics. This is due to the significantly higher stress of
the sensor coupled with UV-cured glue with respect to that observed for
lithium molybdates [44]. Therefore, we characterized 130TeO2 bolometers also
at higher detector plate temperatures and polarized NTD-Ge with stronger
currents reducing the resistances to a few-tens M$\Omega$. The detectors show
a rather high sensitivity at 17 mK, $O$(100 nV/keV), which is reduced by a
factor 10 at subsequent operation at warmer temperatures. It is worth
emphasizing a comparatively small reduction (a factor of 2) in sensitivity of
the metal-coated crystals (irrespective of single or 4-side coating with
either a grid or a layer) with respect to the bare sample, in contrast to
observations with Li2MoO4 crystals [42]. A drastic decrease in the sensitivity
at warmer working temperatures has a minor impact on the measured baseline
noise after optimum filtering, which is around 2–5 keV FWHM, similar to what
was reported for Li2MoO4-based bolometers tested in this set-up in the same /
similar detector structure [44]. It is worth noting the quite long thermal
signals of the 130TeO2 bolometers (figure 9) exhibiting a few seconds long
descending part of the signals, similar to early reported response for such
detector material, which are more than a factor 5 longer than the signals of
the reference Li2100MoO4 bolometers [44]. At the same time, the internal $2\nu
2\beta$-decay activity of 130Te in the enriched crystals is about two orders
of magnitude lower than in the case of 100Mo (due to the difference in the
half-lives, e.g. see recent results in [52] and [53] respectively), thus it
does not represent a notable background in form of pile-ups in contrast with
the behaviour of 100Mo-enriched bolometers [25]. The performance of 130TeO2
thermal detectors are summarized in table 4.

Figure 9: Fragments of average signals of 130TeO2 bolometers operated at a detector plate temperature of 17 mK (top left), 22 mK (top right), and 27 mK (bottom). Table 4: Performance of 130Te-enriched TeO2 bolometers (0.55 kg each) at different temperatures of the detector plate ($T_{plate}$). For each working point of the NTD thermistor (resistance $R_{NTD}$ at a given current $I_{NTD}$), we report the pulse-shape time constants (rise- and decay-time parameters, $\tau_{r}$ and $\tau_{d}$), the detector sensitivity represented by the signal voltage amplitude per unit of the deposited energy ($A_{s}$), and the RMS baseline width (RMSnoise). Channel | $T_{plate}$ | $R_{NTD}$ | $I_{NTD}$ | $\tau_{r}$ | $\tau_{d}$ | $A_{s}$ | RMSnoise
---|---|---|---|---|---|---|---
| (mK) | (M$\Omega$) | (nA) | (ms) | (ms) | (nV/keV) | (keV)
TeO-20 | 17 | 180 | 0.10 | 45 | 370 | 250 | 3.0
| 22 | 26 | 0.55 | 7.0 | 360 | 38 | 2.7
| 27 | 7.0 | 3.0 | 6.0 | 400 | 9 | 2.4
TeO-Al | 17 | 180 | 0.10 | 37 | 420 | 137 | 2.9
| 22 | 25 | 0.55 | 8.0 | 420 | 18 | 3.6
TeO-AlPd | 17 | 170 | 0.10 | 36 | 600 | 101 | 4.7
| 22 | 6.4 | 1.0 | 3.4 | 24 | 12 | 8.7
| 27 | 3.0 | 3.0 | 3.0 | 550 | 7 | 4.6

Detectors calibration and characterization is done with the help of a 232Th
source, inserted inside the lead shielding. Coincidences between 130TeO2
detectors and corresponding Ge bolometers have been established according to
time properties (rise time) of their thermal signals. An example of a two-
dimensional distribution of the amplitudes of the detected events is
illustrated in figure 10. A band of events with a light signal amplitude
clearly depending on the heat energy — appreciable in figure 10 — corresponds
to $\gamma$($\beta$) interactions in the 130TeO2 absorber resulting to the
emission of Cherenkov radiation. Clusters of high-energy particles with no
detectable light signals correspond to events induced by $\alpha$ decays; a
possible tiny scintillation of the TeO2 material is not observed for
$\alpha$-particles as the conditions of light collection are not favorable due
to the absence of a reflective film around the crystal. This is also the main
reason for only partial discrimination between $\alpha$ and $\gamma$($\beta$)
populations. We notice in fact that a much better discrimination was obtained
with a device consisting of an even larger (50-mm cube) crystal placed inside
a holder with a reflective foil [36], characterized by a similar-performance
but smaller light detector. The discrimination will be improved by using an
optimized geometry of the Al electrode with a better Ge-wafer surface
coverage, which will allow us to improve the amplification by a factor 2 [45].
Moreover, a further enhancement in the NTL gain can be achieved by facing
grids to the crystal surface and by applying a twice higher electrode bias
[45].

Figure 10: Distribution of 130TeO2-detected heat signals (TeO-AlPd) in
coincidences with photon signals of the accompanied Ge bolometer (placed above
the crystal) measured over 84-h 232Th calibration in the CROSS set-up at the
LSC (Spain). The light detector was operational in the NTL amplification mode
(50 V bias was applied on the Al electrode) allowing to reduce the baseline
noise to ten(s) eV RMS. The population of $\gamma$($\beta$) events is mainly
induced by $\gamma$s from the external 232Th source, while $\alpha$ particles
detected are originated to either external U $\alpha$ source used (more broad
distributions) or to internal contamination of the crystal (mainly by 210Po).
Figure 11: Energy spectra of events detected by the 0.55-kg 130TeO2 bolometers
TeO-20 (a), TeO-AlPd (b), and TeO-Al (c). Both metal-coated bolometers were
additionally irradiated by a close $\alpha$ source of 238U/234U. The data were
acquired in the CROSS set-up at the Canfranc underground laboratory in the
following conditions: (a) and (b) correspond to measurements at 27 mK over 116
h (reference date is Nov. 2023); (c) represents 168-h-long dataset collected
at 17 mK (reference date is Mar. 2023). The origin of the most intense
$\gamma$ and $\alpha$ peaks observed in the spectra is labeled.

The 130TeO2 bolometers were operated at a higher temperature — regulated at 27
mK — where the performance is expected to be better considering the high NTD
resistances. In these conditions, we performed a relatively long measurement
(116 h) without the external 232Th $\gamma$ source, attempting to probe the
radiopurity level of the samples. In the case of the TeO-Al sample, we used
such type of “background” data collected at 17 mK with a similar duration (168
h), realized 7.5 months before the measurement at the warmest working point.
The energy spectra of the three investigated detectors showing the region of
interest for $\gamma$, $\beta$, and $\alpha$ spectroscopy are depicted in
figure 11. The experiment was not supplied by the air flow with low radon
content, thus the presented data exhibit a significant $\gamma$ background
induced by radioactive daughters (mainly 214Bi) of radon decay, as
demonstrated in [49]. In the $\alpha$ region, we observed only a single peak-
like structure corresponding to the internal 210Po contamination, while a
nearby peak, shifted by 100 keV to lower energies is traced to surface
pollution of detector materials by 210Po (210Pb). In addition, as already
mentioned above, data of TeO-AlPd and TeO-Al detectors exhibit two broad peaks
from an external $\alpha$ source used to irradiate the crystal. It is
interesting to note that the energy of the “bulk” 210Po peak, calibrated by
2614.5 keV $\gamma$ quanta, is found about 3% higher than the nominal value
(5407.5 keV) for all data and crystals, except 17 mK measurements with the
TeO-20 detector showing a $\sim$1% deficit in the estimated energy. The
observed miscalibration of the 210Po $\alpha$ peak is similar to early
reported data of TeO2 bolometers (e.g. a quenching factor of $\sim$1.02 [54,
36] and $\sim$0.99 [55]) and it can be explained by non-linearity of the
bolometric response, as demonstrated in [55], where the quenching factor
1.0076(5) was measured for the 2.3 MeV $\alpha$ peak of 147Sm. This comparison
additionally shows a strong impact of measurements and data analysis on the
estimate of quenching factors for alpha particles and nuclear recoils,
similarly to such studies with other detector materials [56, 57].

The activity of 210Po detected in the 130Te-enriched crystals, as shown in
figure 11, is 1.6(1), 0.8(1), and 2.5(1) mBq/kg for the sample TeO-20, TeO-
AlPd, and TeO-Al respectively. This is a conservative estimate of the 210Po
bulk contamination, because it includes a contribution of the crystal surface
radioactivity, while a comprehensive background model is needed to scrutinize
the bulk and surface components, e.g. as done in [58]. For the same reference
date (Nov. 2023) of the TeO-20 and TeO-AlPd samples, corresponding to the end
of the experiment described here, the 210Po activity in the TeO-Al crystal is
expected at 0.8 mBq/kg. At the beginning of the experiment (Mar. 2023), the
210Po activity in TeO-20 and TeO-AlPd, estimated from the same 168-h-long
dataset, is 4.5(1) and 2.4(1) mBq/kg respectively, in agreement with the
expectations based on 27 mK data (Nov. 2023). Thus, the 210Po content is a few
times lower than that in the precursor samples with natural Te isotopic
abundance. This confirms the efficiency of the additional purification applied
in the process of the 130TeO2 crystal development. Taking into account that
the 228Th (5520 keV) and 220Rn (5590 keV) peaks are below the detection limit,
we can conclude that the activity of 228Th and 226Ra in the crystal bulk is
expected to be at least two orders of magnitude lower than the observed 210Po
content. Altogether, this proves a high internal radiopurity of the 130TeO2
crystals developed for the CROSS experiment, which can also comply demands of
next-generation experiments.

##  4 Conclusions

In the framework of the CROSS project, we developed a production chain of
tellurium dioxide crystals from highly purified 130Te-enriched powder (with
93% of 130Te in the isotopic composition of Te). The directional
solidification method has been used for the additional purification of the
powder, reducing the concentration of impurities by a factor 10, down to a few
ppm level in total. The purest part of the ingot is determined by the
segregation profiles of impurities (surface / bulk of the ingot); the yield is
around 80%. The capability of the selected crystal-growth company (using
Czochralski method) to supply high quality TeO2 samples with the required
internal radiopurity (e.g. below 10 mBq/kg activity of 210Po) has been
validated by the bolometric results obtained above ground with two small-size
samples (20$\times$20$\times$10 mm), produced from a powder with natural Te
isotopic abundance and no additional purification. After that, we proceeded
with the fabrication of six large 130TeO2 crystals with the size of
45$\times$45$\times$45 mm and a mass of 0.55 kg each. The isotopic abundance
of 130Te remains very high in the crystals ($\sim$91%) thanks to the adopted
Czochralski technique (with a small seed), avoiding a significant dilution of
the 130Te-enriched material in the growth seed.

Four over six 130TeO2 samples have been operated as thermal detectors in the
CROSS low-background set-up at the Canfranc underground laboratory (Spain),
using a newly developed detector structure for the CROSS experiment. Half of
the samples tested were coated with either Al film (on 4 sides) or with AlPd
grid (on 1 side). The detectors show high bolometric and spectrometric
performance, as well as a low internal contamination by $\alpha$-active
radionuclides from U/Th chains. The only detected contamination is ascribed to
210Po with an activity of around 1 mBq/kg (reference date is Nov. 2023), which
is a few times lower than in precursors with natural Te composition and no
additional purification. Activities of the most harmful radionuclides for
double-beta decay search experiments, 228Th and 226Ra from 232Th/238U
families, are expected to be at least two orders of magnitude lower than that
of 210Po (i.e. below $\sim$10 $\mu$Bq/kg). A high radiopurity level of newly
developed large-volume 130TeO2 crystals fully satisfies the CROSS demands and,
potentially, requirements of next-generation experiments.

## Acknowledgments

This work is supported by the European Commission (Project CROSS, Grant No.
ERC-2016-ADG, ID 742345), by the Agence Nationale de la Recherche (ANR France;
Project CUPID-1, ANR-21-CE31-0014), and US National Science Foundation (NS
1614611). We acknowledge also the support of the P2IO LabEx (ANR-10-LABX0038)
in the framework “Investissements d’Avenir” (ANR-11-IDEX-0003-01 – Project
“BSM-nu”) managed by ANR, France. This work was also supported by the National
Research Foundation of Ukraine under Grant No. 2020.02/0011 and by the
National Academy of Sciences of Ukraine in the framework of the project
“Development of bolometric experiments for the search for double beta decay”,
the grant number 0121U111684. Russian and Ukrainian scientists have given and
give crucial contributions to CROSS. For this reason, the CROSS collaboration
is particularly sensitive to the current situation in Ukraine. The position of
the collaboration leadership on this matter, approved by majority, is
expressed at https://a2c.ijclab.in2p3.fr/en/a2c-home-en/assd-home-en/assd-
cross/.

## References

  * [1] J.J. Gómez-Cadenas, J. Martín-Albo, J. Menéndez, M. Mezzetto, F. Monrabal and M. Sorel, _The search for neutrinoless double-beta decay_ , _La Rivista del Nuovo Cimento_ 46 (2023) 619.
  * [2] M. Agostini, G. Benato, J.A. Detwiler, J. Menéndez and F. Vissani, _Toward the discovery of matter creation with neutrinoless double-beta decay_ , _Rev. Mod. Phys._ 95 (2023) 025002 [2202.01787].
  * [3] R.L. Workman et al., _The Review of Particle Physics_ , _Prog. Theor. Exp. Phys._ 2022 (2022) 083C01.
  * [4] A. Giuliani, J.J. Gomez Cadenas, S. Pascoli, E. Previtali, R. Saakyan, K. Schaeffner et al., _Double beta decay APPEC committee report_ , 1910.04688v2.
  * [5] V.I. Tretyak and Y.G. Zdesenko, _Tables of double beta decay data — an update_ , _At. Data Nucl. Data Tables_ 80 (2002) 83.
  * [6] A.S. Barabash, _Precise Half-Life Values for Two-Neutrino Double- $\beta$ Decay: 2020 Review_, _Universe_ 6 (2020) 159.
  * [7] M. Wang, W.J. Huang, F.G. Kondev, G. Audi and S. Naim, _The AME2020 atomic mass evaluation (II). Tables, graphs and references_ , _Chinese Phys. C_ 45 (2021) 030003.
  * [8] J. Meija et al., _Isotopic compositions of the elements 2013 (IUPAC Technical Report)_ , _Pure Appl. Chem._ 88 (2016) 293.
  * [9] Y. Chu, Y. Li, Z. Ge, G. Wu and H. Wang, _Growth of the high quality and large size paratellurite single crystals_ , _J. Cryst. Growth_ 295 (2006) 158.
  * [10] M. Barucci, C. Brofferio, A. Giuliani, E. Gottardi, I. Peroni and G. Ventura, _Measurement of low temperature specific heat of crystalline TeO 2 for the optimization of bolometric detectors_, _J. Low Temp. Phys._ 123 (2001) 303.
  * [11] CUORE collaboration, _CUORE crystal validation runs: Results on radioactive contamination and extrapolation to CUORE background_ , _Astropart. Phys._ 35 (2012) 839.
  * [12] CUORE collaboration, _Study of rare nuclear processes with CUORE_ , _Int. J. Mod. Phys. A_ 33 (2018) 1843002.
  * [13] CUORE collaboration, _CUORE opens the door to tonne-scale cryogenics experiments_ , _Prog. Part. Nucl. Phys._ 122 (2022) 103902.
  * [14] C. Brofferio and S. Dell’Oro, _The saga of neutrinoless double beta decay search with TeO 2 thermal detectors_, _Rev. Sci. Instrum._ 89 (2018) 121502.
  * [15] S. Pirro et al., _Present status of MI-BETA cryogenic experiment and preliminary results for CUORICINO_ , _Nucl. Instrum. Meth. A_ 444 (2000) 71.
  * [16] A. Bryant, _A Search for Neutrinoless Double Beta Decay of 130Te_, Ph.D. thesis, University of California, Berkeley, U.S., Spring, 2010.
  * [17] CUORE collaboration, _Exploring the Neutrinoless Double Beta Decay in the Inverted Neutrino Hierarchy with Bolometric Detectors_ , _Eur. Phys. J. C_ 74 (2014) 3096.
  * [18] CUPID collaboration, _CUPID: CUORE (Cryogenic Underground Observatory for Rare Events) Upgrade with Particle IDentification_ , _missing_ (2015) [1504.03599].
  * [19] CUPID collaboration, _R &D towards CUPID (CUORE Upgrade with Particle IDentification)_, _missing_ (2015) [1504.03612].
  * [20] D.R. Artusa et al., _Enriched TeO 2 bolometers with active particle discrimination: towards the CUPID experiment_, _Phys. Lett. B_ 767 (2017) 321.
  * [21] D. Poda and A. Giuliani, _Low background techniques in bolometers for double-beta decay search_ , _Int. J. Mod. Phys. A_ 32 (2017) 1743012.
  * [22] T. Tabarelli de Fatis, _Cerenkov emission as a positive tag of double beta decays in bolometric experiments_ , _Eur. Phys. J. C_ 65 (2010) 359.
  * [23] CUPID collaboration, _CUPID pre-CDR_ , _missing_ (2019) [1907.09376].
  * [24] A. Giuliani, F.A. Danevich and V.I. Tretyak, _A multi-isotope $0\nu 2\beta$ bolometric experiment_, _Eur. Phys. J. C_ 78 (2018) 272.
  * [25] D.M. Chernyak et al., _Random coincidence of $2\nu 2\beta$ decay events as a background source in bolometric $0\nu 2\beta$ decay experiments_, _Eur. Phys. J. C_ 72 (2012) 1989.
  * [26] A. Armatol, _Innovative methods for background rejection in next-generation neutrinoless double beta decay bolometric experiments_ , Ph.D. thesis, Université Paris-Saclay, Orsay, France, September, 2023.
  * [27] CROSS collaboration, _The 0 $\nu$2$\beta$-decay CROSS experiment: preliminary results and prospects_, _JHEP_ 01 (2020) 018.
  * [28] I.C. Bandac et al., _Phonon-mediated crystal detectors with metallic film coating capable of rejecting $\alpha$ and $\beta$ events induced by surface radioactivity_, _Appl. Phys. Lett._ 118 (2021) 184105.
  * [29] A. Armatol et al., _Double beta decay APPEC committee report_ , _Nucl. Instrum. Meth. A_ (2024) [2402.12262].
  * [30] M. Biassoni et al., _An acrylic assembly for low-temperature detectors_ , _Eur. Phys. J. Plus_ 136 (2021) 986.
  * [31] M. Biassoni et al., _Towards the development of a polymer-based assembly for cryogenic detectors for neutrino-less double beta decay_ , _Eur. Phys. J. Plus_ 138 (2023) 384.
  * [32] E.E. Haller, _Advanced far-infrared detectors_ , _Infrared Phys. Techn._ 35 (1994) 127.
  * [33] E. Andreotti et al., _Production, characterization and selection of the heating elements for the response stabilization of the CUORE bolometers_ , _Nucl. Instrum. Meth. A_ 664 (2012) 161.
  * [34] A. Armatol et al., _ZnO-based scintillating bolometers: New prospects to study double beta decay of 64Zn_, _JINST_ 18 (2023) P06026.
  * [35] V. Novati et al., _Charge-to-heat transducers exploiting the Neganov-Trofimov-Luke effect for light detection in rare-event searches_ , _Nucl. Instrum. Meth. A_ 940 (2019) 320.
  * [36] L. Bergé et al., _Complete event-by-event $\alpha$/$\gamma$($\beta$) separation in a full-size TeO2 CUORE bolometer by Neganov-Luke-magnified light detection_, _Phys. Rev. C_ 97 (2018) 032501.
  * [37] M. Mancuso et al., _An experimental study of antireflective coatings in Ge light detectors for scintillating bolometers_ , _EPJ Web Conf._ 65 (2014) 04003.
  * [38] E. Olivieri et al., _Vibrations on pulse tube based Dry Dilution Refrigerators for low noise measurements_ , _Nucl. Instrum. Meth. A_ 858 (2017) 73.
  * [39] E. Gatti and P. Manfredi, _Processing the signals from solid-state detectors in elementary-particle physics_ , _Riv. Nuovo Cim._ 9 (1986) 1.
  * [40] M. Mancuso, _Development and optimization of scintillating bolometers and innovative light detectors for the search for neutrinoless double beta decay_ , Ph.D. thesis, Université Paris-Sud, Orsay, France, February, 2016.
  * [41] T.B. Bekker et al., _Aboveground test of an advanced Li 2MoO4 scintillating bolometer to search for neutrinoless double beta decay of 100Mo_, _Astropart. Phys._ 72 (2016) 38.
  * [42] H. Khalife, _CROSS and CUPID-Mo: future strategies and new results in bolometric search for $0\nu\beta\beta$_, Ph.D. thesis, Université Paris-Saclay, Orsay, France, January, 2021.
  * [43] F.A. Danevich and V.I. Tretyak, _Radioactive contamination of scintillators_ , _Int. J. Mod. Phys. A_ 33 (2018) 1843007.
  * [44] D. Auguste et al., _A novel mechanical design of a bolometric array for the CROSS double-beta decay experiment_ , _submitted to JINST_ (2024) [2405.18980].
  * [45] D. Auguste et al., _Light detectors with thermal signal amplification for CROSS and CUPID $\beta\beta$ decay search experiments_, _to be submitted to JINST_ (2024) .
  * [46] CROSS collaboration, _The new CROSS Cryogenic Underground (C2U) facility: an overview_ , in _XXIX International (online) Conference on Neutrino Physics and Astrophysics (Neutrino 2020), June 22 – July 02, 2020_ , 2020, https://indico.fnal.gov/event/19348/contributions/186315/.
  * [47] A. Armatol et al., _A CUPID Li 2100MoO4 scintillating bolometer tested in the CROSS underground facility_, _JINST_ 16 (2021) P02037.
  * [48] A. Ahmine et al., _Test of 116CdWO4 and Li2MoO4 scintillating bolometers in the CROSS underground facility with upgraded detector suspension_, _JINST_ 18 (2023) P12004.
  * [49] I. Bandac et al., _Li 2100deplMoO4 scintillating bolometers for rare-event search experiments_, _Sensors_ 23 (2023) 5465.
  * [50] P. Carniti, C. Gotti and G. Pessina, _High-Resolution Digitization System for the CROSS Experiment_ , _J. Low Temp. Phys._ 199 (2020) 833.
  * [51] P. Carniti, C. Gotti and G. Pessina, _High resolution filtering and digitization system for cryogenic bolometric detectors_ , _Nucl. Instrum. Meth. A_ 1045 (2023) 167658.
  * [52] CUORE collaboration, _Measurement of the $2\nu\beta\beta$ Decay Half-Life of 130Te with CUORE_, _Phys. Rev. Lett._ 126 (2021) 171801.
  * [53] C. Augier et al., _Measurement of the $2\nu\beta\beta$ Decay Rate and Spectral Shape of 100Mo from the CUPID-Mo Experiment_, _Phys. Rev. Lett._ 131 (2023) 162501.
  * [54] A. Alessandrello et al., _The thermal detection efficiency for recoils induced by low energy nuclear reactions, neutrinos or weakly intercating massive particles_ , _Phys. Lett. B_ 408 (1997) 465.
  * [55] F. Bellini et al., _Response of a TeO 2 bolometer to $\alpha$ particles_, _JINST_ 5 (2010) P12005.
  * [56] V.I. Tretyak, _Semi-empirical calculation of quenching factors for ions in scintillators_ , _Astropart. Phys._ 33 (2010) 40.
  * [57] V.I. Tretyak, _Semi-empirical calculation of quenching factors for scintillators: new results_ , _EPJ Web Conf._ 65 (2014) 02002.
  * [58] CUORE collaboration, _Measurement of the two-neutrino double-beta decay half-life of 130Te with the CUORE-0 experiment_, _Eur. Phys. J. C_ 77 (2017) 13.


# Skyrmion inertia in synthetic antiferromagnets

Sujit Panigrahy  Université Paris-Saclay, CNRS, Laboratoire de Physique des
Solides, 91405, Orsay, France  Sougata Mallick  Université Paris-Saclay, CNRS,
Laboratoire de Physique des Solides, 91405, Orsay, France  João Sampaio
Université Paris-Saclay, CNRS, Laboratoire de Physique des Solides, 91405,
Orsay, France  Stanislas Rohart<EMAIL_ADDRESS>Université Paris-Saclay, CNRS, Laboratoire de Physique des Solides, 91405,
Orsay, France

###### Abstract

We describe the dynamics of magnetic skyrmions in synthetic antiferromagnets
(SAF), with a finite interlayer coupling. Due to the opposite gyrovector of
the skyrmions in the two SAF layers, the coupled skyrmions reach a stationary
regime with a spatial separation, in a direction orthogonal to their velocity.
As a consequence, and contrary to the ferromagnetic situation, a transient
regime necessarily occurs, with a finite acceleration, related to an inertia,
and that limits its time response. The formalism developed here, based on two
coupled Thiele equations, allows a quantitative description of this
phenomenon. The time constant associated to the transient regime scales
inversely with the antiferromagnetic coupling constant. We also show that the
coupling force reaches a maximal value at a finite skyrmion separation. This
sets a maximum velocity limit, beyond which the the coupling force cannot
stabilize the bound state.

Magnetic solitonic textures (domain walls, vortex cores, skyrmions…) can be
described as particles that can be excited and moved by external forces.
Therefore, similarly to real particules, the question of a mass and inertia
naturally rises. It was shown by Döringv. Döring (1948) that a texture
deformation related to the velocity can be responsible for the inertia of
magnetic domain walls. In a Bloch domain wall, the inner magnetization angle
and the velocity appear as conjugated variables, giving rise to inertia.
Magnetic skyrmions have been attracting a large interest for the last ten
years, since the first proposition of skyrmion-based spintronic applications
Fert et al. (2013). Due to their singular topology, they obey to a specific
dynamics which includes a gyrotropic deflectionThiele (1974); Everschor et al.
(2011). The question of a possible mass has been highly debatedBack et al.
(2020). While several theories and experiments have shown that topological
bubbles (stabilized without any Dzyaloshinskii Moriya interaction) display an
inertiaMoutafis et al. (2009); Makhfudz et al. (2012); Büttner et al. (2015a);
Makhfudz et al. (2012), the possibility of a mass in chiral and small
skyrmions remain elusive and debated.Back et al. (2020) Some studies report
the absence of a mass in ferromagnetic skyrmionsKravchuk et al. (2018), while
some cases observe an inertia related to a specific excitationShiino et al.
(2017); Wang et al. (2020) or to the vicinity of magnetic defectsPavlis and
Psaroudaki (2020). In contrast, the recent interest for antiferromagnetic
systemsJungwirth et al. (2016, 2018); Manchon et al. (2019) bring a new type
of skyrmion, where a mass, intrinsic to the texture, naturally arises from the
antiferromagnetic couplingVelkov et al. (2016); Kravchuk et al. (2018); Kim
and Tserkovnyak (2017); Kim et al. (2017).

Among the antiferromagnetic systems, synthetic antiferromagnets (SAF) are
quite appealing to the development of skyrmions, since they take the most of
ferromagnetic systems, where skyrmion are now routinely stabilizedJiang et al.
(2015); Boulle et al. (2016); Moreau-Luchaire et al. (2016); Soumyanarayanan
et al. (2017); Hrabec et al. (2017). In SAF, two ferromagnetic layers are
coupled antiferromagnetically through a non-magnetic spacer, via RKKY
interaction, which allows a tunable antiferromagnetic coupling. Predicted to
display a vanishing skyrmion gyrotropic deflection Barker and Tretiakov
(2016); Zhang et al. (2016) similarly to real antiferromagnetic
skyrmionsBarker and Tretiakov (2016), they enable the stabilization of small
skyrmionsLegrand et al. (2020) that can be moved, when excited by spin-orbit
torquesDohi et al. (2019).

In this paper, we study the skyrmion dynamics in synthetic antiferromagnets
and the consequences of their natural inertia. We derive an analytical
expression of the skyrmion coupling force, including the dependence with the
antiferromagnetic coupling strength and the skyrmion size. The inertia arises
from the opposite gyrovectors in the two layers, due to the antiferromagnetic
alignement, which spatially separates the two skyrmions, in a direction
perpendicular to the applied force. The skyrmion separation therefore appears
as the variable conjugated to the skyrmion pair velocity. The associated mass
is inversely proportional to the coupling strength and decreases when the
skyrmion size increases. We describe the transient regime associated to this
inertia, both for the bound state velocity and for the skyrmion separation,
with transient times of the order of a nanosecond. Our approach allows
investigations beyond the linear approximation used in effective models Velkov
et al. (2016); Kim and Tserkovnyak (2017); Kim et al. (2017), and predicts a
maximum skyrmion velocity beyond which the bound state is destroyed under the
action of the gyrotropic forces in each layer (the two skyrmions forming the
pair behave as independent particles).

##  I Thiele equation modeling

We consider two antiferromagnetically coupled magnetic layers, hosting Néel
skyrmions that form skyrmion pairs or bound states. To investigate their
dynamics, we explicitly describe each skyrmion, rather than using an effective
model that represent the system on a single layer with effective parameters
Velkov et al. (2016); Kim and Tserkovnyak (2017); Kim et al. (2017). We use
the Thiele formalism, initially developed for ferromagnetic systems, that
describes the dynamics of a rigid texture
$\mathbf{m}(\mathbf{r},t)=\mathbf{m}_{0}(\mathbf{r}-\mathbf{v}t)$
($\mathbf{m}$ is the local and time-dependent magnetization, $\mathbf{v}$ its
velocity and $\mathbf{m}_{0}$ is the stationary texture profile) Thiele (1973,
1974). In the two coupled layer system, we use two coupled Thiele equations

$\displaystyle\mathbf{G}\times\mathbf{v}_{1}-\alpha
D\mathbf{v}_{1}+\mathbf{F}_{1}+\mathbf{F}_{2/1}$ $\displaystyle=$
$\displaystyle 0$ (1a) $\displaystyle-\mathbf{G}\times\mathbf{v}_{2}-\alpha
D\mathbf{v}_{2}+\mathbf{F}_{2}+\mathbf{F}_{1/2}$ $\displaystyle=$
$\displaystyle 0.$ (1b)

where $\mathbf{v}_{1}$ and $\mathbf{v}_{2}$ are the skyrmion velocities in
each layers. This extends the Thiele formalism to a coupled systemBarker and
Tretiakov (2016); Büttner et al. (2015b) using the coupling forces
$\mathbf{F}_{2/1}$ and $\mathbf{F}_{1/2}$ that account for the action of the
textures on one another. Note that the action/reaction principle implies
$\mathbf{F}_{2/1}=-\mathbf{F}_{1/2}$. The three other terms in the equations,
similar to those in single ferromagnetic layers, correspond respectively to
the gyrotropic deflection (with $\mathbf{G}=L_{S}t\;n\mathrm{~{}}\mathbf{z}$
the gyrovector, $L_{S}=M_{S}/\gamma$ the angular momentum density, $M_{S}$ the
spontaneous magnetization, $\gamma$ the gyromagnetic ratio and $t$ the film
thickness), the energy dissipation (with $\alpha D=\alpha L_{S}t\;d$ and
$\alpha$ the Gilbert damping parameter) and the force of an external stimulus.
Here, we consider the motion induced by a spin-orbit torque (SOT) due to a
charge current $\mathbf{j}$ which flows in an adjacent layer along the $x$
direction.Sampaio et al. (2013) This charge current is converted to a spin
current $\theta_{H}\mathbf{j}\times\mathbf{y}$, via the spin Hall effect (SHE)
($\theta_{H}$ is the spin Hall angle), polarized along the $y$ direction Miron
et al. (2011); Liu et al. (2012a) and leads to a force
$\mathbf{F}=-\frac{\hbar}{2e}\mathbf{j}\theta_{H}\;f$ ($\hbar$ the reduced
Plank constant and $e$ the elementary charge). The three parameters $n$, $d$
and $f$ are related to the texture profile, as

$\displaystyle n$ $\displaystyle=$
$\displaystyle\iint\left(\frac{\partial\mathbf{m_{0}}}{\partial
x}\times\frac{\partial\mathbf{m_{0}}}{\partial
y}\right)\cdot\mathbf{m_{0}}d^{2}r$ (2a) $\displaystyle d$ $\displaystyle=$
$\displaystyle\iint\left(\frac{\partial\mathbf{m_{0}}}{\partial
x}\right)^{2}d^{2}r$ (2b) $\displaystyle f$ $\displaystyle=$
$\displaystyle\iint\left(m_{0}^{(x)}\frac{\partial m_{0}^{(z)}}{\partial
x}-m_{0}^{(z)}\frac{\partial m_{0}^{(x)}}{\partial x}\right)d^{2}r$ (2c)

with $m_{0}^{(i)}$ the components of $\mathbf{m_{0}}$. The skyrmion profile is
essentially described by four parameters, the skyrmion radius $R$, the
micromagnetic domain wall width parameter $\Delta$ (the lengthscale over which
the magnetization rotates Rohart and Thiaville (2013)), the skyrmion core
polarization $p$ and the chirality $c$ ($c=\pm 1$ respectively for clockwise
and counterclockwise spin rotation), which leads to simple expressions for
$n$, $d$ and $f$. The number $n$ corresponds to the skyrmion topology and only
depends on the skyrmion core polarization, with $n=4\pi p$. The number $d$ is
similar to an exchange integral and involves the magnetization rotation
lengthscale. For $R\gg\Delta$, the skyrmion has a bubble profile Rohart and
Thiaville (2013) and $d=2\pi R/\Delta$ Sampaio et al. (2013). For
$R\ll\Delta$, the integral does not go to zero due to the non-trivial skyrmion
topology Belavin and Polyakov (1975) and $d\rightarrow 4\pi$ for $R\rightarrow
0$ Büttner et al. (2015b). Over the full size range, $d$ can be approximated
by $d\approx 4\pi[\exp\left(-\frac{R}{2\Delta}\right)+\frac{R}{2\Delta}]$ to
account for the two limitsBerges et al. (2022); not (a). The number $f$ is
similar to a Dzyaloshinskii-Moriya interaction (DMI) exchange and involves the
skyrmion chirality, with $f\approx\pi^{2}cR$. In equation 1, $\mathbf{G}$ and
$D$ are calculated from the texture in layer 1. Since layer 2 is aligned in
the opposite direction, the gyrovector is opposite ($n$ is odd in
$\mathbf{m}$) but the dissipation is the same ($d$ is even in $\mathbf{m}$).
Even if $f$ is even in $\mathbf{m}$, we consider the possibility of two
different SOT-induced forces since, in a SAF stack, the torques can hardly be
identical (the spin current source being an adjacent layer, it can hardly be
the same for both magnetic layers since, for one of them, the spin current
source is in the spacer and therefore thinner than for the other one).

##  II Interlayer coupling between two skyrmions

The coupling between the two skyrmions originates from the antiferromagnetic
interaction between the two layers, with the energy
$E_{\mathrm{AF}}=J_{\mathrm{AF}}\iint(\mathbf{m}_{1}.\mathbf{m}_{2}+1)d^{2}r$,
where $\mathbf{m}_{1}$ and $\mathbf{m}_{2}$ are the skyrmion profiles in each
layer and $J_{\mathrm{AF}}$ is the interaction energy constant (positive to
promote the antiparallel alignement). For skyrmions with opposite core
magnetization, this energy is minimized when the skyrmions are aligned in the
$(x,y)$ plane. When separated by $\mathbf{\delta R}$, the energy increases.

For a small separation as compared to $R$ and $\Delta$, i.e. when there is a
significant overlap between the area where the magnetization rotates (skyrmion
periphery), $\mathbf{m}_{2}\approx-\mathbf{m}_{1}-\delta
R\partial_{u}\mathbf{m}_{1}-\frac{1}{2}\delta
R^{2}\partial^{2}_{u}\mathbf{m}_{1}$ (with $\delta R$ the magnitude of the
skyrmion separation and $u$ the coordinate along the separation direction) so
the interaction energy can be developed as

$\displaystyle E_{AF}^{small}$ $\displaystyle\approx$
$\displaystyle\frac{J_{\mathrm{AF}}}{2}\delta
R^{2}\iint\left(\frac{\partial\mathbf{m}_{1}}{\partial
u}\right)^{2}d^{2}r=\frac{J_{\mathrm{AF}}}{2}d\mathrm{~{}}\delta R^{2}$ (3)

where the geometrical parameter $d$ (eq. 2b), also involved in the
dissipation, is naturally found, since it represents an exchange energy. The
coupling force derives from this energy, as

$\mathbf{F}^{small}_{2/1}\approx-J_{\mathrm{AF}}d\mathrm{~{}}\mathbf{\delta
R}.$ (4)

In this linear approach, the antiferromagnetic interaction induces a spring
force between the two skyrmions with $k=J_{\mathrm{AF}}d$ the spring constant.
Note that since $d$ increases with the skyrmion radius, the spring constant
increases with the skyrmion size (in the limit of $R\gg\Delta$, $k$ varies
linearly with the skyrmion size).

For a large separation as compared to $\Delta$, the overlap between skyrmion
peripheries is negligible. For $R>>\Delta$, the antiferromagnetic energy
essentially corresponds to the overlap between the skyrmion cores, considering
rigid skyrmions as circular, homogeneously magnetized domains. The associated
energy

$\displaystyle E^{large}_{AF}$ $\displaystyle\approx$ $\displaystyle
2J_{AF}R^{2}\left[\pi-2\mathrm{~{}}\mathrm{acos}\left(\frac{\delta
R}{2R}\right)+\frac{\delta R}{R}\sqrt{1-\left(\frac{\delta
R}{2R}\right)^{2}}\right]\mathrm{~{}for~{}\delta R<2R}$ (5a)
$\displaystyle\approx$ $\displaystyle 2J_{AF}R^{2}\pi\mathrm{~{}for~{}\delta
R>2R}$ (5b)

leads to the coupling force

$\displaystyle\mathbf{F}^{large}_{2/1}$ $\displaystyle\approx$
$\displaystyle-8RJ_{AF}\sqrt{1-\left(\frac{\delta
R}{2R}\right)^{2}}\mathbf{u}\mathrm{~{}for~{}\delta R<2R}$ (6a)
$\displaystyle\approx$ $\displaystyle 0\mathrm{~{}for~{}\delta R>2R}.$ (6b)

Contrary to the first approximation, this force strength monotonously
decreases with $\delta R$. Note that for $\delta R>2R$ the overlap vanishes
and the coupling energy is constant, the coupling force zero and therefore the
skyrmion pair is broken.

Figure 1:  (a) Skyrmion coupling force vs. skyrmion separation, for different
skyrmion radius. The force is normalized to the antiferromagnetic coupling
constant and the skyrmion radius. The dots corresponds to a numerical
evaluation using two rigid skyrmions described by the profile obtained from
the ansatz in Ref. Tomasello et al., 2018. The full lines correspond to the
small skyrmion separation approximation (eq. 4), and the dashed lines
correspond to the large skyrmion separation approximation (eq. 6). Note that
for $R/\Delta>10$ and for the smallest separations, $F_{2/1}/R$ is independent
on the skyrmion size since, in this range, $d\approx 2\pi R/\Delta$
(therefore, $F_{2/1}\approx 2\pi J_{AF}R\delta R/\Delta$). (b) Maximum
coupling force and (c) corresponding skyrmion separation vs. skyrmion size. In
(c), the full line corresponds to $\delta R_{\mathrm{max}}=R$, valid for small
skyrmions.

The two approximations are compared with a numerical evaluation (see figure
1a). A micromagnetic calculation of the coupling force vs. the skyrmion
separation is difficult since, at a finite separation, the situation is not at
equilibrium. Complex numerical approaches, such as Nudged Elastic Band methods
Mills et al. (1995); Dittrich et al. (2002) could be used, but are too time
consuming in order to scan a broad range of parameters. Instead, under the
approximation of two rigid skyrmions, the force is easily calculated by
translating equilibrium profiles. In the absence of an analytical skyrmion
profile, we consider the following ansatz,Tomasello et al. (2018)

$\tan{\frac{\theta(r)}{2}}=\frac{R}{r}\exp\left(\frac{R-r}{\Delta}\right)$ (7)

with $\theta$ the angle of the magnetization with the sample normal direction,
and $r$ the skyrmion radial coordinate. This formula has been shown to provide
good results, whatever $R/\Delta$. Covering different skyrmion sizes
($R/\Delta$ from 0.1 to 1000), the two different regimes described above are
clearly visible with a good agreement with eqs. 4 and 6. Only for small
skyrmions ($R/\Delta\lesssim 1$), the large skyrmion separation regime
underestimates the coupling strength since, at such small sizes, skyrmions
cannot be approximated by circular homogeneous domains. A maximum coupling
force is observed [slightly larger in the numerical evaluation
($|\mathbf{F}_{2/1,\mathrm{max}|}=9$ to $10RJ_{AF}$) than the one anticipated
from the two approximations ($8RJ_{AF}$)], as plotted in Fig. 1b. The skyrmion
separation $\delta R_{\mathrm{max}}$ that corresponds to this force maximum is
of a few $\Delta$ for the largest skyrmions, and of the skyrmion radius for
the smallest skyrmions (Fig. 1c).

##  III Steady state regime and maximum velocity

###  III.1 Skyrmion separation in the steady state

In a steady state regime, a bound state is formed with
$\mathbf{v}_{1}=\mathbf{v}_{2}=\mathbf{v}$. The solution of eq. 1 is obtained
by summing the two Thiele equations. In the resulting equation, the coupling
force disappears and the bound state velocity is independent on the coupling
strength. Additionally, the bound state displays no gyrotropic deflection in
agreement with established resultsBarker and Tretiakov (2016); Zhang et al.
(2016); Woo et al. (2018); Hirata et al. (2019). Therefore, the skyrmion pair
moves along the total force $\mathbf{F}_{tot}=\mathbf{F}_{1}+\mathbf{F}_{2}$
direction, at a velocity $\mathbf{v}=\mathbf{F}_{tot}/2\alpha D$.

At the steady state velocity, the two skyrmions are separated, so that at the
single layer level, the coupling force compensates the gyrotropic force. In
the large separation regime, the coupling force strength monotonically
decreases with $\delta R$, so a larger separation is unable to further
compensate gyrotropic effects. Therefore, only the small separation regime is
explored, with a separation limit $\delta R_{\mathrm{max}}$ as discussed
before. If the two SOT forces are not identical, the skyrmion with the larger
applied force drags the second one via the interlayer coupling, another source
of skyrmion separation. From the Thiele equations, the steady state skyrmion
separation is

$\mathbf{\delta
R}=\frac{\mathbf{G}\times\mathbf{v}-\frac{1}{2}\Delta\mathbf{F}}{k}$ (8)

with $\Delta\mathbf{F}=\mathbf{F}_{2}-\mathbf{F}_{1}$ the SOT force
difference. The comparison with a micromagnetic simulation on CoFeB based
SAFnot (b) is shown in Fig 2 with a very good agreement. The skyrmion
separation is found to be inversely proportional to the spring constant, and
therefore to the antiferromagnetic coupling constant. In the specific case
where $\mathbf{F}_{1}=\mathbf{F}_{2}$ ($\eta=0$ in Fig. 2a), $\mathbf{\delta
R}$ is orthogonal to the velocity, since the skyrmion separation would only
compensate the gyrotropic effects in each layer. Similarly, neglecting the
gyrotropic vector, the separation would be along $\Delta\mathbf{F}$, therefore
along $\mathbf{v}$, since the skyrmion separation only compensates the force
difference. It is interesting to note that, in the general case, even if the
skyrmion separation depends on the gyrotropic effect and the force difference,
the bound state velocity is independent on the SOT force difference, as
anticipated by summing the two Thiele equations.

The relative importance of the effect due to the force difference decreases at
low damping constant. For a given applied total SOT force $\mathbf{F}_{tot}$,
varying the force repartition between the two layers, eq. 8 becomes

$\mathbf{\delta R}=\frac{\mathbf{G}\times\mathbf{F}_{tot}}{2k\alpha
D}-\frac{\eta\mathbf{F}_{tot}}{2k}.$ (9)

with $\eta$ defined to that $\mathbf{\Delta F}=\eta\mathbf{F}_{tot}$ and where
the result of the steady state velocity $\mathbf{v}=\mathbf{F}_{tot}/2\alpha
D$ bas been used. As a consequence, the smaller the dissipation, the closer
$\mathbf{\delta R}$ from the direction perpendicular to the applied force.

Figure 2:  (a) Micromagnetic simulation of the skyrmion separation as a
function of the SOT induced force repartition. The separation is calculated in
the steady state regime for a given current $j\theta_{H}=3\times 10^{11}$
A/m2, as a function of the repartition between the two layers. The velocity
(-60 m/s, not shown in the graphs) is independent of the repartition of the
force. The skyrmion separation is in good agreement with eq. 8, and only the
separation along the applied force depends on the force repartition. (b)
Scaled skyrmion separation as a function of the skyrmion radius for different
interlayer coupling constants and equal force repartition ($\eta=0$. The line
corresponds to the expected value from the model. The parameters used in the
simulation are those of a CoFeB based SAFnot (b). In (a) the coupling constant
is 0.2 mJ/m2 and the skyrmion radius is 9.2 nm ($D=2.5$ mJ/m2). A rather large
damping constant $\alpha=0.3$ has been used to enhance the influence of the
force difference than at more realistic values. In (b), the skyrmion radius is
varied by changing the DMI from -1.4 to -2.8 mJ/2.

The skyrmion separation is sensitive to the skyrmion size, in particular
through the dependence of $k$. For a given current density, the variation is
not monotonic, since $k$, $D$ and $F$ all vary with the size. This is
described in Fig 2(b), in the case of a current equipartition ($\eta=0$). For
the smallest skyrmions, the separation varies linearly with $R$: in this
regime, $d\approx 4\pi$ so the increase of $\delta R$ is only due to the
increase in velocity (for small skyrmions, $\mathbf{F}_{\mathrm{tot}}$
increases linearly with $R$). The separation reaches a maximum as $d$ starts
to vary significantly ($R>\Delta$). For the large skyrmion size $d\propto R$,
as well as $k$, $D$ and $\mathbf{F}_{\mathrm{tot}}$, which leads to $\delta
R\propto 1/R$.

###  III.2 Maximum velocity

This preceding description, with a steady state regime, is valid as long as
the coupling force can compensate the gyrotropic force. Since, the coupling
force displays a maximum at finite skyrmion separation, a too large velocity
implies a gyrotropic force that cannot be compensated. This leads to a
skyrmion decoupling: their separation increases with time and the skyrmion
dynamics is described by two independent Thiele equations. At the maximum
velocity $v_{\mathrm{max}}$, the strengths of the gyrotropic force reaches the
maximum coupling force so

$v_{\mathrm{max}}=\frac{|\mathbf{F}_{2/1,\mathrm{max}}|}{|\mathbf{G}|}\propto\frac{RJ_{AF}}{|\mathbf{G}|}$
(10)

where the proportionality factor varies slightly only with the skyrmion size
(from 10 for the smallest skyrmions to 9 for the largest skyrmions – see fig.
1b). Note that, a similar formula can be obtained using the small separation
coupling force approximation, and using the maximum skyrmion separation (Fig.
1c) instead of the maximum force ($v_{\mathrm{max}}=k\delta
R_{\mathrm{max}}/|\mathrm{G}|$). However, this approach leads to an
overestimation, since at the separation maximum the coupling force exact
result deviates from the linear approximation, as it can be observed in Fig.
1a. Eq. 10 summarizes the physical origin of the velocity limitation in SAF
skyrmions: the maximum velocity is proportional to the antiferromagnetic
coupling constant which aligns the two skyrmions and is inversely proportional
to the gyrovector strength which pushes the two skyrmions away. Additionally,
$v_{\mathrm{max}}$ is proportional the skyrmion radius, which shows that pairs
of larger skyrmions are better coupled (see Fig 3).

A situation with two different driving forces $\mathbf{F}_{1}$ and
$\mathbf{F}_{2}$ ($\eta\neq 0$) is another source of skyrmion separation.
Therefore, a reduction of the maximum velocity can be anticipated with

$|\mathbf{v}|_{\mathrm{max}}=\frac{|\mathbf{F}_{2/1,\mathrm{max}}|}{\sqrt{\mathbf{G}^{2}+\eta^{2}(\alpha
D)^{2}}}.$ (11)

The maximum velocity is the highest for equal forces ($\eta=0$) and minimum
when the driving force is applied to one layer only ($|\eta|=1$). However, the
correction, proportional to $(\alpha D)^{2}$, is small as compared to the
gyrotropic effect one, except for very large skyrmions (typically, for
$R/\Delta=10$ and $\alpha=0.1$, the maximum velocity is expected to decrease
by 10 % between the two extreme force repartition cases).

Figure 3:  Micromagnetic simulation of the maximum velocity for a skyrmion in
a CoFeB based SAFnot (b), as a function of the antiferromagnetic coupling
constant and for different skyrmion radius. The dots correspond to the
simulations and the lines to the model (eq. 10). The SOT induced force are
identical in both layers.

##  IV Transient regime

Figure 4:  Transient regime for skyrmions in a CoFeB based SAFnot (b). (a)
Skyrmion pair velocity and (b) separation vs. time after application of the
current ($j\theta_{H}=5\times 10^{9}$ A/m2) for different interlayer coupling
constants. The skyrmion radius is 9.2 nm ($D=2.5$ mJ/m2). The dotted line
corresponds to the maximum skyrmion separation according to our model. For
$J_{AF}=10^{-4}$ J/m2, the skyrmion coupling is not sufficient to establish a
steady state regime at this current density. Note that the velocity is
negative since the skyrmion chirality is left-handed.

To reach the stationary regime with spatially separated skyrmions, a transient
regime is expected with independent skyrmion velocities $\mathbf{v}_{1}$ and
$\mathbf{v}_{2}$, as shown in the simulations reported in Fig. 4. Within the
small separation approximation, combining the two Thiele equations provides
two independent equations for the bound state velocity $\mathbf{v}$ and the
skyrmion separation $\mathbf{\delta R}$:

$\displaystyle\mathbf{F}_{tot}-2\alpha D\mathbf{v}$ $\displaystyle=$
$\displaystyle\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{k}\mathbf{\dot{v}}$ (12a)
$\displaystyle\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{2k\alpha
D}\mathbf{\dot{\delta R}}+\mathbf{\delta R}$ $\displaystyle=$
$\displaystyle\frac{\mathbf{G}\times\mathbf{F}_{tot}}{2k\alpha
D}-\frac{\Delta\mathbf{F}}{2k}$ (12b)

The first equation is similar to a single Thiele equation, but with an
additional term which adds an inertia to the system, with a mass

$M=\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{k}.$ (13)

This inertia, which arises from the deformation of the skyrmion pair through
its separation, is similar to the Döring mass of domain walls (in the case of
the Bloch wall, the deformation corresponds to the rotation of the domain wall
magnetization). It is mostly due to the gyrotropic effect that separates the
skyrmions, and a small correction due to the energy dissipation. Note that the
deformation due to the force difference, which is in fact static, does not
participate to the inertia. This shows that a deformation induces inertia only
if it is directly related to the velocity. The second equation describes the
skyrmion separation dynamics. Here, the force difference is involved in the
equation. However it has no influence on the transient regime time scale, and
only affects the steady state value, as described in eq. 8. The two equations
imply a time constant

$\tau=\frac{M}{2\alpha D}=\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{2k\alpha D}.$
(14)

The skyrmion pair mass and the time constant are inversely proportional to the
antiferromagnetic coupling constant, in good agreement with the simulations
shown in Fig. 5. Consequently, the larger the antiferromagnetic coupling
strength, the smaller the inertia and the time constant. Note that, while the
dependence of the mass with the antiferromagnetic coupling has been
established in Ref. Velkov et al., 2016, our approach also shows a significant
dependence with the skyrmion size. Neglecting the dissipation term in the mass
($\alpha D\ll|\mathbf{G}|$), the mass is proportional to $\Delta/R$ for large
skyrmions (for small skyrmions, it does not diverge and tends to
$M=\mathbf{G}^{2}/4\pi J_{A}F$). The time constant is even more affected,
being proportional to $(\Delta/R)^{2}$.

Figure 5:  Transient time vs. the interlayer coupling constant for different
skyrmion radius in a CoFeB based SAFnot (b). The points correspond to the
simulation and the lines to the model (eq. 14). The result is valid for any
current densities, as long as the bound state is stable.

The time dependent equations also enable to understand the mechanics of the
skyrmions in synthetic antiferromagnets. Combining eqs. 4 and 12 in the
absence of a dissipation, we obtain Hamiltonian-like equations

$\displaystyle\dot{X}$ $\displaystyle=$
$\displaystyle\frac{1}{G}\frac{\partial E}{\partial\delta Y}$ (15a)
$\displaystyle\dot{\delta Y}$ $\displaystyle=$
$\displaystyle-\frac{1}{G}\frac{\partial E}{\partial X}$ (15b)

with $E$ the system total energy, $X$ the skyrmion pair position, $\delta Y$
the separation in the direction perpendicular to the direction of motion and
$G=\mathbf{G}.\mathbf{z}$. These equations show that the skyrmion pair
position in one direction is conjugated to the skyrmion pair separation in the
other direction.

The inertia obtained from the transient regime can be compared with the
skyrmion energy variation during the motion. Within the small skyrmion
separation approximation, the skyrmion pair energy in the stationary regime is

$E_{AF}=\frac{1}{2}\left[\frac{\Delta\mathbf{F}^{2}}{4k}+\frac{\mathbf{G}^{2}}{k}\mathbf{v}^{2}\right]$
(16)

The first term is independent on the skyrmion velocity and only depends on the
force difference. This underlines that the force difference does not
participate to the skyrmion inertia. The second term, proportional to
$\mathbf{v}^{2}$, is similar to a kinetic energy. However, the associated
mass, $\mathbf{G}^{2}/k$, is only one part of the inertial mass since the
inertia due to the dissipation does not participate to the kinetic energy. A
similar fact is obtained for the domain wall mass, when comparing the Döring’s
original theory v. Döring (1948), that ignores damping, with more elaborated
approaches.Malozemoff and Slonczewski (1979)

The model is compared with simulations of synthetic antiferromagnetsnot (b)
(see Figure 4). The steady state velocity is in perfect agreement with the
expected one and is independent on the interlayer coupling constant (see
Figure 4a). The transient regime observed in the simulation is more complex
than what is expected from the model, as observed in Fig. 4. While a first
order transition for the velocity and the skyrmion separation is generally
observed, additional oscillations are also visible. These are induced by the
magnetic precession of the regions surrounding the skyrmions when subjected to
SOT. Indeed, SOT not only pushes magnetic textures, but also acts on
homogeneous regions, Miron et al. (2011); Liu et al. (2012a) tilting magnetic
domains.Liu et al. (2012b) This transition, described by the LLG equation,
induces a precession of the magnetization around the new equilibrium
direction. As a consequence, in the simulation, two transient regimes
simultaneously occur, one for the surrounding regions and one for the skyrmion
pair, and interact together. A first order law fit on the velocity or the
skyrmion separation, therefore neglecting the oscillations, yields a transient
time which is close to our model expectation, in particular with a $1/J_{AF}$
dependency (see Figure 5). However, we note that the dependency with the
skyrmion size is not perfectly reproduced by the simulations, except for the
smallest coupling constants. It indicates that the two transient times are of
the same order of magnitude. While the transient time associated with the
domain tilt is obviously independent on the skyrmion size, the one associated
with the skyrmion acceleration is expected to decrease with $R$. Therefore,
for the largest sizes, the transient regime is dominated by the slowest time
constant, independent on the skyrmion size.

##  V Conclusion

We have described the skyrmion dynamics in synthetic antiferromagnets, through
the dynamics of two coupled rigid skyrmions. With a finite interlayer
coupling, the opposite gyrotropic forces tend to spatially separate the
skyrmions, which creates an inertia, inversely proportional to the
antiferromagnetic coupling strength. The skyrmion coupling force displays a
maximum at a finite skyrmion separation and therefore, a bound state maximum
velocity is observed, beyond which the skyrmions are decoupled and behave
independently. This approach can be extended to any antiferromagnetic
situation that can be split into two antiferromagnetically aligned lattices.
Note that in some cases (SAF with different materials and/or thicknesses,
rare-earth/transition metal ferrimagnetic alloys…), the two lattices are not
identical but the equations can still be solved (see appendix) with minor
corrections to the inertia. The main difference between synthetic
antiferromagnets and real antiferromagnets or ferrimagnets is the strength of
the antiferromagnetic coupling constant. In synthetic antiferromagnets, the
interfacial antiferromagnetic coupling constant is of the order of 1 mJ/m2. In
antiferromagnets or ferrimagnets, the volume antiferromagnetic coupling
constant is larger than $10^{8}$ J/m3. This volume interaction should be
integrated along the sample thickness to be used in eq. 4. For a sample
thickness of a few nanometers, it leads to an interaction strength 100 times
larger than in SAF and therefore to a time constant 100 times smaller, which
can generally be neglected.

##  Appendix A Inertia in unbalanced situations

The calculations in the paper concern the inertia for a system where the two
antiferromagnetically-coupled lattices are identical (except for the SOT-
induced forces). For an unbalanced SAF (e.g. different layer thicknesses) or a
rare-earth/transition metal ferrimagnetic alloy, the situation is different
since the gyrovectors and dissipation parameters are not the same in the two
lattices. With minimal approximations, the two Thiele equations can still be
solved.

We consider the following pair of Thiele equations

$\displaystyle\mathbf{G}_{1}\times\mathbf{v}_{1}-\alpha_{1}D_{1}\mathbf{v}_{1}+\mathbf{F}_{1}-k\mathbf{\delta
R}=0$ (17a)
$\displaystyle\mathbf{G}_{2}\times\mathbf{v}_{2}-\alpha_{2}D_{1}\mathbf{v}_{2}+\mathbf{F}_{2}+k\mathbf{\delta
R}=0$ (17b)

where the indices refer to each lattice. Neglecting the terms in
$\alpha_{1}D_{1}-\alpha_{2}D_{2}$, i.e. considering that the dissipation in
the two lattices are close, the velocity is given by

$\displaystyle\mathbf{G}_{tot}\times\mathbf{v}-2\alpha
D\mathbf{v}+\mathbf{F}_{tot}=$
$\displaystyle\frac{(\Delta\mathbf{G})^{2}-\mathbf{G}_{tot}^{2}+(2\alpha
D)^{2}}{4k}\mathbf{\dot{v}}-\frac{\alpha
D}{2k}\mathbf{G}_{tot}\times\mathbf{\dot{v}}.$ (18)

with $\mathbf{G}_{tot}=\mathbf{G}_{1}+\mathbf{G}_{2}$,
$\Delta\mathbf{G}=\mathbf{G}_{1}-\mathbf{G}_{2}$ and
$\mathbf{F}_{tot}=\mathbf{F}_{2}+\mathbf{F}_{1}$. Note that in an
antiferromagnetic situation, $\mathbf{G}_{1}$ and $\mathbf{G}_{2}$ are
antiparallel so $|\mathbf{\Delta G}|\gg|\mathbf{G}_{tot}|$. However, contrary
to the balanced antiferromagnetic situation, the two gyrovectors do not
compensate and a finite deflection remains in the equation. The inertia is
more complex than for the balanced situation, with an additional term. This
term differs from a simple mass, as defined in classical mechanics, since it
involves a cross product between the total gyrovector and the skyrmion pair
acceleration. Since the gyrovector is related to the skyrmion topology, this
additional inertia can be considered as a topology induced inertia. Note
however that this term is small as compared to the first one and can generally
be neglected so that the skyrmion pair time constant is

$\tau\approx\frac{(\Delta\mathbf{G})^{2}-\mathbf{G}_{tot}^{2}+(2\alpha
D)^{2}}{8k\alpha D}$ (19)

close to that of a balanced SAF.

##  Appendix B Spin transfer torque induced skyrmion motion

Another possibility to move skyrmions in magnetic multilayers is the spin-
transfer torque (STT), although it is expected to be less efficient than spin-
orbit torques.Sampaio et al. (2013) The calculation derived from eq. 1 can
also be performed with this other excitation. The coupled coupled Thiele
equations become

$\displaystyle\mathbf{G}\times(\mathbf{v}_{1}-\mathbf{u})-D(\alpha\mathbf{v}_{1}-\beta\mathbf{u})-k\mathbf{\delta
R}$ $\displaystyle=$ $\displaystyle 0$ (20a)
$\displaystyle-\mathbf{G}\times(\mathbf{v}_{2}-\mathbf{u})-D(\alpha\mathbf{v}_{2}-\beta\mathbf{u})+k\mathbf{\delta
R}$ $\displaystyle=$ $\displaystyle 0.$ (20b)

where the effect of a current flowing in the magnetic layers through the
magnetic textures is represented by $\mathbf{u}$ (uniform to a velocity and
proportional to the current density) and a non-adiabaticity parameter
$\beta$.Thiaville, A. et al. (2005); Sampaio et al. (2013) In the limit of the
linear interaction regime, these equations lead to

$\displaystyle-2D(\alpha\mathbf{v}-\beta\mathbf{u})$ $\displaystyle=$
$\displaystyle\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{k}\mathbf{\dot{v}}$ (21a)
$\displaystyle\frac{\mathbf{G}^{2}+(\alpha D)^{2}}{2k\alpha
D}\mathbf{\dot{\delta R}}+\mathbf{\delta R}$ $\displaystyle=$
$\displaystyle\frac{1}{k}\left(\frac{\beta}{\alpha}-1\right)\mathbf{G}\times\mathbf{u}$
(21b)

The first equation for the average velocity corresponds to the motion of a
non-topological magnetic texture under STT, with a mass
$M=\frac{1}{k}[\mathbf{G}^{2}+(\alpha D)^{2}]$. The second equation describes
the dynamics of the skyrmion separation. Both equation display the same time
constant $\tau=M/2\alpha D$. In the stationary regime,
$\mathbf{v}=\frac{\beta}{\alpha}\mathbf{u}$, a typical result for non-
topological textures (including domain walls) moved under STTThiaville, A. et
al. (2005) and $\mathbf{\delta
R}=\frac{1}{k}(\frac{\beta}{\alpha}-1)\mathbf{G}\times\mathbf{u}$. This result
shows that the inertial concepts derived from SOT excitation are general and
not specific to a single excitation: the mass is the same and has its origin
from the opposite gyrotropic effect in both layers due to the
antiferromagnetic coupling.

###### Acknowledgements.

We are grateful to André Thiaville for enlightening discussions and a critical
reading of the manuscript. This work was supported by the French National
Research Agency (ANR) [under contract ANR-17-CE24-0025 (topsky) and a public
grant overseen as part of the “Investissements d’Avenir” program (Labex
NanoSaclay, reference: ANR-10-LABX-0035), SPICY] and by an Indo-French
collaborative project supported by CEFIPRA (IFC/5808-1/2017).

## References

  * v. Döring (1948) W. v. Döring, Z. Naturforsch. 3a, 373 (1948).
  * Fert et al. (2013) A. Fert, V. Cros, and J. Sampaio, Nat. Nano. 8, 152 (2013).
  * Thiele (1974) A. Thiele, J. Appl. Phys. 45, 377 (1974).
  * Everschor et al. (2011) K. Everschor, M. Garst, R. A. Duine, and A. Rosch, Phys. Rev. B 84, 064401 (2011).
  * Back et al. (2020) C. Back, V. Cros, H. Ebert, K. Everschor-Sitte, A. Fert, M. Garst, T. Ma, S. Mankovsky, T. Monchesky, M. Mostovoy, et al., J. Phys. D. 53, 363001 (2020).
  * Moutafis et al. (2009) C. Moutafis, S. Komineas, and J. A. C. Bland, Phys. Rev. B 79, 224429 (2009).
  * Makhfudz et al. (2012) I. Makhfudz, B. Krüger, and O. Tchernyshyov, Phys. Rev. Lett. 109, 217201 (2012).
  * Büttner et al. (2015a) F. Büttner, C. Moutafis, M. Schneider, B. Krüger, C. Günther, J. Geilhufe, C. v. Korff Schmising, J. Mohanty, B. Pfau, S. Schaffert, et al., Nature Phys. 11, 225 (2015a).
  * Kravchuk et al. (2018) V. P. Kravchuk, D. D. Sheka, U. K. Rößler, J. van den Brink, and Y. Gaididei, Phys. Rev. B 97, 064403 (2018).
  * Shiino et al. (2017) T. Shiino, K. Kim, K. Lee, and B. Park, Sci. Rep. 7, 13993 (2017).
  * Wang et al. (2020) D. Wang, H.-B. Braun, and Y. Zhou, _Dynamical mass generation for ferromagnetic skyrmions in two dimensions_ (2020), URL https://arxiv.org/abs/2009.07486.
  * Pavlis and Psaroudaki (2020) A. Pavlis and C. Psaroudaki, Phys. Rev. Research 2, 032058(R) (2020).
  * Jungwirth et al. (2016) T. Jungwirth, X. Marti, P. Wadley, and J. Wunderlich, Nat. Nano. 11, 231 (2016).
  * Jungwirth et al. (2018) T. Jungwirth, J. Sinova, A. Manchon, X. Marti, J. Wunderlich, and C. Felser, Nat. Phys. 14, 200 (2018).
  * Manchon et al. (2019) A. Manchon, J. Železný, I. M. Miron, T. Jungwirth, J. Sinova, A. Thiaville, K. Garello, and P. Gambardella, Rev. Mod. Phys. 91, 035004 (2019).
  * Velkov et al. (2016) H. Velkov, O. Gomonay, M. Beens, G. Schwiete, A. Brataas, J. Sinova, and R. Duine, New. J. Phys. 18, 075016 (2016).
  * Kim and Tserkovnyak (2017) S. K. Kim and Y. Tserkovnyak, Appl. Phys. Lett. 111, 032401 (2017).
  * Kim et al. (2017) S. K. Kim, K.-J. Lee, and Y. Tserkovnyak, Phys. Rev. B 95, 140404(R) (2017).
  * Jiang et al. (2015) W. Jiang, P. Upadhyaya, W. Zhang, G. Yu, M. B. Jungfleisch, F. Y. Fradin, J. E. Pearson, Y. Tserkovnyak, K. L. Wang, O. Heinonen, et al., Science 349, 283 (2015).
  * Boulle et al. (2016) O. Boulle, J. Vogel, H. Yang, S. Pizzini, D. de Souza Chaves, A. Locatelli, T. O. Mentes, A. Sala, L. D. Buda-Prejbeanu, O. Klein, et al., Nat. Nano. 11, 449 (2016).
  * Moreau-Luchaire et al. (2016) C. Moreau-Luchaire, C. Moutafis, N. Reyren, J. Sampaio, C. A. F. Vaz, N. Van Horne, K. Bouzehouane, K. Garcia, C. Deranlot, P. Warnicke, et al., Nat. Nano. 11, 444 (2016).
  * Soumyanarayanan et al. (2017) A. Soumyanarayanan, M. Raju, A. L. Gonzalez Oyarce, A. K. C. Tan, M.-Y. Im, A. Petrović, P. Ho, K. H. Khoo, M. Tran, C. K. Gan, et al., Nat. Mat. 16, 898 (2017).
  * Hrabec et al. (2017) A. Hrabec, J. Sampaio, M. Belmeguenai, I. Gross, R. Weil, S. M. Chérif, A. Stashkevich, V. Jacques, A. Thiaville, and S. Rohart, Nat. Com. 8, 15765 (2017).
  * Barker and Tretiakov (2016) J. Barker and O. A. Tretiakov, Phys. Rev. Lett. 116, 147203 (2016).
  * Zhang et al. (2016) X. Zhang, Y. Zhou, and M. Ezawa, Nat. Com. 7, 10293 (2016).
  * Legrand et al. (2020) W. Legrand, D. Maccariello, F. Ajejas, S. Collin, A. Vecchiola, K. Bouzehouane, N. Reyren, V. Cros, and A. Fert, Nat. Mater. 19, 34 (2020).
  * Dohi et al. (2019) T. Dohi, S. DuttaGupta, S. Fukami, and H. Ohno, Nat. Com. 10, 5153 (2019).
  * Thiele (1973) A. Thiele, Phys. Rev. Lett. 30, 230 (1973).
  * Büttner et al. (2015b) F. Büttner, I. Lemesh, and G. S. D. Beach, Sci. Rep. 8, 4464 (2015b).
  * Sampaio et al. (2013) J. Sampaio, V. Cros, S. Rohart, A. Thiaville, and A. Fert, Nat. Nano. 8, 839 (2013).
  * Miron et al. (2011) I. M. Miron, K. Garello, G. Gaudin, P.-J. Zermatten, M. V. Costache, S. Auffret, S. Bandiera, B. Rodmacq, A. Schuhl, and P. Gambardella, Nature (London) 476, 189 (2011).
  * Liu et al. (2012a) L. Liu, C.-F. Pai, Y. Li, H. W. Tseng, D. C. Ralph, and R. A. Buhrman, Science 336, 555 (2012a).
  * Rohart and Thiaville (2013) S. Rohart and A. Thiaville, Phys. Rev. B 88, 184422 (2013).
  * Belavin and Polyakov (1975) A. A. Belavin and A. M. Polyakov, JETP Lett. 22, 245 (1975).
  * Berges et al. (2022) L. Berges, E. Haltz, S. Panigrahy, S. Mallick, R. Weil, S. Rohart, A. Mougin, and J. a. Sampaio, _Size-dependent mobility of skyrmions beyond pinning in ferrimagnetic gdco thin films_ (2022), URL https://arxiv.org/abs/2208.04384.
  * not (a) To account for the two limit at $R=0$ ($d=4\pi$Büttner et al. (2015b)) and $R\gg\Delta$ ($d=2\pi R/\Delta$Sampaio et al. (2013)), we represent $d$ as $4\pi[\exp(-R/2\Delta)+R/2\Delta)]$. In this formulation, the weight of the $R=0$ limit decreases exponentially with the radius, and $d$ variations are negligible close to $R=0$. A comparison with micromagnetic simulations can be found in Ref. Berges et al., 2022.
  * Tomasello et al. (2018) R. Tomasello, K. Y. Guslienko, M. Ricci, A. Giordano, J. Barker, M. Carpentieri, O. Chubykalo-Fesenko, and G. Finocchio, Phys. Rev. B 97, 060402(R) (2018).
  * Mills et al. (1995) G. Mills, H. Jónsson, and G. K. Schenter, Surf. Sci. 324, 304 (1995).
  * Dittrich et al. (2002) R. Dittrich, T. Schrefl, D. Suess, W. Scholz, H. Forster, and J. Fidler, Journal of Magnetism and Magnetic Materials 250, 12 (2002).
  * Woo et al. (2018) S. Woo, K. M. Song, X. Zhang, Y. Zhou, M. Ezawa, X. Liu, S. Finizio, J. Raabe, N. J. Lee, S.-I. Kim, et al., Nat. Comm. 9, 959 (2018).
  * Hirata et al. (2019) Y. Hirata, D.-H. Kim, S. K. Kim, D.-K. Lee, S.-H. Oh, D.-Y. Kim, T. Nishimura, T. Okuno, Y. Futakawa, H. Yoshikawa, et al., Nat. Nano. 14, 232 (2019).
  * not (b) See Supplemental Material at [URL will be inserted by publisher] for the description of the simulation parameters and the variation of the skyrmion size with the DMI value.
  * Malozemoff and Slonczewski (1979) A. Malozemoff and J. Slonczewski, _Magnetic domain walls in bubble materials_ , vol. 81 (Academic Press New York, 1979).
  * Liu et al. (2012b) L. Liu, O. J. Lee, T. J. Gudmundsen, D. C. Ralph, and R. A. Buhrman, Phys. Rev. Lett. 109, 096602 (2012b).
  * Thiaville, A. et al. (2005) Thiaville, A., Nakatani, Y., Miltat, J., and Suzuki, Y., Europhys. Lett. 69, 990 (2005).