A coordinate Bethe ansatz approach to the calculation of equilibrium and nonequilibrium correlations of the one-dimensional Bose gas

1.2 MB·

by  J. C. Zill

A coordinate Bethe ansatz approach to the calculation of equilibrium and
nonequilibrium correlations of the one-dimensional Bose gas
Jan C Zill,1, Tod M Wright,1 Kar´en V Kheruntsyan,1 Thomas Gasenzer,2,3 and Matthew J Davis1,4
∗
1School of Mathematics and Physics, The University of Queensland, Brisbane QLD 4072, Australia
2Kirchhoff-Institut fu¨r Physik, Universita¨t Heidelberg, Im Neuenheimer Feld 227, 69120 Heidelberg, Germany
3ExtreMe Matter Institute EMMI, GSI Helmholtzzentrum fu¨r Schwerionenforschung, 64291 Darmstadt, Germany
4JILA, University of Colorado, 440 UCB, Boulder, Colorado 80309, USA
(Dated: April 15, 2016)
6 We use the coordinate Bethe ansatz to exactly calculate matrix elements between eigenstates of
1 the Lieb–Liniger model of one-dimensional bosons interacting via a two-body delta-potential. We
0 investigate the static correlation functions of the zero-temperature ground state and their depen-
2 denceoninteractionstrength,andanalyzetheeffectsofsystemsizeinthecrossoverfromfew-body
tomesoscopicregimesforuptosevenparticles. Wealsoobtaintime-dependentnonequilibriumcor-
r
p relation functions for five particles following quenches of the interaction strength from two distinct
A initial states. One quench is from the non-interacting ground state and the other from a correlated
ground state near the strongly interacting Tonks-Girardeau regime. The final interaction strength

4 and conserved energy are chosen to be the same for both quenches. The integrability of the model
1 highly constrains its dynamics, and we demonstrate that the time-averaged correlation functions

