# A two-dimensional magneto-optical trap of dysprosium atoms as a compact
source
for efficient loading of a narrow-line three-dimensional magneto-optical trap

Shuwei Jin1,†, Jianshun Gao1,†, Karthik Chandrashekara1, Christian
Gölzhäuser1, Joschka Schöner1, Lauriane Chomaz1,⋆ 1 Physikalisches Institut,
Universität Heidelberg, Im Neuenheimer Feld 226, 69120, Heidelberg, Germany.

###### Abstract

We report on a scheme for loading dysprosium atoms into a narrow-line three-
dimensional magneto-optical trap (3D MOT). Our innovative approach replaces
the conventional Zeeman slower with a 2D MOT operating on the broad 421-nm
line to create a high-flux beam of slow atoms. Even in the absence of a push
beam, we demonstrate efficient loading of the 3D MOT, which operates on the
narrower 626-nm intercombination line. Adding push beams working at either 421
nm or 626 nm, significant enhancement of the loading rate is achieved. We
reach the best performance, with an enhancement factor of $3.6$, using a push
beam red-detuned to the 626-nm line. With loading rates greater than $10^{8}$
atoms/s achieved at a moderate oven reservoir temperature of $800\,^{\circ}$C,
our method offers similar or greater performance than Zeeman-slower-based
systems. Our 2D-MOT-based approach constitutes a promising first step for
state-of-the-art quantum gas experiments with several advantages over the
Zeeman-slower-based setup and is readily adaptable to other open-shell
lanthanides.

##  I Introduction

Over the last decade, ultracold gases of open-$f$-shell lanthanide atoms, such
as Er and Dy, have become a platform of choice for studying novel quantum
phenomena Chomaz  _et al._ (2022). Their electronic ground state’s structures
provide these atoms with remarkable properties: the presence of a closed outer
$6s$ shell yields electronic transitions with similar properties to those of
Yb or Sr. The open $4f$ shell confers on these atoms an even wider variety of
transitions, a large effective spin, and a large magnetic moment, among the
largest of the periodic table. In particular, the latter feature allows to
explore the quantum effects of long-range and anisotropic interactions using
ultracold gases of open-shell lanthanides Chomaz  _et al._ (2022).

The spectrum complexity brought up by the open-shell character of the magnetic
lanthanides had however refrained the scientific effort to laser-cool such
species for years. In 2006, a breakthrough experiment by McClelland and
Hanssen McClelland and Hanssen (2006) demonstrated the possibility of laser-
cooling Er on the tens-of-MHz-broad transition at 401 nm, despite the
existence of numerous decay channels. This work paved the way for many
experiments bringing open-shell lanthanides to ultracold temperatures Chomaz
_et al._ (2022); Berglund  _et al._ (2008); Youn  _et al._ (2010); Lu  _et
al._ (2010); Sukachev  _et al._ (2010); Frisch  _et al._ (2012); Maier  _et
al._ (2014); Miao  _et al._ (2014); Dreon  _et al._ (2017); Cojocaru  _et al._
(2017); Ilzhöfer  _et al._ (2018); Inoue  _et al._ (2018); Phelps  _et al._
(2020); Lunden  _et al._ (2020) and quantum degeneracy Lu  _et al._ (2011,
2012); Aikawa  _et al._ (2012, 2014); Davletov  _et al._ (2020).

Figure 1:  Experimental Setup. (a) Energy level diagram of Dy in wavelength
$\lambda$ for the levels of total electronic angular-momentum quantum number
$J=8$ and $J=9$. The state’s color indicates its parity (blue, odd; red,
even). The used transitions are marked by the arrows. (b) Sketch of the
experimental setup showing the main vacuum components and optical paths.
inset, typical absorption image of a 3D MOT. (c) Section view of the setup at
the oven and 2D-MOT chamber. The aperture set at the effusion cell output is
illustrated as well as an instance of magnet block. (d) Sketch of the 3D-MOT
capture process within our experimental geometry. The atoms in the jet
entering with $v\leq v_{\rm cap}$ can be stopped within the capture region of
radius $R_{\rm cap}$. This region is limited by the beam extent, and therefore
the viewports’ clear aperture.

A particular advantage of the lanthanides is the existence of closed narrow
transitions (i.e. of linewidth comparable to or smaller than the recoil
frequency, see e.g. Dalibard (2015)). These transitions allow for reaching
ultra-low temperatures even within the Doppler cooling regime and have been
exploited in laser cooling schemes. In particular, magneto-optical traps
(MOTs) of open-shell lanthanides working on ultra-narrow lines have been
loaded from broad-line three-dimensional (3D) MOTs Berglund  _et al._ (2008);
Youn  _et al._ (2010); Lu  _et al._ (2010), similarly to Sr and Yb Katori  _et
al._ (1999); Kuwamoto  _et al._ (1999). Among all the closed narrow
transitions of lanthanides, the intercombination line stands out by its
intermediate linewidth of a few hundreds of kHz, which allows for MOTs with
capture velocities on the order of 10 m/s and Doppler temperatures in the few
$\mu$K range. Following works on Yb Kuwamoto  _et al._ (1999), MOTs working on
the intercombination line have been directly loaded from slow atomic beams,
enabling a simplified cooling scheme Chomaz  _et al._ (2022); Frisch  _et al._
(2012); Maier  _et al._ (2014); Dreon  _et al._ (2017); Cojocaru  _et al._
(2017).

Another particular feature of open-shell lanthanides is their high melting
point $\sim 1000\,\,^{\circ}$C. To produce the slow atomic beam required for
3D-MOT loading, all previous open-shell-lanthanide experiments use a similar
scheme based on a high-temperature oven aligned with an axial Zeeman slower
working on the atom’s broadest cooling transition Chomaz  _et al._ (2022);
Berglund  _et al._ (2008); Youn  _et al._ (2010); Lu  _et al._ (2010);
Sukachev  _et al._ (2010); Frisch  _et al._ (2012); Maier  _et al._ (2014);
Miao  _et al._ (2014); Dreon  _et al._ (2017); Cojocaru  _et al._ (2017);
Ilzhöfer  _et al._ (2018); Inoue  _et al._ (2018); Phelps  _et al._ (2020);
Lunden  _et al._ (2020). Besides Zeeman slowers, two-dimensional (2D) MOTs
have been proven to be convenient high-flux sources of slow atoms for many of
the laser-cooled atomic species as demonstrated in the seminal work of ref.
Dieckmann  _et al._ (1998). This has been applied to a wide range of species,
in particular species using an oven for a first vaporization stage, like Li
Tiecke  _et al._ (2009), Na Lamporesi  _et al._ (2013), Sr Nosske  _et al._
(2017); Barbiero  _et al._ (2020), Yb Vaidya  _et al._ (2015); Saskin  _et
al._ (2019) and more recently even to species needing a refrigerated source
like Hg Guo  _et al._ (2023). Among others, the advantages of 2D-MOT sources
over Zeeman slowers are the high compactness of the setups and the absence of
a direct view between the science chamber and the oven output. This is
particularly beneficial for species like open-shell lanthanides, for which
Zeeman slower light needs to be reflected from a mirror inside the vacuum
setup in order to avoid material deposition on the viewport. In this case, a
2D MOT allows for better optical access around the science chamber and a
simplified integration of a glass cell. Yet despite all these benefits, a 2D
MOT has not been realized for open-shell lanthanides.

In this paper, we report on an apparatus in which a 3D MOT of Dy atoms working
on the intercombination line is loaded from a 2D MOT operating on the broadest
cooling transition. We observe efficient loading of the 3D MOT with comparable
rates $\phi_{\rm 3D}\gtrsim 10^{8}$ atoms/s and saturation atom numbers
$N_{\rm sat}\approx 3\times 10^{8}$ as in Dy experiments based on Zeeman
slowers, but at a lower oven operation temperature. The 3D MOT can be loaded
even in the absence of an additional beam that pushes the atoms from the 2D to
the 3D MOT (hereafter called the push beam). Using a push beam working on
either the broad or the narrow transition, the loading of the 3D-MOT is
significantly enhanced with a maximal enhancement factor of $3.6$ achieved
using a red-detuned push beam close to the intercombination line.

The paper is organized as follows. In Sec. II, we briefly review the relevant
characteristics of Dy atoms and the main components of our experimental setup,
and present the 2D-MOT and 3D-MOT schemes. In Sec. III, the measurement scheme
used for our optimization and characterization of the 2D-/3D-MOT source is
described. In Sec.IV, we report on the optimization procedure applied to our
2D-MOT source and its achievement in the absence of a push beam. In Sec. V, we
further introduce a push beam and describe the observed enhancement of the
3D-MOT loading, comparing broad- and narrow-line push configurations. Finally,
in Sec. VI, the optimal 3D-MOT parameters for its loading in the absence and
presence of push beams are investigated, compared, and comprehended in
relation to the velocity features of the 2D-MOT atomic beam and of the 3D-MOT
capture process.

##  II Experimental setup and 2D/3D-MOT scheme

The Dy electronic levels and transitions of interest for this work are
depicted in Fig. 1 (a) and connect its ground state of total angular-momentum
quantum number $J=8$ to excited states with $J^{\prime}=9$ Martin  _et al._
(1978). The broad transition of wavelength $\lambda_{\rm
421}=2\pi/k_{421}=421.291\,$nm has a linewidth of $\Gamma_{\rm 421}=2\pi\times
32.2\,$MHz and a saturation intensity of $I_{\rm
sat}^{421}=564\,\textrm{W/m}^{2}$. The narrower intercombination line of
wavelength $\lambda_{\rm 626}=2\pi/k_{626}=626.082\,$nm has a linewidth of
$\Gamma_{\rm 626}=2\pi\times 136\,$kHz and a saturation intensity of $I_{\rm
sat}^{626}=720\,\textrm{mW/m}^{2}$. The ground state, 421-nm and 626-nm
excited states have similar $g$ factors of $g_{J}=1.24$,
$g_{J^{\prime}}^{(421)}=1.22$, and $g_{J^{\prime}}^{(626)}=1.29$ respectively.

The 421-nm and 626-nm laser lights used to address these transitions are
generated by two commercial frequency-doubled amplified diode lasers (DLC TA-
SHG PRO) from TOPTICA Photonics AG. The lasers are operated with 800mW and
1.6W output power respectively. Both lasers are frequency-locked to a
commercial ultra-low expansion cavity from Stable Laser Systems using the
offset-sideband Pound-Drever-Hall scheme Black (2001). The cavity has a 1.46
GHz free spectral range and a finesse determined by means of cavity-ring-down
measurements to be 21170(60) at 842 nm and 20760(40) at 626 nm, corresponding
respectively to the wavelengths used to lock the 421-nm and 626-nm lasers.

