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

The holy grail of ion-neutral systems is reaching the s-wave scattering regime. However, most of these systems have a fundamental lower collision energy limit which is higher than this s-wave regime. This limit arises from the time-dependant trapping potential of the ion, the Paul trap. In this work, we studied both theoretically and experimentally, the way the Paul trap parameters affect the energy distribution of an ion that is immersed in a bath of ultra-cold atoms. Heating rates and energy distributions of the ion are calculated for various trap parameters by a molecular dynamics (MD) simulation that takes into account the attractive atom-ion potential. The deviation of the energy distribution from a thermal one is discussed. Using the MD simulation, the heating dynamics for different atom-ion combinations is also investigated. In addition, we performed measurements of the heating rates of a ground-state cooled $\ ^{88}$Sr$^+$ ion that is immersed in an ultra-cold cloud of $\ ^{87}$Rb atoms, over a wide range of trap parameters, and compare our results to the MD simulation. Both the simulation and the experiment reveal no significant change in the heating for different parameters of the trap. However, in the experiment a slightly higher global heating is observed, relative to the simulation.


I. INTRODUCTION
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][8][9][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 fact that the ion is trapped using a Paul trap [14], which is based on a time-dependent potential. This potential can inject energy into the system during a collision. 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 ex-pect 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 powerlaw tail of the distribution depends on the atom-ion mass ratio as well as the specific Paul trap parameters [19][20][21][22][23][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). Excess micromotion 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 nonvanishing 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 II A. The numerical results of the ion distribution dynamics and its dependence on trap parameters and atom-ion combination are given in section II B. In section III A we briefly review the experimental system. The measurements of the heating rates during the first few collisions are presented in section III B.

A. Model
Without atoms, the ion's motion depends only on the static and dynamic confining potentials as described by the Mathieu equation [29], where u i is the position of the ion in the i-th 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], 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 0 and Z 0 are constants arising from the geometry of the dc and rf electrodes, respectively.
In the regime of |a i | , q 2 i 1, the solution to the ion trajectory can be written as, where ω i = Ω 2 a i + q 2 i /2 are the secular frequencies and A i and φ i are the harmonic oscillator amplitude and phase in the i-th direction, respectively.
In the presence of an atom, there is a long-range attractive potential which depends on the relative atom-ion distance, with C 4 the atom polarizability and r the relative atomion 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 Eq. 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 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 1 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 equations of motion (EOM), − Ω 2 4 (a i + 2q i cos (Ωt)) r ion,ï 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 equations of motion (4) using fourthorder Runge-Kutta method until it exits from the interaction sphere. Then, the position of the ion is evaluated by Eq. (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, σ atoms = πR 2 int (the cross section of a rigid sphere) and v th = 8kBTa πmatom 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, −ẑ, 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-Boltzamnn 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 +ẑ), θ, are sampled from the following distributions, 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 = C4 ω 2 mion 1 /6 [27]. If the minimal distance between the ion and the atom neglecting the polarization potential is much larger than R 0 , 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 0 . For the 88 Sr + − 87 Rb system with typical values of q ∼ 0.1 and ω i /2π ∼ 1 MHz this characteristic length-scale is R 0 ≈ 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.
FIG. 1. Langevin rate for atom-ion collisions in a 3D spherical symmetric harmonic trap, normalized to the collision rate of free particles. Each point is an average on 10 4 repetitions. The ion starts at rest and the atoms have a thermal energy distribution with Ta = 6 µK. Red line is the theoretical Langevin rate for two free particles. Yellow line is the limit of an ion with an infinite mass (µ = ma).

Langevin rate for a trapped particle
The cross section for an atom which collides with an ion was calculated in 1905 by Paul Langevin [30], where E col is the collision energy in the center of mass frame. The collision rate is given by Γ = nσ L v 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 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 Fig. 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 Eq. 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 Eq. 5 only with µ = m a , shown by the yellow line.

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 Fig. 2, the energy distributions of the 88 Sr + ion before and after a single collision with a 87 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ω/2π = (0.821, 1.29, 0.583) MHz in the two radial directions and the axial direction, respectively. The corresponding trap parameters areq = (−0.123, 0.123, 0) andā = (−3.7, 1.8, 1.9) · 10 −3 .
We observe substantial heating after a single collision. A characteristic energy scale for the ion energy gain was derived theoretically in Ref. [27], For the parameters of the simulation, W 3D 0 = 472 µ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], 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 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 ion = T n n−2 . This parameter also converges for n → ∞ 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.

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 Fig. 3. Not only the most probable energy is increasing, but also the high-energies part of the distribution develops a powerlaw 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.
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 Fig. 4, we compare our steadystate 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 Ref. [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.

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, a x = a y = − 1 2 a z = −a q x = −q y = q q z = 0.
The average and mode energy gain in the first collision, which we find in the simulation, are shown in Fig.  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 Ion energy distribution at steady state for hardsphere potential (green) and ∼ − 1 r 4 potential (blue) for the same experimental parameters as in Fig. 3. The distribution was calculated for 5 · 10 6 repetitions. Tsallis distribution was fitted using MLE (dashed lines) to all energies after entering steady-state. In the hard-sphere potential simulation, residual EMM equivalent to 50 µK was added, in the absence of the atom-ion polarization potential.
heating. In addition, for a thermal distribution, there is a constant ratio between the mean of the distribution to its mode, E = 3 2 E mode . However, here we can see that this is not satisfied, which is an indication of the non-thermal behavior of the system.
The quantity W 0 (Eq. 6) is indicted in a purple line in Fig. 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 Fig. 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. 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: 174 Yb + − 7 Li, 40 Ca + − 7 Li, 174 Yb + − 87 Rb, 138 Ba + − 87 Rb and 88 Sr + − 87 Rb. All other parameters, the dc and rf confinement and trap rf frequency, were kept constant. In Fig. 7, we can see that choosing a lower atom-ion mass ratio indeed im- For every combination, Tion = T 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 Es is the s-wave energy limit for comparison. proves 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 7 Li atoms is also due to the reduced polarization of the 7 Li relative to the 87 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 I, for systems with 87 Rb atoms the s-wave energy threshold is roughly four or five orders of magnitudes lower than the steady-state temperature whereas systems with 7 Li the threshold is comparable with the obtainable steady-state energies [13,32].

A. Methods
A full review of our experimental system is given elsewhere [33]. Briefly, our system consists of two connected vacuum chambers. In the upper chamber the atoms are collected and cooled to mK temperature by a magnetooptical-trap (MOT) and then evaporativly cooled in a CO 2 quasi-static dipole trap to µK temperature. Subsequently, atoms are loaded into a 1D lattice created by two counter-propagating YAG laser beams. The lattice transfers the atoms to the lower chamber where the ion is trapped. In the lower chamber, the atoms are transferred into a crossed-dipole trap. The position of the crosseddipole trap is controlled by a PZT-controlled mirror. After the ion is spin-polarized and ground-state cooled, the crossed dipole trap is moved so that the atomic cloud overlaps the ion. After a given interaction time the atoms are released from the trap and following a short time-offlight (TOF) are imaged. From the TOF images, the density and the temperature of the atoms are extracted.
The ion is trapped in a segmented linear Paul trap with controlled static (dc) and dynamic (rf) potentials. The ion is initially Doppler cooled, and then ground-state cooled using resolved side-band cooling on the quadruple transition 4d 2 D 5/2 −5s 2 S 1/2 . In order to reduce the heating due to EMM, it is compensated by applying external electric fields and minimizing the coupling to the EMM resolved side-bands of the quadruple transition [33]. This process was performed before each experiment and periodically every ∼40 minutes during the experiment.

B. 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 −3 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 −4 . 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 Fig. 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 IV A). As seen, in all experiments the heating rate is roughly the same within the experimental error, and on the order of 100-200 µK ms . The numerical simulation with the exact experimental parameters predicts minor differences between different parameters, but all within the confidence bound (see Fig.  8). However, the heating rate in the simulation is much lower and around 50 µK ms . 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 II in the supplemental material for estimations). Adding the residual EMM to the simulation changes significantly the heating rates, as can be seen in Fig. 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 Fig. 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 ∼ 1mK 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.

