Dipolar dark solitons

Dipolar dark solitons

Krzysztof Pawłowski Center for Theoretical Physics, Polish Academy of Sciences,
Al. Lotników 32/46, 02-668 Warsaw, Poland
   Kazimierz Rzążewski Center for Theoretical Physics, Polish Academy of Sciences,
Al. Lotników 32/46, 02-668 Warsaw, Poland

We numerically generate, and then study the basic properties of dark soliton-like excitations in a dipolar gas confined in a quasi one dimensional trap. These excitations, although very similar to dark solitons in a gas with contact interaction, interact with each other and can form bound states. During collisions these dipolar solitons emit phonons, loosing energy, but accelerating. Even after thousands of subsequent collisions they survive as gray solitons and finally reach dynamical equilibrium with background quasiparticles. Finally, in the frame of classical field approximation, we verified, that these solitons appear spontaneously in thermal samples, analogously to the type II excitations in a gas of atoms with contact interaction.

03.75.Hh , 03.75.Lm , 05.45.Yv

I Introduction

The main motivation to our work comes from recent progress concerning the role of dark solitons in the physics of the gas with contact interaction only. First, let us remind, that solitons are localized waves propagating through non-linear media without any spreading and also moving through each other without changing the shape. Among equations of broad interest for physicists, there are a few which exhibit solitons. One of them is the non-linear one dimensional Schrödinger equation:


in which one can find bright (the case with attracting forces, ) and dark (the repulsive case ) solitons.

The Eq. (1) applies, among others, to the dynamics of bosons interacting via a -like potential, i.e. the interaction potential between two bosons is equal to . In this context, the Eq. (1) is only an approximation valid in the variously defined weak interaction limit Lieb et al. (2005). Such a gas of bosons, can be described by more fundamental many-body model, the Lieb-Liniger model, which has been used to find rigorously the ground state Lieb and Liniger (1963) as well as the elementary excitations Lieb (1963). It turns out, that there are two branches of elementary excitations, staying on the same footing. The first branch has been identified as the Bogoliubov excitations. The second one, consisting of the so called "type II" excitations, within research performed by the next generations of physicists turned out to be solitonic branch Kulish et al. (1976); Kolomeisky et al. (2000); Jackson and Kavoulakis (2002); Kanamoto et al. (2010); Karpiuk et al. (2015); Syrwid and Sacha (2015). It is now clear that the dark solitons, as the elementary excitations, appear spontaneously in a gas at finite temperature Karpiuk et al. (2012).

The equation (1), is not capturing many important and common situations, even in the weakly interacting limit. For example, it is neglecting possible external traps in which the gas is confined, and it is assuming contact interaction only. In particular, it misses the gas with a dipolar interaction, which is recently studied in more and more laboratories Griesmaier et al. (2005); Aikawa et al. (2012); Lu et al. (2011); Ni et al. (2008). Among many motivations to study the dipolar gas is definitely its unusual excitation spectrum, which in certain situations, close to instabilities, can have local minimum for non-zero quasi momenta, so called roton spectrum Santos et al. (2003). It seems natural to ask if here is a place for the second branch and from what excitations it can be build. Up to our knowledge, there is no study of the elementary excitation spectrum beyond the Bogoliubov approximation, neither by solving the underlying many body problem, nor on the mean field level.

Moreover, although the literature about the bright solitons is very reach (see parts of the review Lahaye et al. (2009)), relatively little is known about dark solitons. For the time being the shapes of dipolar solitons in a gas with both contact and dipolar interaction has been presented Adhikari (2014); Cuevas et al. (2009). Surprisingly, the life time of the dark solitons can be increased, with some tricks even in the higher dimensional space Nath et al. (2008). The dark solitons in non-local nonlinear media were the subject of research in optics. It has been found that two optical dark solitons interact with each other via non-trivial potential, they can even form a bound state Kong et al. (2010). In this context, although the effective forces between photons can be non-local, it is always assumed that they are not long-range.

