Incommensurable matter-wave jets in quasi-1D geometry

Introduction

Incommensurability is encountered in aperiodic crystals with incommensurable phases^1,2,3, charge and spin structures^4,5,6,7,8, twisted moiré bilayer experiments⁹ with Hofstadter butterfly energy structure¹⁰, and surface layers¹¹. In the context of cold atoms, incommensurable optical lattices have been used to study wavefunction localization phenomena in disordered potentials within the Aubry-André model^{12,13,14,15,16,17,18,19,20}.

Matter-wave jets occur when a Bose-Einstein condensate (BEC) is subjected to a modulation of interaction²¹. The modulation leads to the formation of density waves inside the BEC followed by the emission of jets²². The microscopic process behind jet emission is photon-stimulated two-atom collisions^21,23, where a quantum of the magnetic field modulation energy is transferred to the atom pair. First-order jets form from pairs of atoms that each absorb half of the energy of the photon. Further collisions lead to higher-order jets having velocities that can either be rational multiples of first-order jet velocity (commensurate jets) or not (incommensurate jets).

In our previous work we studied the formation of first- and second-order commensurate jets in a quasi-one-dimensional geometry²⁴. We observed the jets experimentally and in numerical simulations, however, no incommensurate jets were seen, because they only occur at higher modulation amplitudes. Previous studies of matter-wave jets in two dimensions (2D) have revealed the formation of density waves inside the condensate prior to the emission of jets and explained the process of formation of higher-order jets^22,23. In ref. ²³ higher-order jets with velocities (sqrt{2}) and 2 times the velocity of first-order jets were directly observed, commensurate with the first order. Additionally, a hint of incommensurate (sqrt{3}) jets was reported. It is important to note that in 2D geometry the (sqrt{2}) jets actually result from a density modulation with the velocity vector at a 45° angle with respect to the primary density modulation, so that each component of the velocity vector is, in fact, commensurate with the primary modulation. In contrast, (sqrt{2}) jets in 1D geometry result from a genuine incommensurability in the BEC.

Here, we experimentally and numerically demonstrate the formation of incommensurable “golden” jets and other higher-order matter-wave jets in a quasi-1D BEC under a single-frequency interaction modulation. The incommensurable golden jets are precursors of the supercontinuum formation, which remains to be experimentally confirmed.

Results

Formation of matter-wave jets

We start by preparing an almost pure BEC of approximately 4000 cesium atoms in a crossed dipole trap. We release it into a channel with radial frequency ω_r = 2π ⋅ 90 Hz and simultaneously change the interaction between the atoms to a value of a_dc that is close to zero via the wide s-wave Feshbach resonance²⁴ with zero-crossing at 17.1 G. This results in the formation of a BEC soliton. For 4000 atoms the BEC forms a soliton at interaction a_dc = − 2.7a₀²⁵, where a₀ is the Bohr radius. To induce the emission of matter-wave jets, we modulate the magnetic field and therefore the interaction between the atoms as (a(t)={a}_{{{{rm{dc}}}}}+{a}_{{{{rm{ac}}}}}sin (omega t)) for a time t_p, where a_ac and ω are the amplitude and the frequency of the modulation, respectively²¹ (see Methods). We let the system evolve for t_e in the channel and take an absorption image after a 15 ms time-of-flight expansion. The results are shown in Fig. 1a. In addition to previously observed jets of first order (J1) and second order (J2)²⁴, we observe two pairs of jets with dimensionless velocities

$${u}_{{{{rm{G2}}}}}=frac{{v}_{{{{rm{G2}}}}}}{{v}_{{{{rm{0}}}}}}=pm frac{1pm sqrt{5}}{2},$$

(1)

where ({v}_{{{{rm{0}}}}}=sqrt{hslash omega /m}) is the velocity of the first-order jets (J1), ℏ is the reduced Planck constant, and m is the mass of an atom. We name these “golden jets” (G2) for their golden ratio dimensionless velocities. The irrational velocity ratio between first order and golden jets implies the coexistence of two incommensurable density waves inside the condensate prior to jet emission. Figure 1a shows three experimental realizations showing golden jets. As can be seen in Fig. 1a the golden jets are not necessarily symmetric around J0 and change in each experimental run despite the supposed symmetry of the involved formation processes. We attribute this to the asymmetry in the initial BEC which is slightly inhomogeneous and fluctuates between experimental runs.

