Dark Solitons in the Unitary Bose Gas

Calzavara, Martino; Salasnich, Luca

doi:10.3390/sym12060957

Open AccessArticle

Dark Solitons in the Unitary Bose Gas

by

Martino Calzavara

and

Luca Salasnich

^*

Dipartimento di Fisica e Astronomia “Galileo Galilei”, Università di Padova, Via Marzolo 8, 35131 Padova, Italy

^*

Author to whom correspondence should be addressed.

Symmetry 2020, 12(6), 957; https://doi.org/10.3390/sym12060957

Submission received: 15 May 2020 / Revised: 22 May 2020 / Accepted: 2 June 2020 / Published: 4 June 2020

(This article belongs to the Section Physics)

Download

Browse Figures

Versions Notes

Abstract

:

We study the dilute and ultracold unitary Bose gas, characterized by a universal equation of state due to the diverging s-wave scattering length, under a transverse harmonic confinement. From the hydrodynamic equations of superfluids we derive an effective one-dimensional nonpolynomial Schrödinger equation (1D NPSE) for the axial wavefunction which, however, also takes into account the transverse wavefunction. By solving the 1D NPSE we obtain meaningful analytical formulas for the dark (gray and black) solitons of the bosonic system.

Keywords:

Bose-Einstein condensate; unitary Bose gas; dark solitons; nonpolynomial Schrödinger equation

1. Introduction

The possibility of tuning the interaction between atoms near Feshbach resonances [1] almost at will, by precision control of external magnetic fields, is allowing the experimenters to explore the physics of ultracold gases in a wide range of values of the coupling constant characterized by the s-wave scattering length. The possibility to control the scattering length via Feshbach resonances makes ultracold atomic gases an excellent setting for studies of strongly correlated behavior. In particular it made possible to observe first a Fermi gas [2,3,4,5,6] and then also a Bose gas at unitarity [7,8,9,10], that is the limit in which the s-wave scattering length is infinitely large. In this so-called unitary regime, the interactions are as strong as allowed by quantum mechanics, and the physics cannot explicitly depend on the scattering length, leading to the possibility of new types of universal behavior [11]. In the case of a Bose gas, the presence of three body recombination processes [12] complicates the experimental realization of the unitary limit, which can be held for a limited amount of time. Nevertheless, the peculiar properties of a unitary gas, which behaves in a universal manner and exhibits scale invariance, encourage us to study the physics of such systems more in depth.

In this paper we investigate solitonic configurations of the unitary Bose gas. More specifically, we analyze axial and transverse density profiles of dark solitons in bosonic alkali-metal atoms at unitarity under the action of a transverse harmonic confinement. A dark soliton is a self-bound solitary wave, whose existence is made possible by the interplay between the repulsive interaction and the presence of a phase gradient. A dark soliton is characterized by a local decrease in density with respect to a uniform background. The depth of the dark soliton crucially depends on the phase of the complex field which describes the Bose gas.

The realizability of black and gray solitons in ultracold gases is well documented in the weak-coupling limit [13,14] and in a unitary Fermi gas [15]. Moreover, the properties of black solitons have been theoretically investigated for superfluid fermions also in the BCS-BEC crossover [16,17,18,19]. It is then a natural question to ask, whether black and gray solitons can be realized also in the case of Bose gases at unitarity, and what their specific properties are. It is important to stress that dark (black and gray) solitons of the unitary Bose gas have not yet been observed experimentally nor theoretically analyzed. Our theoretical paper can help experimental groups because we are giving analytical and semi-analytical formulas for density profiles, widths, and velocities of these dark solitons. All these quantities are universal because, contrary to the dark solitons of weakly-interacting bosonic gases, they do not depend on the s-wave scattering length.

We first derive the equations which govern the dynamics of these systems. The resulting equation is a nonlinear Schrödinger equation for a complex order parameter

ψ (r, t)

of the bosonic system [20]. Then, we find the solutions which are in the form of dark solitons, specifically as stationary objects in a moving frame.

Mimicking the experimental practice of confining ultracold gases in the minimum of a potential, we study the case in which the Bose gas is placed in a harmonic potential with cylindrical symmetry, which is an approximation to a more realistic cigar shaped potential [13,15]. A Gaussian variational approach, whose effectiveness in similar cases has already been proved [21,22,23,24,25,26], is then deployed to obtain an effective one-dimensional nonpolynomial Schrödinger equation (1D NPSE) for the axial wavefunction. The relevance of using the 1D NPSE in the study of solitons is due to the fact that the obtained solutions are a reliable generalization of familiar strictly one-dimensional results. These generalized solutions take properly into account that the transverse width of the cigar-shaped bosonic cloud depends on the axial coordinate.

In particular, we focus on the case of weak transverse confinement, in which more naive attempts to reduce the dimensionality of the problem fail. We integrate analytically the 1D NPSE thanks to the presence of a constant of motion, leaving behind a solution in terms of an integral, which we finally evaluate by means of numerical techniques. We are then able to find the axial density profile, transverse width and phase of the bosonic system. Moreover, we find the relation between density at the minimum and velocity of the soliton. While doing so, we also find that the theory actually depends on just one free parameter. Finally, we comment on the limits of the weak-coupling approximation.

2. From Euler Equations to Nonlinear Schrödinger Equation

Euler equations describe the dynamics of a non-viscous and irrotational fluid, such as a superfluid Bose gas at zero temperature [27]. In presence of an external potential

U (r)

, these equations are given by

\begin{matrix} \frac{\partial n}{\partial t} + \nabla \cdot (n v) = 0 \end{matrix}

(1)

\begin{matrix} m \frac{\partial v}{\partial t} + \nabla (U (r) + A n^{\frac{2}{3}} + \frac{1}{2} m v^{2}) = 0 . \end{matrix}

(2)

The equation of state of the unitary Bose gas at zero temperature, namely the bulk pressure P as a function of the number density n, is given by

P (n) = \frac{2}{5} A n^{\frac{5}{3}},

(3)

where n is the number density of atoms and

A = u ℏ^{2} / m

, with u a universal and adimensional (positive) coefficient [7,28]. Notice that

P (n) = μ (n) n - E (n)

