Nonlocal field theory of quasiparticle scattering in dipolar Bose-Einstein condensates

We consider the propagation of quasiparticle excitations in a dipolar Bose-Einstein condensate, and derive a nonlocal field theory of quasiparticle scattering at a stepwise inhomogeneity of the sound speed, obtained by tuning the contact coupling part of the interaction on one side of the barrier. To solve this problem $ab$ $initio$, i.e., without prior assumptions on the form of the solutions, we reformulate the dipolar Bogoliubov-de Gennes equation as a singular integral equation. The latter is of a $novel$ $hypersingular$ type, in having a kernel which is hypersingular at only two isolated points. Deriving its solution, we show that the integral equation reveals a continuum of evanescent channels at the sound barrier which is absent for a purely contact-interaction condensate. We furthermore demonstrate that by performing a discrete approximation for the kernel, one achieves an excellent solution accuracy for already a moderate number of discretization steps. Finally, we show that the non-monotonic nature of the system dispersion, corresponding to the emergence of a roton minimum in the excitation spectrum, results in peculiar features of the transmission and reflection at the sound barrier which are nonexistent for contact interactions.


Introduction
How do the excitations of a quantum field above a given vacuum, the quasiparticles, propagate in a system which has nonlocal and anisotropic particle interactions?This seemingly simple question can be connected in fact to a plethora of related issues in different branches of physics and mathematics because of the very nature of the particle interactions.Indeed, in early quantum field theory, nonlocal alternatives for field theories were sought by the scientific community within the program of "handling divergences" [1,2], and foundations of such nonlocal quantum field theories were extensively studied, with particular emphasis on S-matrix properties cf., e.g., Refs.[3][4][5][6].However, the success and simplicity in studying lowenergy phenomena afforded by the renormalization program of local quantum field theories eventually won as the primary paradigm.Nevertheless, instances of such nonlocal field theories appear, e.g., as necessary tools in investigating electromagnetic phenomena in material media [7][8][9][10], for studying trans-Planckian physics [11] and the universality of Hawking radiation [12], and to assess the influence of boundaries [13][14][15].
Among the signatures of the dipolar interaction of particular importance for our analysis is the existence of rotonic excitations [29], which are caused by the anisotropy of the interaction, which is partly positive and partly negative.Roton modes, in particular, occur when the dipolar interaction dominates the interaction at high enough densities of the atoms or molecules.They play a pivotal role in the description of the dynamical instability emerging in the dimensional crossover from dynamically stable quasi-1D [30] or quasi-2D [31] condensates to 3D dipoledominated BECs, which are always dynamically unstable.We focus in what follows on quasi-1D trapping.
The appearance of a roton minimum in the dimensional crossover signals a marked departure of the standard Bogoliubov dispersion relation from its contact interaction form, which is what is obtained in a local field theory.Together with the associated maxon maximum, it corresponds to a non-monotonic dispersion.Indeed, examples of intriguing effects related to rotonic excitations include the enhancement of many-body entanglement [32], of density oscillations [33], and the occurrence of roton confinement [34].On the experimental side, rotons in elongated BECs have been observed, e.g., in [35][36][37].
We are however not aware of a solution of the Bogoliubov-de Gennes (BdG) equation describing quasiparticle propagation in dipolar Bose-Einstein condensates in an inhomogeneous setup, e.g., presenting an interface between regions of distinct quasiparticle spectra.Finding such solutions is key to provide a general answer on the apparently basic question posed at the beginning of this Introduction, and we provide in the below an ab initio answer, in which we put the inhomogeneous dipolar BdG equation of a quasi-1D gas into the form of an equivalent singular integral equation, and solve this equation.To the best of our knowledge, this singular integral equation is novel, in that it provides an extension of the well known Cauchy-type singular kernels [38].Specifically, the integral kernel we obtain is a combination of Cauchy-type kernels almost everywhere, with the exception of two isolated points where the singularity is stronger and the kernel becomes hypersingular.Our case is however different from established textbook examples of hypersingular kernels [39], where the set of singular points has nonzero measure.From a more practical perspective, no universally reliable numerical method exists for solving singular integral equations, and each case must be treated differently in order to avoid numerical instabilities [38].We discuss a method suitable for the BdG equation in its hypersingular integral form.We also provide a discretized version of the singular integral equation for the inhomogeneous dipolar BdG equation, and demonstrate its excellent performance for already a moderate number of discretization steps.
To give some intuition why the solution of this problem is nontrivial, note that within the instantaneous approximation for the dipolar interactions, a signal sent towards the barrier will interact with it before and after the signal has reached it.This is in striking contrast to the standard contact interaction case, where the signal interacts with the barrier only locally.Therefore, we expect nontrivial scattering phenomena to emerge.As we will show, these nontrivial phenomena are even more pronounced when roton excitations are involved due to the then increased number of the types of elementary excitations present in the system.Indeed, to assume an inhomogeneous configuration (e.g., a gas containing a sound barrier), as we will show, greatly increases the mathematical complexity of the perturbations in such systems, which in general forbids a fully analytical treatment.
Below, we reveal in detail how quasiparticles propagate in systems containing a sound barrier.Our results represent a major step for constructing a complete nonlocal field theory of dipolar BECs.The particular model we study, which encapsulates all required features, is a trapped BEC at rest with aligned magnetic or electric dipoles, which provides a sound barrier constructed by tuning locally the contact interaction between its particles, which is then separating the system into two regions with distinct sound velocities.Our goal is to make solutions to the nonlocal dipolar BdG equations as analytically amenable as possible.We show how the solutions we find can be used to build the S-matrix, which in the context of wave scattering comprises the reflection and transmission coefficients in such a way that unitarity is manifest.We shall see that when the dipolar interactions are present and the roton minimum exists, the increased number of the types of elementary excitations present in the system implies that the dimension of this matrix is larger than the 2 × 2 S-matrix for the case of contact-only interactions.
We shall discuss two methods of solving the BdG equation, one based on an approximate model, and the second one given in terms of special functions solutions to the novel class of singular integral equations we put forth.The method based on approximating the model has the advantage of allowing for analytic solutions, while the singular integral equation is treated numerically.Our results show that whenever the barrier exists, the dipolar interactions give rise to a continuum of evanescent channels bound to the barrier, which potentially play a role in near-boundary physics like recently explored in the context of the analogue gravity of sonic black holes [40].Furthermore, novel characteristic features include a decrease in the barrier's transmittance when the roton minimum is about to form and the barrier's complete transmittance/reflectance for particular signals when the roton minimum exists even in the limit of "weak" barriers.These findings represent a remarkable departure from the homogeneous system (no barrier), where one may naively expect to see a continuous dependence of sound propagation on barrier height.Yet, complete transmittance/reflectance is observed even for vanishing barriers, near the roton and maxon frequencies, in marked disagreement with the continuous dependence obtained in contact-interaction condensates.In summary, our presentation is organized in three parts as follows.
-The complexity of the scattering problem for our sound barrier model can be traced back to the analytical properties of the dipolar interactions in Fourier space.In section 2 we present the interaction kernel in Fourier space for a quasi-1D dipolar condensate and show that this kernel is always non-analytical.
-In section 3 we present our condensate configuration containing a sound barrier for its sound waves, and solve the scattering problem.Specifically, in subsection 3.1 we write down the BdG equation which models the propagation of small disturbances over our condensate.In subsection 3.2 we classify all the quasiparticle excitations of our condensate, from which the scattering problem is formulated in terms of waves sent towards the barrier and the associated reflected and transmitted parts.Moreover, subsections 3.3 and 3.4 contain the methods we employ to solve the scattering problem.
-Lastly, we discuss in sections 4 and 5 the physical implications of the quasiparticles found in section 3. We show that the scattering process is unitary, and determine in detail how sound waves are scattered by interfaces in dipolar condensate.
Previous studies have dealt with phonon scattering and the associated S-matrix for contact interactions, e.g. in the context of acoustic Hawking radiation [41].Yet, to the best of our knowledge we present the first complete ab initio nonlocal field theory of quasiparticle scattering at an inhomogeneity in the presence of a Bose-Einstein condensate on top of which the quasiparticles reside.We reveal, in particular, the impact of the anisotropy of interactions and the existence of a roton minimum on the scattering matrix.While a recent study explored the scattering properties of quasiparticles in polar dielectrics [42], our results are more general, do not assume any a priori knowledge of boundary conditions imposed by the dipolar metamaterial geometry and constitution and incorporate, in distinction to [42], the existence of a condensate.Considering a dipolar BEC with a stepwise discontinuous contact interaction, the structure of the S-matrix is derived from first principles.Therefore, our study has, as a further application, the potential to describe metamaterials built from dipolar BECs, by establishing a clear recipe of how to predict scattering phenomena in such systems.It thus paves the way towards a plethora of applications obtained by generalizations of our model.For instance, our results can be readily applied for an inhomogeneous extension of the recent experiment reported in [43].In this work, the crossover regime of a dipolar condensate to a supersolid and isolated droplet regimes was obtained by tuning the contact interaction, and studied using Bragg scattering of high energy excitations.By tuning the contact interaction locally, our model predicts the system response at any energy scale.

Dipolar interactions in quasi-1D condensates 2.1 Interaction kernel after dimensional reduction
As explained in the introduction, dipolar condensates can be stabilized by trapping potentials, which avoid the head-to-tail instability of bulk dipolar condensates by dimensional reduction.In this section we discuss the trapping mechanism we adopt, for which our condensate behaves as a quasi-1D stable condensate.We start with a strongly elongated dipolar condensate with its dipoles oriented along a given direction d (|d| = 1), such that its particles interact via the long-range instantaneous interaction energy in terms of the order parameter Φ and dipolar interaction strength C dd . 1 The interaction kernel U d is given by Let us assume the system is subjected to a strong radially symmetric trapping potential in such a way that the order parameter separation ansatz Φ(t, x) = φ ⊥ (|x ⊥ |)φ(t, x) holds, where φ ⊥ is normalized as d 2 x ⊥ |φ ⊥ | 2 = 1 , assuming the geometry presented in Fig. 1.For the particular case of dipolar interactions, under the assumed radially symmetric trapping, it was shown [44] that the only contribution from the interaction kernel in Eq. ( 1) is given by the Fourier transform where θ is the angle between d and the x axis, see for an illustration Fig. 1.It then follows from the order parameter separation ansatz that (∆x 1 C dd = µ 0 d 2 m for magnetic and C dd = d 2 e /ε 0 for electric dipoles, with dipole moments d m and d e , and where µ 0 and ε 0 are permeability and permittivity of the vacuum, respectively.
where G is defined via its Fourier transform G as2 and with ñ⊥ Here, ℓ ⊥ denotes the typical length scale of the transverse trapping.Moreover, we have set as the effective quasi-1D dipole coupling We note at this point that g d > 0 is required for the system to be stable in the thermodynamic limit and for vanishing contact interaction.We also observe that the advantages of presenting the dipolar kernel in Fourier space [Eq.(5)] rather than in configuration space are two-fold.Indeed, in Fourier space the double integral appearing in the energy functional (4) becomes a single integral (by the convolution theorem) which tends to simplify the analysis in particular when the condensate density is constant.Moreover, as we aim to study wave scattering at sound barriers, it is necessary to work with the interaction kernel in Fourier space, i.e., expressed in such a way that the role played by wave vectors become manifest.
For the (commonly employed) particular case of a strong transverse harmonic trapping, one has the Gaussian approximation , where then ℓ ⊥ is the harmonic oscillator length).For harmonic trapping, G(η) = (η 2 /2) exp(η 2 /2)E 1 (η 2 /2), E 1 being the first exponential integral function [44].In our work one particular property of this function plays an important role: it has a discontinuity branch on the imaginary axis of the complex k plane, which (for the Gaussian transverse profile) comes from the function E 1 [45].This feature increases the mathematical complexity of the condensate perturbations when some form of sound barrier exists in comparison to the case of contact-only interactions, and before we proceed to the model, let us pinpoint the origin of such a discontinuity and how it is related to the reduction to the quasi-1D regime.It is, in particular, not a feature of the transverse harmonic trapping per se, but occurs generically for any radial trapping, e.g., also for cylindrical box traps.

Analyticity of the kernel
Assuming an analytical interaction kernel (in Fourier space) is a simplifying hypothesis in nonlocal field theories [12].
for any function f , which leads us to question whether this property is fulfilled by our G.We note, in particular, that if Eq. ( 7) holds, the effect of the long-range interactions on plane waves is "local," and therefore the latter does not prevent the existence of plane wave solutions to the field equations.For the particular case of a dipolar interaction, inspection of Eq. ( 5) reveals the analytical structure of G in the complex plane.It has a discontinuity branch on the imaginary axis as can be seen from the application of the Sokhotski-Plemelj identity 1/(q ± iε) = 1/q ∓ iπδ(q) [46] as k x approaches the imaginary axis.Indeed, if iq for real q is any point on the imaginary axis, then straightforward manipulations lead to the jump magnitude measured by ∆ G(iq The above equation highlights the advantage of writing the interaction kernel in the integral form of Eq. ( 5), as it shows that the branch of G exists for any shape of radial trapping, which in turn determines the discontinuity branch jump through the form factor X in Eq. ( 8).For the Gaussian profile, the latter reads X(q 2 ) = exp(−q 2 /2), while for a cylindrical box trap, for which |φ ⊥ (|x ⊥ |)| 2 = 1/(πℓ 2 ⊥ ) for |x ⊥ | < ℓ ⊥ and zero otherwise, we find X(q 2 ) = 2J 2 1 (|q|)/|q| 2 , where J 1 is a Bessel function [45].
By tracing back from Eq. ( 5) to Eq. ( 3), we see that the existence of this discontinuity branch comes from the poles of the dipolar interaction in Fourier space, and Eq. ( 8) reflects the fact that different radial wavevectors add up to form the quasi-1D system.This is manifest in Eq. ( 3), where the pole at k 2 = k 2 x + k 2 ⊥ = 0 (in Fourier space) is evident.This gives rise to two first order poles at k x = ±ik ⊥ (manifest in Eq. ( 5)), which upon integration produces the discontinuity branch.In conclusion, both the dipolar interaction and the dimensional reduction combine to give rise to the discontinuity branch.
The relevance of Eq. ( 8) to our model is that it greatly modifies the structure and complexity of the perturbations in the system when a sound barrier exists, mainly because the commonly assumed Eq. ( 7) does not hold.Moreover, while noting that we are ultimately interested in the case where the kernel describes dipolar interactions, these same conclusions also hold for a gas of charged bosons (e.g., a Cooper pair gas in a superconductor) whose pairwise interaction follows the Coulomb law.The latter, however, which represents isotropic interactions, does not give rise to a roton minimum.Finally, the particular details of the trapping mechanism enter the analysis only through |φ ⊥ (|x ⊥ |)| 2 , which we assume henceforth to be given by the Gaussian approximation.We shall also omit the subscript x from the momentum k x along the weakly confining direction and denote it as k in what follows.
3 Quasiparticles in the presence of a sound barrier

Formulation of the Bogoliubov-de Gennes problem
We consider a stationary background condensate at rest given by φ = n exp(−iµt) (ħ h = 1), with particle density n and chemical potential µ.In order to model a sound barrier for the phonons in this system, we also allow the particles to interact via the Feshbach-tunable contact term H c = g c dx|φ| 4 /2, for an almost everywhere constant g c with a steplike discontinuity at x = 0.Then, as the local sound velocity is defined as c = n(g c + g d ) (setting the mass of the dipolar atoms or molecules m = 1), we see that this setup corresponds to a system in which the sound velocity has a sudden jump -the sound barrier -at x = 0.Alternatively, from the sound velocity definition, we see that a barrier is also created for an inhomogeneous condensate with a particle density n jump at x = 0, for fixed g c .Bearing in mind that this simplified physical system already requires a complex mathematical treatment, we shall assume at once that the region for x < 0 is dipole-dipole dominated g c /g d ∼ 0 and for x > 0, we have g c /g d > 0.Moreover, we shall assume for the sake of simplicity that the region x > 0 is such that g c prevents the formation of rotonic excitations.The discussion that follows can be easily extended to include the case in which rotonic excitations exist on both sides of the barrier, and also to sound barriers created by a density jump.
Therefore, the dynamics of the condensate is ruled by the system total energy where H d is defined in Eq. ( 4), and U is the 1D external potential, that leads to the 1D nonlocal Gross-Pitaevskii equation where denotes the convolution.In deducing Eq. ( 10), we used that G(∆x) = G(−∆x), as follows from the 3D dipolar kernel (2).Substituting our ansatz φ = n exp(−iµt) in Eq. ( 10) fixes the external potential in terms of the particle density: We note that the term (∂ 2 x n)/2 n vanishes for the case of a homogeneous condensate, and equals a "delta derivative" potential for the case of a sudden varying n at x = 0.The latter was recently explored by [40] in the context of analogue gravity in BECs.
Small disturbances in a stationary condensate are modeled by the Bogoliubov expansion φ = exp(−iµt)( n + ψ), where |ψ| 2 ≪ n, and ψ is a solution of the BdG equation, obtained by linearizing Eq. ( 10): In what follows, we take n to be constant, in which case the BdG equation simplifies to We emphasize, however, that the results below can be straightforwardly extended to the more general case of both n and g c changing at the barrier.We scale from now on lengths with ξ d = 1/ng d , wavevectors with 1/ξ d , and frequencies with 1/ξ 2 d , assuming thereby that g d is always rendered finite and positive.When thus fixing the scale ξ d , it should be kept in mind that g d depends on both the dipole orientation angle θ and the transverse trapping scale ℓ ⊥ via Eq.(6).
Our goal in this work is to study the solutions of Eq. ( 13).They are more easily found in terms of the Nambu field Ψ = (ψ, ψ * ) t , as demonstrated in detail for our type of system in [40].Note that because of stationarity, the quasiparticle modes still assume the general form Ψ(t, x) = exp(−iωt)Ψ ω (x), and Ψ ω satisfies where σ i , i = 1, 2, 3 are the usual Pauli matrices, with σ 4 = + σ 1 .If Ψ ω is a solution for Eq. ( 14), then σ 1 Ψ * ω is also a solution with −ω * .Accordingly, we might focus on the field modes with ω > 0. However, the analytical properties of G prevent Ψ ω from being a finite combination of exponential functions when g c is discontinuous (see for details Appendix A).
Far from the interface the solutions simplify to (a combination) of plane waves of the form Φ k exp(ik x) for constant Φ k , where Here, the dimensionless parameter β = ℓ ⊥ /ξ d = ℓ 2 ⊥ ng d (θ , ℓ ⊥ ) (reinstating here ξ d for clarity) measures the extent to which we are in the quasi-1D regime (in the dipole-dominated case).The proper quasi-1D limit, with all perpendicular motion frozen out, is achieved when β → 0.
The corresponding Bogoliubov dispersion relation is conveniently written as f ω (k) = 0, where as we show in Appendices A and B. We shall present below two routes for obtaining the solutions of Eq. ( 14).The solution for this equation is demonstrated in subsection 3.3, and a model approximation in which we discretize the integral of Eq. ( 5) is put forth in subsection 3.4.

Classification of quasiparticles
We can enumerate all the possible solutions of Eq. ( 14) as presented in Fig. 2 in which each plane wave propagating towards the barrier corresponds to a distinct quasiparticle.
From Fig. 2, we see the characteristic roton minimum formation in the dispersion relation caused by the dipole interactions [29][30][31]44,47], which singles out the spectrum subset defined by Ω (r) < ω < Ω (m) .For the sake of organization, we shall denote by k (respectively p) the wavevector solutions at x < 0 (respectively x > 0).Note that, in contrast to the contactonly interaction case, in each region the rotonic dispersion relation always admits 6 possible wavevector solutions, and for the propagating ones, each corresponding group velocity sign (graph slope) indicates if the solution represents plane waves traveling towards or away from the interface at x = 0. Accordingly, if ω / ∈ (Ω (r) , Ω (m) ), we have only one solution at x ≪ 0 (respectively x ≫ 0) propagating towards the interface, denoted by k in (respectively p in ), and one solution propagating away from it, k 1 (respectively p 1 ).Moreover, we also find at each side of the boundary evanescent channels, i.e., channels that are exponentially suppressed far from the barrier, denoted by k 2 , k 3 , p 2 , and p 3 .These channels have Im k i < 0 and Im p i > 0, i = 2, 3.
The complementary case -ω ∈ (Ω (r) , Ω (m) ) -is the physically richer one.We notice from Fig. 2 that all the six solutions at the left hand side of the barrier represent propagating waves, with three of them, k in1 < k in2 < k in3 propagating towards the barrier.We shall label the channels propagating away from the barrier as k 1 < k 2 < k 3 .To summarize the above discussion, we depict in Fig. 2 lower panel the schematics of all solutions of Eq, ( 14), from which we can state the main result of our work -the solution of the scattering problem -, as follows: For each wave sent towards the barrier at x = 0, we show how to determine the intensity of the reflected and transmitted waves.
The distinct behavior for dipolar interactions when a roton minimum is present comes from the shaded region on the LHS of the top panel Fig. 2, where three modes k in1 , k in2 , and k in3 in the band between Ω (r) and Ω (m) , the roton and maxon frequencies, respectively, can propagate towards the barrier (cf.lower panel on the left).This is in marked distinction to the contact-dominated case on the right (x > 0) of the barrier.
We note that the existence of complex wave vector solutions to the dispersion relation is not a peculiarity of our dipolar condensate as they are also present in inhomogeneous contactdominated configurations [40].Furthermore, evanescent channels admit no physical interpretation in term of quasiparticles, i.e., they, alone, do not represent solutions to the BdG equation per se.The same, however, is not necessarily true for real wave vector solutions, for at least in configurations similar to the one we considered in which an asymptotic regime exists (far from the barrier), propagating channels can be directly linked to the quasiparticles of the asymptotic system.The interpretation of each component appearing in a quasiparticle mode is thus r) , Ω (m) ) x ≪ 0, characterized by a fully dipole dominated interaction; Ω (r) and Ω (m) are the roton and maxon frequencies, respectively.Right: solutions for the region x ≫ 0, where a contribution from contact interaction is present.Lower panel: Asymptotics of all possible quasiparticle modes for the system under study.Zigzagged (respectively curved) lines represent propagating (respectively evanescent) channels.We note that k in in the right panel represents the three possible choices k in1 , k in2 , and k in3 .Also, arrows pointing to the right (respectively left) have V g > 0 (respectively V g < 0), where V g = dω/dk.Lower panel top part (respectively bottom part): schematics of the modes initiating at x < 0 (respectively x > 0).The scattering problem corresponds to a determination of how the transmission and reflection occurs for each possible incoming channel k in and p in .model dependent.In the system we considered the novelty regarding evanescent channels is their increased number, which from a physical perspective is an indication that near the barrier the transient regime caused by the barrier existence is distinct from the contact condensate analogue, as expected from the non-local character of the interactions.

Singular integral equation for the dipolar Bogoliubov-de Gennes problem
The details of how to solve Eq. ( 14) are presented thoroughly in Appendix B, and here we synthesize the main points.Following the asymptotic behavior just presented, we know that the quasiparticle modes are labeled by each incoming signal (k in or p in ), and far from the barrier they reduce to linear combinations of the channels shown in Fig. 2 lower panel.Our strategy here can be understood in terms of the latter channels as follows.For local, contact-only condensates, quasiparticle modes are built by combining the plane waves given by the Bogoliubov dispersion relation, i.e., local solutions to the BdG equation, and imposing matching conditions at the barrier.This procedure fails in the dipolar case, as the plane waves from the dispersion relation are not local solutions to the BdG equation.Accordingly, we shall develop a procedure to build, from the plane waves of Fig. 2 local solutions which can then be combined in the same fashion as for contact-interaction condensates.This "completion" procedure, performed via the BdG equation, gives rise to singular integral equations and associated special functions, but has the benefit of expressing the quasiparticle modes in a manner that leaves manifest their asymptotic properties while allowing for an easier numerical treatment than the direct solution of the BdG equation.
Because each field mode has a continuum of evanescent channels as imposed by the convolution in Eq. ( 14) (Appendix A), solutions can be found with the aid of the ansatz where the k ′ s, p ′ s, φ k , and φ p are given by Eqs. ( 15) and ( 16).The quantities ζ k and ζ p are matrix-valued functions, given by and with Π(q) = q 2 /2 − ωσ 3 σ 4 /(q − q − )(q − q + ).Furthermore, the functions Λ k,q and Λ p,q are solutions of the novel singular integral equation where the integrals are Cauchy principal values, ∆ G(iq) was defined in Eq. ( 8), and the function h(q) is defined by with , the average of G along the discontinuity branch.Finally, the real parameters q − < 0 < q + are the two simple zeros of h: h(q ± ) = 0. We stress that the ansatz (17) was constructed in such a way that each ζ k (x)Φ k , ζ p (x)Φ p is a local solution to the BdG equation, i.e., they are built to satisfy Eq. ( 14) at all points except at the barrier (x = 0).Thus, the continuum of evanescent channels in Eqs.(18) and (19) gives a succinct representation of the fact that the Bogoliubov channels of Fig. 2 fail from being local solutions.Moreover, a few features of the ansatz (17) are revealed by direct inspection of Eq. (20).In general, for analytic interaction kernels, one has ∆ G(iq) := 0, which implies Λ k,q , Λ p,q := 0 for all k's and p's, and the ansatz (17), through Eqs. ( 18) and ( 19), reduce to a finite combination of the exponentials given by the Bogoliubov dispersion relation (16).This is, in particular, the case for contact-type interactions.Furthermore, the decay of Λ k,q , Λ p,q as functions of q depends on the radial trap through ∆ G(iq).For the Gaussian approximation adopted here, Eq. ( 20) shows that Λ k,q , Λ p,q are exponentially suppressed for large q, whereas for the box trap profile, Λ k,q , Λ p,q decay with a power law [cf.Eq. ( 8) and the discussion after it].
In order to find a solution in the form (17), the scattering coefficients S k and S p must be uniquely fixed up to an overall phase and a normalization constant.We note first that substitution of this ansatz into the BdG equation implies (after a lengthy calculation) that a solution in this form exists only if are satisfied, which thus fixes two of the six scattering coefficients.These equations are necessary conditions for Eq. ( 20) to hold (details in Appendix B).The remaining boundary conditions are fixed by standard wave mechanics techniques applied to Eq. ( 14): Ψ ω and ∂ x Ψ ω are continuous at the barrier.
A few more general comments about the ansatz ( 17) are in order.First of all, the solution just found is build from the solutions of the singular integral equation (20), which in fact can be even more difficult to solve than the BdG equation itself.However, this equation can be solved numerically to any precision by means of cubic splines [48].Although this requires a considerable numerical effort, the latter method has an advantage over e.g.collocation schemes [49], as it does not rely on prior assumptions on the form of the solutions (for which explicit examples are contained in Appendix D).The sensitivity towards choosing the appropriate numerical scheme additionally serves to illustrate the mathematical challenge posed by nonlocal field theories.
Also, we see from the definitions in Eqs. ( 18) and ( 19) that in addition to the evanescent channels coming from the dispersion relation (16), the nonlocal dipolar interactions give rise to a continuum of evanescent channels when a sound barrier exists.This can be traced back to the analytical properties of G discussed in Sec. 2, as we see from there that if ∆ G ̸ = 0 (see Eq. ( 8)), then Λ k,q given by Eq. ( 20) is nonvanishing.Furthermore, we know that when the barrier is absent (g c = 0), the solutions to the BdG equation are single propagating exponentials [31,32,44].The numerical implementation of the general solution (17) for this particular case recovers this fact, as we verified when building the numerical solutions presented in Appendix D.
We conclude this presentation of solutions to the BdG equation with a further remark regarding the level of mathematical complexity in this system when compared with sound propagation in contact-only interacting condensates.In the latter, the analogue situation of a sound barrier modeled by a sudden g c jump over a homogeneous condensate results in the quasiparticle modes being combinations of only a few simple exponentials [40].This inspiration, gleaned from contact interactions, motivates us to seek for a different strategy in solving for the scattering coefficients, which consists in an approximation to the kernel (5), such that the corresponding quasiparticle modes are also expressed as a finite combination of exponentials.

Approximate solution
The strategy we follow in this subsection consists in finding an approximate discrete version of the integral equation (20), in such a way that the exact quasiparticles can be recovered via Figure 3: Approximations for the Fourier space kernel G in (5) for two sets of {N , ∆q}.For N = 100 the approximation is indistinguishable on the scale of the figure from the exact result.some limiting process.It turns out that a convenient way of achieving such approximation is to substitute the dipolar kernel integral (5) by a finite sum.Specifically, we take G → G, where and the two parameters ∆q > 0 and the integer N are free.We note that N → ∞ and ∆q → 0 reproduces exactly Eq. ( 5), as depicted in Fig. 3. Furthermore, 2N is the number of (simple) poles in this function, which are located at k = j(i∆q), for −N ≤ j ≤ N , and, naturally, ∆q is the distance between two consecutive poles.The examples of Fig. 3 show the accuracy potential of our approximation.In particular, we see from Fig. 3 that for N = 10 and ∆q = 1/3, i.e., the approximating G containing only 20 poles, a reasonable agreement is already obtained.Let us now build the solutions for Ψ ω for the model of Eq. ( 24).The general asymptotic properties of the possible quasiparticles modes are the same as presented in Fig. 2. In the absence of G, the dispersion relation in Eq. ( 16) is a degree four polynomial equation, whose solutions are explicitly found for each ω > 0, whereas for G given in Eq. ( 24), in addition to these four solutions, another 2N solutions are present, one for each pole in G. Furthermore, in view of the numerator of Eq. ( 24), no dispersion relation solution coincides with the poles of G in the complex k plane.
For the approximate model under study, the solutions for the dispersion relation can be grouped in the single zero level set f ω (k) = 0, as done in Fig. 4. As explained, in each region of the condensate and for each ω > 0, Eq. ( 16) has 4 + 2N solutions, with the real solutions depicted in Fig. 2, and the remaining being necessarily complex.We recall that when G is the exact one given in Eq. ( 5), the wave function presents a continuum of evanescent modes on the imaginary axis.We thus readily see that the approximation implemented here also discretizes this continuum, as indicated by Fig. 4, where in addition to the six channels existing when G is used, we have the extra blue points on the imaginary axis.Furthermore, as N → ∞ and ∆q → 0, these poles clump together, thus recovering the continuum of the solutions (Eq.( 17)) and reinforcing the validity of our approximation.We label the evanescent channels at x < 0 by k j and at x > 0 by p j .Clearly, Im k j ≤ 0 and Im p j ≥ 0.
It turns out that these exponentials, the solutions to the dispersion relation Eq. ( 16) of the approximated model, can be used to construct solutions to the wave equation.Indeed, let us Figure 4: Zero level sets for the real and imaginary parts of Eq. ( 16).Blue circles are solutions to f ω = 0 for the fixed frequency ω = 0.4, i.e., they are the wave vectors solutions to the Bogoliubov dispersion relation in the complex plane.We have chosen N = 10 and ∆q = 1/3 in the discrete representation (24).k r and k i denote the real and imaginary parts of k, respectively.Left: (Right:) dispersion relation at x < 0 (x > 0).The solutions on the imaginary axis tend to a continuum of evanescent channels when N → ∞ and ∆q → 0. look for a solution with the ansatz whose relation to Eq. ( 17) is manifest: the integrals in Eqs. ( 18) and ( 19) are substituted by finite sums of evanescent channels (see Fig. 4).Each incoming signal at the barrier gives rise to a distinct solution, which is associated to reflected and evanescent channels.By carefully counting the number of unknown coefficients S k and S p , we find that for a given N , we have 2 + N S k 's and 2 + N S p 's, leading to a total 4 + 2N unknown coefficients for each field mode.If a solution exists in the form (25), then the coefficients are fixed by the BdG equation ( 14).We note, however, that the application of standard wave mechanics techniques to the field equation only produces 4 boundary conditions, namely, Ψ ω and ∂ x Ψ ω are continuous at the barrier.This is caused because the convolution in Eq. ( 14) is always continuous (see Appendix A).We conclude, therefore, that 2N boundary conditions appear to be "missing," and the solution to this puzzle is found from the form of Eq. ( 25).Indeed, the convolution σ 4 G * Ψ ω of Eq. ( 14) is conveniently written as where Ψω (k) is the Fourier transform of Ψ ω (cf.Appendix A), and the (simple) poles of Ψω are precisely the k's and p's appearing in Eq. ( 25).Thus, as the poles of G(β k) cannot coincide with the poles of Ψω (which are given by the dispersion relation ( 16)), the convolution (26) does not have evanescent channels at the poles of G(β k) if and only if σ 4 Ψω vanishes at these poles.Therefore, as G(β k) has 2N poles at the imaginary axis, a solution of the form (25) exists if and only if for −N ≤ j ≤ N , j ̸ = 0.This gives us an additional set of 2N boundary conditions, as required.This concludes the construction of the quasiparticle modes for the system under study.We emphasize at this point that the approximate solutions described in the above greatly reduce the numerical effort necessary to simulate any observable quantity in the inhomogeneous dipolar BEC.

S-matrix and unitarity
All the quasiparticle modes thus constructed admit the in/out state interpretation: For each field mode, the incoming channel, labeled by k in 's or p in , represents a signal sent towards the boundary at the asymptotic past which emerges at the asymptotic future propagating away from the barrier through the various reflected/transmitted channels.This can be read off directly from the solutions in Eqs. ( 17) and ( 25), as far away from the boundary each elementary excitation reduces to a sum of plane waves.The scattering coefficients involved in all of these processes can be conveniently grouped into an unitary matrix -the S-matrix -which is simply an expression for the conservation of field mode normalization throughout the system's causal development.When only contact interactions are present, this conservation is implied by the BdG equation to be ∂ x J ω = 0, where The important consequence of this equation is that the total flux in the system is conserved: J ω (∞)− J ω (−∞) = 0, giving rise to the S-matrix conveniently defined in terms of the norm-preserving condition where the sums are performed over the propagating channels only.It is to be stressed that this conservation law is not due to the particle number conservation of the full theory implied by its U(1) symmetry.It is a well known fact that the Bogoliubov expansion leads to a theory with a spontaneously broken U(1) symmetry [50].The conservation law we are exploring here is, then, the time-independence of the Bogoliubov scalar product, defined as 〈Ψ, Ψ ′ 〉 := dxΨ † σ 3 Ψ ′ , where Ψ = exp(−iωt)Ψ ω (x), Ψ ′ = exp(−iω ′ t)Ψ ω ′ (x) are any two field modes.Via the BdG equation ( 14), the condition ∂ t 〈Ψ, Ψ〉 = 0 can be shown to be equivalent to ∂ x J ω = 0 for all ω in the absence of dipolar interactions.However, when dipolar interactions are present, J ω is no longer conserved (see for the extended discussion in Appendix C): This should not be read as implying that there is no flux conservation, but that the quantity J ω is no longer a bona fide representation for the system total flux.We show in detail in Appendix C that for both families of solutions found in the previous section the dipolar analogue of Eq. ( 28) acquires the form The meaning of this relation is more readily grasped by constructing the S-matrix.To that end, we take the solution for Φ k to be normalized as where V g = dω/dk is the group velocity.
We have presented here the flux conservation in terms of the approximated model, which in turn implies its validity for the continuum as we take the limit of an infinite number of discretization steps.Also, we note that the particular form of the normalization in Eq. ( 31) is not relevant for the physics of the problem.However, this form is convenient for us because it implies for propagating channels.This form of normalization is one possible choice which bounds the scattering coefficients absolute value to unity, as we shall see now.Let us label the quasiparticle modes according to their incoming channel by adding a superscript to the scattering coefficients.For instance, for p in , we set S k → S p in k in the general solution (25).Moreover, we set to unity the intensity of the incoming channels, i.e., S p in In view of the schematics displayed in Fig. 2 lower panel, Eqs. ( 30) and ( 32) simply mean that the matrix S ω defined by for ω / ∈ (Ω (r) , Ω (m) ) and in case that ω ∈ (Ω (r) , Ω (m) ), is unitary, i.e., S † ω S ω = S ω S † ω = .This "unitarity" enables us to study quasiparticle scattering by the sound barrier as considered, because it ensures that no signal amplitude is lost during quasiparticle propagation and scattering.Finally, a different choice of normalization for Φ k clearly does not spoil the S-matrix as it is given by Eq. ( 30), but the simplified form S ω and its unitarity condition changes.
We finish this section with a disclaimer regarding the S-matrix nomenclature.Naturally, scattering matrices appear generically in distinct branches of physics, and as such it is important that the meaning of the notion of S-matrix be completely clear in each context.We cite, for instance, the two works [51,52] in which S-matrices appear within the same wave mechanics context as in our work.Yet, in [52] the very same mathematical object is not referred to as an S-matrix.We should thus stress that the S-matrix in our context is not the same object as the one occurring in particle scattering in quantum field theory, where it is just representing the Schrödinger equation written in terms of the evolution operator [53].Note that a proof of the unitarity of the S-matrix in such a quantum field theory context for a nonlocal family of scalar quantum field theories was provided in [6].

Transmittance and reflectance: Impact of dipolar interaction
As an application of the solutions constructed we now investigate how the sound barrier transmittance/reflectance is affected by the roton minimum.Once the scattering coefficients are known, we can study how waves sent towards the barrier get reflected and transmitted, the dipolar condensate analogue of the scattering of light rays at interfaces in classical optics.Specifically, far from the barrier (|x| ≫ 1) in general we obtain from Eq. ( 17) or ( 25) that i.e., only the propagating channels (real wave vectors) contribute.For the sake of illustration, if we consider the case where no roton is present, the quasiparticle of frequency ω and indexed by k in far from the barrier assumes the form This is the solution of the BdG equation in our dipolar condensate configuration.Indeed, for contact-only condensates the quasiparticles also assume the form in Eq. ( 36) far from inhomogeneities, and thus we conclude that the scattering coefficients encapsulate all new physical information coming from the dipolar interactions in comparison to contact interactions.We observe that Eq. ( 36) corresponds to a plane wave of "unity intensity" sent towards the barrier from the asymptotic left, namely the term proportional to exp(−iωt + ik in x), which then gives rise to a reflected wave [exp(−iωt | 2 and a transmitted wave Moreover, the unitarity of the process discussed in the previous section shows us that the constraint |S When the roton minimum is present, the same reasoning and interpretation holds, with the difference that more propagating channels might be present in the quasiparticles following the scheme of Fig. 2. We discuss in the below each of the possible cases separately. We use in this section only the approximate model, as it allows for an easier numerical simulation, and leave for the Appendix D some worked out examples of the singular integral BdG equation (20).We recall from the dispersion relation ( 16) that the influence of the dipolar interactions is measured by the coefficient β = ℓ ⊥ : As β increases, the roton minimum emerges and the system eventually becomes dynamically unstable; when β → 0, the system is stable under the proviso that g d > 0 .Therefore, β measures how deep into the quasi-1D regime the system has penetrated.
Furthermore, we note from Eqs. ( 5) and ( 24) that, because of the global factor k 2 , the long-range part of the dipolar interactions is suppressed in the quasi-1D limit β → 0, and the condensate behaves no different than one with only local contact interactions, which is modeled by g = g d for x < 0 and g = g d + g c for x > 0. Thus our model enables us to compare results with the non-dipolar case operating near or within the quasi-1D limit.For the sake of a clear representation, we now treat separately the cases with roton and no roton minimum.

No roton minimum
When rotonic excitations are not present in the system, for each frequency ω in the system spectrum there are two elementary excitations, corresponding to the signals sent towards the barrier at each of its sides, i.e., signals k in and p in .Accordingly, the S-matrix has the form (33) for each frequency subspace.As S ω is unitary, this means that both its row and column vectors form orthonormal bases, which in turn implies that the absolute value of one of its components fixes all the others.Therefore, by calculating |S k in p 1 | 2 -the transmitted intensity through the barrier from left to right -we also determine |S We show in Fig. 5 how the barrier transmittance varies with frequency for several values of β.The dotted-brown curve corresponds to β = 0, and it shows that the barrier with "height" showing how the roton minimum threshold is approached.The brown-dotted curve corresponds to the contact-only regime (β = 0), showing that this barrier is almost transparent for contact-only interactions.Yet, the bending of the dispersion relation leads to a noticeable decay of the barrier transparency even in the absence of a roton minimum, as shown by the continuous, dashed, and dot-dashed curves.Here, discretization parameters are N = 10 and ∆q = 1/3.4.For these parameters, the roton minimum formation threshold is approximately β ∼ 0.689.modeled by g c /g d = 0.2 is almost transparent for contact-only interactions.Yet, we note that a noticeable decay in the transmittance is predicted to occur for an extended range of signal frequencies even in the absence of a roton minimum when β is increased, which thus reveals how wave propagation in this system is sensitive to inhomogeneities.

With roton minimum
When the roton minimum exists, for frequencies not in the grey band shown in the left upper panel of Fig. 2, ω / ∈ (Ω (r) , Ω (m) ), there are only two quasiparticle excitations for each frequency, and the same reasoning used to interpret the case with no roton minimum present can be repeated.We present in Fig. 6 several simulations for this regime, which supplements our findings from Fig. (5) beyond the roton minimum formation, namely, the high energy sector of the theory is sensitive to the bending of the dispersion relation, a feature not present when only contact interactions exist for condensates at rest.Furthermore, within the band ω ∈ (Ω (r) , Ω (m) ) (the shaded area in Fig. 6) the degeneracy of each frequency subspace is 4 for the parameter choice we are investigating, corresponding to the four distinct quasiparticles that can be excited: k in1 , k in2 , k in3 , p in .We stress that this increased degeneracy has no counterpart in condensates whose particles interact only locally.It is instructive to analyze each quasiparticle separately. (r) , Ω (m) ), k in1 excitation We show in Fig. 7 the reflectance and transmittance coefficients of the quasiparticle branch indexed by k in1 , cf.Fig. 2. Reflectance and transmittance measure the fractions of the plane wave signal exp(−iωt + ik in1 x) coming towards the barrier from the left that get reflected and  (r) , Ω (m) ), which cannot be characterized by the coefficient |S k in p 1 | 2 alone.Here, N = 10, ∆q = 1/3.4and β = 0.76.This implies Ω (r) ∼ 0.24 and Ω (m) ∼ 0.45 for the roton and maxon frequencies, respectively.These curves represent supplementary data to the ones presented in Fig. 5: although those sound barriers have weak influence on the phonon sector of the theory, the high-energy sector is sensitive to the bending of the dispersion relation, a feature not present in the contact-only case.transmitted through the various available channels: k 1 , k 2 , k 3 , and p 1 .Therefore, the data in Fig. 7 show us that the barrier is mostly opaque for these waves, as the transmittance in this case satisfies |S k in1 p 1 | 2 ≪ 1.Furthermore, the signal is integrally reflected through the channel k 1 (respectively k 2 ) for frequencies close to the roton (respectively maxon) frequency, and this feature survives even for very small g c /g d = 0.001.It is noteworthy that in the barrier's absence, i.e., g c = 0, any signal of this form is clearly integrally transmitted, and thus the very existence of the barrier leads to the total reflection of k in1 quasiparticles with frequencies near Ω (r)  and Ω (m) . (r) , Ω (m) ), k in2 excitation The analysis of the remaining quasiparticle modes follows the same line of reasoning just presented.Figure 8 presents the scattering coefficients for the k in2 quasiparticle branch.

ω ∈ (Ω
We see from Fig. 8 that, in distinction to the k in1 quasiparticles, plane waves of the form exp(−iωt +ik in2 x) sent towards the barrier from the left are integrally transmitted for frequencies close to the roton minimum, and a transition to a completely opaque barrier is observed as the frequency tends to Ω (m) .We also observe that this opaqueness near maxon frequencies survives even for very small g c /g d = 0.001, in the same way as it happens for the k in1 excitation discussed in the preceding paragraph. (r) , Ω (m) ), k in3 and p in excitations

ω ∈ (Ω
The remaining quasiparticle branches, namely, k in3 and p in are analyzed in the same fashion.They correspond to the plane waves exp(−iωt+ik in3 x), exp(−iωt+i p in x) propagating toward the barrier from the left and right, respectively.We see from Fig. 9 that both waves experience Figure 7: Scattering coefficients for the quasiparticle branch indexed by k in1 in the frequency band ω ∈ (Ω (r) , Ω (m) ).Here, N = 10, ∆q = 1/3.4.This implies Ω (r) ∼ 0.24 and Ω (m) ∼ 0.45 for the roton and maxon frequencies, respectively.We note that the sum |S 1 for all ω, one of the properties of the S-matrix.For the quasiparticles, labeled by k in1 , the barrier is mostly opaque, with the incoming signal being reflected exclusively through the channel k 1 for ω → Ω (r) and through the channel k 2 as ω → Ω (m) .a partially transmitting barrier for frequencies near Ω (m) , in distinction to the k in1 and k in2 excitations discussed above that experience a completely opaque barrier at these frequencies.
Furthermore, for frequencies near the roton minimum, the barrier becomes completely opaque for the k in3 excitations and completely transparent for the p in branch.

Final remarks
The present study provides a systematic route to describe the scattering of quasiparticles in inhomogeneous dipolar Bose-Einstein condensates.To this end, we have studied perturbations in a quasi-1D dipolar condensate in which a sound interface exists, separating the condensate in two regions possessing different sound velocities.The perturbations were built via two distinct methods, one based on the direct simulation of the solutions and one based on a family of approximated models.Both methods explore the fact that the long-range dipolar interaction in Fourier space is not modeled by an analytic kernel, which gives rise to a distinct set of evanescent channels bound to the sound barrier.
As a particular application, we have shown how sound scattering occurs as a function of the perturbation frequency.An intricate pattern of reflectance/transmittance was shown to emerge due to the existence of rotonic excitations, which are due to the interaction anisotropy, which leads to a strong dependence of the barrier transmittance as function of the mode frequency and "polarization."In fact, if rotonic excitations exist, for frequencies within the roton/maxon band, the number of elementary excitations that can be scattered by the barrier is larger than two, as happens for the case where only contact interactions are present.Each of these correspond to a distinct system response to external excitations, and the barrier behaves fully/partially transparent or opaque depending on the excitations.Figure 8: Scattering coefficients for the quasiparticle branch indexed by k in2 in the frequency band ω ∈ (Ω (r) , Ω (m) ).The parameters are the same as in Fig. 7.We observe that these quasiparticles experience a transition from a completely transparent barrier for frequencies near Ω (r) to a completely opaque one, for frequencies near Ω (m) and with the reflection exclusive through k 3 .The channels k 1 and k 2 have negligible participation in the scattering process.
We started from a homogeneous quasi-1D dipolar condensate at rest, and asked how sound propagates in this system if a sound barrier is constructed at which the sound speed experiences a steplike change.As the system admits also rotonic (and maxonic) excitations when dipoleinteraction-dominated, we demonstrated that the S-matrix dimension is enlarged.It should be emphasized that, starting from these two basic ingredients we use-dipolar condensate at rest and a sound barrier -, that only through the exact knowledge of the quasiparticle modes we derived it is possible to unveil the peculiar results presented for the reflectance and transmittance of the increased number of modes available to the system.The properties of the S-matrix we obtain clearly distinguish dipole-interaction-dominated BECs from their their contact-interaction-dominated counterparts.
The importance of the present work thus consists in providing a recipe to build the quasiparticle spectrum in dipolar condensates with inhomogeneous sound velocity.We should also stress that our method is not restricted to condensates with homogeneous densities, as only a few modifications are necessary to study configurations in which the sound interface is caused by density jumps.Indeed, once the knowledge of how to build quasiparticle mode solutions for this type of inhomogeneous dipolar condensates is set, it is straightforward to apply the same technique a number of systems of interest.In particular, the application we have explored here is far from exhausting all the features inherent in this type of system, with examples including a roton almost touching the wavevector axis which could hybridize with low-energy phonons, and the existence of rotonic excitations on both sides of the barrier.As a natural extension of our results, nontrivial effects are expected if the sound barrier is taken to have a finite size.Furthermore, we expect that selective nature of quasiparticle scattering in dipolar BECs through the associated transmission and reflection coefficients will have significant impact on the properties of droplets, which represent naturally emerging interfaces in the condensate due to quantum fluctuations [23,[25][26][27][28].
Our analysis indicates how intricate quasiparticle scattering in inhomogeneous dipolar Figure 9: Scattering coefficients in the frequency band ω ∈ (Ω (r) , Ω (m) ).Left panel: quasiparticle branch indexed by k in3 .Right: quasiparticle branch indexed by p in .The parameters are the same as in Figs.7 and 8.Both families of quasiparticles experience a partially transmitting barrier for frequencies near Ω (m) , whereas for frequencies near the roton minimum, the quasiparticles k in3 (respectively p in ) experience a completely opaque (respectively transmitting) barrier.
BECs with continuous variations of density and coupling will be, where in general no simple solutions to the quasiparticle modes can be found.The extension of our technique to the case of continuous sound barriers therefore promises to reveal further features of quasiparticle scattering in dipolar BECs over those familiar from contact interaction condensates.The quasiparticle modes we constructed are already properly normalized according to the Bogoliubov scalar product.Therefore, quantization of the present nonlocal field theory is a formal step which follows straightforwardly from expanding the bosonic field operator into this complete set of modes.
We finally note that infinitely extended (quasi-)1D Bose-Einstein condensates, according to the Hohenberg theorem [54], do not exist, 3 due to the divergence of nonlocal phase fluctuations in the system.This imposes certain restrictions on the applicability of infinitely extended systems to model experimental realizations, as we expect the phonon part of the spectrum to be sensitive to finite size effects.This, however, should not interfere with the general properties of the scattering processes here investigated -which occur far from the long-wavelength phonon regime -as long as the condensate is still sufficiently large along the weakly confining direction, where the achievable length is subject to a position space version of the Hohenberg theorem [55].
Eq. (B.2 ) is true only if the term inside the outermost curly brackets vanishes, which then yields a Cauchy-type singular integral equation [49].
Equation (B.2 ) can be further simplified by relating Λ q to the scalar functions Λk ′ ,q and Λp,q occurring in Eqs. ( 18)-( 20) of the main text, With this relation, we find that where the function h is defined in Eq. ( 21) of the main text.
A quick inspection of Eq. (B.4 ) suggests that we could define each function Λ to be the solution of the integral equation inside the square brackets, e.g., in such a way that the BdG equation would be satisfied at all points except at x = 0.However, this procedure is "partially" wrong, and the reason for this is that we must understand the space of solutions for these integral equations first.Note that if we remove the integral part in Eq. (B.5 ), then the function Λk ′ ,q necessarily has two poles, at q = q ± , defined by the two zeros of the function h(q ± ) = 0.That happens because the operator which multiplies by h(q) is only invertible when restricted to the space of solutions that diverge at q ± .Because of this feature, we might look for solutions in the form Λk ′ ,q = Λ k ′ ,q (q − q − )(q − q + ) .(B.6) By plugging this ansatz back into Eq.(B.4 ), notice that the integral kernel becomes which is not defined when q → q ± unless Λ k ′ ,q ± = 0.This kind of divergence is well understood and we can trace it back to the groundbreaking work of Ugo Fano [56].Indeed, we have i.e., each double pole contains delta-contributions, leading to the divergence of the integral as q → q ± .Now, in view of Eqs.(B.6 ) and (B.8 ), Eq. (B.4 ) becomes − δ(q − q − )Λ k ′ ,q − − {q − ↔ q + } + q − − q + iq − k ′ + {k ′ ↔ p} = 0 .(B.9) We thus see that the term inside the first curly brackets in Eq. (B.9 ), apart from the delta functions, defines the integral equation ( 20) of the main text.Referring back to Eq. (B.2 ), we conclude that the ansatz (B.1 ) is a solution of the BdG equation ( 14) only if the two delta functions contributions in Eq. (B.9 ) at q = q ± vanish.This leads precisely to the conditions Eqs. ( 22) and ( 23) of the main text.We can then assert that the functions Λ k ′ ,q satisfy the integral equation ( 20), which can be solved numerically to any desired accuracy by means of cubic splines [48], see Appendix D.

C Flux conservation
In this appendix we shall present a proof for the flux conservation stated in Eq. (30).Due to the fact that G tends to G in the limit of a large number of discretization steps, we demonstrate flux conservation for the discrete approximation, from which conservation for the continuum model then follows.We shall start by taking two arbitrary solutions of the BdG equation ( 14), Ψ ω and Ψ ω ′ .By multiplying the equation Ψ ω ′ on the left by Ψ † ω , and the equation for Ψ † ω on the right by Ψ ω ′ , followed by their subtraction and integration with respect to x results in if the solutions are assumed to vanish at infinity.This is just the expression for the Bogoliubov scalar product in the space of classical solutions between any two such solutions.In particular, Eq. (C.1 ) holds as long as is fulfilled.As we show in the main text, the boundary conditions (27) for the approximated model mean exactly that we have G * Ψ ω = G(−iβ∂ x )Ψ ω at the system's solutions.With this, it is straightforward to perform the integral in Eq. (C.2 ): where k = k ω , k ′ = k ω ′ , and similarly for p, p ′ .Equation (C.3 ) includes the result, namely, Eq. ( 30) in the main text, which we are looking for.If we let ω ′ → ω, the off-diagonal terms in Eq. (C.3 ) alongside terms involving evanescent channels vanish in view of Eq. ( 15), and the remaining diagonal terms are precisely the flux conservation presented in Eq. ( 30) of the main text.

D S-matrix examples from solving Eq. (20)
In this appendix, selected explicit scattering coefficients of quasiparticle modes from solving the singular integral BdG equation ( 20) are worked out.As anticipated in the main text, solutions can be built numerically to any accuracy and precision, by calculating the various functions Λ k,q and Λ k,q and following the recipe developed in the main text Subsec.3.3.In Fig. 10, we present some examples for solutions of (20), which are used to construct the Tables 1, 2, and 3.The method we have implemented here is based on the spline technique of [48], which requires a higher computational effort, but on the other hand does not rely on prior assumptions on the form of the solutions besides continuity to hold for its effectiveness.
Figure 10: Examples for solutions obtained from the integral equation ( 20), using the spline method of [48].
We present in Tables 1, 2, and 3 the S-matrix coefficients of three sets of quasiparticle modes in this system.All the matrices have been verified to satisfy the unitarity condition (S † ω S ω = ) numerically.For the indicated parameters, the roton minimum exists and is given by Ω (r) ∼ 0.239701, with Ω (m) ∼ 0.438915.In particular for the quasiparticle modes in Table 3, we set ω = 0.24, which is close to the roton minimum frequency.These various coefficients give us information of how the barrier transmits and reflects for these three frequencies, as explained in the main text Subsec.5.2.

Figure 1 :
Figure1: Schematics of the elongated dipolar condensate under consideration.The system symmetry axis is here taken to be the x axis, and the dipoles are oriented by an external field along the direction d, which defines the angle θ as shown.

Figure 2 :
Figure 2: Upper panel: Bogoliubov dispersion relation.Left: solutions in the regionx ≪ 0, characterized by a fully dipole dominated interaction; Ω(r) and Ω (m) are the roton and maxon frequencies, respectively.Right: solutions for the region x ≫ 0, where a contribution from contact interaction is present.Lower panel: Asymptotics of all possible quasiparticle modes for the system under study.Zigzagged (respectively curved) lines represent propagating (respectively evanescent) channels.We note that k in in the right panel represents the three possible choices k in1 , k in2 , and k in3 .Also, arrows pointing to the right (respectively left) have V g > 0 (respectively V g < 0), where V g = dω/dk.Lower panel top part (respectively bottom part): schematics of the modes initiating at x < 0 (respectively x > 0).The scattering problem corresponds to a determination of how the transmission and reflection occurs for each possible incoming channel k in and p in .

Figure 5 :
Figure 5: Sound barrier transmittance measured by the coefficient |S k in p 1 | 2 for several values of β.Inset: positive branch of the dispersion relation for each curve, showing how the roton minimum threshold is approached.The brown-dotted curve corresponds to the contact-only regime (β = 0), showing that this barrier is almost transparent for contact-only interactions.Yet, the bending of the dispersion relation leads to a noticeable decay of the barrier transparency even in the absence of a roton minimum, as shown by the continuous, dashed, and dot-dashed curves.Here, discretization parameters are N = 10 and ∆q = 1/3.4.For these parameters, the roton minimum formation threshold is approximately β ∼ 0.689.

Table 3 :
S-matrix coefficients for a quasiparticle mode obtained from