Effect of ion-trap parameters on energy distributions of ultra-cold atom–ion mixtures

Co-trapping ultra-cold atoms and cold ions offers new possibilities for exploring low-temperature collisions, that include phenomena such as s-wave scattering, Feshbach resonances [1], shape-resonances [2] and the creation of molecular ions [3]. In addition, it is also a promising platform for performing quantum computations [4], quantum simulations [5] and for studying out-of-equilibrium dynamics [6]. In the last decade these hybrid systems were realized in several experiments, for reviews see [7–10]. The interaction of trapped-ions with ultra-cold thermal clouds [11] and quantum degenerate gases of neutral atoms [12] was studied. Nonetheless, most of the experiments were limited to atom–ion interaction energy which is greater than the energy scale of the quantum phenomena mentioned above. Only recently, collisions at this energy scale were observed in a system with a heavy ion and light atoms [13].

The interaction energy limitation arises from the time-dependent potential of the ion Paul trap [14]. The interruption of the motion, by a collision, in this time-dependent potential can lead to injection of energy into the system. As a result, the system exhibits non-equilibrium dynamics. Already in 1968, it was observed by Major and Dehmelt [15] that collisions of ions with heavy atoms in a Paul trap lead to exponential heating of the ion and subsequently to its loss. Forty years later, DeVoe [16] demonstrated numerically that a single collision cannot cause this enormous heating effect, but it is rather caused by a sequence of collisions which occur at a certain phase of the oscillating trapping potential. He also showed that these consecutive collisions lead to a power-law energy distribution which is not thermal, as one would expect from a thermalization process. The Tsallis distribution [17], originally proposed as a generalization of the Boltzmann–Gibbs statistical mechanics to non-extensive systems, was proposed to describe this energy distribution [16]. In addition to a characteristic temperature parameter, this distribution also has a parameter describing its power-law tail. Only recently it was shown, by the formalism of super-statistics, that this distribution indeed arises in the limit of multiple collisions [18]. The power-law tail of the distribution depends on the atom–ion mass ratio as well as the specific Paul trap parameters [19–24].

The ion's energy distribution has a characteristic energy scale which can have various sources. If the atoms, colliding with the ion, are at a high temperature, their temperature will determine this energy scale. However, if the colliding atoms are at zero temperature, the energy of the ion can be lost in a single collision if the collision occurs at the center of the trap, ideally leading to a zero energy steady-state for the ion as well. However, it was shown [12, 19, 25] that the mean energy of the ion is typically much higher than that of the atoms by more than an order of magnitude. This additional energy scale arises from the fact that static stray electric fields can move the ion equilibrium position such that it will experience non-vanishing rf fields. This effect is called excess micromotion (EMM). EMM can be compensated by applying an external static electric field that moves the ion to the null point of the rf field [26]. It was Cetina et al, in their seminal paper [27], which realized that even without EMM and for an ion in the ground state, a single collision with an atom at zero temperature can still increase significantly the ion's energy. This effect occurs due to the attractive polarization potential between the ion and the atoms that causes the first collision to happen far from the equilibrium position of the ion, and hence at a region with non-vanishing time-dependent electric fields. This heating effect was observed experimentally by Meir et al [28]. There, heating was observed even when the EMM was sufficiently compensated and a cloud of ultra-cold atoms was overlapped with an ion in its ground state. In addition to the heating, a deviation from a thermal distribution was observed.

In this paper we study how the energy distribution and dynamics of an ion, immersed in a bath of ultra-cold atoms, inside a linear Paul trap, and beyond the first collision, depends on different trapping parameters. The heating rates and energy distribution of the ion are calculated by a molecular dynamics (MD) simulation that takes into account the back-action of the polarization of the atom on the ion. This simulation is performed over a wide range of trap parameters as well as for different atom and ion species. In addition, we experimentally measured the heating rates of the ion for different trap parameters and compared our measurements to the results of the MD simulation.

This paper is organized as follows. The MD simulation and its underlying assumptions are described in section 2.1. The numerical results of the ion distribution dynamics and its dependence on trap parameters and atom–ion combination are given in section 2.2. In section 3.1 we briefly review the experimental system. The measurements of the heating rates during the first few collisions are presented in section 3.2.

2.1. Model

Without atoms, the ion's motion depends only on the static and dynamic confining potentials as described by the Mathieu equation [29]

$$\begin{eqnarray}&&{\ddot{u}}_{i}+\displaystyle \frac{{{\rm{\Omega }}}^{2}}{4}\left({a}_{i}+2{q}_{i}\cos \left({\rm{\Omega }}t\right)\right){u}_{i}=0,\end{eqnarray} \tag{ 1 }$$

where u_i is the position of the ion in the ith direction (i = x, y, x), a_i (q_i) is the dc (rf) trap parameter and Ω is the rf drive frequency. In a linear Paul trap, these trap parameters are defined by [29]