Within this work, limited to the mean field description, we look closer at the basic properties of the dark solitons in the dipolar media. The aim of this paper is to study dark dipolar solitons and contrast them with the true solitons existing in the gas with contact interaction only. The analysis is performed within the mean field description, but we include also the finite temperature case. We start with the detailed description of our model in Sec. II. After presentation of the basic features of single soliton (shapes, dynamics) in Sec III, we discuss a motion of two dipolar solitons in Sec. IV. Sec. V is devoted to the analysis of the thermal samples.

Ii The Model

Figure 1: (Color online) Configuration: Trapping potential is given by a steep harmonic trap in the and the directions (called transverse directions), and a long box with periodic boundary conditions in the direction. The dipole moments of the constituent atoms are polarized along axis, such that they mostly repel each other.

We assume that the dipolar gas is confined in the potential


with a tight harmonic trap in the and the directions. The space in the direction is assumed to be finite, with the length and with periodic boundary conditions. These conditions impose a quantization of the momenta in this direction, is a multiple of .

Throughout the paper we use the box units, where , , and are the units of length, time, energy and temperature, respectively. We assume that the dipole moments of the constituent atoms are polarized along the axis , namely when the dipolar potential in the momentum space has the form 111We give the dipolar interaction potentials in the momentum space because this form is more useful in numerical treatment.


where is a dipolar coupling strength. Thus, the atoms move on the circumference of a circle, having the dipole moments perpendicular to the circle-plane. In this configuration dipoles mostly repel each other, what favors the dark solitons. The geometry of our problem is sketched in Fig. 1.

We will investigate the system within the mean field model. The time-dependent mean field wave-function obeys the three dimensional Gross-Pitaevskii equation (GPE)

We assume that the transverse confinement is so tight, that the wave-function stays in the lowest energy level in the and directions, namely when the following single-mode approximation holds:


where is the Gaussian wave function:


where is the transverse oscillator width. The necessary condition underlying this approximation is that both, the thermal energy and the interaction energy, are much below the energy of the first excited state in the transverse direction.


Assuming such a situation we multiply the three dimensional GPE by , integrate the result with respect to and and finally we reach the one dimensional effective GPE


The Fourier transform of the effective one dimensional dipolar potential reads Nath et al. (2008); Cai et al. (2010).


where the constant term with mimics a contact interaction, we introduced one dimensional dipolar coupling strength and is equal to the aspect ratio of the trap. The function which appears in Eq. (8) is equal to , where Ei is the exponential integral. In this paper we aim at studying the role of pure dipolar interactions only, hence from now on, we will assume 222This can be done via Feshbach resonances, setting .

The spatial representation of this potential is visualized in Fig. 2.

Figure 2: (Color online) Effective potential defined in (8) in position representation. For large distances it behaves like , its characteristic range (width at half maximum) is of the order of . Parameters: and .

It is worth stressing that Eq. (7) describes not a movement of single dipoles, but rather the movement of slices of the dipolar gas as marked in Figs. 1 and 2. If such slices are far from each other, then the details of the distribution of dipoles within such slices plays no role, and we can treat them as single large dipoles. In such situations the interaction between them scales as . The situation changes, when the distance between them is of the order of the transverse width . This gives the characteristic range of the effective 1D potential (8). Finally in the limit of slices going to the same place, the effective dipolar potential converges to a constant.

Iii Single solitons

Figure 3: (Color online) Shape of the solitonic-like excitations in the dipolar gas. Results of the imaginary time evolution with the phase constraint given in Eq. (9) are marked with: red solid line for , dashed blue for and dot-dashed green for . In all figures: and . Gray thick line indicates the Shabat-Zakharov dark soliton Zakharov and Shabat (1973) for the gas with contact interaction only, . The blue dashed line corresponds to experimental-like parameters: the mass and the dipole moment of Dysprosium, assuming atoms and the transverse trap equal to kHz.

In the gas with contact interaction only, the phase of a solitonic wave function is changing by around the region with minimal density. This motivates us to look for the solitons in the dipolar gas within the wave function which minimizes the energy but under constraint of a jump of the phase equal to . If we only force the jump of the phase equal to , then the periodic boundary condition would not be satisfied. Hence, except the jump required by the soliton, we add to the constraint on the phase a linearly changing function of the position:


