\documentclass[aps,prl, amsmath, showpacs, preprintnumbers,superscriptaddress, twocolumn,sort&compress,floatfix, amssymb]{revtex4}
\pdfoutput=1
\usepackage{graphicx}
\usepackage{rotating}
\usepackage{dcolumn}
\usepackage{bm}
\usepackage{color}
\usepackage{mathptmx, textcomp}
\usepackage[latin1]{inputenc}
\usepackage{braket}
\usepackage{multirow}
\newcommand{\minitab}[2][l]{\begin{tabular}{#1}#2\end{tabular}}
\bibliographystyle{apsrev}
\hyphenation{Fesh-bach}
\begin{document}
\setcounter{figure}{4}
\section{Supplementary Material: Probing the Excitations of a Lieb-Liniger Gas from Weak to Strong Coupling}
\subsection{Lattice depth calibration and error bars}
The lattice depth $V_{x,y}$ is calibrated via Kapitza-Dirac diffraction. The statistical error for $V_{x,y}$ is 1\%, though the systematic error can reach up to 5\%.
The scattering length $a_s$ is calculated via its dependence on the magnetic field \cite{Mark11MAT} with an estimated uncertainty of $\pm 5 a_0$ arising from systematics in the magnetic field calibration and conversion accuracy. Additionally, the magnetic field gradient for sample levitation \cite{Weber2002MAT,Kraemer2004MAT} leads to a variation of less than $\pm 3 a_0$ across the sample.
\subsection{Atom number distribution across the 1D tubes}
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP1A.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP1B.pdf}
\caption{\label{SuppFIG1} (a) Atom number distribution over the 1D tubes for a total atom number $N=1.1 \times 10^5$ in the BEC and a 3D scattering length during loading $a_{\rm{3D}}=173a_0$. (b) Number of tubes that are filled with $N_{i,j} = N'$ atoms. (c) Total number of atoms in tubes filled with $N_{i,j} = N'$ atoms.}
\end{figure}
In order to calculate a mean value of the interaction parameter $\gamma$ and the Fermi wave vector $k_F$ we need to know the atom number distribution across the ensemble of tubes, as shown in Fig.~\ref{SuppFIG1}. We follow the calculation presented in \cite{Haller11MAT}.
The BEC is adiabatically loaded from a crossed dipole trap into the optical lattice. When the lattice is fully ramped up to its final depth $V_{x,y} = 30 E_R$, the laser beams forming the array of tubes give rise to an additional background harmonic confinement. The total background harmonic confinement of lattice and dipole trap laser beams is measured to $\omega_x = 2 \pi \times 15.4(0.1) \, \rm{Hz}$, $\omega_y = 2 \pi \times 20.1(0.1) \, \rm{Hz}$, and $\omega_z = 2 \pi \times 15.8(0.1) \, \rm{Hz}$. We deduce the atom number $N_{i,j}$ for the tube $(i,j)$ from the global chemical potential $\mu$. Assuming that interactions are sufficiently small during the loading process so that all tubes are in the 1D Thomas-Fermi (TF) regime, the local chemical potential in each tube reads
$$
\mu_{i,j} = \mu - \frac{1}{2} m (\lambda/2)^2 \left( \omega_x^2 i^2 + \omega_y^2 j^2\right) \, ,
$$
with $m$ the mass of the Cs atom and $\lambda=1064.5 \, \rm{nm}$ the wavelength of the lattice light. From $\mu_{i,j}$ the atom number in each tube can be derived via
$$
\mu_{i,j} = \left (\frac{3 N_{i,j}}{4 \sqrt{2}} g_{\rm{1D}} \omega_z \sqrt{m} \right)^{2/3} \, ,
$$
with $g_{\rm{1D}}$ the 1D coupling strength
$$
g_{\rm{1D}} = 2 \hbar \omega_\perp a_s \left( 1- 1.0326 \frac{a_s}{a_\perp} \right)^{-1} \, .
$$
Here, $a_\perp = \sqrt{\hbar / (m \omega_\perp)}$ is the radial oscillator length. The global chemical potential is calculated iteratively from the condition $N = \sum_ {i,j} N_{i,j} ( \mu )$. \\
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP2A.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP2B.pdf}\\
\caption{\label{SuppFIG2} 1D line-density distribution of a tube with $N=30$ atoms at $a_s=150 a_0$ (a) and $a_s=800 a_0$ (b). The blue line shows the result from the numerical solution of the Lieb-Liniger system. The dashed (dotted) gray lines depict the analytic results for the profile in the Thomas-Fermi (Tonks-Girardeau) limit \cite{Dunjko01MAT}. Here, $\omega_\perp=2\pi \times 15$ kHz and $\omega_z=15$ Hz. Finite temperature profiles using the Yang-Yang thermodynamic equations are shown for $T=10$nK (orange line) and $T=30$nK (red line).}
\end{figure}
\subsection{Density profile in the tubes}
For $T=0$ we model the 1D density distribution $n(z)$ in each tube individually by numerically solving the Lieb-Liniger system and making a local density approximation \cite{Lieb63MAT,Dunjko01MAT}. For $T>0$ we solve the Yang-Yang thermodynamic equations of the 1D Bose gas making a local density approximation to calculate $n(z,T)$ \cite{Yang69MAT}. Representative examples for two different values of $a_s$ at zero and finite temperature are shown in Fig.~\ref{SuppFIG2}.
\subsection{Mean $\gamma$ and $k_F$}
From the atom number distribution, we calculate the density profile for each tube individually as described above. The mean 1D density in each tube $n_{i,j}^{\rm{1D}}$ then delivers local and mean values for $\gamma$ and $k_F$
$$
\gamma_{i,j} = \frac{m g_{\rm{1D}}}{\hbar^2 n_{i,j}^{\rm{1D}}} \, , \qquad \gamma = \frac{1}{N} \sum\limits_{i,j} N_{i,j} \gamma_{i,j}
$$
and
$$
k_F^{i,j} = \pi n_{i,j}^{\rm{1D}} \, , \qquad k_F = \frac{1}{N} \sum\limits_{i,j} N_{i,j} k_F^{i,j} \, .
$$
As an example, the distribution of atoms over tubes with local $\gamma_{i,j}$ and $k/k_F^{i,j}$ is depicted in Fig.~\ref{SuppFIG3}(a) and (b) for a specific value of $a_s$. In Fig.~\ref{SuppFIG3}(c) and (d) we show the mean value of $\gamma$ and $k_F$ as a function of $a_s$.
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP3A.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP3B.pdf}\\
\vspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP3C.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP3D.pdf}
\caption{\label{SuppFIG3}(top row) Relative atom number in tubes with local interaction parameter $\gamma_{i,j}$ (a) and a local value for $k/k_F^{i,j}$ (b). As an example, the distributions are calculated for $a_s = 399 a_0$ and $N=1.1 \times 10^5$ resulting in a global $\gamma \, (k_F)$ of 11 $(3.7 \mu \rm{m}^{-1})$. (bottom row) Global values for $\gamma$ (c) and $k_F$ (d) as a function of the scattering length $a_s$. The points show the calculated values with $n_{\rm{1D}}$ computed numerically within the Lieb-Liniger model using a local density approximation. The dashed line in (d) is calculated with the density in the Thomas-Fermi limit while the dotted line shows the result using the density for a Tonks-Girardeau gas. The solid lines are interpolation functions to guide the eye. In (a)-(d) all quantities have been calculated using the mean value for $n_{\rm{1D}}$ in each tube.}
\end{figure}
\subsection{Sampling of the dynamical structure factor over the array of tubes and comparison with the experimental data}
In linear response, the energy $\Delta E (k,\omega,T)$ dumped into a single 1D gas after the Bragg pulse relates to the dynamical structure factor $S(k,\omega,T)$ via \cite{Pitaevskii03MAT}
$$
\Delta E (k,\omega,T) \propto \hbar \omega (1-e^{- \hbar \omega / (k_B T)}) S(k,\omega,T) \, .
$$
Provided that the system is probed at sufficiently large momentum $k$ the DSF of the trapped system can be computed in a local density approximation
$$
S(k,\omega,T) = \frac{1}{2L} \int \limits_{-L}^{L} S_{\rm{hom}}(k,\omega,T ; n(z)) dz \, ,
$$
with $L$ the system length and $S_{\rm{hom}}$ the DSF calculated for a homogeneous system with uniform density. In practice, we find that dividing the profile of each tube in $\approx 10$ homogeneous subsystems approximates the DSF of the trapped gas sufficiently well. The DSF in each tube with $N_{i,j}$ atoms obeys the $f$-sum rule \cite{Pitaevskii03MAT},
$
\int \limits_{-\infty}^{+\infty} \omega S(k, \omega, T) \, d \omega \propto N_{i,j} \, .
$
In combination with detailed balance $S(k,\omega)=e^{\hbar \omega / (k_B T)} S(k,-\omega)$ we find
$$
\int \limits_{0}^{+\infty} \omega S(k, \omega, T) (1- e^{- \hbar \omega / (k_B T)} ) \, d \omega \propto N_{i,j} \, .
$$
This allows us to calculate the dynamical response of the entire ensemble by first normalizing the DSF of each individual tube $(i,j)$ by its $f$-sum, respectively, and then weighting its contribution to the total signal by the number of atoms $N_{i,j}$. For a direct comparison with the experimental signal, we finally normalize the ensemble averaged response to unit area. \\
In the experiment, we extract $\langle p^2 \rangle$ as a function of $\omega$ from the momentum space distribution obtained after time-of-flight. For a direct comparison, each spectrum is normalized to unit area as the overall signal depends on the intensity of the Bragg lasers, which we slightly increase for data taken at stronger interactions. A changing offset due to an overall broadening of the momentum distribution with increasing $\gamma$ is subtracted from the data. This offset does not affect the shape of the excitation spectrum and stems from the broadening of the unperturbed momentum distribution with increasing $\gamma$.
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP4A.pdf}\\
\vspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP4B.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP4C.pdf}\\
\vspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP4D.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP4E.pdf}
\caption{\label{SuppFIG4}(a)-(c) Bragg excitation as a function of the laser pulse length measured in a weakly interacting 3D BEC (a), in 1D at $a_s = 15 a_0$ (b), and in 1D at $a_s = 819 a_0$ (c). The laser power in both Bragg beams is set to $0.5$mW (circles) and $0.2$mW (squares) in (a) and (b), and $0.5$mW (circles) and $0.3$mW (squares) in (c). (d),(e) Bragg excitation as a function of the laser power in both Bragg beams measured in a weakly interacting BEC (circles), in 1D at $a_s = 15 a_0$ (squares), and in 1D at $a_s = 819 a_0$ (diamonds). Here, the laser pulse length is fixed to $5$ms. The solid lines are exponentially damped sinusoids fit to the data. The dashed lines are linear fits in (a)-(c) and quadratic fits in (d) and (e) to the initial increase denoting the regime of linear response.}
\end{figure}
\subsection{The ABACUS algorithm}
The DSF for a homogeneous system is computed using the ABACUS algorithm. The computations are performed for a finite system of length $L$ with a finite particle number $N$ and with periodic boundary conditions. This means that the experimental situation is recovered only in the thermodynamic limit $N,L\rightarrow\infty$ with a density $n=N/L$ fixed through the local density approximation. We have confirmed that the particle numbers $N$ chosen for the computations (up to $N=128$) are large enough for the results to faithfully represent the thermodynamic limit. The only exceptions are the highest temperature curves in Fig.~3 (b)-(d) where some oscillations due to the finite size are still visible at high energies. The correlation function is evaluated numerically by summing contributions according to the Lehmann spectral representation. The spectral sum is infinite and needs to be truncated. The ABACUS algorithm performs this operation in an efficient way capturing the most relevant contributions. The error caused by the truncation is easily tractable by evaluating the $f$-sum rule and, for the presented data, does not exceed 5\% of the total spectral weight.
\subsection{Regime of linear response}
Measuring the DSF via Bragg spectroscopy requires to probe the system in the regime of linear response. This we tested experimentally by validating that the excitations created depend linearly on the pulse length and quadratically on the intensity of the Bragg lasers for the range of parameters used in the experiment, see Fig.~\ref{SuppFIG4} \cite{Brunello01MAT}.
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP5A.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP5B.pdf}
\caption{\label{SuppFIG5} (a) Heating effects during the ramp of the scattering length $a_s$ on the Bragg spectroscopy data. Bragg spectroscopy data taken directly after the ramp of the scattering length to 173$a_0$ within 50 ms (squares) is compared to data taken after ramping to 819$a_0$ and back within 100 ms (circles). (b) Full width at half maximum (FWHM) extracted via the Gaussian fit (see main article) from the Bragg spectroscopy data as a function of $\gamma$. The dashed (solid) arrow denotes an estimate for line broadening due to the finite pulse length (quantum finite size effect) in the non-interacting limit.}
\end{figure}
\subsection{Heating effects on the excitation spectra}
Besides the detailed theoretical analysis of finite temperature effects on the measured DSF reported in the main article, we have checked in an experiment that heating during the ramp of $a_s$ to large values only marginally influences the shape of the excitation spectrum when compared to the effect of interactions. In Fig.~\ref{SuppFIG5}(a) we show Bragg excitation spectra taken at moderate interaction strength, corresponding to the data shown in Fig.~2(b) of the main article. For the two datasets shown, we either ramp directly to 173$a_0$ (squares) or ramp deep into the regime of strong interactions, $a_s$=810$a_0$, and back (circles) prior to applying the Bragg pulse.
\subsection{Interaction-independent spectral width at small $\gamma$}
In the limit of weak interactions we observe a width of the excitation spectra that levels off at a constant value as shown in Fig.~\ref{SuppFIG5}(b). We attribute this to the onset of two effects that prevent the observation of a $\delta$-like peak as expected in a non-trapped homogeneous system for $\gamma \rightarrow 0$. First, the finite length of the Bragg excitation pulse causes Fourier broadening, which results in an estimated residual width of $\approx 120 \rm{Hz}$ in the non-interacting limit \cite{Pitaevskii03MAT}. Second, the finite size of the sample due to the presence of the harmonic trap leads to an uncertainty-limited energy width (not accounted for in the LDA) that can be estimated to $\delta \omega = 2 \hbar k / (m a_z )$ for a sample of non-interacting particles \cite{Golovach09MAT}. Here, $a_z = \sqrt{\hbar/(m \omega_{z})}$ denotes the quantum length scale of the trap.
\begin{figure}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP6A.pdf}
\hspace{2mm}
\includegraphics[width=0.48\columnwidth]{PRL_BRAGG1D_SUPP6B.pdf}
\caption{\label{SuppFIG6} Momentum distribution of two states excited with a single (particle-hole) excitation over the ground state for $\gamma=64$. In the first case (a) the momentum is carried by a hole-like Lieb II excitation, in the second (b) by a particle-like Lieb I excitation. Both states have the same total momentum equal to $k/k_F=1$. The calculation is done for a system of 32 particles. The shaded areas indicate the difference to the ground-state momentum distribution.}
\end{figure}
\subsection{Momentum distribution of Lieb-I and Lieb-II excitations}
In order to illustrate the effects of additional particle-hole excitations on the momentum distribution, we consider two situations in which the total momentum of an excited state (obtained by adding a single particle-hole excitation with total momentum equal to the Fermi momentum on the ground state) is carried either by the hole (Lieb type II) mode (Fig.~\ref{SuppFIG6}(a)) or by the particle (Lieb type I) mode (Fig.~\ref{SuppFIG6}(b)). The computations are performed again using the ABACUS method~\cite{Caux07MAT}. The results show that while in both cases the momentum distribution function is augmented by broadened peaks, the hole excitation leads to an overall broadening of $n(k)$ while the particle excitation appears as a more distinguishable peak. This highlights the collective nature of the Lieb-II hole-like excitations and qualitatively reproduces the features seen in Fig.~4 of the main text.
\bibliographystyle{apsrev}
\begin{references}
\bibitem{Mark11MAT}M. J. Mark, E. Haller, K. Lauber, J. G. Danzl, A. J. Daley, and H.-C. N\"agerl, Phys. Rev. Lett. \textbf{107}, 175301 (2011).
\bibitem{Weber2002MAT} T. Weber, J. Herbig, M. Mark, H.-C. N\"agerl, and R. Grimm, Science \textbf{299}, 232 (2002).
\bibitem{Kraemer2004MAT} T. Kraemer, J. Herbig, M. Mark, T. Weber, C. Chin, H.-C. N\"agerl, and R. Grimm, Appl. Phys. B \textbf{79}, 1013 (2004).
\bibitem{Haller11MAT}E. Haller, M. Rabie, M. J. Mark, J. G. Danzl, R. Hart, K. Lauber, G. Pupillo, and H.-C. N\"agerl, Phys. Rev. Lett. \textbf{107}, 230404 (2011).
\bibitem{Lieb63MAT}E. H. Lieb and W. Liniger, Phys. \textbf{130}, 1605 (1963).
\bibitem{Dunjko01MAT}V. Dunjko, V. Lorent, and M. Olshanii, Phys. Rev. Lett. \textbf{86}, 5413 (2001).
\bibitem{Yang69MAT}C. N. Yang and C. P. Yang, J. Math. Phys. \textbf{10}, 1115 (1969).
\bibitem{Pitaevskii03MAT}L. Pitaevskii and S. Stringari, \textit{Bose-Einstein Condensation}, (Oxford Univ. Press, New York, 2003).
\bibitem{Brunello01MAT}A. Brunello, F. Dalfovo, L. Pitaevskii, S. Stringari, and F. Zambelli, Phys. Rev. A \textbf{64}, 063614 (2001).
\bibitem{Golovach09MAT} V. N. Golovach, A. Minguzzi, and L. I. Glazman, Phys. Rev. A \textbf{80}, 043611 (2009).
\bibitem{Caux07MAT}J.-S. Caux, P. Calabrese, and N. A. Slavnov, J. Stat. Mech. \textbf{P01008} (2007).
\end{references}
\end{document}