Controlling coherence via tuning of the population imbalance in a bipartite optical lattice

The control of transport properties is a key tool at the basis of many technologically relevant effects in condensed matter. The clean and precisely controlled environment of ultracold atoms in optical lattices allows one to prepare simplified but instructive models, which can help to better understand the underlying physical mechanisms. Here we show that by tuning a structural deformation of the unit cell in a bipartite optical lattice, one can induce a phase transition from a superfluid into various Mott insulating phases forming a shell structure in the superimposed harmonic trap. The Mott shells are identified via characteristic features in the visibility of Bragg maxima in momentum spectra. The experimental findings are explained by Gutzwiller mean-field and quantum Monte Carlo calculations. Our system bears similarities with the loss of coherence in cuprate superconductors, known to be associated with the doping-induced buckling of the oxygen octahedra surrounding the copper sites.

We employ an optical lattice with two classes of wells (denoted as A and B) arranged as the black and white fields of a chequerboard. The optical potential is V (x, y) = −V 0 cos 2 (kx) + cos 2 (ky) + 2 cos θ cos(kx) cos(ky) (1) with the tunable well depth parameter V 0 and the lattice distortion angle θ. Adjustment of θ permits controlled tuning of the well depth difference ∆V ≡ 4 V 0 cos(θ) between A and B sites. In the special case θ = π/2 (∆V = 0) both types of sites are equivalent and hence a monopartite lattice arises, while in general the lattice is bipartite. Using the potential of Eq.
(1), we have numerically solved the Schrödinger equation for the single particle problem to obtain the exact band structure, including 14 bands in the plane-wave matrix representation of the Hamiltonian [1]. The single-particle problem is reformulated in terms of a tight-binding model Hamiltonian, where J is the hopping between neighboring sites of different sublattices, J A (J B ) is the hopping coefficient between neighboring sites of sublattice A (B) (see Supplementary Figure 1), and E A (E B ) is the on-site energy of sites belonging to sublattice A (B). We neglect the A → A hopping (henceforth indicated as J A ) along the diagonal lines of the lattice (same for B → B), because for the monopartite lattice (θ = π/2 or ∆V = 0) these hopping coefficients are exactly zero as a consequence of the symmetry of the Wannier functions. For sufficiently small deviations from θ = π/2, we expect that these coefficients are still negligible compared to J A or J B ; this assumption is supported by the full band structure calculation. For θ 0.53 π, this assumption becomes less reliable (see Supplementary Figure 2). By diagonalizing the Hamiltonian in Eq. (2) in momentum space and taking the lattice constant to unity, an analytic expression for the corresponding band structure (depending on the parameters E A , E B , J, J A , J B ) can be derived, When θ is tuned away from zero a gap opens, splitting the lowest band. We denote the two resulting bands by "1" and "2", with the corresponding energies E 1 (k x , k y ) and E 2 (k x , k y ). It is straightforward to verify that In order to determine the parameters of the model Hamiltonian (2), instead of calculating Wannier functions, we use these equations to adjust the tight-binding bands to the exact band structure calculation, finding reasonable agreement up to θ = 0.53 π, as shown in Supplementary Figure 2. The resulting values of the hopping coefficients and the energy difference E A − E B are plotted in Supplementary Figure 3. Since |J A | and |J B | are nearly two orders of magnitude smaller than J, we will neglect them in what follows, as long as J = 0.

SUPPLEMENTARY NOTE 2. MEAN-FIELD PHASE DIAGRAM OF THE BIPARTITE LATTICE MODEL
In this section, we summarize the calculation of the mean-field phase diagram of the bosonic model We restrict ourselves to nearest-neighbor hopping coefficients, which are the only relevant ones, as shown in Supplementary Figure 3. The Hamiltonian in Eq. (5) describes a bipartite lattice, in which one allows for different densities in the two sublattices. A similar problem has been discussed also in Ref. [2]. Extending a standard approach [3,4], we apply a mean-field decoupling of the hopping term [first term in Eq. (5)], with the order parameters ψ A,B ≡ a i(A,B) and the fluctuations δa i (A,B) . Neglecting the second order fluctuations of the fields, one finds where N A denotes the number of sites in the sublattice A. We use H J0 as a perturbation to the interaction part of the Hamiltonian (5), and neglect for the moment the irrelevant constant shift given by the last term in Eq. (7).
Since H J0 is local, the total Hamiltonian contains only local terms and we can apply perturbation theory in each unit cell. The unperturbed Hamiltonian H(J = 0) is diagonal with respect to the number operators and, hence, the eigenstates of H(J = 0) in each unit cell are |n A , n B , where n A and n B are the occupation numbers of the sites A and B, respectively. The energy per unit cell is given by The ground state corresponds to occupations g A and g B determined by the relations U ν (g ν − 1) < µ ν < U ν g ν , with ν = A, B. The first order contribution of the perturbation H J0 vanishes because H J0 does not conserve the number of particles, whereas the second order is found to be Including the previously ignored constant shift and using the fact that at zero temperature the calculated energy is the same as the Helmholtz free energy F , we can write with Here, z = 2d is the coordination number of the lattice; in our case d = 2 and z = 4. According to the (generalized) Landau criterion for continuous phase transitions, the phase boundaries are given by the condition Det[M] = 0. In the phase diagram shown in Supplementary Figure 4(a), one observes a series of lobes corresponding to Mott-insulator phases with occupation numbers that can vary in the two sublattices according to the value of the chemical potentials (see also Supplementary Figure 4(b), where the (g A , g B ) filling of the Mott lobes is explicitly given). Outside the lobes the system is superfluid.

