# Structure, Bose-Einstein condensation and superfluidity

of two-dimensional confined dipolar assemblies

###### Abstract

Low temperature properties of harmonically confined two-dimensional assemblies of dipolar bosons are systematically investigated by Monte Carlo simulations. Calculations carried out for different numbers of particles and strengths of the confining potential yield evidence of a quantum phase transition from a superfluid to a crystal-like phase, consistently with what is observed in the homogeneous system. It is found that the crystal phase nucleates in the center of the trap, as the density increases. Bose-Einstein condensation vanishes at upon entering the crystalline phase, concurrently with the disappearance of the superfluid response.

###### pacs:

## I Introduction

Bose-Einstein Condensation (BEC) and superfluidity (SF) are intimately connected quantum many-body phenomena Leggett (2006), both underlain by long cycles of exchange of identical particles Feynman (1953). In three dimensions BEC and SF appear simultaneously, and become almost synonymous when interactions among particles are weak, a condition normally experimentally attained in trapped ultracold gases Davis et al. (1995); Anderson et al. (1995). In that limit, the fractions of the system that are superfluid and Bose condensed both approach 100% at the temperature . This provides a major justification for theoretical studies based on mean-field techniques Krauth (1996); Heinrichs and Mullin (1998).

On the other hand, the effect of strong inter-particle interactions on SF and BEC is very different. A chief example is provided by liquid He, which approaches the 100% superfluid limit as , while its condensate is largely depleted; indeed, theoretical Moroni and Boninsegni (2004a); Boninsegni et al. (2006a) and experimental Azuah et al. (1997); Glyde et al. (2000) estimates for the condensate fraction in superfluid He yield a value in the neighbourhood of 7-8%, at =0. Furthermore, no theoretical understanding currently exists of the relationship, if any, between the local values of the superfluid and condensate fractions. Broadly speaking, the most quantitative aspects of the relationship between these two quantities remain to be further elucidated.

Impressive advances in the trapping and manipulation of ultracold atoms, appear to offer a pathway to a more in-depth investigation of this relationship, by making it feasible to create in the laboratory “artificial” many-body systems characterized by interactions of variable strength and shape. Thus, one may smoothly interpolate between weakly and strongly interacting regimes, as well as explore the effect on SF and BEC of interactions not (easily) realized in ordinary condensed matter systems.

In recent times, attention has turned to the case of dipolar bosons, further bolstered by the realization of Bose-Einstein condensation of Chromium Griesmaier et al. (2005). The long-range and anisotropic nature of the dipole-dipole interaction leads to many fascinating phenomena Lahaye et al. (2009). The conceptually simplest scenario is that of an assembly of dipolar bosons confined to two dimensions, with their dipole moments aligned by an external (electric or magnetic) field directed perpendicular to the confinement plane. In this case, the interaction between any two particles is purely repulsive, proportional to the inverse cubic power of the distance between them.

Recently, it was shown that a mixture of equal-mass dipolar isotopes, in such a configuration, demixes at finite temperature due to quantum statistical effects Jain and Boninsegni (2011).
A promising direction for implementing trapped dipolar assemblies, is through the use of Rydberg-excited atoms. The large induced dipole moments of such systems may facilitate the realization of solid-like phases Anderson et al. (1998); Pupillo et al. (2010).

The phase diagram of a homogeneous system of purely repulsive dipolar bosons in two dimensions, has been explored by Monte Carlo simulations Büchler et al. (2007); Astrakharchik et al. (2007); Mora et al. (2007). These works have yielded evidence of a quantum phase transition at =0 between a superfluid and a crystal phase, but the quantitative characterization of such a phase transition, e.g., its location and the width of the coexistence region, has proven surprisingly difficult, due to the long-range nature of the interaction. Theoretical arguments Spivak and A. (2004) have been put forth to the effect that a conventional first-order quantum phase transition ought not occur in this system, as the long-ranged nature of the interactions renders the coexistence of two phases of difference densities, separated by a single macroscopic interface, energetically unfavourable. Rather, an “emulsion” should form, consisting of (relatively) large solid domains embedded in the superfluid. Such a scenario has not yet been observed in computer simulations, presumably due to the need of studying systems of size not currently attainable.

