Exchange-induced crystallization of soft-core bosons

The role played by quantum-mechanical exchanges of indistinguishable particles in determining the fluid–solid phase boundary is a subject of fundamental interest in condensed matter and quantum many-body physics. It has long been the conventional wisdom that exchanges should have little or no influence over the freezing–melting phase transition. On its face, this assumption would seem reasonable; after all, in naturally occurring crystals quantum exchanges are strongly suppressed, both by particle localization at lattice sites as well as by the strongly repulsive core at the short distance of any known interatomic or intermolecular potential. Furthermore, the melting of most solids occurs at temperatures at which the average rms excursion of particles away from their lattice sites is dominated by thermal effects, with quantum-mechanical corrections being generally negligible (see, for instance, [1]). For this reason, it has been customary to neglect quantum statistics altogether in theoretical studies of quantum crystals also near the melting line.

Recent work [1] has challenged this assumption, however, by showing that in Bose systems with hard-core type interactions (such as condensed $^{4}$He), quantum exchanges have the effect of greatly expanding the region of stability of the fluid phase, with respect to what it would be if exchanges were not present, i.e. if particles were distinguishable (we henceforth refer to distinguishable quantum particles as 'boltzmannons'). To phrase this more quantitatively, the free energy gain associated with the occurrence of long cycles of permutation of identical particles has the effect of moving the freezing line to a considerably higher density than one would predict based on calculations only including the zero-point motion. For this reason, theoretical studies of the phase diagram of a Bose system neglecting exchanges are likely to incur significant quantitative error in the determination of the solid–fluid phase boundaries, and predict unphysical thermocrystallization (i.e. re-entrance of the solid phase) at finite temperatures [1]. Furthermore, long bosonic exchanges (i.e. comprising a macroscopic fraction of all particles in the system) can underlie and impart significant resilience (i.e. long lifetime) to metastable, glassy superfluid phases. Microscopically, this can be phrased in the language of path integrals, in terms of 'frozen' exchange cycles, in which the paths of many particles become entangled. Because a macroscopic number of single-particle (or, rare multi-particle) tunnelling events are required, in order to disentangle all particles, the system may remain 'stuck' in a metastable disordered, glassy superfluid phase [2].

In this work, we show that the findings of [1] crucially rely on the presence of a 'hard' repulsive core at short distances in the pair-wise interaction v(r). Indeed, if v(r) instead features a 'soft' core (i.e. $v\left( r\to 0 \right)\sim {{\hbar }^{2}}/2m{{d}^{2}}$, where m is the particle mass and d the mean inter-particle distance), the effect of the Bose statistics is in fact the opposite. Specifically, a high-density 'droplet' (or 'cluster') crystal phase, featuring a multiply occupied unit cell [3, 4], is strengthened over the fluid one, again with respect to the physics of a system of boltzmannons (or classical systems featuring the same kind of interaction, e.g. macromolecules [5–7] ⁵ ). Phrased differently, since the energy cost associated with particles laying at a close distance is relatively small, a phase in which each unit cell acts in a sense as a mesoscopic Bose condensate has a lower free energy than the uniform fluid phase ⁶ . We arrived at this conclusion by studying a two-dimensional system of soft-core bosons by means of quantum Monte Carlo simulations. Although the results presented here are for a specific kind of soft-core interaction [8–10], experimentally realizable in an assembly of cold atoms [11], the physics described here are independent of the detailed form of the potential utilized in our study, but only with the presence of a soft core at short inter-particle separation [4].

The remainder of this paper is organized as follows. In the next section, we describe the mathematical model utilized here. We then briefly outline the computational methodology, which is fairly standard and extensively documented in the literature, and illustrate our results in the following two sections. We outline our conclusions in the last section.

We consider here an ensemble of N spin-zero Bose particles of mass m, whose motion is confined to two physical dimensions (a choice made for convenience only, the main physical conclusions being independent of the dimensionality). The system is enclosed in a square cell of area A, with periodic boundary conditions, and is described by the following many-body Hamiltonian:

$$\begin{eqnarray}&&\hat{H}\,=-\frac{{{\hbar }^{2}}}{2m}\underset{i}{\mathop{\mathop{\sum }^{}}}\,\nabla _{i}^{2}+\underset{i<j}{\mathop{\mathop{\sum }^{}}}\,v\left( \left| {{\mathbf{r}}_{i}}-{{\mathbf{r}}_{j}} \right| \right)\end{eqnarray} \tag{ 1 }$$

The specific form of the potential utilized in this study is

$$\begin{eqnarray}&&v\left( r \right)=\frac{{{v}_{0}}}{r_{c}^{6}+{{r}^{6}}}\end{eqnarray} \tag{ 2 }$$

