Bosonic molecules in a lattice: unusual fluid phase from multichannel interactions

TL;DR

This paper demonstrates that multichannel interactions in ultracold bosonic molecules in optical lattices create an unusual fluid phase with suppressed phase coherence, altering the traditional phase diagram and potentially impacting related systems.

Contribution

It introduces the concept of a fluid phase induced by multichannel interactions, characterized by suppressed phase coherence despite number fluctuations, and uses multiple methods to analyze the phase diagram.

Findings

Identification of an unusual fluid phase between superfluid and Mott insulator.

Phase coherence can be arbitrarily suppressed in a two-site system.

Multichannel interactions are significant even without dipolar interactions.

Abstract

We show that multichannel interactions significantly alter the phase diagram of ultracold bosonic molecules in an optical lattice. Most prominently, an unusual fluid region intervenes between the conventional superfluid and the Mott insulator. In it, number fluctuations remain but phase coherence is suppressed by a significant factor. This factor can be made arbitrarily large, at least in a two-site configuration. We calculate the phase diagram using complementary methods, including Gutzwiller mean-field and density matrix renormalization group (DMRG) calculations. Although we focus on bosonic molecules without dipolar interactions, we expect multichannel interactions to remain important for dipolar interacting and fermionic molecules.

Equations

$$H=-t\sum _{⟨i,j⟩}{c}_i^{†}{c}_j+\sum _i\sum _{\alpha }{U}_{\alpha }{d}_{i,\alpha }.$$

$$H=\sum_{\alpha}\big{[}-2t(O_{\alpha}\ket{1,1}\!\bra{+_{\alpha}}+\text{h.c.})+U_{\alpha}\ket{+_{\alpha}}\!\bra{+_{\alpha}}\big{]}.$$

$$H=-2t(\eta c+\eta c^{†})+\sum _{\alpha }{U}_{\alpha }{d}_{\alpha }-\mu n$$

Full text

Introduction.

Over the last several years, multiple species of ultracold nonreactive molecules (NRMs) have been created Takekoshi et al. (2012); Köppinger et al. (2014); Molony et al. (2014); Takekoshi et al. (2014); Molony et al. (2016); Gregory et al. (2016); Park et al. (2015a, b); Will et al. (2016); Park et al. (2016); Wang et al. (2015); Guo et al. (2016), following breakthroughs creating reactive molecules Ni et al. (2008); Heo et al. (2012); Dutta et al. (2014); Deiglmayr et al. (2008); Gadway and Yan (2016); Moses et al. (2017). NRMs allow one to avoid the particle loss that occurs in reactive molecules and thereby study rich many-particle physics under conditions where molecular motion is relevant Carr et al. (2009); Lemeshko et al. (2013); Fedorov et al. (2016); Sundar and Mueller (2013); Gorshkov et al. (2011a, b), without limiting to regimes protected by the quantum Zeno effect Zhu et al. (2014), or where motion is completely frozen Yan et al. (2013); Wall et al. (2015).

However, recent work has shown that the interaction potential between NRMs is more complex than for atoms, with not a few but hundreds of relevant bound states Mayle et al. (2012, 2013); Croft and Bohn (2014); Croft et al. (2017). This significantly modifies the Hamiltonian for NRMs in an optical lattice from the usual Hubbard model Jaksch et al. (1998), even in the absence of dipolar interactions Doçaj et al. (2016); Wall et al. (2017a, b). While these papers derived the new model, they did not attempt to solve it or characterize its features. Hence, in light of this paper, predictions for the many-body physics of NRMs in optical lattices made in prior work that neglected multi-channel interactions (MCIs) must be revised to account for them, or at least justify neglecting them.

In this paper, we show that MCIs qualitatively alter the phase diagram of bosonic NRMs in an optical lattice, using the Hamiltonian previously derived in other work. We use a combination of analytic arguments, Gutzwiller mean-field theory, and density matrix renormalization group (DMRG) calculations to compute and understand the phase diagram of a model that captures the basic features of NRM collisions.

Model.

Bosonic NRMs in an optical lattice are described by Doçaj et al. (2016); Wall et al. (2017a, b)

[TABLE]