$$\begin{eqnarray*}\begin{array}{rcl}{a}_{x,y} & = & -\displaystyle \frac{4e}{{m}_{\mathrm{ion}}{{\rm{\Omega }}}^{2}}\displaystyle \frac{{V}_{\mathrm{DC}}}{{Z}_{0}^{2}}\ {a}_{z}=\displaystyle \frac{8e}{{m}_{\mathrm{ion}}{{\rm{\Omega }}}^{2}}\displaystyle \frac{{V}_{\mathrm{DC}}}{{Z}_{0}^{2}}\\ {q}_{x} & = & -{q}_{y}=\displaystyle \frac{2e}{{m}_{\mathrm{ion}}{{\rm{\Omega }}}^{2}}\displaystyle \frac{{V}_{\mathrm{RF}}}{{R}_{0}^{2}}\ {q}_{z}=0,\end{array}\end{eqnarray*}$$

where m_ion is the ion mass, e is electron charge, V_DC and V_RF are the dc and rf voltages on the corresponding electrodes and R₀ and Z₀ are constants arising from the geometry of the dc and rf electrodes, respectively.

In the regime of $\left|{a}_{i}\right|,{q}_{i}^{2}\ll 1$, the solution to the ion trajectory can be written as

$$\begin{eqnarray}&&{u}_{i}\left(t\right)\approx {A}_{i}\cos \left({\omega }_{i}t+{\phi }_{i}\right)\left[1+\displaystyle \frac{{q}_{i}}{2}\cos \left({\rm{\Omega }}t\right)\right],\end{eqnarray} \tag{ 2 }$$

where ${\omega }_{i}=\tfrac{{\rm{\Omega }}}{2}\sqrt{{a}_{i}+{q}_{i}^{2}/2}$ are the secular frequencies and A_i and ϕ_i are the harmonic oscillator amplitude and phase in the ith direction, respectively.

In the presence of an atom, there is a long-range attractive potential which depends on the relative atom–ion distance,

$$\begin{eqnarray}&&V\left(r\right)=-\displaystyle \frac{{C}_{4}}{2{r}^{4}},\end{eqnarray} \tag{ 3 }$$

with C₄ the atom polarizability and r the relative atom–ion distance. At short distances of few nm, there are electronic exchange interactions which cause a strong repulsion between the particles [3]. In principle, the differential cross section depends on the collision energy and the short-range potential. However, in the discussed situation, since the collision energies are relatively high, the differential cross section is angle independent [19]. Hence, in this work, we model the potential as an infinite barrier at a distance of 5nm, and follows equation (3) otherwise. Since collisions are elastic, when reaching the infinite barrier, the atom and the ion leave in a random direction while conserving the total energy and momentum.

The atom–ion polarization potential is long range. Hence, in order to calculate the ion position, one should integrate the equations of motion (EOM) of the ion and all the atoms, which is a formidable task when there are many degrees of freedom. However, we can reduce the number of calculations by using the fact that the atomic gas is relatively dilute, and hence the inter-particle distance is large. For example, for a gas with density of ${10}^{11}\ \tfrac{1}{{\mathrm{cm}}^{3}}$ there is approximately one atom in a sphere with a radius of ∼1 μm. For that, we can define an 'interaction sphere' with radius R_int = 1.2 μm around the equilibrium position of the ion and the position and velocity of both ion and atom are calculated by solving the EOM

$$\begin{eqnarray}\begin{array}{rcl}{\ddot{r}}_{\mathrm{ion},i} & = & -\displaystyle \frac{1}{{m}_{\mathrm{ion}}}\displaystyle \frac{2{C}_{4}\left({r}_{\mathrm{ion},i}-{r}_{\mathrm{atom},i}\right)}{{\left|{{\boldsymbol{r}}}_{\mathrm{ion}}-{{\boldsymbol{r}}}_{\mathrm{atom}}\right|}^{5}}\\ & & -\ \displaystyle \frac{{{\rm{\Omega }}}^{2}}{4}\left({a}_{i}+2{q}_{i}\cos \left({\rm{\Omega }}t\right)\right){r}_{\mathrm{ion},i}\\ {\ddot{r}}_{\mathrm{atom},i} & = & \displaystyle \frac{1}{{m}_{\mathrm{atom}}}\displaystyle \frac{2{C}_{4}\left({r}_{\mathrm{ion},i}-{r}_{\mathrm{atom},i}\right)}{{\left|{{\boldsymbol{r}}}_{\mathrm{ion}}-{{\boldsymbol{r}}}_{\mathrm{atom}}\right|}^{5}},\end{array}\end{eqnarray} \tag{ 4 }$$

where m_atom is the atom mass.