The superfluid transition at finite temperature has been predicted to be compatible with the Berezinskii-Kosterlitz-Thouless (BKT) theory Filinov et al. (2010). And, while no evidence has been seen of a commensurate supersolid phase Spivak and A. (2004); Mora et al. (2007); Kurbakov et al. (2010), its presence has not been ruled out yet (it is in fact predicted to occur in a system of dipolar bosons in the so-called Rydberg-blockaded regime Cinti et al. (2010)).

It is expected that many of the outstanding issues will be soon elucidated by experiments on spatially confined assemblies of dipolar particles. Two-dimensional confinement can be achieved in the laboratory by means of a harmonic trap in the transverse direction, of frequency . Clearly, the system must also be confined in the remaining two dimensions, as the interactions are purely repulsive, and this is accomplished by means of a second, in-plane harmonic trap, of frequency , realizing the so-called pancake geometry. The number of particles that can be currently trapped is typically of the order of a few thousands, which means that the physics of such a system is strongly influenced by its finite size. Thus, theoretical studies of finite two-dimensional dipolar systems in a harmonic potential seem especially opportune and timely, in order to guide in the interpretation of experiments aimed at inferring bulk properties.

However, a finite system of dipolar bosons is of interest in its own right, chiefly because novel fundamental understanding of the relationship between BEC and SF can be gained.
Indeed, theoretical work on clusters of He Lewart et al. (1988); Sindzingre et al. (1989) has shown that superfluid and condensate fractions can be meaningfully and rigorously defined for systems comprising as few as =70 particles. One should note that, for the two-dimensional systems considered in this work, because the system is finite, reduced dimensionality does not imply absence of BEC at finite temperature Bagnato and Kleppner (1991); Hadzibabic et al. (2008).
Moreover, while superfluid and condensate are uniform in the bulk system, owing to translational invariance, in confinement they can display non-trivial local variations. Of particular import seems the comparison of the local superfluid and condensate fractions, which can lead to original insight into their interplay.

Previously, Lozovik et al. Lozovik et al. (2004), and more recently Pupillo et al. Pupillo et al. (2010), have studied harmonically confined dipolar systems of relatively small size (up to particles).
In this paper, we investigate by means of Monte Carlo simulations the low temperature phase diagram of dipolar Bose assemblies comprising up to particles, a number which appears relevant to current (or planned) experiments.

The purpose of this work is twofold. On the one hand, we compare the ground state phase diagram of the trapped system against the predictions made in Refs. Büchler et al. (2007); Astrakharchik et al. (2007); Mora et al. (2007) for the bulk case. Our aim is to make contact between the theoretical phase diagram of the bulk and the quantitative information that can be provided by experiments on finite-size systems, comprising a relatively small number of particles.
Secondly, we compute both global and local densities of superfluid and Bose-Einstein condensate, based on numerically unbiased estimators, and study how the two quantities evolve as the strength of the interactions in the system increases (as measured by the density in the middle of the trap).
To the best of our knowledge, this is the first system for which both the local superfuid and condensate densities are computed simultaneously for a strongly interacting many-body system, based on a methodology free from approximations.

This paper is organized as follows. In section II we introduce the model Hamiltonian for a two-dimensional system of dipolar Bose particles harmonically confined, and discuss the relevant observables pertinent to structure, Bose-Einstein condensation and superfluidity. We illustrate our results in the following three sections, Section III focusing on the structure of the system with respect to system size and in-trap mean inter-particle distance, Section IV focusing on the superfluidity and condensation, notably both global and local properties, and finally Section V discusses the ground state phase diagram, as well as the momentum distribution. We outline our conclusions in Section VI.

## Ii Model and methodology

We consider a two-dimensional system comprising bosons of mass , confined by a harmonic potential of frequency . Particles then interact via a purely repulsive, pair-wise dipolar potential, given by , being the distance between the particles. The simplest way of realizing such an interaction, is with particles possessing an electric (magnetic) dipole moment, aligned along the transverse direction by an applied electric (magnetic) field Büchler et al. (2007). The quantum-mechanical many-body Hamiltonian is given in dimensionless form by:

(1) |

where is the position of the th particle, and where we have introduced length and energy scales respectively as and . The trap strength parameter is , being the oscillator length. It is also useful to introduce a reduced temperature , with expressed in units of .

