Characteristic Length Scale during the Time Evolution of a Turbulent Bose-Einstein Condensate

: Quantum turbulence is characterized by many degrees of freedom interacting non-linearly to produce disordered states, both in space and in time. In this work, we investigate the decaying regime of quantum turbulence in a trapped Bose-Einstein condensate. We present an alternative way of exploring this phenomenon by deﬁning and computing a characteristic length scale, which possesses relevant characteristics to study the establishment of the quantum turbulent regime. We reconstruct the three-dimensional momentum distributions with the inverse Abel transform, as we have done successfully in other works. We present our analysis with both the two- and three-dimensional momentum distributions, discussing their similarities and differences. We argue that the characteristic length allows us to intuitively visualize the time evolution of the turbulent state.


Introduction
Turbulence is a process that occurs in many types of fluids and a broad range of length scales, and it is characterized by chaotic changes in the flow velocity and pressure. The field of quantum turbulence (QT) investigates turbulence in quantum fluids, mainly liquid helium and trapped Bose-Einstein condensates (BECs) [1,2]. Many features of classical turbulence are not entirely explained, to the extent that Feynman deemed it the most important unsolved problem in classical physics [3]. Consequently, dealing with QT may seem a challenging task [3]. However, turbulence in quantum fluids might be more manageable than its classical equivalent because the vortex circulation is quantized in the former and continuous in the latter. Furthermore, from a technical point of view, the advances in cooling and tuning the interparticle interactions in trapped atomic BECs make them attractive candidates for investigating quantum turbulence and connecting it to related fields [4].
The first observation of turbulence in a trapped BEC, and its self-similar expansion, dates to 2009 [5]. Since then, considerable progress has been made in identifying and characterizing QT. A significant breakthrough in the area was the observation of a particle cascade, which appears as a power law in the momentum distribution [6,7], where δ is a positive constant, and its value depends on the mechanism behind the generation of the turbulent state.
There are some intrinsic obstacles in determining the momentum range, where the power law is observed and its characteristic exponent. The range of momentum scales present in trapped BEC systems is small compared to superfluid helium, for example. The exponent of the power law depends on the mechanism behind the turbulence, and different theoretical models predict exponents to be close together. Unfortunately, experiments do not have the necessary precision to distinguish between them. Hence, strategies other than the power law identification have been employed to identify and characterize QT. Energy and particle fluxes have been used in simulations [8,9] and experiments [10,11] to overcome some of these difficulties. Since turbulence and disorder are intimately related, an approach based on the entropy of turbulent BECs has also been successfully applied [12] as an alternative method to investigate and characterize quantum turbulence.
A typical scale used to study turbulence in liquid helium is the vortex line density. Its time dependence provides evidence of the mechanism behind the turbulent regime [13]. In some 4 He experiments, where visualization techniques are well-developed, the geometry and interactions of vortices can be directly observed [14]. In trapped BECs, where the range of length scales available is much smaller, the visualization techniques have not yet reached the same level of detail.
In this work, we employed a length scale associated with the momentum distribution n(k) to study the onset of turbulence and the turbulent state. It is inspired by the integral length scale, which is a quantity commonly used in classical turbulence. If we assume isotropic flow, then the integral length scale L E can be written in terms of the incompressible kinetic energy spectrum E(k) [13,15,16], . ( Casting Equation (2) in this form also illustrates that it is the length scale that contains most of the energy of the system.
It is known from numerical simulations that in some cases turbulence is mainly in the form of waves, and in sother cases mainly in the form of vortices, depending on the excitation protocol and boundary conditions of the system. Since moving vortices radiate waves and strong waves can create vortices, the relative proportion of waves and vortices depends on the particular experiment. In numerical simulations, one has access to the phase of the wave function [17]. Thus, the circulation can be computed to distinguish vortices from waves. Moreover, in the simulations, one can formally identify the compressible kinetic energy, which comes from waves, and the incompressible kinetic energy related to vortices. Unlike in the simulations, in the experiments, we cannot separate waves and vortices so easily.
However, with current experimental techniques, we can measure the momentum distribution n(k) independently of its origin: waves, vortices, or a combination of both. Hence, with Equation (2) in mind, we define the following length scale, . ( Intuitively, L is associated with the scale where most of the particles reside. In this work, we investigated the behavior of this quantity, and we showed that it is possible to use it to study a turbulent BEC. Besides the difficulties mentioned above, there is also an experimental challenge when studying QT in trapped BECs. The momentum distribution of the cloud is obtained using a two-dimensional (2D) projection of the three-dimensional (3D) condensate. We employed the inverse Abel transform, an integral transform that connects the 2D projection of an axially or spherically symmetric function to its 3D value, to reconstruct the momentum distribution of the three-dimensional cloud. We showed that the results for the characteristic length scale are qualitatively the same if calculated using the two-dimensional projection. This indicates that it is possible to study some aspects of the turbulent states using the experimental data directly, without reconstructing the three-dimensional cloud.
This work is structured as follows. First, we provide the experimental details of how the BECs are produced and excited. Then, we present the momentum distributions in both two and three dimensions. These are used to compute the characteristic length scale and other quantities related to it. We discuss both the implications of our findings regarding the length scale and the projection of the cloud. In the Appendix A, we provide the Abel transforms of momentum distributions relevant to our system.

Experimental Procedure
The first step is the production of a Bose-Einstein condensate in equilibrium. A typical BEC contains ≈ 4 × 10 5 rubidium-87 atoms in the hyperfine state |F, m F = |2, 2 , confined in a Quadrupole-Ioffe configuration (QUIC) magnetic trap of frequencies ω r /2π = 130.7(8) Hz and ω x /2π = 21.8(2) Hz. Hereafter, we adopted the convention of reporting the uncertainties as one standard deviation between parenthesis. Before any excitation is applied, the BEC in equilibrium has a condensate fraction of 70(5)%. The chemical potential at the center of the cloud is µ 0 /k B = 124(5) nK, and the healing length is ξ = 0.15(2) µm.
In Figure 1a, we present schematically the protocol we employed to drive the BEC out of equilibrium. Following the condensate production, an oscillating magnetic field is applied while it is still in the trap. The field is produced by a pair of coils, placed in an anti-Helmholtz configuration, with their axis tilted by a small angle of approximately 5 • with respect to the axis of the trap. The excitation potential is given by V exc (r, t) = A[1 − cos (Ωt)]x /R x , where x is the coordinate in the rotated frame and R x = 42 µm is the in-trap extent of the BEC along the x-axis of the trap. Since the perturbation is not aligned with the axes of the trap, the oscillations generate deformations, displacements, and rotations. Several excitation parameters can be varied, such as the amplitude A, total excitation time, and perturbation frequency. In this work, we performed the parametric excitation of fixed frequency Ω/(2π) = 132.8 Hz, close to the radial trapping frequency, ω r /2π = 130.7(8) Hz. We increased the excitation amplitude until reaching a value where the momentum distribution corresponds to an out-of-equilibrium state, and we look for turbulent characteristics. The energy input to the condensate is related to both the total excitation time and the amplitude of the perturbation. Larger amplitudes need less time to reach similar conditions than it would take for smaller amplitudes. The range of amplitudes and excitation times to obtain a turbulent state was the topic of investigation in previous works [18]. In this work, we chose to apply the excitation protocol during a time t exc = 5τ, where τ = 2π/Ω. The amplitude A is varied, ranging from 0 (no perturbation) to 2.2 µ 0 .
After the excitation is turned off, we hold the cloud for a time t inside the trap, often called holding time, which is varied from 20 to 90 ms. During this time, the temporal evolution of the momentum distribution n(k, t), which we are interested in, occurs. Next, we turn off the trap potential and measure n(k, t) using absorption images taken from the ballistic expansion of the cloud after a time of flight (ToF) of t ToF = 30 ms. In Figure 1b, we show typical absorption images corresponding to an excitation amplitude of A = 1.8 µ 0 for different holding times t.
For each excitation amplitude and holding time, we perform several realizations of the experiment and then average the results. In Figure 2, we show a typical two-dimensional momentum distribution obtained from the absorption images. The distance that an atom has traveled from the center of the cloud, after a time t ToF , is given by r =ht ToF k/m, where m is the mass of a rubidium-87 atom. Thus, in practice, the ToF technique corresponds to a Fourier transform of the spatial distribution, which yields the momentum distribution, n(r) ∝ n(ht ToF k/m). There are known shortcomings of the ToF technique, which do not significantly impact the measurement of our momentum distributions, mainly because the turbulent state is kinetically dominated [19]. For a more detailed discussion, the reader is referred to Reference [12] and references therein. This technique has been used successfully in the past to obtain the momentum distribution of turbulent trapped BECs [6,7].

Momentum Distributions
We performed angular averages on the momentum distributions n(k x , k y , t) obtained from the absorption images, such that the resulting profiles depend only on k = (k 2 x + k 2 y ) 1/2 . The two-dimensional momentum distributions n 2D (k, t) are normalized according to 2π dk k n 2D (k, t) = 1.
As discussed above, an experimental challenge when studying momentum distributions of trapped BECs is that the absorption images correspond to a projection of the cloud. We overcome this difficulty by considering the symmetry of the trapped BEC in momentum space. The inverse Abel transform [20][21][22] has been successfully used in the literature [6,7] to obtain the 3D momentum distribution from its two-dimensional projection. It is an integral transform given by We normalized the distributions according to The signature of a particle cascade is a power law, n(k) ∝ k −δ . Equations (4) and (6), together with dimensional analysis, suggest that if we observe a power law both in the two-dimensional momentum distribution, n 2D (k) ∝ k −δ 2D , and in the three-dimensional one, n 3D (k) ∝ k −δ 3D , then their exponents differ by one, δ 3D − δ 2D = 1. To go beyond simple dimension analysis, in Appendix A, we derive this relation analytically for the case where a power law is present over the whole k-range. Although this toy-model is nonphysical, it sheds light on how we can reconstruct the three-dimensional momentum distribution based on symmetry arguments.
In our system, the low-momenta region of n(k) is dominated by the presence of the condensate, which corresponds to a Gaussian distribution. We show in Appendix A that the inverse Abel transform of a Gaussian function is also a Gaussian with the same width. This symmetry is extremely useful because we can work with the two-dimensional projections for quantities related to the Gaussian shape without the need for 3D reconstruction. This is the case of the temperature, which is related to the width of the Gaussian.
All the arguments presented above indicate that the power-law exponents in an ideal situation would be related through δ 3D − δ 2D ≈ 1. However, the fact that we have the power-law behavior superimposed with the condensate at the low-momenta region of the momentum distribution alters this relation. Hence, we need to verify the exponents with the experimental data.
A power law is simply a line in a log-log plot of the momentum distribution as a function of k. We then look for a time window when, in a certain k-range, n(k) is proportional to k −δ , the particle cascade characteristic of a turbulent cloud. We performed the experiment described in Section 2 employing six different excitation amplitudes. For only the three highest ones, A = 1.8, 2.0, and 2.2 µ 0 , we observed the appearance of a power-law, around t ≈ 35 ms and in the region 10 µm −1 k 15 µm −1 . We found the exponents δ 2D = 3.1(1) and δ 3D = 3.8 (2), which lead to δ 3D − δ 2D = 0.7(3). It is interesting to see that even in our finite-sized non-homogeneous system inside an anisotropic trap, we still have δ 3D − δ 2D close to one.
In Figure 3, we present the time evolution of both the two-and three-dimensional momentum distributions for an excitation amplitude of A = 1.8 µ 0 (which is qualitatively the same for A = 2.0 and 2.2 µ 0 ). As we wait after the external excitation has been turned off, the distribution evolves, promoting the population from low to high momentum values, as can be seen in Figure 3.
Some words regarding the values of the exponents we found are in order. To the best of our knowledge, there is no theoretical work that describes all aspects of the experiments we performed, mainly for two reasons. First, our condensate is produced in an anisotropic trap. From the theoretical perspective, it is much easier to implement periodic boundary conditions and describe, or simulate, bulk systems. Second, the route we take is the inverse of most experiments. We begin with a BEC in equilibrium and then excite it, while, for example, quench experiments usually start with a thermal gas and produce a condensate [23,24]. However, after these considerations, we can compare the momentum distributions and exponent we obtained with other works that share similarities with our experiment.
One of the first predictions for the time evolution of the momentum distribution describing Bose-Einstein condensation in a far-from-equilibrium system is given in Ref-erence [25]. The authors find a plateau in the lower momenta region and a power law at higher momenta, akin to what we observe in our experiments; see Figure 3. Since then, much progress has been made in characterizing turbulent flows in BECs. The advances and state of the art concerning this topic can be found in Reference [26], for example. Here, we will discuss two references that capture the essential physical aspects of our experiment.  The authors of [27] address the topic of turbulence in ultracold Bose gases under the light of the so-called non-thermal fixed points [28]. They consider a variety of scenarios and analyze each region of the momentum distributions. For strong turbulence and a freely (without external energy input or dissipation) decaying initial state, scaling arguments lead to δ 3D = 5 prediction, which was confirmed by accompanying numerical simulations. We did not expect quantitative agreement with their results, since we have dissipation in our system, and there are no indications that we are in the strong turbulent regime. Nonetheless, this is one of the few references that address the decay of turbulence.
Another interesting study that allows comparison, to some extent, to our work is presented in Reference [29]. The authors performed numerical simulations employing the forced-dissipated Gross-Pitaevskii equation to study Bose-Einstein condensation under non-equilibrium conditions. They observed that the momentum distribution for the late time dynamics, that is, after the kinetic stage is over, displays an approximately constant value at low-momenta, in accordance with our findings, and a power law in a higher momentum region with the exponent δ 3D = 4.4 (for some values of the nonlinear term). This value is relatively close to what we observed in this work, δ 3D = 3.8 (2). We should remark that the boundary conditions in their simulations and our experiments are different and that we are exciting a BEC while they are looking at the inverse process. However, the agreement warrants further investigations.

The Characteristic Length Scale
We computed the characteristic length scale given by Equation (3) using the experimental data available, where k D ∼0.05 µm −1 and k d ∼30 µm −1 are the smallest and largest wave vectors we can measure, respectively. It is worth noting that although our definition relies on n(k)/k, the singular behavior when k → 0 will never be reached. The lower limit of the integrals, k D , is inversely proportional to the largest length scale of the system, which can be large concerning other scales, but always finite. Using the momentum distributions obtained with different excitation amplitudes and holding times, we can study the time evolution of the BEC, ranging from a slightly perturbed cloud up to a turbulent state. In Figure 4, we present our results for the characteristic length scale computed with both the two-dimensional projection of the cloud and its three-dimensional reconstruction. Although they differ quantitatively, their qualitative behavior is remarkably the same.
The L(t) value can be interpreted as the evolution of the length scale where most of the particles are located. If we think of 1/k as being a weight in Equation (3), then L(t) is related to the inverse of the momentum value for which n(k)/k is peaked. Nucleation of excitations occurs during the excitation, whether in the form of vortices or waves. Then the interaction of these excitations takes place, leading to different stages of deviation from equilibrium.
For small excitation amplitudes, A = 0.8 µ 0 , the system is only slightly disturbed and removed from equilibrium, but it does not have enough energy to reach what is considered a disordered state. In this case, it evolves differently from the others, and the value of L is approximately constant with time.
For intermediate perturbations, A = 1.4 and 1.6 µ 0 , there is a separation between these results and the smallest amplitude, besides a clear dependence with time. For these excitation amplitudes, we are in a regime best characterized as the onset of turbulence. The characteristic length scale decreases on time, indicating the particle transfer to highermomenta, but slower than the higher excitation amplitudes.
In this work and previous investigations [12], we identified the highest excitation amplitudes with turbulent clouds, A = 1.8, 2.0, and 2.2 µ 0 . The value of L(t) at the end of the processes seems to depend on the amplitude and, more importantly, if we deal with a moderate perturbation, the onset of turbulence, or a state with turbulent characteristics. It is interesting to observe that the turbulent states quickly reach lengths comparable to the healing length (ξ = 0.15(2) µm), where dissipation processes are expected to occur. The behavior of the characteristic length scale can be described by an exponential decay, where L 0 is the extrapolation of the characteristic length scale to the instant when the excitation was introduced, and t 0 is its characteristic time.
In Figure 5, we present the values of L 0 and t 0 fitted to the functional form of Equation (8) for different excitation amplitudes. We did not include the results for A = 0.8 µ 0 , since we obtain a value of t 0 ≈ 750 ms, of the same order as the lifetime of the condensate.  Figure 5a shows that the length containing most particles of the system is approximately the same with respect to the excitation amplitudes. However, if we perform the analysis employing the two-or three-dimensional momentum distributions, we arrive at different values for L 0 . A possible reason for this is that the inverse Abel transform slightly shifts the n(k) profiles toward higher-momena-see Figure 3-which then implies a smaller value of L 0 .
In Figure 5b, we present the characteristic time that the particle transfer takes as a function of the excitation amplitude. It is observed that t 0 decreases with amplitude, as expected, since larger amplitudes lead to a faster formation of excitations and, therefore, speed up the decay. It is possible to see that the results obtained in 2D and 3D are in remarkable agreement. This opens the possibility of studying dynamical processes in 2D without the need to reconstruct the three-dimensional cloud, depending on the quantity of interest.

Discussion and Final Remarks
In this work, we defined and computed a characteristic length scale related to the momentum distribution of a trapped BEC, which allowed us to identify distinct out-ofequilibrium stages: a slightly perturbed cloud, the onset of turbulence, and the turbulent state. This quantity complements the formal analysis of identifying a power-law in the momentum distribution as a hallmark of turbulence.
We also focused our efforts on calculating this characteristic length scale using both two-and three-dimensional momentum distributions. The former is obtained straightforwardly from the experimental data, and the latter is reconstructed based on the symmetry of the cloud. From a technical point of view, it is preferable to work only with the twodimensional distributions, since no assumptions about the symmetry of the cloud have to be made. Although the length scales are affected by the inverse Abel transform, which shifts the momentum distributions to higher momenta, the qualitative behavior calculated in 2D and 3D is remarkably similar. The excellent agreement in the characteristic time of the particle transfer indicates that the two-dimensional analysis may be appropriate to investigate dynamical aspects of these systems.
One important remark is that isotropy is assumed through Equation (3) only for the kinetically dominated regions in Fourier space. Such an assumption can be made even for inhomogeneous cigar-shaped clouds. This is because such large-scale inhomogeneities affect only regions in Fourier space up to the order of k ∼ 2π/L min , where L min is the smaller linear size of the cigar-shaped cloud. In a previous work [11], we studied the impact of anisotropy in the energy transfer during the evolution of turbulence in a trapped BEC. Like the integral length scale, the energy flux can also be computed from the kinetic energy spectrum. We found that the turbulent state can be identified and characterized in terms of the energy flux regardless of whether we employ the whole cloud or just a region close to the major axis of the expanded cloud. We should note that the axial trapping frequencies of this work and of Reference [11] are very close; however, in this work, we employ a radial trapping frequency that is ≈0.7 smaller than the one used in Reference [11]. Thus, the BECs in this work are much less elongated and closer to a spherical shape than the ones in Reference [11]. Therefore the range of validity for Equation (3) is even larger than in previous works.
In future works, we intend to vary the excitation protocol to investigate the changes in the characteristic length. Since there is a compromise between the excitation amplitude and time [18], it may prove insightful to investigate situations where the same amount of energy is introduced in the system but varying the time it takes to be injected from very slow inputs up to abrupt changes.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.
Hence, the Abel transform of a Gaussian is a Gaussian of the same width. This is very convenient in the case of the temperature of the cloud, for example, since it is estimated through the width of the Gaussian profile of the momentum distribution.
The second case we considered is a power-law with exponent δ 3D , with A constant. The standard normalization procedure, Equation (6), is going to fail because this momentum distribution is not valid in the entire domain. In reality, the power-law would be observed over a certain k-range, k i k k f with k i = 0. However, this simplified example will have an interesting result as we will see. The Abel transformation yields where Γ is the gamma function, A is another constant, and δ 2D = δ 3D − 1. The conclusion is that the bidimensional projection of a power-law with an exponent of δ 3D in three dimensions is also a power-law, but with the exponent increased by one. Clearly, the momentum distributions we presented in the main text cannot be fully described by these two simple examples. However, they provide indications of the expected behavior of the projection procedure.