The vacuum apparatus for our experiment is depicted in Fig. 1 (b). It consists
of one high-temperature dual-filament effusion cell (DFC-40-25-WK-2B) from
CreaTec Fischer & Co. GmbH, and of two chambers – the 2D-MOT chamber and the
science chamber – connected via a differential pumping section. Most vacuum
components have been produced by SAES Rial Vacuum out of stainless steel
(grade 316L or 316LN) or titanium. Pressures on the order of $10^{-9}$ mbar
and $10^{-11}$ mbar are achieved in the 2D-MOT and science chambers
respectively.

Solid Dy pieces are placed in the reservoir region of the effusion cell and
are vaporized by heating this region up to $800\,^{\circ}$C. The relatively
low temperature of the reservoir region was chosen to spare material and allow
for a long lifetime of the source. Figure 1 (c) details the setup design from
the effusion-cell output to the 2D-MOT chamber. The oven is inserted into the
2D-MOT chamber, transversely to the 2D-MOT axis ($x$-axis). The distance
between the oven’s last aperture and the center of the 2D-MOT chamber (40.8mm)
is minimized with the aim of maximizing the incoming atomic flux in the
2D-MOT. A Dy atomic vapor jet is formed by a custom-designed set of apertures
located in the hot-lip region of the effusion cell, which is heated up to
$1100\,^{\circ}$C. A last cold aperture, connected to the oven water-cooling
system, enables to filter out the part of the atomic jet exiting the oven at
angles larger than $7.5\,^{\circ}$.

The atoms exiting the oven into the 2D-MOT chamber experience the forces
induced by the cooling beams. The cooling beams are made up of a single 421-nm
laser beam propagating through the chamber in a bow-tie $\sigma_{+}\sigma_{-}$
retro-reflected configuration. The bow-tie plane ($yz$-plane), together with
the magnetic-field configuration (see below), defines the 2D-MOT axis as the
orthogonal $x$ axis (see Fig. 1 (b,c)). The laser beam has a power of
$430\,$mW and waist of $16\,$mm 111In the main-text measurement, we use the
maximal power available while the 2D-MOT power dependence is reported in App.
A.. The laser beam path does not include any active optical components in
order to ensure a high effective power on the atoms. Therefore, the beam
frequency can only be altered through the laser locking point, and its power
via mechanical means.

Eight stacks of permanent magnets are placed symmetrically on both sides of
the 2D-MOT chamber, as partially illustrated in Fig. 1 (c), to provide the
2D-MOT magnetic field. The magnetic field is zero-valued along the 2D-MOT axis
and oriented along the cooling beams on their propagation axes. It has a
roughly uniform gradient in the $yz$ plane over the chamber’s central region
whose magnitude, $b^{\prime}_{\rm 2D}$, can be adjusted by changing the number
of magnets. Hereafter, the values of $b^{\prime}_{\rm 2D}$ correspond to the
theoretical expectations for a perfect arrangement of magnet stacks and no
other magnetic source. Increasing the number of magnets by one in each of the
8 blocks increases the gradient by approximately $4.4$ G/cm. A gradient of up
to $44.4$ G/cm can be generated with our current magnet-holder design. To
optimize the performance of the 2D MOT, we adjust the divergence of the
cooling beam, its detuning from resonance $\delta_{\rm 2D}$, as well as the
number and positions of the magnets. We note that due to the permanent magnets
implemented, adjusting the 2D-MOT magnetic field involves an important manual
aspect.

The atoms trapped in the 2D MOT can travel along the $x$ direction and reach
the center of the science chamber. The distance between the two chamber
centers along $x$ was minimized during our design process. It equals $347$ mm
and includes a $55.7$ mm-long differential pumping section, which is inserted
into the 2D-MOT chamber. The differential pumping section is of conical shape
and has an aperture diameter of 2 mm on one end and 5.9 mm on the other. At
the center of the science chamber, a 3D MOT is formed using three orthogonal
retro-reflected 626 nm laser beams and a pair of magnetic coils. The coils are
aligned along the $z$-axis and connected in anti-Helmholtz configuration to
provide a magnetic gradient $b^{\prime}_{\rm 3D}$ up to 4 G/cm in the current
configuration. Hereafter, the values of $b^{\prime}_{\rm 3D}$ correspond to
the gradient in the $xy$ plane extracted from numerical calculations using our
coil geometry. In the present setup an additional pair of magnetic coils,
aligned along the $z$ axis and connected in Helmholtz configuration can be
used to apply a tunable offset magnetic field along the $z$ direction, aligned
with gravity.

An important constraint related to the narrow-linewidth transition used for
our 3D-MOT scheme is the low capture velocity of the 3D MOT, $v_{\rm cap}$.
The principle of the capturing process is illustrated in Fig. 1 (d): atoms of
the 2D-MOT jet with an initial axial velocity $v_{x}<v_{\rm cap}$ have to be
stopped by the 3D-MOT radiation pressure within the 3D-MOT capture region of
radius $R_{\rm cap}$. An upper bound on $v_{\rm cap}$ can be estimated in the
limit of an infinitely saturated transition, where the atoms scatter photons
at a rate $\Gamma_{626}/2$ independent of the light detuning. In this
approximation, a constant radiation pressure force is exerted onto the atoms
over the 3D-MOT capture region, yielding $v_{\rm cap}\leq\sqrt{2R_{\rm
cap}\hbar\Gamma_{626}k_{626,x}/{m}}$, with $k_{626,x}$ the recoil momentum
transferred by one photon along the $x$ axis ($\hbar$ is the reduced Planck’s
constant), see e.g. Ref Dalibard (2015). With our geometry [see Fig. 1 (d)],
$k_{626,x}=k_{626}/\sqrt{2}$ and even with MOT-beam sizes equal to the
viewports’ clear aperture (yielding $R_{\rm cap}=35/\sqrt{2}$ mm), we estimate
the maximum capture velocity of our 626-nm 3D MOT to be $v_{\rm cap}\lesssim
11\,$m/s. To favor the loading of atoms into the 3D MOT, we thus use
relatively large 3D-MOT beam waists of $w_{\rm 3D}=12\,$mm with a power of
$P_{\rm 3D}\approx 85$ mW per beam, so as to nearly fulfill the above
estimates of the capture radius and velocity. Note that the atoms fall under
gravity when traveling from the 2D-MOT to the 3D-MOT chamber, which might
compromise their capture. With a horizontal velocity of 11 m/s, the falling
distance is 4.9 mm and is smaller than $R_{\rm cap}$. Furthermore, if the
atoms emerge from the 2D-MOT with an 11 m/s velocity oriented 15 mrad upward,
the fall is suppressed. The relevant parameters for the optimization of the
3D-MOT loading are the detuning of the 3D-MOT beams, $\delta_{\rm 3D}$, and
the magnetic-field gradient, $b^{\prime}_{\rm 3D}$.

In the following, we focus on the isotope 164Dy (atom mass $m=164\,$amu with
amu the atomic mass unit), which has the highest natural abundance. To gain
first insights into the expected performance of our 2D-MOT-based source and
the relevant parameter ranges, we perform Monte Carlo simulations of the 2D-
and 3D-MOT capture processes along the lines of ref. Barbiero  _et al._
(2020), see Appendix B and ref. Gao (2022) for details. Using an oven
temperature of $T=1000\,^{\circ}$C, our simulations indicate that a maximal
2D-MOT flux is achieved for $b^{\prime}_{\rm 2D}=31\,$G/cm and $\delta_{\rm
2D}=-2.1\,\Gamma_{421}$ and is estimated to $\phi_{\rm 2D}\approx 3\times
10^{10}$ atoms/s for 164Dy. The 3D-MOT loading process can also be included in
the simulation. However, due to the large fraction of unloaded trajectories
over the full process and our limited number of simulated trajectories (see
Appendix B), the extracted 3D-MOT loading rates show large fluctuations. In
the simulations, loading into the 3D-MOT is detected when adding a push beam
and 3D-MOT loading rates on the order of $9\times 10^{8}$ atoms/s were
extracted for 164Dy. Based on these simulation results and their parameters,
we started our experimental search. In the experiment, we observed a 3D-MOT
loading even without an additional push beam. This experimental observation
serves as the starting point for our optimization process.

##  III Measurement protocol

In the following, we characterize the performance of our 2D-MOT-based atom
source through the achieved 3D-MOT loading rates and atom numbers. Our
experimental procedure is as follows: In a first step we switch on the 2D-MOT,
optional push, and 3D-MOT beams, the offset-field and gradient coils, with
fixed parameter values for a time $t_{\rm load}$. We then switch off the
2D-MOT and push beams, and hold the cloud for 70 ms without changing any other
parameter values except the offset field 222The experiment is performed in the
presence of a typical remnant bias field of 0.4 G mostly directed along $z$.
Applying a $z$-bias field was found to have only a marginal effect on the
3D-MOT loading performances when varied over a few-G range and will not be
further discussed here. It is typically set to 0.4 G at the loading stage..
Finally, we switch off the 3D-MOT light and gradient field, let the cloud fall
and expand for a short time of flight (TOF) $t_{\rm TOF}$ (typically $t_{\rm
TOF}=5\,$ms), and take an absorption image using horizontally linearly
polarized 421-nm light propagating along the $y$ axis. The absorption signal
is recorded on a CMOS camera (Hamamatsu Orca Spark) via an imaging lens
providing a magnification of $0.438$. The imaging pulses last for $25\,\mu$s.
The imaging beam is operated on resonance, with an intensity below
$0.2\,\rm{mW}/\rm{cm}^{2}$ and a waist of about $10\,$mm. An exemplary image
is shown in the inset of Fig. 1 (b).

For the present characterization, we do not compress the gas by decreasing the
light detuning and intensity after the 3D MOT loading stage. Therefore the
cloud has a relatively high temperature of about $500\,\mu{\rm K}$ while the
remnant field is estimated to be around $0.4\,$G. At these temperature and
magnetic-field values, the atomic population is expected to be depolarized,
and all Zeeman substates occupied. Assuming an equal population of all
substates, the light-scattering cross-section, $\sigma$, is identical for any
light polarization and is given by renormalizing the bare cross-section
$\sigma_{0}=3\lambda_{421}^{2}/2\pi$ by the average of the squares of the
Clebsch-Gordan coefficients for our $J=8\rightarrow J^{\prime}=9$ dipole
transition, over the initial Zeeman substates. This yields
$\sigma=0.3725\sigma_{0}$. We observe that the absorption signal does not
depend on the imaging light polarization, which experimentally supports the
use of $\sigma=0.3725\sigma_{0}$.

By integrating over a region of interest in the absorption images, we extract
the atom number $N$ in the 3D MOT at the end of the sequence. The statistical
uncertainty on $N$ is typically of 5% and we have estimated the systematic
uncertainty to be of about 4% with dominant contributions from the imaging
magnification and unequal populations of Zeeman substates at finite
temperature. To optimize our setup, we mostly rely on a simple scheme that
consists in measuring $N$ for a characteristic $t_{\rm load}=4\,$s, hereafter
referred to as $N_{4{\rm s}}$. To further characterize our system, we record
loading curves of $N$ versus $t_{\rm load}$, which we fit to an exponential
growth function, $N(t)=N_{\rm sat}(1-e^{-t/\tau})$, to extract the 3D-MOT
loading rate, $\phi_{\rm 3D}=N_{\rm sat}/\tau$, and the saturation atom number
$N_{\rm sat}$.