We investigate the low temperature properties of this finite system by means of Monte Carlo simulations. We consider systems with particle numbers ranging from 20 to 1000, and trap strengths . We use the Worm Algorithm, in the continuous-space path integral representation Boninsegni
et al. (2006b, a). The methodology is numerically exact, errors being only statistical in character, and reducible to an essentially negligible size with a relatively modest computer time expenditure. Technical details of the calculation are the same as in other studies Mezzacapo and Boninsegni (2006, 2007); Cuervo et al. (2005); the use of the dipolar interaction in (1) entails no particular difficulty.
All of our quoted ground state results, with their uncertainties, are obtained by extrapolating our numerical estimates to the limit of temperature . In general, we have observed that the structural and energetic properties of the assemblies considered here do not change, within the statistical error of our calculations, if the temperature is 0.1, below which the system may be regarded as essentially in its ground state. As customary in this type of numerical study, we have also carried out extrapolation of our estimates to the limit of vanishing time step Cuervo et al. (2005). As a typical example, Fig. 1 shows the total (potential plus kinetic) energy per particle as a function of (imaginary time) step size . As the step size is made sufficiently small, the energy converges to a fixed value.

As mentioned above, one of our aims is to make contact between the ground state properties of the finite assembly considered here, and those of the same system in the bulk limit, where the physics is determined solely by the mean inter-particle distance , being the particle density. For the trapped system, the density is not uniform, but rather is a function of the radial distance from the centre of the trap. Near the centre of the trap, the physical behaviour of the finite system ought resemble most closely that of the bulk; thus, we associate a value of to the trapped system, by extrapolating to the behaviour of the function

(2) |

Details of the extrapolation are given below.

Besides the radial density , we compute the local superfluid (), and condensate () densities. The average superfluid density is computed using the area estimator Sindzingre et al. (1989). The local superfluid density is computed by means of an estimator first proposed by Kwon et al. Kwon et al. (2006):

(3) |

where (we set the Boltzmann constant equal to one), stands for thermal average, is the area swept by the many-particle in the plane, and where and are respectively the contributions to from paths in a shell of radius , and its classical moment of inertia.

The estimation of the condensate fraction is carried out based on the ideas expounded in Refs. Löwdin (1955); Onsager and Penrose (1956). Specifically, we compute by Monte Carlo the angularly averaged one-body density matrix

(4) |

where and are the usual Bose field operators. The matrix has positive-definite eigenvalues, the largest of which is the average number of particles in the maximally occupied natural orbital (condensate), which is the corresponding eigenfunction caz (). The assumption is made here that the condensate should be cylindrically symmetric. As we are only interested in its largest eigenvalue here, and its associated natural orbital , we solve the eigenvalue equation for iteratively (numerically), using the power method. The condensate fraction is then given by with corresponding local condensate density given by .

## Iii Structure

We begin by discussing the structural properties of the system as a function of the density of particles in the trap, in the (ground state) limit. Numerical studies of the ground state phase diagram of a two-dimensional homogeneous system of dipolar bosons have yielded evidence of a superfluid-to-crystal quantum phase transition Büchler et al. (2007); Astrakharchik et al. (2007); Mora et al. (2007). The precise determination of the location and width of the coexistence region has so far proven rather challenging. Büchler et al. Büchler et al. (2007) place it within the relatively large range , and other estimates Astrakharchik et al. (2007); Mora et al. (2007); Filinov et al. (2010) are in general agreement. The question is what remnants, if any, of such a phase transition can be observed on a finite system.

A first, general comment is in order: one must exercise caution when referring to phase transitions in finite-size (mesoscopic) systems. It is well known that, although such systems may exhibit “phase-like behaviour”, the coexistence of two phases, which is typical of a first-order phase transition, is energetically disfavoured by the presence of an interface separating the two phases. Thus, in a finite-size system, especially of the size treatable in any realistic numerical simulation, phase transitions are smoothed out Fisher and Berker (1982). To some extent, therefore, even though we shall make use of this terminology, any reference to “phase” must be understood qualitatively. As well, obviously the fact that the system is confined by a harmonic trap can be expected to influence qualitatively its physics, especially at the border where particle localization is enhanced by the steepness of the confining harmonic potential.

Figure 2 shows instantaneous configurations (i.e., particle world lines) of a system comprising =200 dipolar bosons, for different values of the trap strength , increasing from (a) to (d). These results corresponds to temperatures sufficiently low that the system may be regarded essentially in its ground state. These snapshots offers visual, qualitative insight into the evolution of the system as the density inside the trap is increased. The presence of a featureless fluid phase at low density, and the emergence of crystalline order at high density, with particles arranged on a triangular lattice near the centre of the trap, are clear in Figure 2. Specifically, at high density (i.e., greater trap strength) particle world lines become localized, and exchanges of (indistinguishable) particles suppressed, as shown for the largest trap strength in Figure 2(d).