Multiple collisions are simulated by introducing atoms one after another into the interaction sphere. For each atom, we solve the EOM (4) using fourth-order Runge–Kutta method until it exits from the interaction sphere. Then, the position of the ion is evaluated by equation (2) until the next atom enters the interaction sphere.

The atoms enter the interaction sphere with a rate Γ_atoms = nσ_atoms v_th where n is the atomic density, ${\sigma }_{\mathrm{atoms}}=\pi {R}_{\mathrm{int}}^{2}$ (the cross section of a rigid sphere) and ${v}_{{\rm{th}}}=\sqrt{\tfrac{8{k}_{{\rm{B}}}{T}_{{\rm{a}}}}{\pi {m}_{\mathrm{atom}}}}$ is the thermal velocity of atoms at temperature T_a. The interaction sphere radius must be larger than the amplitude of the ion motion. The initial value of R_int is taken to be 1.2 μm. For a typical secular frequency of ω_i/2π ∼ 1 MHz, this radius corresponds to an ion with energy of ∼300 mK. If the ion has a comparable amplitude to R_int after a collision, the interaction sphere radius is increased and remains at the same size until the end of the calculation for that realization.

Since the atoms density is approximately uniform over the ion's trajectory, we sample the entry point of the atoms uniformly on the sphere. The velocity vector is sampled assuming the atom enters from the south pole, $-\hat{z}$, and then is rotated to the chosen position. The velocity distribution of the atoms is thermal, but the velocity component in the radial direction must be positive, directed only into the sphere. Therefore, the velocity amplitude distribution is not the regular Maxwell–Boltzmann distribution. In order to reproduce the correct velocity distribution of the atoms, the amplitude of the velocity, v, and the azimuth angle, φ, and the polar angle (measured relative to $+\hat{z}$), θ, are sampled from the following distributions

$$\begin{eqnarray*}\begin{array}{rcl}{p}_{T}\left(v\right) & = & \displaystyle \frac{2{m}_{{\rm{atom}}}}{{\left({k}_{{\rm{B}}}T\right)}^{2}}{v}^{3}{{\rm{e}}}^{-\tfrac{{m}_{{\rm{atom}}}{v}^{2}}{2{k}_{{\rm{B}}}T}}\ v\in \left[0,\infty \right]\\ p\left(\theta \right) & = & \sin \left(2\theta \right)\ \theta \in \left[0,\displaystyle \frac{\pi }{2}\right]\\ p\left(\varphi \right) & = & \displaystyle \frac{1}{2\pi }\ \varphi \in \left[0,2\pi \right].\end{array}\end{eqnarray*}$$

The amount of calculations that are needed can be reduced further by taking into account the characteristic length scale at which the polarization force becomes larger than the ion trap force, ${R}_{0}={\left(\tfrac{{C}_{4}}{{\omega }^{2}{m}_{\mathrm{ion}}}\right)}^{\tfrac{1}{6}}$ [27]. If the minimal distance between the ion and the atom neglecting the polarization potential is much larger than R₀, there is no close contact and the change to the ion's energy is negligible. Therefore, a full calculation is performed only for atoms that would approach the ion sufficiently close without taking the polarization potential into account. In addition, the full calculation of the EOM is stopped once the atom moved away from the ion by at least R₀. For the ⁸⁸Sr⁺ − ⁸⁷Rb system with typical values of q ∼ 0.1 and ω_i/2π ∼ 1 MHz this characteristic length-scale is R₀ ≈ 64 nm.

When performing the full integration of the EOM (4) and the atom–ion distance is less than a critical distance of 5 nm, elastic hard-sphere collision is assumed and the atom and ion separate at some random angle. Depending on the rf phase, a temporary bound state can be created as the atom may not have enough energy to escape from the polarization potential. In this case, it collides several more times until it gains enough energy [27].

This process is repeated until the elapsed time from the first collision reaches the total interaction time. In order to get the ion energy distribution, many realizations with random initial conditions of the ion are calculated.

2.2. Results

2.2.1. Langevin rate for a trapped particle

The cross section for an atom which collides with an ion was calculated in 1905 by Langevin [30]

$$\begin{eqnarray*}&&{\sigma }_{{\rm{L}}}=\pi \sqrt{\displaystyle \frac{2{C}_{4}}{{E}_{\mathrm{col}}}},\end{eqnarray*}$$

where E_col is the collision energy in the center-of-mass frame. The collision rate is given by Γ = nσ_Lv and v is the relative atom–ion velocity. The MD simulation does not assume the Langevin collision rate a priori, and hence it can be calculated. The ratio between the rate of atoms that enter the interaction sphere and the Langevin collision rate is given theoretically by