with ${{v}_{0}}>0$. Such a potential describes the interaction between two Rydberg atoms [12] in the so-called Rydberg blockaded regime [13, 14]. The above choice of interaction is motivated by the fact that a quasi-2D Bose assembly with such pair-wise potential can be experimentally realized in an assembly of cold Rydberg atoms [12, 15, 16], which feature strong van der Waals interactions [17]. These are currently utilized in numerous experiments to study long-range interacting effective spin systems [15, 18–21], as well as for applications in quantum optics [22–25] and quantum information science [26–28].

The most important feature of the potential (2) is the soft repulsive core of radius ${{r}_{c}}$, which is the main consequence of the 'Rydberg blockade' mechanism [13], causing a flattening off of the repulsive part at short inter-particle separation ⁷ . The rapidly decaying tail, also repulsive, does not play an important role in the context of this work; indeed, as mentioned in the introduction, the same qualitative behavior shown here can be observed with a broad class of physical potentials displaying a repulsive soft core at a short distance, with the only requirement being the presence of a negative Fourier component [5–7]. Another important feature is that the strength of the repulsive core of the pair-wise interaction can be 'tuned', allowing one to go from the soft to the hard core regime, in which qualitatively different physics arise.

We take ${{\epsilon }_{0}}={{\hbar }^{2}}/mr_{c}^{2}$ as our unit of energy (and temperature, i.e. we set the Boltzmann constant ${{k}_{B}}$ to one), and ${{r}_{c}}$ that of the length. Thus, the density of particles $\rho N/A$ is expressed in units of $r_{c}^{-2}$. The dimensionless parameter ${{V}_{0}}=m{{v}_{0}}/\hbar r_{c}^{4}$ measures the relative strength of the interaction compared with the characteristic kinetic energy ${{\epsilon }_{0}}$.

The phase diagram of this system is similar to that of other soft-core Bose systems [4]. If ${{V}_{0}}\to \infty$, the physics approaches that of an ensemble of hard disks, whose phase diagram features a low-density fluid (gas), turning superfluid at low temperature, and transitioning at sufficiently high density into a crystalline phase with one particle per unit cell (the presence of a weak repulsive tail stabilizes such a crystalline phase at t = 0, even at low density). Multiple occupation crystals occur at density $\rho \gtrsim 1$; in this regime, in which no supersolid phase is observed, the physics of the system in the solid phase mimic those of the Bose Hubbard model [29]. On the other hand, in a range of values of ${{V}_{0}}$ (roughly $1\lesssim {{V}_{0}}\lesssim 20$), at low temperature the system transitions from the fluid phase directly into a crystalline one featuring multiply occupied sites (unit cells) [30]. This is the physical regime of interest here.

We investigated the low-temperature phase diagram of the system described by the Hamiltonian (1) by means of first-principle computer simulations based on the worm algorithm in the continuous-space path integral representation [31]. This is a fairly well-established computational methodology, allowing one to essentially obtain the exact thermodynamics properties of Bose systems at finite temperature, using only the microscopic Hamiltonian as the input. Because the continuous-space worm algorithm is thoroughly illustrated elsewhere [32], we shall not review its implementation here. We only mention that the details of the calculations are standard, as the use of the potential (2) entails no particular difficulty. We utilized the usual fourth-order approximation for the high-temperature density matrix [33]. All of the results reported here are extrapolated to the limit of the zero time step.

The main quantity of interest here is the superfluid density, which we compute by means of the well-known 'winding number' estimator [34]. Most of the calculations for which results are shown here were carried out with a number of particles, of the order of a few hundred, with N = 800 being the largest system size considered. We carried out parallel simulations, at the same thermodynamic conditions, of a system of boltzmannons described by the same Hamiltonian, in order to assess the effect of Bose statistics on the phase diagram. It is worth remembering that the two systems have the same ground state; this is a straightforward consequence of the fact that the ground state wave function of a many-boson system is nodeless [35].

As mentioned above, the regime of interest in this work is that in which the repulsive core of potential (2) is soft, i.e. ${{V}_{0}}\lesssim 20$. In this range of repulsive interaction, at a density $\rho \gtrsim 1$, the system transitions from a fluid to a droplet crystal phase with a site occupation of the order of a few. In particular, considering a density $\rho \approx 1$ and with an incommensurate occupation number per site, it was recently shown as zero-point vacancies cause a superfluid flow of particles through the crystal [36], accordingly with the seminal works of Andreev–Lifshitz–Chester on supersolidity [37, 38]. However, in this study we shall consider a number of particles, 10 or above.

The occurrence of a specific phase, and in particular one that has crystalline long-range order, can be established in a computer simulation by calculating structural quantities like the pair correlation function, which displays marked oscillations in the crystalline phase. Equivalently, its Fourier transform, related to an experimentally measurable quantity known as the static structure factor, features a peak in the correspondence of the wave vector $k=2\pi /a$, with a being the lattice constant (typically $a\sim {{r}_{c}}$). However, the presence of crystalline order can also easily be assessed by visual inspection of particle world lines, an example of which is offered in figure 1, which clearly shows the formation of a droplet crystal for a system of Bose particles.

