Observation of vortex-antivortex pairing in decaying 2D turbulence of a superfluid gas

In a two-dimensional (2D) classical fluid, a large-scale flow structure emerges out of turbulence, which is known as the inverse energy cascade where energy flows from small to large length scales. An interesting question is whether this phenomenon can occur in a superfluid, which is inviscid and irrotational by nature. Atomic Bose-Einstein condensates (BECs) of highly oblate geometry provide an experimental venue for studying 2D superfluid turbulence, but their full investigation has been hindered due to a lack of the circulation sign information of individual quantum vortices in a turbulent sample. Here, we demonstrate a vortex sign detection method by using Bragg scattering, and we investigate decaying turbulence in a highly oblate BEC at low temperatures, with our lowest being ~0.5Tc, where Tc is the superfluid critical temperature. We observe that weak spatial pairing between vortices and antivortices develops in the turbulent BEC, which corresponds to the vortex-dipole gas regime predicted for high dissipation. Our results provide a direct quantitative marker for the survey of various 2D turbulence regimes in the BEC system.

Quantum turbulence (QT) is a state of chaotic flow in a superfluid. Because of its inviscidity and quantized circulation, QT constitutes a unique realm in turbulence research. Decades of study involving superfluid helium have revealed many aspects of QT similar to and different from those of turbulence in classical fluids 1,2 , and atomic Bose-Einstein condensates (BECs) were recently employed to extend the scope of QT studies [3][4][5][6] . One of the experimentally unanswered questions is related to the inverse energy cascade in two-dimensional (2D) QT. It is well known that regarding the 2D turbulence of a classical hydrodynamic fluid, the kinetic energy flows toward large length scales, generating a large-scale flow structure due to small-scale forcing 7 . This phenomenon is qualitatively different from three-dimensional turbulence, where energy is dissipated at small length scales. The key issue regarding 2D QT is whether the inverse energy cascade occurs and consequently leads to the formation of a large superflow structure; this issue has drawn a great deal of recent theoretical attention [8][9][10][11][12][13][14][15][16][17][18][19][20][21][22] . Two-dimensional QT is also relevant to the 2D superfluid phase transition which is associated with free vortex proliferation in the Berezinskii-Kosterlitz-Thouless description 23 .
The turbulent flow of an irrotational superfluid is characterized by the configuration of quantum vortices in the superfluid. In 2D, quantum vortices are topological point defects, and the turbulent superfluid can be depicted as a system of interacting 'vortex' particles. This point-vortex picture was introduced by Onsager in his model, which presented a statistical description of classical 2D turbulence 24,25 . The turbulent state is parameterized with the mean vortex energy, where E v is the incompressible kinetic energy of the system and N v is the total vortex number 16,26 . Figure 1 illustrates two vortex configurations for low and high ε v values. In Fig. 1(a), each vortex is adjoined by an antivortex, i.e., a vortex with opposite circulation and their velocity fields cancel each other out in the far-field, thus lowering ε v . A small dipole of the vortex and antivortex undergoes a linear motion, and the low-ε v states are referred to as a vortex-dipole gas regime 16,17,19 . On the other hand, Fig. 1(b) shows a high-ε v state, where vortices of same circulation signs are clustered, constructively enhancing the superflow velocity. Large vortex clusters are called Onsager vortices, which are anticipated to develop as a result of the inverse energy cascade 14,[16][17][18][19] .
How the mean vortex energy ε v changes as the vortex system evolves underlies the inverse energy cascade problem in 2D QT. The system eventually evolves toward a stationary ground state by decreasing both E v and N v via various dissipation mechanisms, such as sound radiation 27 , mutual friction by coexisting thermal components [28][29][30] , and vortex-antivortex pair annihilation 5,29 . It has been noted that the vortex-antivortex annihilation would facilitate the increase of ε v because the contribution of the annihilated vortex dipole to E v is smaller than that of 2ε v ; thereby, it is called evaporative heating 16 . However, when the system is highly dissipative, E v can decrease quickly, even without decreasing N v , thus lowering ε v . Some theoretical studies raised a question regarding the fundamental possibility of ε v increasing in decaying QT 9,15 .
Atomic BECs with highly oblate geometry provide a suitable system for 2D QT 4, 5 , as the vortex line excitations are strongly suppressed along the tightly confining direction 31,32 . Many numerical studies have been performed using the Gross-Pitaevskii (GP) equation and have indicated that various turbulence regimes can exist in the system parameter space spanned by compressibility 8,9 , dissipation [11][12][13][14][16][17][18][19][20] , and trapping geometry 10,33 . In previous experiments, vortex clustering was examined in a forced annular BEC 4 , and the thermal relaxation of turbulent BECs was investigated 5 . However, full characterization of a turbulent BEC has never been achieved. Such a characterization requires measurements of not only the vortex positions but also their circulation directions. Vortex circulation signs might be determined by tracking the motions of individual vortices 34,35 or by analyzing an interference fringe pattern with a stationary reference sample 36,37 , although this is experimentally challenging using a BEC with a complex vortex configuration. A new imaging technique was proposed in which a BEC is tilted before imaging so that each vortex core shows vortex sign-dependent deformation 38 .
In this study, we conduct spatially resolved Bragg spectroscopy to measure the full 2D vortex configuration of a turbulent BEC. Using this method, we examine the evolution of decaying 2D QT in a BEC at low temperatures, with our lowest being ~0.5T c , where T c is the critical temperature of the trapped sample. We observe the development of weak pair correlations between vortices and antivortices in the turbulent BEC, which corresponds to the vortex-dipole gas regime predicted for high dissipation. This work represents the first full experimental characterization of 2D QT in a BEC system and the results reported herein can be a valuable quantitative reference for theories of atomic superfluid turbulence.