$$\begin{eqnarray}&&\displaystyle \frac{{{\rm{\Gamma }}}_{{\rm{L}}}}{{{\rm{\Gamma }}}_{\mathrm{atoms}}}=\sqrt{\displaystyle \frac{{C}_{4}\pi {m}_{\mathrm{atom}}}{2\,\mu {k}_{{\rm{B}}}{T}_{{\rm{a}}}}}\displaystyle \frac{1}{{R}_{\mathrm{int}}^{2}}.\end{eqnarray} \tag{ 5 }$$

Numerically, this ratio should be equal to the ratio between the number of hard-sphere collision to the total number of events. In order to check if the collisions rate is changed for a trapped particle, we simulate collisions inside a fictitious 3D time-independent harmonic trap for the ion (without rf fields). As can be seen in figure 1, the resulting collision rate is lower than the theoretical value (red line). In the weak trap limit, the rate almost converges to the ratio in equation (5) since the particle is nearly free. In the strong trap limit, the particle can be regarded as a particle with an infinite mass and then the rate is given again by equation (5) only with μ = m_a, shown by the yellow line.

Zoom In Zoom Out Reset image size

2.2.2. Energy distribution after a single collision

We know that the energy distribution of the ion after many collisions is not a thermal distribution. However, already after one collision, we see that the energy distribution deviates from thermal distribution. In figure 2, the energy distributions of the ⁸⁸Sr⁺ ion before and after a single collision with a ⁸⁷Rb atom in a Paul trap are shown (red and yellow lines respectively). The initial energy in each motional mode of the ion is sampled from a thermal distribution of a harmonic oscillator with temperature of 50 μK and the atomic cloud temperature is 6 μK. The Paul trap frequency is Ω/2π = 26.5 MHz and the secular frequencies are $\bar{\omega }/2\pi =(0.821,1.29,0.583)$ MHz in the two radial directions and the axial direction, respectively. The corresponding trap parameters are $\bar{q}=(-0.123,0.123,0)$ and $\bar{a}=(-3.7,1.8,1.9)\times {10}^{-3}$.

Zoom In Zoom Out Reset image size

We observe substantial heating after a single collision. A characteristic energy scale for the ion energy gain was derived theoretically in [27]

$$\begin{eqnarray}&&{W}_{0}^{3{\rm{D}}}=\displaystyle \frac{8}{3\pi }{\left(\displaystyle \frac{{m}_{\mathrm{atom}}}{{m}_{\mathrm{ion}}+{m}_{\mathrm{atom}}}\right)}^{5/3}{\left(\displaystyle \frac{{m}_{\mathrm{ion}}^{2}{\omega }^{4}{C}_{4}}{{q}^{2}}\right)}^{1/3}.\end{eqnarray} \tag{ 6 }$$

For the parameters of the simulation, ${W}_{0}^{3{\rm{D}}}=472\,\mu {\rm{K}}$ (indicated by the black vertical line), which agrees with the most probable energy of the simulation.

In order to quantify the energy distribution by small number of parameters, we fit to a Tsallis distribution [17]

$$\begin{eqnarray*}&&{P}_{\mathrm{Tsallis}}\left(E;T,n\right)=\displaystyle \frac{\left(n-3\right)\left(n-2\right)\left(n-1\right)}{2{\left({{nk}}_{{\rm{B}}}T\right)}^{3}}\displaystyle \frac{{E}^{2}}{{\left(1+\tfrac{E}{{{nk}}_{{\rm{B}}}T}\right)}^{n}},\end{eqnarray*}$$

where E is the ion energy, T gives the energy scale (equivalent to temperature) and n describes the power-law tail. This distribution converges in the limit of $n\to \infty$ to a Maxwell–Boltzmann distribution of a gas in a harmonic trap. In order to compare to a thermal distribution using a single parameter, we define ${T}_{\mathrm{ion}}=T\tfrac{n}{n-2}$. This parameter also converges for $n\to \infty$ to the temperature of a thermal distribution.

Fitting the distribution after a single collision to the Tsallis function using maximum likelihood estimation (MLE) gives T = 129 μK and n = 7.8, comparing to T = 44 μK and n = 19 before the collision (equivalent to a thermal distribution). The simulated energy distribution after a single collision in a harmonic trap with the same secular frequencies, is shown by the green line, for comparison. As seen, here the collision does not give a heating effect but cooling, with T = 36 μK and n = 16. This indicates that both the heating and the deviation from a thermal distribution are due to the presence of the rf fields, and occur even after a single collision.

2.2.3. Energy distribution time dependence

To study the heating process dynamics, we sample the energy distribution of the ion at different times. The system is initialized with the same initial conditions as described in the previous section. The time evolution of the ion energy distribution is shown in figure 3. Not only the most probable energy is increasing, but also the high-energies part of the distribution develops a power-law tail. In addition, as can be also seen in the figure inset, the energy distribution evolution converges to a steady state after ∼10 ms, corresponding roughly to 20 collisions. The values shown in the inset are found by fitting our numerical results to the Tsallis distribution using MLE. Although, as we will show in the following, the Tsallis distribution is only a rough approximation to the ion's energy distribution at steady-state.