IV. DISCUSSION
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 stud- ied 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)

APPENDICES A. Carrier thermometry
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 n-th level of the harmonic oscillator is given by [29] P D (t R ; n) = sin 2 (Ω n,n t R ) , where n is the harmonic oscillator level (in each of the 3 modes), t R is the time of the 674nm pulse and the coupling strength Ω n,n is defined as, Here, L n (x) is the Laguerre polynomial of order n, Ω 0 is the rabi frequency at the ground state and η j = k j x j = 2π λ 2mωj cos θ are the Lamb-Dicke parameters for each mode (λ -laser wavelength, m -ion mass and ω j the harmonic frequency of the j-th mode). Ω 0 , θ 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 P D (t R ) = n P (n) sin 2 (Ω n,n t R ) .
Where P (n) is the energy distribution. However, in our case the probability is given as a function of the total energy P (E). In the classical limit whenn 1 the energy in the i-th mode with n phonons can be approximated as E i (n i ) ≈ ω 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 Eq. (7). The probability for this term is taken from the predefined P (E) probability (Thermal, Tsallis or numerical). Then, for any given pulse time t R the shelving probability is expressed as an integral over the energy 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 sinesquare with faster decay as the temperature is higher.

B. Maximum likelihood estimation
For the experiment, the temperature at a specific interaction time is extracted from the data using Maximum likelihood estimation (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 (t i ; P (E)) and N i number of trails. Therefore, the likelihood of a specific experiment result is In order to simplify the calculations, for the experimental data, the energy distribution for calculating the Rabi flop in Eq. 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 [1] log L(T ± ) − log L(T max ) = − 1 2 , where T max is the maximum of log L. For the MD simulation, at each interaction time, a Rabi flop was calculated by Eq. 9, with the numerical energy distribution given by the simulation. The MLE and the confidence bounds are calculated in the same method as before.

C. EMM estimation
The process of detecting and compensating excess micro-motion (EMM) is described in detail in Ref. [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 ad- 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, where Ω 0 (Ω 1 ) is the carrier (sideband) Rabi frequency and 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, 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 II. 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 side band. 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.