, where

μ (n) = A n^{2 / 3}

is the bulk chemical potential and

E (n) = (3 / 5) A n^{5 / 3}

is the bulk internal energy. In general, in the absence of external confinement, that is,

U (r) = 0

, if the uniform atomic gas is very dilute and at zero temperature, there are only two length scales: the mean distance between atoms

n^{- 1 / 3}

and the s-wave scattering length

a_{s}

of the inter-atomic interaction. For the unitary gas the scattering length

a_{s}

diverges and consequently the remaining length scale is

n^{- 1 / 3}

. In addition, in the Hamiltonian of the many-body problem always appears the ratio

ℏ^{2} / m

. Thus, from dimensional analysis it follows that the energy scales as

ℏ^{2} n^{2 / 3} / m

and the pressure as

ℏ^{2} n^{5 / 3} / m

. The unknown quantity in our problem with the unitary Bose gas is the universal coefficient u. However, the value of

u = ξ {(6 π^{2})}^{\frac{2}{3}} \frac{5}{6}

has been derived in many ways [28,29,30,31,32] and, depending on the procedure, values in the range

ξ \approx 0.4 \div 1.75

have been reported.

We shall now derive a nonlinear Schrödinger equation (NLSE) with a

4 / 3

nonlinearity exponent [20], starting from Equations (1) and (2) and adding a term, which will account for quantum effects. This is done via the mapping

\begin{matrix} n (r, t) = {| ψ (r, t) |}^{2} \end{matrix}

(4)

\begin{matrix} v (r, t) = \frac{ℏ}{m} \nabla α (r, t), \end{matrix}

(5)

where

α (r, t)

is a scalar field. Equation (5) is indeed fully justified by the fact the fluid is irrotational, that is,

\nabla \land v = 0

.

If we substitute Equations (4) and (5) into Equations (1) and (2) we obtain

\begin{matrix} \frac{\partial | ψ |}{\partial t} + \frac{ℏ}{2 m} (\nabla^{2} α) | ψ | + \frac{ℏ}{m} \nabla α \cdot \nabla | ψ | = 0 \end{matrix}

(6)

\begin{matrix} ℏ \nabla \frac{\partial α}{\partial t} + \nabla (U (r) + A n^{\frac{2}{3}} + \frac{ℏ^{2}}{2 m} {| \nabla α |}^{2}) = 0 . \end{matrix}

(7)

Let us integrate Equation (7) and multiply both sides by

| ψ |

- | ψ | ℏ \frac{\partial α}{\partial t} = | ψ | (U (r) + A n^{\frac{2}{3}} + \frac{ℏ^{2}}{2 m} {| \nabla α |}^{2} + C) .

(8)

The integration constant C appeared, but we shall see that it is not going to affect the dynamics. In order to get a mathematically sound Schrödinger equation, instead of the nasty Guerra-Pusterla Schrödinger equation [33,34,35] given by (8), we add the term

- (ℏ^{2} / 2 m) \nabla^{2} | ψ (r, t) |

on the right hand side of Equation (8). Indeed, adopting the same procedure for the weakly-interacting Bose gas, where the scattering length

a_{s}

is small and

P \propto n^{2}

, one obtains the Gross-Pitaevskii equation, whose validity for weak coupling is well established [36,37]. In other words, this additional term implies a quantum correction, explicitly dependent on

ℏ^{2} / m

in the Euler equations.

We can express the system of real equations as a single complex equation by adding Equation (6) multiplied by

i ℏ

to Equation (8). Then we multiply every term by

e^{i α}

, to get

\begin{matrix} [i ℏ \frac{\partial | ψ |}{\partial t} - ℏ | ψ | \frac{\partial α}{\partial t}] e^{i α} & = [- \frac{ℏ^{2}}{2 m} {(- | ψ | | \nabla α |}^{2} + \nabla^{2} | ψ | + i \nabla^{2} α | ψ | + 2 i \nabla α \cdot \nabla | ψ |)] e^{i α} \\ + [(U (r) + A n^{\frac{2}{3}} + C) | ψ |] e^{i α} . \end{matrix}

(9)

Finally, by defining

ψ (r, t) = | ψ (r, t) | e^{i α (r, t)},

(10)

the last equation can be written as a nonlinear Schrödinger equation (NLSE)

i ℏ \frac{\partial ψ}{\partial t} = - \frac{ℏ^{2}}{2 m} \nabla^{2} ψ + (U (r) + C) ψ + A {| ψ |}^{\frac{4}{3}} ψ .

(11)

C results in a potential energy offset. It can be removed by substituting

ψ \to e^{- \frac{i}{ℏ} C t} ψ

. It is important to stress that Equation (10) implies that the scalar field

α (r, t)

is an angle and the circulation of the velocity field

v (r, t) = (ℏ / m) \nabla α (r, t)

is quantized.

3. Nonpolynomial Schrödinger Equation

Equation (11) can also be derived minimizing the action functional

S = \int d t L

S [ψ, ψ^{★}] = \int_{t_{1}}^{t_{2}} d t \int d^{3} r \frac{i ℏ}{2} (ψ^{★} \frac{\partial ψ}{\partial t} - ψ \frac{\partial ψ^{★}}{\partial t}) - \frac{ℏ^{2}}{2 m} {| \nabla ψ |}^{2} - U (r) {| ψ |}^{2} - A {| ψ |}^{\frac{10}{3}}

(12)

with the constraint

\int d^{3} r {| ψ (r, t) |}^{2} = N .

(13)

Let us suppose that the confining potential is harmonic in the

(x, y)

plane with cylindrical symmetry along the z axis

U (r) = \frac{1}{2} m ω_{⊥}^{2} (x^{2} + y^{2}) .

(14)

Thus, we are assuming that there is no confinement along the z axis but there is instead a harmonic trap of frequency

ω_{⊥}

in the two other directions. We consider the following variational ansatz for

ψ

ψ (r, t) = ϕ (z, t) \frac{e x p (- \frac{x^{2} + y^{2}}{2 σ^{2} (z, t)})}{\sqrt{π σ^{2} (z, t)}},

(15)

where the variational parameters are the axial wavefunction

ϕ (z, t)