**Fig. 1: Experimental observation of golden jets and jet formation processes.**

The origin of golden jets can be understood from the conservation of energy (frac{m{v}_{1}^{2}}{2}+frac{m{v}_{2}^{2}}{2}+hslash omega =frac{m{v}_{1}^{{prime} 2}}{2}+frac{m{v}_{2}^{{prime} 2}}{2}), and momentum: (m{v}_{1}+m{v}_{2}=m{v}_{1}^{{prime} }+m{v}_{2}^{{prime} }), where v_i (({v}_{i}^{{prime} })) are the initial (final) velocities of the atoms and the momentum of the photon is neglected. A pair of atoms in J1 is formed from a pair of colliding atoms from the original condensate (J0) by absorbing one quantum of energy from the modulating magnetic field. Similarly, a pair J0-J2 forms from a pair of J1 atoms with the same velocity and a pair G2-G2 from a pair J0-J1. In general the final velocities of the resulting jets are

$${u}_{1,2}^{{prime} }=frac{({u}_{1}+{u}_{2})pm sqrt{{({u}_{1}-{u}_{2})}^{2}+4}}{2},$$

(2)

where u_1,2 are the initial dimensionless velocities of the atom pair. With each absorbed photon the number of possible jets increases. By the number of photons that an atom in a jet has absorbed we define the order of the jets: the initial BEC is zeroth order, J1 are first order, J2, G2 and D2 are second order and so on, see Fig. 1(b) and Table 1. The total number of jets increases as a double exponential with order, going as 1, 3, 11, 109, 11605,…, and can be approximately calculated as (1.794{9}^{{2}^{n}}) (from the fit up to n = 4, where n is the order). Experimentally, we observe J1, J2 and G2, but we have not been able to observe jets with dimensionless velocity (pm sqrt{2}) (D2) or any third order jets. J2 are formed from J1 which move with the same velocity so the process is equivalent to the primary process in a moving frame, while the atoms that form G2 move relative to each other, which reduces the interaction time, making golden jets rare. For D2 the velocity difference is two times higher, explaining why we do not observe them. Additionally, the formation of J2 is more probable due to bosonic stimulation, since they form in pair with J0, which is present in the system beforehand²³.

Table 1 Jet orders, abbreviations, and formation processes

Full size table

Numerical simulation of Gross-Pitaevskii equation

We turn to numerical simulations of 1D Gross-Pitaevskii equation (GPE) to better understand the conditions for the formation of each order of jets. As in our previous work²⁴, we start with a soliton as the initial state and modulate the interaction parameter (k(t)=Na(t)/{a}_{r}={k}_{{{{rm{dc}}}}}+{k}_{{{{rm{ac}}}}}sin (omega t)) for a time t_p = 8π with different frequencies ω and amplitudes k_ac. We leave the wavefunction to evolve for t_e and examine the resulting density profile n(x) = ∣ψ(x)∣². The details about the simulation are provided in Methods and in our previous work²⁴. In the simulations, the time t is given in units 1/ω_r, ω in units ω_r, and distances in units ({a}_{r}=sqrt{hslash /m{omega }_{r}}).

An example of the final density profile can be seen in Fig. 2a where all visible jets are marked. Figure 2b shows the time evolution of the density profile from the start of modulation. In contrast with the experiment, we observe additional jets of third order J3 and G3 which form from J2-J2 and J1-J2 pairs, respectively. Surprisingly, G2 and G3 have comparable amplitude even though G2 jets are second order whereas G3 only arise in third order; this seems to be due to the large amplitude of J2 jets (from which G3 are formed) in the simulation.

**Fig. 2: Numerical simulation of jet formation.**