Here $\sum_{\langle i,j\rangle}$ sums over all nearest neighbors $i$ and $j$, and we define modified creation and annihilation operators $c_{i}^{\dagger}\ket{0}_{i}=\ket{1}_{i}$, $c_{i}^{\dagger}\ket{1}_{i}=\sqrt{2}\sum_{\alpha}O_{\alpha}\ket{2_{\alpha}}_{i}$, and $c_{i}^{\dagger}\ket{2_{\alpha}}_{i}=0$, where $\ket{n}_{i}$ is the $n$ molecule ground state on site $i$, $\ket{2_{\alpha}}_{i}$ are the two-molecule (doublon) eigenstates on site $i$, and $d_{i,\alpha}=\ket{2_{\alpha}}_{i}\bra{2_{\alpha}}_{i}$ counts doublons. There can be as many $\sim 1000$ relevant doublon states Mayle et al. (2013); Doçaj et al. (2016). The $O_{\alpha}$ and $U_{\alpha}$ are molecule-dependent parameters whose values were estimated in Refs. Doçaj et al. (2016); Wall et al. (2017a, b). The $O_{\alpha}$’s satisfy $\sum_{\alpha}O_{\alpha}^{2}=1$, and we order the channels $\alpha$ such that $U_{1}\leq U_{2}\leq\dots$. For simplicity, we neglect dipolar interactions, as is appropriate for zero electric field or for homonuclear molecules, although we expect interesting physics to persist with dipolar interactions. The only assumption for the form of the Hamiltonian that is essential is that the molecular interactions are sufficiently short-ranged compared to the Wannier function length. Under this assumption, an arbitrary multichannel Hamiltonian will have the form we utilized, and we expect this to hold for at least the bialkalis. This is detailed in Refs. Doçaj et al. (2016); Wall et al. (2017a, b).

Figure 1(a) shows a schematic of the processes contained in Eq. (1). A molecule can tunnel from one site to an empty site at rate $t$. Two molecules on two adjacent sites can tunnel onto the same site at a rate $tO_{\alpha}$. Unlike the usual Hubbard model, there are numerous doublon states, indexed by $\alpha$, even in the lowest band.

Two Sites.

We first consider the problem restricted to two sites with two molecules. Despite its simplicity, this calculation reveals the essential characteristic features of the MCI model that will appear in the many-site problem. It is also of direct relevance for understanding small numbers of molecules prepared in tunnel-coupled optical tweezers Kaufman et al. (2014), as has been recently proposed Liu et al. (2017). For this system, the only possible states are $\ket{1,1}=\ket{1}_{1}\!\ket{1}_{2}$ and $\ket{\pm_{\alpha}}=\frac{1}{\sqrt{2}}(\ket{2_{\alpha}}_{1}\!\ket{0}_{2}\pm\ket{0}_{1}\!\ket{2_{\alpha}}_{2})$. The antisymmetric states $\ket{-_{\alpha}}$ decouple and are irrelevant for the ground state. The Hamiltonian in the symmetric sector is

[TABLE]

This is an $(N_{c}+1)\times(N_{c}+1)$ matrix with $N_{c}$ the number of channels. To expose the qualitative phenomena, we choose the $U_{\alpha}$s to be separated by at least an order of magnitude, and the $O_{\alpha}$s to be weighted away from the lowest energy $\ket{2_{\alpha}}$.

Figure 1(b,c,d) shows the coherence $C=C_{1,2}$ where $C_{ij}\equiv\sqrt{\langle c_{i}^{\dagger}c_{j}\rangle}$ and number fluctuations $\sigma=\braket{n_{i}^{2}}-\braket{n_{i}}^{2}$ where $n_{i}=\ket{1}_{i}\!\bra{1}_{i}+2\sum\limits_{\alpha}\ket{2_{\alpha}}_{i}\!\bra{2_{\alpha}}_{i}$, in the ground state of Eq. (2) as a function of $t/U_{1}$ for $N_{c}=1,2,3$. Some features of the MCI model are similar to the one channel case (that is, the standard Bose-Hubbard model): As $t$ is increased, the number fluctuations and coherence both increase from zero to their maximum around $t\sim U_{1}$, a two-site precursor of the Mott insulator to superfluid phase transition.

However, two features qualitatively distinguish the MCI model’s behavior from the single channel’s. The most significant novel feature is the suppressed coherence in a broad intermediate range of the “phase diagram” from $t\sim U_{1}$ to $t\sim U_{2}$ in the shaded portion of Fig. 1(c). Another feature is a dip in the number fluctuations appearing around the crossover at $t\sim U_{2}$. These features are therefore direct consequences of the level structure, as that is the only difference in the models. We will see these survive in the thermodynamic limit.