such that the phases at the edges of the box are equal. Due to this additional term the gas is moving with the speed . When showing the result for the dynamics, we would unrotate this motion.

To find the state with the lowest energy but with the phase equal to the constraint (9) we use the method of the imaginary time evolution, forcing at each step the phase (9). Densities computed with this method for three different aspect ratios , and are shown in Fig. 3. We compare these profiles with the density of a dark soliton from purely contact interacting case Zakharov and Shabat (1973) marked in Fig 3 with the thick gray line. As in both cases, a contact interacting gas and the considered here quasiD dipolar gas, the lowest lying excitations are phonons, then comparing their speed of sound we translate the dipolar coupling coefficient into an effective coupling coefficient . We find . The speed of sound in the dipolar media is equal to . The width of the Shabat-Zakharov soliton is proportional to the healing length, , where is the effective scattering length defined by the relation . The Gross-Pitaevskii equation is supposed to be valid in the weakly interacting limit, namely when , what implies .

One can see that the more elongated is the trap, the closer is the result to the soliton in contact interacting gas, marked in Fig. 3 with the thick gray line. To understand this observation mathematically we look at the system in the momentum space. The range of the effective dipolar potential in the momentum space is equal to , equal in the box units to . Hence, in the limit of a very elongated trap, the range of the effective potential is going to , and the evolution is dominated by its values relatively close to . The derivative of the effective potential Eq. (8) at vanishes, what means that the potential is flat at the origin. Thus locally, where this ’local’ region spans even to high momenta, the interaction potential looks like a constant function, as it would be in the contact interacting case. The ratio between the typical width of the soliton and the range of the potential , written in the momentum space and in the box units is equal to . Thus if the trap is very elongated and the coupling strength is small, we expect that the dipolar soliton will be close to the Shabat-Zakharov solution and moreover that our system will evolve similarly to the gas with contact interactions only.

Figure 4: Density as a function of time presented in the frame rotating with the speed . In the laboratory frame the solitons travel around the circle times. Parameters: , and , as in the Dysprosium like case commented in caption of Fig. 3. The horizontal red line is at time : a time at which a single phonon would pass the whole box.

The presented solutions are only candidates for solitons, having similar jump of the phase as the Shabat-Zakharov solitons. The crucial feature of solitons is their ability to propagate without changing their shape. To verify if this is the case also here, we were setting as initial states the results of the imaginary time evolution, and trace their dynamics obtained via numerical solution of GPE. An example is given in Fig. 4, where we unrotated the constant speed due to the linear dependance of the phase. The evolution time in this figure is about times longer than a time at which a phonon would travel across the whole box. When testing the system in a wide parameter range and for a longer evolution time we always observe the dips which propagate without noticable distortion.

Iv Two solitons

Figure 5: Density as a function of time. The initial distance between solitons is equal to . Parameters as in Fig 3.

The real solitons obey an equivalent of the superposition principle: although the system is non-linear one can construct many solitonic solution based on a single one. As shown in the paper of Shabat and Zakharov Zakharov and Shabat (1973) in the contact interacting gas two solitons, when colliding, do not change the shape. To test if the dipolar solitons can behave similarly, we generate and propagate pairs of them. We follow the method described for a single soliton, but this time we enforce two jumps of the phase, at positions and , and without linear background


Having such constraint on the phase we reach with the help of the imaginary time evolution two solitons and then we propagate them using GPE. An example of such dynamics is shown in Fig 5. The immediate conclusion is that the dipolar solitons interact , in the example shown they attract each other. Hence, the excitations we study are not solitons in the strict mathematical sense. Even though, these excitations resemble the shape for a long evolution times, thus we will use the name soliton-like or, for brevity, even the name soliton.

Figure 6: Acceleration versus distance between solitons for different initial states. Parameters as in Fig 3.

To understand the motion of such a pair we monitor the position of one of the solitons versus time: . Then we numerically compute both the velocity and the acceleration . Combining these quantities we obtain acceleration versus inter-solitonic distance , where . We repeat this procedure for many initial states. One can see that results of many simulations are lying on a common curve, see Fig. 6. It suggests us to use the following classical mechanics picture, in analogy to many cases for dark solitons in the contact interacting gas Theocharis et al. (2005); Busch and Anglin (2000) we treat each soliton as a point-like (negative) mass obeying the Newton law. Assuming , we can numerically integrate to extract inter-solitonic potential divided by a mass .