Results
Vortex sign detection via Bragg scattering. Our vortex sign detection method is based on the velocity sensitivity of Bragg scattering 39 . Let us consider the situation where a BEC with a singly charged vortex is irradiated by a pair of counterpropagating laser beams along the x′ direction [ Fig. 2(a)]. A two-photon process, which imparts momentum → q and energy ε to an atom, occurs resonantly when ε = q 2 /2m + → ⋅ → q v , where m and → v are the atomic mass and velocity, respectively. Here, → = ′ q kx 2 L  and ε = ħδ, where k L is the wavenumber of the two Bragg beams and δ is their frequency difference. For a positive vortex with counterclockwise circulation, the velocity field is given by with → r being the position from the vortex core, and the resonance condition is given by . Because of the Doppler effect, the scattering response is antisymmetric with respect to the Bragg beam axis [ Fig. 2(b)]; thus, the vortex sign can be determined from the position of the scattered atoms relative to the vortex core. The use of Bragg scattering to measure a superfluid velocity field was demonstrated with a rotating BEC 40 . In this work, we probe high-velocity regions near vortex cores to determine the circulation signs of individual vortices.
We conduct experiments using a BEC of 23Na atoms in the |F = 1, m F = −1〉 state in a pancake-shaped hybrid trap composed of optical and magnetic potentials. The trapping frequencies are (ω x , ω y , ω z ) = 2π × (4.3, 3.5, 350) Hz. For an atom number N = 4.0(3) × 10 6 and a condensate fraction of 80%, the Thomas-Fermi radii are (R x , R y , R z ) = (155, 190, 1.9) μm. The condensate chemical potential is μ ≈ h × 510 Hz and the healing length at peak The vortex dynamics is effectively 2D for R z /ξ ≈ 3 31,32 . Two pairs of Bragg beams are employed by retro-reflecting two laser beams with frequencies of ω L and ω L + δ, which are red-detuned by ≈1.7 GHz from the F = 1 to F′ = 2 transition [ Fig. 3(a)]. We apply a Bragg beam pulse for 600 μs after a short time of flight (TOF) of 300 μs which is initiated by releasing the trapping potential. During the short TOF, the condensate rapidly expands along the tightly confined z direction to reduce the optical depth of the sample for the Bragg beams, but the modification of the transverse velocity field is negligible. After an additional TOF of τ = 9 ms, we take an absorption image of the sample [ Fig. 3(b)]. Two atom clouds are scattered out from the condensate in both the ±x′ directions. Since the displacement due to the initial atomic velocity is negligible, i.e., v x′ τ = |δ d |τ/2k L < 5 μm for |δ d |/2π < 2 kHz, the spatial distributions of the two scattered atom clouds reliably reveal the velocity regions that satisfy the Bragg scattering condition in the condensate. We first apply spatially resolved Bragg spectroscopy to a BEC containing a vortex dipole, i.e., one positive and one negative vortex. A vortex dipole is generated by linearly sweeping the center region of the condensate using a repulsive Gaussian laser beam 34,41,42 , and after a period of 2 s, when the two vortices are well separated in the trapped condensate, we probe the sample using the Bragg beams.  Fig. 3(d,g,l and o), it is apparent that the position of the localized scattering region relative to the vortex core becomes inverted with respect to the Bragg beam axis when the vortex sign or the sign of δ d is changed or when the scattering direction is reversed. This is consistent with the aforementioned antisymmetric response of the vortex state to the Bragg scattering. Furthermore, we find that the density profiles of the scattered atom clouds near the vortex cores are quantitatively accounted for by a theoretical estimation including the spectral broadening of the Bragg scattering (see Supplementary Information).
Probing 2D quantum turbulence. Next, we apply the Bragg scattering method to probe the vortex configuration of a turbulent BEC containing a large number of vortices. Turbulence is generated by stirring the condensate using a repulsive laser beam (see the Methods section). The initial vortex number is N v ≈ 26 and the mean intervortex distance is µ We set δ d /2π = −1.1 kHz, which was observed in the previous experiment to yield a localized scattering signal peaking at r ≈ 13 μm from a vortex core [ Fig. 3(d)]. A higher |δ d | generates a more localized signal but the signal-to-noise ratio is poor.
To facilitate the vortex sign determination, we construct a Bragg signal S B (x′, y′) ≡ n + − n − , where n ± (x′, y′) are the density distributions of the ±x′-scattered atom clouds, which are translated to the condensate reference frame [ Fig. 4(d-f)]. Because n + and n − are complementary to each other due to the local mirror symmetry along the Bragg beam line [Fig. 3(d and l)], S B contains vortex-sign information. In S B , the vortex sign is manifested as the sign of the signal derivative along the y′ direction at the vortex position ′ ′ x y ( , ) with a = 13 μm. When many vortices are located in close proximity to each other, the Bragg signal around some vortices might be weak, and it would be necessary to scrutinize the overall vortex configuration to assign the vortex signs (see Supplementary Information). In particular, in the case of a small vortex dipole for which the surrounding velocity field is almost canceled, the scattering signal is absent, and the vortex signs must then be determined based on the crescent shape of their merged density-depleted cores 5 . Thus, S B and the condensate density distribution provide sufficient information to determine the full vortex configuration of the turbulent BEC.
Vortex-antivortex pairing. The complete determination of the vortex configuration enables us to characterize the evolution of the turbulent state of the BEC (Fig. 5). The BEC relaxes as the vortex number N v decreases, respectively. Interestingly, we observe that the polarization variance, δp 2 , decreases during the evolution [ Fig. 5(a)]. If the vortex decay is a vortex-sign-independent process, δp 2 would increase as δp 2 , similar to diffusion by a random walk. The reduction of δp 2 indicates that a polarized turbulent state is forced to decay into a balanced state and it also suggests that vortex-antivortex pair annihilation is the dominant vortex decay mechanism in a turbulent BEC. We measure the mean vortex dipole moment of the BEC, = | ∑ → | = d s r , where s i = ±1 is the sign of the ith vortex and → r i is its position with respect to the condensate center [ Fig. 5(b)]. In our experiment, ≈ . d R 0 1 throughout the decay evolution. The measured value of d is found to be slightly smaller than the mean value obtained for the random vortex distributions sampled for the same N ± , thus excluding the formation of the Onsager vortex state. We note that in the vortex-dipole gas regime, d can be smaller than the mean dipole moment of random distributions for finite N v .
To examine the vortex pair correlations, we evaluate the second-order vortex sign correlation function , where c ij = 1(0) if the i th vortex and its j th nearest neighbor have the same (different) sign 10 . A truly random configuration yields C 2 = 0.5, and like-sign vortex clusters and vortex dipoles are reflected as increases and decreases of C 2 , respectively. Our experimental data show that C 2 ≈ 0.4 [ Fig. 5(d)], indicating that it is more probable to have neighboring vortices with opposite signs.
We perform a further analysis of the measured vortex configurations by applying the vortex classification algorithm introduced by Billam et al. 18,19 : two vortices are assigned as a dipole if they are the nearest neighbors to each other and have opposite signs; a group of same-sign vortices as a cluster if they are closer to each other than to any other opposite-sign vortex; and the remaining vortices as free vortices [ Fig. 5(e-g)]. In recent numerical studies, it was shown that the fractional populations of dipole vortices, clustered vortices, and free vortices according to this classification scheme provide a unique representation of the 2D QT states, suggested as vortex thermometry 43 . We measure the vortex numbers, N d , N c , and N f , of dipoles, clusters, and free vortices, respectively, where The initial turbulence state shows N d ≈ N c ≈ N f , which is a characteristic of the random vortex configuration 43 , and it is observed that as the decay evolution proceeds, the fractional population of dipole vortices increases to N d /N v ≈ 0.45, whereas that of clustered vortices decreases to N c /N v ≈ 0.2. This observation corroborates the vortex-antivortex pairing in the turbulent BEC.