Zoom In Zoom Out Reset image size

While the Tsallis distribution is an exact limit to the energy distribution of the ion under EMM [24], in the case of an energy distribution generated by a simulation that includes the polarization potential interaction, the Tsallis distribution does not faithfully describe the distribution. To see this, in figure 4, we compare our steady-state result to the steady-state distribution obtained by a simulation with only hard-sphere potential and EMM, shown by the solid green line, similar to the one detailed in [19], which is perfectly described by a Tsallis distribution (dashed black line). However, in the case of a polarization potential simulation (solid blue line) the Tsallis function describes the high energy tail reasonably well, but fails to describe the most probable energy and the low-energy part (dashed purple line). Nonetheless, it describes the distribution much better than a thermal distribution (dashed red line) and therefore was used in this work to compare energy distributions in different trap parameters.

Zoom In Zoom Out Reset image size

2.2.4. Trap parameters dependence

We now turn to investigate the dependence of the energy distribution on the Paul trap parameters, a_i and q_i. These parameters are known to have an effect on the energy gain of the first collision [27] and the steady-state power-law [19] in presence of EMM. Here we want to investigate their effect on the steady-state distribution in the absence of EMM. For a linear symmetric Paul trap, the trap parameters are

$$\begin{eqnarray*}&&{a}_{x}={a}_{y}=-\displaystyle \frac{1}{2}{a}_{z}=-a\ {q}_{x}=-{q}_{y}=q\ {q}_{z}=0.\end{eqnarray*}$$

The average and mode energy gain in the first collision, which we find in the simulation, are shown in figure 5. As seen, the energy gain in the first collision depends strongly on the rf voltage, through the parameter q, but shows almost no dependence on the dc confinement, characterized by the parameter a. This indicates that stronger rf fields can transfer more energy to the ion during the collision, whereas the amplitude of harmonic pseudo-potential has less dominant impact on the heating. In addition, for a thermal distribution, there is a constant ratio between the mean of the distribution to its mode, $\langle E\rangle =\tfrac{3}{2}{E}_{\mathrm{mode}}$. However, here we can see that this is not satisfied, which is an indication of the non-thermal behavior of the system.

Zoom In Zoom Out Reset image size

The quantity W₀ (equation (6)) is indicted in a purple line in figure 5. This formula agrees with our simulation for low q values, but deviates from our observations for q larger than ∼0.3. This might be due to the fact that this energy scale was derived in the absence of dc potentials.

In figure 6 the dependence of the steady-state distribution as function of the same trap parameters is shown. The distribution is described here by the Tsallis distribution parameters, T and n, that are extracted from a MLE fit to the simulation results. As before, the rf voltage has a greater influence on the distribution than the dc voltage. Tighter rf confinement leads to higher temperature and lower n. Lower n means heavier power-law tail and a stronger deviation from a thermal distribution. Similarly to the effect on the first collision, the dependence on the dc parameter is weaker. On the limit of weak rf voltage the distribution is also tending to be hotter with heavier power-law tail. This can be attributed to the lower spring constant in the radial directions that leads to increased heating.

Zoom In Zoom Out Reset image size

2.2.5. Atom–ion combinations

So far we have shown the results of the simulation of our own mixture of ⁸⁸Sr⁺–⁸⁷Rb. Different atom–ion combinations are expected to have different heating rates and steady-state distributions [21]. It was shown [27] that choosing lighter atoms and a heavy ion reduces the effect of pulling the ion from the rf null during the collision. We therefore performed the numerical simulation for several atom–ion systems: ¹⁷⁴Yb⁺–⁷Li, ⁴⁰Ca⁺–⁷Li, ¹⁷⁴Yb⁺–⁸⁷Rb, ¹³⁸Ba⁺–⁸⁷Rb and ⁸⁸Sr⁺–⁸⁷Rb. All other parameters, the dc and rf confinement and trap rf frequency, were kept constant. In figure 7, we can see that choosing a lower atom–ion mass ratio indeed improves the energy distribution in two ways: the characteristic temperature, T, is lower, and the power-law exponent, n, is larger, leading to lower probabilities for high energy events. The lower temperature in systems with ⁷Li atoms is also due to the reduced polarization of the ⁷Li relative to the ⁸⁷Rb [31]. The difference in the polarizability and the reduced mass in the center-of-mass frame increase the s-wave energy threshold. As can be seen in table 1, for systems with ⁸⁷Rb atoms the s-wave energy threshold is roughly four or five orders of magnitudes lower than the steady-state temperature whereas systems with ⁷Li the threshold is comparable with the obtainable steady-state energies [13, 32].