Figure 7: Effective interaction potential between solitons divided by mass in the classical picture. The figures in the bottom panel show density versus time for two initial conditions, for the minimum of the potential (left) and for (right). Parameters: , and .

We show the results for two cases: in Fig. 7 for very elongated trap with and weak dipolar interaction and in Fig. 8 for and . In both cases the inter-solitonic potential has a repulsive core and the tail scaling with the distance as . We tested this classical picture starting again from different initial states and watching their evolution. Samples of this investigation are given at the bottom panels in Figs. 7 and 8. Two solitons placed close to each other strongly repel each other. As the inter-solitonic potential has a local minimum, thus there is a possibility for bound states as shown in the panels marked with the letter A. Two far away lying solitons, like in the panel marked with the letter B, they hardly see each other, but after longer time they finally start oscillations in the common potential.

In Figs. 7 and 8 we compare the resulting potential just with the rescaled density of a single soliton , where is the density far from the soliton, approximately equal to 333We compute numerically, it is not fitting parameter.. The agreement between solitonic shape and inter-solitonic potential is surprisingly good. This indicates the origin of the inter-solitonic potential: each soliton is simply moving in the mean field potential created by its partner. As the dark solitons have negative masses, they are climbing towards the local maxima of the density.

Figure 8: Effective interaction potential between solitons divided by mass in the classical picture. The potential consist of repulsive exponentially decaying core and a tail that decays like . The figures in the bottom panel show density versus time for three initial conditions, at the repulsive core (left), close to the minimum of the potential (center) and for (right). Parameters as in Fig 3.

The attentive Reader could notice some peculiarity in Fig. 9 bottom right panel. There, the solitons seem to slightly accelerate. This counter-intuitive effect is not a numerical artefact.

Figure 9: Density as a function of time. The red solid line indicates the speed of sound. Parameters: , and .

To illustrate clearly what is happening we present the evolution for the extreme conditions, , in Fig. 9. Then one can see that the collision is not elastic: colliding solitons radiate phonons. In consequence the energy of solitons decays, and the dark solitons are transformed into the gray solitons. As in the case with contact interacting gas, the gray solitons have some characteristic speed. This extra velocity is thus gained by the solitons due to the inelastic collision.

Hence, the dissipation during the collision leads in fact to acceleration of solitons, as clearly shown in Fig. 9. Actually this effect, although not directly described, seems to exist for the optical solitons with nonlocal nonlinearity studied in Kong et al. (2010) (Fig. 3 therein).

Figure 10: The grayness defined in the text as a function of time. The timescale is so long, that during evolution there occurs hundreds of the collisions between solitons. The superimposed red line corresponds to travel at the speed of sound. The plot of density vs time corresponds to the evolution at short time window, from to . The red line is a fit in which we invert the energy to grayness relation and assume energy with fitting parameters , and , to model the process of the stimulated emission. Parameters as in Fig 3.

Until collision the movement of a pair of solitons can be captured by a classical model. Then, collision goes clearly beyond the classical picture: it is an event when the part of the energy is transferred from solitons into initially unoccupied excitations, i.e. phonons. What is the fate of these solitons? To find the answer we let the system evolve for a long time. During the evolution we monitored minimal density in the system, called here a grayness . In the case of a gas with contact interactions only, this parameters can be translated into the Shabat-Zakharov soliton speed, equal to , and its energy, equal to Jackson and Kavoulakis (2002). The evolution is shown in Fig. 10. As we started with two black solitons, initially . Then, at short times,the grayness increases, first quadratically and then exponentially. At a very long time the quantity starts to saturate. We understand the dynamics in terms of interaction between the subsystem (solitons) with an environment (all phonons and other quasiparticles), taking the idea from the paper Bach et al. (2002). Initially the environment is empty: this is the regime of spontaneous emission, where a coupling between subsystem and the unoccupied modes leads to a slow dissipation. Then, we enter the regime of stimulated emission, where the modes after the stage of spontaneous emission are already occupied, and hence the rate of phonon emission rapidly grows. If the environment would be infinitely large then eventually we would loose the energy from the subsystem completely. However here, we deal with a small environment, in which the phonons always overlap with solitons. This allows for the dynamical equilibration, in which the amount of the energy radiated to phonons is equal to the energy absorbed, as in the qubit interacting with a thermal reservoir. Indeed we see in the inset of Fig. 10, that even after a long evolution time the two solitons survive. The large fluctuations of are due to overlaps between phonons and solitons, but also due to statistical fluctuation of the flow of the energy between solitons and phonos.