Discussion
All the results of our vortex configuration analysis demonstrate that vortex-antivortex pair correlations develop in a turbulent BEC under our experimental conditions. It was anticipated that the characteristics of 2D QT in atomic BECs evolve into the vortex-dipole gas regime as the system's dissipation becomes stronger 16,17,19 , but since the quantitative understanding of dissipation in finite-temperature vortex dynamics is still incomplete 44 , there is no theoretical prediction regarding the critical temperature at which the emergence of the Onsager vortex state can be observed 45 . Thermal damping is typically modeled using a few parameters 46 , but it might be questionable whether the damping effects in various vortex dynamics can be fully captured by the parameters. Our experimental results provide quantitative information regarding the dissipative vortex dynamics in 2D QT. For future reference, our main finding is summarized as follows: a highly oblate turbulent BEC with N v ≈ 20 and ξ ≈ R / 290 evolves at T/T c ≈ 0.5 to a state with C 2 ≈ 0.4, N d /N v ≈ 0. 45, and N c This work can be extended to investigate various 2D QT regimes by changing the system parameters. Although it is highly desirable to reach a low dissipation regime by lowering the sample temperature, this was difficult to achieve in our experiment because the sample was heated during the turbulence generation process. It has been noted that utilizing a steep-wall trap instead of a harmonic trap provides a beneficial condition for the formation of Onsager vortex state 33 . Additionally, it might be conceivable to prepare an Onsager vortex state by merging two oppositely rotating BECs and to investigate its relaxation through the vortex-dipole gas regime at high temperatures.
In summary, we have demonstrated the Bragg scattering method for detecting the quantum vortex circulation sign and have successfully applied it to probe decaying 2D QT in a trapped BEC. Various properties of the turbulent BEC were measured based on its vortex configuration and the development of vortex-antivortex pairing was observed in our experiment at finite temperatures. We expect that the Bragg scattering method presented here will enable a direct experimental study of various 2D QT regimes in the atomic BEC system.

Methods
Vortex state preparation. We generated quantum vortices by stirring the center region of a BEC using a focused repulsive Gaussian laser beam as demonstrated in previous experiments 5,34,41,42 . The 1/e 2 beam width was σ ≈ 10 μm and the potential barrier height was V ≈ h × 8 kHz. When we stirred the condensate, the radial trapping frequencies were ω x,y /2π = 7.5 Hz. A tighter trap is helpful for minimizing the dipole motion of the condensate, which might be induced by the stirring. After the vortex generation, the radial trapping potential was adiabatically ramped down within 2 s to the condition of the main experiment. To generate a vortex dipole, we linearly swept the condensate by translating the laser beam in the −y direction over ≈100 μm with a velocity of v = 0.98 mm/s ≈ 0.25c s , where c s is the speed of sound. For generating turbulence, we stirred the condensate in a sinusoidal manner with an amplitude of 40 μm at 15 Hz for 200 ms.
Data availability. The data that support the findings of this study are available from the corresponding author on reasonable request.