cluding both the repulsive and attractive interspecies in- teraction. We used the many-body version of Rayleigh- Schr¨ odinger perturbation theory in the kinetic energy term with respect to the ground state of the system when the ki- netic energy term is absent, and calculate ground state ener- gies needed to carry out our analysis. This technique was previously used to discuss the phase diagram of the single- species BH model [19–21, 23], extended BH model [24], and of the hardcore BH model with a superlattice [25], and its results showed an excellent agreement with Monte Carlo sim- ulations [23, 25]. Motivated by the success of this technique with these models, here we generalized it to the two-species BH model, hoping to develop an analytical approach which could be as accurate as the numerical ones. 

We derived analytical expressions for the phase boundary between the incompressible Mott insulator and the compress- ible superfluid phase up to third order in the hoppings. We als o proposed a chemical potential extrapolation technique based on the scaling theory to extrapolate our third-order power se- ries expansion into a functional form that is appropriate for the Mott lobes. In particular, when the interspecies interaction is sufficiently large and attractive, we found a re-entrant qua n- tum phase transition from paired superfluid (superfluidity o f composite bosons, i.e. Bose-Bose pairs) to Mott insulator and again to a paired superfluid in all one, two and three dimen- 

[1] M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. Fisher, Phys. Rev. B 40 , 546 (1989). 

[2] M. Greiner, O. Mandel, T. Esslinger, T. W. H¨ansch, and I. Bloch, Nature (London) 415 , 39 (2002). 

[3] T. St¨oferle, H. Moritz, C. Schori, M. K¨ohl, and T. Esslinger, Phys. Rev. Lett. 92 , 130403 (2004). 

[4] I. B. Spielman, W. D. Phillips, and J. V. Porto, Phys. Rev. Lett. 98 , 080404 (2007). 

[5] I. B. Spielman, W. D. Phillips, and J. V. Porto, Phys. Rev. Lett. 100 , 120402 (2008). 

[6] I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80 , 885 (2008). 

[7] A. B. Kuklov and B. V. Svistunov, Phys. Rev. Lett. 90 , 100401 (2003). 

[8] E. Altman, W. Hofstetter, E. Demler, and M. D. Lukin, New J. Phys. 5 , 113 (2003). 

[9] A. Kuklov, N. Prokof’ev, and B. Svistunov, Phys. Rev. Lett. 92 , 050402 (2004). 

[10] A. Isacsson, Min-Chul Cha, K. Sengupta, and S. M. Girvin, Phys. Rev. B 72 , 184507 (2005). 

[11] A. Arg¨uelles and L. Santos, Phys. Rev. A 75 , 053613 (2007). 

[12] C. Trefzger, C. Menotti, and M. Lewenstein, Phys. Rev. Lett. 103 , 035304 (2009). 

[13] Anzi Hu, L. Mathey, I. Danshita, E. Tiesinga, C. J. Williams, and C. W. Clark, Phys. Rev. A 80 , 023619 (2009). 

[14] P. Buonsante, S. M. Giampaolo, F. Illuminati, V. Penna, and A. Vezzani, Eur. Phys. J. B 68 , 427 (2009). 

[15] A. Hubener, M. Snoek, and W. Hofstetter, Phys. Rev. B 80 , 245109 (2009). 

[16] C. Menotti and S. Stringari, Phys. Rev. A 81 , 045604 (2010). 

sions. Since the available Monte Carlo calculations [9, 10] do not provide the Mott insulator to superfluid transition phas e boundary in the experimentally more relevant chemical po- tential versus hopping plane, we could not compare our results with them. This comparison is highly desirable to judge the accuracy of our strong-coupling expansion results. 

A possible direction to extend this work is to consider the limit where hopping of one-species is much larger than the other. In this limit, the two-species BH model reduces to the Bose-Bose version of the Falicov-Kimball model [28], the Fermi-Fermi version of which has been widely discussed in the condensed-matter literature and the Fermi-Bose version has just been studied [29]. It is known for such models that there is a tendency towards both phase separation and density wave order [30], which requires a new calculation partially similar to that of Ref. [24]. One can also examine how the momentum distribution changes with the hopping in the insu- lating phases [23, 31], which has direct relevance to ultracold atomic experiments. 

The author thanks Anzi Hu, L. Mathey and J. K. Freer- icks for discussions, and The Scientific and Technological R e- search Council of Turkey (T ¨ $_{UB}$˙ I TAK) for financial support. 

[17] J. Catani, L. De Sarlo, G. Barontini, F. Minardi, and M. Ingus- cio, Phys. Rev. A 77 , 011603(R) (2008). 

[18] G. Thalhammer, G. Barontini, L. De Sarlo, J. Catani, F. Mi- nardi, and M. Inguscio, Phys. Rev. Lett. 100 , 210402 (2008). 

[19] J. K. Freericks and H. Monien, Phys. Rev. B 53 , 2691 (1996). 

[20] T. D. K¨uhner and H. Monien, Phys. Rev. B, 58 , R14741 (1998). 

[21] P. Buonsante, V. Penna, and A. Vezzani, Phys. Rev. B 70 , 184520 (2004). 

[22] K. Sengupta and N. Dupuis, Phys. Rev. A 71 , 033629 (2005). 

[23] J. K. Freericks, H. R. Krishnamurthy, Y. Kato, N. Kawashima, and N. Trivedi, Phys. Rev. A 79 , 053631 (2009). 

[24] M. Iskin and J. K. Freericks, Phys. Rev. A 79 , 053634 (2009). 

[25] Itay Hen, M. Iskin, and M. Rigol, Phys. Rev. B 81 , 064503 (2010). 2 

[26] Recall that U ↑↓ cannot be greater than or equal to U ↑↑ U $_{↓↓}$, oth- erwise the mixture would be unstable against collapse. In addi- tion, see e.g. Fig. 7 in [13], where TEBD calculations show in one dimension that V ≲ − 0 . 06 U is already sufficient for the Mott insulator to paired superfluid transition. 

[27] Note that, unlike those of single-particle and single-hole exci- tations where dt σ is a constant when d → ∞ , in the case of two-particle and two-hole excitations, dt 2 σ must be kept con- stant when d → ∞ . In this respect, Eqs. (11) and (12) do not contain any finite- d correction at the second order in hopping. 

[28] L. M. Falicov and J. C. Kimball, Phys. Rev. Lett. 22 , 997 (1969). 

[29] M. Iskin and J. K. Freericks, Phys. Rev. A 80 , 053623 (2009); and see references therein. 

[30] S¸. G. S¨oyler, B. Capogrosso-Sansone, N. V. Prokof’ev, and B. V. Svistunov, New J. Phys. 11 , 073036 (2009). 

[31] M. Iskin and J. K. Freericks, Phys. Rev. A 80 , 063610 (2009). 

## VI. ACKNOWLEDGMENTS 

