Bragg spectroscopy of a strongly interacting Bose-Einstein condensate

We study Bragg spectroscopy of a strongly interacting Bose-Einstein condensate using time-dependent Hartree-Fock-Bogoliubov theory. We include approximatively the effect of the momentum dependent scattering amplitude which is shown to be the dominant factor in determining the spectrum for large momentum Bragg scattering. The condensation of the Bragg scattered atoms is shown to significantly alter the observed excitation spectrum by creating a novel pairing channel of mobile pairs.


Introduction
A strongly interacting Bose-Einstein condensate (BEC) has been rich topic to study both theoretically and experimentally. Theoretical problems arise from the need to provide a proper description of the elementary excitations and experimental difficulties from the instability of the BEC to three-body collisions. Research on these systems offers the potential to be highly rewarding by shedding light on other strongly interacting systems, such as superfluid 4 He.
Weakly interacting BECs, on the other hand, are generally well understood. The ground state properties are accurately described by Bogoliubov theory as the measurement of the excitation spectrum in 1999 using Bragg spectroscopy confirmed [1,2]. In Bragg spectroscopy, the condensate is excited using stimulated two-photon Bragg scattering yielding the possibility for very high momentum and energy resolution. The measured spectrum was in very good agreement with the Bogoliubov theory [3,4].
The Bragg spectrum of a strongly interacting BEC was measured in 2008 [5]. The study found significant deviations from the standard Bogoliubov theory, highlighting the need for a more thorough theoretical treatment of the process.
In this work, we study the Bragg specroscopy of a strongly interacting BEC using Hartree-Fock-Bogoliubov (HFB) theory [6]. The theory has been widely applied to various problems, both static and time-dependent. In spite of well known problems in the excitation spectrum, such as the unphysical zero-momentum gap, it has been able to provide correct qualitative features and good agreement with experiments [7,8,9,10]. The present time-dependent theory is an extension of the one used in [8], including both a proper renormalization and allowing for the macroscopic occupation of the Bragg scattered states but neglecting the explicit molecule formation channel. We will show how the appearance of condensate atoms in finite momentum states leads to profound effects in the observed Bragg spectrum.

Hamiltonian
The system is described by the Hamiltonian whereĤ Bragg (t) is the perturbation due to the Bragg field,ĉ k (ĉ † k ) is the annihilation (creation) operator for a boson with momentum k, ǫ k =h 2 k 2 2m , V is the volume, and U k,p is the two-body T-matrix of the atom-atom interaction.
As already observed in [5], the momentum dependence of the scattering amplitude is the dominating effect in the high momentum Bragg scattering and therefore we have included it in the results. A constant interaction strength U k,p = U would correspond to the delta function approximation of the real-space T-matrix (yielding constant U in the momentum space). However, we want to include the energy dependence of the two-body T-matrix in an approximative way. For low momenta k the energy (or the momentum) dependence of the two-body T-matrix for on-shell scatterings (k = p) is where a is the scattering length. In order to satisfy the time reversal symmetry of the scattering process (or the Hermiticity of the Hamiltonian (1)), the energy dependent interaction strength U k,p must be symmetric in the exchange of k and p. To keep the numerical solution tractable, we approximate the interaction strength by two constant interaction strengths U L := U 0,0 and U H := U QB,QB so that U k,p = U L in the low momentum manifold k, p ≪ Q B and U k,p = U H in the high momentum manifold (k and/or p ≈ Q B ). That is, the Bragg scattered condensed atoms will interact through interaction strength U H whereas the atoms in the unscattered condensate will feel U L . While the quantum fluctuations and the Bragg scattering produce excitations with momentum higher than Q B , the typical fraction of such atoms is at most 4% in our calculations. Thus, in order to speed up the calculation, we assume that all fluctuations interact through U L . The operatorĤ Bragg (t) in the rotating wave approximation can be written aŝ where Q B is the momentum of the Bragg field and Ω(t) is the coupling between the atoms and the Bragg field. Notice that the rotating wave approximation is used here to remove the counter rotating terms that violate the energy conservation. The energy of the Bragg field (energy difference between the photons in the two probing laser fields),hω B , is included in the coupling Ω(t). For a pulse of length T we have where Ω is the strength of the coupling and Θ(x) is the Heaviside function. Using a smooth Gaussian pulse instead of the rapid switching on/off would reduce ringing oscillations in many of the results below. However, the step function is more convenient and straight-forward to implement numerically. The Hartree-Fock-Bogoliubov theory for BECs has been widely used and its properties and problems are well known [6,11]. One important problem, especially related to spectroscopy, is the excitation gap at low momenta. This is in violation with the Hugenholtz-Pines theorem [12] and consequently the HFB theory does not yield a proper phonon spectrum. In practice this means that the Bragg spectrum will have some pathological features when the Bragg momentum Q B is small. In that region, we can hope to gain at most a qualitatively correct picture, assuming that the effect of the excitation gap is well understood. On the other hand, for large momenta, as is the case we consider here, the HFB theory is expected to work well.
Another problem in HFB theory relevant to this work is the ultraviolet divergence of the pairing field m 0 when the theory is not correctly renormalized. This is a consequence of the assumption of a contact interaction potential which is unable to properly describe high energy scattering processes. A correct result can still be obtained by employing a contact interaction approximation to the effective two-body T-matrix. However, since the calculation of the many-body T-matrix replicates the same set of diagrams, the bare two-body scattering diagrams need to be removed from the many-body picture. The regularization has been done in several different ways in the literature but effectively all remove the asymptote of the high-momentum scattering terms.
For a Bose-condensed system, the regularization procedure needs to be chosen carefully. Indeed, we have found out that the standard procedure used in several publications [13,14,8,9,10] leads into instabilities by making the energies of low momentum excitations imaginary. This regularization issue will be discussed in more detail below.

