Phase diagrams of the disordered Bose-Hubbard model with cavity-mediated long-range and nearest-neighbor interactions

Recent experiments with ultracold atoms in an optical lattice have realized cavity-mediated long-range interaction and observed the emergence of a supersolid phase and a density wave phase in addition to Mott insulator and superfluid phases. Here we consider theoretically the effect of uncorrelated disorder on the phase diagram of this system and study the two-dimensional Bose-Hubbard model with cavity-mediated long-range interactions and uncorrelated diagonal disorder. We also study the phase diagram of the extended Bose-Hubbard model with nearest-neighbor interactions in the presence of uncorrelated diagonal disorder. The extended Bose-Hubbard model with nearest-neighbor interactions has been realized in the experiment using dipolar interaction recently. With the help of quantum Monte Carlo simulations using the worm algorithm, we determine the phase diagram of those two models. We compare the phase diagrams of cavity-mediated long-range interactions with nearest-neighbor interactions. We show that two kinds of Bose glass phases exist: one with and one without density wave order. We also find that weak disorder enhances the supersolid phase.


I. INTRODUCTION
The interplay between disorder and interaction attracts a lot of attention in condensed matter and statistical physics. A certain degree of disorder is ubiquitous in all condensed matter, but a thorough understanding of these systems is impeded by a poor control over the disorder and competing interactions. On the other hand, ultracold atoms, especially bosons in an optical lattice become an important way to simulate condensed matter systems [1][2][3][4][5][6][7][8][9][10][11][12][13]. In these experiments, interactions and disorder can be tuned independently. The short-range interaction can be realized using Feshbach resonances, while the long-range interactions have been studied using ultracold gases of particles with large magnetic or electronic dipole moments [14][15][16], polar molecules [17,18], atoms in Rydberg states [19][20][21], or cavity-mediated interactions [22,23]. Random potentials are usually produced using speckle patterns [2,24,25], while quasiperiodic potentials can be generated using bichromatic lattices [26]. Other possibilities to engineer disorder include the introduction of localized atomic impurities [4] and holographic techniques which produce point-like disorder [27].
Theoretically, a paradigmatic model to describe the interacting bosonic particles in an optical lattice is the Bose-Hubbard model (BHM). The BHM without disorder and only on-site repulsion features two phases: a superfluid (SF) phase and a Mott insulator (MI) phase. The so-called extended BHM includes nearest-neighbor and/or long-range interactions, as for instance dipolar interactions and cavity-mediated interactions. The phase diagram of the extended BHM with dipolar interactions has been calculated in Ref [28,29]. In contrast to dipolar interactions which decays as 1/r 3 , cavity-mediated long-range interactions are global, which means that the interaction strength between two bosons does not decay with the distance between them. The ground state phase diagram of the extended BHM with cavitymediated long-range interactions has been investigated extensively with the help of mean-field theory [30][31][32][33][34], Gutzwiller ansatz [35,36], quantum Monte Carlo [33,[36][37][38], Variational Monte-Carlo [39], and exact diagonalization [40,41] in 1D, 2D, and 3D. The results show that by adding cavity-mediated long-range interactions, the extended BHM exhibits a richer phase diagram with additional density wave (DW) and supersolid (SS) phases.
Introducing disorder into the standard BHM leads to the emergence of the gapless Bose glass (BG) phase, characterized by finite compressibility and absence of offdiagonal long-range order, always intervenes between the SF phase and MI phase [42,43]. The phase diagram of the disordered extended BHM with nearest-neighbor interactions was calculated for 3D [44,45], and the phase diagram of the disordered BHM with dipolar interactions was calculated for 2D [11]. However, the study of the extended BHM with cavity-mediated long-range interactions in the presence of disorder is still lacking. Whether the disordered potential enhances or suppresses the DW and SS phases here is still unknown.
In this paper, we use quantum Monte Carlo simulations based on the worm algorithm [46] to study the phase diagram of the two-dimensional Bose-Hubbard model with cavity-mediated long-range interactions and uncorrelated disorder. The paper is organized as follows: in section II, we introduce the Hamiltonian of the system with cavity-mediated long-range and nearestneighbor interactions. In section III, we discuss various phases and the corresponding order parameters. In section IV A, we present the phase diagrams of the 2D extended BHM with cavity-mediated long-range interactions and uncorrelated disorder. On the mean-field level the BHM with cavity-mediated interactions is identical to the BHM with nearest-neighbor repulsion -up to a renormalization of the chemical potential and the on-site potential [33]. Therefore, for comparison, we study in section IV B the phase diagram of the extended BHM with nearest-neighbor interactions in the presence of uncorrelated disorder. The extended BHM with nearestneighbor interactions was experimentally realized in [47]. Finally, section V concludes this paper.

II. HAMILTONIAN
In the following, we consider bosons trapped in an optical lattice with both short-range on-site and cavitymediated long-range interactions in the presence of disordered potential. The system is described by the Hamiltonian [32,37,38]: where the first term is the kinetic energy characterized by the hopping amplitude t. Here · · · denotes nearest neighboring sites on an underlying square lattice of linear size L with periodic boundary conditions, a † i (a i ) are bosonic creation (annihilation) operators satisfying the bosonic commutation relations. The second term is the short-range on-site repulsive interaction with interaction strength U s . Here, n i = a † i a i is the particle number operator. The third term is the cavity-mediated long-range interaction with interaction strength U l , the summations i ∈ e and j ∈ o denote summing over even and odd lattice sites respectively [32]. The fourth term is the chemical potential term with chemical potential µ shifted by the on-site random disordered potential ε i , where ε i is uniformly distributed within the range [−∆, ∆]. ∆ is the disorder strength. We set the unit of energy and length to be the hopping amplitude t. For each U l /t, U s /t, and ∆/t, we average over 100-200 realizations of disorder. On the mean-field level the BHM with cavity-mediated interactions is identical to the BHM with nearestneighbor repulsion -up to a renormalization of the chemical potential and the on-site potential [33]. In order to check this, we also consider the disordered 2D BHM with nearest-neighbor interactions defined by the Hamiltonian: Here, the first term is the kinetic energy with hopping amplitude t. The second term is the short-range on-site interaction with the interaction strength U s . The third term is the repulsive interaction with interaction strength U nn between bosons on nearest neighboring sites. The fourth term is the disordered potential term coupled with the chemical potential term. For each U nn /t, U s /t, and ∆/t, we average over 500-1000 realizations of disorder.

III. PHASES AND ORDER PARAMETERS
In this section, we list the phases we find in model 1 and 2 and the corresponding order parameters in table I. Each of the phases listed in table I corresponding to a unique combination of the order parameters. Here, three order parameters are needed to separate those quantum phases: superfluid stiffness ρ, structure factor S(π, π), and compressibility κ.  A non-vanishing superfluid stiffness ρ signifies offdiagonal long-range order, and it is easily accessible in QMC simulations using world line algorithms by calculation the winding number [48]. The superfluid stiffness is then given by: Here, W is the winding number. d is the dimension of the system and here, d = 2. L is the linear system size and β is the inverse temperature. The structure factor characterizes diagonal long-range order and is defined as: Here, k is the reciprocal lattice vector with k = (π, π) for the density wave with a checker board pattern and N = L × L is the system size.
The compressibility measures the density fluctuations and it is defined as:

IV. GROUND STATE PHASE DIAGRAMS
In this section, we present the ground state phase diagram for fixed particle density n i = 1 (note that in    this case the chemical potential in model (1) and (2) is superfluous) for cavity-mediated long-range interactions ( Fig. 1 (a)-(c)) and nearest-neighbor interactions ( Fig. 1 (e)-(g)). The x-axis is the on-site interaction U s /t and the y-axis is the disorder strength ∆/t, here we set the hopping amplitude t = 1. Table I summarizes the quantum phases in Fig. 1 and the corresponding order parameters: superfluid stiffness ρ, structure factor S(π, π), and compressibility κ. The phase boundary is determined by considering cuts through the x-axis (U s /t) and calculating the above three order parameters as a function of disorder strength ∆/t, as illustrated in Fig. 2

and 4.
A. Long-range interaction Figure 1 (a)-(c) shows the phase diagrams of the disordered BHM with cavity-mediated long-range interactions at interaction strength U l /t = 5, 10, and 16 at filling factor n i = 1. Without long-range interactions, the phase diagram of the disordered BHM at filling factor n i = 1 contains three phases: an SF phase, a MI phase, and a BG phase [8]. We use system size L = 16 and measure the three order parameters as a function of disorder strength ∆/t for various on-site interactions U s /t to determine the phase diagrams. Other system sizes have also been used to make sure the transition points are within the error bars. Figure 2 (a) shows the superfluid stiffness ρ and compressibility κ as a function of disorder strength ∆/t at U s /t = 20 for U l /t = 5. The structure factor S(π, π) is zero at fixed U s /t = 20 for different disorder strength. When the disorder strength ∆/t < 4, the system is in the MI phase with zero superfluid stiffness and zero compressibility. At disorder strength 4 < ∆/t < 6, the system is in the BG phase with finite compressibility but no superfluidity. As the disorder strength increases, the system goes to the SF phase at disorder strength 6 < ∆/t < 50. Finally, at large disorder strength ∆/t > 50, the superfluidity is destroyed and the system enters the BG phase. Figure 1 (a) shows the phase diagram at the interaction strength U l /t = 5. Compared with the phase diagram of the disordered BHM without long-range interactions, the shape of the phase boundaries at U l /t = 5 does not change but the region of the SF phase shrinks. For example, at U s /t = 20, the SF phase disappears around disor-  superfluid stiffness ρ (red circles), structure factor S(π, π) (orange triangles), and compressibility κ (purple rectangles) as a function of disorder strength ∆/t at Us/t = 20 for U l /t = 5 (a), Us/t = 20 for U l /t = 10 (b), and Us/t = 30 for U l /t = 16 (c).
der strength ∆/t ∼ 50, while for the disordered BHM without long-range interactions, the SF phase exists up to ∆/t ∼ 70 [8]. This is because the cavity-mediated long-range interaction tends to localize the particles in a 'checkerboard' pattern which suppresses superfluidity. Figure 1 (b) shows the phase diagram at the interaction strength U l /t = 10. Compared with Fig. 1 (a), the shape of the phase diagram boundaries does not change but with the SS phase emerges inside the SF phase at lower disorder strength. The SS phase has both diagonal long-range order and off-diagonal long-range order and is characterized by a finite superfluid stiffness ρ and a finite structure factor S(π, π). Figure 2 (b) shows the three order parameters as a function of disorder strength ∆/t at U s /t = 20 for U l /t = 10. At disorder strength 9 < ∆/t < 20, the system has a finite superfluid stiffness ρ and a finite structure factor S(π, π), implying that the system is in the SS phase. In the absence of disorder, at U l /t = 10, the DW to SF phase transition happens around U s /t ∼ 15.5 [39]. Interestingly, by adding disorder to the system, the DW phase is transformed to the SS phase at weak disorder strength and the SS phase exists even around U s /t ∼ 22, implying that weak disorder enhances the SS order. Here, disorder transfers a solid into a percolating supersolid [44,45] which is a percolating superfluid coexists with a solid.   Figure 2 (c) shows the three order parameters as a function of disorder strength ∆/t at U s /t = 30 for U l /t = 16. At disorder strength 16 < ∆/t < 32, the system has finite compressibility κ and finite structure factor S(π, π), which shows that the system is in the DS phase. We find that, at lower U s /t, we have the DW solid to DS phase transition first and then the DS to BG phase transition. This can be explained by the theory of inclusions [42,43], which states that a compressible glassy phase, the DS phase is surrounded incompressible phase, the DW phase. The glassy phase here is the DS since both of the DS phase and DW phase have finite structure factor. Figure 3 shows the density maps at fixed U l /t = 16 and ∆/t = 26 for on-site interaction U s /t = 25, 30, and 35, respectively. The radius of a red circle at a given site is proportional to the density at that site. At on-site interaction U s /t = 25, the density map shows clearly the density wave pattern. As the on-site interaction increases, the system losses the density wave pattern and the structure factor decreases. By further increases the disorder strength, the DS phase is destroyed in favor of the BG phase.
For model 1, the region of parameter space corresponding to U s /t < 18 has not been explored extensively. The reason is that the finite-size effects are much more pronounced at small on-site interactions, and at the same time, the long-range interaction term in the Hamiltonian slows down the worm update in the algorithm more and more when the system size increases such that systems sizes beyond L=16 are computationally inaccessible to us.

B. Nearest-neighbor interaction
In this subsection, we study the ground state phase diagram of model 2. On the mean-field level the BHM with cavity-mediated interactions is identical to the BHM with nearest-neighbor repulsion -up to a renormalization of the chemical potential and the on-site potential -and the mean-field phase diagrams are identical [33]. Here we demonstrate that the true phase diagram, for 2D in the presence of disorder, are only vaguely similar and actually shows significant differences.     shows the phase diagrams of the disordered BHM with nearest-neighbor interaction at the interaction strength U nn /t = 5, 7, and 10 at filling factor n i = 1, respectively. For the clean system without disorder, in the classical limit U s /t → ∞, the ground states are known [49][50][51]. At filling factor n i = 1, when the nearest-neighbor interaction and on-site interaction satisfy zU nn /U s < 1, the ground state is the MI state. While for zU nn /U s > 1, the ground state is the DW state. Here z is the coordination number and z = 4 in 2D.
Figure 4 (a) shows the superfluid stiffness ρ, structure factor S(π, π), and compressibility κ as a function of disorder strength ∆/t at U s /t = 16 for U l /t = 5 and L = 16. The structure factor decreases as the disorder strength ∆/t increases. At lower disorder strength ∆/t < 6, the system is in the DW phase with zero superfluid stiffness, zero compressibility, but finite structure factor. At disorder strength 6 < ∆/t < 12, the system is in the SS phase with finite superfluid stiffness and finite structure factor. Large disorder tends to destroy the DW order and the system goes in the SF phase at disorder strength 12 < ∆/t < 32. Further increasing the disordered potential results in the destruction of the SF phase in favor of the BG phase. Figure 1 (d) shows the phase diagram at interaction strength U nn /t = 5. When the on-site interaction U s /t > 20, there exists SF, MI, and BG phases, while when U s /t < 20, there exists DW, SF, BG, and DS phases. For the clean system, the SS phase exists when the on-site interaction strength U s /t < 5 [39]. At weak disorder strength, we find that the SS phase exists around 10 < U s /t < 18. Here, weak disordered potential enhances the SS phase [44,45]. Figure 4 (b) shows the three order parameters as a function of disorder strength ∆/t at U s /t = 28 for U nn /t = 7. At lower disorder strength 0 < ∆/t < 10, the system is in a region displaying phase separation (PS). At larger disorder strength ∆/t > 10, the SF phase appears. The emergence of the SF phase for increasing disorder is due to the formation of a percolating SF cluster as described in [45]. Finally, strong disorder destroys the SF phase in favor of the BG phase. Figure 4 (d) shows the density map of PS at U s /t = 28 and ∆/t = 6, where the DW phase is separated from the MI phase. Figure 1 (e) shows the phase diagram at the interaction strength U nn /t = 7. Compared with the phase diagram in FIG. 1  (d), as we increase the nearest-neighbor interaction, the SF phase shrinks and we find a region displaying PS at lower disorder strength around U s /t ∼ 28. Figure 4 (c) shows the three order parameters as a function of disorder strength ∆/t at U s /t = 30 for U nn /t = 10. At disorder strength 16 < ∆/t < 24, the system has finite compressibility κ and finite structure factor S(π, π), which shows that the system is in the disordered solid (DS) phase. Figure 1 (f) shows the phase diagram at the interaction strength U nn /t = 10. At this interaction strength, there is no SF phase anymore and all bosons are localized. Phase separation occurs around U s /t ∼ 40 at lower disorder strength. Interestingly, in addition to the DW and BG phases, a disordered solid phase emerges. At lower U s /t, the DW goes to DS and then BG phase as the disorder increases. The DS phase intervenes between the DW and BG phases since both the DW and DS phases have a finite structure factor.

V. CONCLUSION
Comparing the phase diagrams of model 1 and 2, we can see that the phase diagrams of the extended BHM with cavity-mediated long-range interactions and nearest-neighbor interactions only vaguely similar and actually display many significant differences. The main difference is that weak disorder leads to the PS in the nearest-neighbor interaction case. And in the weak disorder region, the phase diagram for the nearest-neighbor interaction changes around zU s /U nn ∼ 1, where the DW dominates for zU s /U nn < 1 while MI dominates for zU s /U nn > 1. There is no such change for the extended BHM with cavity-mediated long-range interaction. Here one always finds the DW phase for small U s and the MI phase for large U s . In conclusion, the phase diagrams for the disordered BHM with cavity-mediated long-range and nearest-neighbor interactions are vaguely similar but with significant differences in the size of phases and the existence of the region of phase separation.