More quantitative information is offered by plots such as those of Figure 3, which shows radial density profiles for a trapped system comprising =200 particles, for different values of the confining strength of the trapping potential. The upper panel of the figure shows the local density of particles . At low density (i.e., low values of ), is essentially constant near the centre, whereas at the edge of the trap, where the confining potential dominates the energetics, a shell structure forms. On the other hand, at high density the dipole-dipole interaction energy dominates in the centre of the trap, where displays oscillations that are consistent with a regular arrangement of particles.

On the expectation that the physical behaviour of the system should mimic that of the bulk near the centre of the trap, we may compare results for a trapped system with those in the thermodynamic
limit by attributing an effective value of to the confined system. As explained above, we do so by taking the value of at (i.e., the centre of the trap), or, in case of oscillatory behaviour (e.g., triangles in Figure 3(b)), we take the value of averaged over few periods of oscillation near the centre. The lower panel of Figure 3 shows typical results for , as defined as in (2). We estimate the uncertainty in the values of determined in this way, quoted in the remainder of the paper, to be of the order of 5%.

For a system with =200 particles, the results shown in Figure 3 suggest that a transition from a uniform liquid phase to a crystalline one, takes place for roughly below . A more precise assessment of the density region within which the transition takes place, can be obtained by studying superfluidity and Bose-Einstein condensation; we come back to this later on.

In order to establish a more quantitative connection between the physics of the trapped assembly and that of a bulk system, we compute the pair correlation function in a circular region of radius , centered at . Such a “local pair correlation function” inside this region, , has the normalization

(5) |

where is the average number of particles in the bounded region, and the
corresponding average density is given by .

Figure 4 shows the computed for trapped systems with different number of particles, specifically =100 (diamonds), 200 (squares), and 1000 (circles), in the limit. Results are shown for different trap strengths , adjusted to yield a value of in the middle of the trap around 0.03 (top panel) and 0.085 (bottom panel). Here, we took and 0.3 respectively, i.e., values smaller than the radius at which the density goes to zero in the trap.
For , the ground state of the bulk system is crystalline, a superfluid for . Both figures also show (solid lines) the corresponding pair correlation functions for the bulk system with the same values of , also computed in the 0 limit (a separate calculation was performed for this case, using standard bulk methodology, with periodic boundary conditions).

The first observation is that in both cases, on increasing the number of particles approaches that of the bulk system, as expected. Moreover, the finite size of the system, combined with the fact that the trap is harmonic, enhances particle localization and therefore strengthens crystalline order. Results shown in Figure 4(a) correspond to in the middle of the trap; displays marked oscillations for both =200 and 1000 particles. However, for the smaller system size the height of the peaks is lower, and the distance between successive peaks is also different from what observed in the bulk. This suggests that, although the physics of the system is solid-like in both cases, for the smaller system it is substantially influenced by the trap, whereas with five times more particles it is already rather close to that of the bulk. Indeed, the pair correlation function in the center of the trap is virtually indistinguishable from that of the bulk system beyond the first peak, within our statistical errors. It displays the pronounced oscillations that one expects for a crystalline phase.

Figure 4(b) shows results for the case . Here too, the local pair correlation function has a noticeable dependence on the size of the system. However, in this case the distance between the peaks is roughly the same in all cases, i.e., the influence of the edge on the physics of the system at the center of the trap is less significant than in the case previously discussed (not surprisingly, as the trap is weaker). As the size of the system is increased from to particles, the most significant change is the height of the first peak of , which decreases by roughly 30%. Beyond the first peak, the pair correlation functions for =200 and =1000 are essentially indistinguishable from one another, as well as from that of the bulk system (again within the statistical uncertainties of our calculation), and are consistent with liquid-like behavior.

Figure 5 shows snapshots as in Figure 2 for at low temperature, in correspondence of two values of shown in Fig. 4. For , an orderly arrangement of atoms in a triangular lattice occurs in the center of the trap, confirming the picture given by consideration of the local pair correlation function. It is important to note, however, that exchanges still take place, primarily in the peripheral disordered region, where the system goes from bulk crystal (center of the trap) to the arrangement on concentric rings characteristic of the surface region.