We note that the resonance frequencies for the 421-nm and 626-nm lights are
extracted from the 3D-MOT images. The 421-nm resonance frequency is directly
extracted from the maximum in the absorption signal when scanning the
frequency of the imaging light. The resonance frequency for the 626-nm light
is determined by shining an additional 626-nm beam onto the atoms during the
first milliseconds of their time of flight and monitoring the atom-number
depletion versus the beam frequency via absorption imaging after time of
flight.

##  IV Optimization of the 2D MOT, without push beam

In our setup, the optimization of the 2D-MOT parameters has been performed,
following the first observation of a 3D-MOT loading, in the absence of a push
beam. Here we report on the protocol followed in the optimization process and
the performance achieved in this configuration. We start the optimization
process by setting the 2D-MOT parameters close to the simulated optimum (see
Sec. II), where a 3D-MOT loading is detected ($b^{\prime}_{\rm
2D}=35.4\,$G/cm, $\delta_{\rm 2D}=-2\Gamma_{421}$). We scan the 3D-MOT
gradient and detuning to maximize the 3D-MOT loading, yielding
$b^{\prime}_{\rm 3D}=0.9$ G/cm and $\delta_{\rm 3D}=-55.8(3)\Gamma_{626}$.
Using this 3D-MOT setting, we go on by changing the 2D-MOT parameters as
follows: choosing a given number of magnets per stack, we start by setting the
magnet stacks at their design positions. In this configuration, we optimize
the divergence of the 2D-MOT cooling beam 333We note that the beam divergence
is particularly sensitive because the power balancing between counter-
propagating beams in the 2D MOT is relevant not only at the center but all
along its axis and in particular close to the differential pumping stage. and
its detuning $\delta_{\rm 2D}$ by maximizing $N_{4{\rm s}}$. We then adjust
the position of each magnet stack and iterate on the beam parameters.

For each number of magnets per stack implemented, we identify the best values
of the 2D-MOT parameters and we finally record the full 3D-MOT loading curve
in the optimized configuration. Examples of such loading curves and their fits
(see Sec. III) are shown in Fig. 2 (a) inset. Figure 2 (a) displays the
dependence of the fitted loading rate $\phi_{\rm 3D}$ with the magnetic
gradient $b^{\prime}_{\rm 2D}$. A maximum loading rate of $\phi_{\rm
3D}=2.7(2)\times 10^{7}$ atoms/s and a saturation atom number of $N_{\rm
sat}=9.9(2)\times 10^{7}$ is found for the configuration with 6 magnets per
stack ($b^{\prime}_{\rm 2D}\approx 26.7\,$G/cm) and $\delta_{\rm
2D}=-1.95(1)\Gamma_{421}$. We compare these experimental observations to the
expectations from our Monte Carlo simulations, see Sec. II, and Appendix B for
details. So as to obtain reliable results despite limited sampling, we do not
include the 3D-MOT loading step and simulate only the 2D-MOT loading process.
Figure 2 (a) shows the simulated 2D-MOT flux as a function of $b^{\prime}_{\rm
2D}$ with an oven at $T=800\,^{\circ}$C and at the value of $\delta_{\rm 2D}$
for which this rate is maximal, like in the experimental procedure. We observe
that the optimal $b^{\prime}_{\rm 2D}$ and $\delta_{\rm 2D}$ are comparable
yet respectively slightly smaller and larger in magnitude for the experimental
$\phi_{\rm 3D}$ compared to the simulated $\phi_{\rm 2D}$ (optimum at
$b^{\prime}_{\rm 2D}=32\,$G/cm and $\delta_{\rm 2D}=-2.1\Gamma_{421}$). Both
in experiment and theory, we find that the optimum value of $|\delta_{\rm
2D}|$ slightly varies over the explored $b^{\prime}_{\rm 2D}$-range, by about
$0.2\,\Gamma_{421}$ and $0.7\,\Gamma_{421}$ respectively. The similarity
between the experiment and simulation results is remarkable given the
approximations made in the simulations, see Appendix B. The simulations do not
necessarily provide the quantitative optima for different parameters, but are
a suitable tool to identify the relevant range of parameters, see Sec. II.

Figure 2:  2D-MOT optimization without push beam.(a) Experimental 3D-MOT
loading rate $\phi_{\rm 3D}$ (circles, left axis) and simulated 2D-MOT flux
$\phi_{\rm 2D}$ (green line and dots, right axis) as a function of the 2D-MOT
magnetic gradient $b^{\prime}_{\rm 2D}$. The 3D-MOT parameters were fixed to
$b^{\prime}_{\rm 3D}=0.9$ G/cm, $\delta_{\rm 3D}=-55.8(3)\Gamma_{626}$ and the
2D-MOT configurations (including $\delta_{\rm 2D}$) were individually
optimized, see text. The inset shows the experimental 3D-MOT loading curves
and exponential fits from which the rates are extracted, and the errorbars
show the standard deviation of three experimental runs. (b) Experimental
$\phi_{\rm 3D}$ (blue circles) and simulated $\phi_{\rm 2D}$ (green line and
dots) as a function of the 2D-MOT detuning $\delta_{\rm 2D}$ at the
experimental optimal gradient, $b^{\prime}_{\rm 2D}=26.7\,$G/cm. The 3D-MOT
parameters were fixed to $b^{\prime}_{\rm 3D}=0.42$ G/cm, $\delta_{\rm
3D}=-42.6(3)\Gamma_{626}$. In (a) and (b), the errorbars are the $63\%$
confidence interval from the fit. The shaded area shows the standard deviation
of three simulation runs.

In the optimized magnetic configuration identified above, we further
investigate the influence of the 2D-MOT detuning and record the full 3D-MOT
loading curve for different values of $\delta_{\rm 2D}$. The extracted loading
rates are shown in Fig. 2 (b). Over the investigated $\delta_{\rm 2D}$ range
($1\Gamma_{421}$), the variations of $\phi_{\rm 3D}$ are rather symmetric
around its maximum, and its value changes by less than $50\%$. In Fig. 2 (b),
we also show the simulated 2D-MOT flux at $b^{\prime}_{\rm 2D}=26.7\,$G/cm.
The variations of the simulated $\phi_{\rm 2D}$ match those of the
experimental $\phi_{\rm 3D}$ well. Experimentally, the maximum is found at
$\delta_{\rm 2D}=-1.95(1)\Gamma_{421}$ with a loading rate of $\phi_{\rm
3D}=2.7(1)\times 10^{7}$ atoms/s and a saturation atom number of $N_{\rm
sat}=12.2(1)\times 10^{7}$. We note that different 3D-MOT parameters were used
compared to Fig. 2 (a), which explains the slightly different loading
performances (see also Sec.VI).

##  V Push-beam enhancement

Figure 3:  Push-beam Effect. Enhancement factor in the atom number in 3D-MOT
at $t_{\rm load}=4s$, $N_{4{\rm s}}$, obtained by the addition of a 421-nm (a)
and 626-nm (b) push beam, as a function of the push-beam power $P_{\rm push}$
and detuning $\delta_{\rm push}$. The 3D-MOT parameters were set to
$b^{\prime}_{\rm 3D}=0.42$ G/cm, $\delta_{\rm 3D}=-42.6(3)\Gamma_{626}$. The
reference without push beam is $N_{4{\rm s}}=7.5(3)\times 10^{7}$. $P_{\rm
push}$ is tuned through the amplitude of the driving signal of an acousto-
optic modulator. The $P_{\rm push}$ values are extracted via a rescaling of
the measured (for (a), frequency-dependent) power in front of the entrance
viewport at fixed diffraction efficiency. (c) Velocity-dependent force (1)
exerted by the 421-nm (blue lines) and 626-nm (red lines) push beams for the
characteristic parameters highlighted by triangle (solid lines) and star
(dashed lines) symbols in (a) and (b), respectively. The forces are calculated
assuming no magnetic field. The inset shows a zoom-in in the relevant low-
velocity region. The grey-shaded area illustrates a Gaussian approximation of
the expected velocity distribution of the 2D-MOT atomic beam.

To further increase the loading rate of our 3D MOT, we implement a push-beam
scheme, as typically done in 2D-MOT setups, see e.g. refs. Dieckmann  _et al._
(1998); Tiecke  _et al._ (2009); Lamporesi  _et al._ (2013); Nosske  _et al._
(2017); Barbiero  _et al._ (2020); Dörscher  _et al._ (2013); Vaidya  _et al._
(2015); Saskin  _et al._ (2019); Guo  _et al._ (2023); Kwon  _et al._ (2023):
we additionally shine a beam propagating through the apparatus along $+x$, see
Fig.1 (b). This beam has a frequency close to one of the transitions of Dy
and, via radiation pressure, pushes the atoms from the 2D-MOT to the science
chamber in a velocity-selective way. In the case of Dy, similarly to Yb
Dörscher  _et al._ (2013); Vaidya  _et al._ (2015); Saskin  _et al._ (2019),
two convenient choices are possible: the push beam can be near-resonant either
with the 626-nm intercombination line or with the broad 421-nm one. To
determine which choice is more beneficial, we implemented both sequentially
and compared their experimental achievements. We note that, in both cases, a
beam-walking optimization yields a configuration in which the push beam goes
through the differential pumping stage and is detected on the opposite side of
the science chamber. The push-beam waist is $w_{\rm push}\approx 0.8$ mm, and,
in order to control the effect of the push beam, we vary its power $P_{\rm
push}$, and detuning $\delta_{\rm push}$. In either case (421-nm or 626-nm),
an optimization on the push-beam parameters clearly improves the 3D-MOT
loading.

We perform a systematic experimental study of the 3D-MOT atom number $N_{4{\rm
s}}$ while varying the push-beam parameters in a range where improved 3D-MOT
loading is found. Figures 3 (a) and (b) show the enhancement factor defined by
the ratio of $N_{4{\rm s}}$ obtained in a configuration with and without a
push beam operating at [Fig. 3 (a)] 421 nm and [Fig. 3 (b)] 626 nm. Here we
use $b^{\prime}_{\rm 3D}=0.42$ G/cm, $\delta_{\rm 3D}=-42.6(3)\Gamma_{626}$,
and the reference without push beam is $N_{4{\rm s}}=7.5(3)\times 10^{7}$ [see
also Fig. 4 (a),(b),(c)]. In both Fig. 3 (a) and (b), the enhancement factor
at fixed $P_{\rm push}$ shows a maximum for varying $\delta_{\rm push}$. We
now describe the power dependence of this maximum. In either case (421-nm and
626-nm), the $\delta_{\rm push}$ value at which the maximum is found increases
in absolute value for increasing $P_{\rm push}$. Yet, the power dependence
shows different features in the two cases. In particular, the value of the
maximum enhancement factor itself shows distinct variations with $P_{\rm
push}$. In the case of the 421-nm push beam, this maximum enhancement factor
is nearly independent of $P_{\rm push}$. In contrast, in the 626-nm case, it
shows an overall optimum in power, corresponding to $P_{\rm push}=18$ mW.
Overall, the push beam on the narrow transition is found to outperform the one
on the broad transition and yields the maximal enhancement factor observed of
$3.6(2)$.