and the width

σ (z, t)

of the transverse Gaussian wavefunction. The choice of a Gaussian function in the

(x, y)

plane is justified by the presence of a harmonic potential in the

(x, y)

plane. In fact, when the nonlinear term of the Schrödinger equation is negligible, the ground state wavefunction of the problem in the

(x, y)

plane is indeed a Gaussian, where

σ = a_{⊥}

with

a_{⊥} = \sqrt{ℏ / (m ω_{⊥})}

the characteristic length of the transverse harmonic confinement.

By plugging Equation (15) into the action in Equation (12) we obtain

\begin{matrix} S [ϕ, σ] = \int_{t_{1}}^{t_{2}} d t \int d z \frac{i ℏ}{2} (ϕ^{★} \frac{\partial ϕ}{\partial t} - ϕ \frac{\partial ϕ^{★}}{\partial t}) - \frac{ℏ^{2}}{2 m} | \frac{\partial ϕ}{\partial z} |^{2} \\ + \int_{t_{1}}^{t_{2}} d t \int d z - [\frac{ℏ^{2}}{2 m σ^{2}} + \frac{1}{2} m ω_{⊥}^{2} σ^{2}] {| ϕ |}^{2} - \frac{3}{5} \frac{A}{{(\sqrt[]{π} σ)}^{\frac{4}{3}}} {| ϕ |}^{\frac{10}{3}}, \end{matrix}

(16)

under the approximation of neglecting the space derivative of

σ (z, t)

. In the case of weakly-interacting bosons, which are very well described in three dimensions by the Gross-Pitaevskii equation (i.e., the 3D Schrödinger equation with cubic nonlinearity) this approximation has also been found to be quite good when

σ (z, t)

depends strongly but monotonically on z [21,26].

Minimizing the action with respect to

ϕ

with the constraint (13), we obtain

i ℏ \frac{\partial ϕ}{\partial t} = [- \frac{ℏ^{2}}{2 m} \frac{\partial^{2}}{\partial z^{2}} + \frac{ℏ^{2}}{2 m σ^{2}} + \frac{1}{2} m ω_{⊥}^{2} σ^{2} + \frac{A}{{(\sqrt[]{π} σ)}^{\frac{4}{3}}} {| ϕ |}^{\frac{4}{3}} - μ] ϕ

(17)

and by doing the same with respect to

σ

\frac{σ^{4}}{a_{⊥}^{4}} = 1 + \frac{4}{5} \frac{A}{π^{\frac{2}{3}}} \frac{m}{ℏ^{2}} {| ϕ |}^{\frac{4}{3}} σ^{\frac{2}{3}},

(18)

where

a_{⊥} = \sqrt{ℏ / (m ω_{⊥})}

. We shall refer to Equation (17) endowed with (18) as the nonpolynomial Schrödinger equation [21,22,23,24,25,26].

Let us rewrite the equations in terms of adimensional quantities. We can do that with the substitutions

\begin{matrix} z \to z a_{⊥} \\ t \to \frac{t}{ω_{⊥}} \\ σ (z, t) \to σ (z, t) a_{⊥} \\ ϕ (z, t) \to \sqrt[]{\frac{N}{a_{⊥}}} ϕ (z, t) \\ μ \to ℏ ω_{⊥} (μ + 1) \\ A \to A \frac{ℏ^{2}}{2 m} . \end{matrix}

Equation (17) then becomes

i \frac{\partial ϕ}{\partial t} = [- \frac{1}{2} \frac{\partial^{2}}{\partial z^{2}} + \frac{1}{2 σ^{2}} + \frac{1}{2} σ^{2} - 1 + \frac{A}{2} \frac{N^{\frac{2}{3}}}{π^{\frac{2}{3}}} \frac{{| ϕ |}^{\frac{4}{3}}}{σ^{\frac{4}{3}}} - μ] ϕ,

(19)

while Equation (18) becomes

σ^{4} = 1 + \frac{2}{5} \frac{A}{π^{\frac{2}{3}}} N^{\frac{2}{3}} {| ϕ |}^{\frac{4}{3}} σ^{\frac{2}{3}} .

(20)

Remember that in this equation

σ

is adimensional—it is in units of the characteristic length

a_{⊥}

of the transverse harmonic confinement. In Equation (20) one can identify two regimes: a regime of strong transverse confinement where

σ ≃ 1

and consequently the system is effectively mean-field one-dimensional, but also a regime of weak transverse confinement where

σ ≫ 1

. Thus, the regime of weak transverse confinement is obtained by neglecting the first addend on the right hand side of Equation (20), supposing that

A N^{\frac{2}{3}} {| ϕ |}^{\frac{4}{3}} ≫ 1

. This gives rise to

σ = {(\frac{2}{5} A)}^{\frac{3}{10}} {(\frac{N}{π})}^{\frac{1}{5}} {| ϕ |}^{\frac{2}{5}},

(21)

which plugged back in Equation (19) yields the following equation for

ϕ

i \frac{\partial ϕ}{\partial t} = [- \frac{1}{2} \frac{\partial^{2}}{\partial z^{2}} + \frac{1}{2} {(\frac{5}{2 A})}^{\frac{3}{5}} {(\frac{π}{N})}^{\frac{2}{5}} {| ϕ |}^{- \frac{4}{5}} + \frac{7}{10} {(A)}^{\frac{3}{5}} {(\frac{5}{2} \frac{N}{π})}^{\frac{2}{5}} {| ϕ |}^{\frac{4}{5}} - μ - 1] ϕ .

(22)

In the same approximation we can neglect the term

\propto (A N^{\frac{2}{3}} {| ϕ |}^{\frac{4}{3}})^{- \frac{3}{5}}

, and obtain

i \frac{\partial ϕ}{\partial t} = [- \frac{1}{2} \frac{\partial^{2}}{\partial z^{2}} + γ {| ϕ |}^{\frac{4}{5}} - μ] ϕ,

(23)

where we defined

γ = \frac{7}{10} A^{\frac{3}{5}} {(\frac{5}{2} \frac{N}{π})}^{\frac{2}{5}} .

(24)

Finally, let us also express

σ

in terms of

γ