SUPPLEMENTARY NOTE 3. EFFECT OF THE TRAP
We now discuss the effect of the additional harmonic trap potential. We set U A = U B = U , which is a very good approximation for θ 0.53 π. In Supplementary Figure 5, horizontal sections through the mean-field phase diagram are plotted for fixed values of V 0 . The lobes for µ B < 0 correspond to Mott phases with occupations (g A , g B ) = (g, 0), with g integer. For different values of θ, we also plot the lines L(θ) given by where ∆µ(θ) = E A − E B is the difference of the local energies E A and E B determined through Eq. (4). According to the local density approximation, one can define a local chemical potential with a maximal value in the center of the trap fixed by the total particle number, which decreases towards the edge of the trap. Hence, the phases encountered locally along a radial path pointing outwards from the trap center are given by the homogeneous phase diagram, when following the lines L(θ) towards decreasing values of µ A /U . The lines L(θ) shift to large, negative values of µ B /U as θ increases. As discussed in the main text, this means that the population of the B sites decreases and eventually vanishes. Hence, the density profile evolves into a wedding cake structure where only the A sites are populated, i.e., most atoms contribute to pure A-site Mott shells (g, 0) separated by narrow superfluid films, also with negligible B population (see Supplementary Figure 6 for an example of density profiles calculated with the Gutzwiller ansatz). The plot also shows that for increasing V 0 the Mott lobes (g, 0) cover an increasing area in the phase diagram, while at the same time the lines L(θ) shift towards lower values of µ B /U . This explains why the value of ∆V c , at which one observes a sudden loss of the visibility, reduces when V 0 is increased.

SUPPLEMENTARY NOTE 4. PERTURBATIVE RESULTS FOR THE VISIBILITY IN THE ASYMPTOTIC LIMIT.
The regime where the imbalance between A and B sites is large can be studied using perturbation theory up to second order [5] when the filling is chosen to be integer in the homogeneous case. In the limit where the hopping term is neglected (which is also the mean-field ground state), the ground state is given by a perfect Mott insulator of the form (g A , g B ) = (g, 0) Let us start from the first order term. The only non-vanishing terms are the ones for which a particle is removed from a site A and moved to one of the nearest-neighbor B sites. The energy difference is ∆ = U (g − 1) + ∆µ and the first order correction has thus the form The quadratic correction is such that a particle is removed from an A site, moved to a nearest-neighbor B site and from there it is transferred again to an A site which is different from the original one. The final A site can be a nearest-neighbor A site or a next-nearest-neighbor one. The correction becomes The ground state is therefore where the first term is simply the unperturbed term with a wave function renormalilzation.
We can now calculate the momentum distribution where N s is the number of unit cells in the system. The visibility V is calculated at momenta k max = (0, 0) and k min = ( √ 2π, √ 2π). Therefore, where r 1 ≡ cos( √ 2π) ≈ −0.266 and r 2 ≡ cos( √ 8π) ≈ −0.858, and we eventually find The visibility obtained in Eq. (21) is of the order 10 −1 in the highly imbalanced regime for filling between 2 and 3 (see Supplementary Figure 7(a)). The second order processes contribute significantly, as can be observed in Supplementary  Figure 7(b). In the theory just discussed, we did not include the contributions given by the bare next-nearest-neighbor hopping processes (J A ), despite the fact that the ground state (17) effectively includes this type of hopping contributions through virtual transitions. The reason is that the effective hopping processes contribute more substantially to the visibility than the bare ones (not displayed here). In Supplementary Figure 8, the experimental data for the visibility (extracted from Supplementary Figure 9) are plotted for V 0 = 10.8 E rec and V 0 = 11.44 E rec . The behavior of the visibility at large imbalance, where the system is deeply in a Mott insulator phase, can be described by Eq. (21), where the average fillingḡ has been used as a fitting parameter.