These experimental findings can be comprehended from a description of the
push-beam effect, which is rooted in the radiation pressure force,

$\displaystyle\boldsymbol{F}_{\text{push}}({v}_{x})$ $\displaystyle=$
$\displaystyle\hbar
k\boldsymbol{e}_{x}\frac{\Gamma}{2}\frac{s_{0}}{(1+s_{0})}\frac{1}{1+4\frac{\left({\delta_{\rm
push}}-{k}v_{x}\right)^{2}}{(1+s_{0})\Gamma^{2}}}.$ (1)

Here $k$ is the push-beam wavenumber, $s_{0}=2P_{\rm push}/(\pi w_{\rm
push}^{2}I_{\rm sat})$ is the saturation parameter, $\Gamma$ and $I_{\rm sat}$
are the associated transition’s linewidth and saturation intensity. The force
(1) is directed along the push-beam propagation direction of unit vector
$\boldsymbol{e}_{x}$, and its magnitude depends on the atom’s velocity along
$x$, $v_{x}$, through the Doppler effect. More precisely it is a Lorentzian
function of $v_{x}$, of center $\delta_{\rm push}/k$, of width
$\sqrt{1+s_{0}}\,\Gamma/k$, and of amplitude on resonance $\hbar
k\frac{\Gamma}{2}\frac{s_{0}}{(1+s_{0})}$. The effects of the push-beam
parameters are as follows: Varying $\delta_{\rm push}$ changes the velocity
class with which the force is resonant. Changing $P_{\rm push}$ alters $s_{0}$
and has a twofold effect: (i) it scales the resonant amplitude of the force up
by the factor $\frac{s_{0}}{1+s_{0}}$ and (ii) it broadens the range of
velocities addressed by the force by $\sqrt{1+s_{0}}$. The effects (i) and
(ii) dominate at low and high saturation respectively. Note that Eq. (1)
assumes no effect of the magnetic field, which theoretically cancels along the
propagation axis of the push beam.

The drastically different linewidths of the 421-nm and 626-nm transition
result in a different operation of the corresponding push beams. The 421-nm
push beam may generate much larger forces than the 626-nm beam. In particular,
a force magnitude corresponding to the saturated 626-nm case is generated with
a 421-nm beam of saturation parameter as low as $s_{0}=0.0021$. Furthermore,
the bare velocity widths of the force $\Gamma/k$ are also different, equal to
$0.085$ m/s and $13.6\,$m/s for the 626-nm and 421-nm transition respectively.
This is to be compared with the expected spread of velocity distribution of
the 2D MOT of about $4\,$m/s defined by its half width at half maximum.

These different linewidths, together with the powers used in experiment 444For
the 421-nm beam, the range of power extends from $P_{\rm push}\approx 20$ to
$\approx 400\,\mu$W. Lower powers were not used due to the difficulty of
obtaining a stable power value in this range., result in broadly different
force profiles in the optimal push configurations. Figure 3 (c) sketches such
force profiles as expected from Eq. (1). In both the 421-nm and 626-nm cases,
improved loading conditions are achieved for red detunings, $\delta_{\rm
push}<0$, which corresponds to a resonant pushing of atoms traveling against
the push beam propagation ($v_{x}<0$). Yet, the resonant velocity classes are
vastly different with $v_{x}\sim-10$ m/s ($v_{x}\sim-100$ m/s) for the 626-nm
(421-nm) case. Therefore, in the 421-nm case, mostly the detuned ”tail” of the
radiative force is involved in pushing the atoms. On the contrary, the much
weaker force generated by the 626-nm light is used on and close to its
resonance.

The 421-nm and 626-nm push beams also operate at considerably different
saturation parameters of $s_{0}\leq 0.5$ and $3\times 10^{3}<s_{0}<60\times
10^{3}$, respectively. Following the discussion on the effects (i)-(ii), a
change in power thus affects the force profiles of the 421-nm and 626-nm
lights differently. As illustrated in Fig. 3 (c) through two typical
situations of different $P_{\rm push}$, an increase in power for the 421-nm
beam mostly yields a scaling up of the resonant force magnitude. In contrast,
the 626-nm push operates at roughly constant (saturated) resonant force
magnitude and the change in power mainly yields a broadening of the resonance.

Relevant to the push-beam effect on our 2D-MOT are the force profiles in the
small $|v_{x}|$ range encompassing the 2D-MOT velocity distribution, see Fig.
3 (c) inset. The power dependences of the force in this range explain the
behaviors of the maximal enhancement observed in Figs. 3 (a) and 3(b). The
overall shift to larger $|\delta_{\rm push}|$ values for increasing $P_{\rm
push}$ is justified as follows: increasing $P_{\rm push}$ at fixed
$\delta_{\rm push}$ yields a detrimental effect of pushing the atoms with
positive $v_{x}$ too much (via power scaling or broadening of the force) such
that the final $v_{x}$ may exceed the 3D-MOT capture velocity $v_{\rm cap}$.
Instead, shifting $\delta_{\rm push}$ to larger negative values prevents this
effect and additionally yields the benefit of pushing atoms with larger
negative $v_{x}$ back towards the 3D MOT.

This shift has however different impacts in the 421-nm and 626-nm cases due to
their different operation regimes. In the 421-nm case, the small $|v_{x}|$
range corresponds to a far-detuned regime in which the amplitude scaling is
well compensated by a shift of the resonance. Therefore, when increasing
$P_{\rm push}$, $\delta_{\rm push}$ is adjusted such that the force profile in
this $|v_{x}|$ range is kept basically unchanged. The push-beam effect is thus
nearly power-independent and so is the enhancement factor. In the 626-nm case,
instead, the force profile in the optimal-enhancement condition changes with
$P_{\rm push}$ in the small $|v_{x}|$ range, affecting the push efficiency. At
low power, the range of velocity classes addressed by the force is small
compared to the velocity distribution itself, yielding a low push efficiency.
As described above, increasing $P_{\rm push}$, the 626-nm force broadens with
a saturated resonant amplitude. Therefore if one tries to keep a constant push
effect (i.e. force profile) in the small negative $v_{x}$ range while
increasing $P_{\rm push}$, the saturation and power broadening effects imply
an increasing detrimental push effect on the positive $v_{x}$ range. Hence, at
too large power, either the push effect on $v_{x}>0$ or on $v_{x}<0$ is not
optimal, and an optimum efficiency is expected at an intermediate power. At
the optimum $P_{\rm push}=18$ mW, the 626-nm push force has an expected
velocity width of about $7$ m/s, comparable to the expected width of the
2D-MOT distribution. The fact that the 626-nm force profile can be adapted to
strongly push the $v_{x}\lesssim 0$ class with a reduced impact on the
$v_{x}\sim 10\,$m/s one may explain the observed better performance of the
626-nm push beam. We note that a foreseeable way to further improve the
loading of the 3D MOT, instead of relying on the simple power broadening of
the push beam, is to introduce sidebands in the push-beam frequency, see e.g.
Refs Vaidya  _et al._ (2015); Saskin  _et al._ (2019); Maier  _et al._ (2014);
Dreon  _et al._ (2017). This approach could allow tailoring the velocity-
dependent lineshape of the force to effectively address different velocity
classes.

For the case of the 626-nm push beam, we additionally observe a rather
unexpected behavior in Fig. 3 (b): $N_{4{\rm s}}$ shows multiple local maxima
when varying $\delta_{\rm push}$ at fixed $P_{\rm push}$. This effect may be
explained by the possible presence of a remnant magnetic field along the push-
beam propagation axis. Such a magnetic field makes the vectorial nature of the
push-beam transition ($J=8\rightarrow J^{\prime}=9$) become relevant and
modifies the simple picture of Eq. (1) and Fig. 3 (c). In particular, this
yields different resonant conditions for the push-beam light components
driving the $\sigma_{+}$, $\pi$, and $\sigma_{-}$ transitions respectively.
Therefore, the total force profile, given by the sum of these three
contributions, would then present three distinct peaks with different resonant
velocities, and amplitudes and widths set by the light-polarization
composition. The relevance of this effect is supported by an observed change
in the relative strengths of the maxima in $N_{4{\rm s}}$ when changing the
push-beam polarization. We note that the best performance is found with a push
beam of horizontal linear polarization.

Figure 4:  3D-MOT with and without push beam. (a-c) Experimental 3D-MOT atom
number, $N_{4{\rm s}}$, as a function of the 3D-MOT detuning $\delta_{\rm 3D}$
and gradient $b^{\prime}_{\rm 3D}$, in (a) the absence, or the presence of a
(b) 421-nm and (c) 626-nm push beam. The other parameters are set to the
optimal values identified in Figs. 2 and 3. In particular, we set (b)
$\delta_{\rm push}=-8.3\Gamma_{421}$, $P_{\rm push}=0.26$ mW and (c)
$\delta_{\rm push}=-82.3(3)\Gamma_{626}$, $P_{\rm push}=12$ mW. (d) 3D-MOT
detuning, $\delta_{\rm 3D}^{*}$, at which the maximal $N_{4{\rm s}}$ is
achieved at fixed $b^{\prime}_{\rm 3D}$. The values are extracted from the
scans of (a, black square), (b, blue triangle), (c, red circle) and shown as a
function of $b^{\prime}_{\rm 3D}$. (e) Optimal values of $N_{4{\rm s}}$,
$N^{*}$, as a function of $b^{\prime}_{\rm 3D}$. Same code as (d). (f) 3D-MOT
loading rate as a function of the 3D-MOT beam power, $P_{\rm 3D}$, with a
626-nm push beam of $\delta_{\rm push}=-82.3(3)\Gamma_{626}$, $P_{\rm
push}=18$mW, and $b^{\prime}_{\rm 3D}=0.42$ G/cm, $\delta_{\rm
3D}=-42.6(3)\Gamma_{626}$. The inset shows the corresponding full 3D-MOT
loading curve (circles) and their exponential fits (lines). The different
colors correspond to the different $P_{\rm 3D}$ values. In (d-f), the error
bars show the standard deviation of three experimental repetitions.

##  VI 3D-MOT parameters and capture with and without push beam