σ = \sqrt{\frac{4}{7}} \sqrt{γ} {| ϕ |}^{\frac{2}{5}} .

(25)

It will prove useful to study also the limit of strong confinement, that is instead given by

σ = 1

:

i \frac{\partial ϕ}{\partial t} = [- \frac{1}{2} \frac{\partial^{2}}{\partial z^{2}} + γ_{S C} {| ϕ |}^{\frac{4}{3}} - μ] ϕ .

(26)

γ_{S C} = \frac{A}{2} {(\frac{N}{π})}^{\frac{2}{3}} = \frac{5}{4} {(\frac{4}{7} γ)}^{\frac{5}{3}} .

(27)

4. Dark Solitons in the 1D NPSE

A soliton is a solution to Equation (23) of the form

ϕ (z, t) = f (z - v t) e^{i α (z, t)},

(28)

where

f (ζ) \geq 0

and

α (z, t)

are real functions, and v is the velocity of the soliton. Once we substitute (28) in (23) we obtain the equations for f and

α

\begin{matrix} f^{″} & = 2 [f (\frac{\partial α}{\partial t} + \frac{1}{2} {(\frac{\partial α}{\partial z})}^{2} - μ) + γ f^{\frac{9}{5}}], \end{matrix}

(29)

\begin{matrix} f^{'} v & = (f^{'} \frac{\partial α}{\partial z} + \frac{1}{2} f \frac{\partial^{2} α}{\partial z^{2}}) . \end{matrix}

(30)

Upon inspection of Equation (30), we find that the right hand side should depend only on

ζ

, so

\frac{\partial α}{\partial z}

is only function of

ζ

. Let us assume the following

α (z, t) = θ (z - v t) + β (t) .

(31)

If we plug this back into Equations (29) and (30) we get

\begin{matrix} f^{″} & = 2 [f (\frac{d β}{d t} - v θ^{'} + \frac{1}{2} {(θ^{'})}^{2} - μ) + γ f^{\frac{9}{5}}], \end{matrix}

(32)