Gutzwiller Mean Field Theory.

To understand Eq. (1)’s ground state beyond the two-site case, we employ two approaches. The first is Gutzwiller mean field theory (MFT), which approximates the ground state as a spatial product state Fisher et al. (1989); Rokhsar and Kotliar (1991). Its accuracy for the similar Bose-Hubbard model is discussed in Carr et al. (2010). The second is DMRG, which provides accurate results for one-dimensional systems Schollwöck (2011). We show the qualitative features found with the two-site calculation remain in the more accurate calculations.

The Gutzwiller MFT Hamiltonian is (adding a chemical potential $\mu$ to control the particle number)

[TABLE]

acting on the single site Hilbert space $\ket{0}$, $\ket{1}$, $\ket{2_{\alpha}}$. Here, $\eta$ is defined to be $C_{i,i+1}=\braket{c}$. For other lattices, including higher dimensions, replace $t$ with $zt/2$, where $z$ is the lattice coordination number. To find the MFT ground state of Eq. (3), we iteratively compute $\eta$ until self-consistency with $\langle c\rangle$ in the ground state is achieved 111Our convergence criterion is $||\eta_{m}|-|\eta_{m-1}||<p$ for $p=0.001$, where $\eta_{m}$ is the coherence for the $m$’th step of the iteration. The MFT captures the expected, sharp phase transition between the Mott insulator (small $t$) and superfluid (large $t$), which was smeared out in the two-site case due to its finite size. The calculated Mott insulator-superfluid phase boundary is similar to the one-channel case. Note that, because of our on-site particle number restriction, only the vacuum, $n=1$, and $n=2$ Mott lobes appear.

Figure 2 (left three columns) shows that the key features of the two-site MCI model persist to an infinite system approximated in the MFT. Namely, the intermediate region U_{1}\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{<}}t\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{<}}U_{2} displays a suppressed coherence relative to the single channel model (by $\sim 31\%$), whereas the number fluctuations are barely affected (\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{<}}4\%). The suppression in coherence can be made even larger by other choices of $O_{\alpha}$. There is a small dip in number fluctuations at the crossover $t$ and large $\mu$.

MFT reveals that the MCIs alter the phase boundary’s shape. Because the bandwidth of doublons (hopping on a background of singly occupied sites) decreases as $O_{1}$ decreases, the Mott lobe tilts up and to the right as $O_{1}$ decreases. This is similar to the renormalization of doublon tunneling rate that occurs in strongly interacting atoms, as discussed in Refs. Duan (2005); Hazzard and Mueller (2010); von Stecher et al. (2011); Büchler (2010); Wall and Carr (2012a, 2013), although there the matrix elements are affected through the band mixing rather than the MCIs.

Numerical solution in one dimension.

We use DMRG to calculate the density, $C_{L/3,2L/3}$ (chosen to avoid finite size effects), and site-averaged number fluctuations across the phase diagram in a one dimensional (1D) chain Wall and Carr (2012b); Jaschke et al. (2017). This provides highly accurate quantitative results, although some care was required in order to obtain converged results with the MCIs. We performed calculations for $L=30$ sites and 15 sweeps with discarded weight decreasing from $10^{-3}$ to $10^{-15}$, resulting in maximum bond dimensions of $\sim 350$. We estimate that these give $C$ and $\sigma$ to within $\sim 2$% in a $O_{1}=O_{2}=1/\sqrt{2}$ calculation. For the calculation in the figure, jaggedness is the result of the opacity imbalance being at the limit of what our numerical routines can reasonably handle. However, we are confident in the broad, qualitative features displayed due to its agreement with our MFT calculations within constant factors and accounting for distortion from the Legendre transform in $\mu$ at high $t/U_{1}$, and due to qualitative similarity to the cited case of known good convergence.

Fig. 2 reveals that the features of the MCI Hamiltonian found using the MFT and two-site approximations survive in the DMRG and are thus true properties of Eq. (1), at least in 1D. Specifically, there is an intermediate regime of suppressed coherence, a dip in number fluctuations, and a tilt of the Mott insulator/superfluid phase boundary.

Dependence on $O_{\alpha}$s and $U_{\alpha}$s.