Equations of motion
Using the Hamiltonian (1), we derive the Heisenberg equation of motion For condensed states, k = l Q B , where l is an integer, this operator is allowed to have a nonvanishing expectation value ĉ n QB =: ψ n but otherwise the first non-zero terms include higher order correlators such as the normal Green's function ĉ † kĉ k and the anomalous Green's function ĉ kĉ− k . These correlators describe both thermal and quantum fluctuations and give corrections to the Gross-Pitaevskii equation for the condensate. The equation of motion (4) produces an infinite series of higher order correlators. We truncate this series at the two operator correlator level.
Due to the Bragg field, the number of relevant mean-fields is large and the resulting equation of motion is potentially involved. However, for sufficiently weak Bragg pulses, multiphoton scattering, in which the atom gets kicked twice by the Bragg field into momentum state ±2 Q B , is very unlikely. In practice, we include only the condensate states with momenta − Q B , 0, and Q B .
The equations of motion for the condensates are and for the fluctuations The Hartree shift of low momentum atoms is is the fraction of atoms in the excitations, and h H = 2U H n, where n is the total density of the gas. The off-diagonal fluctuation density where m n = k ĉ kĉ− k+n QB . Notice that the first three equations of motion (5), (6), and (7) are mean-field condensates but the last one (8) is an equation of motion for a fluctuation operator. From the fluctuation operator we form equations of motion for the fluctuation fields ĉ † kĉ p and ĉ kĉ − p . All anomalous fluctuation fields m 0 and m ±1 are ultraviolet divergent and need to be regularized.
Instead of coupling the condensate ψ 0 directly into the excitations ĉ † kĉ k , the Bragg field Ω provides the coupling into mobile condensate states ψ ±1 . Initially these mobile condensate states are empty but the Bragg pulse will break the translational symmetry of the condensate by rotating the condensed atoms into a superposition of the zero-momentum condensate state and the mobile condensate states. The twobody scattering processes couple these atoms into excitations, eventually leading into dephasing of the single-particle density matrix and, in principle, to fragmentation into separate condensates ψ 0 and ψ ±1 . As will be shown below, this coupling is relatively strong as it leads into rapid decay of the excited condensates. However, the present theory is unable to describe the dephasing process.
The set of equations of motion above shows that the Bragg field creates also new excitation fields, such as ĉ † kĉ k+ QB and ĉ kĉ− k± QB . In the calculations below, we have included the equations of motion for all these fields in addition to the standard excitation fields ĉ † kĉ k and ĉ kĉ− k . It is these mean-fields and the condensate fields ψ 0 , ψ ±1 that we propagate in real time in our theory. As our analysis below shows, all these fields are needed for a full description of even relatively weak Bragg pulses in which only a small fraction of atoms is excited by the field. In particular, the backward scattered condensate ψ −1 is important in order to obtain the proper Bogoliubov spectrum in the weakly interacting limit.
In all the numerical calculations below, we have considered a uniform 85 Rb condensate of density 10 14 cm −3 .