As mentioned above, the crystal to superfluid transition in the bulk is expected to take place roughly in the interval. The results obtained in this work show that a harmonically trapped system comprising at least =1000 particles ought allow one to observe experimentally some of the relevant physics. On the other hand, the behavior of systems comprising just few hundred particles is still largely influenced by the confining potential, in the same range of . A more quantitative discussion of the physics of the trapped system, and it connection to the bulk superfluid to insulator quantum phase transition, is offered in section V, based on the study of the superfluid properties of the system.

## Iv Superfluidity and condensation

Figure 6 shows a typical result for the superfluid and condensate fractions of a dipolar assembly. In this particular case, the number of particles is 100, and the trap strength yields a value of in the middle of the trap approximately equal to 0.22, as a function of temperature. As the temperature is lowered below , both the superfluid and condensate fractions become non-zero, saturating to very different values in the limit. Specifically, the superfluid fraction approaches unity, whereas the condensate fraction is less than 20% at . On the other hand, in the large limit (weaker trap strength) superfluid and condensate fractions are close in magnitude, both close to unity at low temperature, as expected for a more dilute system. That inter-particle interactions have a stronger (depleting) effect on the condensate than the superfluid fraction, is of course not a novel feature of this many-body system. The same trend is predicted theoretically Moroni and Boninsegni (2004b) and observed experimentally Glyde et al. (2011) in condensed helium.

It is well understood theoretically, that SF and BEC are both underlain by quantum-mechanical exchanges of identical particles Ceperley (1995). In order to explore this notion in greater detail, in Figure 7 we show the relative frequency of occurrence of exchange cycles involving a variable number of particles (here, ). These results are for four different values of the inter-particle distance in the center of the trap, namely at low temperature. In all four cases, the superfluid fraction is 100%, (within statistical uncertainties), whereas the condensate fraction goes from 54% for the the most dilute, to just 6% (i.e., a mere six particles) for the densest assembly. This difference in the fraction of the system condensed to the same single-particle quantum-mechanical state is directly related to the much greater frequency with which permutation cycles (of any length) occur at lower density, i.e., with weaker interactions. At the same time, these results clearly show that there is no obvious, direct numerical connection between the length of the longest observed exchange cycle, or the statistics of occurrence of cycles comprising specific numbers of particles, and the superfluid fraction. Remarkably, the superfluid signal in the low temperature limit is unaffected by the significant differences in the statistics of exchange cycles at the four thermodynamic conditions considered.

As explained in section II, we can compute local superfluid and condensate density, gaining greater insight into the relationship between SF and BEC. Figure 8 shows the total, superfluid and condensate radial density profiles, for a system of =100 particles and two different values of inter-particle separation in the center of the trap. The system is essentially in the ground state, i.e., the superfluid fraction in both cases. As one can see, the condensate density becomes suppressed relative to the superfluid one, as the inter-particle separation is decreased. Significantly, however, both local quantities are almost uniformly distributed throughout the system. For the superfluid density, this is consistent with recent findings for parahydrogen clusters Mezzacapo and Boninsegni (2008). The condensate density displays a markedly different beahvior in the vicinity of the surface, depending on the degree of particle localization. At the largest average inter-particle separation (lowest density), the system near the surface is highly dilute, and the local condensate fraction approaches 100%, as observed in helium droplets as well Lewart et al. (1988). As the density of the system is increased, particles near the surface become increasingly localized, and BEC is locally suppressed, while remaining finite near the center of the trap, where the physical behavior is that of a superfluid.

It is worth noting that the behaviour of the system of dipolar bosons studied here deviates significantly from that of hardcore bosons studied previously by DuBois and Glyde DuBois and Glyde (2001). In that study, it was found that the local condensate density becomes severely depleted at the center of the trap, as the density is increased. In particular, the condensate forms a shell at the edge of the trap for strongly interacting systems. This effect is absent here, a fact that can presumably be attributed to the long-range nature of the dipole-dipole interaction, which in turn renders the use of any local density approximation questionable.

As mentioned above and as will be quantitatively discussed in the next section, on decreasing the average inter-particle separation a crystalline phase eventually nucleates in the center of the trap, at which point both superfluidity and condensation become suppressed in the system. The absence of a sharp phase transition in a finite system results in a value of the superfluid fraction less than one in the ground state, as superfluid and crystalline phases coexist. In the next section, we determine the ground state phase diagram of the system by calculating the superfluid and condensate fractions as a function of system size and inter-particle separation.

## V Ground state phase diagram

As discussed in Sec. I, a =0 (quantum) phase transition between a crystalline, insulating phase and a superfluid phase is observed in the bulk system, driven by the chemical potential (i.e., the density). One of the issues of interest here is, whether any remnant of such a transition can be observed, either by simulation or possible experiments, in a finite system.
In Figure 9 (top), we show values of the superfluid fraction corresponding to different interparticle distances in the middle of the trap. These values are computed for various system sizes (between and ), at a temperature sufficiently low to be regarded as ground state estimates.

Although no sharp phase transition can be observed in a finite system, the data clearly show a crossover between two regimes, a superfluid () and crystal (), taking place roughly in the same region of in which it is observed in numerical studies of the ground state bulk phase diagram Büchler et al. (2007) (shaded region of Fig. 9). It is worth mentioning that, in the range , during a sufficiently long simulation the system switches back and forth between two different regimes, characterized by zero and finite superfluid density. This behavior, already reported in previous works Böning et al. (2008), is not unlike what observed in small clusters of parahydrogen Mezzacapo and Boninsegni (2006, 2007); it suggests that the two phases are energetically very close, and in practice renders lengthy simulations necessary, in order to assign a reliable value to the superfluid fraction. This is the reason for the (relatively) large statistical uncertainties of some of the estimates reported in Figure 9.

The results in Figure 9 (top) suggests that melting and freezing values for at =0 should lie in the range 0.055 0.07, which is consistent with bulk studies ^{1}^{1}1In particular, in two separate quantum Monte Carlo studies of the bulk system the transition was predicted to occur either in the region Büchler et al. (2007) or in the region Astrakharchik et al. (2007), in the latter case by extrapolating the system energy to the thermodynamic limit.. A precise determination of the coexistence region is obviously beyond the scope of what can be achieved in a study of a finite system.

Figure 9 (bottom) shows the condensate fraction for the same set of simulation parameters for which the superfluid fraction is given in the top part. Here, we define “condensate fraction” , i.e., the largest eigenvalue of the cylindrically averaged one-particle density matrix, divided by the total number of particles. In order to facilitate the comparison of results obtained for systems of different size, we have subtracted from all estimates of , as 1 is the lowest possible value of . Such a correction is non-trivial for the smallest systems.

In the thermodynamic limit, at any finite temperature, as thermal fluctuations destroy off-diagonal long-range order. In a finite system, however, can be finite, and indeed it is found to be finite for values of above 0.07, i.e., in the superfluid region. Remarkably, the data in Figure 9 (bottom) are consistent with the vanishing of concurrently with that of the superfluid fraction , as the system is compressed to the crystalline phase.

In an infinite system, the single-particle state in which bosons condense is a plane wave of momentum , simply based on considerations of translational invariance. This results in a (theoretical) -like peak at in . Experimental evidence of phase coherence throughout the system can be detected in the of a spatially confined system as well, through the appearance of a peak at at low temperature, even though harmonic confinement has the effect of broadening it, its width being of order 1/.

In experiment with trapped cold atoms, the momentum distribution is routinely measured via time-of-flight absorption imaging. As a guide for possible experiments, it is therefore instructive to compute the momentum distribution in the limit , for two systems which have a superfluid and a crystalline ground state.

For a two-dimensional cylindrically symmetric system, the momentum distribution is given by

(6) |

where is the spatially averaged one-body density matrix, computed by simulation, and is the zeroth order Bessel function Boninsegni (2009).
Figure 10 shows the computed for a system of =200 dipolar bosons at three different temperatures, for a trap strength yielding a value of in the middle of the trap, i.e., well into the superfluid region of the phase diagram. The development of a sharp peak at =0 as the temperature is lowered is clear; the value of the superfluid fraction at low approaches 100% in this case. No evidence of such a peak is observed in our simulations for values of for which the assembly of dipolar particles does not display a measurable superfluid response.

Thus, the results of this study suggest that much of the physics of the crystal-superfluid transition in dipolar bosons can be probed experimentally on a system comprising as few as several hundred particles.

## Vi Discussion and Conclusion

We have investigated by extensive Quantum Monte Carlo simulations the properties of dipolar bosons confined in two-dimensional traps. In particular we have analyzed the structure of the system, superfluidity and Bose-Einstein condensation, as well as the ground state phase diagram for the Hamiltonian (1) for different system sizes and as a function of in-trap density. We have considered system sizes up to , a number relevant to current and planned experiments. It is therefore our expectation that a direct comparison between theory and experiment should be possible.

We find that the physics of the system in the middle of the trap closely mimics that expected in an infinite system, at least in the region of density wherein the superfluid-insulator quantum phase transition is expected, if the number of particles in the trap is of order of a few hundreds. In particular, a transition between superfluid and crystal like phases is evident, as shown by our results for the density profiles, local pair correlation function as well as the superfluid and condensate fractions. By considering system sizes up to , we verify that such a transition occurs for an inter-particle separation of 0.055 0.07, which is consistent with the most reliable theoretical predictions for the bulk case Büchler et al. (2007).

In order to gain deeper understanding of the elusive relationship between superfluidity and Bose-Einstein condensation, a study of the global and local superfluid properties has also been performed. The general conclusion is consistent with that obtained for systems with different interactions, namely that strong inter-particle interactions have a much more pronounced quantitative effect on the condensate than on the superfluid fraction. However, our results are consistent with the concurrent disappearance of both as the system freezes at high density. With the exception of a narrow region near the edge of the trap, where BEC is either enhanced or suppressed depending on particle localization, superfluid and condensate density are largely uniform throughout the system.

In striking contrast from the case of hardcore bosons DuBois and Glyde (2001), BEC is not depressed at the center of the trap in the presence of strong interactions (i.e., at high density). Such a qualitative difference points to the importance of the long range nature of the interactions in (1), with the inference that a local density approximation cannot be straightforwardly applied here.

We have also investigated the momentum distribution as a function of temperature. The presence of coherence at low temperature, if the system features a superfluid ground state, is rendered evident by the appearance of a pronounced peak at zero momentum, in the low temperature limit. This seems to indicate that important signatures of a quantum phase transition predicted in the bulk system, should be readily observed in current or planned experiments.

## Acknowledgements

The authors thank Fabio Mezzacapo for useful discussions. This work was supported in part by the Natural Science and Engineering Research Council of Canada under research grant 121210893, and by the Alberta Informatics Circle of Research Excellence (iCore).

## References

- Leggett (2006) A. J. Leggett, Quantum Liquids (Oxford University Press, 2006).
- Feynman (1953) R. P. Feynman, Phys. Rev. 91, 1291 (1953).
- Davis et al. (1995) K. B. Davis, M. O. Mewes, M. R. Andrews, N. J. van Druten, D. S. Durfee, D. M. Kurn, and W. Ketterle, Phys. Rev. Lett. 75, 3969 (1995).
- Anderson et al. (1995) M. Anderson, J. Ensher, M. Matthews, C. Wieman, and E. Cornell, Science 269, 198 (1995).
- Krauth (1996) W. Krauth, Phys. Rev. Lett. 77, 3695 (1996).
- Heinrichs and Mullin (1998) S. Heinrichs and W. Mullin, J. Low Temp. Phys. 113, 231 (1998).
- Moroni and Boninsegni (2004a) S. Moroni and M. Boninsegni, J. Low Temp. Phys. 136, 129 (2004a).
- Boninsegni et al. (2006a) M. Boninsegni, N. Prokof’ev, and B. Svistunov, Phys. Rev. E 74, 036701 (2006a).
- Azuah et al. (1997) R. T. Azuah, W. G. Stirling, H. R. Glyde, M. Boninsegni, P. E. Sokol, and S. M. Bennington, Phys. Rev. B 56, 14620 (1997).
- Glyde et al. (2000) H. R. Glyde, R. T. Azuah, and W. G. Stirling, Phys. Rev. B 62, 14337 (2000).
- Griesmaier et al. (2005) A. Griesmaier, J. Werner, S. Hensler, J. Stuhler, and T. Pfau, Phys. Rev. Lett. 94, 160401 (2005).
- Lahaye et al. (2009) T. Lahaye, C. Menotti, L. Santos, M. Lewenstein, and T. Pfau, Reports on Progress in Physics 72, 126401 (2009).
- Jain and Boninsegni (2011) P. Jain and M. Boninsegni, Phys. Rev. A 83, 023602 (2011).
- Anderson et al. (1998) W. R. Anderson, J. R. Veale, and T. F. Gallagher, Phys. Rev. Lett. 80, 249 (1998).
- Pupillo et al. (2010) G. Pupillo, A. Micheli, M. Boninsegni, I. Lesanovsky, and P. Zoller, Phys. Rev. Lett. 104, 223002 (2010).
- Büchler et al. (2007) H. P. Büchler, E. Demler, M. Lukin, A. Micheli, N. Prokof’ev, G. Pupillo, and P. Zoller, Phys. Rev. Lett. 98, 060404 (2007).
- Astrakharchik et al. (2007) G. E. Astrakharchik, J. Boronat, I. L. Kurbakov, and Y. E. Lozovik, Phys. Rev. Lett. 98, 060405 (2007).
- Mora et al. (2007) C. Mora, O. Parcollet, and X. Waintal, Phys. Rev. B 76, 064511 (2007).
- Spivak and A. (2004) B. Spivak and K. S. A., Phys. Rev. B 70, 155114 (2004).
- Filinov et al. (2010) A. Filinov, N. V. Prokof’ev, and M. Bonitz, Phys. Rev. Lett. 105, 070401 (2010).
- Kurbakov et al. (2010) I. L. Kurbakov, Y. E. Lozovik, G. E. Astrakharchik, and J. Boronat, Phys. Rev. B 82, 014508 (2010).
- Cinti et al. (2010) F. Cinti, P. Jain, M. Boninsegni, A. Micheli, P. Zoller, and G. Pupillo, Phys. Rev. Lett. (2010).
- Lewart et al. (1988) D. S. Lewart, V. R. Pandharipande, and S. C. Pieper, Phys. Rev. B 37, 4950 (1988).
- Sindzingre et al. (1989) P. Sindzingre, M. L. Klein, and D. M. Ceperley, Phys. Rev. Lett. 63, 1601 (1989).
- Bagnato and Kleppner (1991) V. Bagnato and D. Kleppner, Phys. Rev. A 44, 7439 (1991).
- Hadzibabic et al. (2008) Z. Hadzibabic, P. Kruger, M. Cheneau, S. P. Rath, and J. Dalibard, New Journal of Physics 10, 045006 (2008).
- Lozovik et al. (2004) Y. Lozovik, S. Volkov, and M. Willander, JETP letters 79, 473 (2004).
- Boninsegni et al. (2006b) M. Boninsegni, N. Prokof’ev, and B. Svistunov, Phys. Rev. Lett. 96, 070601 (2006b).
- Mezzacapo and Boninsegni (2006) F. Mezzacapo and M. Boninsegni, Phys. Rev. Lett. 97, 045301 (2006).
- Mezzacapo and Boninsegni (2007) F. Mezzacapo and M. Boninsegni, Phys. Rev. A 75, 033201 (2007).
- Cuervo et al. (2005) J. E. Cuervo, P.-N. Roy, , and M. Boninsegni, J. Chem. Phys. 122, 114504 (2005).
- Kwon et al. (2006) Y. Kwon, F. Paesani, and K. B. Whaley, Phys. Rev. B 74, 174522 (2006).
- Löwdin (1955) P.-O. Löwdin, Phys. Rev. 97, 1474 (1955).
- Onsager and Penrose (1956) L. Onsager and O. Penrose, Phys. Rev 104, 576 (1956).
- (35) The statement made in Ref. Krauth (1996), that the number of particles in the condensate corresponds to the length of the longest observed exchange cycle, is incorrect.
- Fisher and Berker (1982) M. E. Fisher and A. N. Berker, Phys. Rev. B 26, 2507 (1982).
- Moroni and Boninsegni (2004b) S. Moroni and M. Boninsegni, J. Low Temp. Phys. 136, 129 (2004b).
- Glyde et al. (2011) H. R. Glyde, S. O. Diallo, R. T. Azuah, O. Kiricheck, and J. W. Taylor, Phys. Rev. B 83, 100507 (2011).
- Ceperley (1995) D. M. Ceperley, Rev. Mod. Phys. 67, 279 (1995).
- Mezzacapo and Boninsegni (2008) F. Mezzacapo and M. Boninsegni, Phys. Rev. Lett. 100, 145301 (2008).
- DuBois and Glyde (2001) J. L. DuBois and H. R. Glyde, Phys. Rev. A 63, 023602 (2001).
- Böning et al. (2008) J. Böning, A. Filinov, P. Ludwig, H. Baumgartner, M. Bonitz, and Y. E. Lozovik, Phys. Rev. Lett. 100, 113401 (2008).
- Boninsegni (2009) M. Boninsegni, Phys. Rev. B 79, 174203 (2009).