Zoom In Zoom Out Reset image size

Table 1.  Comparison of steady-state energy for various atom–ion species. For every combination, ${T}_{\mathrm{ion}}=T\tfrac{n}{n-2}$ from a fit to the steady-state distribution of the simulation, E_col is the collision energy in the center-of-mass frame and E_s is the s-wave energy limit for comparison.

Ion-atom	$\mu \left(\mathrm{amu}\right)$	${T}_{\mathrm{ion}}\left(\mu {\rm{K}}\right)$	${E}_{\mathrm{col}}\left(\mu {\rm{K}}\right)$	Es(μK)
174Yb+−7Li	6.7	38	1.5	6.4
40Ca+−7Li	6	78	11.6	8.18
174Yb+–87Rb	58	680	227	0.0443
138Ba+−87Rb	53.4	740	286	0.052
88Sr+−87Rb	43.7	1200	596	0.0778

3.1. Methods

3.2. Results

Due to experimental limitations, we cannot test the full a–q space which was simulated. First, we are limited with the maximal values of a ≲ 2 × 10⁻³ and q ≲ 0.15 due to possible voltage breakdown between adjacent electrodes. On the other hand, in order to perform a ground-state cooling and a carrier thermometry, the ion should be in the Lamb–Dicke regime. This sets lower bounds on the secular frequencies, and hence also on the trap parameters. In our case, a lower limit of ∼100 kHz for the secular frequencies, implies q ≳ 0.045 and a ≳ 10⁻⁴. Second, due to a single Doppler-cooling beam used in the experiment, we need to break the radial symmetry in our trap. We do so by applying an additional dc voltage which creates a frequency difference of ∼100 kHz between the two radial modes.

The observed temperatures for different interaction times in different trap parameters are shown in figure 8(a). At each interaction time, a Tsallis temperature was extracted from the Rabi nutation data by MLE assuming a constant power-law (see appendix A). As seen, in all experiments the heating rate is roughly the same within the experimental error, and on the order of $100-200\,\mu {\rm{K}}\,{\mathrm{ms}}^{-1}$. The numerical simulation with the exact experimental parameters predicts minor differences between different parameters, but all within the confidence bound (see figure 8). However, the heating rate in the simulation is much lower and around $50\,\mu {\rm{K}}\,{\mathrm{ms}}^{-1}$.

Zoom In Zoom Out Reset image size

This discrepancy can be explained by several reasons. First, it can arise from a systematic error in estimating the number of atoms through absorption imaging. Another systematic error can arise from the fact that the numerical simulation did not take into account the remaining EMM after compensation. This systematic may vary between different experiments (see table 2 in the supplemental material for estimations). Adding the residual EMM to the simulation changes significantly the heating rates, as can be seen in figure 9 for one of the experiments. The observed results can be explained by an additional uncompensated EMM of ∼500 μK. However, this heating cannot be explained by the effect of uncompensated EMM alone, without the polarization potential, which would give much lower temperatures (dashed line in figure 9). The effect of heating due to the residual EMM is added to the dominant heating caused by the polarization potential pulling. Although we cannot bound EMM in the system below ∼1 mK (equivalent to ∼2.5 V m⁻¹ residual electric field) from our compensation calibrations, our actual EMM is probably at most few 100's μK, since otherwise the observed temperatures would be an order of magnitude higher.

Zoom In Zoom Out Reset image size

Table 2.  EMM energy estimation and residual electric field amplitude in the radial directions for each experiment. The energy and electric field were estimated by comparing the shelving probability on the EMM-sideband transition relatively to the carrier transition. The EMM that enters from the stability of the electrodes is negligible in all experiments.

