Sounds waves and fluctuations in one-dimensional supersolids

Our first system is a Bose gas of magnetic atoms confined in a tube shaped potential [see Fig. 1(a)]. The single particle Hamiltonian is

$$H_{\mathrm{sp}}=-\frac{\hbar^{2}\nabla^{2}}{2m}+\frac{1}{2}m\omega^{2}(y^{2}+z^{2}),$$

(1)

where $\omega$ is the transverse harmonic confinement angular frequency, and the atoms are free to move along the $x$-axis. The atomic dipole moments are polarized along $z$ by a bias field and the interactions are described by the potential

$$U(\mathbf{r})=\frac{4\pi a_{s}\hbar^{2}}{m}\delta(\mathbf{r})+\frac{3a_{\mathrm{dd}}\hbar^{2}}{mr^{3}}\left(1-3\frac{z^{2}}{r^{2}}\right),$$

(2)

where $a_{s}$ is the $s$-wave scattering length, $a_{\mathrm{dd}}=m\mu_{0}\mu_{m}^{2}/12\pi\hbar^{2}$ is the dipole length, and $\mu_{m}$ is the atomic magnetic moment. The ratio $\epsilon_{\mathrm{dd}}=a_{\mathrm{dd}}/a_{s}$ characterises the relative strength of the dipole-dipole to s-wave interactions. When this parameter is sufficiently large the ground state undergoes a transition to a crystalline state with modulation along $x$. Quantum fluctuation effects are necessary to stabilize this state and further details about this model are given in Appendix A.1. The value of $\epsilon_{\mathrm{dd}}$ where the transition occurs depends on $\omega$ and $\rho$, and the transition can be continuous or discontinuous. In the results we present here are for ¹⁶⁴Dy atoms with $a_{\mathrm{dd}}=130.8\,a_{0}$ and $\omega/2\pi=150\,$Hz.

Our second system is a BEC of atoms free to move in 1D, but interacting via a soft-core interaction potential $U_{\mathrm{sc}}(x)=U_{0}\theta(a_{\mathrm{sc}}-|x|)$, where $a_{\mathrm{sc}}$ is the core radius and $U_{0}$ is the potential strength [see Fig. 1(b)]. It conventional to define the dimensionless interaction parameter

$\displaystyle\Lambda=\frac{2ma_{\mathrm{sc}}^{3}U_{0}\rho}{\hbar^{2}},$

(3)

and continuous transition occurs from a uniform to a modulated state at the critical value $\Lambda_{c}=21.05$. Further details about this model are given in Appendix A.2.

By solving the (generalized) Gross-Pitaevskii equations for these models we determine the energy minimising wavefunction $\psi_{0}$. These solutions can be found in a single unit cell of length $a$ along $x$. If the state is crystalline, then $a$ is the lattice constant (for uniform states the energy is independent of $a$). Take $\mathcal{E}$ to be the energy per unit length of this state for density $\rho$ and lattice constant $a$ (as defined in Appendix A). In the thermodynamic limit where the system is free to choose the microscopic length $a$ to minimise the energy, the ground state is

$\displaystyle\mathcal{E}_{0}(\rho)=\min_{a}\mathcal{E}(\rho,a),$

(4)

which can be used to implicitly determine $a(\rho)$ as the value of $a$ that minimises $\mathcal{E}$ at density $\rho$.

Here our interest is in a Lagrangian description of the low-energy, long-wavelength features of a system where 1D crystalline order can occur along the $x$-direction. We work in a thermodynamic limit where the ground state is determined by the average linear number density $\rho$ in addition to the interaction and relevant transverse confinement parameters. Long-wavelength perturbations of the ground state can be characterized by three fields: $\{\delta\rho(x,t),u(x,t),\phi(x,t)\}$. Here $\delta\rho(x,t)$ denotes variations in the average density, $u$ denotes the deformation field of the crystal lattice along the $x$ direction (relative to the unstrained ground state), and $\phi$ denotes the superfluid phase field. For treating small departures from equilibrium the Lagrangian density can be expanded to quadratic order in these fields as

	$\displaystyle\mathcal{L}=$	$\displaystyle-\hbar\,\delta\rho\,\partial_{t}\phi-\frac{\alpha_{\rho\rho}}{2}(\delta\rho)^{2}-\frac{\alpha_{uu}}{2}(\partial_{x}u)^{2}-\alpha_{\rho u}\delta\rho\,\partial_{x}u$
		$\displaystyle+\frac{1}{2}m\rho_{n}\Bigl{(}\partial_{t}u-\frac{\hbar}{m}\partial_{x}\phi\Bigr{)}^{2}-\rho\frac{\hbar^{2}}{2m}(\partial_{x}\phi)^{2},$		(5)