So far, we chose $U_{\alpha}$s and $O_{\alpha}$s to illustrate the qualitative effects that emerge from the MCIs. However, the values chosen were not very realistic. Now we examine the dependence on these parameters and incorporate values consistent with expectations Doçaj et al. (2016); Wall et al. (2017b). The consequences of changing $U_{\alpha}$ are straightforward. We previously found a series of regions separated at values of tunneling $t\sim U_{\alpha}$. For example, in the two channel case, $C$ and $\sigma$ increase from zero at small $t$ to a plateau at $t\sim U_{1}$, and then at $t\sim U_{2}$ the coherence again increases to the value at which it saturates. Changing the $U_{\alpha}$ merely changes the locations of these crossovers, and if the $U_{\alpha}$s are not well-separated, the crossovers blend together.

The consequences of changing $O_{\alpha}$ are more intricate, as shown in Fig. 3. Increasing $O_{2}/O_{1}$ increases the $t/U_{1}$ at which $\sigma$ and $C$ initially increase to the intermediate plateau. The value of $C$ in the intermediate plateau decreases with increasing $O_{2}/O_{1}$, and $C\rightarrow 0$ as $O_{1}\rightarrow 0$. Increasing $O_{2}/O_{1}$ also causes the “turn-off” of $\sigma$ to occur at smaller $t/U_{1}$.

These effects can be understood by considering the following limit: $U_{1}\ll U_{2}$, $O_{1}\ll O_{2}$, and $U_{1}/O_{1}\ll U_{2}/O_{2}$. This limit is illuminating because it is contrived to separate the energy scales of the various features, but the features persist beyond this limit. We consider three regimes. (1) First consider the regime tO_{1}\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{<}}U_{1}. Here $tO_{2}\ll U_{2}$ by our assumptions, and the state $\ket{+_{2}}$ can be neglected. Thus, the Hamiltonian reduces to a single channel model with $U=U_{1}$. For $tO_{1}\ll U_{1}$, the state is $\ket{1,1}$ – a Mott insulator with $\sigma=0$ and $C=0$ – and beyond tO_{1}\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{>}}U_{1} it crosses over to $(\ket{1,1}+\ket{+_{1}})/\sqrt{2}$ with $\sigma=1/2$ and $C=O_{1}$. Note the $\sigma$ takes its maximum value, but the coherence is suppressed to $O_{1}$. (2) As $t$ increases further, $\ket{+_{2}}$ becomes relevant and it repels the $\ket{1,1}$ level. For $tO_{2}\ll U_{2}$ we can treat the term in Eq. (2) that couples $\ket{1,1}$ to $\ket{+_{2}}$ in second order perturbation theory to find $\Delta E_{1,1}=\frac{-(2tO_{2})^{2}}{U_{2}}$. Then the energy gap between $\ket{1,1}$ and $\ket{+_{1}}$ becomes $U_{1}+(2tO_{2})^{2}/U_{2}$. When $(tO_{2})^{2}/U_{2}$ becomes on the order of $tO_{1}$, i.e. $tO_{2}\sim O_{1}U_{2}/O_{2}$, the transition between $\ket{1,1}$ and $\ket{+_{1}}$ will be sufficiently off-resonant that the ground state is again simply $\ket{1,1}$. Therefore, $C$ and $\sigma$ return to zero. (3) Finally, for very large $t$ such that tO_{2}\mathrel{\hbox to0.0pt{\lower 4.0pt\hbox{\hskip 1.0pt\sim}\hss}\raise 1.0pt\hbox{>}}U_{2}, the coupling of $\ket{+_{2}}$ to $\ket{1,1}$ dominates any coupling to the $\ket{+_{1}}$ state, and we can work with the Hamiltonian projected to $\ket{1,1}$ and $\ket{+_{2}}$. Consequently, increasing $t$ beyond $tO_{2}\sim U_{2}$, the ground state becomes $(\ket{1,1}+\ket{+,2})/\sqrt{2}$ and thus $\sigma\rightarrow 1/2$ and $C\rightarrow 1$.

Now that we understand the effects of general $O_{\alpha}$ and $U_{\alpha}$, we consider typical values for them. These depend on as yet unknown (either theoretically or experimentally) species-specific molecular properties, as well as other parameters of the experiment. Hence, we will take a statistical approach, in which we sample parameters from an appropriate probability distribution informed by a combination of random matrix theory, transition state theory, and quantum defect theory Doçaj et al. (2016); Wall et al. (2017a, b). As a distribution that roughly captures the features of this approach, we take the energies from a uniform distribution from between $U_{\alpha}=0$ and $U_{\alpha}=1000$, and the $O_{\alpha}$ are sampled from a normal distribution with zero mean then normalized so that $\sum_{\alpha}O_{\alpha}^{2}=1$. The most artificial aspect of this choice is that we consider only positive energy states, although both negative and positive states are equally likely. This choice is mainly a convenience to preserve stability at zero temperature, but it may be qualitatively appropriate in some regimes and regardless is a first step towards a fuller understanding.