SUPPLEMENTARY NOTE 5. QUANTUM MONTE CARLO RESULTS IN 2D
Here, we describe the QMC procedure used to obtain the results for the critical values of the interaction at the tip of the g = 1, 2, 3 lobes in the phase diagram of the homogeneous Bose-Hubbard model for the monopartite lattice. We make use of the worm algorithm, as implemented in the ALPS libraries [6,7]. By measuring the superfluid stiffness, we are able to distinguish between the two phases of the homogeneous system on a square lattice. We use a finite-size scaling to determine the position of the critical point, keeping the product of the temperature T and the linear size L of the lattice constant [8]: T × L = 0.1 U . The comparison with a precise QMC calculation at zero temperature for filling g = 1 [9] allows us to conclude that the choice of temperature is adequate to describe the zero-temperature system. We study the system for different values of the ratio J/U , while keeping the chemical potential constant and equal to µ/U = 0.371, 1.427, 2.448, for g = 1, 2, 3, respectively. The choice of µ for g = 1 is comparable with the results in Ref. [9], while the choices of µ for g = 2 and g = 3 are taken from Ref. [10]. We set the maximum on-site occupation number to be g + 2 (and g + 3 for g = 1), thus allowing for more processes than just particle-hole excitations.
Our goal is to evaluate the lobe positions using a more reliable method than the usual mean-field approach [4]. We stress that a higher precision in the determination of the critical points, that could be obtained by lowering the temperature, increasing the system size and using a finer scan of the area around the lobe tip, is beyond the scope of this work, as it would not be relevant in the comparison with experimental results. For g = 1, 2, 3, we find the following values of (J/U ) c : 0.0597 ± 0.0001, 0.0352 ± 0.0001, 0.0250 ± 0.0001, where the errors are due to the use of a finite grid for J/U . These values are in agreement with a high-precision T = 0 result for g = 1 [9], and with estimates given in Ref. [10], based on the use of the effective potential method and Kato's perturbation theory (we observe that the latter values of (J/U ) c are systematically smaller than the ones we find). In Supplementary Figure 10, we show the finite-size scaling, as done in Ref. [8].

SUPPLEMENTARY NOTE 6. ANALOGIES TO HIGH-Tc SUPERCONDUCTORS
Our work may shed some light also on the behavior of similar condensed-matter systems, where loss of phase coherence occurs due to a structural modification of the lattice. One possible example are high-T c cuprates. Although the phenomenon of superconductivity occurs due to paired electrons, and here we are studying bosons, our system could bear some similarities with the cuprates, if one considers the scenario of pre-formed Cooper pairs at a higher temperature scale, as suggested by several theoretical and experimental works [11][12][13][14]. In this case, the onset of superconductivity at T c would correspond simply to phase coherence of the pre-formed "bosons".
The first discovered high-T c cuprate, La 2−x Ba x CuO 4 (see Supplementary Figure 11) was found to exhibit a dip in the critical temperature at the doping value x = 1/8. Later, the same phenomenon was shown to occur for La 2−x Sr x CuO 4 when La was partially substituted by some rare earth elements, like Eu or Nd [15]. This feature was long known as the 1/8 mystery, but further investigations of the materials have shown that it is connected to a structural transition from a low-temperature orthorhombic (LTO) into a low-temperature tetragonal (LTT) phase [16], see also Supplementary Figure 11. This structural transition corresponds to a buckling of the oxygen octahedra surrounding the copper sites, which changes the nature of the copper-oxygen lattice unit cell [16]. By increasing the concentration y of Nd in La 2−x−y Nd y Sr x CuO 4 , superconductivity is eventually destroyed. The onset for the disappearance of superconductivity depends also on the Sr doping x, but actually there is a universal critical angle θ c = 3.6 deg for the buckling of the oxygen octahedra, after which superconductivity cannot survive [17].
Until now, most of the theoretical studies of high-T c cuprates have concentrated on the 2D square copper lattice, but it is well known that the actual superconducting plane is composed of copper and oxygen forming a Lieb lattice, and that the dopants sit on the oxygen (see Supplementary Figure 11). The role of the LTO/LTT structural transition is mostly to shift two of the four in-plane oxygen atoms, which were slightly out of the plane, back into it. Although essentially more complicated than the problem studied here, the critical buckling angle θ c = 3.6 deg for the destruction of superconductivity [17] bears similarities with the critical deformation angle θ c (or equivalently ∆V c ) that we found in this work. We hope that our results will foster further investigations of the specific role played by the oxygen lattice in high-T c superconductors, and its importance in preserving phase coherence.