where $\rho_{n}$ is the normal density, with $\rho=\rho_{s}+\rho_{n}$, and $\rho_{s}$ being the superfluid density. Similar forms of Lagrangian densities have been obtained by other authors (e.g. see [21, 22, 18, 19, 20]), although here we follow the approach of Yoo and Dorsey [18]. A discussion about the relationship between these treatments is given in Refs. [18, 19]. The phase gradient relates to the superfluid velocity as $v_{s}=\frac{\hbar}{m}\partial_{x}\phi$ and the normal velocity relates to the time-derivative of the lattice deformation field $v_{n}=\partial_{t}u$ [17, 18].

In addition to the superfluid density, which determines the rigidity of the state to twists in the phase $\phi$ [31], three other generalized elastic parameters of the system $\{\alpha_{\rho\rho},\alpha_{\rho u},\alpha_{uu}\}$ also appear in Eq. (5). Here we discuss how these can be obtained from ground state calculations.

For 1D crystals the superfluid fraction $f_{s}=\rho_{s}/\rho$ has

	$\displaystyle f_{s}^{+}$	$\displaystyle=\frac{a}{\rho}\left[\int_{\mathrm{uc}}\frac{dx}{\int dy\,dz\,|\psi_{0}|^{2}}\right]^{-1},$		(6)
	$\displaystyle f_{s}^{-}$	$\displaystyle=\frac{a}{\rho}\int dy\,dz\left[\int_{\mathrm{uc}}\frac{dx}{|\psi_{0}|^{2}}\right]^{-1},$		(7)

as upper and lower bounds, respectively [32, 33]. Hereon, we neglect the transverse coordinates $(y,z)$ when applying results to the soft-core case. Both bounds are identical for the soft-core case, and provide the exact superfluid fraction (e.g. see [27]). For the dipolar case the bounds are close to each other and taking $f_{s}=\frac{1}{2}(f_{s}^{+}+f_{s}^{-})$ provides an accurate estimate of the superfluid fraction (see [34]).

The other generalized elastic parameters can be determined from the energy density $\mathcal{E}(\rho,a)$, for the ground state under the constraint of lattice constant being $a$ and the mean density being $\rho$. The first elastic parameter is defined by

$\displaystyle\alpha_{\rho\rho}=\left(\frac{\partial^{2}\mathcal{E}}{\partial\rho^{2}}\right)_{a},$

(8)

and relates to the isothermal compressibility at constant strain:

$\displaystyle\tilde{\kappa}=\frac{1}{\rho^{2}\alpha_{\rho\rho}}.$

(9)

The second parameter describes the effects of straining a site at constant density

$\displaystyle\alpha_{uu}=a^{2}\left(\frac{\partial^{2}\mathcal{E}}{\partial a^{2}}\right)_{\rho}.$

(10)

Noting that lattice strain, given by $\partial_{x}u$, relates to changes in the lattice constant according to $\delta a/a=\partial_{x}u$. In higher dimensional crystals the $\alpha_{uu}$ term generalizes to the elastic tensor. For the 1D crystal $\alpha_{uu}$ can also be identified as the layer-compressibility used to characterize smectic materials [35, 19]. Finally, we define the density-strain coupling parameter by the mixed partial derivative

$\displaystyle\alpha_{\rho u}=a\left(\frac{\partial^{2}\mathcal{E}}{\partial\rho\partial a}\right).$

(11)