To quantify the size of jets, we calculate the integrals over an interval around each peak in the density profile in momentum space ∣ψ(q)∣², which gives the fraction of the atoms in each type of jets. For J0 we integrate from −0.56 to 0.56 and for other jets we integrate an area ± 0.16 around each jet. The total integral is normalized, ∫∣ψ(q)∣²dq = 1, and ψ(q) is the Fourier transform of ψ(x). The fraction in the central peak J0 only starts decreasing after the threshold amplitude for J1 formation is reached, as can be seen in Fig. 2c. The J1 first increases steeply with amplitude and then levels out due to J0 depletion and exitations to the higher-order jets. It can be seen that only the first order has a sharp threshold (Fig. 2c), while the higher orders increase smoothly from zero once the first order exists (Fig. 2d). D2 are not observed even in simulations due to the onset of a wide background at high amplitudes needed to create them (discussed below).

Based on calculated atom fractions, we construct a phase diagram of visibility of different sets of jets, shown in Fig. 3a. A new region is declared to exist when the fraction for the next jet exceeds 0.006. We choose this value so that it is much larger than the noise in the simulation, but smaller than the atom fraction of the jets at high modulation amplitudes. As the amplitude k_ac increases, additional jets appear, represented by new regions in the phase diagram. As discussed above, the transitions are actually continuous and not discrete, except from the “no jets” to the “J1” region. G3 jets occur in parts of the region marked with “J1+J2+J3+G2″ in Fig. 3a.

**Fig. 3: Comparison of theoretical and experimental phase diagrams of jet formation.**

Figure 3b shows the final density profiles in momentum space for a fixed frequency ω = 33. The onset of J1, J2 and other jets can be clearly seen as the modulation amplitude increases. For very high amplitudes the discrete jets start disappearing in the increasingly uniform background of a multitude of higher-order jets. The observed spread from the single momentum mode to a wide momentum distribution resembles the generation of supercontinua in nonlinear optical media, where various instabilities cause an initially narrow frequency spectrum to disperse into a wide continuum²⁶. The boundary for the supercontinuum in Fig. 3a is determined with an entropy-like measure (S=-int| psi (q){| }^{2}log left(| psi (q){| }^{2}right)dq), where we integrate over the whole interval. We scale the threshold with the square root of the modulation frequency, to take into account the scaling of the jet velocities with modulation frequency. When (S > 0.38sqrt{omega }) we say that it has transitioned into the supercontinuum region (the transition is actually gradual, the threshold value was chosen to be above numerical noise and visually without separate jets). The threshold boundary for J1 formation is proportional to ω^1/2, as previously reported^21,24.

Experimental phase diagram of jet formation

We explore the phase diagram of jet formation experimentally as well. Figure 3c, e show the experimental phase diagrams of jet formation for modulation times t_p = 10 ms and 30 ms, respectively. We determined the regions of the phase diagram where we observe J1, J1 and J2 (J1+J2), and J1, J2 and G2 (J1+J2+G2).

In contrast to simulation, in the experiments, the number of atoms in jets varies from shot to shot; for example, see Fig. 1a. This was previously observed in 2D²⁷, where the distribution of population in angular modes resembles the distribution of thermal radiation. More information on the statistics of the number of atoms in each type of jet is available in the Supplementary Note 2 and Supplementary Fig. 3. Because of this variance, gathering statistics is necessary to determine the presence of the phase at each point in the phase diagram. We choose the boundary to be at the point where at least one image out of ten shows the corresponding jets. Examples of averages of 50 absorption images of jets from different parts of the phase diagram are shown in Fig. 3d, f and in the Supplementary Fig. 1. Both Fig. 3d, f show experimental images for modulation frequency 3 kHz, but (f) is for the shorter modulation time t_p = 10 ms and (d) is for t_p = 30 ms. Note that in Fig. 3d the distance between the jets and the initial condensate varies for different modulation amplitudes. This is because jets are emitted sooner if the modulation amplitude is higher (see Supplementary Fig. 6). This effect is not as pronounced in Fig. 3f due to the shorter modulation time.