\begin{matrix} f^{'} v & = (f^{'} θ^{'} + \frac{1}{2} f θ^{″}) . \end{matrix}

(33)

4.1. Phase

Let us solve Equation (33) in order to obtain

θ^{'}

in terms of f. We can rewrite it as

v \frac{{(f^{2})}^{'}}{f} = \frac{{(f^{2} θ^{'})}^{'}}{f},

(34)

for

f (ζ) \neq 0

. After integration with respect to

ζ

we have

f^{2} (ζ) (v - θ^{'} (ζ)) = D,

(35)

where D is an integration constant. If we assume that

\begin{matrix} lim_{ζ \to \pm \infty} | f (ζ) | & = & f_{\infty}, \end{matrix}

(36)

\begin{matrix} lim_{ζ \to \pm \infty} θ^{'} (ζ) & = & 0, \end{matrix}

(37)

and we get

D = v f_{\infty}^{2}

. The expression for

θ^{'}

we are looking for is then

θ^{'} (ζ) = v (1 - \frac{f_{\infty}^{2}}{f^{2} (ζ)}),

(38)

valid for

f (ζ) \neq 0

.

4.2. Modulus

Let us plug Equation (38) back into Equation (32)

f^{″} = 2 \frac{d β}{d t} f - v^{2} f + v^{2} \frac{f_{\infty}^{4}}{f^{3}} - 2 μ f + 2 γ f^{\frac{9}{5}} .

(39)

We can choose

β

in order to get simpler equations. Let us suppose

β (t) = \frac{v^{2}}{2} t,

(40)

so that Equation (39) becomes

f^{″} = v^{2} \frac{f_{\infty}^{4}}{f^{3}} - 2 μ f + 2 γ f^{\frac{9}{5}} .

(41)

Equation (41) can be written in the form of a Newton equation for a conservative force field

\begin{matrix} f^{″} & = & - \frac{d W}{d f}, \end{matrix}

(42)

\begin{matrix} W & = & \frac{1}{2} v^{2} \frac{f_{\infty}^{4}}{f^{2}} + μ f^{2} - \frac{5}{7} γ f^{\frac{14}{5}} . \end{matrix}

(43)

The total energy K is then a conserved quantity for all

ζ

K = \frac{1}{2} f^{' 2} + W,

(44)

and energy balance gives us the following first order equation

K = \frac{1}{2} f^{' 2} + \frac{1}{2} v^{2} \frac{f_{\infty}^{4}}{f^{2}} + μ f^{2} - \frac{5}{7} γ f^{\frac{14}{5}},

(45)

which is satisfied by the two branches

f^{'} = \pm \frac{1}{| f |} \sqrt{2 K f^{2} - v^{2} f_{\infty}^{4} - 2 μ f^{4} + \frac{10}{7} γ f^{\frac{24}{5}}} .

(46)

This equation can be put in the integral form with the separation of variables, giving

| ζ - ζ_{0} | = \int_{f (ζ_{0})}^{f (ζ)} d f \frac{f}{\sqrt{2 K f^{2} - v^{2} f_{\infty}^{4} - 2 μ f^{4} + \frac{10}{7} γ f^{\frac{24}{5}}}},

(47)

where, since the integral is positive, we have chosen the + sign when

ζ - ζ_{0} \geq 0

and the − sign otherwise. Therefore,

f (ζ) = f (| ζ - ζ_{0} |)

has even parity with respect to reflections around

ζ_{0}

, and we just need to study the case

ζ - ζ_{0} \geq 0

.

We now have to find the values of the parameters

μ, K, f_{\infty}

. This can be done by choosing boundary conditions for f.

4.3. Black Solitons

Let us find an odd parity solution (which features a node in the origin for the density n) with a positive horizontal asymptote at

+ \infty

. We set

ζ_{0} = 0

, and since f is even, the phase

θ

will account for the change of sign. Then, the following boundary conditions are needed

f (0) = 0, f (+ \infty) = f_{\infty}, f^{'} (+ \infty) = f^{″} (+ \infty) = 0 .

(48)

This is possible only if

v = 0

(ζ = z)

in order to get rid of the term in Equation (45) which diverges as

f \to 0

. The condition

f^{″} (+ \infty) = 0

applied to Equation (41) gives

μ = γ f_{\infty}^{\frac{4}{5}},

(49)

while

f^{'} (+ \infty) = 0

applied to (45) gives

K = \frac{2}{7} γ f_{\infty}^{\frac{14}{5}} .

(50)

Equation (47) then becomes

z = \int_{0}^{f (z)} \frac{d f}{\sqrt{\frac{4}{7} γ f_{\infty}^{\frac{14}{5}} - 2 γ f_{\infty}^{\frac{4}{5}} f^{2} + \frac{10}{7} γ f^{\frac{14}{5}}}} .

(51)

It is possible to get rid of

f_{\infty}

with the rescaling

h (z) = \frac{f (z)}{f_{\infty}}

(52)

and the parameter

δ = γ f_{\infty}^{4 / 5},

(53)

while

σ

does not change:

σ (z) = \sqrt{\frac{4}{7}} \sqrt{γ} f {(z)}^{\frac{2}{5}} = \sqrt{\frac{4}{7}} \sqrt{δ} h {(z)}^{\frac{2}{5}}

. We obtain

z = \frac{1}{\sqrt{δ}} \int_{0}^{h (z)} \frac{d h}{\sqrt{\frac{4}{7} - 2 h^{2} + \frac{10}{7} h^{\frac{14}{5}}}} .

(54)

Plots of the numerical integration of Equation (54) are shown in Figure 1, for three values of

δ

. Accordingly with those integrations, also

σ (z)

is computed and shown in the plot.

We may define

\tilde{z} = \sqrt{δ} z

and

\tilde{σ} = \frac{σ}{\sqrt{δ}}

, so that the theory has no free parameters. This way there is only one plot to make, but we want to explicitly show in the plots the effects of changing the interaction, so for the moment we are keeping z and

σ

.

4.4. Gray Solitons

Let us find an even parity solution with a positive horizontal asymptote. We set

ζ_{0} = 0

once again, and require the following boundary conditions

f^{'} (0) = 0, f (+ \infty) = f_{\infty}, f^{'} (+ \infty) = f^{″} (+ \infty) = 0 .

(55)

The condition

f^{″} (+ \infty) = 0

applied to Equation (41) gives

μ = \frac{1}{2} v^{2} + γ f_{\infty}^{\frac{4}{5}},

(56)

f^{'} (+ \infty) = 0

applied to Equation (45) gives

K = f_{\infty}^{2} [v^{2} + \frac{2}{7} γ f_{\infty}^{\frac{4}{5}}],

(57)

while by equating K inside Equation (45) evaluated at 0 and at ∞, applying

f^{'} (0) = 0

we get a relation between

v, γ

and

f_{0}

\frac{v^{2}}{γ} = \frac{2}{7} f_{0}^{2} \frac{7 f_{\infty}^{\frac{4}{5}} f_{0}^{2} - 5 f_{0}^{\frac{14}{5}} - 2 f_{\infty}^{\frac{14}{5}}}{2 f_{\infty}^{2} f_{0}^{2} - f_{\infty}^{4} - f_{0}^{4}} .

(58)

Equation (47) then becomes

ζ = \frac{1}{\sqrt[]{γ}} \int_{f_{0}}^{f (ζ)} \frac{d f f}{\sqrt{\frac{10}{7} f^{\frac{24}{5}} - \frac{v^{2}}{γ} f^{4} - 2 f_{\infty}^{\frac{4}{5}} f^{4} + \frac{v^{2}}{γ} 2 f_{\infty}^{2} f^{2} + \frac{4}{7} f_{\infty}^{\frac{14}{5}} f^{2} - \frac{v^{2}}{γ} f_{\infty}^{4}}} .

(59)

Also in this case we can get rid of

f_{\infty}

in the same way as before, obtaining

\begin{matrix} {(\frac{v}{v_{s}})}^{2} = \frac{5}{7} h_{0}^{2} \frac{5 h_{0}^{\frac{14}{5}} - 7 h_{0}^{2} + 2}{{(1 - h_{0}^{2})}^{2}} . \end{matrix}

(60)

\begin{matrix} ζ = \frac{1}{\sqrt{δ}} \sqrt{\frac{7}{2}} \int_{h_{0}}^{h (ζ)} \frac{d h h}{\sqrt{2 h^{2} - 7 h^{4} + 5 h^{\frac{24}{5}} - (2 h_{0}^{2} - 7 h_{0}^{4} + 5 h_{0}^{\frac{24}{5}}) \frac{{(1 - h^{2})}^{2}}{{(1 - h_{0}^{2})}^{2}}}}, \end{matrix}

(61)

where of course

h_{0} = \frac{f_{0}}{f_{\infty}}

and

v_{s} = \sqrt{\frac{2}{5} δ}

is the speed of sound, that is, the propagation velocity of an infinitesimal perturbation in the fluid density.

Now the theory has two free parameters, namely

h_{0}

and

δ

, while v is given in terms of them. In order to be more concise, we now define

\tilde{ζ} = \sqrt{δ} ζ

,

\tilde{σ} = \frac{σ}{\sqrt{δ}}

and

\tilde{v} = \frac{v}{v_{s}},

(62)

so that we are left with a theory that depends only on

h_{0}

(or

\tilde{v}

). The scaled velocity

\tilde{v}

ranges from 0 (black soliton) to 1 (sound wave).

Plots of the numerical integration of Equation (61) are shown in Figure 2, for three values of

\tilde{v}

. Accordingly with those integrations,

\tilde{σ} (\tilde{ζ})

is also computed and shown in the plot.

4.5. Phase

In the case of the black soliton, Equation (38) becomes

θ^{'} = 0

. Since we want the solution to be odd, we choose

θ (z) = - s g n (z) \frac{π}{2} .

(63)

In the general case, we can now compute h and substitute inside Equation (38) to obtain

θ^{'}

. Let us express also

θ

in terms of

\tilde{ζ}

:

θ (ζ) = θ (\frac{\tilde{ζ}}{\sqrt{δ}}) = : \tilde{θ} (\tilde{ζ})

(64)

it follows that

θ^{'} (ζ) = \sqrt{δ} {\tilde{θ}}^{'} (\tilde{ζ}),

(65)

therefore

{\tilde{θ}}^{'} (\tilde{ζ}) = \sqrt{\frac{2}{5}} \tilde{v} (1 - \frac{1}{h^{2} (\tilde{ζ})}) .

(66)

Last equation upon integration yields

\tilde{θ}

. A plot of the numerical integration of Equation (66) is shown in Figure 3 for four values of

\tilde{v}

, including

\tilde{v} = 0

, which is the black soliton case.

5. Weak vs. Strong Transverse Confinement

In the previous sections we have mainly investigated the regime of weak transverse confinement where

σ ≫ 1

. See the discussion of Equation (20). Let us now briefly analyze some properties of the regime of strong transverse confinement where

σ ≃ 1

. In this case the system is truly one dimensional and the equations are much simpler.

By performing computations which are very close to the ones we have already shown in the weak confinement case, one can derive from (26) that in the case of a black soliton the following relation holds

z = \frac{1}{\sqrt{δ_{S C}}} \int_{0}^{h (z)} \frac{d h}{\sqrt{\frac{4}{5} - 2 h^{2} + \frac{6}{5} h^{\frac{10}{3}}}},

(67)

where

δ_{S C} = γ_{S C} f_{\infty}^{\frac{4}{3}}

. Instead, for a gray soliton we obtain

\begin{matrix} {(\frac{v}{v_{s}^{S C}})}^{2} = \frac{3}{5} h_{0}^{2} \frac{3 h_{0}^{\frac{10}{3}} - 5 h_{0}^{2} + 2}{{(1 - h_{0}^{2})}^{2}} . \end{matrix}

(68)

\begin{matrix} ζ = \frac{1}{\sqrt{δ_{S C}}} \sqrt{\frac{5}{2}} \int_{h_{0}}^{h (ζ)} \frac{d h h}{\sqrt{2 h^{2} - 5 h^{4} + 3 h^{\frac{16}{3}} - (2 h_{0}^{2} - 5 h_{0}^{4} + 3 h_{0}^{\frac{16}{3}}) \frac{{(1 - h^{2})}^{2}}{{(1 - h_{0}^{2})}^{2}}}}, \end{matrix}

(69)

where

v_{s}^{S C} = \sqrt{\frac{2}{3} δ_{S C}}

. In dimensional units, taking into account Equation (18), denoting

ρ_{0} = ρ (0)

the minimal axial density of the dark soliton and

ρ_{\infty} = ρ (\pm \infty)

the bulk axial density of the dark soliton, one finds that under the condition

ρ_{\infty} a_{⊥} ≪ 1

the dark soliton is surely in the strictly 1D regime of strong transverse confinement at any point of the axial coordinate. Instead, under the condition

ρ_{0} a_{⊥} ≫ 1

the dark soliton is surely in the 3D regime of weak transverse confinement at any point of the axial coordinate. Clearly, in the case of the black soliton, where

ρ_{0} = 0

, near the minimum the bosonic cloud cannot be in the regime of weak transverse confinement.

It is important to stress that the comparison of our theoretical results with future experiments is constrained not only by the short lifetime of unitary bosonic gas due to fast three-body recombinations but also by the lifetime due to the snake instability [16,17,18,19]. A dark soliton, that is not strictly 1D, has a snaking transverse oscillation which breaks the axial symmetry and eventually destroys the soliton. Since the NPSE preserves axial symmetry, it is not suitable for investigating snaking oscillations and the dynamics of the dissolving dark soliton. However, it is possible to write a triaxial NPSE by using two transverse width

σ_{x} (z, t)

and

σ_{y} (z, t)

in the ansatz of Equation (15); see Reference [38] for details in the case of weakly-interacting bosons.

6. Conclusions

We have derived the dynamical equations that govern the motion of a unitary Bose gas at zero temperature. Using these equations, we have studied dark solitons in a unitary Bose gas, confined in a cylindrically symmetric potential. We have reduced the dimensionality of the 3D nonlinear Schrödinger problem by means of the 1D nonpolynomial Schrödinger equation approach, which also takes dynamically into account the transverse width. We have solved the equations employing for the most analytical techniques, in order to find axial density, transverse width, phase and velocity of both black and gray solitons, in the weak and strong confinement regimes. We have found that the weak confinement approximation breaks down in the case of a black soliton.

Our theoretical predictions could be tested experimentally employing Bose gases of alkali-metal atoms at ultra-low temperatures, whose scattering length can be tuned by means of Feshbach resonance techniques to reach the unitary limit. Despite the difficulties which presently limit the time during which such systems can be observed due to three-body losses but also to the snake instability, we are confident that the recent developments in the experimental techniques are going to allow in the future to put our predictions at test.

Author Contributions

Conceptualization, L.S.; Investigation, L.S. and M.C.; analytical calculations, M.C.; software M.C.; writing–original draft preparation, M.C.; writing–review and editing, L.S. and M.C.; supervision, L.S.; funding acquisition, L.S. All authors have read and agreed to the published version of the manuscript.

Funding

Partial support from the BIRD project “Time-dependent density functional theory of quantum atomic mixtures” of the University of Padova.

Acknowledgments

L.S. thanks Alberto Cappellaro and Andrea Tononi for useful discussions.

Conflicts of Interest

The authors declare no conflict of interest.

References

Inouye, S.; Andrews, M.R.; Stenger, J.; Miesner, H.-J.; Stamper-Kurn, D.M.; Ketterle, W. Observation of Feshbach resonances in a Bose-Einstein condensate. Nature 1998, 392, 151–154. [Google Scholar] [CrossRef]
O’Hara, K.M.; Hemmer, S.L.; Gehm, M.E.; Granade, S.R.; Thomas, J.E. Observation of a Strongly Interacting Degenerate Fermi Gas of Atoms. Science 2002, 298, 2179–2182. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Regal, C.A.; Ticknor, C.; Bohn, J.L.; Jin, D.S. Creation of ultracold molecules from a Fermi gas of atoms. Nature 2003, 424, 47–50. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Strecker, K.E.; Partridge, G.B.; Hulet, R.G. Conversion of an Atomic Fermi Gas to a Long-Lived Molecular Bose Gas. Phys. Rev. Lett. 2003, 91, 080406. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Jochim, S.; Bartenstein, M.; Altmeyer, A.; Hendl, G.; Riedl, S.; Chin, C.; Hecker Denschlag, J.; Grimm, R. Bose-Einstein Condensation of Molecules. Science 2003, 302, 2101–2103. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Gehm, M.E.; Hemmer, S.L.; Granade, R.S.; O’Hara, K.M.; Thomas, J.E. Mechanical stability of a strongly interacting Fermi gas of atoms. Phys. Rev. A 2003, 68, 011401. [Google Scholar] [CrossRef] [Green Version]
Li, W.; Ho, T.-L. Bose Gases near Unitarity. Phys. Rev. Lett. 2012, 108, 195301. [Google Scholar] [CrossRef] [Green Version]
Rem, B.S.; Grier, A.T.; Ferrier-Barbut, I.; Eismann, U.; Langen, T.; Navon, N.; Khaykovich, L.; Werner, F.; Petrov, D.S.; Chevy, F.; et al. Lifetime of the Bose Gas with Resonant Interactions. Phys. Rev. Lett. 2013, 110, 163202. [Google Scholar] [CrossRef]
Makotyn, P.; Klauss, C.E.; Goldberger, D.L.; Cornell, E.A.; Jin, D.S. Universal dynamics of a degenerate unitary Bose gas. Nat. Phys. 2014, 10, 116–119. [Google Scholar] [CrossRef] [Green Version]
Fletcher, R.J.; Gaunt, A.L.; Navon, N.; Smith, R.P.; Hadzibabic, Z. Stability of a Unitary Bose Gas. Phys. Rev. Lett. 2013, 111, 125303. [Google Scholar] [CrossRef] [Green Version]
Giorgini, S.; Pitaevskii, L.P.; Stringari, S. Theory of ultracold atomic Fermi gases. Rev. Mod. Phys. 2008, 80, 1215. [Google Scholar] [CrossRef] [Green Version]
Roberts, J.L.; Claussen, N.R.; Cornish, S.L.; Wieman, C.E. Magnetic Field Dependence of Ultracold Inelastic Collisions near a Feshbach Resonance. Phys. Rev. Lett. 2000, 85, 728. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Burger, S.; Bongs, K.; Dettmer, S.; Ertmer, W.; Sengstock, K.; Sanpera, A.; Shlyapnikov, G.V.; Lewenstein, M. Dark Solitons in Bose-Einstein Condensates. Phys. Rev. Lett. 1999, 83, 5198. [Google Scholar] [CrossRef] [Green Version]
Denschlag, J.; Simsarian, J.E.; Feder, D.L.; Clark, C.W.; Collins, L.A.; Cubizolles, J.; Deng, L.; Hagley, E.W.; Helmerson, K.; Reinhardt, W.P.; et al. Generating Solitons by Phase Engineering of a Bose-Einstein Condensate. Science 2000, 287, 97–101. [Google Scholar] [CrossRef] [PubMed]
Ku, M.J.H.; Mukherjee, B.; Yefsah, T.; Zwierlein, M.W. Cascade of Solitonic Excitations in a Superfluid Fermi gas: From Planar Solitons to Vortex Rings and Lines. Phys. Rev. Lett. 2016, 116, 045304. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Antezza, M.; Dalfovo, F.; Pitaevskii, L.P.; Stringari, S. Dark solitons in a superfluid Fermi gas. Phys. Rev. A 2007, 76, 043610. [Google Scholar] [CrossRef] [Green Version]
Cetoli, A.; Brand, J.; Scott, R.G.; Dalfovo, F.; Pitaevskii, L.P. Snake instability of dark solitons in fermionic superfluids. Phys. Rev. A 2013, 88, 043639. [Google Scholar] [CrossRef] [Green Version]
Lombardi, G.; Van Alphen, W.; Klimin, S.N.; Tempere, J. Soliton core filling in superfluid Fermi gases with spin-imbalance. Phys. Rev. A 2017, 96, 033609. [Google Scholar] [CrossRef] [Green Version]
Van Alphen, W.; Lombardi, G.; Klimin, S.N.; Tempere, J. Dark soliton collisions in superfluid Fermi gases. New J. Phys. 2018, 20, 053052. [Google Scholar] [CrossRef] [Green Version]
Adhikari, S.K.; Salasnich, L. Nonlinear Schrödinger equation for a superfluid Bose gas from weak coupling to unitarity: Study of vortices. Phys. Rev. A 2008, 77, 033618. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L.; Parola, A.; Reatto, L. Effective wave equations for the dynamics of cigar-shaped and disk-shaped Bose condensates. Phys. Rev. A 2002, 65, 043614. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L.; Parola, A.; Reatto, L. Transition from three dimensions to one dimension in Bose gases at zero temperature. Phys. Rev. A 2004, 70, 013606. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L.; Parola, A.; Reatto, L. Quasi-one-dimensional bosons in three-dimensional traps: From strong-coupling to weak-coupling regime. Phys. Rev. A 2005, 72, 025602. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L.; Malomed, B.A. Vector solitons in nearly one-dimensional Bose-Einstein condensates. Phys. Rev. A 2006, 74, 053610. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L.; Cetoli, A.; Malomed, B.A.; Toigo, F.; Reatto, L. Bose-Einstein condensates under a spatially modulated transverse confinement. Phys. Rev. A 2007, 76, 013623. [Google Scholar] [CrossRef]
Salasnich, L.; Malomed, B.A.; Toigo, F. Matter-wave vortices in cigar-shaped and toroidal waveguides. Phys. Rev. A 2007, 76, 063614. [Google Scholar] [CrossRef] [Green Version]
Pitaevskii, L.P.; Stringari, S. Bose-Einstein Condensation; Oxford Univ. Press: Oxford, UK, 2003. [Google Scholar]
Rossi, M.; Salasnich, L.; Ancilotto, F.; Toigo, F. Monte Carlo simulations of the unitary Bose gas. Phys. Rev. A 2014, 89, 041602. [Google Scholar] [CrossRef] [Green Version]
Cowell, S.; Heiselberg, H.; Mazets, I.; Morales, J.; Pandharipande, V.R.; Pethick, C.J. Cold Bose Gases with Large Scattering Lengths. Phys. Rev. Lett. 2002, 88, 210403. [Google Scholar] [CrossRef] [Green Version]
Song, J.L.; Zhou, F. Ground State Properties of Cold Bosonic Atoms at Large Scattering Lengths. Phys. Rev. Lett. 2009, 103, 025302. [Google Scholar] [CrossRef] [Green Version]
Lee, Y.-L.; Lee, Y.-W. Universality and stability for a dilute Bose gas with a Feshbach resonance. Phys. Rev. A 2010, 81, 063613. [Google Scholar] [CrossRef] [Green Version]
Diederix, J.M.; van Heijst, T.C.F.; Stoof, H.T.C. Universality and stability for a dilute Bose gas with a Feshbach resonance. Phys. Rev. A 2011, 84, 033618. [Google Scholar] [CrossRef] [Green Version]
Guerra, F.; Pusterla, M. A nonlinear schrödinger equation and its relativistic generalization from basic principles. Lett. Nuovo Cim. 1982, 34, 351–356. [Google Scholar] [CrossRef]
Salasnich, L.; Toigo, F. Extended Thomas-Fermi density functional for the unitary Fermi gas. Phys. Rev. A 2008, 78, 053626. [Google Scholar] [CrossRef] [Green Version]
Salasnich, L. Hydrodynamics of Bose and Fermi superfluids at zero temperature: The superfluid nonlinear Schrödinger equation. Laser Phys. 2009, 19, 642–646. [Google Scholar] [CrossRef]
Gross, E.P. Structure of a quantized vortex in boson systems. Nuovo Cim. 1961, 20, 454–477. [Google Scholar] [CrossRef] [Green Version]
Pitaevskii, L.P. Vortex lines in an imperfect Bose gas. Sov. Phys. JETP 1961, 13, 451–454. [Google Scholar]
Salasnich, L. Generalized nonpolynomial Schrödinger equations for matter waves under anisotropic transverse confinement. J. Phys. A: Math. Theor. 2009, 42, 335205. [Google Scholar] [CrossRef] [Green Version]

Figure 1. Black soliton. Normalized axial density profile

ρ (z) = h^{2} (z) = {(f (z) / f_{\infty})}^{2}

and transverse width profile

σ (z)

vs. axial coordinate z, for three values of the parameter

δ = γ | f_{\infty} |^{4 / 5}

, with

γ

given by Equation (24) and

f_{\infty}

the bulk value of the soliton wavefunction.

Figure 1. Black soliton. Normalized axial density profile

ρ (z) = h^{2} (z) = {(f (z) / f_{\infty})}^{2}

and transverse width profile

σ (z)

vs. axial coordinate z, for three values of the parameter

δ = γ | f_{\infty} |^{4 / 5}

, with

γ

given by Equation (24) and

f_{\infty}

the bulk value of the soliton wavefunction.

Figure 2. Gray soliton. Upper panel: Scaled axial density profile

ρ (\tilde{ζ}) = h {(\tilde{ζ})}^{2} = {(f (z) / f_{\infty})}^{2}

vs. scaled comoving axial coordinate

\tilde{ζ} = (z - v t) \sqrt{δ}

. Lower panel: Scaled transverse width

\tilde{σ} (\tilde{ζ})

vs.

\tilde{ζ}

. Notice that

δ = γ | f_{\infty} |^{4 / 5}

, with

γ

given by Equation (24) and

f_{\infty}

the bulk value of the soliton wavefunction.

\tilde{v} = v / v_{s}

is the velocity rescaled by the sound velocity

v_{s}

.

Figure 2. Gray soliton. Upper panel: Scaled axial density profile

ρ (\tilde{ζ}) = h {(\tilde{ζ})}^{2} = {(f (z) / f_{\infty})}^{2}

vs. scaled comoving axial coordinate

\tilde{ζ} = (z - v t) \sqrt{δ}

. Lower panel: Scaled transverse width

\tilde{σ} (\tilde{ζ})

vs.

\tilde{ζ}

. Notice that

δ = γ | f_{\infty} |^{4 / 5}

, with

γ

given by Equation (24) and

f_{\infty}

the bulk value of the soliton wavefunction.

\tilde{v} = v / v_{s}

is the velocity rescaled by the sound velocity

v_{s}

.

Figure 3. Scaled phase

\tilde{θ} (\tilde{ζ})

of gray (

\tilde{v} > 0

) and black (

\tilde{v} = 0

) solitons vs. the scaled comoving axial coordinate

\tilde{ζ}

.

Figure 3. Scaled phase

\tilde{θ} (\tilde{ζ})

of gray (

\tilde{v} > 0

) and black (

\tilde{v} = 0

) solitons vs. the scaled comoving axial coordinate

\tilde{ζ}

.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Calzavara, M.; Salasnich, L. Dark Solitons in the Unitary Bose Gas. Symmetry 2020, 12, 957. https://doi.org/10.3390/sym12060957

AMA Style

Calzavara M, Salasnich L. Dark Solitons in the Unitary Bose Gas. Symmetry. 2020; 12(6):957. https://doi.org/10.3390/sym12060957

Chicago/Turabian Style

Calzavara, Martino, and Luca Salasnich. 2020. "Dark Solitons in the Unitary Bose Gas" Symmetry 12, no. 6: 957. https://doi.org/10.3390/sym12060957

APA Style

Calzavara, M., & Salasnich, L. (2020). Dark Solitons in the Unitary Bose Gas. Symmetry, 12(6), 957. https://doi.org/10.3390/sym12060957

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

Dark Solitons in the Unitary Bose Gas

Abstract

1. Introduction

2. From Euler Equations to Nonlinear Schrödinger Equation

3. Nonpolynomial Schrödinger Equation

4. Dark Solitons in the 1D NPSE

4.1. Phase

4.2. Modulus

4.3. Black Solitons

4.4. Gray Solitons

4.5. Phase

5. Weak vs. Strong Transverse Confinement

6. Conclusions

Author Contributions

Funding

Acknowledgments

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI