Bound states of Dipolar Bosons in One-dimensional Systems

We consider one-dimensional tubes containing bosonic polar molecules. The long-range dipole-dipole interactions act both within a single tube and between different tubes. We consider arbitrary values of the externally aligned dipole moments with respect to the symmetry axis of the tubes. The few-body structures in this geometry are determined as function of polarization angles and dipole strength by using both essentially exact stochastic variational methods and the harmonic approximation. The main focus is on the three, four, and five-body problems in two or more tubes. Our results indicate that in the weakly-coupled limit the inter-tube interaction is similar to a zero-range term with a suitable rescaled strength. This allows us to address the corresponding many-body physics of the system by constructing a model where bound chains with one molecule in each tube are the effective degrees of freedom. This model can be mapped onto one-dimensional Hamiltonians for which exact solutions are known.

In this paper, we will address the few-body structures that can be expected in an array of several tubes. A schematic setup for three tubes is shown in figure 1. The orientation of the dipolar molecules in the tubes can be controlled by an external electric or magnetic field depending on whether magnetic atoms or heteronuclear molecules are used. This opens up a host of interesting phenomena since the potentials of two dipoles in a single tube and the potential between two dipoles in different tubes change magnitude and sign as one changes the angles shown in figure 1.
Here we are concerned with the important question of existence of N -body bound states in multi-tube geometries. In order to address this, we need to carefully consider the limits of small dipole moment and we consequently develop a perturbative approach to weakly bound states in one dimension that can handle the dipolar potential. At stronger coupling, we consider a harmonic approximation. Both are compared with numerical results from the stochastic variational method that allow us to test the accuracy of the analytical formalism. A comparison is also made to exact results for N -body bosonic systems interacting through zero-range attractive forces, and we ask whether dipolar systems in one dimension approach the exact results in the given limits and can in turn be effectively described by (properly renormalized) zero-range interactions. This is particularly important for many-body studies since zero-range interactions are not only convenient to work with but also often allow analytical results to be obtained. We also consider the two-body scattering dynamics of the dipoles and compare the results with the Schematic representation of the setup for the case of three equidistant 1D tubes along the x-direction in the x y-plane with a distance between adjacent tubes of d. The dipolar moment, D, points along the direction specified by the two angles φ and θ as defined in the figure. Dark filled circles indicate the dipolar atoms in the tubes. The tubes have a thickness given by l. In typical experimental setups the tubes also have an in-tube confinement. This is indicated in the figure by the shrinking of the tubes at both ends. This in-tube confinement will, however, be neglected in this work. zero-range results. Finally, we discuss the impact of our results for the many-body physics of dipolar bosons in 1D geometries.

Basic setup
We consider the setup depicted schematically in figure 1, i.e. an array of equidistant 1D tubes containing dipolar particles with dipole moments aligned by an external field that does not interfere with the tubular geometry. The dipolar particles are identical bosons with mass m and dipole moment D. The potential between two dipoles with coordinates (0, 0, 0) and (x, nd, 0), V (n, x), can be written in the form where d is the distance between two adjacent tubes and n is an integer such that nd is the intertube distance (n = 1 for nearest-neighbour (adjacent) tubes, n = 2 for next-nearest neighbours and so on). Here we have defined the dimensionless dipole coupling strength, where D is the dipole moment (absolute value of the vector D in figure 1). For n = 0, we obtain the intratube interaction in the limit where the tubes are strictly 1D. In experimental setups, arrays of 1D tubes are constructed by applying optical lattices to the dipolar gas [7,8]. In this case the tubes are not strictly 1D but will have some width along the transverse direction that is determined by the laser intensity. This can be translated into a Gaussian wave packet in the transverse direction that will be increasingly localized in space as the laser intensity increases. An effective interaction that takes this into account by integrating out a Gaussian wave packet can then be obtained [25], which will in the limit of zero Gaussian width reduce to the expression in (1). We will assume that the lattice is very strong so that the strict 1D expression above is valid for the interaction of particles in different tubes, which is accurate when the transverse width, l, is much smaller than the intertube distance, d. Corrections to this picture have been discussed in [24,25]. Below we will return to the question of finite transverse width when we treat bound states with more than one particle per tube. From here on we will adopt units¯h 2 md 2 and d for all energies and lengths.
The potential above has the interesting property that for the case of n 1 and as a function of the angles we thus see that we can obtain both positive and negative values of the integrated interaction. The criterion for the existence of a two-body bound dimer in the limit of small U is a negative integral [27,28]. We thus see that one can control the presence of a weakly bound dimer by changing φ and θ. Below we will mostly address the situation where φ = π/2 and θ = 0, i.e. dipoles that are oriented perpendicular to the tubes. In this case, the integral in (2) is − 2U n 2 and this means that the system can produce bound dimers between any two particles that are located in different tubes for any n > 0.
However, if a given few-body state has more than one particle in a single tube, then the interlayer interaction will be either attractive or repulsive depending on the sign of 1 − 3 cos 2 (θ ) cos 2 (φ) (see (1) with n = 0). Here we will only consider the regime where this quantity is positive, i.e. the case for which two dipoles in a single tube repel each other in order to avoid any collapsing states within the tubes.