Figure 4 shows the effects of including realistic $U_{\alpha}$s and $O_{\alpha}$s. Rather than plotting the phase diagram, we parametrically plot $C$ versus $\sigma$ for varying $t$ at fixed $n=1$. We plot this way because the phase diagram – for example the Mott insulator-superfluid boundary – can fluctuate wildly due to trivial rescalings of the $U_{\alpha}$. By plotting parametrically we avoid the trivial rescalings: for example, in the single channel case this plot would be independent of $U$.

Figure 4 reveals that the MCIs have a strong effect on the phase diagram with remnants of the qualitative features that we have already identified. First consider the two-site, $N_{c}=5$, results shown in Fig. 4(a). The bottom-most curve is the single channel result, with the Mott insulator in the bottom left and the superfluid in the upper right. All of the MCI curves lie above and to the left of this, demonstrating the suppression of coherence at fixed $\sigma$ in the intermediate phase, just as in the earlier examples. For a typical sample, $\sigma$ can be enhanced by several tens of percent, and up to a factor of 4 in $\sim 10$% of our samples. The same behavior is seen in the DMRG calculations [Fig. 4(b)] with $L=15$. Finite-size effects suppress the sharp feature in $\sigma$ near $C\sim 0$, which is recovered for $C_{ij}$, $j\gg i$, and $L\to\infty$. Finally, Fig. 4(c) shows the two-site case for physically relevant $N_{c}$, where the bulk of the curves pull away from the $N_{c}=1$ curve. We expect this to persist in the many-body case, but it is not technically feasible to perform accurate DMRG calculations with this many channels.

Conclusions.

We demonstrated that multichannel interactions lead to characteristic features in the phase diagram of NRMs in an optical lattice, focusing on the simplest case of bosonic molecules without an applied electric field (including homonuclear ground state molecules Herbig et al. (2003); Danzl et al. (2008, 2010); Reinaudi et al. (2012); Stellmer et al. (2012)). The most striking difference with the usual one-channel Hubbard phase diagram is a large regime at intermediate $t/U$ between the usual Mott insulator and superfluid where coherence is suppressed while number fluctuations remain large. In this sense, the system acts as a normal fluid: particles move, but there is arbitrarily small off-diagonal order.

This work opens many questions. One major question is whether the steady state of Eq. (1) from any realistic initial condition can even support a superfluid. Another arises from our omission of channels with $U_{\alpha}<0$. The approximate probability distribution of $U_{\alpha}$, as given in Refs. Doçaj et al. (2016); Wall et al. (2017a, b), has equal likelihood of positive and negative $U_{\alpha}$. Negative values of $U_{\alpha}$ will cause the bosons to trivially clump together on a single site in the ground state in the thermodynamic limit Oelkers and Links (2007), although energy or kinetics suppressing triple occupancy can prevent this, as pointed out Diehl et al. (2010a, b) for the two-channel case Radzihovsky et al. (2004); Romans et al. (2004). Nevertheless, our calculations are relevant. First, experimental situations are possible where all the $U_{\alpha}$ are positive, especially when only a few channels are present. Second, at finite temperature, the system may be prevented from collapse, and the ground state phase diagram may point to features of this finite temperature phase diagram. In the future, it will be crucial to understand the phase diagram and dynamics of other systems with multichannel collisions, including fermionic, spinful, or dipolar molecules Lahaye et al. (2009); Baranov et al. (2012).

Acknowledgements.

We acknowledge Rick Mukherjee and Ian White for conversations. K.R.A.H. thanks the Aspen Center for Physics, which is supported by National Science Foundation grant PHY-1066293, for its hospitality while part of this work was performed. This work was supported with funds from the Welch foundation, Grant No. C-1872, and in part by the Data Analysis and Visualization Cyberinfrastructure funded by NSF under grant OCI-0959097 and Rice University. K.D.E. acknowledges support by a scholarship from the Physics and Astronomy Department at Rice University.