Experiment name	${T}_{\mathrm{EMM}}\left(\mu {\rm{K}}\right)$	E $\left({\rm{V}}\,{{\rm{m}}}^{-1}\right)$
RF 18 (#1)	860 ± 250	2.8
RF 18 (#2)	700 ± 80	2.5
RF 16	1300 ± 270	2.7
RF 21	330 ± 50	2.5
DC 390 V	850 ± 160	2.8
DC 600 V	1600 ± 230	1.4

Previous studies have shown different heating mechanisms in ultra-cold atom–ion collisions. After eliminating those mechanisms, the inherent heating effect of ion pulling from the zero rf point was shown to have a dominant heating effect on the first collision. Here we studied how this process depends on the different choices of Paul trap parameters, and how it affects the ion energy distribution after few to many collisions. We have used a numerical simulation to gain and compare the energy distribution in steady state and for different atom–ion species, a task that can be experimentally hard to perform. We have shown that the heating effect is weakly dependant on the trap parameters, and for experimental purposes is practically independent. We have also shown that the ion's energy distribution clearly deviates from a thermal one. The distribution features a power-law tail, tough it is not described adequately by the previously proposed Tsallis distribution. The experimental measurement suggests that the heating due to residual EMM adds up to the heating due to the polarization potential. In order to considerably decrease this heating effect, apart from working with light atoms and heavy ion, a trapping method without oscillating field for the ion is required, for example, optical trapping of the ion [34].

This work was supported by the Crown photonics center, the Israeli Science foundation, the Israeli ministry of Science Technology and Space and the European Research council (Consolidator grant 616919-Ionology).

After the interaction with the atoms, the temperature distribution is extracted from electron shelving on the quadruple transition. For a carrier transition on the quadruple transition, the shelving probability of the excited state ('dark state') for the nth level of the harmonic oscillator is given by [29]

$$\begin{eqnarray*}&&{P}_{D}\left({t}_{R};{\boldsymbol{n}}\right)={\sin }^{2}\left({{\rm{\Omega }}}_{{\boldsymbol{n}},{\boldsymbol{n}}}{t}_{R}\right),\end{eqnarray*}$$

where ${\boldsymbol{n}}$ is the harmonic oscillator level (in each of the 3 modes), t_R is the time of the 674 nm pulse and the coupling strength ${{\rm{\Omega }}}_{{\boldsymbol{n}},{\boldsymbol{n}}}$ is defined as

$$\begin{eqnarray}\begin{array}{rcl}{{\rm{\Omega }}}_{{\boldsymbol{n}},{\boldsymbol{n}}} & \equiv & {{\rm{\Omega }}}_{0}\left|\left\langle n| \exp \left\{{\rm{i}}\displaystyle \sum _{j}{\eta }_{j}\left({\hat{a}}_{j}{{\rm{e}}}^{-{\rm{i}}\nu t}+{\hat{a}}_{j}^{\dagger }{{\rm{e}}}^{{\rm{i}}\nu t}\right)\right\}| n\right\rangle \right|\\ & = & \ {{\rm{\Omega }}}_{0}{{\rm{\Pi }}}_{j}{{\rm{e}}}^{-\tfrac{{\eta }_{j}^{2}}{2}}{L}_{{n}_{j}}\left({\eta }_{j}^{2}\right).\end{array}\end{eqnarray} \tag{ 7 }$$

Here, ${L}_{n}\left(x\right)$ is the Laguerre polynomial of order n, Ω₀ is the Rabi frequency at the ground state and ${\eta }_{j}={k}_{j}{x}_{j}=\tfrac{2\pi }{\lambda }\sqrt{\tfrac{{\hslash }}{2m{\omega }_{j}}}\cos \theta$ are the Lamb–Dicke parameters for each mode (λ—laser wavelength, m—ion mass and ω_j the harmonic frequency of the jth mode). Ω₀, θ and ω_j are measured independently.

Due to the ion's energy distribution, the population is divided over many harmonic oscillator levels. Then the probability to be in the excited state is

$$\begin{eqnarray}&&{P}_{D}\left({t}_{R}\right)=\displaystyle \sum _{{\boldsymbol{n}}}P\left({\boldsymbol{n}}\right){\sin }^{2}\left({{\rm{\Omega }}}_{{\boldsymbol{n}},{\boldsymbol{n}}}{t}_{R}\right),\end{eqnarray} \tag{ 8 }$$

where $P\left({\boldsymbol{n}}\right)$ is the energy distribution. However, in our case the probability is given as a function of the total energy $P\left(E\right)$. In the classical limit when $\bar{n}\gg 1$ the energy in the ith mode with n phonons can be approximated as ${E}_{i}\left({n}_{i}\right)\approx {\hslash }{\omega }_{i}{n}_{i}$. The summation over the n's is preformed by taking logarithmic spaced n's for each mode up to some cutoff value, then the corresponding Rabi frequency is calculated by equation (7). The probability for this term is taken from the predefined $P\left(E\right)$ probability (Thermal, Tsallis or numerical). Then, for any given pulse time t_R the shelving probability is expressed as an integral over the energy

$$\begin{eqnarray}&&{P}_{D}\left({t}_{R};P\left(E\right)\right)={\int }_{E}P\left(E\right){\sin }^{2}\left({\rm{\Omega }}\left(E\right){t}_{R}\right){\rm{d}}{E}.\end{eqnarray} \tag{ 9 }$$

Qualitatively, a 'cold' ion will give a high contrast sine-square function, since only few spectral components are dominant, whereas 'hot' ion will give dephasing sine-square with faster decay as the temperature is higher.

For the experiment, the temperature at a specific interaction time is extracted from the data using MLE. For each interaction time with the atoms, a Rabi flop is taken at five different pulse time t_i. For each time, the experiment is repeated N_i times and x_i dark events were observed. The number of dark events has a binomial distribution with probability ${P}_{D}\left({t}_{i};P\left(E\right)\right)$ and N_i number of trails. Therefore, the likelihood of a specific experiment result is

$$\begin{eqnarray}&&L\left(P\left(E\right)| \left\{{x}_{i},{N}_{i}\right\}\right)=\displaystyle \prod _{i}\left(\displaystyle \genfrac{}{}{0em}{}{{N}_{i}}{{x}_{i}}\right){P}_{D}{\left({t}_{i};P\left(E\right)\right)}^{{x}_{i}}\times {\left(1-{P}_{D}\left({t}_{i};P\left(E\right)\right)\right)}^{{N}_{i}-{x}_{i}},\end{eqnarray} \tag{ 10 }$$

In order to simplify the calculations, for the experimental data, the energy distribution for calculating the Rabi flop in equation (8) is assumed to be a Tsallis distribution with a constant power-law n = 4. Hence, the fitting problem is reduced to finding of a single parameter, the temperature T, which maximizes the log-likelihood function. Since the likelihood function has a Gaussian shape, confidence bounds, T_±, of 1σ are found from the condition [35]

$$\begin{eqnarray*}&&\mathrm{log}L({T}_{\pm })-\mathrm{log}L({T}_{\max })=-\displaystyle \frac{1}{2},\end{eqnarray*}$$

where T_max is the maximum of $\mathrm{log}L$. For the MD simulation, at each interaction time, a Rabi flop was calculated by equation (9), with the numerical energy distribution given by the simulation. The MLE and the confidence bounds are calculated in the same method as before.

The process of detecting and compensating EMM is described in detail in [33] and will discussed here only briefly. The EMM arises as a result of non-vanishing rf field in the minimum of the pseudo-potential, for example, due to a uniform dc field. For a static electric field E_dc, the ion motion has an additional term oscillating in the rf frequency. Adding this term to equation (2) [26]

$$\begin{eqnarray}&&{u}_{i}\left(t\right)\approx \left[\displaystyle \frac{{{eE}}_{{\rm{dc}},i}}{m{\omega }_{i}^{2}}+{A}_{i}\cos \left({\omega }_{i}t+\phi \right)\right]\left[1+\displaystyle \frac{{q}_{i}}{2}\cos \left({\rm{\Omega }}t\right)\right].\end{eqnarray} \tag{ 11 }$$

The amplitude of the EMM

$$\begin{eqnarray}&&{u}_{{\rm{EMM}},i}=\displaystyle \frac{{q}_{i}{{eE}}_{{\rm{dc}},i}}{2m{\omega }_{i}^{2}}.\end{eqnarray} \tag{ 12 }$$

In the resolved sideband spectroscopy, this motion causes additional sidebands in the rf frequency. The relative coupling between the carrier and the first EMM sideband is given by

$$\begin{eqnarray}&&\displaystyle \frac{{{\rm{\Omega }}}_{1}}{{{\rm{\Omega }}}_{0}}\approx \displaystyle \frac{{\boldsymbol{k}}\cdot {{\boldsymbol{u}}}_{{\rm{EMM}}}}{2},\end{eqnarray} \tag{ 13 }$$

where Ω₀ (Ω₁) is the carrier (sideband) Rabi frequency and ${\boldsymbol{k}}$ is the laser wave-vector. Minimization of the EMM is done by applying an external electric field and minimizing the Rabi frequency. The EMM is mainly in the radial plane of the trap, and therefore the minimization is done with two orthogonal laser beams and two orthogonal electrodes giving constant electric field in the trap center. The EMM in the axial direction due to static electric field is negligible by trap design. However, rf fields can still persist in the axial direction and they are compensated using a rf fields injected along the axial direction with a controllable phase, relative to the rf of the trap. The EMM in terms of energy is given by

$$\begin{eqnarray}&&{k}_{{\rm{B}}}{T}_{\mathrm{EMM}}=\displaystyle \frac{1}{4}m{{\rm{\Omega }}}^{2}{u}_{\mathrm{EMM}}^{2}=\displaystyle \frac{1}{4}m{{\rm{\Omega }}}^{2}{\left(\displaystyle \frac{2{{\rm{\Omega }}}_{1}}{k\cos {\theta }_{i}{{\rm{\Omega }}}_{0}}\right)}^{2},\end{eqnarray} \tag{ 14 }$$

where θ_i is the angle between the laser beam to the radial plane (measured independently).

Our residual EMM energies after compensation for each experiment are summarized in table 2. The confidence bounds are derived from the projection noise of the shelving probability after compensation. The extracted energies are considerably large than expected. This can be attributed to different sources which are not EMM. First, the pulse time is relatively long, in order to detect weak coupling to the sideband. But, in this timescale the decoherence (for example, due to magnetic field noise), can be dominant. For very low shelving probability, decoherence will increase the population in the excited state, and therefore it will look as EMM.