Regularization of the ultraviolet divergence
Before studying the Bragg spectroscopy any further, we will address some issues related to the regularization of the anomalous fluctuation fields. Using Matsubara Green's functions, the anomalous fluctuation field m 0 can be written as where G 0 ( K) is the bare Bose Green's function and G( K) is the dressed Green's function obtained for some self-energy Σ. The four-vector K = (iω, k) consists of a Matsubara frequency iω and a three-dimensional momentum k. The regular Hartree-Fock-Bogoliubov theory gives, and is given by, where However, the important point is that the bare Green's function needs to include the Hartree shifts, i.e.
Regularization can now be done self-consistently by removing the free particle diagrams where the free particle propagator is given by the same energy shifted Green's function G 0 . Thus we remove from m 0 the term This regularization differs from the standard scheme used in the literature [13,14,8,9,10] by the energy shift 2U n − µ in the bare Green's function G 0 ( K). Neglecting the energy shift, i.e. replacing ǫ ′ k by ǫ k in the regularizing term m 2B 0 will make the regularizing part larger than the initial anomalous term |m 0 | < m 2B 0 . This has the unfortunate effect of turning the regularized anomalous pairing field positive. The low momentum energy gap in the HFB spectrum √ −4U n 0 m 0 becomes then imaginary and the low energy excitations turn unstable. However, the imaginary energies are small and the corresponding lifetimes are long, so that the instability can easily go unnoticed. In addition, instabilities apply only to the very lowest energy excitations requiring a dense grid in the momentum space in order to play a role.
Including a constant energy shift δ in the bare Green's function does not affect the two-body scatterings but guarantees m 0 < 0 if δ is large enough. Indeed, choosing δ = 2U n − µ, as in (11), is enough to guarantee positivity of m 0 . However, this would require an iterative solution of δ as the chemical potential µ depends on δ. Here we approximate δ = U n but the results are relatively insensitive to the actual choice of δ as long as it is large enough. Regularization of the mobile pairing fields m ±1 proceeds in an identical manner.

Bragg spectroscopy of a weakly interacting BEC
We will first study the Bragg spectroscopy of a weakly interacting Bose gas. Since the vast majority of the atoms are in the initial condensate state ψ 0 , we are interested only in the energy of the excited condensate state with momentum Q B . In the weakly interacting regime, we can safely ignore the fluctuation parts of the finite momentum pairing and density fields m ±1 and n 1 . The equations of motion of the excited condensates are now (assuming that the populations of the excited condensate states ψ 1 and ψ −1 are small) and ih d dt In the absence of the Bragg coupling Ω(t) = 0, these equations can be solved by the Bogoliubov transformation, showing that the energies of the excited condensates at momenta ± Q B do indeed match the energies of the corresponding excitations in the Hartree-Fock-Bogoliubov specturm. Notice that, in addition to the Bragg field, the zero-momentum pairing field ∆ 0 acts as a source for the two condensates ψ ±1 even though it cannot level out the difference in the populations due to momentum conservation (i.e. a zero momentum pair of atoms can be turned into a pair ψ 1 ψ −1 which conserves momentum but also the population difference in the two states). This stability of the population difference between the ψ 1 and ψ −1 condensates underlines the importance of the off-diagonal fields ∆ ±1 and n 1 for the decay of the multiply condensed state. Indeed, it is only through these fields that the system finally decays into an equilibrium state. However, this decay process is slow in a weakly interacting gas, and can therefore be neglected in most cases.  allows the atoms to lower their energies even further by pair formation and this is reflected also in the Bragg spectrum. Fig. 1 shows the decay of the excited condensates after (and during) the Bragg pulse. The decay lifetime scales as h/U n, showing that the effect of the mobile pairing field is large even though the anomalous fields m ±1 are small. Notice that the total momentum is not quite conserved after the pulse. This is due to a finite momentum cutoff in the summations, but the error can be made arbitrarily small by increasing the momentum grid size. Fig. 2 shows the Bragg spectrum for a Bragg pulse exciting roughly 10% of atoms. Here and in the rest of the figures Bragg detuning hν B refers to the Bragg field energy offset from the free particle resonance hν =hω B −h 2 Q 2 B 2m . The figure shows two peaks corresponding to the forward and backward scattering and the shift from the free particle line with increasing interaction strength. The line shape also changes from a Lorentzian (for weakly interacting gas) into an asymmetric peak (for strongly interacting gas) because also quantum fluctuations are affected by the Bragg field and the corresponding atoms have the Bragg resonance at higher detuning than the condensate atoms (because the excitation spectrum is a concave function of momentum k). For finite temperatures, the asymmetricity of the peak would be even more pronounced as the condensate fraction is reduced.

