Excitation modes of bright matter wave solitons
Abstract
We experimentally study the excitation modes of bright matter wave solitons in a quasionedimensional geometry. The solitons are created by quenching the interactions of a BoseEinstein condensate of cesium atoms from repulsive to attractive in combination with a rapid reduction of the longitudinal confinement. A deliberate mismatch of quench parameters allows for the excitation of breathing modes of the emerging soliton and for the determination of its breathing frequency as a function of atom number and confinement. In addition, we observe signatures of higherorder solitons by adjusting the size of the wave packet before the quench, and by following the subsequent time evolution and shape oscillations. Our experimental results are compared to analytical predictions and to numerical simulations of the onedimensional GrossPitaevskii equation.
The dispersionless propagation of solitary waves is one of the most striking features of nonlinear dynamics, with multiple applications in hydrodynamics, nonlinear optics and broadband longdistance communications T. Dauxois (2006). In particular, solitons with onedimensional (1D) propagation, such as in optical fibers, have multiple technical applications Haus and Wong (1996). In fiber optics, 1D “bright” solitons, i.e. solitons presenting a local electric field maximum, have been observed Mollenauer et al. (1980). They exhibit a dispersionless flow and excitation modes such as breathing or higherorder modes, which influence the propagation of the wave packets Satsuma and Yajima (1974); Mollenauer et al. (1980); Stolen et al. (1983). Matter waves can also display solitary dispersion properties. Typically, bright matter wave solitons are created in quasi1D systems by quenching the particle interaction in a BoseEinstein condensate (BEC) from repulsive to attractive Khaykovich et al. (2002). Recent experiments demonstrated the collapse Donley et al. (2001), collisions Nguyen et al. (2014), reflection from a barrier Marchant et al. (2013), and the formation of trains Strecker et al. (2002); Nguyen et al. (2017) of bright matter wave solitons.
In this letter, we experimentally study the excitation modes of a single bright matter wave soliton. In previous studies, other dynamical properties have been observed, such as centerofmass oscillations of solitons in an external trap Nguyen et al. (2014), excitations following the collapse of attractive BECs Donley et al. (2001); Cornish et al. (2006), and quadrupole oscillations of attractive BECs in three dimensions (3D) Everitt et al. (2015). Here, we probe the fundamental breathing mode of a single soliton by measuring its oscillation frequency and the time evolution of its density profile. In addition, we observe signatures of higherorder matter wave solitons. Higherorder solitons can be interpreted as stable excitations with a periodic oscillations of the density profile and phase, or as a bound state of overlapping modes Satsuma and Yajima (1974); Carr and Castin (2002). Unlike trains of separate solitons, a higherorder soliton is a single, connected wave packet with a coherent evolution and with a constant total size. We generate higherorder solitons by matching the initial size of the BEC before the interaction quench to values close to those required for second and thirdorder solitons.
The shapepreserving evolution of a matter wave is due to a balancing of dispersive and attractive terms in the underlying 3D GrossPitaevskii equation (GPE) Dalfovo et al. (1999). For quasi1D systems with tight radial confinement, we can approximate the matter wave in the 3DGPE by the product of a Gaussian wave function for the radial direction and a function for the longitudinal direction (see Supplementary Material Sup ()). Depending on the ansatz for the Gaussian with either constant or varying radial sizes, satisfies either the 1DGPE or the nonpolynomial Schrödinger equation (NPSE) Salasnich et al. (2002). We reference to the analytical solutions of the 1DGPE in the manuscript, but use both equations in our numerical simulations Sup ().
For the 1DGPE, solutions for the normalized longitudinal wave function are of the form
(1) 
with a single parameter, , that determines both the longitudinal size and the amplitude of the soliton. Solitons form with a value of that minimizes the total energy and that provides a compromise between the kinetic and the interaction energies. This is illustrated in Fig. 1b, which shows the energy of the wave packet for varying sizes Parker et al. (2007). The kinetic energy provides a potential barrier for small that prevents the collapse of the soliton, while its spreading is inhibited by the interaction energy, which increases for large .
Even without an external longitudinal potential, the soliton is stable against small perturbations of . In a way, a bright matter wave soliton creates its own trapping potential, which defines its size and excitation modes. Variational methods provide accurate predictions of its size at the energy minimum which can be calculated analytically Carr and Castin (2002) or numerically Parker et al. (2007). For the fundamental solution (order ) of the 1DGPE with an atom number , swave scattering length and radial trapping frequency , the size corresponds to the healing length at the peak density of the soliton, i.e. Carr and Castin (2002); Parker et al. (2007). Here, is the radial harmonic oscillator length. Small deviations of close to the energy minimum lead to oscillations of the soliton size. We use those oscillations resulting from an initial mismatch of to experimentally measure the selftrapping frequency of the soliton potential.
Our experimental sequence starts with a BoseEinstein condensate of cesium (Cs) atoms in the state with a condensate fraction of approximately 70% and a scattering length of , where is Bohr’s radius. The BEC is generated in a crossed optical dipole trap formed by the horizontal and vertical laser beams and as illustrated in Fig. 1a. To compensate the strong gravitational force for Cs atoms, we levitate the atoms with a magnetic field gradient that pushes the atoms upwards and balances gravity. An additional magnetic offset field allows us to tune the scattering length by means of a broad magnetic Feshbach resonance with a zero crossing at 17.119 G Berninger et al. (2013). Details about our experimental setup and the levitation scheme can be found in reference Di Carli et al. (2019).
The solitons in our setup are confined to a quasi1D geometry with almost free propagation along the horizontal direction and strong radial confinement of Hz provided by laser beam . In quasi1D geometry, bright matter wave solitons collapse for large densities and interactions Donley et al. (2001), which for our typical experimental scattering length of approximately corresponds to a critical atom number of 2500 Carr and Castin (2002). As a result, we need to strongly reduce the atom number to avoid collapse, modulation instabilities Strecker et al. (2002) and threebody loss Everitt et al. (2017) for a deterministic and reproducible creation of the soliton. We remove atoms with a small additional magnetic field gradient, which pushes the atoms over the edge of the optical dipole trap. Our precise control of magnetic field strengths allows us to reduce the atom number down to atoms, with a reproducibility of for 600 atoms and for 4500 atoms, measured as the standard deviation of the atom number in 50 consecutive runs. A removal period of 4 s and smooth ramps of the magnetic field strength are necessary to minimize excitations of the BEC. Following the removal procedure we measure residual fluctuations of the width of the BEC below 3.5%.
We generate the matter wave soliton with a quench of the scattering length towards attractive interaction (), and by a reduction of the longitudinal trap frequencies (). When changing and independently, the quenches excite inward and outward motion respectively. Usually, it is desirable to minimize the excitations of the soliton by matching the initial ThomasFermi density profile of the BEC closely to the density profile of the soliton (see inset Fig. 1a) and by adjusting the quench parameters. However, we deliberately mismatch the quench parameters to create breathing oscillations of the soliton in order to study its selftrapping potential. The cloud size at a hold time after the quench is measured by fitting the function , with fitparameters and , to the density profiles. The density profiles are determined from absorption images after a free expansion time of 16 ms (Fig. 1c).
The response of the atomic cloud to the different steps of the quench protocol is presented in Fig. 1d. We first quench only the longitudinal confinement for atoms from Hz to Hz in 4 ms while keeping the interaction strength constant (blue diamonds in Fig. 1d). The BEC starts an outwards motion with an oscillation frequency of , which is well expected for a BEC in the ThomasFermi regime Menotti and Stringari (2002); Haller et al. (2009). In a second quench protocol, we additionally quench the interaction strength to (green squares in Fig. 1d). In order to minimize excitations of the soliton, we set Hz to match the width of the BEC to the expected size of the soliton. As a result, we observe almost dispersionless solitons with a linear increase of of (green line in Fig. 1d). Finally, we deliberately mismatch the size of the BEC with Hz, and generate small amplitude oscillations of the soliton with a frequency of Hz (red circles in Fig. 1d). This breathing frequency of the soliton is significantly larger than any possible breathing frequency of a BEC, and larger than the breathing frequency of noninteracting atoms of Hz. We observe no discernible oscillation in the radial direction after the quenches.
In a second experiment, we demonstrate that the breathing frequency, , depends on the interaction strength, a property typical of the nonlinear character of the soliton. The interaction term in the 1DGPE depends on the product , and we expect that the variation of scattering length and atom number have equivalent effects on . We choose to change , since the initial removal process is independent of the interaction quench, and it allows us to study the breathing frequency without changing the quench protocol and without excitations of additional modes. We measure the breathing frequency, , for varying with fixed parameters Hz, Hz, , , Hz (see red circles in Fig. 2). The values of decrease for lower , and they approach the breathing frequency for noninteracting atoms in a harmonic trap (dashed gray line).
We compare our experimental data points to two theoretical models. For a numerical simulation of the 1DGPE, we use the ansatz in Eq. 1 to calculate the starting conditions, and we determine the breathing frequency from a spectral analysis of the time evolution of the wave function (blue triangles in Fig. 2a) Sup (). An analytical approximation for the breathing frequency (red line) can be calculated with a Lagrangian variational analysis at the energy minimum of the 3DGPE Carr and Castin (2002); Sup ()). We find that both models agree well with the trend of the measurements of , although our experimental data points are systematically lower for large than our theoretical predictions. We speculate that this is due to nonharmonic contributions to the energy of the soliton on the breathing oscillations for finite oscillation amplitudes (see Fig. 1b).
To determine the influence of the trapping potential, we measure the variation of as we reduce the longitudinal trapping frequency . Smaller values of result in larger equilibrium sizes of the soliton, and we need to reduce the initial trap frequencies to keep the oscillation amplitudes comparable during the measurements. The typical difference between and is approximately Hz. Two regimes can be identified in Fig. 2b where the soliton frequency changes with varying the trap frequencies . For large values of the trap dominates the breathing of the soliton and approaches twice the trap frequency . For small values of , interactions dominate the breathing of the soliton and reaches a constant value. This offset of the breathing frequency is a result of the “selftrapping” potential of a free soliton.
Again, we compare the experimental results with our theoretical model (see red line in Fig. 2b) and the numerical simulations of the 1DGPE (see blue band in Fig. 2b for simulations between and ). The simulations predicts a lower breathing frequency for the free soliton than the analytical approximation, but all curves are within the uncertainly range of the experimental data.
Breathing oscillations of close to the equilibrium size are not the only possible excitation modes of solitons. The existence of higherorder solitons has been predicted in the nonlinear Schrödinger equation (NSE) Satsuma and Yajima (1974), and observed for optical solitons in silicaglass fibers Mollenauer et al. (1980); Stolen et al. (1983). A soliton of order can be interpreted as a bound state of strongly overlapping solitons Carr and Castin (2002). By exploiting the equivalence of the NSE and 1DGPE, similar effects were later proposed for bright matter wave solitons Carr and Castin (2002); Golde et al. (2018), where it was suggested that higherorder solitons can be generated by a rapid increase of the attractive interaction strength. However, it is experimentally easier to avoid a second interaction quench by instead increasing the initial size of the soliton. The wave function in Eq. 1 describes a 1DGPE soliton of order with a constant size . Higherorder solitons form for the same initial shape of the wave function, but with an initial size that is a specific multiple of the healing length , i.e. for .
Figure 3 shows numerical simulations of the 1DGPE for the time evolution of second and thirdorder solitons with initial sizes and (see Sup ()). Large initial soliton sizes lead to the periodic formation of local maxima and minima of the density profile. Striking characteristics of the time evolution are the periodic development of a sharp central peak with side wings for the secondorder soliton (Fig. 3a,b), and the periodic formation of a broad doublepeak structure for the thirdorder soliton (Fig. 3c,d). Other initial sizes, different from , lead to a “shedding” of the atomic density in the direction. The wave packet loses particles until its size and shape match the next (lower ) higherorder soliton Satsuma and Yajima (1974).
The experimental observation of higherorder matter wave solitons is challenging for several reasons. The available parameter regime for , , and is limited by avoiding collapse or splitting of the soliton, by our imaging resolution of m, and by the observation time. For the creation of secondorder (thirdorder) solitons, we select , (), Hz, Hz, and we adjust the initial size of the BEC before the quench by tuning the trap frequency and the scattering length. In addition, we raise the quench duration to ms to increase reproducibility and to further reduce residual excitations. As a result, there is an uncertainty in the initial size when we generate solitons of secondorder (see Fig. 4a,b) and thirdorder (see Fig. 4c,d). We measure an initial size of m (m) directly after the quench with a freeexpansion time of ms before imaging. The values are close to the expected values of m (m).
The absorption images show the periodic formation of a central peak (Fig. 4a), and the additional formation of a doublepeak structure (Fig. 4c) for the two cases close to the second and thirdorder soliton. In agreement with the theory, we find that the evolution of the density profile of higherorder solitons is susceptible to small perturbations and to small changes of the starting conditions. Our measured density profiles are typically well reproducible up to the first minimum of the soliton width, but show increasing fluctuations for longer observation times. In Fig. 4a,c we postselected the images of the density profiles to illustrate the shape oscillations. For a more reliable analysis that does not rely on shape recognition and on image selection, we choose a parameter to quantify the width of the wave packet and we compare the evolution of the parameter to 1DGPE simulations.
The periodic change of the soliton shape can be well identified by evaluating the second moment of the density distribution , i.e.
(2) 
where is the mean position of the wave packet. The time evolution of shows a minimum for the sharp peak of the secondorder soliton and a local maximum for the double peak structure of the thirdorder soliton (Fig. 3b,d). We clearly identify the features of the secondorder soliton for our experimental data in Fig. 4b, with minima of close to 150 ms and 425 ms and a maximum close to 270 ms. The experimental data agree well with the numerical simulation of the 1DGPE using the initial size m as the only adjustable parameter.
Similarly, the time evolution of in Fig. 4d shows maxima due to the double peak structures of the density profiles. We find that the shape of the oscillations and the period match best to a simulation with an initial size of m, which is smaller than expected, and which results in simulated values of that are smaller than those of the observations. The deviation may be caused by the finite quench duration, by a residual trapping potential, or by the excitation of additional breathing oscillations which superimpose the shapeoscillation of the higherorder solitons. In addition, the excitation of a pure thirdorder soliton is difficult to achieve, and we speculate that the measurements reported in Fig. 4d correspond to an intermediate state with a possible final relaxation to a lower order solitary wave.
In conclusion, we experimentally studied the creation and the excitation of breathing modes of bright matter waves solitons in quasionedimensional geometries. We formed the solitons by quenching a BEC of cesium atoms from repulsive to attractive interactions and by a rapid reduction of the longitudinal confinement. A deliberate mismatch of quench parameters allowed us to excite breathing modes of the emerging solitons, and to measure the “selftrapping” frequency . In addition, we observed signatures of second and thirdorder solitons by adjusting the initial size of the wave packets before the quench. A detailed understanding of the excitation modes of bright solitons is essential to apply matter wave solitons for interferometric precision measurements Cuevas et al. (2013); Polo and Ahufinger (2013); McDonald et al. (2014); Helm et al. (2015).
The authors would like to thank L. D. Carr for helpful initial discussions during his visit. We acknowledge financial support by the EU through the Collaborative Project QuProCS (Grant Agreement 641277). GLO acknowledges support from the European Training network ColOpt (European Union Horizon 2020 program). AdC acknowledges financial support by EPSRC and SFC via the International MaxPlanck Partnership. The work was also supported in part by the EPSRC Programme Grant DesOEQ (Grant No EP/P009565/1).
30
 T. Dauxois (2006)