The range of frequencies in the experimental and the theoretical phase diagram are similar (the frequencies shown in Fig. 3a range from 0 to 7.2 kHz (from 0 to 80 in units of ω_r = 2π ⋅ 90 Hz), whereas the experimental phase diagrams in Fig. 3c, e were measured for frequencies between 1 kHz and 7 kHz). On the other hand, the values of k_ac required for excitation of jets are approximately two times lower in the experiment than in simulation. The difference is most probably because the simulation describes a true 1D system, while in the experiment, the BEC is a 3D object with jets confined into quasi-1D system by a dipole trap beam.

For t_p = 30 ms, we are not able to observe golden jets for modulation frequencies higher than 3 kHz, even though they are observable using the shorter t_p = 10 ms. In the simulation, the modulation time is t_p = 8π, which corresponds to 44 ms. In the experiment, we used shorter modulation times, because longer modulation times cause more atom losses and increase the incoherent background, formed by the excited atoms that are not part of the jets²⁴, making the occurrence of jets rarer and more difficult to observe (see Supplementary Figs. 1 and 2). However, we see from the experimental phase diagrams that the shape of the threshold does depend on t_p. For a short t_p there is no strong dependence on the modulation frequency, whereas for the longer t_p the threshold boundary for the formation of jets is proportional to ω^1/2, as expected from the simulation and previous experiments²⁴. Additionally, the longer modulation times sometimes cause the jets to be emitted twice, making the interpretation of the absorption images less reliable.

As we increase the amplitude of modulation in the experiment, the total number of atoms starts decreasing and the incoherent background between the jets becomes stronger. This can be seen in Fig. 3c and e, where we plot the thresholds at which the total number of atoms is 25% and 50% smaller than the number of atoms in the initial BEC. The decrease in the number of atoms is likely due to partial collapses of the BEC, caused by the strong attractive interactions the BEC is subjected to during the modulation, leading to three-body losses. The incoherent background is formed by the excited atoms that are not part of the jets. Their velocities have random directions and that appears as a homogeneous background between J0 and J1, and to a lesser extent between J1 and J2 (best seen in Supplementary Fig. 2). The momentum distribution of the incoherent background remains of the same width for increasing amplitude and does not increase beyond the position of J2 jets as the supercontinuum would. Atom losses and increased fraction of incoherent background, both not accounted for in the simulation, decrease the effective modulation amplitude and prevent the supercontinuum from forming.

To summarize, the experimental phenomena that the 1D GPE simulation does not account for are the atom losses into the incoherent background and three-body losses. For a complete description of the experiment, one would not only need to take into account the three-dimensionality of the experiment, but also the higher order interaction between atoms (e. g. three-particle interaction) and incoherent processes, that are not described by the mean field picture. Despite all of these differences, the discrepancy between the simulation and the experiment only becomes significant at higher modulation amplitudes.

Discussion

In this work, we show the formation of matter-wave jets with golden ratio dimensionless velocity ({u}_{{{{rm{G2}}}}}=(1+sqrt{5})/2), which form when a BEC in a quasi-1D geometry is subjected to a single-frequency interaction modulation. These jet velocities and the corresponding density waves are not commensurate with the previously observed jets with dimensionless velocities u_J1 = 1 and u_J2 = 2. A hint of incommensurate (sqrt{3}) jets has previously been observed in similar experiments with 2D geometry²³, but, to the best of our knowledge, not in 1D. It needs to be noted that the occurrence of incommensurate density waves can easily be achieved with two-frequency driving (like in ref. ²⁸, but with the frequencies that produce first-order jets that are not commensurable, see Supplementary Note 3 and Supplementary Fig. 4), but in the case of golden jets incommensurability emerges from single-frequency driving.

Numerically and experimentally, we explore the formation of incommensurable jets for a wide range of modulation frequencies and amplitudes, resulting in a phase diagram of possible jets, while also uncovering an unexpected supercontinuum region. The supercontinuum likely arises due to the double-exponential growth of the number of possible modes, while the wavenumber only grows linearly, resulting in a densely populated momentum space. The supercontinuum is potentially experimentally observable if the incoherent background could be reduced by increasing the radial confinement and suppressing three-body losses.