Figure 11: Density as a function of time at three different temperatures. The thick black line depicts a travel with the speed of sound. Parameters: , and and .

V Thermal solitons

As mentioned in the introduction, the ultracold gas with contact interaction possesses not only phonons but also dark solitons Karpiuk et al. (2012). These solitons are the elementary excitations forming the second branch of excitations derived in 1963 by Elliot Lieb Lieb (1963). We would like to check if non-Bogoliubov excitations exist also in a quasi 1D dipolar gas. To fully verify it, one should solve a many-body problem, an equivalent of the Lieb-Liniger model for a dipolar gas. This task is definitely beyond this paper. Instead, still relying on the mean field model, we can at least check if the solitons, as candidates for "type II" excitations in a dipolar gas, appear spontaneously at finite temperatures and if their number obeys the thermal statistics. At the first glance, the presence of such solitons seem unlikely: the system is not integrable and we have already shown that the dipolar soliton-like excitations are not robust, as they collide inelastically. On the other hand, we also demonstrated that two solitons immersed in the gas full of phonons can stay in the dynamical equilibrium. The latter observation encourages us to look also at the thermal equilibrium.

As in Karpiuk et al. (2012), we use the classical field approximation (CFA), in which one replaces the quantum field operator with the classical field . This method is commonly used at zero-temperature, where is interpreted as the macroscopically occupied (by all particles) single particle orbital, called the GPE wave function. The idea of using the same replacement to describe a gas at non-zero temperature relies on the fact that for the degenerate gas there are many substantially occupied single-particle orbitals, which sum up to a classical field .

To compute the values of target variables averaged over a statistical ensemble one performs averaging over the appropriately chosen set of the classical fields or/and time-averaging over the fields propagated with the help of the GPE equation or its generalization.

There is variety of different classical field methods devised to describe ultracold gases publication (2013). Here we use a simple version described in Brewczyk et al. (2007). First, we generate many initial states using a method based on the Bogoliubov description. As our method of generating thermal samples is not sufficiently good, we evolve the field with the help of GPE, we trace the system for sufficiently long time and compute time averages.

More precisely, in each realization, the initial wave-function reads


where the wave vector is with an integer. The sum over wave-vectors is limited by the with , where is the temperature and is the chemical potential. The cut-off is chosen such that with the classical field method one would recover the quantum results at least for the ideal gas Witkowska et al. (2009, 2010). The expansion coefficients should be drawn from the thermal classical distribution. Even on the classical level this drawing can be a cumbersome task. The details of our procedure of drawing initial set of are given in the Appendix. As the procedure gives unsatisfactory effect (especially for higher temperatures), we have to propagate each guessed classical field over a long time to allow for thermalization. Finally, after this auxiliary stage of thermalization, in the next stage we further propagate the state to perform averaging over time. The short-time windows of the evolution during the last stages are illustrated in Fig 11. We see there three examples of the evolution of the density of the gas, with the fraction of the -momentum mode from (top panel), through (middle panel), (bottom panel).

The speed of sound is depicted with the solid black line. In the middle and the bottom panels one can see dips in the density, propagating much slower than the speed of sound. We checked that these dips were stable even at long timescales. Thus we expect that they are dark solitons.