M. P. T. Dauxois, Physics of Solitons (Cambridge University Press, Cambridge, 2006).
 Haus and Wong (1996) H. A. Haus and W. S. Wong, Rev. Mod. Phys. 68, 423 (1996).
 Mollenauer et al. (1980) L. F. Mollenauer, R. H. Stolen, and J. P. Gordon, Phys. Rev. Lett. 45, 1095 (1980).
 Satsuma and Yajima (1974) J. Satsuma and N. Yajima, Prog. Theor. Phys. Suppl. 55, 284 (1974).
 Stolen et al. (1983) R. H. Stolen, W. J. Tomlinson, and L. F. Mollenauer, Opt. Lett. 8, 186 (1983).
 Khaykovich et al. (2002) L. Khaykovich, F. Schreck, G. Ferrari, T. Bourdel, J. Cubizolles, L. D. Carr, Y. Castin, and C. Salomon, Science 296, 1290 (2002).
 Donley et al. (2001) E. A. Donley, N. R. Claussen, S. L. Cornish, J. L. Roberts, E. A. Cornell, and C. E. Wieman, Nature 412, 295 (2001).
 Nguyen et al. (2014) J. H. V. Nguyen, P. Dyke, D. Luo, B. A. Malomed, and R. G. Hulet, Nat. Phys. 10, 918 (2014).
 Marchant et al. (2013) A. L. Marchant, T. P. Billam, T. P. Wiles, M. M. H. Yu, S. A. Gardiner, and S. L. Cornish, Nat. Commun. 4, 1865 (2013).
 Strecker et al. (2002) K. E. Strecker, G. B. Partridge, A. G. Truscott, and R. G. Hulet, Nature 417, 150 (2002).
 Nguyen et al. (2017) J. H. V. Nguyen, D. Luo, and R. G. Hulet, Science 356, 422 (2017).
 Cornish et al. (2006) S. L. Cornish, S. T. Thompson, and C. E. Wieman, Phys. Rev. Lett. 96, 170401 (2006).
 Everitt et al. (2015) P. J. Everitt, M. A. Sooriyabandara, G. D. McDonald, K. S. Hardman, C. Quinlivan, M. Perumbil, P. Wigley, J. E. Debs, J. D. Close, C. C. N. Kuhn, and N. P. Robins, (2015), arXiv:1509.06844 .
 Carr and Castin (2002) L. D. Carr and Y. Castin, Phys. Rev. A 66, 063602 (2002).
 Dalfovo et al. (1999) F. Dalfovo, S. Giorgini, L. P. Pitaevskii, and S. Stringari, Rev. Mod. Phys. 71, 463 (1999).
 (16) See Supplemental Material for additional information .
 Salasnich et al. (2002) L. Salasnich, A. Parola, and L. Reatto, Phys. Rev. A 65, 043614 (2002).
 Parker et al. (2007) N. G. Parker, S. L. Cornish, C. S. Adams, and A. M. Martin, J. Phys. B At. Mol. Opt. Phys. 40, 3127 (2007).
 Berninger et al. (2013) M. Berninger, A. Zenesini, B. Huang, W. Harm, H.C. Nägerl, F. Ferlaino, R. Grimm, P. S. Julienne, and J. M. Hutson, Phys. Rev. A 87, 032517 (2013).
 Di Carli et al. (2019) A. Di Carli, C. Colquhoun, S. Kuhr, and E. Haller, New J. Phys. (2019), 10.1088/13672630/ab1bbd.
 Everitt et al. (2017) P. J. Everitt, M. A. Sooriyabandara, M. Guasoni, P. B. Wigley, C. H. Wei, G. D. McDonald, K. S. Hardman, P. Manju, J. D. Close, C. C. N. Kuhn, S. S. Szigeti, Y. S. Kivshar, and N. P. Robins, Phys. Rev. A 96, 041601 (2017).
 Menotti and Stringari (2002) C. Menotti and S. Stringari, Phys. Rev. A 66, 043610 (2002).
 Haller et al. (2009) E. Haller, M. Gustavsson, M. J. Mark, J. G. Danzl, R. Hart, G. Pupillo, and H.C. Nägerl, Science 325, 1224 (2009).
 Golde et al. (2018) J. Golde, J. Ruhl, M. Olshanii, V. Dunjko, S. Datta, and B. A. Malomed, Phys. Rev. A 97, 1 (2018).
 Cuevas et al. (2013) J. Cuevas, P. G. Kevrekidis, B. A. Malomed, P. Dyke, and R. G. Hulet, New J. Phys. 15, 063006 (2013).
 Polo and Ahufinger (2013) J. Polo and V. Ahufinger, Phys. Rev. A  At. Mol. Opt. Phys. 88, 1 (2013).
 McDonald et al. (2014) G. D. McDonald, C. C. N. Kuhn, K. S. Hardman, S. Bennetts, P. J. Everitt, P. A. Altin, J. E. Debs, J. D. Close, and N. P. Robins, Phys. Rev. Lett. 113, 013002 (2014).
 Helm et al. (2015) J. L. Helm, S. L. Cornish, and S. A. Gardiner, Phys. Rev. Lett. 114, 1 (2015).
 Cowell et al. (2002) S. Cowell, H. Heiselberg, I. E. Mazets, J. Morales, V. R. Pandharipande, and C. J. Pethick, Phys. Rev. Lett. 88, 210403 (2002).
 Suarez (2015) P. Suarez, “An introduction to the split step fourier method using matlab.” (2015).