In a final study, we investigate the optimal 3D-MOT settings and their
variations without, with a 421-nm, and with a 626-nm push beam. We set the
push-beam parameters to the values providing the optimal enhancement factors
in Fig. 3. Figures 4 (a), (b), and (c) show the atom number $N_{4{\rm s}}$,
while scanning the 3D-MOT detuning $\delta_{\rm 3D}$ and magnetic gradient
$b^{\prime}_{\rm 3D}$ in these three push-beam configurations. The three plots
differ in their overall magnitude but show qualitatively similar variations
with $\delta_{\rm 3D}$ and $b^{\prime}_{\rm 3D}$. In all three cases, for each
$b^{\prime}_{\rm 3D}$ value, $N_{4{\rm s}}$ shows a maximum (noted $N^{*}$)
when varying $\delta_{\rm 3D}$. The optimum is found at a negative detuning,
$\delta_{\rm 3D}^{*}$, whose magnitude $|\delta_{\rm 3D}^{*}|$ increases with
$b^{\prime}_{\rm 3D}$. An overall optimum in $N_{4{\rm s}}$ is found at a
finite value of the $(b^{\prime}_{\rm 3D},\delta_{\rm 3D})$ set.

In the following, we compare the variations of $\delta_{\rm 3D}^{*}$ and
$N^{*}$ with $b^{\prime}_{\rm 3D}$ quantitatively in the three aforementioned
configurations in order to understand the capture process. We extract
$\delta_{\rm 3D}^{*}$ from the three sets of experimental data of Fig. 4 (a),
(b), and (c) and display them as a function of $b^{\prime}_{\rm 3D}$ in Fig. 4
(d). We observe that the variations of $\delta_{\rm 3D}^{*}$ versus
$b^{\prime}_{\rm 3D}$ are similar in the presence of 421-nm and 626-nm push
beams but slightly differ from the case without push beam. Furthermore, in all
three cases, $\delta_{\rm 3D}^{*}$ appears to obey a roughly linear dependence
on $b^{\prime}_{\rm 3D}$: $\delta_{\rm 3D}^{*}=\mu_{R}b^{\prime}_{\rm
3D}+\delta_{0}$. Based on a simple theory of the 3D-MOT-capture process
inspired from ref. Dalibard (2015) that we develop in Appendix C (see also
Fig. 1 (d)), we can interpret the linear-dependence parameters in terms of
effective capture parameters (see Eq. (9)). The offset detuning relates to an
effective capture velocity in the limit of small gradients
($b^{\prime}\rightarrow 0$), $v_{0}$, via $\delta_{0}=-k_{626}v_{0}/\sqrt{8}$,
matching the relation expected in optical molasses. The slope relates to an
effective capture radius with large gradients ($b^{\prime}\rightarrow\infty$):
$\mu_{R}\propto R_{\infty}$. A linear fit yields
$v_{0}=[7.8(1),7.5(1),5.6(2)]\,$m/s and $R_{\infty}=[17(1),19.1(7),32(2)]\,$mm
for [with 626-nm push, with 421-nm push, without push]. The increase in
$v_{0}$ and decrease of $R_{\infty}$ found when adding a push beam reveal a
change in the velocity distribution of the atomic beam emerging from the 2D
MOT. It evidences the boost in velocities and the decrease in spreading
provided by the push beams.

Our theory of the 3D-MOT capture process implies variations of the capture
radius and velocity with $b^{\prime}_{\rm 3D}$. In particular, $v_{\rm cap}$
increases linearly with $b^{\prime}_{\rm 3D}$. Ultimately, this limits the
validity of our description to an intermediate gradient range as the capture
parameters are bounded by the geometry and maximal radiative force, see Sec.
II and Appendix C. These variations also enable us to comprehend the changes
with $b^{\prime}_{\rm 3D}$ in the loading efficiency and therefore in $N^{*}$.
Let us first describe the experimental observations. Figure 4 (e) depicts the
variations of $N^{*}$ with $b^{\prime}_{\rm 3D}$, as extracted from Figs. 4
(a), (b), and (c). We observe that $N^{*}$ varies with $b^{\prime}_{\rm 3D}$
following a similar trend between the three configurations: Starting from
small $b^{\prime}_{\rm 3D}$, $N^{*}$ sharply increases and then slowly
decreases when increasing $b^{\prime}_{\rm 3D}$. A maximum of $N^{*}$ is found
at an intermediate $b^{\prime}_{\rm 3D}$ whose value depends on the push-beam
configuration. The optimum $b^{\prime}_{\rm 3D}$ is the lowest when no push
beam is used and takes the value $b^{\prime}_{\rm 3D}=0.31\,$G/cm. The optimum
is shifted to larger values when using a push beam, namely $b^{\prime}_{\rm
3D}=0.42$ G/cm ($0.51$ G/cm) with the 626-nm (421-nm) push beam. We also
observe a steeper decrease of $N^{*}$ at large $b^{\prime}_{\rm 3D}$ with a
push beam compared to the case without push beam. Let us now understand this
behavior based on the capture-process theory. The increase of $N^{*}$ at small
$b^{\prime}_{\rm 3D}$, is justified by the corresponding increase of $v_{\rm
cap}$, enabling to capture a larger fraction of the atomic distribution. The
gain earned by increasing $v_{\rm cap}$ saturates once it encompasses the full
velocity distribution of the atomic beam or reaches the upper bound of $v_{\rm
cap}\lesssim 11\,$m/s imposed by our 3D-MOT configuration as introduced in
Sec. II. For larger values of $b^{\prime}_{\rm 3D}$, a decrease of the capture
efficiency is foreseen since the radiation pressure profile is not anymore
optimized for the low velocities of the atomic jet. The weaker dependence of
$N^{*}$ on $b^{\prime}_{\rm 3D}$ in the absence of a push beam at large
$b^{\prime}_{\rm 3D}$ shall relate to the different velocity distributions in
the atomic jet between the configurations. Using the relation of Appendix C,
we can estimate the capture velocities for the optimal $b^{\prime}_{\rm 3D}$
to $v_{\rm cap}=[10.1(2),10.7(2),8.8(3)]$ m/s for [626-nm push, 421-nm push,
without push]. The different optimal $b^{\prime}_{\rm 3D}$ can therefore be
interpreted as a requirement to increase the capture velocities when
introducing the push beams.

Overall, the largest $N^{*}$ is found with the 626-nm push beam. The optimal
values of the parameters determined within our optimization process are
reported in Table 1. Finally, we study the full loading curve of the 3D MOT in
this identified optimal configuration. Additionally, we investigate the effect
of changing the power of the 3D-MOT beams $P_{\rm 3D}$. The loading curves and
the extracted loading rates are displayed in Fig. 4 (f). At the previously set
value of $P_{\rm 3D}\approx 85$ mW, we measure a loading rate of $\phi_{\rm
3D}=1.10(2)\times 10^{8}$ atoms/s and a saturation number of $N_{\rm
sat}=2.80(2)\times 10^{8}$. By increasing $P_{\rm 3D}$, both the saturation
atom numbers and loading rates first sharply increase and then continue
increasing at a slower rate. With the maximal power accessible in the present
setup of $\approx 130\,$mW per beam, we find a maximal loading rate of
$\phi_{\rm 3D}=1.28(4)\times 10^{8}$ atoms/s and a saturation atom number of
$N_{\rm sat}=3.76(9)\times 10^{8}$. We note that a gain in the loading rate
and saturation number seems still possible by increasing the power yet
further. In this work, we have simply relied on the power broadening of the
3D-MOT radiative force to increase the capture of our narrow-line 3D MOT.
Additional schemes, such as spectral broadening Vaidya  _et al._ (2015);
Saskin  _et al._ (2019); Maier  _et al._ (2014); Dreon  _et al._ (2017) or
angled-slowing beams Seo  _et al._ (2020); Lunden  _et al._ (2020), could be
considered for potential further performance enhancement. We also note that
the available power for the 2D-MOT cooling beam is presently technically
limited, and we expect that a larger available power will enhance the present
performance. The 3D-MOT loading achieved while reducing $P_{\rm 2D}$ is
reported in App. A. Finally, we stress that increasing the oven temperature
drastically increases the loading rate. As declared earlier, we decided to
proceed with the optimization of our setup at a relatively low oven reservoir
temperature for enhancing the lifetime of our source.

|  Parameter |  Value | Unit
---|---|---|---
2D MOT |  $\lambda_{\rm 2D}$ |  421 | nm
$P_{\rm 2D}$ |  430 | mW
$w_{\rm 2D}$ |  16 | mm
$b^{\prime}_{\rm 2D}$ |  26.7 | G/cm
$\delta_{\rm 2D}$ |  -1.95 | $\Gamma_{421}$
|  $s_{\rm 2D}$ |  1.9 |
Push |  $\lambda_{\rm push}$ |  626 | nm
$P_{\rm push}$ |  18 | mW
$w_{\rm push}$ |  0.8 | mm
$\delta_{\rm push}$ |  -82.3 | $\Gamma_{626}$
|  $s_{\rm push}$ |  25000 |
3D MOT (loading) |  $\lambda_{\rm 3D}$ |  626 | nm
$P_{\rm 3D}$ |  85 | mW
$w_{\rm 3D}$ |  12 | mm
$b^{\prime}_{\rm 3D}$ |  0.42 | G/cm
$\delta_{\rm 3D}$ |  -42.6 | $\Gamma_{626}$
|  $s_{\rm 3D}$ |  520 |
Table 1: Values of the different relevant parameters for optimal operation of
our 2D-/3D-MOT scheme. The waists $w_{\rm 2D}$, $w_{\rm push}$, $w_{\rm 3D}$,
and the powers $P_{\rm 2D}$, $P_{\rm 3D}$ were simply set to their values and
not optimized. Other values are the results of our optimization process
described in this paper.

##  VII Conclusions

We have demonstrated the successful operation of a Dy intercombination line 3D
MOT loaded from a 2D MOT working on the broad 421-nm transition. The addition
of a push beam operating close to the intercombination line allows for the
best loading performance, with a more-than-three-fold increase in the loading
rate. We observe loading rates of $\phi_{\rm 3D}>1\times 10^{8}$ atoms/s and a
saturation number of $N_{\rm sat}\approx 3\times 10^{8}$. This is similar or
better compared to other cold-atom Dy setups based on Zeeman slowers despite
the lower oven temperatures at which we operate, see e.g. Youn  _et al._
(2010); Lu  _et al._ (2010); Maier  _et al._ (2014); Dreon  _et al._ (2017);
Ilzhöfer  _et al._ (2018); Lunden  _et al._ (2020). We note that the loading
of the intercombination-line 3D MOT is a promising first step for quantum gas
experiments. In particular, following the loading, a compression step can be
applied in which the power and the detuning absolute value of the 3D-MOT beams
are decreased Maier  _et al._ (2014); Dreon  _et al._ (2017); Ilzhöfer  _et
al._ (2018). By applying such a step to our samples, temperatures of
$15\,\mu$K have been achieved with negligible atom loss.

We note that our setup has several advantages compared to Zeeman-slower-based
ones. These include its compactness (our system is less than 1 m long), its
lower energetic consumption thanks to the use of permanent magnets and lower
oven temperature, as well as the absence of a direct view between the oven and
the center of the science chamber which reduces collisions with hot atoms
without the need of a mechanical shutter and allows for a full optical
switching of the atomic beam compatible with fast-cycling experiments. In
future developments, a glass cell could directly be substituted for the
metallic chamber and therefore allow for even greater optical access and
faster magnetic-field control without e.g. the need for additional transport
of the atomic cloud.