Zoom In Zoom Out Reset image size

A mean-field treatment [30] shows that the physics of the system in the ground state are governed by the single dimensionless parameter $\alpha {{V}_{0}}\rho$. We have verified, by direct numerical simulation, that this assertion holds quantitatively for both Bose and Boltzmann statistics at low temperature. To illustrate this point, in figure 2 we show the results for the superfluid fraction ${{f}_{S}}$ computed for the Bose system in the ground state limit (i.e. temperature $T\to 0$) for three different densities, namely 4.53, 6.78 and 11.33 (roughly corresponding to 10, 15 and 25 particles per unit cell, respectively). We show the results as a function of the renormalized interaction parameter α, in a range corresponding to the values of the interaction strength $3\lesssim {{V}_{0}}\lesssim 8$. Within the statistical errors of our calculations, the values of the superfluid density all fall on the same curve, with three different regimes clearly identifiable. Specifically, at low α the superfluid fraction is $\sim 1$, as the system is in the fluid phase; as α is increased to a value close to 28, ${{f}_{S}}$ abruptly drops to a lower (but finite) value, as expected for a superfluid system breaking transational invariance [39]. Finally, as α is increased even further (close to 37), the system transitions into an insulating droplet crystal phase, with negligible particle tunnelling across adjacent unit cells.

Zoom In Zoom Out Reset image size

Figure 3 shows the superfluid fraction ${{f}_{S}}$ for the Bose system at finite temperature, computed for $\alpha =28$, at which value the system displays a supersolid phase at low temperature. We consider here two values of the density, namely $\rho =9.02$ and 11.33; the simulated system comprises N = 800 particles. We plot the values of ${{f}_{S}}$ as a function of the reduced temperature $tT/\rho$, and observe the collapse of the data. As one can see, ${{f}_{S}}$ starts off at a value slightly less than 0.6 at t = 0, as the system is in the supersolid phase, and jumps up to a higher value at $t\approx 0.4$, in correspondence with the melting of the crystal into a uniform superfluid. The numerical data for ${{f}_{S}}$ in the superfluid regime are fitted in the usual way, based on Berezinskii–Kosterlitz–Thouless theory [40, 41], to obtain an estimate of the superfluid transition temperature.

Zoom In Zoom Out Reset image size

We construct a schematic phase diagram of both the Bose and the distinguishable particle systems in the $\left( \alpha ,t \right)$ plane through a number of vertical 'cuts' at different values of α. The resulting phase diagrams are shown in figure 4. The first obvious observation, aside from the fact that no superfluid phase can exist in a system of boltzmannons at any finite temperature, and that at exactly t = 0 the phase diagram is of course the same, is the much greater region of stability of the droplet crystal phase in the Bose system. This is the most remarkable outcome of this study, as it runs counter to the notion that quantum-mechanical exchanges should in principle favor a uniform phase in Bose systems. This is indeed what is observed in numerical studies of hard-core Bose systems [1], and it is what takes place in the system studied here as well, in the ${{V}_{0}}\gg 1$ (hard core) limit. However, if the repulsive core of the interaction is soft enough that a droplet crystal is present in the phase diagram, then quantum-mechanical exchanges of Bose particles actually strengthen the crystalline phase. Thus, Bose statistics stabilize the crystal phase at a higher temperature than in the system of distinguishable particles.

Zoom In Zoom Out Reset image size

The greater stability of the crystal in the Bose system can be understood in entropic terms. In a system of distinguishable particles, the crystal melts into a fluid due to the greater entropy of the phase with higher symmetry. On the other hand, in the Bose system at low temperature, the entropy of a normal fluid phase in which exchanges are only local in character is comparable to that of a crystal in which exchanges occur among particles confined to within the same unit cell (droplet). Thus, the thermodynamically preferred phase is that of lower energy, i.e. the crystal. As the strength of the repulsion increases, the number of particles in a droplet decreases and the crystalline phase becomes entropically less competitive with the fluid.

Summarizing, we have shown that in a Bose system characterized by pair-wise interactions with a repulsive core at short distances, quantum-mechanical exchanges can act to stabilize either the fluid or solid phase, depending on the strength of the repulsive interaction. While in the hard-core limit, exchanges strengthen the fluid phase, the opposite is true in a system in which the core is soft enough to allow the formation of a cluster (droplet) crystal phase at low temperature. We have shown this effect for a two-dimensional system of Rydberg atoms, but the result is quite general, and in particular is independent of the long-range part of the potential.

We thank Tommaso Macrì, Guido Pupillo and Nikolai Prokof'ev for valuable discussions. This work was supported in part by the Natural Science and Engineering Research Council of Canada. MB gratefully acknowledges the hospitality of the Max-Planck Institute for the Physics of Complex Systems, where the work was carried out.