Bragg spectroscopy of a strongly interacting BEC
From spectra such as shown in Fig. 2 we determine the position of the peak maximum and the Fig. 3 shows how this evolves as a function of pulse length. For very short pulses the linewidth is very large due to the lifetime broadening and the two peaks in the spectra overlap. Thus, the position of the peak maximum depends strongly on the pulse length. Once the pulse is long enough, so that the linewidth of the pulse is less than the distance between the backward and forward scattering peaks, the Bragg resonance can be resolved. There is still some slow drift in the    peak position due to the depletion of the initial condensate, causing some error in the determination of the resonance energy. As long as the total fraction of excited atoms is small enough (less than 10%), the Bragg spectrum is insensitive to the strength of the Bragg coupling. Therefore, in order to be able to resolve the Bragg resonance well, one can use very weak pulses. Furthermore, a Gaussian Bragg pulse would reduce the initial oscillations in the spectrum peak position. In an actual experimental setup, the maximum length of the pulse is limited due to the finite size of the system and the inhomogenous trapping potential. The Fig. 3 shows also how well the time-dependent Hartree-Fock-Bogoliubov theory reproduces the resonance energy of the corresponding static theory when the mobile pairing fields are neglected. Fig. 4 shows the Bragg resonance energy as a function of the interaction strength. At strong interactions, the resonance drops far from the mean-field shift line. Notice that the Bogoliubov line includes the effect of k-dependent scattering length and it agrees surprisingly well with the HFB theory. Without the k-dependent scattering length, the Bogoliubov line would be close to the mean-field shift line. The mean-field shift and the Bogoliubov line are calculated using a static theory but the Hartree-Fock-Bogoliubov lines are from the present time-dependent theory. The presence of the mobile pairing fields drops the resonance energy even further.
Neglecting the effect due to the mobile pairing fields ∆ ±1 , the effect of the na 3 corrections introduced by the HFB theory (as compared to the Bogoliubov result) is very small despite much larger effect on the chemical potential (roughly 10% for na 3 = 0.015). Interestingly, in the na 3 expansion of the corrections to the Bogoliubov theory, the leading order term in the chemical potential (the Lee-Huang-Yang (LHY) correction [15,16]) is 74%. However, the next higher order correction [17] is roughly 150% showing that the standard na 3 expansion is breaking down. It is also interesting to notice that the Popov approximation of the HFB theory (HFB-Popov theory) does reproduce correctly the LHY results in the leading order of the na 3 expansion. However, there are also higher order terms included. The HFB theory does actually agree very well with the gapless HFB-Popov theory for large Bragg momenta, and therefore we expect that the LHY terms are properly included in the present theory. The different treatment of the beyond LHY terms may be the reason for the discrepancy between our HFB results and the Beliaev theory [18] based results in [5].
The qualitative features of Fig. 4 agree with the experimental data in [5]. However, since the present study has been done for a uniform gas, this work should be extended to a trapped gas using local density approximation for proper comparison with the experimental results.

Summary
To summarize, we have studied the Bragg spectroscopy of a strongly interacting BEC using time-dependent Hartree-Fock-Bogoliubov theory. Taking into account momentum dependence of the scattering amplitude, the theory is in qualitative agreement with the experiment done on strongly interacting Rb-85 condensate. The most surprising effect comes from the creation of mobile pairing fields and the subsequent change in the excitation spectrum. While the present experimental results cannot confirm this effect, it should be more visible in the low-momentum excitation spectrum. Another interesting question would be the relation to fragmented condensates [19]. Because of the coherence of the Bragg field, the different condensates at momenta n Q B have a well defined phase-difference. The coupling to the excited states should, in principle, dephase the condensate, leading into a fragmented state. The Bragg scattering induced mobile pairing fields have also an interesting connection to the FFLO-type pairing [20,21] in polarized Fermi gases that would be worth further research.