Another interesting advantage of 3D MOTs of heavy atoms working on narrow
transitions is the possibility to remove the MOT beam coming from the top, see
e.g. Ilzhöfer  _et al._ (2018). This setting allows for greater optical
access, and e.g. easier integration of a high-resolution objective. In our
setup, we observed that a 3D-MOT could also be loaded in such a 5-beam
configuration but we did not perform an optimization in this setting yet.

Moreover, we note that we achieved and observed 2D MOTs of other isotopes of
Dy, namely ${}^{162}{\rm Dy}$ and ${}^{163}{\rm Dy}$. We also tried and
successfully loaded a 3D MOT of ${}^{163}{\rm Dy}$. ${}^{161}{\rm Dy}$ could
not be loaded on a first try and we suspect that a repumping frequency should
be added to the 2D-MOT light.

Our novel scheme thus constitutes a very favorable platform on which to build
more complex experiments. To cite only two examples, based on such a 3D MOT,
one could proceed with (i) loading a dipole trap and performing evaporative
cooling to quantum degeneracy, or (ii) loading single atoms in arrays of
optical tweezers. Both platforms are highly promising candidates for quantum
simulation or quantum computation purposes Bloch  _et al._ (2008, 2012); Gross
and Bloch (2017); Kaufman and Ni (2021). Finally, we also note that this
scheme should be readily adaptable to other open-shell lanthanide species such
as Er.

We note that another setup based on a similar 421-nm 2D-MOT loading a 626-nm
3D-MOT of Dy atoms has been developed in the group of I. Ferrier-Barbut (Bloch
et al., in prep). We have greatly benefited from exchanges between our groups.

###### Acknowledgements.

First and foremost, we thank the Heidelberg quantum-gas community for their
constant technical and scientific support along with the numerous fruitful
discussions. This includes S. Jochim, M. Weidemüller, M. Oberthaler, F.
Jendrzejewski, and their groups with a special mention to the HQA for sharing
many of their design thoughts. We thank I. Ferrier-Barbut and his group for
open exchanges and discussions especially during the design process. We thank
M. Barbiero for providing his simulation code. We thank J. Wilson and J.
Thompson for enlightening discussions based on their Yb 2D/3D-MOT setup. We
further thank A. Patscheider, D. Petter, G. Natale, P. Ilzhöfer, E. Kirilov,
J. Beugnon, J. Dalibard, and C. Weitenberg for numerous early-stage
discussions and technical advice. We thank T. Yefsah and M. Rabinovic for
sharing technical designs. We thank L. Hoenen, P. Holzenkamp, V. Salazar
Silva, B. Bader for their technical assistance. This work is funded by the
European Research Council (ERC) under the European Union’s Horizon Europe
research and innovation program under grant number 101040688 (project 2DDip),
and by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)
through project-ID 273811115 (SFB1225 ISOQUANT) and under Germany’s Excellence
Strategy EXC2181/1-390900948 (the Heidelberg Excellence Cluster STRUCTURES).
Views and opinions expressed are however those of the authors only and do not
necessarily reflect those of the European Union or the European Research
Council. Neither the European Union nor the granting authority can be held
responsible for them. J. G. acknowledges support from the International Max
Planck Research School for Quantum Dynamics (IMPRS-QD).

$\dagger$ These authors contributed equally.

$\star$ Correspondence and requests for materials should be addressed to
<EMAIL_ADDRESS>
## References

  * Chomaz  _et al._ (2022) Lauriane Chomaz, Igor Ferrier-Barbut, Francesca Ferlaino, Bruno Laburthe-Tolra, Benjamin L Lev, and Tilman Pfau, “Dipolar physics: a review of experiments with magnetic quantum gases,” Reports on Progress in Physics 86, 026401 (2022).
  * McClelland and Hanssen (2006) J. J. McClelland and J. L. Hanssen, “Laser cooling without repumping: A magneto-optical trap for erbium atoms,” Phys. Rev. Lett. 96, 143005 (2006).
  * Berglund  _et al._ (2008) Andrew J. Berglund, James L. Hanssen, and Jabez J. McClelland, “Narrow-line magneto-optical cooling and trapping of strongly magnetic atoms,” Phys. Rev. Lett. 100, 113002 (2008).
  * Youn  _et al._ (2010) Seo Ho Youn, Mingwu Lu, Ushnish Ray, and Benjamin L. Lev, “Dysprosium magneto-optical traps,” Phys. Rev. A 82, 043425 (2010).
  * Lu  _et al._ (2010) Mingwu Lu, Seo Ho Youn, and Benjamin L. Lev, “Trapping ultracold dysprosium: A highly magnetic gas for dipolar physics,” Phys. Rev. Lett. 104, 063001 (2010).
  * Sukachev  _et al._ (2010) D. Sukachev, A. Sokolov, K. Chebakov, A. Akimov, S. Kanorsky, N. Kolachevsky, and V. Sorokin, “Magneto-optical trap for thulium atoms,” Phys. Rev. A 82, 011405(R) (2010).
  * Frisch  _et al._ (2012) A. Frisch, K. Aikawa, M. Mark, A. Rietzler, J. Schindler, E. Zupanič, R. Grimm, and F. Ferlaino, “Narrow-line magneto-optical trap for erbium,” Phys. Rev. A 85, 051401(R) (2012).
  * Maier  _et al._ (2014) T. Maier, H. Kadau, M. Schmitt, A. Griesmaier, and T. Pfau, “Narrow-line magneto-optical trap for dysprosium atoms,” Opt. Lett. 39, 3138–3141 (2014).
  * Miao  _et al._ (2014) J. Miao, J. Hostetter, G. Stratis, and M. Saffman, “Magneto-optical trapping of holmium atoms,” Phys. Rev. A 89, 041401(R) (2014).
  * Dreon  _et al._ (2017) Davide Dreon, Leonid A Sidorenkov, Chayma Bouazza, Wilfried Maineult, Jean Dalibard, and Sylvain Nascimbene, “Optical cooling and trapping of highly magnetic atoms: the benefits of a spontaneous spin polarization,” Journal of Physics B: Atomic, Molecular and Optical Physics 50, 065005 (2017).
  * Cojocaru  _et al._ (2017) I. S. Cojocaru, S. V. Pyatchenkov, S. A. Snigirev, I. A. Luchnikov, E. S. Kalganova, G. A. Vishnyakova, D. N. Kublikova, V. S. Bushmakin, E. T. Davletov, V. V. Tsyganok, O. V. Belyaeva, A. Khoroshilov, V. N. Sorokin, D. D. Sukachev, and A. V. Akimov, “Light-assisted collisions in ultracold tm atoms,” Phys. Rev. A 95, 012706 (2017).
  * Ilzhöfer  _et al._ (2018) P. Ilzhöfer, G. Durastante, A. Patscheider, A. Trautmann, M. J. Mark, and F. Ferlaino, “Two-species five-beam magneto-optical trap for erbium and dysprosium,” Phys. Rev. A 97, 023633 (2018).
  * Inoue  _et al._ (2018) Ryotaro Inoue, Yuki Miyazawa, and Mikio Kozuma, “Magneto-optical trapping of optically pumped metastable europium,” Phys. Rev. A 97, 061607(R) (2018).
  * Phelps  _et al._ (2020) Gregory A Phelps, Anne Hébert, Aaron Krahn, Susannah Dickerson, Furkan Öztürk, Sepehr Ebadi, Lin Su, and Markus Greiner, “Sub-second production of a quantum degenerate gas,” arXiv preprint 2007.10807 (2020).
  * Lunden  _et al._ (2020) William Lunden, Li Du, Michael Cantara, Pierre Barral, Alan O. Jamison, and Wolfgang Ketterle, “Enhancing the capture velocity of a dy magneto-optical trap with two-stage slowing,” Phys. Rev. A 101, 063403 (2020).
  * Lu  _et al._ (2011) Mingwu Lu, Nathaniel Q. Burdick, Seo Ho Youn, and Benjamin L. Lev, “Strongly dipolar bose-einstein condensate of dysprosium,” Phys. Rev. Lett. 107, 190401 (2011).
  * Lu  _et al._ (2012) Mingwu Lu, Nathaniel Q. Burdick, and Benjamin L. Lev, “Quantum degenerate dipolar fermi gas,” Phys. Rev. Lett. 108, 215301 (2012).
  * Aikawa  _et al._ (2012) K. Aikawa, A. Frisch, M. Mark, S. Baier, A. Rietzler, R. Grimm, and F. Ferlaino, “Bose-einstein condensation of erbium,” Phys. Rev. Lett. 108, 210401 (2012).
  * Aikawa  _et al._ (2014) K. Aikawa, A. Frisch, M. Mark, S. Baier, R. Grimm, and F. Ferlaino, “Reaching fermi degeneracy via universal dipolar scattering,” Phys. Rev. Lett. 112, 010404 (2014).
  * Davletov  _et al._ (2020) E. T. Davletov, V. V. Tsyganok, V. A. Khlebnikov, D. A. Pershin, D. V. Shaykin, and A. V. Akimov, “Machine learning for achieving bose-einstein condensation of thulium atoms,” Phys. Rev. A 102, 011302(R) (2020).
  * Dalibard (2015) Jean Dalibard, “Une brève histoire des atomes froids,”  Collège de France (2015).
  * Katori  _et al._ (1999) Hidetoshi Katori, Tetsuya Ido, Yoshitomo Isoya, and Makoto Kuwata-Gonokami, “Magneto-optical trapping and cooling of strontium atoms down to the photon recoil temperature,” Phys. Rev. Lett. 82, 1116–1119 (1999).
  * Kuwamoto  _et al._ (1999) T. Kuwamoto, K. Honda, Y. Takahashi, and T. Yabuzaki, “Magneto-optical trapping of yb atoms using an intercombination transition,” Phys. Rev. A 60, R745–R748 (1999).
  * Dieckmann  _et al._ (1998) K. Dieckmann, R. J. C. Spreeuw, M. Weidemüller, and J. T. M. Walraven, “Two-dimensional magneto-optical trap as a source of slow atoms,” Phys. Rev. A 58, 3891–3895 (1998).
  * Tiecke  _et al._ (2009) T. G. Tiecke, S. D. Gensemer, A. Ludewig, and J. T. M. Walraven, “High-flux two-dimensional magneto-optical-trap source for cold lithium atoms,” Phys. Rev. A 80, 013409 (2009).
  * Lamporesi  _et al._ (2013) G. Lamporesi, S. Donadello, S. Serafini, and G. Ferrari, “Compact high-flux source of cold sodium atoms,” Review of Scientific Instruments 84 (2013), 063102.
  * Nosske  _et al._ (2017) Ingo Nosske, Luc Couturier, Fachao Hu, Canzhu Tan, Chang Qiao, Jan Blume, Y. H. Jiang, Peng Chen, and Matthias Weidemüller, “Two-dimensional magneto-optical trap as a source for cold strontium atoms,” Phys. Rev. A 96, 053415 (2017).
  * Barbiero  _et al._ (2020) Matteo Barbiero, Marco G. Tarallo, Davide Calonico, Filippo Levi, Giacomo Lamporesi, and Gabriele Ferrari, “Sideband-enhanced cold atomic source for optical clocks,” Phys. Rev. Applied 13, 014013 (2020).
  * Vaidya  _et al._ (2015) V. D. Vaidya, J. Tiamsuphat, S. L. Rolston, and J. V. Porto, “Degenerate bose-fermi mixtures of rubidium and ytterbium,” Phys. Rev. A 92, 043604 (2015).
  * Saskin  _et al._ (2019) S. Saskin, J. T. Wilson, B. Grinkemeyer, and J. D. Thompson, “Narrow-line cooling and imaging of ytterbium atoms in an optical tweezer array,” Phys. Rev. Lett. 122, 143002 (2019).
  * Guo  _et al._ (2023) Changlei Guo, Valentin Cambier, James Calvert, Maxime Favier, Manuel Andia, Luigi de Sarlo, and Sébastien Bize, “Exploiting the two-dimensional magneto-optical trapping of ${}^{199}\mathrm{Hg}$ for a mercury optical lattice clock,” Phys. Rev. A 107, 033116 (2023).
  * Martin  _et al._ (1978) W C Martin, R Zalubas, and L Hagan, “Atomic energy levels - the rare earth elements. (the spectra of lanthanum, cerium, praseodymium, neodymium, promethium, samarium, europium, gadolinium, terbium, dysprosium, holmium, erbium, thulium, ytterbium, and lutetium). [66 atoms and ions],” Technical Report (United States) 01, (1978).
  * Black (2001) Eric D. Black, “An introduction to Pound–Drever–Hall laser frequency stabilization,” American Journal of Physics 69, 79–87 (2001).
  * Note (1) In the main-text measurement, we use the maximal power available while the 2D-MOT power dependence is reported in App. A.
  * Gao (2022) Jianshun Gao, “A first two-dimensional magneto-optical trap for dysprosium,” PhD Admission Thesis, Heidelberg University (2022).
  * Note (2) The experiment is performed in the presence of a typical remnant bias field of 0.4 G mostly directed along $z$. Applying a $z$-bias field was found to have only a marginal effect on the 3D-MOT loading performances when varied over a few-G range and will not be further discussed here. It is typically set to 0.4 G at the loading stage.
  * Note (3) We note that the beam divergence is particularly sensitive because the power balancing between counter-propagating beams in the 2D MOT is relevant not only at the center but all along its axis and in particular close to the differential pumping stage.
  * Dörscher  _et al._ (2013) Sören Dörscher, Alexander Thobe, Bastian Hundt, André Kochanke, Rodolphe Le Targat, Patrick Windpassinger, Christoph Becker, and Klaus Sengstock, “Creation of quantum-degenerate gases of ytterbium in a compact 2D-/3D-magneto-optical trap setup,” Review of Scientific Instruments 84 (2013), 043109.
  * Kwon  _et al._ (2023) Minho Kwon, Aaron Holman, Quan Gan, Chun-Wei Liu, Matthew Molinelli, Ian Stevenson, and Sebastian Will, “Jet-loaded cold atomic beam source for strontium,” Review of Scientific Instruments 94, 013202 (2023).
  * Note (4) For the 421-nm beam, the range of power extends from $P_{\rm push}\approx 20$ to $\approx 400\,\mu$W. Lower powers were not used due to the difficulty of obtaining a stable power value in this range.
  * Seo  _et al._ (2020) Bojeong Seo, Peng Chen, Ziting Chen, Weijun Yuan, Mingchen Huang, Shengwang Du, and Gyu-Boong Jo, “Efficient production of a narrow-line erbium magneto-optical trap with two-stage slowing,” Phys. Rev. A 102, 013319 (2020).
  * Bloch  _et al._ (2008) Immanuel Bloch, Jean Dalibard, and Wilhelm Zwerger, “Many-body physics with ultracold gases,” Rev. Mod. Phys. 80, 885–964 (2008).
  * Bloch  _et al._ (2012) Immanuel Bloch, Jean Dalibard, and Sylvain Nascimbene, “Quantum simulations with ultracold quantum gases,” Nature Physics 8, 267–276 (2012).
  * Gross and Bloch (2017) Christian Gross and Immanuel Bloch, “Quantum simulations with ultracold atoms in optical lattices,” Science 357, 995–1001 (2017).
  * Kaufman and Ni (2021) Adam M Kaufman and Kang-Kuen Ni, “Quantum science with optical tweezer arrays of ultracold atoms and molecules,” Nature Physics 17, 1324–1333 (2021).
  * Giordmaine and Wang (1960) J. A. Giordmaine and T. C. Wang, “Molecular beam formation by long parallel tubes,” Journal of Applied Physics 31, 463–471 (1960).
  * Olander and Kruger (2003) Donald R. Olander and Valerie Kruger, “Molecular Beam Sources Fabricated from Multichannel Arrays. III. The Exit Density Problem,” Journal of Applied Physics 41, 2769–2776 (2003).
  * Wohlleben  _et al._ (2001) W Wohlleben, Frédéric Chevy, K Madison, and Jean Dalibard, “An atom faucet,” The European Physical Journal D-Atomic, Molecular, Optical and Plasma Physics 15, 237–244 (2001).
  * Scherer (2010) Philipp OJ Scherer,  _Computational Physics_ (Springer, 2010) pp. 139–140.
  * Lett  _et al._ (1989) P. D. Lett, W. D. Phillips, S. L. Rolston, C. E. Tanner, R. N. Watts, and C. I. Westbrook, “Optical molasses,” J. Opt. Soc. Am. B 6, 2084–2107 (1989).
  * Kohel  _et al._ (2003) James M. Kohel, Jaime Ramirez-Serrano, Robert J. Thompson, Lute Maleki, Joshua L. Bliss, and Kenneth G. Libbrecht, “Generation of an intense cold-atom beam from a pyramidal magneto-optical trap: experiment and simulation,” J. Opt. Soc. Am. B 20, 1161–1168 (2003).

