Equilibrium Phases of Tilted Dipolar Lattice Bosons

The recent advances in creating nearly degenerate quantum dipolar gases in optical lattices are opening the doors for the exploration of equilibrium physics of quantum systems with anisotropic and long-range dipolar interactions. In this paper we study the zero- and finite-temperature phase diagrams of a system of hard-core dipolar bosons at half-filling, trapped in a two-dimensional optical lattice. The dipoles are aligned parallel to one another and tilted out of the optical lattice plane by means of an external electric field. At zero-temperature, the system is a superfluid at all tilt angles $\theta$ provided that the strength of dipolar interaction is below a critical value $V_c(\theta)$. Upon increasing the interaction strength while keeping $\theta$ fixed, the superfluid phase is destabilized in favor of a checkerboard or a stripe solid depending on the tilt angle. We explore the nature of the phase transition between the two solid phases and find evidence of a micro-emulsion phase, following the Spivak-Kivelson scenario, separating these two solid phases. Additionally, we study the stability of these quantum phases against thermal fluctuations and find that the stripe solid is the most robust, making it the best candidate for experimental observation.


I. INTRODUCTION
Experimental progress in trapping and controlling ultra cold atoms and molecules has led to the observation of magnetic and electric dipolar interactions in a variety of systems [1][2][3][4][5][6][7][8][9][10][11][12][13]. These systems can be used as quantum simulators to study quantum transport and dynamical and equilibrium properties of models featuring long-range and anisotropic dipolar interactions. Longrange dipolar interactions have been predicted to stabilize a plethora of exotic quantum phases such as pwave superfluids, superfluids of multimers, charge density waves, stripe solids, and supersolids. Moreover, dipolar interactions play an important role in many models of strongly correlated systems, excitons with spatially separated electrons and holes [14][15][16], and frustrated quantum magnets [17][18][19].
The recent success in creating a gas of polar molecules in an optical lattice is the first step towards the realization of a low-entropy initial state which is a key ingredient for any quantum simulator [7,20]. The ability to change the alignment of the dipole moment by applying an external field has sparked numerous studies on the many-body phases of dipolar gases. Many-body bosonic and fermionic dipoles with dipoles aligned perpendicular to the lattice or parallel to the lattice have been studied extensively [21][22][23][24][25][26][27][28][29]. However, the theoretical and numerical effort to study dipolar gases with dipoles aligned at arbitrary tilt angles has either focused on continuous systems [30][31][32][33][34][35][36] or has employed non-exact methods such as functional renormalization group [37], mean field theory and variational approaches [38][39][40].
In the following we present a systematic study of a system of hard-core, dipolar bosons trapped in a twodimensional lattice. We treat the alignment of the dipoles as a parameter which can be adjusted via the application of an external field. We use the worm algorithm [41], a type of Path Integral Quantum Monte Carlo method, and study the zero-and finite-temperature phase diagrams of the system at half-filling as a function of the tilt angle and the strength of dipolar interaction. We find that at low interaction strength the system is in the superfluid state for any value of the tilt angle. Upon increasing the interaction strength the superfluid phase is destabilized in favor of either a checkerboard or a stripe solid depending on the tilt angle. We find evidence of the Spivak-Kivelson scenario [23,[42][43][44] in the transition region between the two solid phases where an emulsion phase of the checkerboard and stripe solids is observed. Additionally, we study the robustness of these quantum phases against thermal fluctuations and show that the solid phases survive at temperatures higher than the critical temperature for the disappearance of the superfluid phase. In particular, due to the anisotropy of the dipolar interaction, the stripe solid turns out to be the most robust phase, making it the best candidate for experimental observation. We give predictions for actual experimental setups and temperatures required to observe solid phases. This paper is structured as follows: in section I we discuss the Hamiltonian describing the system. In section II we present the zero and finite-temperature phase diagrams. In section III we discuss how harmonic confinement affects our results. In section IV we explore possible experimental realizations and provide temperature estimates. Finally, section V concludes the paper.

II. HAMILTONIAN
We study a system of hardcore, dipolar bosons with induced dipole moment d, confined by a two-dimensional optical lattice with lattice constant a and by an external harmonic trap. A schematic representation of this setup, in the absence of the external confinement, is shown in figure 1. In this system, the dipoles are aligned parallel to each other along the direction of polarization which forms a tilt angle θ with an axis perpendicular to the plane of the optical lattice.
where the first term describes the kinetic energy of the system with hopping energy J and the second term is the dipole-dipole interaction with strength V = d 2 /a 3 and r ij = | − → r i − − → r j |. α ij characterizes the angle between the direction of the polarization and the relative position of the two particles given by r ij . Here, a † i (a i ) are the bosonic creation (annihilation) operators with the usual commutation relations and a †2 i = 0 , n i = a † i a i . We use . . . to denote nearest neighboring sites. Finally, µ i = µ− ξ=x,y w ξ ξ 2 i , where w ξ and ξ i are the strength of harmonic confinement and the coordinate of site i along axis ξ, respectively, and µ is the chemical potential which sets the total number of particles.

III. ZERO-AND FINITE-TEMPERATURE PHASE DIAGRAMS
In this section, we present the zero-and finitetemperature phase diagrams of the system described by Eq. 1 at half-filling and in the absence of an external harmonic confinement. The chemical potential is set to ensure the half filling condition, n = N/N sites = 0.5. The effect of the harmonic confinement is discussed in section III. Our results are based on Path Integral Quantum Monte Carlo (QMC) simulations using the worm algorithm (WA). Unless otherwise noted, we simulate systems with spacial dimensions L×L, where L = 24, 30, 36, and 42, with N sites = L × L. To extract the ground state phase diagram, we work at inverse temperature β = L/J which ensures that the system is effectively at zero temperature .
The main panel in Figure 2 shows the zero-temperature phase diagram in the V /J − θ plane which features three phases: a superfluid phase (SF), a checkerboard solid (CB) phase, and a stripe solid (SS) phase. The superfluid phase possesses off-diagonal long-range order and is characterized by finite superfluid stiffness ρ s , which can be extracted from simulations by measuring the winding number in space . The diagonal order in the solid phases is characterized by a finite value of structure factor S(k) = r,r exp[ik(r − r )] n r n r /N , where k is the reciprocal lattice vector. We use k = (π, π) and k = (0, π) to identify the CB and SS phases, respectively.
At low interaction strength, the system is in a superfluid phase for any value of the tilt angle θ. The superfluid phase is destabilized towards a solid phase as the interaction is increased above a critical value V c (θ)/J. Filled squares mark the onset of the CB solid, which forms at lower θ, while filled triangles mark the onset of the stripe phase which appears at larger θ.
Using energy considerations, it is easy to see that the anisotropic nature of the dipolar interaction leads to the stabilization of the CB solid at θ = 0 and the SS phase at θ = π/2. Moreover we observe that at larger values of the tilt angle the superfluid phase is less stable against increasing the dipolar interaction strength. Indeed, at a tilt angle θ = π/2, the dipolar interaction strength needed to destroy superfluidity in favor of a solid is about a factor of two smaller than what needed at a tilt angle θ = 0. At intermediate tilt angles, θ ∼ π/6, the competition between the two solid orders renders both unstable. This adds to the robustness of the superfluid phase, where a dipolar interaction V /J ∼ 6 is required in order to destroy the off-diagonal long range order.
We find the solid-SF transition points by studying finite size effects on ρ s and S(k). We determine the error bars for the solid to superfluid transition by monitoring the order parameters characterizing the two phases, namely ρ s and S(k). We assure that the order parameter characterizing the solid phase vanishes with increasing system size as we cross the phase boundary to the superfluid phase and vice versa. Given the nature of the interaction and the dimensionality of the system, a first order phase transition between two ordered states is not be allowed [23,[42][43][44] and a second order phase transition is highly unlikely. First we investigate the nature of the phase transition between the CB and SS phases V /J = 8. To this end we perform scans along the cut marked I on the main panel of figure 2. We find that the CB-SS phase transition at T = 0 is neither first-order nor second order and follows the Spivak-Kivelson micro-emulsion scenario as predicted in [23,[42][43][44] which is indicated by the shaded region in Fig. 2. We provide examples of the density map of the system in this region in figure 3(a). Here, each circle corresponds to a different site, and its radius is proportional to the local density. In this region the system features a micro-emulsion phase, where stripes of CB phase are embedded in the SS phase or vice versa.
We have also investigated the nature of the CB-SF phase transition by performing a similar scan along V /J = 4 line labeled as II, for a system with L = 100. In figure 3(b) we plot the structure factor S(π, π) as a function of θ using filled circles. The superfluid density ρ s is shown using filled triangles. We observe hysteresis in both S(π, π) and ρ s with a width of the hysteresis loop of only ∆θ ∼ 0.015. While this observation suggests that the system undergoes a weak first-order phase transition, this is not consistent with theoretical predictions stated  figure 2 at two values of θ shows evidence of micro emulsion phase. Here, each circle corresponds to a different site, and its radius is proportional to the local density. (b) The CB structure factor S(π, π) and superfluid density ρs as a function of θ along cut II. We find evidence of hysteresis which indicates that this transition behaves as a weak first-order phase transition for finite systems.
above. Moreover our inability to observe a micro emulsion phase intervening between the superfluid and the CB phase for system sizes up to L = 100 suggests that we may be far away from the thermodynamic limit. Hence if the micro-emulsion phase exists, the size of emulsion would be larger than this lattice size, and correspondingly large systems should be used in order to observe it experimentally.
Next we present an investigation of the robustness of the quantum phases described above against thermal fluctuations. As expected, we find the solid phases to be the most robust against thermal fluctuations. We have performed scans over V /J at two tilt angles θ = 0.271 and θ = 0.84 corresponding to SF-CB and SF-SS transitions at zero temperature. Our results for θ = 0.271 are summarized in figure 4. The phase boundaries are extracted from finite size scaling analysis of ρ s and S(π, π). We were not able to resolve the phase boundaries within the shaded region with the system sizes considered in this paper. Within this region we expect the system to undergo either a weak first order phase transition or feature a micro-emulsion phase at zero temperature. Upon increasing the temperature, thermal fluctuations destroy the SF phase in favor of a normal fluid via a Kosterlitz-Thouless (KT) transition [46]. In figure 4(b) we show ρ s as a function of T /J for L = 18, 24, 36, 42, 48, and 66 at V /J = 2, indicated as cut I in figure 4(a). In the thermodynamic limit, a universal jump is observed at a the critical temperature given by ρ s = 2mk B T c /π 2 . In a finite size system this jump is smeared out as seen in figure 4(b). We extract the critical temperature using the finite size scaling procedure described in [47]. The dashed line in figure 4(b) corresponds to ρ s = T /π and its intersection points with each ρ s vs. T /J curve are used to find T c as shown in figure 4(c). We find T c /J ≈ 0.66 at V /J = 2. The CB solid melts in favor of a normal fluid via a two-dimensional Ising transition. We use standard finite size scaling as shown in the inset of figure 4 where we plot the scaled structure factor S(π, π)L 2β/ν , with 2β/ν = 0.25 as a function of T /J for L = 18, 24, 36, 42, 48 along cut II, at V /J = 4. The crossing indicates a critical temperature T c /J = 0.812 ± 0.002.
We have performed a similar analysis at fixed tilt angle θ = 0.84 where the zero-temperature phase diagram features the SF and SS phases. The results are shown in figure 5(a). Critical temperatures are found using the same methods as described for the SF and CB phases. It is worth noting that while both the SS and CB phases are more robust against thermal fluctuations compared to the SF phase, at any given V /J the SS phase has critical temperature roughly twice that of the CB phase.
Finally, we look at the finite temperature properties of the CB-SS micro-emulsion phase at V /J = 8. In figure 5(b) we show S(k) for the CB and SS phases as a function of T /J at θ ≈ 0.48 and L = 100. We use squares and circles to show S(π, π) and S(π, 0), respectively. The filled and open symbols correspond to two different initial conditions chosen for the simulations at each T /J. Since the equilibrium density distribution of the microemulsion phase is affected by the choice of initial conditions we observe large fluctuations in S(k). These results are in agreement with what reported in [23]. We find that at T /J 1.2 the system is in the normal fluid phase, with both SS and CB phases disappearing. If T /J < 1.2 the micro-emulsion phase can be dominated by either solid order, depending on the initial conditions. Finite temperature properties of the CB-SS micro-emulsion phase at V /J = 8, θ ≈ 0.48, and L = 100. Squares and circles correspond to S(π, π) and S(π, 0), respectively. Filled and empty symbols refer to different initial conditions. At T /J 1.2 the system is in the normal fluid phase, while at T /J < 1.2, the micro-emulsion phase is dominated by either solid order, depending on the initial conditions.

IV. HARMONIC CONFINEMENT
In a typical experimental system particles are subject to the optical lattice potential as well as an external harmonic confinement. The effect of the harmonic confinement can be taken into account in the QMC simulations through a site-dependent chemical potential as shown in Hamiltonian (1).
Experimentally, using a sufficiently deep optical lattice suppresses double occupancy. In this limit, the hard-core condition is valid for all tilt angles. The variation of the chemical potential provides a scan over density hence resulting in coexistence of phases realized at half-filling (as Each circle corresponds to a different site, and its radius is proportional to the local density. Right panels: Statistics of the off-diagonal correlator a † ( ri)a( rj) (see text for detail). A brighter color corresponds to regions of the trap where offdiagonal long-range order is present. Top left and right panels: V /J = 2, θ = 1.34, N ∼ 511. A Mott insulator phase with filling factor one is observed at the center of the trap, while a stripe solid is stabilized in most of the outer shell. Bottom left and right panels: V /J = 6, θ = 0.27, N ∼ 1217. A superfluid phase is observed at the center of the trap, followed by a shell of CB solid which gives way to a superfluid outer shell once the harmonic potential is strong enough destroy the solid order. described above) with phases stabilized at other fillings. We use a lattice depth V 0 = 50E R , where E R is the recoil energy of the molecule and a harmonic confinement of ω x,y = 2π 5 Hz. Our parameters reflect typical experimental choices. The top and bottom left panels in figure 6 present the equilibrium density distribution of the dipoles, where each circle corresponds to a different site, and its radius is proportional to the local density. In the presence of the harmonic trapping potential we cannot use the usual winding number relations to measure the superfluid density in the different regions of the trap. However we can use the off-diagonal correlator a † ( r i )a( r j ) , to map the superfluid regions within the trap. The panels on the right-hand side of figure 6 show the regions in the trap where the correlator is longranged. Each panel is normalized such that the brighter points are those exhibiting long-ranged correlations and thus robust superfluidity.
The top left and right panels correspond to V /J = 2, tilt angle θ = 1.34 and N ∼ 511. The particles at the center of the trap form a Mott insulator (MI) with unit filling while a stripe solid is stabilized in most of the outer shell. The dark regions in the top right panel indicate an absence of off-diagonal long-range order, corresponding to the SS and MI phases. Given the steepness of the harmonic potential, we do not observe a superfluid shell separating the two insulating phases, however offdiagonal long-range order is present adjacent to the MI shell. This is due to the mismatch between the symmetry of the SS and the trap, which destabilizes the SS and allows for the build up of phase coherence in this region.
The bottom left and right panels of figure 6 corresponds to V /J = 6 and tilt angle θ = 0.27 for which the CB phase is stabilized at half filling. The total number of particles in the trap is N ∼ 1217. We observe the presence of superfluid phase at the center of the trap, followed by a shell of CB solid which gives way to a superfluid outer shell once the harmonic potential is strong enough destroy the solid order. The superfluid regions are characterized by local density different than 1 and 0.5 as the left panel shows, while off-diagonal long-range order is present as seen in the right panel. On the other hand, the left panel clearly shows a CB structure at density 0.5 in the inner shell with corresponding absence of off-diagonal long-range order (dark color in right panel).

V. EXPERIMENTAL REALIZATION
The system described in section I can be experimentally realized using bosonic polar molecules trapped in a two-dimensional lattice, provided that the lattice is deep enough to suppress double occupancy hence ensuring the hard-core condition. Once ground-state molecules are formed in the lattice their dipole moments can be aligned through the application of an external static electric field, realizing an interaction term of the form described in Eq. (1). While equilibrium states such as the ones described above have not been prepared yet, rapid experimental advances will allow for the realization of these states in the future. Our calculations suggest that, unless dipoles are tilted at intermediate angles θ ∼ π/6, for which superfluidity survives even for suppressed hopping strength J/V ∼ 0.17, solid phases would be the best candidate for experimental observation since they are more robust against thermal fluctuations. For example, for a system of 39 K 87 Rb polar molecules with permanent dipole moment d = 0.57D, trapped in an optical lattice with lattice constant a = 532 nm and depth V 0 = 40E r , with tunneling rate J/h = 104 Hz [7], one can achieve J/V ∼ 0.3. Hence, observation of the stripe solid would be possible at temperatures ∼ 10 nK. With the same experimental setup the emulsion phase can be observed at temperatures ∼ 6 nK. If we consider 23 Na 40 K with permanent dipole moment d = 2.7D [48] and same trapping parameters as above, one can achieve J/V ∼ 0.0136 which allow for observation of the CB solid phase and the SS phase at temperatures ∼ 109 nK and ∼ 250 nK respectively. It should be noted that similar behavior is expected in experimental systems featuring atoms and molecules with magnetic dipole moments, such as Er 2 , Cr, and Dy [1,9,13].

VI. CONCLUSIONS
In summary, we presented the zero-and finitetemperature phase diagrams of a system of hard-core, dipolar lattice bosons at half-filling as a function of the alignment of the dipole moments, characterized by the tilt-angle θ, and the strength of the dipolar interaction. At zero-temperature, the system features three phases: superfluid, checkerboard solid and stripe solid. The superfluid phase is present at all tilt angles, provided that the interaction strength is below V c (θ)/J, and upon increasing the interaction strength the system enters one of the two solid phases depending on the value of θ. We ob-served signatures of a micro-emulsion phase transition, observing an emulsion phase of the checkerboard and stripe solids in the region separating the two solid phases. We expect solid phases to be stabilized by lattice dipolar fermions as well as shown in [37]. Finally, we studied the robustness of these phases against thermal fluctuations and showed that the solid phases can be observed experimentally at temperatures up to ∼ 100 nK.

VII. AKNOWLEDGMENTS
This work was supported by the NSF (PIF-1415561, PIF-1211914 and PFC-1125844), AFOSR, AFOSR-MURI, NIST and ARO individual investigator awards. The computation for this project was performed at the OU Supercomputing Center for Education and Research (OSCER) at the University of Oklahoma (OU).