The Euler-Lagrange equations obtained from Eq. (5) are

	$\displaystyle\hbar\partial_{t}\phi$	$\displaystyle=-\alpha_{\rho\rho}\delta\rho-\alpha_{\rho u}\partial_{x}u,$		(12)
	$\displaystyle\rho_{n}(m\partial^{2}_{t}u-\hbar\partial_{tx}\phi)$	$\displaystyle=\alpha_{uu}\partial^{2}_{x}u+\alpha_{\rho u}\partial_{x}\delta\rho,$		(13)
	$\displaystyle\partial_{t}\delta\rho$	$\displaystyle=-\rho_{n}\partial_{tx}u-\rho_{s}\frac{\hbar}{m}\partial^{2}_{x}\phi.$		(14)

These describe the hydrodynamic evolution for small departures from equilibrium.

We look for normal mode solutions of Eqs. (12) - (14) of the form $X(x,t)=X_{w}e^{i(qx-\omega t)}$, where $X$ denotes our three fields appearing in the Lagrangian and $q$ is a quasimomentum. The resulting linear system is

$\displaystyle\mathsf{M}\begin{pmatrix}\delta\rho_{w}\\ \phi_{w}\\ u_{w}\end{pmatrix}=\mathbf{0},$

(15)

where

$\displaystyle\mathsf{M}$

$\displaystyle=\begin{pmatrix}\alpha_{\rho\rho}&-i\hbar\omega&iq\alpha_{\rho u}\\ i\hbar\omega&\hbar^{2}q^{2}\rho_{s}/m&-\hbar\omega q\rho_{n}\\ iq\alpha_{\rho u}&\hbar\omega q\rho_{n}&\omega^{2}\rho_{n}m-q^{2}\alpha_{uu}\end{pmatrix}.$

(16)

Nontrivial solutions require the matrix $\mathsf{M}$ to be singular, i.e. its determinant $\Delta_{\mathsf{M}}$ must be zero, where

$\displaystyle\Delta_{\mathsf{M}}=m\hbar\rho_{n}(\omega^{2}-c_{+}^{2}q^{2})(\omega^{2}-c_{-}^{2}q^{2}).$

(17)

Here we have introduced

$\displaystyle mc_{\pm}^{2}$

$\displaystyle=\frac{1}{2}(a_{\Delta}\pm\sqrt{a_{\Delta}^{2}-4b_{\Delta}}),$

(18)

where

	$\displaystyle a_{\Delta}$	$\displaystyle=\rho\alpha_{\rho\rho}-2\alpha_{\rho u}+\frac{\alpha_{uu}}{\rho_{n}},$		(19)
	$\displaystyle b_{\Delta}$	$\displaystyle=\frac{\rho_{s}}{\rho_{n}}(\alpha_{\rho\rho}\alpha_{uu}-\alpha_{\rho u}^{2}).$		(20)

Thus the hydrodynamic excitation frequencies of the two solutions are gapless and given by

$\displaystyle\omega_{\pm}=c_{\pm}q,$

(21)

with $c_{\pm}$ being the speeds of sound.

We compare the speeds of sound determined from direct calculation of the excitations to the hydrodynamic result of Eq. (18) in Figs. 2 and 3 for the dipolar and soft-core systems, respectively. The direct calculation of the excitations involves solving the Bogoliubov-de Gennes (BdG) equations (e.g. see [36, 7, 8, 9, 37]). We make this comparison across the uniform superfluid to supersolid transition. In Figs. 2(a) and 3(a) we show the density contrast $\mathcal{C}$ defined as

$\displaystyle\mathcal{C}=\frac{\max\varrho_{0}-\min\varrho_{0}}{\max\varrho_{0}+\min\varrho_{0}},$

(22)

where $\varrho_{0}(x)=\int dydz|\psi_{0}(\mathbf{x})|^{2}$ is the line density along $x$. The contrast reveals the appearance of density modulations, and hence acts as an order parameter for the formation of crystalline order. In the dipolar system, at the densities we consider here [38, 37], the crystalline order develops continuously as $\epsilon_{\mathrm{dd}}$ increases past the critical value of $\epsilon_{\mathrm{dd},c}$. For the soft-core system the transition is also continuous, and occurs when the interaction parameter $\Lambda$ exceeds $\Lambda_{c}$. In both systems, the superfluid fraction decreases with increasing density contrast.

The numerical solution of the BdG equations [37] yields the energies $\hbar\omega_{\nu q}$ and amplitude functions $\{\mathrm{u}_{\nu q},\mathrm{v}_{\nu q}\}$, where $\nu$ is the band index and $q$ is the quasimomentum along $x$. We then obtain the speeds of sound $c_{\nu}$ of gapless branches from the BdG results by fitting $\omega_{\nu q}$ to $c_{\nu}q$, for small values of $q$. In the uniform superfluid phase only a single branch is gapless, and we identify the speed of sound of this branch as $\nu=+$. In the crystalline state there are two gapless branches, and we assign the $\nu$ of the lower (upper) branches as $-$ ($+$). Our results for the speeds of sound are shown in Figs. 2(b) and 3(b) for the dipolar and soft-core systems, respectively. There are significant differences in behavior revealed by these results, particularly in the nature of the discontinuity of the speeds of sound at the transition point. For the dipolar system $c_{+}$ jumps downwards as $\epsilon_{\mathrm{dd}}$ crosses the transition¹¹1For the $\rho$ considered the jump is small and hard to see. Ref. [37] shows results at other densities where the discontinuity is more clearly revealed.. In contrast, for the soft-core system the $c_{+}$ speed of sound jumps upwards at the transition. Other aspects of the behavior near the transition will be discussed further in Sec. IV.3.

Since the BdG speeds of sound and the hydrodynamic results are in good agreement, in Figs. 2(c) and 3(c) we show the difference between the results. The relative error is $\lesssim 1\%$ for the dipolar gas and significantly smaller for the soft-core model, except very close to the transition. The error for the dipolar model is close to the accuracy that we can determine the speeds of sound from the BdG calculations. We also show the difference arising from neglecting the density-strain term (i.e. setting $\alpha_{\rho u}=0$, which makes our theory equivalent to that in Refs. [19, 25]), noting that over the parameter range considered this causes a $\sim 5\%$ relative upward-shift in the $c_{+}$ speed of sound and a smaller change in $c_{-}$. Furthermore, the role of the density-strain term vanishes as we approach the transition.

The hydrodynamic results presented here are determined by generalized elastic parameters (i.e. $\{f_{s},\alpha_{\rho\rho},\alpha_{\rho u},\alpha_{uu}\}$), which are determined from the ground state calculations as discussed in Sec. III.2. The superfluid fraction results are shown in Figs. 2(a) and 3(a), and the other elastic parameters are shown in Figs. 2(d) and 3(d). These results reveal that the crystal-dependent elastic terms $\alpha_{\rho u}$ and $\alpha_{uu}$ both vanish as we approach the transition from the modulated side. For both systems we find that $\rho\alpha_{\rho u}$ is typically smaller than $\alpha_{uu}$, suggesting that neglecting this parameter can be a reasonable first approximation. The parameter $\alpha_{\rho\rho}$ is non-zero in both phases, but changes discontinuously at the transition (also see [37]).

From the hydrodynamic results for the speed of sound in Eq. (18) we derive limiting results to describe their behavior approaching the transition and deep into the crystal regime. It is useful to introduce two energy scales as

	$\displaystyle mc_{\kappa}^{2}$	$\displaystyle\equiv\rho\frac{\alpha_{\rho\rho}\alpha_{uu}-\alpha_{\rho u}^{2}}{\alpha_{uu}}=\frac{1}{\rho\kappa},$		(23)
	$\displaystyle mc_{u}^{2}$	$\displaystyle\equiv\frac{\alpha_{uu}}{\rho_{n}},$		(24)

being the bulk compressibility and lattice compressibility energies, respectively, defining the associated characteristic speeds $c_{\kappa}$ and $c_{u}$. Here

$\displaystyle\kappa=\frac{\alpha_{uu}}{\rho^{2}(\alpha_{\rho\rho}\alpha_{uu}-\alpha_{\rho u}^{2})},$

(25)

is the isothermal compressibility. Because $\alpha_{\rho u}$ is relatively small in both our models $\kappa\approx\tilde{\kappa}$ [see Eq. (9), Figs. 2(d) and 3(d)]. The lattice compressibility is only non-zero in the modulated state. While both $\alpha_{uu}$ and $\rho_{n}$ go to zero as we approach the transition from the modulated side, $mc_{u}^{2}$ has a non-zero value [see Figs. 2(d) and 3(d)].

For a system deep in the modulated regime the ground state consists of atoms relatively localized in each unit cell (see insets to Fig. 5). This regime should be well-approximated by the linear-chain model, in which a 1D solid is approximated by a set of masses joined by springs, yielding a simple theory for the phonons [40]. This idea was recently explored and found to be a good description of the dynamics of few droplets of a dipolar BEC confined by a three-dimensional cigar-shaped harmonic trap [41].

To apply this model to the systems we consider here, we denote the center-of mass position of the droplet at the $j$-th site by the coordinate $x_{j}$, with the equilibrium position being $x_{j}^{0}=ja$. In the linear-chain model the equation of motion for the droplet positions is

$\displaystyle mN_{\mathrm{uc}}\ddot{x}_{j}=-k_{\mathrm{sp}}(2x_{j}-x_{j-1}-x_{j+1}),$

(39)

where $N_{\mathrm{uc}}=\rho a$ is the number of atoms in each unit cell (localized droplet) and $k_{\mathrm{sp}}$ is the spring constant. We identify the spring constant from how the energy of the droplet²²2Taking the droplet energy to be the energy in the unit cell, i.e. $E=\mathcal{E}a$. $E$ varies with the inter droplet distance $a$:

$\displaystyle k_{\mathrm{sp}}=\left(\frac{\partial^{2}E}{\partial a^{2}}\right)_{N_{\mathrm{uc}}},$

(40)

notably with $N_{\mathrm{uc}}$ held constant. This linear-chain model has a well-known phonon dispersion relation of the form

$\displaystyle\omega_{\mathrm{sp}}(q)=\frac{2c_{\mathrm{sp}}}{a}\left|\sin\left(\frac{qa}{2}\right)\right|,$

(41)

where $c_{\mathrm{sp}}=a\sqrt{{k_{\mathrm{sp}}}/{mN_{\mathrm{uc}}}}$ is the speed of sound and $k$ is the wavevector. Evaluating this result in terms of the generalized elastic parameters yields Eq. (34) for $c_{\mathrm{sp}}$.

In Fig. 5 we compare the dispersion relations for the lowest two gapless bands obtained from the BdG calculations to the linear-chain dispersion relation (41), noting that the latter only depends on the elastic parameters $k_{\mathrm{sp}}$ obtained from ground state solutions. The cases considered for these results are chosen to be deep in the crystal regime with $f_{s}\approx 8\%$ (3%) for the dipolar (soft-core) case. The results show that the low $q$ behavior of the upper gapless band is reasonably well predicted by $\omega_{\mathrm{sp}}(q)$ [also see $c_{\mathrm{s}}$ and $c_{+}$ in Figs. 2(b) and 3(b)]. At higher $q$ values, approaching the band edge, the linear chain dispersion slightly overestimates (underestimates) the frequency of the upper excitation band for dipolar (soft-core) model. These results confirm that the excitations in this band are closely related to the classical crystal phonon modes.

The total density fluctuation

$\displaystyle\!\!\delta\varrho(x)$

$\displaystyle=2\mathrm{Re}\left\{c_{\nu q}\!\int dy\,dz\,\psi_{0}(\mathbf{x})[\mathrm{u}_{\nu q}(\mathbf{x})-\mathrm{v}_{\nu q}(\mathbf{x})]\right\},$

(44)

is the leading order expression in $c_{\nu q}$ for $\varrho(x)-\varrho_{0}(x)$, and we have taken $\psi_{0}$ to be real. In the modulated phase the density varies rapidly with $x$, with a dominant periodicity set by the lattice constant $a$. Here we will focus on the density fluctuations arising on length scales much larger than the lattice constant, i.e. with wave vectors in the first Brillouin zone. For a single quasiparticle ($c_{\nu q}=1$) the only non-zero fluctuation is at wavevector $q$ and $-q$, with

$\displaystyle\delta\tilde{\rho}_{\nu q}=\int d\mathbf{x}\,e^{-iqx}\psi_{0}(\mathbf{x})[\mathrm{u}_{\nu q}(\mathbf{x})-\mathrm{v}_{\nu q}(\mathbf{x})].$

(45)

We show an example of the density change to a modulated ground state with the addition of a quasiparticle in Fig. 6(a), and the associated long-wavelength total density fluctuation is visualized in Fig. 6(c) [also see Appendix B].

With crystalline order we can also consider the effect of the excitation (42) on the crystal structure. This is characterized by the displacement field $u(x)$ which describes the distance a point originally at $x$ displaces. We identify the lattice site locations as the local maxima of the linear density. The unperturbed locations being $x_{j}^{0}=ja$ (i.e. maxima of $\varrho_{0}$), and the locations for the perturbed case being $x_{j}=x_{j}^{0}+u(x_{j}^{0})$ (i.e. maxima of $\varrho$). The displacement field can be expressed as

$\displaystyle u(x)=2\mathrm{Re}\{c_{\nu q}u_{\nu q}e^{iqx}\}.$

(46)

where the amplitude $u_{\nu q}$ is a measure of how the quasiparticle displaces the site at the origin³³3In calculations we evaluate $u_{\nu q}=-\delta\rho_{\nu q}^{\prime}(0)/\varrho_{0}^{\prime\prime}(0)$, where the prime is derivative with respect to $x$ and $\delta\rho_{\nu q}(x)=\int dy\,dz\,\psi_{0}[\mathrm{u}_{\nu q}-\mathrm{v}_{\nu q}]$. In Fig. 6(b) we show (46) for the perturbed density of Fig. 6(a), and give an example of its relationship to the site displacements [inset to Fig. 6(b)].

To understand the effect of lattice deformation on the density we consider a weak displacement field of the form (46) on a supersolid state with equilibrium density profile $\varrho_{0}(x)$. Assuming that the atom number in each cell remains constant then the purely-deformed density profile is given by

$\displaystyle\varrho_{u}(x)={\varrho_{0}\left(x-u(x)\right)}{\left(1-\partial_{x}u(x)\right)},$

(47)

where the first factor describes the displaced density profile and second factor adjusts the envelope of density profile to keep the atom number per cell constant. The change in cell width, i.e. strain, leads to a change in density across the system, giving rise to the lattice strain density fluctuation. Fourier transforming $\delta\varrho_{u}(x)\equiv\varrho_{u}(x)-\varrho_{0}(x)$ for a single quasiparticle yields [cf. Eq. (45)]

$\displaystyle\delta\tilde{\rho}_{u,\nu q}=-iq\rho u_{\nu q}L,$

(48)

and $\delta\tilde{\rho}_{u,\nu-q}=\delta\tilde{\rho}_{u,\nu q}^{*}$ as the nonzero contributions in the first Brillouin zone. In Eq. (48) $L$ is the length of the $x$ region integrated over.

We illustrate the lattice strain density fluctuation construction in Fig. 6. Using the identified deformation field $u(x)$ [Fig. 6(b)], in Fig. 6(a) we show the purely-deformed density profile $\varrho_{u}(x)$, with associated density fluctuation shown in Fig. 6(c).

In Fig. 6(a) we observe that while the maxima and minima of the density $\varrho(x)$ (i.e. lattice distortion) coincide with those of $\varrho_{u}(x)$, these functions are not identical. The difference arises from the tunnelling of particles between sites, leading to a change in the atom number per site. This is generally referred to as the defect density [14, 17, 18], and can be defined as $\varrho_{\Delta}(x)=\varrho(x)-\varrho_{u}(x)$, or equivalently

$\displaystyle\delta\tilde{\rho}_{\Delta,\nu q}=\delta\tilde{\varrho}_{\nu q}-\delta\tilde{\varrho}_{u,\nu q}.$

(49)

We illustrate the defect density fluctuation in Fig. 6(c).

The decomposition of the total density fluctuation into strain and defect contributions is compared for an example case in Fig. 6(c). Here we see that the strain and defect densities have a relatively large amplitude and are out of phase, such that the total density fluctuation is small. This behavior of the lattice straining in the opposite sense to the flow of atoms in the lowest gapless band reveals a general property of the lowest energy band we have used as an example in Fig. 6(a)-(c). Similar behavior has been observed in experiments with dipolar supersolids as a low frequency Nambu-Goldstone excitation [10, 11, 5].

It is useful to define the quantities

	$\displaystyle S_{\rho,\nu}(q)$	$\displaystyle=\frac{|\delta\tilde{\rho}_{\nu q}|^{2}}{N},$		(50)
	$\displaystyle S_{u,\nu}(q)$	$\displaystyle=\frac{|\delta\tilde{\rho}_{u,\nu q}|^{2}}{N},$		(51)
	$\displaystyle S_{\Delta,\nu}(q)$	$\displaystyle=\frac{|\delta\tilde{\rho}_{\Delta,\nu q}|^{2}}{N},$		(52)

to characterise the strength of the total density, strain density and defect density fluctuations, arising from the excitation with quantum numbers $\nu$ and $q$, where $N=L\rho$. We note that $S_{\rho}(q)=\sum_{\nu}S_{\rho,\nu}(q)$ is the usual static structure factor, and that aspects of the static and dynamical structure factors for the dipolar supersolid have been discussed in Refs. [5, 12, 37].

In Fig. 6(c) we show the dispersion relations, and in Fig. 6(e)-(g) the related fluctuation strengths introduced in Eq. (50)-(52), for the lowest two gapless bands of a dipolar supersolid. The lowest band ($\nu=-$, blue) has strong strain and defect density fluctuations as $q\to 0$, but these are out-of-phase, leading to a small total density fluctuation. The second band ($\nu=+$, red) has weak defect density fluctuations, and the strain density fluctuation makes the dominant contribution to the total density fluctuation. Hence the second band can be considered a crystal-like excitation, consistent with the conclusions from the linear-chain model presented in Sec. IV.4.

The hydrodynamic description (5) also allows the calculation of the density correlation functions (see [18]). In particular, our interest is in the imaginary part of the density response function, which has the form

$\displaystyle\chi_{o}^{\prime\prime}(q,\omega)$

$\displaystyle=\frac{1}{2}Nq\pi\sum_{\nu}\xi^{\nu}_{o}[\delta(\omega-c_{\nu}q)-\delta(\omega+c_{\nu}q)],$

(53)

where the subscript $o$ can take the values $\{\rho,u,\Delta\}$ to represent the total density, strain density and defect density, respectively. This result is limited to the lowest two bands (i.e. $\nu=\{-,+\}$) in the modulated regime and to a single band (i.e. $\nu=\{+\}$) in the uniform regime. The hydrodynamic results for the speeds of sound $c_{\nu}$ have already been introduced, and $\xi_{o}^{\nu}$ is a correlations length that represents the contribution of band $\nu$ to the response function. For the modulated case these are given by

	$\displaystyle\xi_{\rho}^{\pm}$	$\displaystyle=\frac{\hbar}{m\rho}\frac{\rho mc_{\pm}^{2}-\frac{\rho_{s}}{\rho_{n}}\alpha_{uu}}{mc_{\pm}(c_{\pm}^{2}-c_{\mp}^{2})},$		(54)
	$\displaystyle\xi_{u}^{\pm}$	$\displaystyle=\frac{\hbar\rho}{m\rho_{n}}\frac{mc_{\pm}^{2}-\rho_{s}\alpha_{\rho\rho}}{mc_{\pm}(c_{\pm}^{2}-c_{\mp}^{2})},$		(55)
	$\displaystyle\xi_{\Delta}^{\pm}$	$\displaystyle=\frac{\hbar\rho_{s}}{m\rho\rho_{n}}\frac{\rho mc_{\pm}^{2}-(\alpha_{uu}-2\rho\alpha_{\rho u}+\rho^{2}\alpha_{\rho\rho})}{mc_{\pm}(c_{\pm}^{2}-c_{\mp}^{2})}.$		(56)

In the uniform case the only non-zero quantity is $\xi_{\rho}^{+}\to\xi_{\kappa}$, with $\xi_{\kappa}=\hbar/mc_{\kappa}$ being the usual expression for the healing length of a Bose-Einstein condensate.

At zero temperature and for $\omega>0$, the response function relates to the dynamic structure factor as $\chi_{o}^{\prime\prime}(q,\omega)=\pi S_{o}(q,\omega)$ (e.g. see [39]), and from this we can obtain the static structure factors $NS_{o}(q)=\int d\omega S_{o}(q,\omega)$. The static structure factors can be decomposed as $S_{o}(q)=\sum_{\nu}S_{o,\nu}(q)$ where $S_{o,\nu}(q)$ is the contribution from band $\nu$, which were defined in Eqs. (50)-(52). From hydrodynamic result (53) we obtain

$\displaystyle\lim_{q\to 0}S_{o,\nu}(q)=\frac{1}{\pi N}\int_{0}^{\infty}d\omega\,\chi_{o}^{\prime\prime}(q,\omega)=\frac{1}{2}q\xi^{\nu}_{o}.$

(57)

In Figs. 6(e)-(g) we compare the $S_{o,\nu}(q)$ functions obtained from the BdG calculations Eq. (50)-(52) to the hydrodynamic result for a dipolar supersolid. All the density fluctuations vanish proportional to $q$ in the $q\to 0$ limit, and the hydrodynamic result is seen to agree with the BdG results in this limit. A survey over a broader parameter regime is shown in Figs. 7 and 8 for the dipolar and soft-core systems, respectively. In these plots we only consider the hydrodynamic character of the fluctuations and compare $\xi_{o}^{\nu}$ to the $q\to 0$ limit of $2S_{o,\nu}(q)/q$ calculated from the excitations⁴⁴4We use a BdG excitation with $q\sim 10^{-2}/a$ to extract the limiting behavior..

Our results show that strain and defect density fluctuations are reasonably strong in the $\nu=-$ band of the dipolar supersolid, but these contributions cancel out to suppress the total density functions [Fig. 7(a)]. In contrast for the same band of the soft-core supersolid the lattice strain density fluctuations are heavily suppressed [Fig. 8(a)]. As noted in Sec. IV, these two systems differ in the relative size of their characteristic energies $mc_{\kappa}^{2}$ and $mc_{u}^{2}$. For the dipolar system the bulk compressibility dominates and it is favorable to reduce the total density fluctuations. In contrast for the soft-core system the lattice compressibility dominates and it is favourable to reduce lattice motion. Interestingly $\nu=+$ band dominates the lower band for total density fluctuations in the dipolar supersolid, whereas for the soft-core case the $\nu=-$ band dominates close to the transition.

For the $\nu=+$ band the role of defect fluctuations in the modulated state is more important for the soft-core system close to the transition. However, for both models deep into the crystalline regime, the defect contribution is negligible and we find $\xi_{u}^{+}\approx\xi_{\rho}^{+}$. This is consistent with this band becoming a purely crystal excitation.

Finally, we note that the response functions relate to important sum rules. Notably, the compressibility sum rule

$\displaystyle\int_{-\infty}^{\infty}\frac{d\omega}{\pi}\frac{\chi_{\rho}^{\prime\prime}(q,\omega)}{\omega}=\frac{N\alpha_{uu}}{\rho(\alpha_{\rho\rho}\alpha_{uu}-\alpha_{\rho u}^{2})}=N\hbar\rho\kappa,$

(58)

[cf. Eq.(25)] and the $f$-sum rule

$\displaystyle\int_{-\infty}^{\infty}\frac{d\omega}{\pi}\omega\chi_{\rho}^{\prime\prime}(q,\omega)=N\frac{\hbar q^{2}}{m}.$

(59)