University of Amsterdam
Institute for Theoretical Physics
Master’s project
Classical Solitons in the
Quantum Nonlinear Schr¨
odinger
Equation
Author:
Yuri van Nieuwkerk
Student ID:
6179703
Supervisor:
Prof. Dr. Jean-S´
ebastien Caux
Examiner:
Dr. Vladimir Gritsev
September 1, 2015 - August 6, 2016
60 ECTS
Abstract
In this thesis, classical solitons are constructed in the context of the quantum Nonlinear Schr¨odinger equation (NLS). Solitons arise as stable, localized solutions to the classical NLS equation, which can be seen as a mean field version of its fully quantum mechanical counterpart. This thesis aims to identify exact quantum states underlying such mean field solitons.
For the attractive r´egime of the quantum model, Wadachi et al. [1] have identified such quantum objects: when bound states from the Bethe Ansatz framework, known as strings, are superposed to form a wave packet, the density profile of a classical soliton is recovered. Under time evolution, however, such profiles are subject to decay. This thesis introduces a wave packet for which the decay process is described as an asymptotic series. For late times, the resulting density profile behaves as a flat Gaus-sian which broadens with time. For early times, the classical soliton is stable. The characteristic lifetime separating these r´egimes diverges in the large-N limit for weak interactions. This behavior can be reconciled with results for the mean field soliton by treating its center of mass as a quantum observable.
For the repulsive scenario, Sato et al. [2] have superposed holes in the ground state Fermi sea to reproduce the soliton’s density profile. Although this approach fails to produce infinite lifetimes in general, Sato et al. have presented numerical evidence of divergent lifetimes in the weakly interacting limit. This thesis supports these claims with analytical means. In addition to this, a number of new candidates for localization are proposed. A density notch, constructed by Girardeau et al [3] at infinite repulsion, is reproduced in an interaction-independent description. The corresponding lifetimes scale as the square of the interaction strength. A second type of localized state is produced using split Fermi seas at weak coupling, leading to a lifetime proportional to the number of particles, which signifies stability in the thermodynamic limit. The shape of these distributions is studied with bosonization techniques, leading to results which correspond strikingly well with numerical data.
Contents
1 Introduction 7
2 The nonlinear Schr¨odinger equation in cold atomic gases 9
2.1 Bose-Einstein condensation . . . 9
2.2 The Gross-Pitaevskii equation . . . 12
2.3 Mean field solitons . . . 16
2.4 Summary and look ahead . . . 20
3 The quantum nonlinear Schr¨odinger equation 23 3.1 Exact eigenstates . . . 23
3.2 Exact spectrum . . . 28
3.3 Density form factors . . . 32
3.4 The thermodynamic limit . . . 34
4 Bright solitons 40 4.1 N -strings . . . 40
4.2 Result: time-evolving bright solitons . . . 42
4.3 Summary of results . . . 45
5 Dark solitons 46 5.1 Single hole states . . . 46
5.2 Result: degenerate states at strong interaction . . . 50
5.3 Result: split Fermi seas at weak interaction . . . 58
6 Conclusions 65
A Decay time of a mean field soliton 67
B Split Fermi seas: extension of Gaudin’s expansion for weak coupling 69 C Split Fermi seas: bosonization identities and finite-size corrections 72
Chapter 1
Introduction
“I was observing the motion of a boat which was rapidly drawn along a narrow channel by a pair of horses, when the boat suddenly stopped - not so the mass of water in the channel which it had put in motion; it accumulated round the prow of the vessel in a state of violent agitation, then suddenly leaving it behind rolled forward with great velocity, assuming the form of a large solitary elevation, a rounded, smooth and well defined heap of water, which continued its course along the channel apparently without change of form or diminution of speed. I followed it on horseback, and overtook it still rolling on at a rate of some eight or nine miles an hour, preserving its original figure some thirty feet long and a foot to a foot and a half in height. Its height gradually diminished, and after a chase of one or two miles I lost it in the windings of the channel. Such, in the month of August 1834, was my first chance interview with that singular and beautiful phenomenon which I have called the Wave of Translation, . . . ”
-John Scott Russel, Report on Waves. York, September 1844. This thesis investigates the occurrence of solitons in one-dimensional bosonic gases. Solitons, or Waves of Translation, as John Scott Russel called them, are a classical phenomenon. They occur in systems where Planck’s constant can be neglected in one way or another, be it because temperatures are high, because particle numbers are large, or simply because a measurer like John Scott Russel has to exhaust his horse in order to cover them. One-dimensional bose gases, on the other hand, are an intrinsically quantum mechanical environment. Their exact wave functions, along with important observables, can be obtained in the Bethe Ansatz framework. Nevertheless, the behavior of the one-dimensional bose gas at low energies is approximated by a classical equation of motion called the nonlinear Schr¨odinger equation, whenever the number of particles is large enough. Strikingly, this equation predicts the existence of solitons, offering a fine example of a classical object emerging from a quantum mechanical context. This raises the main question in this project: which exact quantum mechanical wave functions start to behave as classical solitons whenever the nonlinear Schr¨odinger equation becomes valid?
Ultracold bosonic gases have been subjected to intense research, ever since Albert Einstein’s prediction [4] of a phase transition to a collective state of matter in 1924. This phase transition can only be described in a statistical framework which is fundamentally quantum mechanical, as proposed by Satyendra Nath Bose [5]. The resulting process of Bose-Einstein condensation, as it came to be called, has turned out to contain a wealth of interesting physics. It was applied to superfluidity by Fritz London in 1938, after which Lev Landau and Nikolay Bogolyubov devel-oped ways to describe its elementary excitations. The concept of Bose-Einstein condensation was mathematically formalized by Oliver Penrose and Lars Onsager in the 1950’s. In the same decade, Eugene Gross and Lev Pitaevskii independently developed a mean field description of interactions in a Bose-Einstein condensate [6, 7], which reduces to the nonlinear Schr¨odinger equation for free systems in one dimension.
Although Bose-Einstein condensation in a strict sense does not occur in one-dimensional sys-tems, the nonlinear Schr¨odinger equation, or one-dimensional Gross-Pitaevskii equation in the absence of a potential, has turned out to display interesting features. In 1971, Toshio Tsuzuki showed the existence of single soliton solutions to the equation [8]. This result was extended
greatly by Vladimir Zakharov and Aleksey Shabat in 1972, who subjected the problem to the mathematical rigor of the inverse scattering method [9]. This method allows for a systematic description of solitons of higher order.
In the mean time, an exact theory for the one-dimensional bose gas had been developed, based on the seminal Ansatz [10] by Hans Bethe in 1931. In 1963, Lieb and Liniger used this Ansatz to find exact eigenstates for the system, which has since been referred to as the Lieb-Liniger model [11,12]. The corresponding thermodynamics were described in 1969 by Yang and Yang [13]. A new direction was laid bare at the end of the 1970’s, when Fadeev and Sklyanin developed the quantum inverse scattering method [14]. This method offers an algebraic framework in which Bethe Ansatz results can be reproduced, with the advantage that it enables the computation of correlation functions [15–17].
The algebraic similarities between the inverse scattering method and its quantum mechanical counterpart soon raised the following question: is there a similar correspondence between classical solitons in the mean field model and an exact quantum mechanical object? Fadeev and cowork-ers established such a connnection in 1976 [18]. Based on a resemblance between the objects’ dispersion relations, solitons in the repulsive mean field model were associated with holes in the Lieb-Liniger ground state. Wadachi and colleagues [1] established a similar connection between solitons in the attractive mean field model and objects called strings, which characterize states in the attractive Lieb-Liniger model.
A problem with the above correspondences, however, is that the quantum mechanical objects under consideration are translationally invariant eigenstates of the total momentum operator, whereas solitons are highly localized. To reproduce mean field solitons using exact quantum states, one is thus forced to superpose such translationally invariant momentum modes to form wave packets. This has recently been attempted for the repuslive case by Sato et al. [2, 19, 20]. Although such attempts can reproduce the correct density profiles at t = 0, these start to loose their shape as soon as they are time-evolved. Similar results were observed by Lai and Haus [21] when constructing wave packets in the attractive scenario.
In the 1990’s, the advent of experimental cooling techniques further increased the interest in cold bosonic gases. Walraven and Silvera have produced the first instance of a dilute ultracold gas using evaporative cooling [22] as early as 1979. Combining this technique with laser cooling, a BEC was created for the first time by Cornell and his coworkers at JILA in 1995 [23]. With very similar techniques, quasi one-dimensional bose gases have since been created in laboratories around the world [24–27], leading to a great revival of theoretical research, both into mean field descriptions and exact techniques. A similar revival has occurred in the investigation of solitons, thanks to their observation in both the repulsive [28] and attractive [29] scenarios.
These experimental results increase the urgency behind this project’s main question: how can classical solitons in the one-dimensional bose gas, known from mean field theory and observed in experiment, be described in an exact way? To answer this question, this thesis is organized as follows: in Chapter 2, cold bosonic gases are introduced. The evidence for the existence of solitons is summarized, using the mean field description by the nonlinear Schr¨odinger equation, and two examples of experimental realizations. Via an argument due to Castin [30], a prediction is given for soliton lifetimes, thus leading to a set of requirements which quantum objects need to fulfill, if they are to describe a soliton. In Chapter, 3 an exact description of the quantum model is given by the Bethe Ansatz method. This leads to a basis of exact eigenfunctions, a framework for thermodynamics and a description of the relevant correlation functions. Chapter 4 offers the first results: for the attractive scenario, a superposition of string solutions is found which reproduces the single soliton solution. Importantly, the proposed superposition also reproduces the correct timescales of decay, as predicted by Castin [30]. Chapter 5. attempts to do the same for solitons in the repulsive r´egime. Decay times observed by Sato et al. [20] are reproduced with analytical means. In addition to this, two new localized states are proposed. The first one coincides with a result by Girardeau et al. [3] in the infite coupling limit, but is derived in a form which generalizes it to arbitrary coupling strength. Lifetimes of this object are determined numerically. The second object seems to be undescribed in the literature: it is a density notch built up from split Fermi seas [31] at weak interaction. Bosonization techniques are used to obtain the notch depth as a function of system size, in an approximation which coincides strikingly well with numerical results. In Chapter 6, the main question is revisited and conclusions are presented.
Chapter 2
The nonlinear Schr¨
odinger
equation in cold atomic gases
2.1
Bose-Einstein condensation
Albert Einstein’s last great discovery was inspired by a young physicist named Satyendra Nath Bose. In his 1924 article [5], Bose showed how Planck’s formula for the expected number of photons nk occupying mode k with energy 0(k),
hnki = 1
e0(k)/T− 1, (2.1)
can be rederived from a purely quantum mechanical starting point. Einstein was inspired by Bose’s reasoning, and generalized it to cover all possible identical particles with integer spin, a class which Paul Dirac would later designate as bosons. In a recently rediscovered manuscript [4], Einstein made the prophetic prediction that Equation 2.1 gives rise to a phase transition at a critical temperature Tc. Below this temperature, Einstein expected a gas of free bosons to occupy a single quantum state with a macroscopic number of particles. These particles form a collective state of matter which cannot be described by classical physics, and which would later be referred to as a Bose-Einstein Condensate.
2.1.1
Bose-Einstein Condensation in an ideal three-dimensional bose
gas
The most straightforward example of Bose-Einstein Condensation, and the one described by Ein-stein in his original article [4], occurs in a gas of N free massive bosons placed in a three-dimensional box of volume V . Switching from the case of photons, where the chemical potential is zero, to massive bosons with negative chemical potential µ, the total number of particles is obtained by summing the distribution 2.1 over all modes k, using energies 0(k) = k2− µ:
N V = 1 V X k 1 e(k2−µ)/T − 1. (2.2)
Obviously, Equation 2.2 must be satisfied for all temperatures. This means that for any decrease in temperature, µ needs to increase, getting closer to zero. Is there a finite temperature Tc at which the chemical potential µ vanishes? In three dimensions the answer is yes. As long as Equation 2.1 defines a smooth distribution, the sum can be turned into an integral, after which solid angles are integrated and a geometric series can be recognized in the denominator. This series then sums to a Riemann zeta function ζ, to give
N V (hnki smooth) = 1 (2π)3 ˆ V d3k e(k2−µ)/T − 1 = e 2µ T T 4π 3/2 ζ 3 2 . (2.3)
2.1. BOSE-EINSTEIN CONDENSATION
From this expression it becomes clear that the chemical potential vanishes below a critical tem-perature given by Tc= 4π N V 2/3 ζ 3 2 −2/3 . (2.4)
Below Tca paradox presents itself: since the chemical potential is bounded from above by zero, it can not increase any further, as T drops below Tc. This means that there are no solutions to Equation 2.3 whenever T < Tc, since the right hand side of this equation will always be smaller thanN
V in that case. Nevertheless, Equation 2.2 must still be fulfilled. The only difference between Equations 2.2 and 2.3 is the fact that a sum has been replaced by an integral, under the assumption that the distribution 2.1 is sufficiently smooth everywhere. It is this assumption of smoothness which apparently breaks down for T < Tc. Since the thermal cloud can no longer accommodate N particles for these temperatures, a certain number of particles N0will occupy the lowest energy mode at k = 0, leading to a singularity in the distribution 2.1. The difference between Equation 2.3 and Equation 2.2 is then accounted for by these N0 particles occupying the zero momentum mode, forming a singular peak in the distribution which is not picked up when moving from a sum to an integral in Equation 2.3. This anomalously large collection of particles in the ground state is referred to as the condensate.
The number of particles in the zero momentum mode can be determined using Equation 2.2. The right hand side describes the number of particles outside of the condensate, N − N0, which is equal to N for T = Tc. Writing N = 4πTc
3/2
ζ 32 and N0= N − 4πT 3/2
ζ 32, one then finds N0 N = 1 − T Tc 3/2 (2.5) for T ≤ Tc. Clearly, as T decreases beyond Tc, a macroscopic number of particles starts to occupy the same quantum state of vanishing momentum. A Bose-Einstein Condensate (BEC) has formed. “The theory is pretty, but is there truth in it?” Einstein posed this question in a letter [32] to Paul Ehrenfest, whilst discussing his proposed condensation process. The final and affirmative answer would take 70 years to arrive: it was in 1995 that Eric Cornell and his coworkers at JILA reported the first realization of a BEC in87Rb atoms [23]. Wolfgang Ketterle and his MIT group would soon reach a similar result using sodium [33].
2.1.2
Broken symmetry and low dimensions
The experiments at MIT established an important feature of BECs: the coherent nature of the condensate. Since wave functions in quantum mechanics are defined up to an overall phase factor, they are said to possess a U (1) symmetry. But as particles start to condense into a single quantum state, this symmetry is spontaneously broken: all single particle wave functions obtain the same, well-defined phase. Ketterle et al. were able to implicitly detect this phase using the interference of two freely expanding BECs [34].
As is known from Noether’s theorem, the existence of a symmetry implies the conservation of an associated quantity. In the case of U (1) symmetry, it turns out that the total number of particles is conserved. Apparently, this conservation of particles is no longer true for a BEC, since the associated symmetry is broken. In fact, a BEC containing N 1 particles is essentially equivalent to a BEC containing N −1 particles, that is, the two states have an overlap approaching 1 as N → ∞.
The coherent nature of a BEC also manifests itself in its spatial correlations. Following the treatment by Pitaevskii and Stringari, [35], the enhancement of correlations becomes clear from the expectation value of the one-body density operator
n(1)(x, x0) =D ˆΨ†(x) ˆΨ(x0)E, (2.6) where the quantum field ˆΨ(x) ˆΨ†(x)annihilates (creates) a particle at position x, whilst obeying bosonic commutation relations:
2.1. BOSE-EINSTEIN CONDENSATION
Figure 2.1: Page two of Einstein’s 1924 manuscript “Quantentheorie des einatomigen idealen Gases” [4]. The second paragraph closes with Einstein’s seminal prediction of a phase transition, inspired by Satyendra Nath Bose: “Es tritt eine Scheidung ein; ein Teil ‘kondensiert’, der Rest bleibt ein ‘ges¨attiges ideales Gas’.”
2.2. THE GROSS-PITAEVSKII EQUATION
For a uniform and isotropic system, the expectation value 2.6 only depends on the relative coor-dinate s = x − x0 and it can be written in terms of Fourier modes as
n(1)(s) = 1 V
ˆ
ddk n(1)(k) e−ik·s, (2.8)
where n(1)(k) is the bosonic distribution over modes k. Now as long as this distribution is smooth, Equation 2.8 vanishes as s → ∞. However, we have seen that the existence of a BEC means that n(1)(k) develops a singularity near k = 0 and behaves like
n(1)(k) ≈ N0δ(k) + ˜n(1)(k), (2.9)
with ˜n(1)(k) regular at k = 0. In this case, the spatial correlator 2.8 obtains a finite value even at an infinite separation,
n(1)(s → ∞) = N0
V . (2.10)
These long-range correlations in the off-diagonal elements (that is, x 6= x0) of the one-body density operator are referred to as off-diagonal long-range order.
Both broken U (1) symmetry and off-diagonal long-range order are essential features of BECs. They are manifestations of a continuous symmetry that is broken spontaneously. This fact has far-reaching consequences for the occurrence of BECs in dimensions one and two. It has been proven by Mermin and Wagner [36], Hohenberg [37] and Coleman [38] that continuous symmetry can not be spontaneously broken at a finite temperature in infinite systems of dimension smaller than three. As a result of this fact, Bose-Einstein condensation in two dimensions can only occur at T = 0. In the one-dimensional case, which is the subject of this thesis, the situation is even more constrained and a true BEC will never form. Nevertheless, the theory of BECs offers a number of interesting predictions about one-dimensional physics, which will be treated in Section 2.3. To enable such predictions, an effective description of the dynamics at T = 0 will be given in terms of a classical differential equation, known as the Gross-Pitaevskii equation.
2.2
The Gross-Pitaevskii equation
In Einstein’s initial treatment of BECs in three dimensions, interactions between particles were neglected altogether. This is not as restrictive as it may seem: in spite of significant quantita-tive differences, the qualitaquantita-tive aspects of BECs described in Section 2.1, such as broken U (1) symmetry and off-diagonal long-range order, are still observed for sufficiently weak interactions. The experiments at JILA and MIT [23], [33] have confirmed this fact from the outset. In one-dimensional systems, however, interactions play a far more crucial role, as will become clear in Section 2.3.
2.2.1
Pseudopotential for weak interactions
To describe the effect of interactions, consider the following Hamiltonian: ˆ H = ˆ ddxh− ˆΨ†(x) ∆xΨ(x) + ˆˆ Ψ†(x)V (x) ˆΨ(x) + c ˆΨ†(x) ˆΨ†(x) ˆΨ(x) ˆΨ(x) i , (2.11)
with Heisenberg equation of motion i∂ ˆΨ(x, t)
∂t =
h ˆΨ(x, t), ˆHi
=−∆x+ V (x) + 2c ˆΨ†(x, t) ˆΨ(x, t) ˆΨ(x, t). (2.12) The effect of two-body interactions is included in the last term, which is as a result of the more general two-body interaction term with phenomenological potential U (x1− x2),
ˆ Hint.= ˆ ddx1 ˆ ddx2h ˆΨ†(x1) ˆΨ†(x2)U (x1− x2) ˆΨ(x1) ˆΨ(x2) i , (2.13)
2.2. THE GROSS-PITAEVSKII EQUATION
being treated in the Born approximation. In this approximation of the scattering problem, which is treated exhaustively in the review paper by Castin [30], the relevant two-particle Schr¨odinger equation is written as an integral equation in relative coordinates, for which the iterative series solution is truncated at first order. The only point where the interaction potential U (x1−x2) plays a role at this order is under an integral sign. Whenever the potential is short-ranged compared to the average spacing between particles and collision energies are low, this integral reduces to a number a, referred to as the scattering length. This means that under these conditions, any pseudopotential reproducing the correct scattering length will have the same effect as U (x1− x2). In the one-dimensional case, it is then sufficient to replace the phenomenological potential U (x), which is often quite elaborate, by a much simpler contact pseudopotential,
U (x) → c δ(x). (2.14)
The coupling constant c is related to the scattering length as c = 2π~
2a
m , (2.15)
for which units have been temporarily restored. In higher dimensions, pseudopotentials need to contain certain regularizing operators which are of no immediate consequence to this thesis and for which the reader is referred to the explicit treatment in [30]. What matters here is that the above approximation is valid whenever the average interparticle separation, given in terms of the average density n =N
V, exceeds the scattering length, or 1
d √
n a. (2.16)
This condition defines what is known as the weakly interacting r´egime.
2.2.2
Mean field description
In the presence of a BEC, the equation of motion 2.12 can be approximated by a mean field description. Following Pitaevskii and Stringari [35], this becomes clear by expanding ˆΨ(x) in terms of bosonic annihilation operators ˆai pertaining to a basis ϕi(x):
ˆ
Ψ(x) = ϕ0(x)ˆa0+ X
i
ϕi(x)ˆai. (2.17)
The operators ˆai satisfy bosonic commutation relations: h
ˆ ai, ˆa†j
i
= δi,j, [ˆai, ˆaj] = 0. (2.18)
Their matrix elements are given by the square roots of occupation numbers Ni for the associated modes ϕi(x):
ˆ
ai|. . . , Ni, . . .i = p
Ni|. . . , Ni− 1, . . .i . (2.19) For the case of a BEC state, the zero-momentum mode ϕ0(x) becomes macroscopically occupied, meaning that the corresponding matrix element of ˆa0is much larger than its commutator, namely √
N0 1. This suggests the approximation of neglecting the commutator altogether, by replacing ˆ
a0by its relevant matrix element of size √
N0. In doing so, any quantum fluctuations caused by the nonvanishing commutator are neglected. The object ϕ0(x)ˆa0is treated as a classical, commuting field Ψ0,
ϕ0(x)ˆa0→ ϕ0(x) p
N0≡ Ψ0(x), (2.20)
which is often referred to as the condensate wave function. Indeed, the condensate wave function is the expectation value of ϕ0(x)ˆa0 on a state containing a Bose-Einstein Condensate, which will be made precise in Equation 2.24. A quantum field is thus replaced by its expectation value, which is why 2.20 is called a mean field approximation.
2.2. THE GROSS-PITAEVSKII EQUATION
After the above replacement, the collection of noncondensed particles is still described by a quantum field, given by ˆΨ0(x) ≡P
iϕi(x)ˆai. This leads to the final substitution ˆ
Ψ(x) → Ψ0(x) + ˆΨ0(x). (2.21)
If BEC in the strict sense is absent, the mean field replacement 2.21 might still be useful: whenever a quantum field can be described as a mean value with small quantum fluctuations around it, 2.21 applies. In those cases, Ψ0(x) and ˆΨ0(x) no longer serve as condensate and non-condensate fields, respectively, but acquire a new interpretation. An example of this is the one-dimensional case described in Section 2.3, in which the lowest momentum mode is not occupied macroscopically, but where particles still display collective behavior for which the mean field description leads to valuable predictions.
Having identified a classical function Ψ0(x) to approximate the field operator ˆΨ(x), the dy-namics of the mean field Ψ0(x, t) can be described by neglecting the quantum fluctuations due to the noncondensed quantum field ˆΨ0(x) altogether. Substituting definition 2.21 into Equation 2.12 and setting ˆΨ0(x) to zero yields the Gross-Pitaevskii equation, as first published in [6] and [7]:
i∂Ψ0(x, t) ∂t = −∆x+ V (x) + 2c |Ψ0(x, t)| 2 Ψ0(x, t). (2.22)
Stationary solutions to this equation are time-evolved using the chemical potential µ,
Ψ0(x, t) = Ψ0(x)e−iµt. (2.23)
The easiest way to understand this is by interpreting the condensate wave function as the expecta-tion value of ϕ0(x)ˆa0on a state containing a BEC. As observed in Section 2.1, the overlap between a BEC containing N0particles and one containing N0− 1 particles approaches 1 as N0→ ∞. For that reason, hN0| ϕ0(x)ˆa0|N0i = ϕ0(x) p N0hN0, . . . |N0− 1, . . .i | {z } ≈1 ≈ Ψ0(x). (2.24)
Then, using the Schr¨odinger picture to describe time evolution, ˆ a0|N0, ti = p N0|N0− 1, ti = p N0e−iE(N0−1)t|N0− 1i |N0, ti = eiE(N0)t|N0i . (2.25)
Recognizing the fact that E(N0) − E(N0− 1) ≈ ∂N∂E
0 = µ, one obtains Ψ0(x, t) = hN0, t| ϕi(x)ˆai|N0, ti = ϕi(x)
p
N0e−iµt≡ Ψ0(x)e−iµt. (2.26) The stationary solution 2.23 can be inserted into 2.22 to obtain the time-independent Gross-Pitaevskii equation,
−∆x+ V (x) + 2c |Ψ0(x)| 2
− µΨ0(x) = 0. (2.27)
2.2.3
Uniform ground state and collapse
In this thesis, situations will be studied where the external potential V (x) can effectively be neglected. For repulsive interactions, where g > 0, this leads to a uniform ground state solution to the Gross-Pitaevskii equation 2.22,
Ψ0(x, t) = √
n0e−iµt. (2.28)
Inserting the solution 2.28 into 2.22 shows that
2.2. THE GROSS-PITAEVSKII EQUATION
in this case. For attractive interactions, g < 0, the uniform ground state 2.28 turns out to be unstable with respect to small perturbations, which cause it to collapse to a droplet. An intuitive argument for this claim can be made by describing the condensate wave function in a density-phase representation,
Ψ(x, t) =pn(x, t)eiφ(x,t). (2.30)
Following the treatment in the review [39], one then considers small fluctuations around the mean value of the ground state density 2.28, via n(x, t) = n0+ δn(x, t). Inserting Equations 2.30 and 2.29 into the Gross-Pitaevskii equation 2.22 and retaining only terms that are linear in δn and derivatives of φ, a little algebra yields separate equations for real and imaginary parts of the linearized equation: δ ˙n = −2∆xφ (2.31) − ˙φ = ∆x δn 2n0 + gδn. (2.32)
In order to consider small fluctuations of the density, it is instructive to investigate what happens to any low momentum modes in the expansion of δn,
δn(x, t) ∼ cos (k · x − ωt) , (2.33)
with |k| 2n0|g|. For these modes, the first term on the right-hand side of Equation 2.32 can be neglected. Taking the time derivative of Equation 2.31 and inserting ˙φ from Equation 2.32 then leads to
δ¨n = 2n0g∆xδn, (2.34)
or
ω2= 2n0g k2. (2.35)
The case of g < 0 clearly leads to imaginary frequencies ω = ±ikp2n0|g| for which the lowest mode will start to grow exponentially over time:
δn ∼ ekt √
2n0|g|. (2.36)
Clearly, a negative value of g causes small perturbations in density, with Fourier modes satisfying |k| 2n0|g|, to blow up exponentially over time, leading to a spatial collapse of the initially uniform density profile.
A similar result can be obtained by considering elementary excitations in the Bugoliubov-de Gennes formalism, for which we refer to the treatment in [39]. In this formalism, the noncondensed, or excited, quantum field ˆΨ0 is treated to second order, yielding the excitation spectrum
k = p
k2+ 2n
0gk. (2.37)
For negative interactions, these energies become imaginary whenever k < 2n0|g|, in precise accor-dance with the result mentioned above. For finite systems, where momenta experience a spacing of order ∼ 1/L, a sufficiently small negative interaction can still guarantee positive energies in Equations 2.35 and 2.37 for any nonzero momentum, as long as |g| 1
2n0L. In an infinite system however, the allowed momenta form a continuum, and can be arbitrarily close to zero. This guar-antees the existence of a collapsing spatial mode satisfying k < 2n0|g| for any negative interaction strength in an infinite system.
To conclude, the Gross-Pitaevskii equation has been introduced as an effective description of an interacting BEC. This description neglects noncondensed particles and replaces the BEC field operator by its mean value, which is a commuting object. This implies that any quantum fluctu-ations, encoded in the commutation relations of the field operator, are excluded in the mean field
2.3. MEAN FIELD SOLITONS
framework of the Gross-Pitaevskii equation. For repulsive interactions and in the absence of a potential, the uniform solution has been identified as the ground state. For attractive interactions, this solution has turned out to be unstable under large-wavelength perturbations, ruling it out as a viable ground state candidate. In the next Section, the true ground state for the one-dimensional mean field scenario will be identified. It turns out that this ground state can be described in terms of a bright soliton.
2.3
Mean field solitons
For a Bose gas constrained to move in a single direction along which a potential is absent, the Gross-Pitaevskii equation 2.22 reduces to the nonlinear Schr¨odinger (NLS) equation,
i∂Ψ(x, t) ∂t = −∂ 2 ∂x2 + 2c|Ψ(x, t)| 2 Ψ(x, t). (2.38)
This is a completely integrable differential equation, which arises from an infinite-dimensional Hamiltonian system. This implies that it has an infinite number of conserved quantities, meaning that they all commute with the Hamiltonian under a poisson bracket. The first three of these quantities are the number of particles, the total momentum and the total energy itself, that is,
N = ˆ ∞ −∞ dx |Ψ(x)|2 (2.39) P = ˆ ∞ −∞ dx Ψ(x) ∂ ∂xΨ ∗(x) (2.40) E = ˆ ∞ −∞ dx ∂ ∂xΨ(x) 2 + c|Ψ(x)|4 . (2.41)
The NLS equation has been solved by the inverse scattering method. This method was origi-nally developed by Lax [40] in order to solve the Korteweg-deVries (KdV) equation, but was shown to apply to the NLS equation by Zakharov and Shabat in 1972. Their article [9] is devoted to solutions which are know as solitons. Solitons have been observed for the first time in 1834 by John Scott Russel in a canal near Edinburgh [41]. They can nowadays [42] be defined as localized waves which retain their shape over time and which, importantly, can engage in interactions with other solitons, from which they emerge without losing their shape, up to an acquired phase shift. One of the achievements of Zakharov and Shabat was to show how such interactions between mul-tiple solitons within the NLS equation can always be reduced to pairwise interactions. Although such multiple-soliton states are of great interest, this thesis will focus on a better understanding of quantum mechanical counterparts of the single soliton, which is sometimes called a solitary wave.
2.3.1
Bright solitons - mean field and experiment
In contrast to the three-dimensional case, where an attractive Bose gas tends to collapse to a droplet, the one-dimensional equation has solutions where the gas clusters into solitons, with some characteristic width ξb. For an infinite system, the single-soliton solution to the attractive NLS-equation reads as follows:
Ψbright(x, t) = √
n0sech √
n0¯c(x − vt − x0) eivx/2, (2.42) where ¯c = −c has been defined, v is the velocity of the soliton and x0is the position of its center. The density n0 at this central point can be determined by the normalization condition
N = ˆ ∞
−∞
|Ψbright(x, t)|2dx, (2.43)
which yields the final form of the single bright soliton, Ψbright(x, t) = r ¯ cN2 4 sech x − vt − x0 ξb eivx/2, (2.44)
2.3. MEAN FIELD SOLITONS
with characteristic width
ξb= 2
N ¯c. (2.45)
In the case of a finite system with periodic boundary conditions, solitons are expressed as Jacobi theta functions displaying the correct periodicity properties. The solution 2.42 can be considered as a limiting case of these finite-size expressions, for which the reader is referred to the exhaustive treatment in [43].
The energy of a bright soliton solution is negative, and can be directly found from the Hamil-tonian 2.41, Ebright= − ¯ c2N3 12 + N v2 4. (2.46)
This shows that soliton solutions are energetically favoured over the unstable uniform solution mentioned in Section 2.2.
In fact, setting v = 0 for the soliton in Equation 2.42 yields the ground state of the attractive bose gas in its mean field description. This ground state is highly degenerate: any value of the soliton center, x0, leads to an allowable ground state configuration. This raises the following question: at which position x0 in a system does the soliton materialize? Castin [44] answers this question by treating the soliton center x0 as an observable subject to a quantum mechanical uncertainty relation. He then describes the ground state as a superposition of such soliton states with center x,
|Ψbright; N i = ˆ
dx |Ψbright; N ; xi , (2.47)
such that a measurement of the system leads to a symmetry-breaking collapse to a single soliton position x0.
Since x0 is now treated as a quantum mechanical observable ˆx0, it is subject to unitary time evolution and a quantum uncertainty relation. Castin [30] shows that the uncertainty in the soliton position x0 increases with time. In his treatment, the variance of x0 becomes comparable to the soliton width, ξ, after a decay time
τdec.= √
N π ξ 2
4 . (2.48)
For completeness, a derivation of Equation 2.48 is reproduced in Appendix A. Apparently, even bright solitons arising from an integrable differential equation can give rise to finite decay times, whenever the soliton center x0 is treated as a quantum operator.
In an experimental context, a natural form of symmetry breaking is inherent to the creation pro-cess of a bright soliton. The first realization of a bright soliton in a BEC, reported by Khaykovich et al. [29], offers a clear example of this. In these experiments, a gas of7Li atoms was cooled in a series of steps and Bose condensed in an optical trap, which can be described as a superposition of harmonic oscillators in three spatial directions. Following [29], the Gross-Pitaevskii energy functional corresponding to Equation 2.22 is given by
EGP= ˆ d3x |∇Ψ(x)|2+ c|Ψ(x)|4+1 2ω 2 ⊥(x2+ y2) + ω2zz 2 |Ψ(x)|2 . (2.49)
As indicated in Equation 2.15, the coupling constant c is proportional to the scattering length a, which can be magnetically tuned using Freshbach resonance. A plot of the scattering length as a function of the applied magnetic field is given in the left panel of Figure 2.2, taken from [29]. As can be seen, a small region exists in which scattering lengths are negative.
The trapping frequencies in the x- and y- direction are chosen large enough for the gas to be essentially one-dimensional, thus preventing a collapse, whilst at the same time staying small enough to ensure the validity of the Born approximation. After Bose-Einstein Condensation has occurred, the trapping frequency in the z-direction is adiabatically turned to zero, allowing for
2.3. MEAN FIELD SOLITONS
the formation of the bright soliton solution from Equation 2.42. Clearly, translational symmetry has been broken at this point: the minimum-energy position for the bright soliton to form is determined by the center of the trap in the z-direction, here set to z0= 0.
As a consequence, Khaykovich et al. [29] propose the variational Ansatz
Ψ(x) ∼ 1 pσ2 ⊥lz sech(z/lz) exp −x 2+ y2 2σ2 ⊥ , (2.50)
which is a product of the ground state of the harmonic oscillator with the bright soliton. The latter is now centered at z = 0, due to the symmetry-breaking process imposed by the original trap. The corresponding energy graphs as a function of characteristic length lz have been plotted in [29] for different values of N |a|, as shown in the central panel of Figure 2.2. As can be seen, specific values of N |a| lead to a stability region for finite lz. It is in these regions that a bright soliton can form. Its width ξb is then associated with the characteristic length lz for which the energy attains a local minimum.
Indeed, upon release into a one-dimensional waveguide, the Bose Condensed7Li-gas of char-acteristic width lz retains its shape over the full course of the experiment, for the values of a and N |a| indicated pictorially in Figure 2.2. A bright soliton seems to have formed. This behavior is contrasted with the expansion of an ideal BEC, which disperses over time, in the absence of the nonlinear term in the nonlinear Schr¨odinger equation 2.38.
Figure 2.2: Figures taken from [29]. Left : scattering length a due to Freshbach resonance as a function of magnetic field. The black square corresponds to an ideal BEC, the white circle to an attractive bose gas and the black circle to an attractive bose gas displaying a bright soliton. Center : energy curves for the Ansatz 2.50 as a function of characteristic length, for different values of N |a|. The solid line represents a value of N |a| where the existence of a local minimum allows for the formation of a stable soliton. Right : absorption images of the 7Li BEC after its release into a one-dimensional waveguide. The top image corresponds to an ideal gas (black square in the left panel) and the bottom image shows a soliton-like BEC (black circle in the left panel). No dispersion has been observed over the course of the experiment.
2.3.2
Grey solitons - mean field and experiment
For the nonlinear Schr¨odinger equation 2.38 with repulsive interactions, that is, c > 0, the single-soliton solution was found by Tsuzuki in 1971 [8]. A year later, Zakharov and Shabat [9] obtained the same solution, along with its generalization to more complex solitons, via the inverse scattering method. The single soliton solution to Equation 2.38 reads
Ψgrey(x, t) = √ n iv cs + s 1 − v 2 c2 s tanh " (x − vt − x0)) ξg s 1 − v 2 c2 s #! , (2.51) where cs = √
2nc is the velocity of sound in the system and ξg = 1/ √
nc is the grey soliton’s characteristic width. The term grey soliton derives from the decrease in density around the
2.3. MEAN FIELD SOLITONS
soliton’s center: the soliton forms a density notch with depth n1 − v2 c2 s
. In fact, a stationary soliton with v = 0 has zero density at its central coordinate x0, which is why it is often called a dark soliton. As the soliton velocity v increases from zero to the speed of sound, the soliton’s density notch fills up completely, until it dissolves into a flat profile at v = cs. Again, it should be stressed that the soliton 2.51 is an infinite-size result. Finite-size solitons with the correct periodicity are expressed in terms of Jacobi theta functions, which have been treated in detail by Carr et al. in [45].
Using the NLS momentum 2.40 and Hamiltonian 2.41, the soliton’s energy can be calculated:
Egrey = 2nξc2 s 3 1 − v 2 c2 s 3/2 + ncN. (2.52)
Apparently, solitons with higher velocities, and therefore smaller notch depths, are energetically favored. Following [39], this image becomes more pronounced when considering the soliton’s momentum, p = ∂Egrey ∂v = −2nξ s 1 −v 2 c2 s v. (2.53) The quantity −2nξq1 −v2 c2 s
can be interpreted as a negative mass, with a modulus which decreases for increasing velocities. Amusingly, a dark soliton is in some respects the opposite of a relativistic particle: it is an object which can decrease its energy by increasing its velocity, thereby lowering its mass, until it disappears when it reaches the speed of sound. Such dissipative dynamics, including interactions with the noncondensed fraction of the gas, have been described by Fedichev et al. in [46]. It will be clear by now that in contrast with bright solitons, grey solitons do not constitute the ground state of the system.
A way towards experimental realization of grey solitons is suggested by considering a density-phase representation, where ΨGrey(x, t) =p|Ψgrey(x, t)|2eiφ(x,t). In this notation, the density and phase become |Ψgrey(x, t)|2= n 1 − 1 − v 2 c2 s sech2 " (x − vt − x0)) ξg s 1 − v 2 c2 s #! , (2.54) φ(x, t) = π 2 − arctan cs v s 1 − v 2 c2 s tanh " (x − vt − x0)) ξg s 1 − v 2 c2 s #! . (2.55)
As can be seen, grey solitons experience a phase jump in a limited region around x = x0+ vt, which grows monotonously with v. For this reason, grey solitons are sometimes referred to as phase kinks.
The first experimental realization of grey solitons, reported by Burger et al. [28], was achieved by enforcing such a phase kink on a BEC of 87Rb atoms. This method, often called phase imprinting, was first proposed in [47] and relies on the fact that the local phase of any wave function evolves by application of the exponentiated Hamiltonian,
Ψ(x, t) = Tt h e−i ´t 0dt 0H(x,t0)i Ψ(x, 0), (2.56)
where Ttis the time-ordering operator. This means that subjecting different regions of a BEC to different Hamiltonians for a time τ can lead to a relative change in phase. Explicitly, one applies
H(x, t) = (
H0(x, t), in region A,
H0(x, t) + V (x, t)ϑ(τ − t)ϑ(t), in region B,
(2.57)
where ϑ denotes the Heavyside step function. If V (x, t) acts in the same way on all relevant eigenstates of H0 and the duration of the perturbation τ is smaller than the other timescales in the system governed by H0, the effect of V (x, t) can be seen as addition of a phase φ ≈
´τ
0 V (x, t)dt to the wave function of region B.
2.4. SUMMARY AND LOOK AHEAD
In the experiments performed by Burger et al. [28], the perturbing potential V was supplied by a far detuned laser beam shining on half of the condensate. Its exposure time was chosen as τ =cs
ξ, which is the time it will take the fastest possible density wave to traverse the correlation length of the system. As can be seen in Figure 2.3, taken from [28], phase imprinting the condensate causes a grey soliton-like object to emerge. The object disappears after a time of approximately ∼ 10 ms. The authors explain this as a form of dissipation, caused by interactions with the thermal cloud, which breaks the integrability of the NLS equation. As mentioned before, the details of such disipative dynamics have been described by Fedichev et al. in [46].
In our experimental setup (see [16]), condensates con-thermal cloud, are produced every 20 s. The fundamen-Hz along the axial and radial
Figure 2.3: Images taken from [28]. Left : absorption images of a87Rb BEC after phase imprinting and a short time-of-flight period. Grey soliton structures can be seen to traverse the condensate, which vanish due to dissipation after a period of ∼ 10 ms. Right : numerical simulations performed by Burger et al. [28] showing the emergence and time evolution of a grey soliton-like object after phase imprinting.
2.4
Summary and look ahead
One-dimensional systems do not allow for Bose-Einstein Condensation in the strict sense of a macroscopic occupation number at zero momentum. Nevertheless, both experiments and mean field treatments in the absence of a potential point to a macroscopic occupation of more compli-cated states, which behave like solitons.
2.4.1
Overview
In the attractive case, the mean field ground state is formed by a bright soliton, with density
|Ψbright(x, t)|2= ¯ cN2 4 sech 2 ¯cN 2 (x − vt − x0) . (2.58)
Here, the soliton center at t = 0 is given by x0. In a translationally-invariant context, each measurement will correspond to a different value of x0, breaking the translational symmetry of the system. Treating x0as a quantum mechanical observable in this way gives rise to finite lifetimes of the bright soliton’s expectation value after the translational symmetry has been broken: the variance of ˆx0will become comparable to the width of the soliton after a time
τdec.= √
N π ξ 2
4 . (2.59)
This time is designated as a decay time, because the spreading of hˆx0i over time will manifest itself in similar spreading of the expectation value|Ψbright(x, t)|2 , as the central coordinate x0 in Equation 2.58 becomes more and more uncertain over time.
2.4. SUMMARY AND LOOK AHEAD
The repulsive case gives rise to excited states in the form of grey solitons, with density
|Ψgrey(x, t)|2= n 1 − 1 − v 2 c2 s sech2 " (x − vt − x0)) ξg s 1 −v 2 c2 s #! . (2.60)
Apart from the above density profile, this grey soliton is characterized by a phase jump from 0 to π across its central coordinate x0.
2.4.2
Goal
With the above results in mind, the question from Section 1 can be made more precise. It should read: which exact quantum mechanical wave functions reproduce |Ψbright(x, t)|2 or |Ψgrey(x, t)|2 whenever the NLS equation becomes valid? To answer this question, this project focuses on the diagonal part of the one-body density operator, introduced in Equation 2.6:
hˆρ(x, t)i ≡D ˆΨ†(x, t) ˆΨ(x, t)E. (2.61) The NLS description becomes valid whenever the mean field expansion ˆΨ(x) → Ψ0(x) + ˆΨ0(x) can be made. To zeroth order in this expansion, Equation 2.61 becomes
D ˆΨ†(x, t) ˆΨ(x, t)E≈ |Ψ(x, t)|2. (2.62) This project’s goal can therefore be phrased as follows:
Goal. This project aims to find exact quantum states |Ψi for which hΨ| ˆΨ†(x, t) ˆΨ(x, t)|Ψi behaves like |Ψbright(x, t)|2 or |Ψgrey(x, t)|2.
In order to gather the necessary ingredients to reach this goal, it is useful to zoom in on the expectation value 2.61. Expanding a generic state in a specific basis, labeled by µ and λ,
|Ψi =X λ
Cλ|λi , (2.63)
the expectation value 2.61 can be expanded as hΨ| ˆΨ†(x, t) ˆΨ(x, t)|Ψi =X
µ,λ
eix(Pλ−Pµ)e−it(Eλ−Eµ)C
µCλ∗hλ| ˆΨ†(0, 0) ˆΨ(0, 0)|µi , (2.64)
where translation operators ei ˆP xand unitary time evolution e−i ˆHthave been used. To find states |Ψi for which Equation 2.64 takes on a solitonic density profile, a number of observations should be made in advance:
• Equation 2.64 can be read as a time- and space-dependent Fourier transform of the quantity hµ| ˆΨ†(0, 0) ˆΨ(0, 0)|λi ≡ hµ| ˆρ(0, 0)|λi , (2.65) which is called the density form factor. It is therefore natural to look for a form factor which behaves like the Fourier transform of a soliton profile. Alternatively, one could look for coefficients Cµ,λ so that the same can be said about the product CµCλhµ|ˆρ(0, 0)|λi. • The energies Eλ,µand momenta Pλ,µare not independent, in general. A nontrivial dispersion
relation Eλ(Pλ) leads to a decay of the density profile 2.64, whenever different spatial modes move at different phase velocities. To control such processes, expressions for these momenta and energies need to be found in an appropriate basis of states.
2.4. SUMMARY AND LOOK AHEAD
• A way to achieve stability is by looking for states which all share the same energy, Eλ = Eµ for all µ and λ. Time evolution is then completely trivial. A stable, moving density distribution can be obtained by restricting oneself to states with a linear dispersion, for which Eλ= vPλ. In that case,
eix(Pλ−Pµ)e−it(Eλ−Eµ)= ei(x−vt)(Pλ−Pµ), (2.66) meaning that hΨ| ˆΨ†(x, t) ˆΨ(x, t)|Ψi moves coherently.
As it turns out, there is a basis in which the above energies, momenta and form factors can be obtained explicitly. This basis is obtained in the next chapter, using the Bethe Ansatz.
Chapter 3
The quantum nonlinear
Schr¨
odinger equation
The one-dimensional bose gas can be described by the field theoretical Lieb-Liniger Hamiltonian, ˆ HLL= ˆ dx ∂ ∂x ˆ Ψ†(x) ∂ ∂x ˆ Ψ(x) + c ˆΨ†(x) ˆΨ†(x) ˆΨ(x) ˆΨ(x) , (3.1)
which is the one-dimensional version of Equation 2.11 in the absence of an external potential. The Hamiltonian 3.1 describes a single species of bosons constrained to move along a line, interacting through a pseudopotential of strength c that has zero range. As explained in Section 2.2, this pseudopotential is realistic as long as the Born approximation is valid. The bosonic quantum fields ˆΨ(x) in Equation 3.1 obey canonical equal-time commutation relations
h ˆΨ(x), ˆΨ†(y)i= δ(x − y), h ˆΨ(x), ˆΨ(y)i=h ˆΨ†(x), ˆΨ†(y)i= 0, (3.2) and they act on a Fock vacuum |0i and its dual h0| which are defined as
ˆ
Ψ(x) |0i = 0 = h0| ˆΨ†(x), h0|0i = 1. (3.3)
The Heisenberg equation of motion corresponding to the field theoretical Hamiltonian 3.1 is i∂ ∂t ˆ Ψ(x, t) =h ˆΨ(x, t), ˆHLL i = − ∂ 2 ∂x2 + 2c ˆΨ †(x, t) ˆΨ(x, t) ˆ Ψ(x, t). (3.4)
This equation is often referred to as the quantum nonlinear Schr¨odinger equation. In contrast to the classical nonlinear Schr¨odinger equation 2.38, which describes the mean value of ˆΨ as a classical function, Equation 3.4 is an operator equation. It respects the bose field’s commutation relations, including the resulting uncertainty relations between observables. In this sense, Equation 3.4 is an exact equation, within the range of validity of the Born approximation. It governs the full quantum field theory of bosons in one dimension, beyond any mean field description.
3.1
Exact eigenstates
Any eigenstate of the field theoretical Hamiltonian 3.1 simultaneously diagonalizes the particle number operator ˆN and the total momentum operator ˆP ,
ˆ N = ˆ dx ˆΨ†(x) ˆΨ(x), (3.5) ˆ P = −i 2 ˆ dx ˆ Ψ†(x) ∂ ∂x ˆ Ψ(x) − ∂ ∂x ˆ Ψ†(x) ˆ Ψ(x) , (3.6)
3.1. EXACT EIGENSTATES
as can be seen by a direct computation of the commutatorsh ˆHLL, ˆN i
=h ˆHLL, ˆP i
= 0. A general eigenstate |Ψi of ˆHLLis therefore characterized by a fixed number of particles N . Parametrizing |ΨNi by of a generic set of N complex numbers, {λ1, . . . , λN}, it can be expressed in terms of a complex-valued amplitude χN(x1, . . . , xN|λ1, . . . , λN), called the many-body wave function:
|ΨN(λ1, . . . , λN)i = ˆ
dx1. . . dxNχN(x1, . . . , xN|λ1, . . . , λN) ˆΨ†(x1) . . . ˆΨ†(xN) |0i . (3.7) By applying the commutation relations 3.2, the time-independent Schr¨odinger equation
ˆ
HLL|ΨN(λ1, . . . , λN)i = EN|ΨN(λ1, . . . , λN)i (3.8) is seen to be equivalent to the quantum mechanical problem for the many-body wave function,
HLLχN(x1, . . . , xN|λ1, . . . , λN) = ENχN(x1, . . . , xN|λ1, . . . , λN), (3.9) where the quantum mechanical version of ˆHLLis given by the Lieb-Liniger Hamiltonian
HLL= N X j=1 − ∂ 2 ∂x2 j + 2cX k<l δ(xk− xl). (3.10)
The total momentum operator is expressed in this quantum mechanical notation as P = −iX
j ∂ ∂xj
. (3.11)
3.1.1
The Bethe Ansatz
Strikingly, the Lieb-Liniger Hamiltonian 3.10 can be diagonalized exactly by means of the Bethe Ansatz. This method was developed by Hans Bethe in his famous 1931 article on spin chains [10]. In 1963, Elliott Lieb and Werner Liniger showed [11, 12] how Bethe’s approach can be applied to the one-dimensional bose gas Hamiltonian 3.10, which since carries their name. Their method is reviewed here, in a form which strongly owes to the subsequent work of Gaudin [48] and Korepin, Izergin and Bogolyubov [15].
Consider the Lieb-Liniger Hamiltonian defined on an indefinite axis, which will later be given a length L with periodic boundary conditions. An essential ingredient in solving the system is the realization that away from the points xj = xk with j 6= k where particle positions coincide, HLL reduces to a system of free bosons, solved by plane waves. In fact, the system defined by equations 3.9 and 3.10 is equivalent to the system defined by the free particle Schr¨odinger equation
N X j=1 −∂ 2 ∂x2 j χN = ENχN, (3.12)
on the open domain
DN : x1< x2< . . . < xN, (3.13) subject to the N − 1 boundary conditions
∂ ∂xj+1 − ∂ ∂xj χN x j+1−xj=0+ = c χN ∀ j < N. (3.14)
Due to bosonic exchange symmerty, any solution valid on DN can be straightforwardly ex-tended to any other open region with a different particle ordering. The boundary conditions 3.14 are obtained for each neighbouring pair of particles (j, j + 1) separately, by transforming to center of mass and relative coordinates, xcm=
xj+1+xj
3.1. EXACT EIGENSTATES
equation 3.9 with respect to xrelover a vanishing interval [−, ] around xrel= 0 then leads to the condition −1 2 ∂ ∂xj+1 − ∂ ∂xj χN xj+1−xj=0+ xj+1−xj=0− + c χN = 0. (3.15)
Here, the first term stems from the fundamental theorem of calculus, the second term is a result of the delta function in the Lieb-Liniger Hamiltonian and the vanishing right-hand side is a consequence of the continuity of the wave function at xj+1= xj. Since the bosonic wave function χN is symmetric under the exchange of xj+1and xj and the differential operator in Equation 3.15 is antisymmetric under the same operation, Equation 3.14 is obtained.
Now that the Lieb-Liniger model has been reduced to a free model on N ! open domains subjected to N − 1 boundary conditions 3.14, the Bethe Ansatz can be formulated. It describes a wave function as a sum over may-body plane waves, weighed by amplitues A which depend on elements π of the symmetric group SN with N ! elements. In short, the Bethe wave function reads
χN(x1, . . . , xN|λ1, . . . , λN) = X
π∈SN
A(π)eiPjxjλπ(j). (3.16)
The Schr¨odinger equation 3.12 on the open domain 3.13 is satisfied by the Bethe wave function 3.16 with energy E = N X j=1 λ2j, (3.17)
whilst the total momentum of the wave function is simply given by
P = N X
j=1
λj. (3.18)
Owing to this relation between the set of spectral parameters {λi} and the total momentum in the system, the complex numbers λi are often called rapidities, or quasimomenta.
3.1.2
Scattering phases
The boundary conditions 3.14 have to be satisfied by an appropriate choice for the coefficients A(π). To find these coefficients, one observes that any two elements π and ˜π of the symmetric group SN can be related to each other by a finite application of nearest-neighbour exchanges σj ∈ SN: ˜ π = Y j∈Σ σj π. (3.19)
Here, σj is used as a shorthand for σ(j)↔(j+1), the element of SN exchanging positions j and j + 1. This subscript notation should not be confused with the parentheses in the notation π(j), which denotes some element π ∈ SN acting on index j. Σ is some finite sequence of integers smaller than N , describing a sequence of nearest-neighbour exchanges bringing about the transformation from π to ˜π. Such a sequence is not unique.
The existence of a decomposition 3.19 means that any ratio between coefficients can be ex-pressed as a product of ratios of coefficients differing by a nearest neighbour exchange:
A(π) A(˜π) = A(π) A(σj1π) A(σj1π) A(σj2σj1π) . . .A(σjk−1. . . σj2σj1π) A(˜π) . (3.20)
Altough a decomposition such as Equation 3.20 is non-unique, any choice of a series of nearest neighbour exchanges Σ should lead to the same ratio A(π)A(˜π). This nontrivial requirement is an incar-nation of the Yang-Baxter equation, described in Section 3.1.3. As a consequence of Equation 3.20,
3.1. EXACT EIGENSTATES
any ratio between coefficients is known whenever an expression can be found for A(π)/A(σjπ), the ratio between a generic coefficient A(π) and a coefficient A(σjπ) differing by a nearest neighbour exchange at a generic particle index j.
It turns out that the required generic ratio, A(π)/A(σjπ), is exactly the information encoded in the jthboundary condition from Equation 3.14. Applying this boundary condition to the Bethe wave function 3.16 yields
i λπ(j+1)− λπ(j) − c A(π) + i λπ(j)− λπ(j+1) − c A(σjπ) = 0, (3.21) as can be seen by equating those terms for which π and ˜π differ by a single exchange (j) ↔ (j + 1) in Equation 3.16, after application of 3.14. One thus finds the ratio
A(π) A(σjπ) = −c − i λπ(j)− λπ(j+1) c + i λπ(j)− λπ(j+1) = −e −iφ(λπ(j)−λπ(j+1)), (3.22) where φ(λ) = 1 iln c + iλ c − iλ. (3.23)
The function φ equals
φ(λ) = 2 arctanλ
c, (3.24)
if λ is real-valued. It describes the phase acquired by the wave function whenever two particles pass one another. For this reason, it is referred to as a scattering phase.
3.1.3
The Yang-Baxter relation
Now that the ratios between coefficients differing by a nearest-neighbour exchange have been found, decompositions like 3.20 fix the coefficients A(π) up to an overall phase factor. This is a nontrivial statement: apparently, any scattering process in the Lieb-Liniger model can be decomposed into subsequent two-body scattering processes. As mentioned, such decompositions are non-unique. Using expression 3.22, one can check that different decompositions lead to the same total scattering phase. To this end, it suffices to consider the simplest nontrivial case, namely the three-body ratio
A(λ1λ2λ3) A(λ3λ2λ1) =A(λ1λ2λ3) A(λ1λ3λ2) A(λ1λ3λ2) A(λ3λ1λ2) A(λ3λ1λ2) A(λ3λ2λ1) =A(λ1λ2λ3) A(λ2λ1λ3) A(λ2λ1λ3) A(λ3λ2λ1) A(λ2λ3λ1) A(λ3λ2λ1) , (3.25) where the particle ordering corresponding to each amplitude is indicated explicitly. Equation 3.25 is a manifestiation of the Yang-Baxter relation, which is a central relation in exactly solvable mod-els. It implies that any scattering process can be decomposed into two-body scattering processes and that any order of such two-body processes leads to the same final scattering phase.
Due to the conservation laws 3.5 and 3.6, any two-body scattering process conserves the set of spectral parameters {λi} characterizing the wave function. The Yang-Baxter relation ensures that the same is true for any scattering process, as can be seen by decomposing it into subsequent two-body processes. This property was described by Bill Sutherland [49] as nondiffractive scattering: a scattering process which merely leads to a redistribution of a fixed set of M rapidities {λi}
M i=1 amongst the M particles involved. In such a scattering process, any symmetric or antisymmetric function of the rapidities is conserved. Nondiffractive scattering is therefore closely connected to the existence of an infinite set of conserved charges, thus strongly constraining the dynamics of the model.
The Yang-Baxter relation, nondiffractive scattering and the existence of an infinite tower of conserved charges are often quoted as requirements for a quantum many-body problem to be quantum integrable. This term is inspired by the close analogy with classical integrability, which was mentioned in Section 2.3, and which is similarly characterized by a large number of conser-vation laws and a Yang-Baxter relation. However, the notion of quantum integrability is not as well-defined as its classical counterpart [50]. The relation between the algebraic structures in the two realms is an interesting one (see, for instance [51] and [52]), which will not be pursued in this thesis.
3.1. EXACT EIGENSTATES
3.1.4
Bethe Ansatz wave function for the Lieb-Liniger model
The form of the wave function on DN, including the correct boundary conditions, is now entirely fixed by the nearest neighbour relation 3.20, up to normalization and an overall U (1) symmetry. Using the antisymmetry of the scattering phase φ(−λ) = −φ(λ), one can check that the following form of the coefficients satisfies Equation 3.20 on DN:
A(π) = (−1)[π]e2i P
N ≥j>k≥1φ(λπ(j)−λπ(k)) for x ∈ D
N, (3.26)
where (−1)[π] is the parity of π. Next, this result is extended to all of RN. Here one should take into account that the left-hand side of the boundary condition 3.14 changes sign whenever xj and xj+1are exchanged. This can be accounted for by multiplying each scattering phase by a factor sgn(xj− xk) in 3.26. At the same time, bosonic exchange statistics can be guaranteed by multiplying the entire wave function by a factorQ
N ≥j>k≥1sgn(xj− xk), leading to χN(x1, . . . , xN|λ1, . . . , λN) = N Y N ≥j>k≥1 sgn(xj− xk) X π∈SN (−1)[π]e2i P N ≥j>k≥1sgn(xj−xk)φ(λπ(j)−λπ(k))eiPjxjλπ(j). (3.27)
The norm N can be computed using a Jacobian determinant due to Gaudin [48], in a conjecture which was proven by Korepin et al. [15]
For positive values of the interaction strength in a finite system, the coefficients A(π) offer an interpolation between a permanent at c = 0,
χN(x1, . . . , xN|λ1, . . . , λN) = N X
π∈SN
eiPjxjλπ(j), (3.28)
and a determinant for c−1= 0, symmetrized to bosonic statistics: χN(x1, . . . , xN|λ1, . . . , λN) = N Y N ≥j>k≥1 sgn(xj− xk) X π∈SN (−1)[π]eiPjxjλπ(j). (3.29)
In a system of length L, the permanent in Equation 3.28 is the wave function for free bosons, where momenta will be given by bosonic quantum numbers Jj:
λj= 2π
LJj, Jj ∈ Z, c = 0. (3.30)
The determinant in equation 3.29 is the wave function of the Tonks-Girardeau gas, proposed in a different notation by Marvin Girardeau in 1960 [53], and described by the Lieb-Liniger model in its c → ∞ limit. As will become clear in Section 3.2, rapidities are linearly related to integers (half-odd integers) Ij for an odd (even) number of particles in this r´egime:
λj = 2π
L Ij, c → ∞. (3.31)
The quantum numbers Ij are distinct for different values of the particle index j. In fact, j 6= k implies λj 6= λk for any value of c 6= 0: away from the free boson limit given by Equation 3.28, the wave function 3.27 is antisymmetric under the exchange of two momenta. As a consequence, one obtains the following theorem:
Theorem 1. For c 6= 0, the Bethe wave function χN, given by Equation 3.27, vanishes whenever two rapidities coincide.
The Bethe wave function thus displays a type of fermionic exclusion principle, which exem-plifies a broader tendency of one-dimensional systems to mix bosonic and fermionic behavior. Heuristically, this can be seen as a result of the one-dimensional geometry, in which particles cannot be exchanged without passing through each other, meaning that scattering phases and statistical phases cannot be separately considered.
3.2. EXACT SPECTRUM
3.2
Exact spectrum
3.2.1
Bethe Equations
A description of the behavior of rapidities away from the free boson limit 3.30 and the Tonks-Girardeau limit 3.31 can be given by placing the system on a finite interval L and applying periodic boundary conditions,
χN(x1, . . . , xj−1, xj+ L, xj+1, . . . , xN|λ) = χN(x1, . . . , xj−1, xj, xj+1, . . . , xN|λ), (3.32) for each value of j. Using Equation 3.27, these conditions immediately yield
eiλjL= (−1)N −1e−iPNk=1φ(λj−λk)=Y k6=j
λj− λk+ ic λj− λk− ic
, for j = 1, . . . , N. (3.33)
The N equations above are called the Bethe equations. They offer nonlinear and fully coupled relations between all rapidities in the system. As such, the Bethe equations are highly nontrivial: although the Hamiltonian 3.10 only contains interactions at zero range, the spectral parameters in the Bethe Ansatz framework are all sensitive to each other in a nonlinear way. Heursitically, the boundary condition 3.32 can be seen as the following requirement: moving particle j over a full circle around the periodic system, thus scattering it across all other particles, should lead to the same result as when it had remained at its initial position, xj. This means that the additional factor of eiλjL in the Bethe wave function 3.27 has to be canceled exactly by the product of N − 1 scattering phases 3.22, acquired in the process. This statement made precise by Equation 3.33.
Taking the logarithm of the Bethe Equations 3.33 leads to the illuminating form
λj= 2π LIj− 1 L N X k=1 φ(λj− λk). (3.34)
The numbers Ij result from the periodicity of the complex exponential function. They are integers (half-odd integers) Ijfor an odd (even) number of particles, due to the sign (−1)N −1on the right-hand side of Equation 3.33. From Equation 3.34, it becomes clear that the Bethe equations constitute a nonlinear mapping from the set of numbers {Ij} to the set of rapidities {λj}. A number of theorems have been proven about this mapping, which will be described in the next chapter.
3.2.2
Bethe roots for general, positive coupling
In the work by Korepin et al. [15], the following theorems about the mapping between a set of quantum numbers {Ij} and the set of rapidities {λj} are proven:
Theorem 2. All solutions λj of the Bethe equations 3.34 are real numbers for c > 0.
Theorem 3. For c > 0, the solutions to the Bethe equations 3.34 exists and are uniquely deter-mined by a proper set of quantum numbers {Ij}.
Because of Theorem 3, the set {Ij} will be referred to as the quantum numbers describing a Bethe wave function.
Theorem 4. For c > 0, Ij> Ik implies λj > λk and Ij = Ik implies λj= λk.
Theorem 5. For c > 0, any difference between rapidities corresponding to Ij< Ik is bounded as follows: 2π L Ij− Ik 1 +2nc ≤ |λj− λk| ≤ 2π L(Ij− Ik). (3.35)
3.2. EXACT SPECTRUM
Theorem 6. For c > 0, the energy functional given by Equation 3.17 is minimized by the choice of quantum numbers
IjG.S.= − N − 1 2
+ j − 1, with j = 1, . . . , N. (3.36)
The corresponding set of rapidities {λj} constitutes the ground state of the Lieb Liniger model. Due to its similarity to the ground state of fermionic systems, this state is referred to as a Fermi sea.
For the proofs of these theorems, the reader is referred to the explicit treatment in [15]. Heuristically, the above theorems state that for c > 0, the Bethe equations 3.34 offer an injective, monotonously increasing mapping from a proper set of quantum numbers {Ij} to a set of rapidities {λj}. In Figure 3.1, an example of such a mapping is given. Theorem 6 is immediate for the Tonks-Girardeau gas, via Equation 3.31. The Theorem then follows for generic c > 0 from the non-intersecting nature of the roots λj, which is proven by Theorem 5. Furthermore, Theorem 5 states that the minimum spacing between consecutive rapidities vanishes as c tends to zero. This limit will be considered in the next Section.
It is appropriate to close this Section with the following theorem:
Theorem 7. The collection of Bethe wave functions corresponding to all possible proper sets of integer (half-odd integer) quantum numbers for an odd (even) number of particles constitutes a complete basis for the Hilbert space pertaining to the c > 0 Lieb-Liniger Hamiltonian.
This theorem was proven in [13], starting from the Tonks-Girardeau limit and extending to arbitrary coupling via a continuity argument.
B
B
C
C
A
A
- 20 -10 10 20 30 40 50c - 5 5 10 Re@ ΛDFigure 3.1: Example of the mapping, defined by the Bethe Equations 3.34, between some set of quantum numbers {Ij} (left ) and the corresponding set of Bethe roots {λj} for different values of the coupling c (right ). At c → ∞ (region C ), a regular spacing at positions λj = 2πLIj can be observed, whereas roots are seen to accumulate around bosonic quantum numbers nj via λj = 2πLnj as c tends to zero (region B ). For negative coupling (region A), nothing is depicted, since rapidities attain complex values there, as described in Section 3.2.4.
3.2.3
Bethe roots for small coupling
From the logarithmic Bethe equations 3.34, it is clear that the Tonks-Girardeau result for rapidities given by Equation 3.31 follows from the vanishing of the scattering phase as c tends to infinity. For the case where c tends to zero, the Bethe equations can be reformulated to arrive at a similar result. To see this, one rewrites the scattering phase for c > 0 as
φ(λ) = 2 arctanλ
c = sgn(λ)π − 2 arctan c