Few-body bound states
We now present our results for up to five dipolar bosonic particles in multi-tube geometries. The discussion will involve analytical tools for addressing the limits of weak and strong dipolar interactions which will be compared with numerics using the exact stochastic variational approach [29,30]. This allows us to determine the range of validity of the analytical methods that we employ. Along the way we will make a detailed comparison between the 1D case with tubes and the two-dimensional (2D) case of multiple layers [7,[31][32][33][34][35][36][37]. We will also compare the results with the the results of McGuire [38] for the exact ground state energy of an N -boson system with pairwise zero-range interactions in 1D. Finally, we will investigate few-body states with more than one particle in a single tube for the case of perpendicular dipoles (θ = 0 and φ = π/2) and with tilted angles (θ = 0 and φ < π/2).

Two-body states of two dipoles in two tubes
The first case is one dipolar particle in each of two adjacent tubes with dipole moments oriented perpendicular to the tubes (φ = π/2 and θ = 0). This configuration has the Schrödinger equation where x is the relative distance of the two dipoles along the tube axis, is the eigenenergy in units of¯h 2 md 2 and is the wave function. We start with the analytically accessible limits of strong /U ∼ 1 and weak /U 1 binding.
3.1.1. Weak binding. From [28], we know that a bound state exists for potentials that fulfil For this bound state, we write the solution of (3) in the form where = −κ 2 , κ > 0 [39]. Notice that we are 'normalizing' the wave function at x = ∞ and then integrating our way to any other value of x. For x < 0 the term e −κ x appears to diverge, but the second term will also contribute and the combination yields the correct answer. We now obtain the equation for the binding energy For weak binding this can be solved iteratively for small κ, i.e.
Here we have used that (κ, y) ≈ 1 + O(κ). Upon inserting the potential in (1) with n = 1 and carrying out the integrals, we find that which implies that to leading order in U , = −κ 2 ∝ −U 2 . These results were initially discussed by Landau [27] and Simon [28] but the procedure we present here allows simpler access to the higher-order terms. We note that the first term in (6) corresponds to a delta-function potential −g¯h 2 md 2 δ( y d ) with g = 2U (where we have reintroduced explicit units for clarity). Numerically we find that this potential gives a good approximation for the real potential with U < 0.2. This is a valuable conclusion for studies that depend on a delta-function representation of the dipolar interaction in one dimension in order to describe many-body physics [21].
A clear difference from the 2D case is seen at this point. The weak-binding result in 1D is ∝ −U 2 , whereas in 2D it is ∝ − exp(−8/U 2 ) [29,[33][34][35]. This is consistent with our expectation that weakly bound 1D dipolar few-body bound states are generally more stable than similar systems in 2D, and any sort of external perturbation would have a less severe effect in the 1D setup.
We can also consider the case where U < 0 corresponds to one of the dipoles in the two tubes pointing in opposite directions, which should be possible to achieve using ac fields [40]. In this case the 2D setup with two adjacent planes will always have a bound dimer [30], whereas in 1D we find numerically that it takes a finite strength U < −4.98 to bind the two-body system. We can understand the latter if we consider κ = 0 in (5). This transforms the condition for the appearance of the bound state. For small U we have U + U 2 π 16 = 0, which has two solutions: one corresponds to bound states with U > 0 and another to bound states with U < − 16 π ∼ − 5.09. The latter is within 2% of the numerical value. However, notice that the U < 0 case opens up the possibility of binding and unbinding the dimer by tuning U as originally envisioned in [31]. The inset is a magnification of the small-U regime. Note the huge difference at small U between 1D and 2D due to the exponential decrease of the binding energy at small coupling strength in 2D.

Strong binding 1D.
For this case, we may assume that our wave function is strongly localized near the origin, i.e. x 2 / → 0, such that the probability to be in the classically forbidden region is small. Under these conditions, we can use standard techniques to obtain the binding energy and the wave function. First we decompose the potential near the origin, Using the assumption of strong localization, we see that up to terms O( 1 √ U ) we can solve the harmonic oscillator problem for the first two terms and include the third term via perturbation theory. This yields A comparison of the numerical results for all values of U to the weak and strong binding expansions is shown in figure 2. Both the limits of weak and strong binding are well reproduced by the different approximation schemes. We see differences between the numerical and analytical results only around U = 1. For completeness, we also summarize all the analytical results and their regimes of validity in table 1.

N -body chains
We now proceed to discuss the case where we have N -body chains that consist of one dipolar particle in each of N adjacent layers or tubes. In the 2D case, these sorts of structures have Table 1. Comparison between 1D and 2D binding energies of two dipoles in two different layers or tubes with perpendicularly oriented dipole moments. The regimes of validity of the different analytical formulae are indicated in parentheses. 1D generated a lot of recent interest since one expects highly non-trivial few-and many-body dynamics [31,33,[41][42][43][44][45][46]. There is a similar interest in the 1D multi-tube configurations [15,19,21,22]. The Hamiltonian for the N -body 1D system is where x i j is the relative distance of the ith and jth particles along the tube direction and n i j is the integer that signifies how far apart the tubes are that hold the ith and jth particles (n = 1 for adjacent tubes and so on). For the case of pairwise potentials of equal strength and zero range, this Hamiltonian was studied in the classical papers of Lieb and Liniger [47,48]. McGuire has shown that for this delta-function case, the bound state energy, E N , of the N -body problem can be written very elegantly as E N /E 2 = N (N 2 − 1)/6, where E 2 is the two-body bound state energy of the same delta-function potential [38].
Here we consider a system of N dipoles placed in N different tubes with the potential in (10) given by (1). This is intrinsically a case with long-range interactions that has also generated some classical studies. In the particular case where the interaction behaves as r −2 , this is the Sutherland-Calogero model [49,50]. In the present case with dipoles we have a long-range tail that behaves as r −3 . We would like to build a connection to the zero-range studies for the dipolar systems in 1D. For general polarization angles, however, the interactions are anisotropic and more difficult to compare with the isotropic zero-range models. In this section we therefore study the case of perpendicular polarization (φ = π/2 and θ = 0) where the anisotropy is absent.
As we have demonstrated above, for U < 0.2 we get universal behaviour in a sense that the potential md 2 h 2 V = −2U δ( x d ) accurately describes our system. Now we want to establish how well this zero-range approximation works on the few-body level with N > 2. To do so we first consider a toy model where all N particles interact with the same interaction corresponding to n = 1 in (1). This is precisely the case where it makes sense to compare to the analytical result E N /E 2 = N (N 2 − 1)/6. In this way, a feeling for when the zero-range approximation works. For realistic chains n is not always equal to one. This means that these are easier to describe using a zero-range approximation since particle pairs that are located several tubes apart in a chain have a potential that effectively corresponds to a smaller value of U (as compared to the n = 1 case). The numerically calculated chain binding energies relative to the two-body energy, E N /E 2 , for both the toy model and the real dipolar chains in 1D are presented in figure 3. We observe the delta-function behaviour in all cases for U < 0.15, and the strong coupling behaviour is slowly approached with increasing interaction strength. In current experiments with 40 K-87 Rb molecules [7,8], we have U 0.1 so this is weak coupling in the current context. In the strong coupling limit, the leading term for the toy model is , which is simply the number of pairs that all contribute an energy E 2 . For the realistic chain, the expression is instead where is the harmonic number of order k. We will elaborate more on the strongcoupling limit below.
We see a large difference between the toy model and the real chain in figure 3 in the weakcoupling limit U 1. This is an important conclusion of our work that has implications for many-body studies on dipolar systems that work with delta-function approximations for the dipolar interactions. The toy model is relevant for an equidistant triangular tube configuration as studied in [21]. Our results for this N = 3 case show a very flat profile for all values of U and from the point of view of energetics it is thus a reasonable approximation to use the delta- The N = 4 and 5 would be applicable to a tube configuration with four and five nearest neighbours, respectively (the latter would, of course, not be possible with standard crystal lattices). However, in these cases the variation of the energies with U is much more drastic, indicating that a zero-range description is only good for very small U .
For the real chain systems, the curves are rather flat for all N . This is connected to the fact that the terms with n > 1 in (10) will give a contribution that is suppressed. Due to the flatness of the curves in figure 3, we can obtain a good approximation for the N -body binding energy of the realistic chains by using an effective model for N -body systems with The blessing of working with a zero-range interaction comes at the price of having an N -dependent coupling. This ensures that the N -body bound states present in this geometry are properly described. We will use this effective dipolar strength to study the corresponding manybody system below. We stress that the delta-function approximation with U eff (N ) coupling is only accurate for weak interaction, i.e. small U . In the case of N = 5 from U → 0 and up to U = 10, we find a difference of about 20% on the value while for N = 3 and 4 the difference is slightly less. So the flatness and thus the reliability of the delta-function approach is on the level of 20% for U 10 and N 5. It is very important to note that we reproduce the strong-coupling limit, i.e. large U , with this choice of U eff . This means that if we use the potential in (12) as an effective potential between N bosonic particles in 1D, we obtain the correct strong-coupling binding energy (11). This all builds on the flatness of the energy in figure 3 and is only as good as the uncertainty quoted above. This is important to note as below we will look at low-energy many-body physics using this effective interaction.

Chains in the oscillator approximation.
As we have seen for small U , the zero-range approximation gives reliable results. However, we would also like to discuss the strong-binding limit where U 1 and provide a simple analytical procedure for obtaining the energies and wave functions. This can be done using exact harmonic models [51,52] where the full twobody interaction (the dipolar interaction in this case) is represented by a harmonic oscillator with parameters that are carefully chosen to reproduce energetics and structure at the two-body level [45]. We expect this to be a good approximation for U 1 since the potential in (1) has a deep pocket in this case. However, it turns out that this approximation can be very accurate even for moderate U in the case of 2D multi-layered systems with chains [53].
In figure 4, we show a comparison of the oscillator approximation to exact numerical calculations for the case of perpendicularly oriented dipoles (φ = π/2 and θ = 0). Note that the horizontal axis begins at U = 2 (for U < 2 the oscillator does not reproduce the exact results well). We get very good agreement for N = 3, while for higher particle numbers there is a clearly visible deviation. This stems from issues with the two-body interactions coming from the outermost tubes where the oscillator approximation starts to fail since the contributions are not well described by an oscillator (the effective U for these terms is small for any U value considered here). However, the overall agreement is still within a few per cent for U > 5. One way to improve the agreement for higher particle numbers could be to introduce an effective three-body force (more generally an N -body force) in addition to the two-body potential terms. Such an approach is very common in nuclear physics and in effective field theory studies of universal three-body physics [54,55].  We now change the angle φ to investigate the effects of anisotropy. Figure 5 shows the case with φ = φ m and θ = 0, where φ m satisfies cos 2 φ m = 1 3 . Recalling the discussion after (1), this implies that the potential between two dipoles in the same tube vanishes. First we observe that the overall energy scale becomes smaller, because the depth of the potential and the net volume of the potential are smaller than for the perpendicular polarization case in figure 4. For N = 3 and 4, we again see very good agreement between the oscillator approximation and the numerics. The largest chain with N = 5 shows what appears to be an even better agreement than for the φ = π/2 and θ = 0 case. What is also striking is that the oscillator results and the numerics cross around U ∼ 5. This happens since our method ensures that the oscillator reproduces the slope of the numerical results in the limit U 1 (as discussed for the 2D case in [53]). The optimization of the parameters for the two-body oscillator potential for 2 U 5, however, produces a slightly different slope and a crossing thus takes place.

Non-chain bound complexes
We now consider the case where there are two particles in one tube. The simplest configuration of this sort is two tubes with three particles in total, i.e. one particle in one tube and two particles in the adjacent tube. At this point, we need to take the width of the tube that we discussed in section 2 explicitly into account in order to work with an intratube interaction that best describes the experimental situation. This can be done by folding the potential with a Gaussian wave packet of width l that represents the wave function of the system in the transverse direction (along the y-axis in figure 1). This yields the intratube interaction [25]  where and λ = d/l and U rep is the strength of the intratube interaction. This means that the intratube potential has been regularized by the finite width of the transverse confinement and does not have the strict 1D form x −3 expected in the limit of λ → ∞. In a system where the dipole moment is the same for all particles in all tubes, we have U rep = U , where U is the strength of the intertube interaction. The dipole orientation and the induced dipole moment are typically controlled by applying an external electric or magnetic field that is constant across the whole system. However, if one applies also a field gradient then it should be possible to obtain a system where the induced dipole moment varies from tube to tube (or layer to layer in the 2D multilayer case). In this case we can have a situation where U rep = U and we thus vary both quantities independently in the current study. As discussed previously, we will consider only the repulsive case where 1 − 3 cos 2 φ cos 2 θ > 0 in order to avoid possible intratube collapse of the system. Unless explicitly stated, we use λ = 5 everywhere (see [25] for a discussion of the effects of changing λ). We start by considering the case where all particles are perpendicularly polarized, i.e. φ = π/2 and θ = 0. Observe that if U rep = 0, we have a trimer energy that is below the dimer energy, and if U rep → ∞ we cannot bind three particles. We now want to find the value U cr rep at which the three-body system is bound for a given value of U , so we keep U fixed and vary U rep . In [36] it has been shown that U cr rep /U is an increasing function of U by an argument that does not depend on the dimensionality and thus applies equally well to layers in 2D and tubes in 1D. In turn we only need to consider U → ∞ to obtain an upper bound for U cr rep /U . This upper bound is a decreasing function of λ, and for λ = 5 we obtain U cr rep /U < 0.41. In order to find the lower bound for U cr rep /U , we will use the variational principle with a normalized three-body wave function of the form (x, y) = (x)ϕ(y), where (x) is the normalized dimer wave function and we have used Jacobi coordinates √ 2x = x 1 − x 2 and √ 3/2y = x 3 − (x 1 + x 2 )/2. From this we obtain the upper bound for the three-body energy where the kinetic term for the y-coordinate is T y = −¯h Here V rep is given in (13) and V is given by (1) with n = 1. Equation (5), in principle, allows us to estimate where E var 3 < E 2 holds, which gives us the lower bound for U cr rep . However, in general this equation is not easy to analyse so we obtain the lower bound using the simpler condition for the effective potential acting on ϕ in (15), which is that it must have negative net volume. The limit of zero net volume requires which can be simplified to This lower bound is a decreasing function of λ, and for λ = 5 we obtain U cr rep /U > 2/25. Collecting our results, we have for λ = 5 that 0.41 > U cr rep /U > 2/25. This means that if we could manipulate the values of the intratube and intertube interactions independently, we could create or destroy the three-body bound state in two tubes with perpendicular polarization. This effect is completely due to the finite range of the potential, because for zero-range potentials any value U rep /U will produce a bound state (the system becomes equivalent to the fermionic states studied by McQuire in [56] for which a bound state always exists).
In figure 6 we show the critical value, U cr rep /U , as a function of U for both the cases of 1D tubes, and for comparison we also plot the results of the same study for a bilayer in 2D [36]. Both these calculations were performed with λ = 5. We see a clear difference between the critical values for three-body bound state formation in the two cases with the 1D situation allowing for larger repulsive interactions. This is connected to the price one pays for localizing particles in a bound states complex. In 2D we need to localize in two independent directions which is more expensive than 1D. As expected, 1D tubes are more favourable for formation of these larger complexes.
Having discussed how a tuning of the intratube repulsion can bind a three-body state, we now consider an alternative route to bind this state, which is by tilting the angle, φ. We assume that U r ep = U , and are looking for values of φ and U for which a bound state exists. For a given angle φ we introduce the notation U cr (φ), which is the value at which the two-body dimer (one particle in each of two adjacent tubes) energy is equal to the three-body energy which defines the threshold. We start with three particles as before but tilted by angles θ = 0, φ ∈ [φ m , π/2]. We already proved that for φ = π/2 we do not have a bound state if U rep = U . However, we know that for φ = φ m we always have a three-body bound state, because the intralayer interaction vanishes. This implies that there is a range of angles where the three-body state exists. To approach this problem analytically, we can calculate the net volume as above. This yields Putting it to zero, we obtain the smallest angle below which the three-body bound state exists for any value of U . The condition for a three-body bound state is φ cos −1 λ 2 −2 3λ 2 −2 , which in the case of λ = 5 yields φ 0.975 or 55.9 • . In the strict 1D limit where λ → ∞, we find the result φ φ m . The function λ 2 −2 3λ 2 −2 is decreasing with λ, which is a result of the behaviour of the function in (14) that has a finite value for zero argument for any finite λ. In strict 1D we expect that the repulsion will not allow this three-body bound state due to the cubic singularity x −3 of the repulsive potential. For the more realistic case with a finite λ, we can have bound states for the angle in the regime φ > φ m as well. An additional bound can be found by using classical arguments, i.e. neglecting kinetic energy in the Schrödinger equation and studying just the potential landscape [36]. This yields a bound of φ 1.15 or 66 • so that for any U we do not have three-body bound states in the case where U rep = U . In summary, we have found that for φ m < φ 1.15 there exists a value of U cr (φ) ∈ [0, ∞), such that for U > U cr (φ) the bi-tube trimer is bound.
In figure 7 we plot U cr as a function of φ for three particles in two tubes in the case where U rep = U , i.e. the case of current experimental interest where the dipole strength is controlled by one external field that is constant over the entire samples. The points indicate the values at which the dissociation points have been obtained numerically (connecting lines are a guide to the eye). The three-body bound states are present on the left-hand side in the diagram since here the intratube repulsion is not strong enough to break the bound state apart. The results obtained , where the trimer in two tubes (one and two particles) and the tetramer in three tubes (one, two and one particles) dissociate as a function of the tilting angle φ for θ = 0. The dissociation thresholds were calculated at the dots and squares with connecting lines to guide the eye. The vertical lines mark the angles beyond which a classical calculation (neglecting kinetic energy for all particles) predicts that no bound states can occur because the system is overall repulsive.
here are consistent with results [25] for U < 15 where a comparison can be made. For U > 15, we find a very sharp increase in U cr (φ) that takes place at φ ∼ 1.1, a few degrees before a classical calculation renders the system unbound at φ = 1.15.
In order to address the question of additional tubes and their effect on bound states with more than one particle per tube, we have also considered a four-body state in three tubes (as in figure 1) with two dipoles in the middle tube and one in each of two outer tubes. The results are also presented in figure 7. We again see a sharp increase for large φ before we get to the classically unbound value of φ = 1.34 or 76.5 • . However, there is a significant increase of region in which bound states with two particles in a single tube are present due to the extra attraction from adding a particle in an adjacent tube on the other side of the complex as compared to the three-body system with two and one in two tubes. Other four-body states can occur in two tubes as discussed in [25], but the parameter regime is severely suppressed due to the extra repulsion in the three-plus-one and two-plus-two four-body states in two tubes. We therefore expect that the case of three tubes could have interesting 'quasi'-Wigner crystal states for large U beyond the clustered states studied in [18] where the clusters are not of the same size in each tube (here with cluster size one in the outer tubes and two in the inner tube). We can continue this kind of game for five tubes in a one-two-one-two-one configuration and so on.

Scattering
While we are mainly concerned with negative energy bound state structures in this paper, we now present a short overview of some scattering results that can be deduced directly based on the formalism we introduced above. One basically only needs to make an analytical continuation of the bound state discussions above in order to arrive at results for scattering of dipolar bosons in 1D tubes. Since current experiments with cold molecules aim for the low-temperature regime where quantum degeneracy can be reached, we will focus exclusively on the low-energy scattering regime here.
The formal solution of the scattering problem, , satisfies the equation where the boundary condition is an out-going flux to the right with coefficient one of the form e ikx [39]. Here we define the reflection, R(k), and transmission, T (k), coefficients through (k, x → −∞) = T (k) −1 e ikx + R(k) T (k) e −ikx to be consistent with standard 1D potential scattering. By taking the limit x → −∞, we obtain For k = iκ, κ > 0, poles of T (k) define bound states with the binding energy defined through the equation which is consistent with (5) and allows us to use the same techniques as before to obtain iterative improvements of the solution.
Consider now the case of dipoles oriented perpendicular to the layers, i.e. φ = π/2 and θ = 0. In the limit of weak binding and low-energy scattering, U 1 and kd U , we find for the transmission coefficient of two particles in the two different tubes and for the case of two particles in the same tube we find These results are deceptively simple, yet they contain all the information about the dipolar physics in 1D tubes at low energy in the small-U limit. As expected they vanish with k as k → 0, which is a manifestation of the Wigner threshold law. Also, as κ b ∝ U in the weakcoupling limit (see table 1), we see that both expressions become large at very small U . This makes sense since the potential becomes very weak (and likewise the derivative of the potential).
We also see that in the strict 1D limit, λ → ∞, the transmission is completely suppressed and total reflection is expected. This is the emergence of the impenetrable boson regime called the Tonks-Girardeau limit [57,58]. This limit can be studied with dipolar bosons in a 1D setup with an in-tube harmonic trap as discussed in [59].
Since we would like to obtain an effective model for dipolar systems of N -chains using zero-range interactions, we need to relate the transmission coefficients above to the zero-range equivalent. If we write the potential as V (x) = gδ(x), we find for low-energy scattering. We see that we recover the g = −2U result for the case of two dipoles in different tubes (using our units of energy and length given in section 2), and we find g = U λ 2 for dipoles in the same tube.

Many-body physics
We now want to address the case of a many-body system of N adjacent tubes with M particles in each tube. We will focus on the case of perpendicular orientation of the dipoles with respect to the tubes (θ = 0 and φ = π/2) and comment on the case of general angles in the end. In the perpendicular case, we expect no bound complexes beyond the chain structures that are bound states of N bosons with one in each of N adjacent tubes. Recent Monte Carlo studies of the many-body system indicate that the longest chains possible, the N -chain for N tubes, are the relevant degrees of freedom in these systems [22,44], and we therefore focus our attention on these chains. The results that we have presented thus far show that one can make an effective model of chains using zero-range interactions with suitably adjusted coupling strength as in (12) which reproduces the chain energy in the strong-coupling limit. This is extremely useful for manybody physics since a number of 1D N -body models with zero-range interactions are exactly solvable as discussed above. For the case of bosons with attractive zero-range interactions, one in fact finds a gas-like state similar to the Tonks-Girardeau gas of repulsively interacting bosons, denoted the super-Tonks-Girardeau gas [60][61][62]. For both attractive and repulsive zerorange interactions, it is also possible to describe these systems using bosonic Luttinger liquids (see [63] for a recent review of 1D bosonic systems). These studies have been boosted by recent experimental success in realizing these 1D quantum gases [64,65]. For dipolar bosons this Luttinger liquid behaviour was studied by several authors [9,10,13,14].
Here we want to address the many-body behaviour by using the few-body physics information that we have derived above. For N adjacent tubes with M particles in each, we can build a system of M chains, and then we may study the effective interchain interactions. Using a zero-range interaction model allows us to map the problem into the realm of exactly solvable 1D systems and/or Luttinger liquid studies and thus infer the many-body dynamics. Here we will consider the low-energy physics only, i.e. kd U . This is important for our use of (23) for the intratube repulsion that we will return to below.
The effective zero-range coupling in (12) has the property that it (approximately) reproduces the N -body binding energy of a single chain for large U as discussed after (12). However, if one considers two chains then one realizes that this effective interaction also describes the attractive force between two different chains. This follows readily from the fact that the intertube interactions that are responsible for the chain-chain interactions are the same as those that provide the intertube attraction in a single chain. If we ignore for a moment any intratube repulsion, we can estimate the ground state energy of a system with M chains of length N to be simply where E cc (N ) is the two-body energy between two different N -chains (considered as effective degrees of freedom in the system) calculated using the potential md 2 The first term is obtained by mapping to the result obtained by McGuire [38] for M bosonic chains (valid for any U eff and thus any U ) and the second term comes from the internal energy of each of the M chains of length N . The N -dependence is given implicitly through U eff .
The repulsive intratube interaction has thus far been ignored. We would now like to include it in a simple way through the use of (23). This was derived under the assumptions U 1 and kd U . The latter is not a problem, while at first glance U 1 is not consistent with the fact that we treat the attractive interactions of the chains in a large U limit. However, we note that (23) has the correct qualitative behaviour when U increases; the repulsive potential suppresses transmission as T ∝ U −1 . We thus assume that we can still use the expression in (23) and thus also a delta-function potential with g = U λ 2 . This is effectively a weak-coupling expression and thus our uncertainty lies in pushing the attractive part, (12), to the weak-coupling limit. This has a 20% uncertainty for U 10 for N 5 as noted below (12). So if we stick to systems with N 5, the argument we have here is only expected to hold to within that uncertainty.
Assuming that we can push the attractive part to weak coupling and thus use a deltafunction for both attraction and repulsion, it is not difficult to take into account through a zerorange interaction term of the form md 2 h 2 V r (x) = NU λ 2 δ( x d ). The factor of N is obtained by noting that there will be one repulsive contribution from each tube separately. In our effective model of chains, we thus have a competition between two interactions of opposite sign, V a (x) and V r (x). We can use this to estimate when the system of chains constitutes an attractive or a repulsive Bose gas by determining when V a (x) = V r (x). Since both are proportional to U , this will give a purely geometric condition. If we use effective delta-function interactions between chains of length N and define then V a (x) = V r (x) corresponds to F = 0, while F > 0 gives an attractive and F < 0 a repulsive chain system. Note that this is not dependent on the number of chains, M, since we are considering the balance of attraction and repulsion by looking at the effective delta-function two-body interactions between two chains that one would put in the M-body Hamiltonian for M chains as the two-body interaction term. What is important is then the sign of the effective two-body term which is determined by (26) and depends only on the chain length N . We thus see immediately that except for very small λ, the system will be in the repulsive regime for any N 2. Small λ corresponds to a soft transverse confinement that requires careful treatment beyond our quasi-1D formalism. We therefore conclude that the system of dipolar chains in the intermediate to strong-coupling limit (U 10) with dipoles oriented perpendicular to the tubes should behave like a repulsive Tonks-Girardeau gas. The effective parameters introduced here can be used to estimate ground state properties using the Lieb-Liniger model [47,48] and/or Luttinger liquid formalism can be used to obtain correlation functions, etc.
However, we caution that we mix the limits of strong and weak coupling and that carries the uncertainty discussed above. A means of improving this is to carry out the derivation of (23) to higher order. That would generate a different function of U that could then be inserted in (26). The crucial point is that going to higher order in (23) would give a function with higher powers of U in the denominator of (23). Consequently, the second term in (26) would have a factor U if we go to second order in (23). Since it corresponds to a repulsion, it will also not change the sign of the second term in (26). Since the first term in (26) does not carry factors of U , this means that at higher order the repulsive part will be even stronger and we expect this to lead to the same conclusion: F(N , λ) < 0 and effectively a repulsive system of chains for larger values of U where the chains are well bound and thus well defined.
If we tilt the dipolar orientation away from perpendicular, we expect a smaller effective coupling due to the angular factor in (1). The repulsive part will, however, remain dominant and we expect the same conclusion as above to hold for small tilting. If we go to the value φ m and θ = 0 where the intratube repulsion is zero, we naturally obtain an attractively interacting system. However, here we need to worry about more complicated complexes than bound states of chains as we have discussed above (this was first pointed out in [24]). One could imagine that if we take a long chain (large N ), then it may be possible to attract one more particle to bind to this chain around the middle tube. This may open up the possibility for scattering events between chains where this extra particle is exchanged, akin to transfer reactions known for instance from nuclear physics. The nature of potentially more complicated bound structures constitutes an interesting direction for future work.

Conclusions and outlook
We have studied dipolar bosonic particles confined to a setup consisting of a number of 1D tubes. For the case of equidistant adjacent tubes, we calculated the few-body bound state structures for up to five particles using both analytical and numerical tools in the weak-and strong-binding limits. For the case of perpendicularly oriented dipoles, we find that chains with one particle in each of a number of adjacent tubes are the most stable structures, and that more complicated bound structures with multiple particles in a single tube are unlikely unless one tilts the dipole orientation away from the perpendicular direction. The bound states in our study should be directly observable using experimental techniques such as lattice shaking [66], RF spectroscopy [67,68] or in situ optical detection [24,25].
Our results for the weak-coupling limit can be used to calculate the low-energy scattering properties of the system and make an effective zero-range model for bound chains of dipolar bosons that takes both attractive and repulsive terms into account. This allows us to use our knowledge of the few-body physics to address the many-body physics of these systems that can be mapped onto an exactly solvable model in 1D and Luttinger liquid dynamics. Estimation of the magnitudes of the competing attractive and repulsive terms in a system of chains indicates that for perpendicular orientation the repulsion is always dominant and a Tonks-Girardeau type of system should emerge. This could be changed by tilting the angle of the dipoles with respect to the tubes, which changes the balance between repulsive intratube and attractive intertube terms. However, the potential presence of bound state complexes with two or more particles per tube needs to be considered in this case. The latter will be the focus of future investigations.
In the limit of strong binding, we demonstrated that models based on replacing the true dipolar potential by a harmonic interaction are quite accurate for describing the energetics. This yields exactly solvable N -body models that could, for instance, be used to study thermodynamic properties in 1D dipolar systems [69,70]. Another intriguing possibility is to study a system where one or more of the tubes have perpendicular dipoles that point in the opposite direction. These flips can be obtained by using ac fields [71]. For two dipoles in two adjacent tubes with opposite orientation of the dipoles, we proved here that in the limit of weak dipole moment there will be no bound two-body state. However, this does not, in general, rule out the possibility of a bound three-body state, a so-called Borromean state [72]. We could thus imagine a regime where three-body states are the lowest non-trivial bound states and thus gives rise to an effective gas of trimers. While all of the questions addressed in this paper have been concerned with bosonic dipolar particles, we expect that similar bound state structures should arise when using fermionic dipoles. The fermionic case should correspond closely to studies on 1D multi-component systems [73,74] if one maps the component index onto the tube index.