The demonstrated formation of incommensurable matter-wave jets could, for example, be used in multispecies experiments. The incommensurable density waves inside one species could be used as a disordered potential for the other species (provided there exist appropriate Feshbach resonances). This would eliminate the need for an external disordered potential or a regular lattice that is usually needed to capture impurities, which act as a source of disorder^{12,13,14,15,16,17,18,19,20}.

Methods

Preparation of the BEC

We start the preparation of the BEC with a collimated beam of atoms slowed down using a Zeeman slower. The slowed atoms are captured by a magneto-optical trap in the main experimental chamber. We usually load ~ 60 million atoms with temperature ~ 70 μK into the magneto-optical trap over 10 s. After loading we compress the atoms and use optical molasses to cool the atoms to ~13 μK and transfer them to the internal ground state F = 3. The atoms are then transferred into the Raman lattice. With degenerate Raman sideband cooling we cool about 20 million atoms in the absolute ground state (leftvert F=3,{m}_{F}=3rightrangle) to 1 μK, which we then transfer into the large crossed dipole trap. After 0.5 s of plain evaporation in the dipole trap we start ramping up the dimple trap. The dimple beams are linearly ramped to their maximum powers in 1.5 s and held for another 400 ms to let the atoms thermalize. Then the large dipole trap is turned off and evaporation in the dimple trap begins. The powers of dimple beams are exponentially ramped down in 6 s and after 2 s we start with additional exponential rampdown of the magnetic field gradient. The combination of the two techniques leads to efficient evaporation and at the end of it we produce a BEC of approximately 4000 atoms. More details can be found in refs. ^25,29.

Numerical simulation

We use the split-step spectral method to solve the Gross-Pitaevskii equation

$$ifrac{partial psi }{partial t}=-frac{1}{2}frac{{partial }^{2}psi }{partial {x}^{2}}+2k(t)| psi {| }^{2}psi ,$$

(3)

in the dimensionless form³⁰, where time t is given in units 1/ω_r, ω in units ω_r, and distances in units ({a}_{r}=sqrt{hslash /m{omega }_{r}}). k is the interaction parameter k(t) = Na(t)/a_r. We solve the equation on an uniform 1D grid of 65536 points in the interval [−1200:1200]. The initial state of the BEC is determined by an imaginary-time evolution, while the calculation of the real-time dynamics is performed with a time-step of τ = 0.001. The code is written in Wolfram Language (Mathematica) and it can solve for the ground state through imaginary-time propagation of an initial approximation, and determine the time-dependence in real-time. Most of the computation time is spent in performing the (inverse) Fourier transform, which is performed highly efficiently in Mathematica, thus the performance of the code is comparable to that of compiled codes while it is significantly more flexible. Further details and examples can be found in ref. ³¹.

Calibration of modulation amplitude

We measured the amplitude of the modulation of interaction by measuring the splitting of a Feshbach resonance at 19.9 G³². In a typical measurement of Feshbach resonance, we measure the number of atoms remaining in the trap after switching the magnetic field to a chosen value for 100 ms. Around the magnetic field corresponding to the Feshbach resonance we observe enhanced atom losses. We repeat the same measurement with the addition of an oscillating magnetic field at the various frequencies used for jet formation. In this case, we observe that the dip in the number of atoms is widened (see Supplementary Note 4 and Supplementary Fig. 5). The width of the dip corresponds to the peak-to-peak amplitude of the magnetic field oscillations. For jet formation we modulate the magnetic field around 17.1 G, where the wide s-wave Feshbach resonance of cesium has the zero crossing³². The interaction as a function of the magnetic field around the zero crossing is known³², and we use it to determine the amplitude of the interaction modulation. By measuring the effect on the atoms, we directly measure the magnetic field at the position of the atoms. The measurements are done for the same amplitudes and frequencies as for jet formation, meaning that all the non-nonlinearities and other effect of the electrical circuit producing the oscillating magnetic field are taken into account, including the amplifier (Lab Gruppen IPD 1200).