following quenches from these two distinct initial conditions are both nonthermal and moreover
]
s distinct from one another.
a
g PACSnumbers: 02.30.Ik,67.85.-d,05.30.Jp
-
t
n
I. INTRODUCTION correlations at large distances, with an exponent given
a
u in terms of the equation of state that is known exactly
q TheLieb–Linigermodelofaone-dimensionalBosegas from the thermodynamic Bethe ansatz [16]. The alge-
t. with repulsive delta-function interactions is a paradig- braic Bethe ansatz provides a determinantal representa-
a matic example of an exactly solvable continuous, inte- tion of correlations, from which their asymptotic behav-
m
grable many-body quantum system [1]. In particular, it ior can be extracted [17]. More recently, exact expres-
d- has served as the context for the development of the- sions for local second- and third-order correlations [18–
n oretical tools that have subsequently been widely ap- 20], together with exact results for the one-body corre-
o plied in the study of integrable systems, such as the lation function at asymptotically short distances [21] in
c so-called “thermodynamic Bethe ansatz” functional rep- terms of the equation of state have been derived.
[
resentation, which provides the exact equation of state, Away from the asymptotic short- and long-range

2 excitation spectrum [1], and bulk parameters [2] of the regimes, the behavior of correlation functions is less well
v system in the thermodynamic limit. However, the cal- known. For intermediate interaction strengths and arbi-
4 culationofcorrelationfunctionsfromtheexactsolutions trary length scales one must resort to numerics to deter-
3
provided by the Bethe ansatz is notoriously difficult. minethecorrelationfunctions. Resultsforthelatterhave
4
At zero temperature, exact closed-form solutions for been obtained using numerical methodologies including
0
0 some equilibrium correlation functions are known in the quantumMonteCarlo[22,23],anddensitymatrixrenor-
. Tonks–Girardeaulimitofinfiniteinteractionstrength[3– malization group approaches [24]. A recently developed,
1
7]. This comparatively tractable limit also allows for integrability-basedapproachcombinesthedecomposition
0
6 some strong-coupling expansion results for large but fi- of correlation functions into sums over matrix elements
1 nite interactions [7–10]. In the opposite weakly inter- (form factors)of certainsimple operatorsbetweenBethe
: acting quasi-condensate regime, a mean-field approach ansatz eigenstates [25, 26]. This approach has generated
v
can be used to describe the system [11] and a Bogoli- results,forexample,forstaticanddynamicalequilibrium
i
X ubov method can be used to determine the low-lying ex- correlationsatzeroandfinitetemperatureforsystemsof
r citation spectrum [12], relying on small density fluctua- up to N 100 particles [27]. Other finite temperature
a tions. Fewer results are available for intermediate inter- results for≈correlation functions have been obtained us-
action strengths, away from the strongly-interacting and ing imaginary time stochastic gauge methods [28, 29],
weakly-interactingregimes. ThedevelopmentoftheLut- taking the non-relativistic limit of a relativistic field the-
tinger liquid description of quantum fluids [13] and the ory [30], utilizing Fermi–Bose mapping for the strongly
related formalism of conformal field theory [14, 15] have interactinggas[9,31,32],employingperturbativeexpan-
lead to the prediction of power-law scaling for first-order sionsintemperatureandinteractionstrength[33],aswell
as combining the thermodynamic Bethe ansatz with the
Hellmann–Feynman theorem [34].
∗ [email protected] Experiments with ultracold quantum gases are able
2
to realize effectively one-dimensional systems by tightly systems of only a few particles. However for small par-
confining the gas in two of the three spatial dimen- ticle numbers N 7 we obtain numerically exact re-
≤
sions, either using optical lattice potentials or atom-chip sults for ground-state correlations, which are described
traps [35–50]. These experiments are now probing the in Sec. III. Our results demonstrate that local correla-
predictions of the Lieb–Liniger model. The configura- tions in the strongly interacting regime are already close
bility of quantum-gas experiments allows for so-called to their thermodynamic-limit values for these few-body
quenches of the system, in which Hamiltonian param- to mesoscopic systems.
eters of the system are abruptly changed, and thus for An additional advantage of our methodology is that it
the study of the Lieb-Linger model out of equilibrium, can also calculate overlaps between Lieb–Liniger eigen-
providing even greater challenges for theory. states corresponding to any two interaction strengths,
The dynamically evolving correlations of the Lieb– whichallowsustostudythedynamicsofquenchesofthe
Liniger gas in nonequilibrium scenarios are currently a interactionstrengthbetweenarbitraryvalues. InSec.IV
topic of significant interest, and a number of theoretical we utilize this property to study the effects of integrabil-
approacheshavebeenapplied. Notableexamplesinclude ity on the relaxation of the Lieb-Liniger model following
exactdiagonalizationunderalowmomentumcutoff[51– such a quench. In particular, we compare two nonequi-
55], mapping of the hard-core Tonks–Girardeau gas to librium quench scenarios with the same final Hamilto-
free spinless fermions [56–63], phase-space methods [64], nian and state energy, but beginning from starkly dif-
dynamicBogoliubov-likeapproximations[65]andtensor- ferent initial states. Statistical mechanics would predict
network methods [66, 67]. References [68–71] employed thatthesystemwouldrelaxtothesamethermalstatein
nonperturbativeapproximativefunctional-integralmeth- bothcases,butduetotheintegrabilityoftheLieb-Linger
ods, while in Ref. [72] a dynamical Luttinger-liquid ap- modelnotonlyarethetime-averagedstatesfollowingthe
proachwastaken. Othercalculationsmakeexplicituseof two quenches non-thermal, they are also distinct. Af-
the integrability of the system. These are based on vari- tercharacterizingandcomparingthenonequilibriumdy-
ousBetheansatzapproaches,andincludeutilizingFermi- namics following both quenches, we conclude in Sec. V.
Bose mapping [73, 74] and strong coupling expansions of
the coordinate Bethe ansatz wave function [75–77], com-
II. COORDINATE BETHE-ANSATZ
bining the algebraic Bethe ansatz with other numerical
METHODOLOGY
methods [78–80], and using the Yudson contour-integral
representation for infinite-length systems [81, 82]. Re-
cently,itwasconjecturedthatthedynamicsfollowingan A. Lieb–Liniger model eigenstates
interaction strength quench are captured by a thermo-
dynamic Bethe ansatz saddle point state and excitations The Lieb–Liniger model [1] describes a system of N
around it — the so-called quench action approach [83– indistinguishable bosons subject to a delta-function in-
88]. In the spirit of the methodology of Refs. [25, 89], teraction potential in a periodic one-dimensional (1D)
Gritsevet al.[78]investigatedaquenchfromγ =0 geometry of length L. We work in units such that (cid:126)=1
→∞
bycombiningalgebraicBetheansatzexpressionsforform and the particle mass m = 1/2, and so the Hamiltonian
factors with truncated sums over states, and employing of this system reads
Monte Carlo summation over the eigenstate components
of the initial state. N ∂2 N
Hˆ = +2c δ(x x ), (1)
In this paper we take a different approach, and calcu- − ∂x2 i− j
i=1 i i<j
late correlation functions of the Lieb-Linger model, both (cid:88) (cid:88)
inandoutofequilibrium,bycalculatingmatrixelements wherecistheinteractionstrength. ThecoordinateBethe
between Lieb–Liniger eigenstates directly within the co- ansatz yields eigenstates λ of Hamiltonian (1) with
j
ordinate Bethe ansatz formalism. Given the known ex- spatial representation [17]|{ }(cid:105)
pressions for the coordinate-space forms of Lieb–Liniger
eigenstates, we generate symbolic expressions for matrix ζ ( x ) x λ
elements of operators between these states in terms of {λj} { i} ≡(cid:104){ i}|{ j}(cid:105)
N
the Bethe rapidities. The numerically obtained values of
= A exp i x λ
therapiditiescanthenbesubstitutedtoyieldessentially {λj} m σ(m)
numerically exact values for the matrix elements. (cid:88)σ (cid:104) m(cid:88)=1 (cid:105)
icsgn(x x )
k l
In our previous work we applied this methodology to 1 − , (2)
× − λ λ
quenches from the ideal gas ground state to positive γ k>l(cid:18) σ(k)− σ(l) (cid:19)
(cid:89)
for up to N =5 particles [90]. In Sec. II we provide the
where the rapidities λ (or quasimomenta) are solutions
details of the methodology, and describe how it can be j
of the Bethe equations
usedtocalculatethematrixelementsoftheLieb–Liniger
eigenstates. These symbolic expressions, and thus the
N
2π 2 λ λ
computational cost of evaluating them, grow combina- λ = m arctan j − k . (3)
j j
torially with particle number, restricting the method to L − L c
k=1 (cid:18) (cid:19)
(cid:88)
3
The quantum numbers m are any N distinct integers ψ(t) = C (t) λ we have
(half-integers) in the casejthat N is odd (even) [2], and | (cid:105) {λj} {λj} |{ j}(cid:105)
(cid:80)
ofσ1d,e2n,o..te.s,Na s.umThoevneroramllaNliz!aptieornmcuotnatsitoannstσre=ad{sσ[1(j7)]} (cid:104)Oˆ(cid:105)t ≡(cid:104)ψ(t)|Oˆ|ψ(t)(cid:105)
(cid:80){ } = C (t)C (t) λ Oˆ λ , (7)
{∗λ(cid:48)j} {λj} (cid:104){ (cid:48)j}| |{ j}(cid:105)
A = k>l(λk−λl) , {(cid:88)λj}{(cid:88)λ(cid:48)j}
{λj} N!det{M{λj}(cid:81)} k>l[(λk−λl)2+c2] 1/2(4) wρˆher=eas in aρSstEatistλical enλsem,bwlee fiwnitdh density matrix
where M (cid:2) is the N N m(cid:81)atrix with elements(cid:3) SE {λj} {λj}|{ j}(cid:105)(cid:104){ j}|
{λj} × (cid:80)
Oˆ = ρSE λ Oˆ λ . (8)
N 2c (cid:104) (cid:105) {(cid:88)λj} {λj}(cid:104){ j}| |{ j}(cid:105)
[M ] =δ L+
{λj} kl kl c2+(λk λm)2 In this article, we focus in particular on the normalized
(cid:16) m(cid:88)=1 − (cid:17) mth-order equal-time correlation functions
2c
. (5)
− c2+(λk−λl)2 g(m)(x1,...,xm,x(cid:48)1,...,x(cid:48)m;t)
Ψˆ (x ) Ψˆ (x )Ψˆ(x ) Ψˆ(x )
The rapidities determine the total momentum P = † 1 ··· † m (cid:48)1 ··· (cid:48)m
, (9)
Nj=1λj and energy E = Nj=1λ2j of the system in each ≡ [(cid:104)nˆ(cid:68)(x1)(cid:105)···(cid:104)nˆ(xm)(cid:105)(cid:104)nˆ(x(cid:48)1)(cid:105)···(cid:104)nˆ(x(cid:48)m)(cid:105)(cid:69)]1/2
eigenstate. The ground state of the system corresponds
t(cid:80)o the set of N rapidities(cid:80)that minimize E and consti- where Ψˆ( )(x) is the annihilation (creation) operator for
†
tutethe(pseudo-)Fermiseaofthe1DBosegas[17]. The the Bose field and nˆ(x) Ψˆ (x)Ψˆ(x). Here and in the
†
Fermi momentum following we drop the tim≡e index t of the state vectors.
Since the Hamiltonian we consider in this article is
2πN 1 translationally invariant along the periodic volume of
k = − (6)
F
L 2 lengthL,themeandensity nˆ(x) nisconstantinboth
(cid:104) (cid:105)≡
time and space, and g(m)(x ,...,x ,x ,...,x ;t) =
1 m (cid:48)1 (cid:48)m
is the magnitude of the largest rapidity occurring in the Ψˆ (x ) Ψˆ (x )Ψˆ(x ) Ψˆ(x ) /nm. The cor-
ground state in the Tonks–Girardeau limit of strong in- (cid:104)rela†tio1n··· f†uncmtions (cid:48)1 ··g·(m)(x(cid:48)m,(cid:105)...,x ,x ,...,x ;t)
1 m (cid:48)1 (cid:48)m
teractions [3]. The only parameter of the Lieb–Liniger can therefore be expressed as the expectation val-
model in the thermodynamic limit is the dimensionless ues of the operators gˆ(m)(x ,...,x ,x ,...,x )
cfionuitpelisnygstγem≡s,cp/hny,swichaelrqeunan≡titNie/sLalsisotdheepe1nDddonentshiteyp.aIrn- Ψˆ†(x1)···Ψˆ†(xm)Ψˆ(x(cid:48)1)···Ψˆ(x1(cid:48)m)/nm.m W(cid:48)1e note(cid:48)m tha≡t
for the same reasons as above the matrix elements
ticle number N (see, e.g., Sec. IIIC), whereas the length λ gˆ(m)(x ,...,x ,x ,...,x ) λ are invariant
Lofoursystem, andthereforealsothedensityn, arear- (cid:104){ (cid:48)j}| 1 m (cid:48)1 (cid:48)m |{ j}(cid:105)
under global coordinate shifts x x + d and thus,
bitrary. Consequently,inthisarticlewewillspecifyboth →
without loss of generality, we can set one of the spatial
N andγ. Unlessspecifiedotherwise, wemeasuretimein
variablestozero. Forthefirst-ordercorrelationfunction,
units of kF−2, energy in units of kF2, and length in units the matrix elements are
of kF−1.
λ gˆ(1)(0,x) λ λ Ψˆ (0)Ψˆ(x) λ
(cid:104){ (cid:48)j}| |{ j}(cid:105)≡(cid:104){ (cid:48)j}| † |{ j}(cid:105)
N
= dx dx ζ (0,x ,...,x )
B. Calculation of correlation functions and overlaps n (cid:90) 1··· N−1 {∗λ(cid:48)j} 1 N−1
ζ (x,x ,...,x ). (10)
× {λj} 1 N−1
As the eigenstates λj form a complete basis [91] The evaluation of the integral in Eq. (10) is complicated
|{ }(cid:105)
for the state space of the Lieb–Liniger model, the ex- bythesignfunctioninEq.(2)andtheassociatednonan-
pereacttoartioOˆn ivnalaueS(cid:104)cOhˆr(cid:105)¨otd=ingTerr-{pρˆi(ctt)uOˆre} doefnasnityarmbiattrrairxyρˆo(pt)- anlaytteiscixtiesanindζx{λjc}o(i{nxciid}e).wHhoewreevaenry, wtweocapnarutsicelethceooBrodsie-
k l
can be expressed as a sum of matrix elements of Oˆ be- symmetry of the wave function ζ ( x ) to reexpress
tween the states λ . In particular, in a pure state this matrix element as a sum of in{λtej}gr{alsi}
j
|{ }(cid:105)
N 1
N! −
λ gˆ(1)(0,x) λ = dx dx ζ (0,x ,...,x )ζ (x ,...,x ,x,x ,...,x ),
(cid:104){ (cid:48)j}| |{ j}(cid:105) n (cid:88)(cid:96)=0 (cid:90)RN−1,(cid:96)(x) 1··· N−1 {∗λ(cid:48)j} 1 N−1 {λj} 1 (cid:96) (cid:96)+1 N−1
(11)
4
over the ordered domains [10]
(x): 0 x < <x <x<x < <x L. (12)
M,j 1 j j+1 M
R ≤ ··· ··· ≤
Substituting the coordinate-space form [Eq. (2)] of the Lieb-Liniger eigenfunctions, we obtain
N! ic ic
λ gˆ(1)(0,x) λ = A A 1 1+
(cid:104){ (cid:48)j}| |{ j}(cid:105) n {λj} ∗{λ(cid:48)j}(cid:88)σ (cid:88)σ(cid:48) j(cid:89)>k(cid:16) − λσ(j)−λσ(k)(cid:17)j(cid:89)(cid:48)>k(cid:48)(cid:16) λ(cid:48)σ(cid:48)(j(cid:48))−λ(cid:48)σ(cid:48)(k(cid:48))(cid:17)
N 1 N 1
− −
exp(iλ x) dx dx exp i (λ λ )x , (13)
× (cid:88)(cid:96)=0 σ((cid:96)+1) (cid:90)RN−1,(cid:96)(x) 1··· N−1 (cid:32) m(cid:88)=1 σ((cid:96)+1)(m)− (cid:48)σ(cid:48)(m+1) m(cid:33)
where σ((cid:96)+1) = (σ(1),...,σ((cid:96)),σ((cid:96)+2),...,σ(N)). The matrix elements of the second-order correlation operator
gˆ(2)(0,x) Ψˆ (0)Ψˆ (x)Ψˆ(x)Ψˆ(0)/n2 are similarly given by
† †
≡
N! ic ic
λ gˆ(2)(0,x) λ = A A 1 1+
(cid:104){ (cid:48)j}| |{ j}(cid:105) n2 {λj} ∗{λ(cid:48)j}(cid:88)σ (cid:88)σ(cid:48) j(cid:89)>k(cid:16) − λσ(j)−λσ(k)(cid:17)j(cid:89)(cid:48)>k(cid:48)(cid:16) λ(cid:48)σ(cid:48)(j(cid:48))−λ(cid:48)σ(cid:48)(k(cid:48))(cid:17)
N 2 N 2
− −
exp i(λ λ )x dx dx exp i (λ λ )x , (14)
× (cid:88)(cid:96)=0 (cid:16) σ((cid:96)+2)− (cid:48)σ(cid:48)((cid:96)+2) (cid:17)(cid:90)RN−2,(cid:96)(x) 1··· N−2 (cid:32) m(cid:88)=1 σ(1,(cid:96)+2)(m)− (cid:48)σ(cid:48)(1,(cid:96)+2)(m) m(cid:33)
whereσ(1,(cid:96)+2) =(σ(2),...,σ((cid:96)+1),σ((cid:96)+3),...,σ(N))andσ(1,(cid:96)+2) isdefinedanalogouslyintermsoftheelementsof
(cid:48)
σ . In the limit x 0 this expression simplifies somewhat, and in general the matrix elements of the local mth-order
(cid:48)
correlation operat→or gˆ(m)(0) [Ψˆ (0)]m[Ψˆ(0)]m/nm are given by the expression
†
≡
N! ic ic
λ gˆ(m)(0) λ = A A 1 1+
(cid:104){ (cid:48)j}| |{ j}(cid:105) nm {λj} ∗{λ(cid:48)j}(cid:88)σ (cid:88)σ(cid:48) j(cid:89)>k(cid:16) − λσ(j)−λσ(k)(cid:17)j(cid:89)(cid:48)>k(cid:48)(cid:16) λ(cid:48)σ(cid:48)(j(cid:48))−λ(cid:48)σ(cid:48)(k(cid:48))(cid:17)
N m N m
− −
dx dx exp i (λ λ )x , (15)
× 1··· N−m (cid:32) σ(m+n)− (cid:48)σ(cid:48)(m+n) n(cid:33)
(cid:88)(cid:96)=0 (cid:90)RN−m n(cid:88)=1
where the domain :0 x <x < <x L. We note, moreover, that Eqs. (13)–(15) include as degenerate
M 1 2 M
R ≤ ··· ≤
cases the diagonal matrix elements (cf. Ref. [10]) appropriate to the calculation of correlations in the ground state
(Sec. III) and in statistical ensembles (Sec. IV).
The calculation of correlation functions from Eqs. (13)–(15) involves the evaluation of integrals of the general form
M L xM x(cid:96)+2
dx dx exp i κ x = dx eiκMxM dx eiκM 1xM 1 dx eiκ(cid:96)+1x(cid:96)+1
(cid:90)RM,(cid:96)(x) 1··· M (cid:32) m(cid:88)=1 m m(cid:33) (cid:90)x M (cid:90)x M−1 − − ···(cid:90)x (cid:96)+1
x x(cid:96) x2
dx eiκ(cid:96)x(cid:96) dx eiκ(cid:96) 1x(cid:96) 1 dx eiκ1x1, (16)
×(cid:90)0 (cid:96) (cid:90)0 (cid:96)−1 − − ···(cid:90)0 1
where (for the repulsive interactions c>0 considered in x ) arising from the primitive integral [92]
m+1
this article) the κ are real numbers. A single closed
m
form for this integral does not exist, as in general one dxxpeikx = (i/k)p+1Γ(p+1, ikx)
or more κ may vanish, and this must be handled sep- − −
m (cid:90)
aedragteeloyfftrhome ptahreticcauslaerofseκtms o(cid:54)=f 0ra.pHidoiwtieevser,λgiveanndknoλwl- = p!(i/k)p+1eikx p (−ikx)s, (17)
{ j} { (cid:48)j} − s!
(and permutations σ and σ(cid:48)), and thus of the locations (cid:88)s=0
ofzeroexponentsκ =0inEq.(16), eachindividualin-
m
tegralofthisformcanbereducedtoanalgebraicexpres- in the case that κm is nonzero, or from dxxp other-
sion in terms of κ . More specifically, each successive wise. In our calculations, the construction of algebraic
integration dx{ myi}elds a term (involving, in general, expressions for the integrals occurring in E(cid:82)qs. (13)–(15)
m
in terms of the rapidities λ is efficiently performed by
j
(cid:82) a simple computer algorithm that accounts for and com-
bines the symbolic terms that arise from these succes-
5
sive reductions. We note that, e.g., each matrix element gressing to smaller values of γ using initial guesses given
λ gˆ(1)(0,x) λ isasumofN integralsover(N 1)- bylinearextrapolationofthesolutionsatstrongerinter-
(cid:104){ (cid:48)j}| |{ j}(cid:105) −
dimensionaldomainsandthattheintegrandineachcase action strengths.
comprises (N!)2 terms [10], illustrating the dramatically We note that this algorithmic approach also provides
increasing computational cost of evaluating correlation for the efficient and accurate calculation of the overlaps
functions with increasing N. Nevertheless, the explicit λ µ between eigenstates of Hamiltonian (1) cor-
j j
(cid:104){ }|{ }(cid:105)
closed-form expression for the integral produced by our respondingtodifferentvaluesofγ,whichwemakeuseof
algorithm can be evaluated to obtain a numerically ex- inouranalysisofnonequilibriumdynamicsinSec.IV.In
act result by substituting in the values of the rapidities. particular, the overlap between an arbitrary eigenstate
ThelatterareobtainedbysolvingEq.(3)numericallyus- λ of Hˆ at a finite interaction strength γ > 0 and
j
|{ }(cid:105)
ing Newton’s method, starting in the Tonks–Girardeau the noninteracting ground state 0 , with constant spa-
| (cid:105)
regime of strong interactions γ 1 and iteratively pro- tial representation x 0 =L N/2, is simply given by
i −
(cid:29) (cid:104){ }| (cid:105)
N
N! ic
λ 0 = A 1+ dx dx exp i λ x , (18)
(cid:104){ j}| (cid:105) LN/2 {λj}(cid:88)σ j(cid:89)>k(cid:16) λσ(j)−λσ(k)(cid:17)(cid:90)RN 1··· N (cid:32)− n(cid:88)=1 σ(n) n(cid:33)
whichcaneasilybeevaluatedsemi-analyticallyusingour herence with increasing γ [16]. As is well known, true
algorithm. In practice we find that the results we ob- long-range order, lim g(1)(x) = n > 0 [98, 99], is
x 0
tain for the overlaps from our evaluation of Eq. (18) prohibited in an intera→c∞ting homogeneous 1D Bose gas
agree with the recently derived closed-form expressions in the thermodynamic limit, even at zero temperature
for these quantities [84, 93–95], which imply in particu- (see Ref. [97] and references therein). Indeed the Lieb–
lar that λ 0 1/λ2 as any λ . Liniger system is quantum critical at zero temperature,
(cid:104){ j}| (cid:105)∝ j j →∞
and the asymptotic long-range behavior of g(1)(x) is a
power-law decay (so-called quasi-long-range order) [17].
III. GROUND-STATE CORRELATION Thispower-lawscalingofg(1)(x)isonlyexpectedtobe
FUNCTIONS realized at separations x large compared to the healing
length ξ =1/√γ and, in a finite periodic geometry such
As a first application of our methodology we calcu- as we consider here, is curtailed by the finite extent L of
late the correlation functions of the Lieb–Liniger model the system (see, e.g., Ref. [16]). Indeed, for γ =0.1, the
in the ground state for up to N = 7 particles. In this power-law decay is not visible in our finite-sized calcu-
case, we need to evaluate only the diagonal elements of lation, although as the interaction strength γ increases
Eqs.(13)–(15)intheground-statewavefunction,thereby g(1)(x) exhibits behavior consistent with power-law de-
obtainingexactalgebraicexpressionsforcorrelationfunc- cay over an increasingly large range of x, see Fig. 1(a).
tions in terms of the ground-state rapidities, which are In particular, for γ (cid:38) 10, our results for g(1)(x) seem
themselves determined to machine precision (Sec. IIB). to converge toward the asymptotic scaling of the Tonks–
The ground-state correlations of the Lieb–Liniger model Girardeau limit (black dot-dashed line) with increasing
have been considered extensively in previous works (see γ.
Refs.[96,97]andreferencestherein),andwecompareour Due to the translational invariance of our system, the
exact mesoscopic results to those obtained with various first-order correlations of the Lieb–Liniger ground state
othermethodsandapproximations,forfinitesystemsizes are encoded in the momentum distribution
as well as in the thermodynamic limit. This allows us to L
clarify the utility and limitations of calculations, such as n(kj)=n dxe−ikjxg(1)(x), (19)
ourshereandinRef.[90],thatinvolveonlysmallparticle (cid:90)0
numbers. which, inour(cid:101)finite periodicgeometry, isonlydefinedfor
discrete momenta k = 2πj/L, with j an integer. In
j
Fig.1(b)weplotthemomentumdistributionsn(k )cor-
j
A. First-order correlations respondingtothefirst-ordercorrelationfunctionsg(1)(x)
shown in Fig. 1(a). The first feature that w(cid:101)e note in
We begin by considering the first-order correlation Fig. 1(b) is that for all interaction strengths, n(k) ex-
function g(1)(x) g(1)(0,x) in the ground state of the hibits a power-law decay n(k) k 4 (dot-dashed black
−
≡ ∝
Lieb–Liniger model. In Fig. 1(a) we plot g(1)(x) for line) at high momenta. This is a universal result for
(cid:101)
N = 7 particles for a range of interaction strengths γ, delta-functioninteractionsinonedimension[21,89,100]
(cid:101)
whichexhibitstheexpecteddecreaseinspatialphaseco- (and indeed also in higher dimensions [101]). The effects
6
(a) 1 1 (c) coherence g(2) (x)=1 N 1 (horizontal dashed line) is
γ=0 − −
therefore independent of x. As the interaction strength
0.8
γ =0 γ is increased, the second-order coherence is increasingly
1)(x) 0.5 γγγ===110.01 γγγ ===110.01 00..46(2)g(x) sinugplpyreesnshedanacetdzeartosesppaatriaatliosnespaxra(cid:38)tio2nkF−an1.d cOosrcrielslaptoionnds-
(g γ=100 γ =100 0.2 in g(2)(x) develop at finite x as the system enters the
0.3 γ=∞ x−21 γ =∞ 0 strongly interacting regime γ (cid:29) 1 [9, 17] and, in par-
ticular, for γ = 100 (dashed cyan line), our numerical
0.2 0.5 1 2 5 0 2 4 6 8
x(unitsof kF−1) x(unitsof kF−1) rTeosnukltss–Garirearpdreaacuticlaimllyitirnedsuislttin(sgouliisdhabblaleckfrlionme)t[h3e]. exact
(b) 1 (d)
An alternative representation of the second-order cor-
10 1
− 0.8 relations of the ground state is given by the static struc-
γ=0 ture factor S(k), which is related to g(2)(x) by [11]
)10−4 γ=0.1 0.6 k)
k γ=1 (
n(10−e7 γ=10 0.4 S S(kj)=1+n Ldxe−ikjx g(2)(x) 1 . (20)
γ=100 0.2 −
k−4 γ= (cid:90)0 (cid:104) (cid:105)
10−10 ∞ 0 In Fig. 1(d) we present the structure factors S(k) cor-
100 101 0 1 2 3 4 5 responding to the correlation functions g(2)(x) shown in
k(unitsof kF) k(unitsof kF) Fig. 1(c). For all values of γ, S(0) = 0 due to particle-
numberconservationandtranslationalinvariance. Inthe
FIG. 1. (Color online) One- and two-body correlations
ideal-gas limit (red circles) S(k ) =1 for all nonzero k .
j j
in the Lieb–Liniger ground state, for N = 7 particles. (a)
In the opposite limit of a Tonks–Girardeau gas
Non-localfirst-ordercoherenceg(1)(x). Theblackdot-dashed
|(lixbn|)e−C1i/no2droricefasapteTosnodnthiknesg–aGzseyirrmaor-ptdteeomatuipcegrlaoasntiugnr-ertahmnegotemhebernemhtuoamdviyondraisgmt(r1iic)b(luximt)ioi∝tn. Sγ=∞(kj)=(cid:40) |kj|(12−k1FN−1) |kkjj|>≤22kkFF, (21)
| |
n(cid:101)(kj). Theblackdot-dashedlineindicatestheuniversalhigh-
momentumpower-lawscalingn(cid:101)(k) k−4 commontoallpos- which tends, in the thermodynamic limit, to the well-
∝
itive interaction strengths [21]. (c) Non-local second-order known result (see, e.g., Ref. [9]) S(k) = k /2k for
F
| |
coherence g(2)(x). (d) Corresponding static structure factor k 2k , and S(k) = 1 for k > 2k . Just as for
F F
| | ≤ | |
S(k). g(2)(x), we observe that for γ = 100 (cyan plus sym-
bols), our numerical results for S(k) are almost iden-
tical to the known exact expression [Eq. (21)] for the
ofthefiniteextentLofthesystemonthefirst-ordercor- Tonks–Girardeau limit (black crosses). For smaller val-
relations are again evident in this momentum-space rep- ues of γ our mesoscopic results for S(k) appear consis-
resentation. For γ = 0.1, no deviation from the k−4 tentwiththoseofRefs.[22,25],obtainedusingquantum
∝
scaling is observed for the smallest (nonzero) momenta Monte Carlo and algebraic-Bethe ansatz techniques, re-
kj that can be resolved in the periodic geometry. For spectively.
larger values of the interaction strength, n(k) departs We now focus in more detail on local correlation func-
from the k−4 scaling at increasingly large values of k tions. We note that the local second-order coherence
∝
withincreasingγ,anddevelopsahumpatmomentanear has recently been proposed as a measure of quantum
(cid:101)
kF for γ (cid:38) 10 [89]. We note that although the small-k criticality in the 1D boson system [102], while the local
behaviorofn(k)tendstowardsthe k−1/2scalingexhib- third-order correlations have received increasing atten-
∝
ited by the Tonks–Girardeau gas in the thermodynamic tion both in theory [103] and experiment [47, 104–106].
limit, the rounding off of the power-law decay of g(1)(x) The local fourth-order correlations for the Lieb–Liniger
(cid:101)
as x L/2 precludes n(k) from reaching the known modelhavealsobeeninvestigated[107]. InFig.2,weplot
→
asymptotic k 0 behavior in our finite geometry. the local second-order coherence g(2)(0) (solid red line),
→
together with the local third-order coherence g(3)(0) =
(cid:101)
[Ψˆ (0)]3[Ψˆ(0)]3 /n3 (dotted green line), and the local
†
(cid:104) (cid:105)
B. Second-, third-, and fourth-order correlations fourth-order coherence g(4)(0) = [Ψˆ (0)]4[Ψˆ(0)]4 /n4
†
(cid:104) (cid:105)
(dashed blue line) for N =7 particles and a broad range
In Fig. 1(c), we present the nonlocal second-order co- of interaction strengths γ. For comparison, we also plot
herence g(2)(x) g(2)(0,x,x,0), which provides a mea- the asymptotic results obtained in the Bogoliubov limit
≡
sure of density-density correlations, for N = 7 parti- of weak interactions (γ 0) in the thermodynamic
→
cles at a range of interaction strengths γ. In the lim- limit [12, 18] (left-hand dot-dashed lines). The numer-
iting case of an ideal gas (γ = 0), the ground state of ical results for small γ are broadly comparable to these
the system is a Fock state of N particles in the zero- thermodynamic-limitresults. However,forthesmallpar-
momentum single-particle mode, and the second-order ticlenumbersconsideredhere,thesuppressionofg(2)(0),
7
100 (a) 1 1 (c)
0.8
10 4
(n)(0) 10−−8 g(2)(0)×101 (1)g(x) NNN ===345 NNN ===345 00..46(2)g(x)
g 10 12 g(3)(0) N =6 N =6 0.2
− g(4)(0) 10−1 0.5 N =7x−1/2K N =7 0
10−16 × 10−1 100 1010 2 4 6 8 10
10−3 10−1 101 103 x(unitsof (πn)−1) x(unitsof (πn)−1)
γ
(b) 100 1 (d)
k 1+1/2K
− 0.8
FIG. 2. (Color online) Interaction-strength dependence of
10 2
tLhieebl–oLcianligseercognrdo-u,ntdhisrtda-tea,nfdorfoNur=th-7orpdaerrticcolehse.reTnoceaiidn vtihse- n(k) −e N =3 NN ==34 0.6 (k)
ibility, we plot g(2)(0) scaled by a factor of 101, and g(4)(0) 10−4 N =4 N =5 0.4 S
stocatilcedwbeayka- (faγctor1o)fa1n0d−1st.rDonogt--cdoausphleidnglin(γes in1d)iceaxtepraessysimonps- NN ==56 NN ==67 0.2
for g(2)(0), g(3(cid:28))(0) and g(4)(0) in the therm(cid:29)odynamic limit 10−6 N =7 0
(see text). 100 101 0 1 2 3 4 5
k(unitsof πn) k(unitsof πn)
g(3)(0), and g(4)(0) due to interactions in the limit of FIG. 3. (Color online) Dependence of first- and second-
ordercorrelationsintheLieb–Linigergroundstateonparticle
small γ is overshadowed by the suppression due to the
number N for γ = 10. (a) First-order correlation function
finite population of the system [20]. At larger γ, the ef-
g(1)(x). (b) Corresponding momentum distribution function
fects of interactions dominate, and the numerical results
n(cid:101)(kj). Black dot-dashed lines in (a) and (b) indicate the
converge closely to the appropriate strong-coupling ex-
asymptotic infrared scaling of g(1)(x) and n(k), respectively,
(cid:101)
pressions [18] (right-hand dot-dashed lines). We note,
with Luttinger parameter K = 1.40 (see text). (c) Second-
therefore, that the local correlations of the Lieb–Liniger order correlation function g(2)(x). (d) Corresponding static
ground state, and particularly their scaling with γ, ap-
structure factor S(k). The black dot-dashed lines in (c) and
peartobequiteinsensitivetotheinfraredcutoffimposed (d) represent the phenomenological expressions of Ref. [108]
by the finite extent of our system in the strongly inter- forg(2)(x)andS(k)inthethermodynamiclimit,respectively.
acting regime γ 1.
(cid:29)
of πn [(πn) 1], so as to avoid a potentially misleading
−
C. System-size dependence dependence of the unit of length on N [cf. Eq. (6)].
In Fig. 3(a) we plot g(1)(x) for particle numbers N =
Theresultswehaveobtainedsofarindicatethat,asex- 3,4,5,6, and 7. For small x, the curves fall nearly per-
pected,thesmallsizeofoursystemleadstocorrectionsto fectlyononeline. Thesamebehaviorcanbeobservedfor
correlationfunctionsascomparedtotheirknownasymp- thelarge-ktailofthecorrespondingmomentumdistribu-
totic forms in the thermodynamic limit. However, our tion n(k), which we plot in Fig. 3(b). Indeed, at larger
results also suggest that the effects of finite system size momenta k (cid:38) 2πn, n(k) appears to exhibit a rapid col-
are comparatively less important for local correlations, lapse to a single curve with increasing N [21, 109]. How-
(cid:101)
particularly in the limit of large interaction strengths ever, the differences in n(k) are so small that they can
(cid:101)
γ 1. To further elucidate the potential significance not be seen in Fig. 3(b). For small momenta, our choice
(cid:29)
of finite-size effects in our calculations of nonequilibrium of units implies an increasing resolution with increasing
(cid:101)
dynamics[90],herewegiveabriefcharacterizationofthe particlenumber,specificallyk1 =2π/L (πn)−1 =2/N.
×
dependence of correlation functions of the Lieb–Liniger However, this lowest resolvable momentum seems to fall
ground state on the particle number N at a fixed value ononelineforincreasingparticlenumber,indicatingthat
of the interaction strength γ. the infrared behavior of large systems can be at least
Specifically we consider the case for γ = 10, as this partly accessed by our mesoscopic system sizes.
valueplacesthesysteminthestronglyinteractingregime Luttinger-liquid theory predicts a long-range power-
γ 1 (which appears less sensitive to finite-size effects law decay g(1)(x) x 1/2K, where the Luttinger pa-
−
th(cid:29)antheweaklyinteractingregimeγ (cid:46)1), whilestillex- rameter K can be∝ca|lc|ulated from the thermodynamic
hibitingsignificantdeviationsfromtheTonks–Girardeau limitoftheBetheansatzsolution(see,e.g.,Refs.[16,17]
limit (see, e.g., Ref. [9]). Whereas elsewhere in this pa- andreferencestherein). ForourparameterswehaveK =
per we quote momenta (lengths) in units of kF (kF−1), 1.40, implying an asymptotic scaling g(1)(x) ∝ |x|−0.357
in comparing results between systems with different par- [black dot-dashed line in Fig. 3(a)]. This corresponds to
ticle numbers N we quote momenta (lengths) in units a power-law behavior n(k) k 1+1/2K = k 0.643 [16]
− −
∝| | | |
(cid:101)
8
[dot-dashed line in Fig. 3(b)] for small momenta. We interaction-strength quenches, from different initial in-
note that this infrared scaling is a true many-body ef- teraction strengths γ , to a common final value of the
0
fect and as such does not show up for N = 2 particles. coupling γ. Specifically, we consider a quench from the
Indeed, one can show analytically that, for N = 2, the non-interacting limit γ = 0 (similar to those previously
0
momentumdistributionn(k) (λ2 k2) 2 andthusk 4 studied in Refs. [63, 78, 84, 90, 112–115]) and a quench
∝ 1− − −
is the highest power in the series expansion of n(k). from the correlated ground state obtained for a strong
In Fig. 3(c) we plot the nonlocal second-order coher- interaction strength γ =100. As Hˆ(γ) is time indepen-
(cid:101) 0
ence g(2)(x) for γ = 10 and N = 3,4,5,6, and 7. The dentfollowingthequench,energyisconservedduringthe
(cid:101)
corresponding static structure factor S(k) is shown in dynamics. We choose the final interaction strength after
Fig. 3(d). In Fig. 3(d) we also plot (black dot-dashed the two quenches such that the postquench energy is the
line)theformofS(k)resultingfromthephenomenologi- same in both cases.
calexpressionproposedinRef.[108](seealsoRef.[110]). The statistical description of the dynamics of suffi-
Thisexpressioninvolvesthelimitingdispersionsandedge ciently ergodic systems is usually based on the assump-
exponents of the Lieb–Liniger model, which we obtain tion that the energy is the sole integral of motion, such
by numerically solving the appropriate integral equa- that the equilibrium system is entirely determined by its
tions [1, 111]. We also plot the corresponding predic- energy. If this would be the case for our system, the
tion for g(2)(x) (black dot-dashed line) in Fig. 3(c). We two quenches would lead to the same equilibrium state.
note that the numerical results for our mesoscopic sys- However,thedynamicsaccordingtotheintegrableLieb–
temsare,ingeneral,ratherclosetothephenomenological LinigerHamiltonianarestronglyconstrainedbythecon-
thermodynamic-limit expressions even for the relatively served quantities other than the total energy. By per-
small particle numbers considered here. formingtwodifferentquenchestothesamefinalHamilto-
nianandenergy,weinvestigatetheeffectsofintegrability
on the postquench evolution of the Lieb–Liniger system.
IV. APPLICATION TO NONEQUILIBRIUM The conserved energy following the quench is the en-
DYNAMICS ergy of the system at time t=0+,
E ψ(0+)Hˆ(γ)ψ(0+)
We now apply our methodology to the nonequilibrium γ0→γ ≡(cid:104) | | (cid:105)
dynamicsoftheLieb–Linigermodel. Specifically,wecon- dEG(γ)
=E (γ )+(γ γ ) (24)
sider the evolution of a system, initially prepared in the G 0 − 0 dγ γ0
groundstateofHamiltonian(1)withinteractionstrength (cid:12)
(cid:12)
γ , following an abrupt change, at time t = 0, of the whereEG(γ0)istheenergyofthegroundstate(cid:12) ψ0 ofthe
0 | (cid:105)
interaction strength to a distinct value γ = γ0 — a so- initial Hamiltonian Hˆ(γ0) and we used the well-known
called“interactionquench”. Theevolution(cid:54)ofthesystem result g(2)(0) = n 2N 1dE (γ)/dγ [18], which implies
γ − − G
following such a quench is generated by Hamiltonian (1) thatE isgivenbyfollowingthetangenttothecurve
with interaction strength γ, which we denote by Hˆ(γ) E (γ)γa0→t γγ out to γ. Here, g(2)(0) ψ gˆ(2)(0)ψ
hereafter. The time-evolving state is given at all times isGthe local0second-order cohereγn0ce in≡th(cid:104)e i0n|itial st|at0e(cid:105).
t>0 by In the case of a quench from the noninteracting ground
state (γ =0), Eq. (24) reduces to the simple expression
0
|ψ(t)(cid:105)={(cid:88)λj} C{λj}e−iE{λj}t|{λj}(cid:105), (22) Eim0p→aγrt=ed(Nto−th1e)nsy2sγte[m66,d9u0r]i,ngimtphleyiqnugetnhcahtdtihveeregneesrgays
γ [63]. By contrast, in a quench from the Tonks–
where λ are the eigenstates of Hˆ(γ) with energies Gi→rar∞deaulimitγ toafiniteinteractionstrengthγ
|{ j}(cid:105) 0 →∞
E , and C λ ψ are the overlaps of the thefinalenergyisboundedfromabove,E E ( ),
{λλj} with th{eλjin}it≡ial(cid:104)s{tajt}e| ψ0(cid:105) . The expectation value by the ground-state energy of the Tonks–∞G→irγa≤rdeaGu g∞as.
j 0
|{ }(cid:105) | (cid:105)
of an arbitrary operator Oˆ in the state ψ(t) is given by Nevertheless, according to Eq. (24), a final interaction
| (cid:105) strength 0<γ <100 such that E =E does
Oˆ ψ(t)Oˆ ψ(t) exist. ∗ 100→γ∗ 0→γ∗
t
(cid:104) (cid:105) ≡(cid:104) | | (cid:105)
Here, we consider quenches of N = 5 particles, and
= C{∗λ(cid:48)j}C{λj}ei(E{λ(cid:48)j}−E{λj})t(cid:104){λ(cid:48)j}|Oˆ|{λj}(cid:105). determine this final interaction strength to machine pre-
{(cid:88)λj}{(cid:88)λ(cid:48)j} cision, inferring a value γ∗ =3.7660... from numerical
(23) solutions for the energy and local second-order coher-
ence of the ground state at finite γ (Sec. IIIB). We note
We use the methodology described in Sec. II to eval- that although the overlaps C of the initial state ψ
uaλte bOˆothλthetohvaetrlaapppseCar{λinj}Eaqn.d(2t3h)e. matrix elements wcailtlhytihnetheiegecnassteatoefstohfeHˆqu(γen∗)c{hλcja}frnombeγcal=cul0at[e9d3–a9n5a],l|yfto0ir(cid:105)-
(cid:104){ (cid:48)j}| |{ j}(cid:105) 0
One of the features of our methodology is that it the quench from γ = 100 no closed-form expressions
0
allows us to describe quenches between arbitrary in- for these quantities are known, and thus their numeri-
teraction strengths. In this paper we consider two cal values must be determined using the semi-analytical
9
100 noninteracting initial state (γ =0), as time evolves the
0
localsecond-ordercoherencedecaysfromitsinitialvalue
g(2)(0,t=0)=1 N 1 beforesettlingdowntofluctuate
−
t) 10−1 about the diagon−al-ensemble expectation value g(2)(0)
0, DE
(2)g( 10−2 γγ ==100→0→γ∗γ∗ (whiothrizroensutaltlsdoobtt-adiansehdedfolrinseim).ilTarhiqsubeenhcahveisoorfisthceoninstisetreanct-
DE,γ =0 γ tion strength from zero to a positive value in Ref. [90].
∗
→
DE,γ =100 γ For the quench from γ =100, the value of g(2)(0) in the
∗ 0
→
10 3 initial “fermionized” state is g(2)(0) 10 3. In this case
− −
≈
10 3 10 2 10 1 100 101 g(2)(0,t)risesastimeprogresses,andthenexhibitssome-
− − −
t(units of kF−2) whatirregularoscillationsaboutgD(2E)(0)(horizontalsolid
line). We observe that the decay (growth) of g(2)(0,t)
FIG. 4. (Color online) Time evolution of local second-order to its diagonal-ensemble value and the onset of irregular
correlations for N = 5 particles following quenches of the oscillations about this value occur on comparable time
interactionstrengthtoafinalvalueγ∗=3.7660... frominitial scales in the two quenches.
values γ = 0 (red dotted line) and γ = 100 (blue dashed
0 0 We note that the predictions of the diagonal ensem-
line). The horizontal solid (dot-dashed) line indicates the
predictionofthediagonalensembleforg(2)(0)forthequench ble for the local second-order coherence gD(2E)(0) are very
similar for the two quenches, despite the significant dif-
from γ =100 (γ =0).
0 0
ference between the values of g(2)(0) in the two initial
states. However, they are clearly distinct — g(2)(0) for
DE
methodology described in Sec. IIB. thequenchfromthenoninteractingstateisinfactlarger
Animportantsummaryofthepostquenchexpectation than that for the quench from the correlated state by an
value of an operator [Eq. (23)] is provided by the time- amount 0.0125, demonstrating that the system retains
≈
averaged value some memory of its initial state in the long time limit
as is expected for an integrable system. We analyze this
O = lim 1 τdt ψ(t)Oˆ ψ(t) . (25) difference in more detail in Sec. IVC.
τ τ (cid:104) | | (cid:105) Wenowturnourattentiontothetimeevolutionofthe
→∞ (cid:90)0
fullnon-localsecond-ordercorrelationfunctiong(2)(x,t).
Neglecting degeneracies in the spectrum of Hˆ(γ ) (see In Fig. 5(a) we show the dependence of g(2)(x,t) on sep-
∗
discussion in Appendix B), such averages are given by aration x for the quench from the noninteracting initial
the expectation values Oˆ = Tr ρˆ Oˆ of operators state at four representative times. [Note that the upper
DE DE
Oˆ in the density matrix(cid:104) (cid:105) { } limit x = 2πkF−1 of the x axis in Fig. 5(a) corresponds
to x = L/2 in the present case of N = 5 particles.]
ρˆ = C 2 λ λ (26) At t = 0 (horizontal solid line), the second-order coher-
DE | {λj}| |{ j}(cid:105)(cid:104){ j}| ence has the constant form of the noninteracting ground
{(cid:88)λj} state. At short times (e.g., t = 0.01kF−2, red dashed
of the diagonal ensemble [116, 117]. line) a minimum in g(2)(x) develops at zero separation,
Formally, the sums in Eq. (22),(23), and (26) range together with the corresponding maximum required by
over an infinite number of eigenstates λ , and thus the conservation of Ldx g(2)(x,t) [66]. As time pro-
|{ j}(cid:105) 0
the basis over which ψ(t) is expanded must be trun- gresses a wave pattern of maxima and minima develops
cated in our numerica|l cal(cid:105)culations. By only including and propagates away(cid:82)from the origin (e.g., t = 0.1kF−2,
elaigrgenersttahtaens swoimtheatnhraebshsoolludt,ewiencitoinals-issttaetnetloyvneerlgalpec|tCs{mλja}ll| lgirneee)n, tdhoettdeidstliinncet).mBaxyimtiamaentd=m1inkimF−2a(obfluge(2)d(oxt,-td)ahshaevde
contributions to correlation functions from weakly occu- broadened in such a way that they are no longer clearly
pied eigenstates and minimize the truncation error for a distinguishable and the correlation function agrees rea-
givenbasissize. Wequantifythistruncationerrorbythe sonably well with its diagonal-ensemble form (black dot-
vwieoldaitsicounsssoifntAheppneonrmdiaxliAza.tion and energy sum rules, as dFaigs.he5d(bl)inwe)efoshroswmatlhlesefpuallrastpiaocnes axn(cid:46)d 0ti.m25e×de2pπeknF−d1e.nIcne
of g(2)(x,t) following a quench from γ = 0, which gives
0
a more complete picture of the development of a corre-
A. Evolution of two-body correlations lation wave at short length scales and its propagation
to larger values of x as time progresses. The correlation
InFig.4weplotthetimeevolutionofthelocalsecond- waveweobservehereisconsistentwiththeresultsofpre-
order coherence g(2)(0,t) for N = 5 particles following viousinvestigationsofthedynamicsfollowingthesudden
quenches of the interaction strength from initial values introduction of repulsive interactions in an initially non-
γ =0 (red dotted line) and γ =100 (blue dashed line) interacting gas [63, 64, 66, 78, 118].
0 0
to the common final value γ . For the quench from the In Fig. 5(d) we plot the spatial form of g(2)(x,t) for
∗
10
FIG. 5. (Color online) Time evolution of the nonlocal second-order coherence function g(2)(x,t) following quenches of the
interaction strength to γ∗ from initial values [(a),(b),(c)] γ0 = 0 and [(d),(e),(f)] γ0 = 100. All data is for N = 5 particles.
[(a),(d)] Correlation function g(2)(x,t) at four representative times t. Black dot-dashed lines indicate the predictions of the
diagonal ensemble for the equilibrium form of this function. [(b),(e)] Evolution of coherence g(2)(x,t) and [(c),(f)] change in
cxo(ht)eroefntcheegfi(2r)s(txe,xt)tr−emg(u2m)(xo,ftt=he0)coforrreslahtoirotntwimaevse,tw≤h0ic.h5kyF−i2el.dBxlacktl0i.n5e16s±i0n.0(1c2)aannddx(f)int0d.4ic9a6±te0.p00o5wfeorr-laqwuefintcshetsoftrhoempγos0it=ion0
∝ ∝
and γ =100, respectively.
0
thequenchfromγ =100atthesamefourrepresentative with the “fermionic” depression around x = 0 lessening
0
times considered in Fig. 5(a). Despite the obvious dis- rather than growing in magnitude, a similar pattern of
tinctionthattheinitial(t=0,solidgreyline)correlation propagating correlation waves in g(2)(x,t) can again be
function is in the fermionized regime with g(2)(0) 1, seen.
(cid:28)
the behavior of g(2)(x,t) for this quench is qualitatively
The correlation-wave pattern common to both
similar to that observed for the quench from γ = 0, in
0 quenches is more clearly exhibited by the change
that at early times (e.g., t = 0.01kF−2, red dashed line), g(2)(x,t) g(2)(x,0) in the correlation function following
deviationsfromg(2)(x,t=0)occuronlyatsmallsepara- −
the quench, which we plot in Figs. 5(c) and 5(f). This
tionsx(cid:28)2πkF−1. Moreover,astimeevolvesandg(2)(0,t) representation of the postquench second-order coherence
increasestowardsg(2)(0),largermodulationsofg(2)(x,t) of the system reveals a remarkably similar pattern of
DE
aboutitsinitialfunctionalformdevelop(e.g.,t=0.1kF−2, propagating waves in both cases, although the maxima
green dotted line). At later times (e.g., t = 1kF−2, blue and minima of the two wave patterns are inverted rela-
dot-dashed line), g(2)(x,t) is close to g(2)(x) at small tive to one another. Fitting a power law to the position
sceuprsairoantsioanwsaxy f(cid:46)rom0.2i5t a×t2laπrkgF−er1,xb.uItneFxhiDgiE.bi5t(se)lawrgeepeloxt- xco(rt)reolafttihoenfiwrastvepsr,owpeagfiantdinxge∝xttr0e.5m16u±m0.0o1f2efaocrhthofetqhueetnwcho
from γ = 0 and x t0.496 0.005 for the quench from
the full space and time dependence of g(2)(x,t) follow- 0 ∝ ±
γ = 100, which we indicate by the solid black lines in
ing the quench from γ = 100. Although the behavior 0
0 Figs.5(c)and5(f). Thesepower-lawtrajectoriesarecon-
of g(2)(x,t) here obviously differs from that following a
sistentwiththe“telescoping”x t1/2 behaviorobtained
quench from the noninteracting initial state [Fig. 5(b)], ∝
for a quench γ = 0 in Ref. [63], and for quenches
→ ∞