Random Matrix Theory Analysis of a Temperature-Related Transformation in Statistics of Fano–Feshbach Resonances in Thulium Atoms

Recently, the transformation from random to chaotic behavior in the statistics of Fano–Feshbach resonances was observed in thulium atoms with rising ensemble temperature. We performed random matrix theory simulations of such spectra and analyzed the resulting statistics in an attempt to understand the mechanism of the transformation. Our simulations show that, when evaluated in terms of the Brody parameter, resonance statistics do not change or change insignificantly when higher temperature resonances are appended to the statistics. In the experiments evaluated, temperature was changed simultaneously with optical dipole trap depth. Thus, simulations included the Stark shift based on the known polarizability of the free atoms and assuming their polarizability remains the same in the bound state. Somewhat surprisingly, we found that, while including the Stark shift does lead to minor statistical changes, it does not change the resonance statistics and, therefore, is not responsible for the experimentally observed statistic transformation. This observation suggests that either our assumption regarding the polarizability of Feshbach molecules is poor or that an additional mechanism changes the statistics and leads to more chaotic statistical behavior.


Introduction
Fano-Feshbach resonances play an important role in the control of interatomic interactions [1,2]. These resonances enable the scattering lengths of the elastic binary collisions [3][4][5] to be changed, turning on and off specific interactions, making Fano-Feshbach resonances a key instrument in quantum simulations [6][7][8]. In the case of lanthanide atoms, Fano-Feshbach resonances have recently attracted a large amount of attention due to the theoretical and experimental demonstration of chaotic statistics in the distribution of these resonances [9][10][11]. We note that, prior to the lanthanide, chaotic spectra were considered in highly excited (Rydberg) states H and He [12] as well as complex atoms and ions [13]. The origin of chaos in Rydberg systems was attributed to the Coulomb interaction [14] or interaction with external fields [15] and is likely not related to the case of lanthanide atoms.
Contrary to the case of erbium and dysprosium, in which chaotic behavior is an intrinsic property of the atomic system, thulium atoms demonstrate both random and chaotic statistics, depending on the temperature of the atomic ensemble [10], in experiments performed in an optical dipole trap (ODT) when temperature was changed together with the trap depth.
Changing the temperature results in two effects: (1) on one side, resonance positions shift with temperature [16], and (2) on the other side, new resonances associated with non-zero mutual angular momentum in the open channel appear as the temperature rises [17][18][19]. Increasing resonance densities could be responsible for the observed change in resonance statistics. Moreover, resonance shifts can change the mutual spacing of the resonances, and therefore change the statistics. Finally, the fact that the temperature changes simultaneously with ODT depth can cause an additional shift in the resonance position, which can also affect the statistics. Such an additional shift in resonance position was indeed detected in previous work [10].
In the present study, we performed random matrix theory (RMT) simulations of the Fano-Feshbach resonances in atomic thulium and investigated the resulting transformations of resonant statistics with changing temperature and ODT power. We show that the temperature-assisted emergence of independent D resonances does not cause random to chaotic statistical changes with rising temperature. Somewhat surprisingly, while demonstrating minor statistical effects within the ODT power range investigated, the resonance shift associated with the Stark effect does not introduce a significant trend, and therefore cannot be responsible for the experimentally observed transformation toward more chaotic statistics.