Figure 12: Number of dark solitons versus fraction of atoms in the -momentum mode. Histogram corresponds to the numerically treated dipolar case with parameters , and . The green points are computed analogously to the dipolar case, but for and . Finally the blue dashed line is the expected number of solitons assuming the equipartition of energy and assuming that dark solitons as quasi-particles, namely number of solitons is approximated with , where is the estimated temperature and is the estimated energy of the dark soliton. The details of the estimation are in the text. In each case we assumed atoms. The results for the dipolar gas and the gas with contact interaction comes from time averaging performed on realizations of the classical field, each drawn as described in the Appendix.

To convince ourselves we prepare the following test: for each thermal sample we performed a long evolution, during which we average over time the fraction of the least energetic -momentum mode , the fluctuation of this fraction and the number of the deepest dips, with the density below . The last quantity, called the number of dark solitons, is presented in Fig. 12. We estimated the temperature of each sample using formulas for average number of motionless particles, i.e. , for the ideal gas Wilkens and Weiss (1997). As the crosscheck, we compute what would be fluctuation of the latter average in the ideal gas and compare it with the previously numerically found . The agreement was within a few percent, hence we concluded that: 1) the interactions in the gas are still so small that the statistics seems close to the ideal gas case and 2) our estimation for temperature is reasonable within a few percent. We also estimated the energy of dark solitons, relying on the formulas for the gas with contact interactions only Jackson and Kavoulakis (2002) as done in the previous sections. This gives the estimation for the energy of the dipolar dark soliton (in box units). Finally, having the estimated energy of the dark soliton, and estimated temperature of each sample we computed the estimated number of dark solitons , assuming their equipartition as it should be for quasiparticles. This expected number of solitons is given by the blue dashed line in Fig. 12. We also repeat all the steps for the gas with the contact interaction only.

Vi Conclusions

We investigated the dipolar gas confined in the potential having the shape of a long tube, with the dipole moments polarized perpendicular to the long axis. In this geometry we find robust localized excitations, very similar to the Shabat - Zakharov solitons appearing in the gas with contact interaction only. In contrast to the contact interacting case, pairs of these soliton-like dipolar excitations mutually interact, and are able to form two-soliton bound states. Moreover, the collisions between such pairs are not elastic: during collision they emit phonons, hence they loose energy and because of this they accelerate. After thousands of subsequent collisions they can reach dynamical equilibrium with the already emited quasiparticles. Finally we show that these excitation appear spontaneously at finite temperatures and their number obeys thermal statistics. The latter observation is closely related to the contact interacting case, where the appearance of the dark solitons at the non-zero temperatures is an evidence that they are the non-Bogoliubov elementary excitations, predicted by Elliot Lieb in 1963.

There are many open questions concerning these excitations. One should check their stability in the real three dimensional calculation. It is interesting if they exist also in the roton regime. There is the question of the effect of a harmonic trap along the long axis (here the box with periodic boundary was assumed). One needs also benchmarks from other methods used for the finite-temperature simulations, as they are quite developed for the dipolar gases Ticknor (2012); Bisset et al. (2012); Blakie et al. (2009); Bisset et al. (2011). Probably the most important task would be a verification if these excitations are really elementary and if they are present in the many-body model.

We are grateful to T. Karpiuk, T. Sowiński, P. Deuar and M. Gajda for valuable discussions. The work was supported by the (Polish) National Science Center Grant No. DEC-2012/04/A/ST2/00090.

Vii Appendix: Generating thermal states within classical field approximation

In the Bogoliubov method one approximates the grand canonical Hamiltonian with , where is the Bogoliubov spectrum and the operators and are the annihilation and creation bosonic ladder operators for quasiparticles, respectively. The operator , annihilating a quasiparticle with quasi-momenta is related to the bosonic operators annihilating a particle with a given momenta via relation:


where and . In the classical version, we replace , , and with the complex numbers denoted as , , and , respectively. We introduce the classical probability distribution with the partition function defined such that . Within this approximation the number of quasiparticles with momenta is given by the exponential distribution with the parameter . To generate the initial state foe a given temperature and for a given number of atoms we follow the steps:

  • For we draw the number of quasiparticles using the exponential distribution.

  • For we draw a phase with the uniform distribution on the interval .

  • We construct the amplitude .

  • We compute the amplitudes .

  • We compute the amplitude of the -momentum mode, .

  • We construct the classical field according to Eq. (12).


Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description