##  Appendix A 3D-MOT loading as a function of the 2D-MOT power

In the main text, we use the maximal power available of 421-nm light for our
2D-MOT cooling beam. This is done in order to maximize the loading performance
while we expect the saturated regime not to be achieved for the 2D-MOT
performances on the broad transition, see also Tiecke  _et al._ (2009);
Lamporesi  _et al._ (2013); Nosske  _et al._ (2017); Barbiero  _et al._
(2020). In this appendix, we report on the change in the loading performances
with 2D-MOT cooling-beam power, $P_{\rm 2D}$. Fig. AS1 shows the relative
variations in the loading rate $\phi_{\rm 3D}$ and saturation atom number
$N_{\rm sat}$ with power together with their fits to exponential-growth
function $A\left(1-e^{-(P_{\rm 2D}-P_{0})/P_{\rm sat}}\right)$. The fit
results show that the saturation regime is not achieved for the available
power and increased performances are expected for increased power.

Figure AS1:  2D-MOT power dependence. Relative values of the experimental
3D-MOT loading rate $\phi_{\rm 3D}$ and saturation atom numbers $N_{\rm sat}$
as a function of the 2D-MOT power $P_{\rm 2D}$. The values at $P_{\rm
2D}=400$mW are used as references. The parameters are $\delta_{\rm
push}=-46.54(13)\,\Gamma_{626}$, $P_{\rm push}=12.5$ mW, and $b^{\prime}_{\rm
3D}=0.57$ G/cm, $\delta_{\rm 3D}=-39.19(13)\Gamma_{626}$ and $\phi_{\rm
3D}^{(400)}=8.1(2)\times 10^{7}$, $N_{\rm sat}^{(400)}=2.28(1)\times 10^{8}$.
The errorbars are the $63\%$ confidence interval from the fit. The lines are
exponential growth fits yielding $A=[3.84,1.42]$ $P_{0}=[117,106]$ mW and
$P_{\rm sat}=[923,234]$ mW, respectively.

##  Appendix B 2D-MOT Simulation

To estimate the performances of our 2D-MOT scheme as a function of the values
of its parameters we performed Monte Carlo simulations of the trajectories of
the atoms within the radiative-pressure force field, closely following the
line of ref. Barbiero  _et al._ (2020). Here we review the essential elements
of the simulation, for more details see Gao (2022).

The simulation starts at the oven output at time $t=0$. Here, the positions
and velocities of $N_{\rm sim}$ (typically $N_{\rm sim}=50000$) atoms are
sampled assuming a uniform spatial distribution over the aperture and a
thermal ($T=800\,^{\circ}$C) velocity distribution emerging from a tubed
aperture. We use the transition regime approximation to describe its angular
distribution Giordmaine and Wang (1960); Olander and Kruger (2003). The
distribution is cut to $7.5^{\circ}$ and only velocity $v<200\,$m/s are
simulated.

At times $t>0$, the atoms are submitted to the force field that results from
the four passes of the 2D-MOT cooling beam via:

$\boldsymbol{F}_{\text{2D}}(\boldsymbol{r},\boldsymbol{v})=\sum_{i=1}^{4}\hbar\boldsymbol{k}_{i}\frac{\Gamma_{421}}{2}\frac{s_{i}(\boldsymbol{r},\boldsymbol{v})}{1+\sum_{i}s_{i}(\boldsymbol{r},\boldsymbol{v})},$
(2)

where $\boldsymbol{k}_{i}$ is wavevector of the beam at pass $i$
($\boldsymbol{k}_{i}=k_{421}\frac{\pm\boldsymbol{e}_{y}\pm\boldsymbol{e}_{z}}{\sqrt{2}}$
with $\boldsymbol{e}_{x,y,z}$ the unit vector in direction $x,y,z$) and
$s_{i}$ is given by

$s_{i}(\boldsymbol{r},\boldsymbol{v})=\frac{s_{i,0}(\boldsymbol{r})}{1+4\left(\delta_{\boldsymbol{r},\boldsymbol{v}}^{(i)}/\Gamma_{421}\right)^{2}},$
(3)

with $s_{i,0}(\boldsymbol{r})=I_{i}(\boldsymbol{r})/I_{\rm sat}^{421}$ the
local saturation parameter, $I_{i}(\boldsymbol{r})=2P_{\rm 2D}/\pi w_{\rm
2D}^{2}\cdot\exp{(-2r^{2}/w_{\rm 2D}^{2})}$ the local intensity Wohlleben  _et
al._ (2001). The local detuning $\delta_{\boldsymbol{r},\boldsymbol{v}}^{(i)}$
in beam-pass $i$ is position- and velocity-dependent due to the Zeeman effect
in the 2D-MOT magnetic field $\boldsymbol{B}_{\rm 2D}(\boldsymbol{r})$ and the
Doppler effect:

$\delta_{\boldsymbol{r},\boldsymbol{v}}^{(i)}=\delta_{\rm
2D}-\boldsymbol{k}_{i}\cdot\boldsymbol{v}+\frac{\mu_{B}}{\hbar}g_{J^{\prime}}^{421}\Delta
m_{i}|\boldsymbol{B}(\boldsymbol{r})|,$ (4)

where $\Delta m_{i}=\pm 1$ depending if the beam pass $i$ drives the
$\sigma_{\pm}$ transition (i.e. defined by the sign of
$\boldsymbol{B}\cdot\boldsymbol{k}_{i}$). We use an ideal magnetic environment
with $\boldsymbol{B}_{\rm 2D}(\boldsymbol{r})=b^{\prime}_{\rm
2D}(y\boldsymbol{e}_{z}+z\boldsymbol{e}_{y})$.

The position $\boldsymbol{r}$ and velocity $\boldsymbol{v}$ of each atom are
computed in an iterative manner using a four-step Runge-Kutta method with a
time step $\Delta t=50$ $\rm{\mu s}$ Scherer (2010). At each time step,
Newton’s second law is applied with Eq. (2) providing the instantaneous
acceleration. Besides Eq. (2), we also include the effect of spontaneous
emission approximately by an additional random force over the time step
$\Delta t$:

$\boldsymbol{F}_{\text{spont}}(\boldsymbol{r},\boldsymbol{v})=\hbar
k_{421}\sqrt{\frac{\Gamma_{421}}{2}\frac{\sum_{i}s_{i}(\boldsymbol{r},\boldsymbol{v})}{(1+\sum_{i}s_{i}(\boldsymbol{r},\boldsymbol{v}))}\Delta
t}\hat{\boldsymbol{e}}_{\text{spont}},$ (5)

where $\hat{\boldsymbol{e}}_{\text{spont}}$ is a random unit vector, which we
randomly draw at each time step of the simulation Lett  _et al._ (1989); Kohel
_et al._ (2003). We note that Eqs. (2-5) approximate the complex
$J=8\rightarrow J^{\prime}=9$ structure of Dy to a simpler $J=0\rightarrow
J^{\prime}=1$ structure. This accounts for the light polarization effect but
neglects the spinful character of the ground state. Furthermore, the non-
closed character of the cooling transition is omitted Lu  _et al._ (2010).

We stop the simulation after a total time of $10$ ms and count the number of
atoms $N_{\rm cap}$ whose trajectories are within $25\,$mrad solid angle from
the $x>0$-axis. We then extract the 2D-MOT flux via $\phi_{\rm 2D}=\phi_{\rm
sim}(T)N_{\rm cap}/N_{\rm sim}$ where $\phi_{\rm sim}(T)$ is the flux of atoms
exiting the oven at the temperature $T$ accounting for the limits set in angle
and magnitude for the velocities of the simulated atoms. For
$T=800\,^{\circ}$C ($T=1000\,^{\circ}$C), $\phi_{\rm sim}(T)=1.5\times
10^{10}$ atoms/s ($\phi_{\rm sim}(T)=1.9\times 10^{12}$ atoms/s) for 164Dy.
The actual oven temperature certainly lies in between these two values due to
the higher hot lip temperature also partly heating the reservoir. The maximum
flux $\phi_{\rm 2D}$ found in Fig. 2, using $T=800\,^{\circ}$C, is estimated
in absolute value to $\phi_{\rm 2D}=3\times 10^{8}$ 164Dy atoms/s while using
$T=1000\,^{\circ}$C, we find a maximum $\phi_{\rm 2D}=3.5\times 10^{10}$ 164Dy
atoms/s. Therefore the fraction of loaded trajectories in our simulation is on
the order of $2\%$ and very slightly larger at lower oven temperatures. We
also note that the total oven flux is about three orders of magnitude larger
than $\phi_{\rm sim}$.

Simulations further including the push beam and 3D-MOT beams and counting the
number of atoms loaded in the 3D MOT were also performed. In this case, the
force (1) is added to the 2D-MOT force, and another force similar to (2) but
summing over the six 3D-MOT beam passes and centered at $x_{0}=+347$ mm is
included. In this case, we count the number of atoms trapped within a 10 mm
radius around $x_{0}$ while extending the simulation time to $260$ ms. We then
calculate the simulated $\phi_{\rm 3D}$ similarly to $\phi_{\rm 2D}$. Using a
push beam we find that the fraction of loaded trajectories is again reduced by
a factor of $\sim 20$ to $30$ compared to the 2D-MOT results, yielding very
weak signals, $N_{\rm cap}$ of the order of ten with the used $N_{\rm
sim}=50000$. This weak signal precludes reliable simulation results for the
3D-MOT loading as highlighted in Secs. II and IV.

##  Appendix C Basic theory of the 3D-MOT capture process

We consider the process in which a 3D MOT with a 6-beam geometry similar to
the sketch of Fig. 1 (d) captures atoms from an atomic beam propagating along
$+x$ axis. A simple understanding can be drawn, following the lines of ref.
Dalibard (2015), based on the phase-space-dependent detuning of the cooling
beams. Accounting for our beam geometry, the detuning of an atom of position
$x^{\prime}$ and velocity $v_{x}$ writes:

$\delta_{x^{\prime},v_{x}}^{\pm}=\delta_{\rm
3D}\pm\left(\frac{k_{626}v_{x}}{\sqrt{2}}+\frac{\mu_{626}^{\prime}b^{\prime}_{\rm
3D}x^{\prime}}{\hbar}\right),$ (6)

where $+$ ($-$) holds for the two cooling beams propagating against (along)
the atomic beam. Here, $x^{\prime}$ denotes the position of the atoms with
respect to the science chamber center ($x^{\prime}=x-x0$),
$\mu^{\prime}_{626}=g_{J^{\prime}}^{626}\mu_{B}$ is the 626nm-transition
excited state differential magnetic moment and $\mu_{B}$ the Bohr magneton.
The capturing process mostly relies on the effect of the two beams propagating
against the atomic beam since it can be seen from Eq. (6) that they are
brought closer to resonance for $v_{x}>0$ (i.e. as in the atomic beam), or for
$x^{\prime}>0$ (i.e. for atoms that would go beyond the trap center).
Hereafter, our capturing model neglects the effect of the two other beams and
simply uses $\delta_{x^{\prime},v_{x}}^{+}$.

We define $v_{\rm cap}$ as the maximal velocity of the atoms that can be
stopped within the capture region of radius $R_{\rm cap}$. Here, we assume
that both $v_{\rm cap}$ and $R_{\rm cap}$ depend on $b^{\prime}_{\rm 3D}$.
Following ref. Dalibard (2015), we then formulate that the optimal capture
configuration is achieved by adjusting $\delta_{\rm 3D}$ and $b^{\prime}_{\rm
3D}$ such that the atoms with a velocity $v_{x}=v_{\rm cap}(b^{\prime}_{\rm
3D})$ are on resonance at the entrance and exit of the capture region (see
highlighted points in Fig. 1 (d)). This yields $\delta_{x^{\prime}=-R_{\rm
cap},v_{x}=v_{\rm cap}(b^{\prime}_{\rm 3D})}^{+}=0$ at the entrance, and
$\delta_{x^{\prime}=R_{\rm cap},v_{x}=0}^{+}=0$ at the exit. This ensures a
maximal radiative force for this class of atoms. The entrance and exit
conditions respectively yield:

$\displaystyle\delta_{\rm 3D}^{*}(b^{\prime}_{\rm 3D})$ $\displaystyle\approx$
$\displaystyle\frac{\mu_{626}^{\prime}b^{\prime}_{\rm 3D}R_{\rm
cap}(b^{\prime}_{\rm 3D})}{\hbar}-\frac{k_{626}v_{\rm cap}(b^{\prime}_{\rm
3D})}{\sqrt{2}},$ (7) $\displaystyle\delta_{\rm 3D}^{*}(b^{\prime}_{\rm 3D})$
$\displaystyle\approx$ $\displaystyle-\frac{\mu_{626}^{\prime}b^{\prime}_{\rm
3D}R_{\rm cap}(b^{\prime}_{\rm 3D})}{\hbar}.$ (8)

where $\delta_{\rm 3D}^{*}(b^{\prime}_{\rm 3D})$ is the value of the detuning
providing the optimal capture configuration for the gradient $b^{\prime}_{\rm
3D}$. Combining these relations yields

$\delta_{\rm 3D}^{*}(b^{\prime}_{\rm
3D})\approx-\frac{\mu_{626}^{\prime}b^{\prime}_{\rm 3D}R_{\rm
cap}}{\hbar}\approx-\frac{k_{626}v_{\rm cap}}{\sqrt{8}}.$ (9)

The approximate relation (9) can be combined with the linear dependence
observed experimentally, $\delta_{\rm 3D}^{*}(b^{\prime}_{\rm
3D})=\mu_{R}b^{\prime}_{\rm 3D}+\delta_{0}$, and allows giving a simple
interpretation of the parameters $\mu_{R}$ and $\delta_{0}$. The offset value
at $b^{\prime}_{\rm 3D}=0$, $\delta_{0}=\delta_{\rm 3D}^{*}(b^{\prime}_{\rm
3D}=0)$ can be written as $\delta_{0}=-k_{626}v_{0}/\sqrt{8}$ with $v_{0}$ the
effective capture velocity at small gradients $b^{\prime}_{\rm 3D}\rightarrow
0$, matching the formula expected in optical molasses Dalibard (2015). The
scaling factor $\mu_{R}$ can instead be written as
$\mu_{R}=-\mu_{626}^{\prime}R_{\infty}/\hbar$ with $R_{\infty}$ the effective
capture radius of the 3D-MOT with large gradient, $b^{\prime}_{\rm
3D}\rightarrow\infty$.

The linear relation can also be combined with (9) to extract the dependence of
the capture velocity and radius on $b^{\prime}_{\rm 3D}$:

$\displaystyle v_{\rm cap}(b^{\prime}_{\rm 3D})$ $\displaystyle=$
$\displaystyle v_{0}+\frac{\sqrt{8}\mu_{626}^{\prime}b^{\prime}_{\rm
3D}R_{\infty}}{\hbar k_{626}},$ (10) $\displaystyle R_{\rm
cap}(b^{\prime}_{\rm 3D})$ $\displaystyle=$ $\displaystyle
R_{\infty}+\frac{\hbar k_{626}v_{0}}{\sqrt{8}\mu_{626}^{\prime}b^{\prime}_{\rm
3D}}.$ (11)

We find that $v_{\rm cap}$ increases and $R_{\rm cap}$ decreases with the MOT
gradient, which matches intuitively expected behaviors. We note that both
relations present nonphysical divergences, either at large gradients for
$v_{\rm cap}$ or small gradients for $R_{\rm cap}$, which evidence the
limitations of our simple model. More precisely, the value of $R_{\rm cap}$
and $v_{\rm cap}$ are ultimately bounded by our beam geometry and the maximal
force that can be applied, as discussed in Sec. II.