Materials and Methods
The collisions of two atoms were considered in the center of a mass reference frame, where the motion of these two atoms could be presented as a reduced mass that is falling at the origin. This way the collision can naturally be considered from the point of view of the diatomic molecular potential. Due to internal degrees of freedom, namely components of the total angular momentum M F , the molecule is subject to integration within the external magnetic and electrical fields. In our experiments, each atom was prepared in the lowest energetic state, m F = −4. The Feshbach resonance is observed when the energy of the colliding atoms is equal to the energy of the molecular bound state, corresponding to a different atomic m F value [20].
Such resonances can be categorized by the relative orbital momentum of the colliding particles, or molecular orbital momentum l. Traditionally, for lanthanide atoms [9][10][11], the resonances are characterized by open channels. Thus, S-resonances correspond to l = 0, and D-resonances correspond to the next possible state of the bosons partial wave l = 2. Due to the so-called "centrifugal barrier" [21], that value in general depends on the collision direction and is of the order from 200 to 300 µK, while D-resonances appear only at relatively high temperatures, and S-resonances can be observed at effectively low temperatures [10].
With this in mind, we have considered two experimentally observed [10] sets of resonances: (1) the S-resonances obtained at 2 µK and (2) the S+D resonances observed at 12 µK. We note here that, in the experiment, the temperature was changed together with ODT depth, and therefore the Stack shift of the atomic levels should be considered along with other parameters.
RMT is a quite general approach that allows evaluation of the statistical properties of complex system eigenstates [22]. It was initially developed to study complex scattering properties in nuclear physics [23] and was recently used to examine the collisional properties of ultracold lanthanide atoms [9,11,24]. The RMT approach uses an ensemble of random Hamiltonians, instead of an exact Hamiltonian that is unknown, and preserves some of the statistical properties of the energy levels.
We use an adaptation of the approach presented in references [9,24] to model the molecular bound states and corresponding Feshbach resonances of thulium atoms.
The molecular Hamiltonian of diatomic Tm 169 molecule is modeled in the atomic basis with full atomic moments of each atom |F 1 , m F1 |F 2 , m F2 together with partial wave angular momentum and its projections |l, m l : where H b , H z . and H S have only diagonal matrix entries. H b and H cpl correspond to the molecular Hamiltonian of two interacting atoms in the free space, while H z and H S stand for the Zeeman and Stark shifts, respectively. In conjunction, H z and H S shift the specific energy state of the two atoms. H cpl determines the off-diagonal entries and is responsible for anisotropic coupling originating from the dipole-dipole and anisotropic part of the van der Waals dispersion interaction. Following reference [24], entries of the matrix H b were randomly sampled so that nearest-neighbor spacings obey Exponential distribution (i.e., distribution of distances between adjacent elements in sequence generated by Poisson process) and have mean value of . This corresponds to independent resonance positions distributed according to Poisson distribution and assumes no chaos in the initial distribution of energy levels. Taken as the perturbation term to atomic states, H S entries are given by: where F is the total atom angular momentum quantum number and m F is its projection quantum number, I is the peak intensity of the trapping laser beam, and is the Stark energy shift, which depends on laser beam intensity and geometry as well as atom polarizability [25]. This assumes that the experimentally determined [25] polarizability of free atoms remains the same in the weakly bound molecular state. For each hyperfine component (i.e., That leads to 16 states with equal probabilities corresponding to all possible total angular momentum projections in both hyperfine sublevels (see detail in Supplementary Materials). H Z entries are given by where M F = m F1 + m F2 , g F is the Lande g-factor depending on the hyperfine component, µ B is the Bohr magneton, and B is the magnetic field value that causes shift of initial and bound states energy level through the Zeeman term. H cpl entries were sampled according to normal distribution with zero mean and v cpl variance, where v cpl sets the strength of coupling. With the matrix entries distributed according to the Gaussian orthogonal ensemble, this off-diagonal coupling Hamiltonian introduces chaotic component into the random (Poissonian) distribution of the molecular bound states. In other words, the chaos is arising from the anisotropic coupling. For each realization of the Hamiltonian H, we calculated the corresponding eigenvalues on the magnetic field from 0 to 20 G, the same range used in the experiment [10], and obtain a set of bound-state energies and magnetic field values. The Feshbach resonance position was found by determining the magnetic field at which the energy of the bound state is equal to the entrance channel energy.

Results
The initial state of thulium gas in our trap is F = 4, m F = −4 [26], meaning that the entrance channel energy corresponds to two free m F = −4 atoms. Thus, we follow [9] and subtract the entrance channel energy at each magnetic field value, so that the position of the eigenvalue crossing zero corresponds to the Feshbach resonance position (see Figure 1).  Table 1.
We characterized the chaotic behavior of the system by analyzing nearest-neighbor spacings in the Feshbach resonance (NNSFR) spectra. In previous analyses of distributions of resonance positions and transformation from random to chaotic statistics, the Berry-Robnik distribution [22] was used. However, this distribution employs a complicated analytical expression and, therefore, is not suitable for the analysis below. An alternative approach is to utilize the fully analytical Brody distribution [9,11,24]. The Brody distribution is an empirical distribution allowing to quantify between exponential and Wigner-Dyson distributions smoothly changing from purely Exponential at 0 η = to purely Wigner-Dyson at 1 η = .
The Brody distribution is given by [22,27]: Here, b is the normalization constant given by: To gain insight into random to chaos transformations [10] in the statistics of Feshbach spectrum resonances, the S and D resonances were modeled as independent spectra having their own mean energy spacings,  Table 1. These parameters were chosen so that the individual spectra generated using the RMT model would reproduce the experimental S and D resonance spectra statistics with regard to nearest neighbor spacings of Feshbach resonances' distribution parameters, such as resonance densities,   Table 1. We characterized the chaotic behavior of the system by analyzing nearest-neighbor spacings in the Feshbach resonance (NNSFR) spectra. In previous analyses of distributions of resonance positions and transformation from random to chaotic statistics, the Berry-Robnik distribution [22] was used. However, this distribution employs a complicated analytical expression and, therefore, is not suitable for the analysis below. An alternative approach is to utilize the fully analytical Brody distribution [9,11,24]. The Brody distribution is an empirical distribution allowing to quantify between exponential and Wigner-Dyson distributions smoothly changing from purely Exponential at η = 0 to purely Wigner-Dyson at η = 1.
The Brody distribution is given by [22,27]: Here, b is the normalization constant given by: , where Γ is the Gamma function. To gain insight into random to chaos transformations [10] in the statistics of Feshbach spectrum resonances, the S and D resonances were modeled as independent spectra having their own mean energy spacings, S /h = 5.6 MHz for S-resonances and D /h = 10.7 MHz for D-resonaces, respectively, and coupling strength parameters, ν S cpl = 0.7 MHz and η RMT S = 0.07 ± 0.03ν D cpl /h = 2.9 MHz. Details of the simulation, such as the number of included levels, discretization, etc., are summarized in Table 1. These parameters were chosen so that the individual spectra generated using the RMT model would reproduce the experimental S and D resonance spectra statistics with regard to nearest neighbor spacings of Feshbach resonances' distribution parameters, such as resonance densities, ρ S 2.2µK = 2 and ρ D 12µK = 1, and Brody parameters, η S = 0.08 and η D = 0.21. The corresponding experimental NNSFR empirical cumulative distribution function (ECDF) obtained for S and D resonances, alone, are shown in Figure 2A  Mixed spectrum analysis shows that the situation becomes more random, resulting in a Brody parameter of 0.07 0.03 Figure 3B). In contrast, experimental mixed spectra display significantly more chaotic behavior, with a Brody parameter of Figure 3A). Mixed spectrum analysis shows that the situation becomes more random, resulting in a Brody parameter of η RMT S+D = 0.07 ± 0.03 (see Figure 3B). In contrast, experimental mixed spectra display significantly more chaotic behavior, with a Brody parameter of η S+D = 0.63 (see Figure 3A). Mixed spectrum analysis shows that the situation becomes more random, resulting in a Brody parameter of 0.07 0.03 Figure 3B). In contrast, experimental mixed spectra display significantly more chaotic behavior, with a Brody parameter of 0.63 S D + η = (see Figure 3A). Previous research has proposed that the reason for statistics transformation may be a change in the resonances caused by the Stark shift [10]. The tensor polarizability of the thulium atom is quite substantial [25] and, therefore, can naturally repulse resonances. Thus, one could expect some correlations between S-and D-resonance spectra due to the presence of the Stark shift in the optical dipole trap. To address this possibility, we simulated spectra statistics using trapping beam power values ranging from 0 to 5 W for both S-and S+D resonances and used the experimentally determined polarizability of free atoms to calculate shifts in both opened and closed channels. The Brody parameter was fitted for each trapping beam power level. The results (the Brody parameter η versus ODT beam power) are presented in Figure 4. Simulation parameters and results corresponding to experimental conditions are summarized in Table 2. Previous research has proposed that the reason for statistics transformation may be a change in the resonances caused by the Stark shift [10]. The tensor polarizability of the thulium atom is quite substantial [25] and, therefore, can naturally repulse resonances. Thus, one could expect some correlations between S-and D-resonance spectra due to the presence of the Stark shift in the optical dipole trap. To address this possibility, we simulated spectra statistics using trapping beam power values ranging from 0 to 5 W for both S-and S+D resonances and used the experimentally determined polarizability of free atoms to calculate shifts in both opened and closed channels. The Brody parameter was fitted for each trapping beam power level. The results (the Brody parameter η versus ODT beam power) are presented in Figure 4. Simulation parameters and results corresponding to experimental conditions are summarized in Table 2.     Figure 4 shows that, while a noticeable correlation exists with the experimental range of intensities, the Brody parameter for each distribution does not change significantly. Thus, the experimentally observed transformation from random to chaotic behavior could not be explained by modeling of S-and D-resonances independently. A potential weakness in the above calculation is the use of the free atom polarizabilities for both opened and closed channels. It is possible that the molecular polarizability is considerably larger than the atomic polarizability, which would enhance the impact of the Stark effect. However, this scenario is not very likely. In the simulations, the power was varied across a much larger range than that utilized experimentally. So, if the polarizability were underestimated, the effects would be revealed at larger power values. The Brody parameter seems, instead, to decrease, if it changes at all, at high ODT power values (see Figure 4). Therefore, we surmise that there is some other mechanism responsible for the correlation between and/or within the S-and D-resonance spectra.

Conclusions
The RMT analysis of the Feshbach spectra of thulium atoms was performed in two cases: (1) changes in resonance density due to the appearance of new resonances with rising temperature and (2) changes in resonance density due to the Stark shift caused by changes in ODT power. Assuming independent S-and D-resonances, our calculations were unable to explain the transformation from random to chaotic statistics observed in the previously reported experiments. In all the scenarios considered herein, the simulated Brody parameter η ≈ 0.1 is significantly smaller than the experimentally observed value of η = 0.63. This suggests the existence of a significant correlation between the S-and D-resonance spectra. Alternatively, it is possible, that some details of molecular interactions not included in the rather common model appear to be important in the formation of the chaos. Another possibility is that RMT does not provide a sufficient representation for anisotropic interaction between Tm atoms in an ultracold regime.