Quantum Field Theory for the Three-Body Constrained Lattice Bose Gas -- Part II: Application to the Many-Body Problem
QQuantum Field Theory for the Three-Body Constrained Lattice Bose GasPart II: Application to the Many-Body Problem
S. Diehl,
1, 2
M. A. Baranov,
1, 2, 3
A. J. Daley,
1, 2 and P. Zoller
1, 2 Institute for Quantum Optics and Quantum Information of the Austrian Academy of Sciences, A-6020 Innsbruck, Austria Institute for Theoretical Physics, University of Innsbruck, A-6020 Innsbruck, Austria RRC “Kurchatov Institute”, Kurchatov Square 1, 123182 Moscow, Russia
We analyze the ground state phase diagram of attractive lattice bosons, which are stabilized by athree-body onsite hardcore constraint. A salient feature of this model is an Ising type transition froma conventional atomic superfluid to a dimer superfluid with vanishing atomic condensate. The studybuilds on an exact mapping of the constrained model to a theory of coupled bosons with polynomialinteractions, proposed in a related paper [11]. In this framework, we focus by analytical meanson aspects of the phase diagram which are intimately connected to interactions, and are thus notaccessible in a mean field plus spin wave approach. First, we determine shifts in the mean field phaseborder, which are most pronounced in the low density regime. Second, the investigation of the strongcoupling limit reveals the existence of a “continuous supersolid”, which emerges as a consequence ofenhanced symmetries in this regime. We discuss its experimental signatures. Third, we show thatthe Ising type phase transition, driven first order via the competition of long wavelength modes atgeneric fillings, terminates into a true Ising quantum critical point in the vicinity of half filling.
PACS numbers: 03.75.Hh,03.75.Kk,11.15.Me,67.85.Hj,64.70.Tg
I. INTRODUCTION
It was recently recognised that two-body and three-body loss processes for bosons in an optical lattice couldgive rise to effective models involving two-body andthree-body hardcore constraints, respectively. The two-body case was observed in an experiment with Feshbachmolecules [1, 2], while it has been proposed theoreticallyto take advantage of strong three-body loss to createa three-body hardcore constraint in bosonic [3, 4] andfermionic [5] lattice systems. The mechanism behind theconstraint is that the dissipative process suppresses co-herent tunnelling processes that would create double ortriple occupation and lead to loss.A salient feature of a bosonic lattice gas with three-body onsite constraint is the possibility to tune it toattractive two-body interactions. The associated dimerbound state formation has a profound effect on the many-body system, resulting in an Ising-type quantum phasetransition from a conventional atomic superfluid to adimer superfluid with vanishing atomic order parame-ter but nonzero pairing correlation. The possibility ofobserving Ising type behavior in cold atomic gases hasbeen uncovered earlier by Radzihovsky et al. [6, 7] andRomans et al. [8] in the context of resonant Bose gasesin the continuum, i.e. at low densities. This, however,turns out to be challenging due to the poor stability ofthe molecular Bose gas close to the resonance [9]. Here,we encounter a weak coupling analog of this scenario onthe lattice, in which the stabilization of the system is pro-vided by the blockade mechanism leading to the 3-bodyhardcore constraint. Besides this feature, the presenceof the lattice leads to intriguing enrichments comparedto the continuum physics, as we will demonstrate in thispaper.The qualitative picture for the Ising transition can be obtained within a simple Gutzwiller approach, in whichthe three-body constraint is easily built in via choice ofthe ansatz wave function [3]. However, this treatmentleaves a number of questions unanswered, which arise onvarious length scales in the problem, and misses out im-portant – even qualitative – aspects of the phase diagramas we will show. On the microscopic scale, this concernsthe bound state formation, as well as the correct formof the effective theory for dimers in the strong couplinglimit. On the intermediate scales, relevant to the ther-modynamics, one may wonder to what extent the phaseborder obtained within the mean field is quantitativelyaccurate. Finally, a thorough analysis of the competitionof the long range low energy degrees of freedom is nec-essary to answer the question of the true nature of thephase transition. We note that all these effects are tiedto interactions, thus not available in a simple spin waveextension of the mean field theory.This paper is the second one of a sequence of two re-lated papers. In Ref. [11], we have developed a quan-tum field theoretical framework which makes it possibleto analytically address the above questions in two andthree spatial dimensions. It is based on an exact map-ping of the constrained lattice boson model to a coupledtheory of two unconstrained bosonic degrees with poly-nomial interactions. In the related paper [11], we haveconcentrated on the formal development of this mapping,and performed calculations in the “vacuum limits” cor-responding to zero and maximum filling n = 0 ,
2, whichare characterized by the absence of spontaneous symme-try breaking. In the present paper we apply this for-malism to the many-body problem. We concentrate onthe three interaction-related aspects of the many-bodyproblem mentioned above: First, we address the quan-titative question of shifts in the phase border, with theresult that they are pronounced at low densities, while a r X i v : . [ c ond - m a t . qu a n t - g a s ] A ug FIG. 1. Phase diagram for the attractive 3-body hardcoreconstrained attractive Bose-Hubbard model. The black curvebelow is the mean field result, while red and blue curve corre-spond to d = 2 , n = 0) and di-holes ( n = 2) takesplace at the red and blue crosses for d = 2 ,
3, determining theendpoints of the critical lines (Sec. III). A bicritical point,characterized by energetically degenerate but different orders(superfluid and charge density wave) is reached asymptoti-cally at half filling. It can be detected experimentally ramp-ing a superlattice (Sec. IV). An Ising quantum critical point,connecting the two ordered phases, is predicted in the vicinityof half filling, while the correlation length is large but finiteaway from this point (Sec. V). basically absent as the filling increases to its maximum n = 2. Second, making use of the perturbative resultsobtained in [11], we consider the many-body physics inthe strong coupling regime, and predict the existence of anew collective mode at half filling n = 1, whose presenceresults from a symmetry enhancement from the conven-tional phase rotation symmetry U (1) (cid:39) SO (2) exhibitedby bosonic systems to an SO (3) symmetry. We proposean experiment to test this scenario, exploring its conse-quences both analytically as well as using exact numericalmethods in one dimension. These studies lead us to callthe system in this regime a “continuous supersolid” – asupersolid with a tunable ratio between the superfluidand the charge density wave order parameters (cf. ananalogous phenomenon in magnetic [12, 13] and attrac-tive fermion [14] systems). Third, in a long wavelengthanalysis the phase transition turns out to be first orderfor generic fillings due to the Coleman-Weinberg mecha-nism [15]. This is in line with the low density continuumanalysis, which has been carried out in detail in [7]. Inour constrained lattice system, however, we find that theradiatively induced first order transition terminates intoa true Ising quantum critical point in the vicinity of halffilling, which connects the two ordered phases of atomicand dimer superfluid. Its origin may be traced back to azero crossing of the dimer compressibility together witha sequence of Ward identities, thus being protected bysymmetry. An estimate of the correlation length suggestsa broad domain of intermediate fillings 1 / (cid:46) n (cid:46) / II. QUANTUM FIELD THEORY FOR THEMANY-BODY PROBLEM
In this section we address two aspects which are partic-ularly relevant for the many-body physics and have notbeen discussed in [11]: The realization of Goldstone’stheorem in our constrained model, and the equation ofstate. To prepare for this discussion and set the notation,we review the construction of the quantum field theoryin Sec. II A, also making the paper rather self-contained.The reader familiar with the construction, and the readerwho is more interested in the physics results of this work,may jump this section.
A. Review of the Construction
The starting point for our analysis is the Bose-Hubbardmodel with a three-body onsite hardcore constraint, H = − J (cid:88) (cid:104) i,j (cid:105) a † i a j − µ (cid:88) i ˆ n i + U (cid:88) i ˆ n i (ˆ n i − , a † ≡ , (1)Here, a i , a † i are the bosonic creation and annihilation op-erators, J is the hopping matrix element µ the chemi-cal potential and U the onsite interaction energy. Thesummation in the first term is performed over near-est neighbors. Because of the constraint, the originalbosonic onsite Hilbert space is reduced to the three states | α (cid:105) , α = 0 , , | α (cid:105) = t † α,i | vac (cid:105) = ( α !) − / (cid:0) a † (cid:1) α | vac (cid:105) , (cid:88) α t † α,i t α,i = (2)from some auxiliary “vacuum” state | vac (cid:105) . The operatorsare not independent but obey a holonomic constraint asindicated above. The Hamiltonian in terms of operators t α reads H = − J (cid:88) (cid:104) i,j (cid:105) (cid:2) K (10) i K (10) † j + 2 K (21) † i K (21) j (3)+ √ K (21) i K (10) † j + K (10) i K (21) † j ) (cid:105) − µ (cid:88) i (ˆ n ,i + 2ˆ n ,i ) + U (cid:88) i ˆ n ,i , where K (10) i = t † ,i t ,i , K (21) i = t † ,i t ,i , ˆ n α,i = t † α,i t α,i . Note that in this representation of the constrained Hamil-tonian, the conventional roles of interaction and hoppingare reversed: while the interaction enters the quadraticpart of the Hamiltonian, the hopping term gives rise toeffective kinematic interactions. The representation istherefore ideally suited in a strong coupling limit.In a naive Gross-Pitaevski treatment of the Hamilto-nian, achieved by replacing the operators with complexvalued amplitudes t α,i → f α,i in Eq. (3), reproduces pre-cisely the Gutzwiller mean-field energy, i.e. a classicalHamiltonian field theory for spatially varying amplitudes f α,i , where the holonomic constraint (cid:80) i f ∗ α,i f α,i = 1 isthe normalization of the onsite wave function. We nowshow how one can introduce a convenient descriptionof the theory on the quantum level. To illustrate themethod we consider the case of vanishing density n → ≤ n ≤ t ,i operators in terms of t ,i and t ,i oper-ators using the constraint. Writing t ,i = | t ,i | exp i ϕ i ,we observe that the phase ϕ i is unphysical: it can beeliminated via a local redefinition of the remaining oper-ators, t α,i = t α,i e i ϕ i ( α = 1 , t as real, and replace t ,i = X / i , X i = 1 − ˆ n ,i − ˆ n ,i in K (10) i . Obviously, the square roots are impracticablefor any quantum field theory because they give rise tovertices of arbitrarily high order. To eliminate this prob-lem, we use the fact that the matrix elements of X / i and X i on our subspace are the same: either 1 or 0.Consequently, on the subspace we may replace K (10) i = t † ,i (cid:112) X i → t † ,i X i , X i = (1 − ˆ n ,i − ˆ n ,i ) , (4)and analogous for the hermitian conjugate. More for-mally, the replacement can be justified by noting thatthe constraint operator is a projection, X i = X i , andthat the Taylor representation for a function of such anoperator is f ( X ) = f (0)(1 − X )+ Xf (1) . With this im-plementation of the constraint, the remaining operators See [11] for a subtlety in deriving this formula. t , t can be treated as standard bosonic operators actingin a complete Hilbert space H = (cid:81) i H i , where H i = {| n i (cid:105)| m i (cid:105)} , n i , m i = 0 , , ... is a bosonic Hilbert spacefor “atoms” t and “dimers” t at each site i : | n i (cid:105) =( n i !) − / ( t † ) n i | (cid:105) i and | m i (cid:105) = ( m i !) − / ( t † ) m i | (cid:105) i . Theonsite Hilbert spaces H i can naturally be splitted intoa physical subspace P i with n i + m i = 0 or 1, andan orthogonal unphysical one U i with n i + m i > H i = P i ⊕ U i . Important for our construction is that theHamiltonian H has no matrix elements between physicaland unphysical subspaces, (cid:104) u | H | p (cid:105) = (cid:104) p | H | u (cid:105) = 0, where | p (cid:105) ∈ P = (cid:81) i P i and | u (cid:105) ∈ U = (cid:81) i U i , and, therefore, isblock diagonal, H = H P + H U . As a result, these sub-spaces do not mix during evolution, and all quantities,both dynamical and statistical, factorize. For example,for the partition function one has Z = Tr exp( − βH ) = Z P + Z U (5)= (cid:88) { p } (cid:104) p | exp( − βH P ) | p (cid:105) + (cid:88) { u } (cid:104) u | exp( − βH U ) | u (cid:105) . Consequently, if we find a way to discriminate betweenthe physical and unphysical contributions, we may indeedconceive the operators t , as conventional bosonic ones.Such a setting is provided by using the effective action to encode the physical information of the theory, see e.g.[17]. It is defined as the Legendre transform of the freeenergy W [ j ] = log Z [ j ] (we introduce a source term j =( j , j † , j , j † ) and use ˆ ξ = ( t † , t , t † , t )):Γ[ ξ ] = − W [ j ] + (cid:90) j T ξ, ξ ≡ δW [ j ] δj , (6)where the new variable ξ = (cid:104) ˆ ξ (cid:105) is the field expectationvalue or the “classical” field. The effective action has thefollowing representation in terms of a functional integral,exp − Γ[ ξ ] = (cid:90) D δξ exp − S [ ξ + δξ ]+ (cid:90) j T δξ, j = δ Γ[ ξ ] δξ . (7)where δξ ≡ ˆ ξ − ξ , and the Euclidean action S = (cid:82) τ (cid:80) i t † ∂ τ t + t † ∂ τ t + H [ t , t ]. The Hamiltonian nowis to be interpreted as a function for classical thoughfluctuating, time dependent fields. The last identity inEq. (7) is the full quantum equation of motion, and theequilibrium situation we are interested in is specified by j = 0, where no mixing between the physical and theunphysical sector occurs. Usually, the most general formof the effective action is only restricted by the symme-tries of the microscopic theory. Since, as shown above, nocouplings mapping from U ↔ P are generated, we haveidentified a means to distinguish physical vs. unphysi-cal contributions by writing down the most general formfor the effective action for the physical sector by directlyexcluding couplings which would violate this constraint.Now we generalize the procedure to arbitrary density.We first follow [18] but then apply our exact procedure.While we have so far replaced the t operator, whichgenerates the mean field vacuum state | Ω (cid:105) = (cid:81) i t † ,i | vac (cid:105) ,we now consider a more general mean field vacuum, | Ω (cid:105) = (cid:89) i (cid:0) (cid:88) α r α exp( iαφ ) | α (cid:105) i (cid:1) (8)= (cid:89) i (cid:0) (cid:88) α r α exp( iαφ ) t † α,i (cid:1) | vac (cid:105) ! = (cid:89) i b † ,i | vac (cid:105) . For site independent amplitude moduli r α , these statesallow for the description of homogeneous ground stateswith spontaneous phase symmetry breaking: If, e.g., all r α (cid:54) = 0, the requirement of a fixed spontaneously chosenoverall condensate phase φ requires the phase relation θ α = αφ ; this fixed phase relation is the manifestation ofspontaneous symmetry breaking in the Fock space. Wecan now introduce a new set of operators b † α ( α = 0 , , b † creates the mean field vacuum and will beeliminated. Such a transformation is performed via atwo-parameter unitary rotation, whose rotation anglesare chosen such that the new operators fluctuate aroundthe new vacuum state and do not feature expectationvalues, b † α,i = ( R θ R χ ) αβ t † β,i (9)with the explicit form of the rotation matrices R θ = cos θ/ θ/ φ − sin θ/ − φ θ/ , (10) R χ = χ/ − sin χ/ i φ χ/ − i φ cos χ/ . A finite θ ( χ ) corresponds to a finite amplitude in | (cid:105) ( | (cid:105) )and we will see below how these quantities are fixed viathe Goldstone theorem. The precise relation is r = cos θ/ , r = sin θ/ χ/ , r = sin θ/ χ/ . (11)At this point we can repeat the steps described above forthe case n = 0 in complete analogy. The constraint isimplemented via the replacement b ,i → X i ≡ − b † ,i b ,i − b † ,i b ,i , b † ,i b ,i → X i . (12)The second line is simply a rearrangement of the holo-nomic constraint. The resulting bosonic Hamiltonian,which is then quantized by means of a functional inte-gral, is rather complex, and we will analyze it below.However, it exhibits a simple structure, H = E GW + H SW + H int . (13) E GW is the Gutzwiller mean field energy and H SW de-scribes the quadratic spin wave theory . The corrections A linear contribution, as naively expected in the expansion aboutthe condensate, does not occur due to Goldstone’s theorem, seebelow. to the mean field phase diagram, as well as nontrivial ef-fects in the deep infrared physics which we are interestedin here, are not captured at this quadratic level. Theyare all encoded in the interaction part H int .Choosing the qualitative form of the ground state pre-requisites a certain knowledge about the physics of thesystem. Equipped with the right qualitative groundstate, we can then perform quantitative calculations be-yond the mean field level based on our mapping. Indeed,Eq. (13) suggests an interpretation of our constructionas an exact requantization procedure of the Gutzwillermean field theory. This is in complete analogy to the con-ventional treatment of e.g. bosonic continuum systemswith broken symmetries, where in a first step a certainorder parameter is chosen and the theory is expandedaround it. However, on the lattice the right choice ofthe qualitative features of the ground state might be lessobvious. For example, spatial modulations of the orderparameter are possible, such as exhibited by charge den-sity waves. This is easily incorporated in the formalism,and such a situation will be indeed encountered in Sec.IV.In sum, we have obtained the following simple result:supplying the most general form of the effective actionwith a constraint principle, the evaluation can proceed asin a standard polynomial boson theory. Similar to sym-metries, the restrictions on the full theory leverage overfrom the microscopic theory. Unlike symmetries, the rel-evance of the constraint depends on scale, being restric-tive on short distances, while on long distances powercounting arguments lead to an effectively unconstrainedthough interacting spin wave theory with two degrees offreedom (see below). In practical computations, we canevaluate a theory of standard coupled bosonic fields. Thisopens up the powerful toolbox of modern quantum fieldtheoretical methods for calculations in onsite constrainedmodels. B. Goldstone’s Theorem and the ConstraintPrinciple
We will derive Goldstone’s theorem from a compari-son of the full quantum equation of motion and the effec-tive potential. The latter is defined as the homogeneouspart of the effective action, obtained by inserting tem-porally and spatially homogeneous field configurations U = Γ[ b † , b , b † , b ] /V d +1 ( V d +1 = M d /T is the quantiza-tion volume, M the number of lattice sites in each latticedirection). The possible dependences of the effective ac-tion and potential on the fields are strongly restrictedby both symmetry and constraint principle; the effectivepotential is further limited by the requirement of homo-geneity. We will show here that the constraint leads to anadditional U (1) invariant on which the effective potentialmay depend with no counterpart in unconstrained theo-ries, but we will also demonstrate that it does not breakthe validity of Goldstone’s theorem – in line with the in-tuition that the microscopic constraint would not affectthe long wavelength physics too strongly. Note, however,that the constraint has an impact on the long wavelengthphysics, as it is indirectly responsible for the presenceof the Ising quantum critical point close to unit filling n = 1 (cf. Sect. V B). Therefore, a thorough discussionof Goldstone’s theorem seem adequate.Let us construct the most general dependence of theeffective potential on the variables b , b . For simplic-ity of the presentation, we focus on a spontaneouslybroken symmetry for the dimers ( θ (cid:54) = 0 , π ), while theatoms are in the normal phase ( χ = 0 , π ). The latterfield can therefore be excluded from the following con-siderations. There are two possible terms associated tothe original t degree of freedom that might appear inthe effective action: Either it appears as a local com-bination ˆ n ,i = t † ,i t ,i , or as a bilocal (in general, n -local) combination, such that the constraint has to betaken into account via proper combination with t , e.g. t † ,i t ,i . While the local combination respects the U (1)symmetry, the second term must appear with a con-jugate partner as t † ,i t ,i t † ,j t ,j . In order to implementthe finite density, we now apply the rotation prescription t ,i = sb ,i + cb ,i , t ,i = cb ,i − s ∗ b ,i and subsequentlyimpose the constraint b ,i → X i . (Here and in the follow-ing, we abbreviate s = sin θ/ − φ , c = cos θ/ U ( ρ, λ † λ ) , (14) ρ = ( sX + cb † )( sX + cb )= s + cs ( Xb + b † X ) + ( c − s ) b † b − s b † b ,λ = ( sX + cb † )( cX − sb )= cs + c b † X − s Xb − csb † b − csb † b , where X, b denote the zero momentum and frequencycomponents of these field expressions, and without loss ofgenerality we have chosen s real. Note that neglecting theconstraint by setting X →
1, and considering low densi-ties, s → θ/ , c →
1, we recover the standard quadraticform for the condensate density from the local combina-tion, ρ = ˆ b † ˆ b , with ˆ b = s + b . However, the constraintprinciple requires a more complicated form of ρ , as well asthe account for a second invariant λ † λ . In the following,we will be concerned with first and second derivatives ofthe effective potential with respect to b , b † , which areevaluated at the physical point b = b † = b = b † = 0.Thus, we may set X → , b , b † → U can be further restricted. For b ,i ( τ ) , b † ,i ( τ ), we intro- duce the basis of hermitian fields σ i ( τ ) , π i ( τ ), b ,i ( τ ) = 1 √ (cid:16) σ i ( τ ) + i π i ( τ ) (cid:17) , (15) b ( q ) = 1 √ (cid:16) σ ( q ) + i π ( − q ) (cid:17) , which as the original fields do not carry expectation val-ues. Here we have used the Fourier conventions b ,i ( τ ) = (cid:90) q e i qx i b ( q ) , b † ,i ( τ ) = (cid:90) q e − i qx i b † ( q ) , (16) x i = ( τ, x i ) , q = ( ω, q ) , (cid:90) q = (cid:90) dω π (cid:88) q . We calculate the local combination in terms of these op-erators, ρ = s + ( c − s ) ( σ + π ) + √ csσ, (17)( σ = σ ( q = 0) , π = π ( q = 0)) and we observe that λ † λ ,to the relevant quadratic order, can be written as λ † λ = s + ( c − s ) ρ − cs ) σ . (18)Thus, the most general dependence of the effective po-tential on the homogeneous fields σ, π is given by U ( ρ, σ ) . (19)Now we study the mass matrix, which can be calcu-lated from the effective potential as the second deriva-tive with respect to σ, π . In particular, the form of theeffective potential implies for the π mass or gap ∂ U ( ρ, σ ) ∂π∂π (cid:12)(cid:12)(cid:12) σ = π =0 = (cid:16) ∂ ρ∂π U (cid:48) + (cid:16) ∂ρ∂π (cid:17) U (cid:48)(cid:48) (cid:17)(cid:12)(cid:12)(cid:12) σ = π =0 (20)(primes denote derivatives w.r.t. the invariant ρ ) with ∂ρ∂π = 0 , ∂ ρ∂π = c − s . (21)To complete the derivation of Goldstone’s theorem, wecalculate the equation of motion for σ from the effectiveaction, but immediately specialize to the case of homo-geneous fields, δ Γ δσ i ( τ ) (cid:12)(cid:12)(cid:12) hom = ∂ U ( ρ, σ ) ∂σ = 2 σ ∂ U ( ρ, σ ) ∂σ + √ cs U (cid:48) ! = 0 . (22)By construction we have σ = 0 as the solution of theequation of motion, and furthermore in the presence ofspontaneous symmetry breaking cs (cid:54) = 0. Thus U (cid:48) = 0 . (23)This simple relation indicates the presence of the gap-less Goldstone mode: The π gap calculated in Eq. (20)vanishes due to Eq. (23), and since ∂ρ/∂π = 0, cf. Eq.(21). This property is protected by the U (1) symmetryof the problem. Though the form of the effective poten-tial is more complicated than in the continuum at lowdensities, where the effective potential depends only onthe low density limit of the invariant ρ , we can explicitlyprove Goldstone’s theorem.Note, that the equation of motion (22) for σ also ex-cludes any homogeneous linear term in this field. Thesame is true for the π mode. Such terms would not becompatible with the equilibrium condition (22). This ex-cludes nonzero couplings from the homogeneous terms ∼ σ, π or b , b † in the effective action. Furthermore, viaEq. (14) the linear terms b , b † are connected by the con-straint principle to cubic terms: only the combinations ∼ b X, Xb † occur in the effective potential. Thus, com-bining Goldstone’s theorem and the constraint principle,we see that the cofficients of the terms b X, Xb † mustvanish, i.e.: δ Γ δb † ,i δb ,i δb ,i (cid:12)(cid:12)(cid:12) hom = δ Γ δb † ,i δb ,i δb ,i (cid:12)(cid:12)(cid:12) hom (24)= δ Γ δb † ,i δb ,i δb † ,i (cid:12)(cid:12)(cid:12) hom = δ Γ δb † ,i δb ,i δb † ,i (cid:12)(cid:12)(cid:12) hom = 0 . In the presence of an atomic condensate χ (cid:54) = 0, anal-ogous equilibrium conditions can be derived for the b field.In the symmetric phases n = 0 , θ = 0 , π ), no distinc-tion between the phase and the amplitude mode appearsand the mass matrix is degenerate. In this case, Gold-stone’s theorem reduces to the condition for the existenceof a dimer/di-hole bound state. C. Equation of State
The equation of state is obtained as the average overthe particle number operator ˆ N = (cid:80) i ˆ n i , ˆ n i = ˆ n ,i +2ˆ n ,i . Thus, after rotation we have for the particle den-sity n = (cid:104) ˆ N (cid:105) /M d = 2 | s | + ( c − | s | )[ (cid:104) b † b (cid:105) + 2 (cid:104) b † b (cid:105) ]+2 c [ s ∗ (cid:104) Xb (cid:105) + s (cid:104) b † X (cid:105) ] = − ∂ U ∂µ . (25)The second equality results from the path integral rep-resentation of the effective action, and is due to the cou-pling − µ ˆ N in the microscopic action. The connected two-point functions are given by the traces of the full Green’sfunctions for b and b . A convenient shorthand to re-late the connected Green’s functions to its one-particleirreducible counterpart is ξ = ( b † , b , b † , b ), (cid:104) b † b (cid:105) = (cid:104) ξ ξ (cid:105) = Tr G , (cid:104) b † b (cid:105) = (cid:104) ξ ξ (cid:105) = Tr G ,G ab = (Γ (2) − ) ab , Γ (2) ab = δ Γ δξ a δξ b . (26)where we suppress spatial or momentum indices. Tr runsover these as well as over the internal (field space) indices. More explicit formulae will be discussed in the next sec-tion. Furthermore, the three-point correlation is relatedto the one-particle irreducible three-point vertex via Eq.(25) (cf. e.g. [17]) (cid:104) Xb (cid:105) = −(cid:104) b † b b (cid:105) − (cid:104) b † b b (cid:105) (27)= (cid:88) a,b,c Tr[ G a G b + G a G b ] G c Γ (3) abc , Γ (3) abc = δ Γ δξ a δξ b δξ c . At this point, we stress that the parameter | s | in theequation of state (25) must not be interpreted as the con-densate fraction, though the formal appearance naivelysuggests such an interpretation. Instead, 2 | s | should beseen as the classical or mean field contribution to the to-tal particle density, and the rest of the equation is due tofluctuations on top of this mean field state. A standardinterpretation of the above equation is only possible inthe low density limits n → ,
2. Omitting the three-pointcorrelations, Eq. (25) reduces to leading order to the fa-miliar form from thermodynamics in the continuum for θ →
0, while taking a similar structure for θ → π , θ ≈ n = 2( δθ/ + (cid:104) b † b (cid:105) + 2 (cid:104) b † b (cid:105) , (28) θ ≈ π : n = 2 − [2( δθ/ + (cid:104) b † b (cid:105) + 2 (cid:104) b † b (cid:105) ] . In these cases, δθ/ θ = π/ (cid:104) b i (cid:105) = sc , however with the valueof θ determined from the implicit condition Eq. (23).More generally, we emphasize that Eqs. (23,25) providethe two exact, but implicit conditions that determine thetwo parameters θ, µ . A further nonzero expectation valuefor the the single atom degree of freedom, described by χ (cid:54) = 0, adds a further such condition analogous to Eq.(23).Our effective action formalism is capable to describethe system at any finite temperature. Calculations inthis regime are beyond the scope of this paper, but let ussketch how the high temperature disordered phase is de-scribed within our theory. Increasing the temperature inthe system will populate the connected parts of Eq. (25)increasingly such that at the phase transition to the sym-metric phase without symmetry breaking f , f →
0. Inother words, the condensate angles vanish, θ, χ →
0. Wemay interpret this scenario as the complete population ofthe “vacuum amplitude” f , which is needed to fulfill theholonomic constraint but does not enter the equation ofstate. The effect of destruction of the order parameterscan also be seen from the condition m π = ∂ U /∂π = 0.A finite temperature will act to generate a positive ther-mal mass or gap contribution, such that at some temper-ature there exists no finite θ, χ and a gapless mode ceasesto exist. At this point, where Goldstones’s theorem canno longer be satisfied, the symmetry broken phases be-come unstable and the system enters the disordered hightemperature phase. III. ASF–DSF PHASE BORDER
In this section we embark the calculation of the phaseborder. We will study the phase border by approach-ing it from the dimer superfluid side where there is noatomic condensate, and calculate at which interactionstrength the atoms become unstable towards an atomicsuperfluid. Thus, we first provide the explicit form ofthe Hamiltonian in the presence of a dimer superfluid,but for atoms in the normal phase. We then consider thelow density limits n → ,
2. In these limits, we can estab-lish a controlled small density expansion describing thedeviation from exactly n = 0 ,
2. The central objects forthe discussion are the atomic and dimer (di-hole) Greenfunctions, which we know exactly in the limits n = 0 , n = 0 ,
2) Green functions due to the condensatemean field. The dominant fluctuations in these limits arethus vacuum fluctuations renormalizing the Green func-tions, while the many-body effects can be captured interms of a Bogoliubov or spin wave theory. More specif-ically, we find that vacuum fluctuations strongly modifythe relation µ ( U ) compared to the mean field relation µ ( U ) = − U/
2, while the role of many-body effects con-sists mainly in depletion effects in the equation of state. We find that the high energy vacuum fluctuations have amuch more pronounced quantitative effect on the phaseborder than the condensate depletion in the limit n → n → n = 1, which takes place ateven stronger coupling, and thus more deeply bound two-particle states. We therefore propose an extrapolation ofthe scheme from the controlled limits n ≈ , n ≈ A. Rotated Hamiltonian for the dimer superfluidphase
Let us now focus on the phase border to the dimer su-perfluid state. As anticipated above, we address it fromthe DSF side where | (cid:105) is not macroscopically populatedand thus χ = 0. The kinetic and potential energy oper-ators read, in the new field coordinates, K (10) i = b † ,i ( cb ,i − s ∗ b ,i ) , K (21) i = ( s ∗ b † ,i + cb † ,i ) b ,i ,P (11) i = b † ,i b ,i , P (22) i = ( s ∗ b † ,i + cb † ,i )( sb ,i + cb ,i ) , (29)with c ≡ cos θ/ , s ≡ sin θ/ φ as above. We can nowwrite the Hamiltonian operator in terms of the new vari-ables, and implement the constraint via b ,i → X i =(1 − ˆ n ,i − ˆ n ,i ), absorbing the phase of b into the re-maining two degrees of freedom as discussed in Sec. II A.Further making use of the projective property X i = X i we find H [ b , b ] = H (10)kin + H (21)kin + H (split)kin + H pot , (30) H (10)kin + H (21)kin = − J (cid:88) (cid:104) i,j (cid:105) (cid:104) ( c + 2 | s | ) b † ,i X i X j b ,j + (2 | s | + c ) b † ,i b ,i b † ,j b ,j − c (cid:0) sb † ,i X i b † ,j b ,j + s ∗ b † ,i b ,i X j b ,j (cid:1)(cid:105) ,H (split)kin = −√ J (cid:88) (cid:104) i,j (cid:105) (cid:104) ( c − | s | ) b † ,i b ,i X j b ,j + c (cid:0) s ∗ X i b ,i X j b ,j − sb † ,i b ,i b † ,j b ,j (cid:1) + h . c . (cid:105) ,H pot = ( − µ + U ) M | s | + ( − µ + U )( c − | s | ) (cid:88) i b † ,i b ,i + [ − µ − ( − µ + U ) | s | )] (cid:88) i b † ,i b ,i +( − µ + U ) c (cid:88) i (cid:104) sb † ,i X i + s ∗ X i b ,i (cid:105) . We remind the reader that, as shown in Ref. [11] andbriefly discussed in Sect. II A, these operators b , maybe interpreted as standard bosonic operators. The cubicterm in the second line of H pot can be omitted from theoutset, and we will do this in the following: As arguedabove, the coefficient of the linear part has to vanish dueto the equation of motion, i.e. the equilibrium condition, and the cubic parts are connected to the linear ones viaEq. (24), such that their coefficient has to vanish as well.However, this does not exclude the possibility of nonlocal cubic terms as they appear in the kinetic terms of theHamiltonian. The total Euclidean action in the presenceof condensation reads S [ b , b ] = (cid:90) dτ (cid:16) (cid:88) i b † ,i ∂ τ b ,i + b † ,i ∂ τ b ,i + H [ b , b ] (cid:17) , (31)where the Hamiltonian is to be interpreted in the Heisen-berg picture and as a function of classical field variables.Quantizing this theory with the functional integral leadsprecisely to the representation of the effective action Eq.(7). B. Low density limits n ≈ , In the following, we will analyze the theory in the vicin-ity of the physical vacua where n ≈ ,
2, described by θ ≈ , π . The limits n = 0 , H n = − (cid:88) (cid:104) i,j (cid:105) (cid:2) g n, b † ,i X i X j b ,j + g n, b † ,j b ,j b † ,i b ,i + √ J ( b † ,i b ,i X j t ,j + b † ,j X j b † ,i b ,i ) (cid:105) + (cid:88) i [( U − µ n )ˆ n ,i − µ n ˆ n ,i ] . (32)The operators b , represent the bosonic single and two-particle excitations, corresponding to atoms and dimersresp. holes and di-holes. Here, for n = 0 we have g , = J, g , = 2 J, µ = µ and for n = 2, g , = 2 J, g , = J, µ = − µ + U . At these points the exact solution ofthe (two-body) scattering problem, and thereby an exactcalculation of the atomic and dimer Green’s function, isavailable as shown in [11]. While the case of the atomicGreen function is trivial as there are no renormalizationeffects in the vacuum, for the dimers/diholes we find theresults G − d ( ω ; µ , k ) = U + (cid:104) (cid:90) d d q (2 π ) d − (cid:15) q + (cid:15) q − k ) + i ω − µ (cid:105) − ,G − h ( ω ; µ , k ) = U + 34 (i ω − µ ) + 14 (cid:104) (cid:90) d d q (2 π ) d − (cid:15) q + (cid:15) q − k ) + i ω − µ (cid:105) − . (33)We will now perform a controlled expansion in the con-densate angle deviation δθ (cid:28) θ = 0 , π . It corresponds to a Bogoliubov approximationfor the condensation physics, but with coupling constantsobtained from the exact solution of the ”vacuum” scatter-ing problem. The procedure amounts to a resummationof ladder diagrams. These vacuum fluctuations are re-sponsible for strong shifts in the phase border as we willsee.Our expansion is defined with Hamiltonians of the form H = H n + δθ δH n + O ( δθ ) . (34)The additional Hamiltonian δH n generates new scatter-ing vertices which are O ( δθ ). Diagrams with more thanone of the new vertices may thus be discarded, and wemay restrict our attention to diagrams with at most oneof them. They will be discussed in a moment.In the low density cases n ≈ , θ = 0, we replace c = 1 , s = δθ/ δH = − J (cid:88) (cid:104) i,j (cid:105) √ (cid:2) X i b ,i X j b ,j − b † ,i b ,i b † ,j b ,j (cid:3) + X i b ,i b † ,j b ,j + h . c . (35) Similarly, around θ = π , setting s = 1 , c = − δθ/ δH = − δH . (36)Note, that the zero order Hamiltonians H n are relatedby a more complicated transformation of parameters, re-flecting the absence of a particle-hole symmetry.Let us now discuss the impact of the additional Hamil-tonian. The stability of the ASF phase is encoded in thefull atomic mass matrix, i.e. the inverse Green’s func-tion G − ( ω ; µ, q ) at zero frequency and momentum: Ifthe eigenvalues of the mass matrix are all positive, thephase without atomic condensate (i.e. the condenseddimer phase) is stable. The instability towards a statewith atomic condensate can thus been inferred from thevanishing of an eigenvalue of this matrix, ordet G − ( ω = 0; µ, q = 0) = 0 . (37)Thus we discuss the beyond mean field effects modifyingthe inverse atomic Green’s function. From the exact so-lutions of the vacuum problems at θ = 0 , π we know thatin these limits the inverse atom Green’s function is not di-rectly renormalized: there are no diagrams in the vacuumlimits which cause renormalization, but clearly, the func-tion µ ( U ) entering the atom propagator changes whentaking the exact dimer or di-hole Green’s function intoaccount. We now concentrate on the effects of δH n . Atlinear order in δθ , we find a direct condensate contribu-tion to the atom inverse propagator on the off-diagonal.This is the contribution familiar from Bogoliubov theoryin the low density limit, and we see that our generaliza-tion to arbitrary density produces such a structure alsoat high density. Now we have to consider the effect of thenew vertices. As argued above, we can restrict ourselvesto diagrams carrying a single one of them. We focus ondiagrams which renormalize the inverse atom propaga-tor. These are tadpole diagrams. The diagrams renor-malizing the diagonal entries must be O ( δθ ) in order to ensure particle number conservation. The diagramsrenormalizing the off-diagonal entries must involve oneof the new vertices O ( δθ ), and the trace over the innerline scales with a function f ( δθ ) with f (0) = 0. Con-sequently, the fluctuation contributions are higher thanlinear order for both diagonal and off-diagonal entriesand can be discarded. Thus, the full atomic mass matrixat O ( δθ ) reads (we separate true potential (binding) en-ergy from kinetic energy, µ ( U ) = E b ( U ) / − Jz at n ≈ µ ( U ) = − E b ( U ) / U + 2 Jz at n ≈ z = 2 d thelattice coordination number) G − ,θ =0 ( ω = 0; µ, q = 0) = (cid:18) − E b ( U ) / √ Jzδθ/ √ Jzδθ/ − E b ( U ) / (cid:19) , (38) G − ,θ = π ( ω = 0; µ, q = 0) = (cid:18) − E b ( U ) / − √ Jzδθ/ − √ Jzδθ/ − E b ( U ) / (cid:19) . Hence we conclude that the dominant effect beyond meanfield theory, which implies the simple linear relation µ ( U ) = U/ µ ( U ) (or, equivalently, E b ( U )) as a function of the in-teraction strength U . These high energy fluctuations areresponsible for shifts of the critical point. Note that thesign on the off-diagonal of the second equation is unphys-ical as it can be absorbed by a phase rotation of the orderparameter. Thus, the equations have the same form.The critical interaction strength can now be extractedfrom the characteristic equation (37) which reads explic-itly, in both limits,( − E b ( U c ) / − Jzδθ ) = 0 . (39)The physical solution is given by E b ( U c ) = − √ Jzδθ ,i.e. in the vacuum limit δθ → d = 3, the binding energy startsquadratically due to the non-analyticity in the fluctua-tion integral. In contrast, in d = 2, the fluctuation inte-gral features the well known logarithmic behavior. Thisyields, in dimensionless units, the physical solutions n → d = 3 : ˜ U c = ˜ U (cid:0) / σ √ δθ (cid:1) , (40) d = 2 : ˜ U c = 2 π (cid:0) log √ δθ Λ (cid:1) − n → d = 3 : ˜ U c = ˜ U (cid:0) − / σ | ˜ U |− | ˜ U | √ δθ (cid:1) ,d = 2 : ˜ U c = − π (cid:0) log √ δθ Λ (cid:1) − FIG. 2. ASF-DSF phase boundary for the attractive 3-bodyhardcore constrained Bose-Hubbard model. The black lowercurve is the mean field result, while red upper and blue mid-dle curve correspond to d = 2 ,
3. The crosses denote theendpoints of the critical lines. For high densities, the d = 2 , with numbers ˜ U ≈ − / , σ ≈ . , Λ ≈ .
50, and ˜ U ≈− / n → non-analytic dependence of thecritical interaction strength on the condensate density δθ which is due to the strong fluctuation effects in thevicinity of the bound state formation. This is in contrastto the mean field result, which shows a linear dependenceon the condensate angle δθ = √ n (see [3] and Eq. (45)below).0 C. Phase Diagram
The analysis of the low density limits in the last para-graph reveals that the strongest beyond mean field effectcomes from the high energy vacuum fluctuations renor-malizing the relation µ ( U ) away from its mean field value µ = U/
2. In our practical implementation of the calcula-tion of the phase boundary, we will rely on this separationof vacuum and many-body effects. For the calculation ofthe phase diagram at a fixed total density n , we now alsotake modifications of the mean field equation of state n = 2 | s | into account. We discuss the equation of state (25) in an approximation where we omit the three-pointcorrelations from the outset. Furthermore, in our con-crete computations, we restrict ourselves to the calcula-tion of the atomic depletion (cid:104) b † b (cid:105) . We will find that thiscontribution is small compared to the condensate partand thus has a small influence on the phase border only– the dominant effect shifting the border stems from thevacuum fluctuations leading to a strong modification of µ ( U ). Based on this observation, we do not expect thatthe dimer depletion strongly modifies the phase border.In order to calculate the atomic depletion, we first con-sider the quadratic spin wave theory, S [ b † , b ] = − J (cid:88) (cid:104) i,j (cid:105) (cid:104) ( c + 2 | s | ) b † ,i b ,j + √ c (cid:0) s ∗ b ,i b ,j + sb † ,i b † ,j (cid:1)(cid:105) + [ − µ − ( − µ + U ) | s | )] (cid:88) i b † ,i b ,i (41)= 12 (cid:90) q ( b † ( q ) , b ( − q )) (cid:18) i ω + p q ∆ q ∆ q − i ω + p q (cid:19) (cid:18) b ( q ) b † ( − q ) (cid:19) ,p q = p + δp q , p = − µ (1 − s ) − U s − Jz (1 + s ) , δp q = 2(1 + s ) δ(cid:15) q , ∆ q = ∆ + δ ∆ q , ∆ = 2 √ Jzcs, δ ∆ q = 4 √ csδ(cid:15) q ,δ(cid:15) q = J (cid:88) λ (1 − cos q λ ) . The matrix in the second line is the inverse atomic Greenfunction in frequency and momentum space G − ( ω ; µ, q ).With these preparations, the approximate equation ofstate and the atomic depletion is found to be n = 2 | s | + ( c − | s | ) (cid:104) b † b (cid:105) , (42) (cid:104) b † b (cid:105) = Tr G ( ω ; µ, q ) = 12 (cid:90) d d q (2 π ) d (cid:16) p q (cid:113) p q − ∆ q − (cid:17) , where we have performed the frequency integral by clos-ing the contour in the upper half plane.As stated above, we consider the renormalization ef-fects on the inverse atom propagator which are presentalready in the vacuum problem. These are encoded inthe value of the chemical potential µ ( U ) and depend ondimension. The condition determining µ reads G − d/h ( ω = 0; µ, q = 0) = 0 , (43)where in the vicinity of n = 0 we use the full dimerGreen’s function G d , and close to n = 2 the di-hole ex-pression G h as given in Eq. (33). In contrast, in themean field approximation which uses the “bare” inversedimer propagator i ω − µ + U , the above condition eval-uates to µ ( U ) = −| U | / G − ( ω = 0; µ, q = 0) = p − ∆ = 0 . (44) We solve the system of equations (42,43,44) in two andthree dimensions numerically. In particular, Eq. (43) de-couples from the other two equations within our approx-imation scheme yielding the renormalized relation µ ( U ),such that we merely need to solve (42,44) with µ ( U ) asan input. The result for the phase diagram is plotted inFig. 1. We compare these results to the mean field ap-proximation, which uses µ = U/ n = 2 | s | for theequation of state, resulting in the critical interaction U c Jz = − (cid:112) − n/ √ n ) . (45)As anticipated above, the beyond mean field effects aremainly due to the fluctuations accompanying the boundstate formation, which strongly modify the relation µ ( U )as compared to the mean field value, while we find a sub-dominant role of the many-body depletion effects.Theshape of the phase boundary directly reflects the non-analytic, dimension dependent behavior associated to thebound state formation. The quantitative effect is morepronounced below half filling than above. This may betraced back to the fact that the domain where nonpertur-bative fluctuation effects play a role is smaller in the highdensity regime than at small densities, cf. [11]. A simplepicture can be given as follows: In the limits n → , G − ) and the microscopic bound state formation (zeroeigenvalue of G − d/h ) coincide and fluctuation effects onthe phase boundary are substantial. Moving away fromthese limits, the absolute value of the critical interac-tion increases, and the microscopic bound state is al-1 FIG. 3. Condensate depletion due to the atomic two-pointfunction in two dimensions. The overall shape of the curvewith zero crossing and sign change is determined by the func-tion c − | s | premultiplying the two-point function, cf. Eq.(25). The small overall size shows that the depletion effectsproduce only very tiny corrections to the phase border. ready tightly bound at the point where atom criticalityis reached. The critical line then approaches the meanfield phase boundary up to small perturbative correctionsThus, though our approximation is lacking a strict order-ing principle when moving away from the limits n → , IV. MANY-BODY PHYSICS IN THE STRONGCOUPLING REGIME AND A CONTINUOUSSUPERSOLID
In this section we investigate the system in the stronglycorrelated limit J/ | U | →
0. In particular, we identify abicritical point at half filling of atoms ( n = 1), at whichhomogeneous superfluid order (spontaneous phase sym-metry breaking) and charge density wave order (spon-taneous translation symmetry breaking) are degenerate.The bicritical point is due to a symmetry enhancementfrom the conventional U (1) (cid:39) SO (2) to SO (3), whichis seen to be intimately connected to the 3-body con-straint. We term the system in this regime a “contin-uous supersolid” due to the degeneracy of phase andtranslation symmetry breaking orders, where the orderparameter may be rotated continuously from one to theother without energetic cost. This behavior is in contrastto other occurrences of supersolidity in bosonic systems[19]. Though this state is only reached asymptotically, itgoverns the physics in strong coupling and close to halffilling, and we work out the observable consequences ofthis situation. We also propose a simple experiment toverify this scenario. A. Analytical Approach
Before embarking the calculation, let us stress thatour beyond mean field approach is indispensable to set-tle these issues. Indeed, a straightforward compari-son of the simple Gutzwiller mean field energies of thedimer superfluid and a charge density wave state (CDW)yields degenerate energies for these two states for all fill-ings: The Gutzwiller mean field CDW state is given by | CDW (cid:105) = (cid:81) i even | (cid:105) i | (cid:105) i +1 , and for a fixed average den-sity n has energy density E CDW /M d = ( U/ n , whichprecisely equals the mean field energy density E DSF /M d of the dimer superfluid for all particle densities. In con-sequence, the question of the correct ground state cannotbe decided within the simple Gutzwiller mean field theory(though a superfluid is clearly more plausible for incom-mensurate fillings). It is necessary to first integrate outthe high energy single particle degrees of freedom, mak-ing the dimers true propagating and interacting physicalexcitations. Moreover, even the second order perturba-tion theory is not fully conclusive as we will see. Thedeviation from the second order result, calculated in [11],plays a key role in the following discussion.In [11] we have calculated the effective theory in thestrong coupling limit: Perturbatively integrating out thesingle particle excitations up to fourth order, and tak-ing the constraint principle for the effective action intoaccount, we found the low energy effective Hamiltonian H eff = − t (cid:88) (cid:104) i,j (cid:105) (cid:0) t † ,i X i X j t ,j − λ ˆ n ,i ˆ n ,j (cid:1) − µ d (cid:88) i ˆ n ,i (46)with t = 2 J / | U | , µ d the effective dimer chemical po-tential and the dimensionless ratio of nearest-neighbourinteraction to hopping λ = v/ (2 t ) discussed below. Sincein the perturbative limit there are only virtual single par-ticle excitations, we may replace the constraint operator X i = 1 − ˆ n ,i − ˆ n ,i → − ˆ n ,i . In this case we have thefollowing mapping to effective spin 1 / s + j = ( s xj + is yj ) ≡ ( − ) j t † ,j X j , (47) s − j = ( s xj − is yj ) ≡ ( − ) j X j t ,j ,s zj = ˆ n ,j − / , where on the bipartite lattice with sublattices A and B we use ( − ) j = + for j ∈ A and ( − ) j = − for j ∈ B . Upto a constant the Hamiltonian then takes the form H eff = 2 t (cid:88) (cid:104) i,j (cid:105) (cid:0) s xi s xj + s yi s yj + λs zi s zj (cid:1) − µ d (cid:88) i s zi . (48)The anisotropy parameter λ = v/ (2 t ) evaluates to λ = 1in the second order perturbation theory, corresponding toan isotropic antiferromagnetic Heisenberg model – notethe sign change in t due to the sublattice dependent sign2in (47). The fourth order calculation yields λ = 1 − z − (cid:0) JU (cid:1) < . (49)This result has been derived in [11], cf. Sect. V.D.2,Eq. (60). It is obtained as the ratio of dimer-dimerinteraction and dimer hopping coefficient calculated atfourth order. Next-to-nearest neighbour terms also ap-pear at fourth order, but their numerical coefficient ismuch smaller than for the nearest neighbours due to therestricted pathways contributing to these terms, and arethus neglected.For λ = 1 and half filling of atoms n = 1, where theterm involving the chemical potential vanishes, the sys-tem exhibits a symmetry enhancement from SO (2) (cid:39) U (1) (corresponding to rotations in the x − y − plane, orphase rotations, generated by exp i θ z S z ∝ exp i θ z ˆ N , withglobal operators S α = (cid:80) i s αi , ˆ N = (cid:80) i ˆ n ,i ) to SO (3)(corresponding to arbitrary rotations on the Bloch sphereexp i (cid:80) α θ α S α ).The SO (3)-invariance of the quadratic part in (48)also implies a simple transformation behavior of the to-tal Hamiltonian under a discrete particle-hole or chargeconjugation transformation: The special choice U =exp(i πS x ) implements the mapping s ± j → U − s ± j U = s ∓ j , s zj → U − s zj U = − s zj . (50)Under such a transformation N → M d − N and, hence, S z → M d / − S z . The particle-hole symmetry makesthe phase diagram of deeply bound dimers symmetricunder the replacement ˆ n → − ˆ n . In general, such asymmetry is absent. Moreover, it is also not present inthe opposite limit of strong repulsive interactions, whichis asymmetric when ˆ n is replaced with 2 − ˆ n .The SO (3)-invariance is a peculiar feature of the lead-ing order perturbation theory. Its physical origin is wellunderstood in terms of the geometric argument whichrelates hopping and interaction paths, cf. Sec. IV D inthe companion paper [11]. At second order, no other in-teraction processes can occur, and thus the hopping andinteraction constants must be equal. However, as seenin [11], at fourth order perturbation theory additionalinteraction processes yield λ <
1, thus reducing the sym-metry to SO (2), and also spoiling the particle-hole sym-metry. In addition, several other terms are generated,which describe next-to-nearest neighbour hopping and One may wonder about implicit density effects for the perturba-tive calculation, that is, if the perturbative calculation at zerodensity is sufficient for the calculation of the effective theoryfor all densities. This may be discussed by studying the limits n = 0 ,
2. At second order, the emergent particle-hole symme-try ensures that the perturbative results coincide for both cases,while at fourth order, differences occur, pointing at the abovementioned implicit density effects. However, due to the identicaldiagrammatic structure one still finds λ < n = 2, which isthe crucial ingredient for the argument presented here. interaction, or three- and four-spin interactions. Never-theless, the proximity to the Heisenberg point λ = 1 hasan impact on the phase diagram, and we will use its wellknow properties to understand the phase diagram andthe nature of the low lying excitations in the perturba-tive regime.The order parameter for this model is given by the ex-pectation value of the N´eel vector ˆ N α = (cid:80) j ( − ) j s αi . Itsvector character is under SO (3), [ S α , ˆ N β ] = i (cid:15) αβγ ˆ N γ ,i.e. global spin rotations transform the N´eel vector com-ponents into each other. Translating back to the origi-nal boson language, a finite (cid:104) ˆ N z (cid:105) corresponds to chargedensity wave order, while finite values of (cid:104) ˆ N x (cid:105) , (cid:104) ˆ N y (cid:105) in-dicate dimer superfluid (DSF) order. For the isotropicHeisenberg model without magnetic field (or at half fill-ing), [ H, S α ] = 0 for all α , and thus CDW and DSF or-der are degenerate. The perturbative limit of our modelthus realizes a bicritical point [20] with two competing or-ders. Such an enhancement of internal symmetries is wellknown in magnetic systems [12, 13], but less common andintuitive in physically realizable bosonic models, whichusually only exhibit phase symmetry. Due to the degen-eracy of phase and translation symmetry broken states,both order parameters are generically nonzero, and wemay term the state a continuous supersolid, whose ex-perimental implications are studied below.We can make this discussion even more explicit whenchanging from the spin to a hardcore boson language s − j = ( − j h j , s + j = ( − j h † j , s zj = h † j h j − , (51)where the hardcore bosons obey h † j ≡
0. We con-sider infinitesimal SO (3) transformation with the pa-rameters ε α (cid:28) α = ( x, y, z )), δs αj = i (cid:2) S, s αj (cid:3) , where S = (cid:80) iβ ε β s βi . The explicit form of the above transfor-mation reads ( s ± j = s xj ± is yj ) δs + j = − iεs zj − iε z s + j , (52) δs − j = iε ∗ s zj + iε z s − j ,δs zj = − i ε ∗ s + j + i εs − j , with ε = ε x + iε y . In terms of bosonic hardcore operatorsone thus obtains δh † j = ( − j (cid:20) − iε ( h † j h j −
12 ) (cid:21) − iε z h † j ,δh j = ( − j (cid:20) iε ∗ ( h † j h j −
12 ) (cid:21) + iε z h j ,δ ( h † j h j −
12 ) = ( − j (cid:104) − i ε ∗ h † j + i εh j (cid:105) . (53)Note that the last terms in the first and second lines arejust usual gauge transformations.The change in the operators results in the change oftheir mean-field values. If one introduces the usual super-fluid order parameter ψ = (cid:80) j (cid:104) h j (cid:105) and the CDW orderparameter N ≡ (cid:104) ˆ N z (cid:105) = (cid:80) j ( − j (cid:104) ( h † j h j − ) (cid:105) , then the3corresponding change in the order parameters is δψ = iε z ψ + iε ∗ N , (54) δ N = − i ε ∗ ψ ∗ + i εψ. It is easy to check that the above transformation leavesthe combination N + | ψ | invariant. We therefore canconclude that the SO (3) symmetry corresponds to canon-ical transformations of the dimer operators, which changeboth superfluid and CDW order parameters, but leavethe combination N + | ψ | invariant.Another important point about the SO (3) symmetryis that it is broken for n (cid:54) = 1 not on the Hamiltonian level(the Hamiltonian with λ = 1 is always SO (3) symmetric)but on the level of the subset of states (with a fixed S z ),on which it has to be minimized. In a generic case n (cid:54) = 1(and, hence, S z (cid:54) = 0) the subspace reduces the symme-try down to SO (2) (cid:39) U (1) gauge group. In the case n = 1 with S z = 0, however, the subspace contains themanifold of spin-singlet (and, therefore, SO (3) symmet-ric) states, which have the lowest energy. The symmetrytransformation corresponds simply to the motion on thismanifold.A similar scenario (an enhancement to a pseudo SU (2)symmetry) is actually observed in attractive latticefermion systems [14]. Similar to the fermion system, thesymmetry enhancement is thus a unique consequence ofthe 3-body hardcore constraint. Indeed, attractive lat-tice bosons without such constraint, analyzed in detailby Petrosyan, Schmidt et al. [21], show a different be-havior: Due to the possibility of virtually occupying alattice site with three atoms, it is found λ = 4. Thisplaces the unconstrained attractive bosons in the “Isinglimit”, which was analyzed further in the latter reference.As we find λ < SO (2) → SO (3) implies the emergence of a secondgapless, and therefore collective, Goldstone mode. For aweakly explicitly broken SO (3) symmetry, one still has anear gapless collective mode with experimentally observ- able consequences. We propose an experiment, which isbased on the idea of explicitly rotating the macroscopicN´eel vector from the x − y − plane representing superfluid-ity to the z − axis, realizing a CDW ordered state. We willalso show that this experiment allows to quantitativelycharacterize the pseudo Goldstone mode.To favor CDW ordering, we explicitly break the latticetranslation symmetry via introduction of a superlatticeshifting the single particle energies on adjacent sites: − µ d (cid:88) i ˆ n ,i → − µ A (cid:88) i ∈ A ˆ n ,i − µ B (cid:88) j ∈ B ˆ n ,j (55)= − µ (cid:88) i ˆ n ,i + ¯ ν (cid:16) (cid:88) i ∈ A ˆ n ,i − (cid:88) j ∈ B ˆ n ,j (cid:17) ,µ = µ A + µ B , ¯ ν = µ A − µ B . Here, µ is and average chemical potenial and ¯ ν an imbal-ance parameter. Now we calculate the mean field groundstate as well as the spectrum of excitations of the effectivelow energy Hamiltonian (46), using the rotation formal-ism (cf. Sec. II A). The approach is fully equivalent tothe leading order 1 /S expansion, which is not a well con-trolled expansion for S = 1 /
2, but is known to yield themain features of the Heisenberg model in external fields.At the Heisenberg point, the SO (3) spin symmetryrequires an enlarged parameter space for the order pa-rameter describing the ground state of the system. Weconsider an ansatz parameterized by two angles for therotation of the t degree of freedom (in contrast to the U (1) (cid:39) SO (2) case with a single rotation angle for t ), R ( θ, ϕ l ) = (cid:18) cos( θ + ϕ l ) / θ + ϕ l ) / φ − sin( θ + ϕ l ) / − φ cos( θ + ϕ l ) / (cid:19) . (56)Here, l = A, B is an index which depends on the sublat-tice A or B , thus enabeling the description of a spatiallymodulated phase. For example, the choice θ = 0 , ϕ A =0 , ϕ B = π describes a charge density wave. The homoge-neous choice θ (cid:54) = 0 , ϕ A = ϕ B = 0 describes a superfluidground state. The rotation matrix is only 2 × b ,i → X i = 1 − ˆ n ,i , to second order we obtain4 E mf /M d = − tz ( s A c A s B c B − λs A s B ) − µ s A + s B ) + ¯ ν s A − s B ) , (57) H lin = (cid:2) tz (cid:0) − c B s B ( c A − s A ) + 2 λs B s A c A (cid:1) + ( − µ + ¯ ν ) c A s A (cid:3) (cid:88) i ∈ A ( b † ,i + b ,i )+ (cid:2) tz (cid:0) − c A s A ( c B − s B ) + 2 λs A s B c B (cid:1) + ( − µ − ¯ ν ) c B s B (cid:3) (cid:88) j ∈ B ( b † ,j + b ,j ) ,H SW = t (cid:88) (cid:104) i,j (cid:105) (cid:2)(cid:0) − ( c A c B + s A s B ) − λc A s A c B s B (cid:1) ( b † ,j b ,i + b † ,i b ,j )+ (cid:0) c A s B + c B s A − λc A s A c B s B (cid:1) ( b ,j b ,i + b † ,i b † ,j ) (cid:3) + (cid:0) tz (2 c A s A c B s B − λs B ( c A − s A )) + ( − µ + ¯ ν )( c A − s A ) (cid:1) (cid:88) i ∈ A b † ,i b ,i + (cid:0) tz (2 c A s A c B s B − λs A ( c B − s B )) + ( − µ − ¯ ν )( c B − s B ) (cid:1) (cid:88) j ∈ B b † ,j b ,j . Here s l = sin( θ + ϕ l ) etc., and we have set the sponta-neously chosen phase φ = 0 without loss of generality. M is the total number of sites in each lattice direction.Since we break the lattice translation symmetry via ourchoice of the ansatz for the ground state, it is importantto be careful with the position indices – b ,i are locatedon the sublattice A , b ,j on the sublattice B . Eventuallywe are interested in a situation at fixed density. The localdensity operator to quadratic order takes the form t † ,i t ,i = s l + c l s l ( b † ,i + b ,i ) + ( c l − s l )ˆ n ,i , (58)where l = A ( B ) for i ∈ A ( B ). Thus, in the mean fieldapproximation the equation of state reads n = N/M d = ( s A + s B ) . (59)At half filling n = 1, this implies s A = c B , s B = c A .Together with the relations s l + c l = 1, within this ap-proximation everything may be expressed in terms of e.g. s A alone. In particular, the mean field energy determin-ing the ground state takes the form E mf /M d = − tzs A (1 − s A )(1 − λ ) + ¯ ν (2 s A − . (60)The ground state is found from identifying the stableminima with respect to variation in θ , and thus we have dropped the contribution from the chemical potential,as it contributes a rotation angle independent constantonly for effectively fixed density. For ¯ ν = 0, the groundstate for λ < s A =1 /
2. For λ >
1, the charge density wave with s A =0 , s B = 1 is favored. At the Heisenberg point λ = 1, bothstates are degenerate in accord with the exact symmetryargument. Now we consider the relevant case λ < ν away from zero by ramping the superlattice,the superfluid acquires a spatial modulation, s A = (1 + ν ) / , s B = (1 − ν ) / , ν = ¯ ν/tz ( λ − | ν c | = 1 , | ¯ ν c | = tz (1 − λ ) ≈ tz ( z − J/U ) (61)the SF is destroyed in favor of the CDW. Thus, rampingthe superlattice corresponds to rotating the N´eel orderparameter. As we will see below, the critical value cor-responds precisely to the gap of the pseudo-Goldstonemode. Hence, via measurement of the SF correlations[22], which cease to exist at ¯ ν c , one can quantitativelydetermine the characteristic property of the additionalcollective mode.The chemical potential µ is determined from the equi-librium condition that the linear terms vanish, evaluat-ing to µ = λtz independent of ¯ ν . Inserting this and theabove expression for s A , and switching to the Lagrangeformalism, we obtain the Gaussian action S = 12 (cid:90) q,q (cid:48) δ ( q − q (cid:48) )( b ( − q ) , b † ( q )) (cid:18) g q h q ( − ω ) h q ( ω ) g q (cid:19) (cid:18) b ( q (cid:48) ) b † ( − q (cid:48) ) (cid:19) , (62) g q = tz ( λ − (1 − ν ) + 1)˜ (cid:15) q , h q ( ω ) = i ω + tz (1 + λ − (1 − ν )˜ (cid:15) q ) , ˜ (cid:15) q = d (cid:88) λ cos q λ . The spectrum of excitations may be computed from thecondition that the determinant of the fluctuation matrix vanish. We obtain ω ( q ) = ± tz (cid:0) ((1 + [ λ − − ν ])˜ (cid:15) q + 1)(1 − ˜ (cid:15) q ) (cid:1) / . (63)5For ¯ ν = 0 and λ <
1, the dispersion simplifies to ω ( q ) = ± tz (cid:0) ( λ ˜ (cid:15) q + 1)(1 − ˜ (cid:15) q ) (cid:1) / , and there is a singleGoldstone mode at q = 0, corresponding to the spon-taneously broken U (1) (cid:39) SO (2) symmetry in the dimersuperfluid. For λ (cid:46)
1, there is a second near gapless modewith gap tz (1 − λ ) located at the edges of the Brillouinzone q = π . At the Heisenberg point λ = 1, this gapcloses, and the system features the two Goldstone modescorresponding to the spontaneously broken SO (3) sym-metry. The system is then at the bicritical point wherethe order parameter can be freely rotated on the Blochsphere. In the general case, the gap of the second neargapless mode is given by∆ = tz (1 − λ )(1 − ν ) , (64)and we observe that we reach a point where there aretwo exactly gapless modes by tuning ¯ ν → ¯ ν c , ν →
1. Inthis case, the two gapless modes correspond to the char-acteristic excitations on an antiferromagnetic, or CDW,ground state, and there is no superfluid order as theBloch vector is confined to the z − axis. In Tab. I, wesummarize the dispersions found in the different densityregimes in the leading order perturbation theory limit λ = 1, and ¯ ν = 0.In summary, we propose a conceptually simple exper-iment that allows to rotate the macroscopic N´eel vectororder parameter via ramping a superlattice. The mea-surement of superfluid and density-density correlations[22] allows to monitor this rotation, as well as to mea-sure the gap of the collective pseudo-Goldstone mode,which is the hallmark of the proximity of the system tothe bicritical point with enhanced symmetry. Alterna-tive experiments for the investigation of this proximityinclude a direct measurement of the dispersion relationvia Bragg spectroscopy on the lattice [23], or analyzingthe system subject to slow rotation, which also acts as acurrent defavoring SF against CDW order [24].We further comment on the relation of our spatiallymodulated superfluid for nonvanishing ¯ ν to a supersolid.The latter is defined as a state with simultaneously andspontaneously broken phase and translation symmetry.In our case, both symmetries are broken, but the transla-tion symmetry breaking is explicit and not spontaneous.Though the correlations are those of a supersolid, wewould not term the state as such.Finally, we note that the evolution of the system fromrepulsive to attractive coupling may be viewed as a tran-sition from a spin 1 model (3 onsite states) to a spin 1/2model. The x-y ordered phases of these two models areseparated by the Ising transition discussed in more detailin the next section. B. Complementary exact numerical study in onedimension
We now investigate how the key features of these re-sults manifest themselves in a 1D system. This can be
TABLE I. Low energy dispersions in the perturbative regime(second order). z is the dynamic exponent, ω ∼ | q | z . n = 0 0 < n < n = 1 1 < n < n = 2 z − − − K CDW K DSF − − − x N=60N=50 CDW(x)DSF(x)
FIG. 4. Computations of the ground state on a 1D latticewith 60 sites for
U/J = −
20 using TEBD methods. (a) Cor-relation functions characterizing the CDW and DSF phasesin a system with open boundary conditions ploted on a log-log scale as a function of distance x . Values are shown for N = 60 particles, where the correlation functions are almostidentical, and N = 50 particles, where the decay of the CDWcorrelations are significantly more rapid than the DSF cor-relations (b) Algebraic decay exponents K DSF and K CDW that are fitted to the envelope of the correlation functions forvarying mean density n . Error bars show typical errors in thefitted decay. done by computing the ground state of the constrainedBose-Hubbard model using the Time Evolving BlockDecimation (TEBD) algorithm [25]. Note that we op-timise our algorithm for the conserved total number ofparticles [26], analogously to the optimisation for goodquantum numbers in Density Matrix RenormalisationGroup methods [27]. In Ref. [3] we already observedquasi off-diagonal long range order in the Dimer Super-fluid (DSF) correlation function (cid:104) b † i b † i b i + x b i + x (cid:105) , togetherwith exponential decay of off-diagonal elements in thesingle-particle density matrix (cid:104) b † i b j (cid:105) . This indicated thetransition between the ASF and DSF phases in the 1Dsystem.Here we particularly investigate the situation near half-filling n = 1, paying attention to the interplay betweenDSF order and CDW order, characterized by the density-density correlation function CDW ( x ) = (cid:104) n i n i + x (cid:105) −(cid:104) n i (cid:105)(cid:104) n i + x (cid:105) . In Fig. 4 we compare the DSF and CDWcorrelation functions for the ground state on a 60 site lat-tice with U/J = −
20 and open boundary conditions. InFig. 4a we plot the correlation functions both for N = 60(half filling of dimers) and N = 50. At half filling the al-gebraic decay of these correlation functions is essentiallythe same, indeed the correlation functions are essentiallyequal, indicating coincidence of CDW and DSF orders inthis state. Whilst reducing the total number of particleson the lattice to N = 50 does not significantly changethe DSF, the density-density correlation function decays6 − − −
10 0 − − − − − U/J K i K CDW K DSF
FIG. 5. Numerical validation of symmetry enhancement SO (2) → SO (3). The plot shows the exponents K DSF and K CDW describing the algebraic decay of DSF and CDW or-der in the ground state as a function of the ratio
U/J at unitfilling n = 1. For sufficiently large interactions, the coinci-dence of the decay exponents signals the degeneracy of thetwo kinds of order. These calculations were performed withTEBD methods for 60 particles on 60 sites. Open bound-ary conditions were used, but the decay exponents were fittedwithin the central 30 sites. The fitting errors are similar tothose in Fig. 4 (b). Number of Schmidt coefficients retainedin TEBD calculations χ = 200. much more rapidly in the ground state, in addition tolarge superimposed oscillations. This relative sensitivityof the correlation functions is characterised in Fig. 4b,where we show the result of fitting an algebraic decay x K i to the envelope of each of the correlation functions.Again, we see that the decay of CDW and DSF correla-tions is identical within fitting errors at unit filling, butthe CDW is very sensitive to deviations from unit filling,and it is dominated away from n = 1 by the DSF.In Fig. 5, we study the approach of the bicritical pointat fixed half filling n = 1 as a function of the ratio ofhopping and interaction J/U . The plot clearly shows thesymmetry enhancement from the conventional U (1) (cid:39) SO (2) to an SO (3) symmetry: The two decay exponentsdescribing DSF and CDW order approach each other forsufficiently strong attractive onsite interaction U , thusindicating the degeneracy of the two different kinds oforder.In an experiment it would be difficult to produce asetup with an exactly commensurate number of parti-cles and lattice sites. One way to observe emergence ofthe CDW order, though would be to prepare the sys-tem in a harmonic trapping potential, where the densitywould vary across the trap. In Fig. 6 we investigate theground state for 30 particles on 60 lattice sites in thepresence of such an external harmonic trap. In Figs. 6a,bwe show the correlation functions for CDW and DSF or-der as they vary across the trap. We note that the DSForder is significant throughout the occupied region. We x n ( x ) −3 −2 −1 x D S F ( x ) V SL =0V SL =0.01J x xyy x FIG. 6. Computations of the ground state with 30 particleson 60 lattice sites with
U/J = −
20, in the presence of a har-monic trap with on-site potential V ( x ). (a) Shaded plot ofthe CDW correlation function (cid:104) n x n y (cid:105) − (cid:104) n x (cid:105)(cid:104) n y (cid:105) (with inter-polated colours, and diagonal elements not shown), showingsubstantial order near x = y = 20 and x = y = 40, where themean filling factor n ∼ V ( x ) = V tr ( x − . , V tr = 1 . × − . (b) Shaded plot of the DSF correla-tion function DSF ( x − y ) = (cid:104) b † x b † x b y b y (cid:105) (with interpolatedcolours), showing substantial order across the occupied re-gion of the lattice, for the same trap parameters as parta. (c) Density n ( x ) for the same parameters as parts a,b.(d) Plots of the DSF correlation function DSF ( x ) for a trap V ( x ) = V tr ( x − , V tr = 0 . / ≈ . × − , with andwithout an additional superlattice potential V SL / (cid:80) i ( − ) i . have checked in addition that across this region, the off-diagonal elements decay algebraically as a function ofdistance. On the other hand, the CDW correlations aremost significant in regions near unit filling. For the trapparameters chosen here, this occurs near sites 20 and 40,as shown in Fig. 6c, where we also see significant oscilla-tions in the density, which are also characteristic of theappearance of CDW order. In Fig. 6d we then investi-gate how the order can be manipulated by the additionof a weak superlattice. We see that the addition of an al-ternating potential on the order 0 . J is sufficient to sig-nificantly increase the algebraic decay exponent for DSForder. Because the system size is small, it was difficult toobtain reliable results for the algebraic decay exponent ofCDW order, but our calculations indicate that applyingsuch a superlattice can indeed be used to select the dom-inant order for a system in the presence of a harmonictrap.Using t-DMRG methods we can also investigated pos-sible time-dependent preparation of the continuous su-persolid beginning from a Mott insulating state in thepresence of a superlattice, analogously to the studies per-formed in Ref. [3]. Beginning in an insulating state withtwo atoms in the lowest wells of a period two superlat-tice, it is possible to prepare a state with n = 1 and U/J = −
20 in a timescale of the order of 100 J − , with7good fidelity of the DSF correlation functions, providedthat a sufficiently strong constraint can be imposed sothat no loss events occur on the timescale of the ramp. V. LONG WAVELENGTH LIMIT: NATURE OFTHE PHASE TRANSITION
Even at low energies, non-linearities in the effectiveaction may in principle have an impact on the physi-cal observables, such as the nature of the phase tran-sition. Such a scenario is known as Coleman-Weinbergphenomenon [15]: Two near gapless degrees of freedomare coupled to each other, in a way that a phase tran-sition which one of them undergoes is driven first orderdue to the long wavelength fluctuations of the other: thefirst order transition is radiatively induced.In our problem, indeed we face competing low energydegrees of freedom at the ASF-DSF transition: First,there is the gapless Goldstone mode present in the dimercondensate, which does not undergo qualitative changesat the ASF-DSF transition point. Second, at the Isingtype transition one expects a Z degree of freedom toemerge in the low energy sector for the atom degreesof freedom. A possible coupling between those low en-ergy degrees of freedom may or may not give rise to aColeman-Weinberg mechanism.Here we study this question by means of a systematicderivative expansion of the effective action. At “low”densities n ≈ , n = 1. Thus, we identify a true Ising quantumcritical point in our system, connecting the two orderedASF and DSF phases.Note that the scenario crucially hinges on the controlover a coupling of the two near gapless modes close tothe transition. It is evident that the discussion cannot be lead based on a simple quadratic spin wave theory. A. Low Energy Derivative Expansion
Our strategy is as follows: We will approach the phasetransition from the DSF side, where there is not yet anatomic condensate, and tune the atomic mass parame-ters to criticality from there. For this purpose, we drawthe low energy, continuum limit of the effective actioncorresponding to Eq. (30). We then identify the relevantlow energy fluctuations and integrate out the massive de-grees of freedom. We arrive at an action that describesthe dimer Goldstone physics, the Ising degree of freedomas well as a cubic coupling of Goldstone mode to Isingdensity. The derivation is similar to the one presentedin [7] in the continuum, differs however in the crucialaspect that the model discussed there features alreadymicroscopic propagating dimer degrees of freedom. Herewe show how such terms are generated via successive in-tegration of the massive degrees of freedom.At low energies, the action corresponding to the Hamil-tonian (30) encounters two immediate simplifications:First, we consider the constraint X i = 1 − ˆ n ,i − ˆ n ,i :The density operators are less relevant than the number1 at low energy. Consequently we replace X i →
1. Bythis replacement, we effectively drop the local constraintfor the atoms and dimers. Physically, this is justifiedfrom the fact that infrared fluctuations with wavelengthsmuch larger than the lattice spacing do not resolve singlesites – as stated above, while symmetries provide scaleindependent restrictions on the form of the effective ac-tion, the relevance of the constraint principle depends onscale. Second, we draw the continuum limit. Our origi-nal Hamiltonian (30) often contains bilocal terms. In thequadratic sector, the resulting spatial derivative termsare kept: they describe spatial propagation and may beleading in the infrared for zero mass terms encounteredclose to the phase transition. However, in the interac-tion terms we drop the gradient couplings if they appearin combination with a local one, which in comparison isalways more relevant in the sense of the renormalizationgroup. Finally, we drop the local quartic terms couplingatoms with dimers, which are subleading in comparisonwith the local cubic ones. The corresponding action reads S = S [ b † , b ] + S [ b † , b ] + S int [ b † , b , σ, π ] , (65) S [ b ] = (cid:90) x b † [ ∂ τ − µ − ( − µ + U ) s ) − J (1 + | s | )( z + (cid:52) )] b − √ Jcs (cid:0) b ( z + (cid:52) ) b + c.c. (cid:1) ,S [ b ] = (cid:90) x b † [ ∂ τ + ( − µ + U )( c − s )] b ,S int [ b , b ] = (cid:90) x Jz (cid:2) cs (cid:0) b † + b (cid:1) ˆ n − √ c − s ) (cid:0) b † b + c.c. (cid:1)(cid:3) . ( (cid:82) x = (cid:82) dτ d d x, x = ( τ, (cid:126)x ) , (cid:52) the Laplace operator. We omit the ( τ, (cid:126)x ) dependence of the field for brevity.) Here8and in the following we have chosen s real without lossof generality.In the next step we identify the relevant phase fluc-tuations. The terms in the action (65) are seen to be intwo classes: The first one is made up of field combinationswhich transform according to U (1) × U (1) (the first phaserotation acts on b and the second on b ), i.e. they do notlock the phases. In contrast, the cubic interaction termsin the last line of Eq. (65) lock the phases such that theresidual symmetry is a single U (1). (Such a mechanismbreaking U (1) × U (1) → U (1) is a consistency check forour theory, which emerges from a constrained version ofthe Bose-Hubbard model, in turn only possessing a sin-gle U (1) symmetry.) The dominant temporal and spatialphase fluctuations thus originate from the vicinity of thephase constraint emerging from the phase locking of theatomic to the dimer phase, θ ( x ) = 2 θ ( x ). To bring out the physics of these fluctuations, it is convenient toperform a local gauge transformation on the b field suchas to absorb the θ fluctuations [7]. Here we work incartesian coordinates for the fluctuating fields, and con-sequently the gauge transformation is realized linearly.The gapless phase fluctuations of the dimer field are rep-resented by its imaginary part, b ( x ) = ( σ ( x )+i π ( x )) / √ b ,we introduce dressed fields according to b ( x ) → ˜ b ( x ) = b ( x )(1 − i κ π ( x )) , (66) κ = ( c − s ) / (2 √ cs ) . Now the gauge transformed action can be calculated. Inthis expression, we only keep leading terms which areaffected at linear order in the infinitesimal rotation. Theresult is S [˜ b † , ˜ b ] = 12 (cid:90) x (˜ b † , ˜ b ) (cid:32) ∂ τ + m − J (1 + s ) (cid:52) − √ Jcs ( z + (cid:52) ) − √ Jcs ( z + (cid:52) ) − ∂ τ + m − J (1 + s ) (cid:52) (cid:33) (cid:32) ˜ b ˜ b † (cid:33) , (67) S int [˜ b † , ˜ b , σ, π ] = (cid:90) x √ iκ∂ τ π ˜ b † ˜ b + Jzσ (cid:104) √ cs ˜ b † ˜ b − ( c − s ) (cid:0) ˜ b ˜ b + ˜ b † ˜ b † (cid:1)(cid:105) with m ≈ | U | / − Jz (1 + s ), using − µ ≈ | U | / n ≈
0, cf. Sec. III). As a preparation for the eliminationof the massive modes, we further introduce hermitianfields for the single particle excitations ˜ b = ϕ + i ψ , suchthat the action reads S [ ϕ, ψ ] = 12 (cid:90) x ( ϕ, ψ ) (cid:32) m − ξ (cid:52) i ∂ τ − i ∂ τ m − − ξ − (cid:52) (cid:33)(cid:32) ϕψ (cid:33) ,S int [ ϕ, ψ, σ, π ]= (cid:90) x i κ √ ∂ τ π ( ϕ + ψ ) + σ ( λ − ϕ + λ + ψ ) , (68)where ξ ± = J ( c ± √ s ) , m ± = −| U | / − zξ ± , λ ± = Jz (3 cs/ √ ± ( c − s )), and again we keep only leadingterms. As appropriate for the phase mode, the field π interacts with the atomic fields ϕ and ψ only throughits time derivative, while the field σ interacts directly.Note that m > m − , and upon approaching the phasetransition, m hits zero prior to m − [6, 7]. Indeed thecondition m = 0 coincides with Eq. (45) if we also usethe mean field equation of state n = 2 s . For vanishing m + , we then find m − ≈ √ csJz = 4 (cid:112) n (1 − n/ Jz within the mean field approximation for the high energyphysics. Hence, the field ψ (the imaginary part of theatomic field ˜ b ) remains massive for any density 0 < n < ϕ becomes soft and plays the role ofan Ising field. The resulting effective action for the fields ϕ , π , and σ reads S [ ϕ, π, σ ]= (cid:90) x (cid:110) σ ( M − ξ σ ∆) σ + i σ∂ τ π + i κ √ ∂ τ πϕ + ζ ( ∂ τ π ) − ξ π ∆ π + 12 ϕ ( m − Z ϕ ∂ τ − ξ ∆) ϕ (cid:111) , (69)where M ∼ λ − /m − , ζ ∼ κ /m − , ξ ∼ ξ σ ∼ λ − ξ − /m − ,and Z ϕ ∼ m − − . (Note that in the limit | U | (cid:29) Jz bothfields ϕ and ψ are massive and, after integrating them outperturbatively, we get Eq. (46) for the effective dimerHamiltonian of the b field.) The field σ now becomesmassive and can be integrated out as well. The finaleffective action for the fields ϕ and π is S eff [ ϕ, π ] = (cid:82) x (cid:110) ϕ ( − Z ϕ ∂ τ − ξ ∆ + m ) ϕ + λϕ + π ( − Z∂ τ − ξ ∆) π + i κ √ ∂ τ πϕ (cid:111) , (70)with Z ∼ M − and λ ∼ λ /M . This action describes acoupled theory for the Goldstone mode π and the Isingmode ϕ . Note that here we also keep a fourth order Isingcoupling. Its presence being rooted in the tree-level σ exchange, this coupling is positive. Thus, the low energytheory contains an Ising part, i.e. a real field with quarticpotential which exhibits Z symmetry breaking when m turns negative. If this part of the action were isolated,9the transition would be in the Ising universality class,and therefore of second order. In the presence of theGoldstone-Ising coupling, more care needs to be taken:In general, a coupling of two (near) gapless real bosonicdegrees of freedom can lead to a fluctuation induced firstorder phase transition, known as the Coleman-Weinbergphenomenon [15]. The Ising self-interaction λ and theIsing-Goldstone coupling κ can be compared via naivepower counting : the canonical dimension of λ is 3 − d ,and that of κ is (3 − d ) /
2. Thus, in any dimension thecorresponding terms have the same degree of relevanceand therefore compete with each other.The form of the action (70) coincides with the one ob-tained in [7] from the continuum Feshbach model. Therenormalization group analysis of the action (70) fornonzero κ has been performed in d = 3 by Frey andBalents [28] at T = 0, and extended to nonzero tempera-ture by Lee and Lee [29], revealing a Coleman-Weinbergphenomenon. Thus, for a generic κ (cid:54) = 0 the phase transi-tion will be driven first order. This scenario is realized inthe low density limits n ≈ ,
2, where our conclusion thusmatches the expectations from the continuum, which wasanticipated in [6, 8] and discussed in detail in [7].However, the lattice offers the possibility to penetratethe regime where n ≈
1. Here, an intriguing situationappears: There exists a point in the phase diagram atwhich the coefficient of the cubic terms vanishes exactly ,which happens due to the zero crossing of the coupling κ . From Eq. (66) we have κ ∼ c − s . Working withthe mean field equation of state n = 2 s , one concludesthat this takes place at n = 1. In reality, renormaliza-tion effects will add contributions to the naive value of κ . Furthermore, inspection of the full equation of state(25) suggests further shifts from the naive expectation,but we have seen in Sec. III B that close to n = 1 theseare small. Thus, we expect the decoupling point to belocated in the close vicinity of the commensurate point n = 1. We provide further evidence for this expectationfrom a symmetry argument in the next section. B. Symmetry argument for the Ising quantumcritical point
The decoupling of Goldstone and Ising mode at a spe-cial point in the phase diagram can also be obtained froma symmetry argument. Being based on a combination ofthe phase locking symmetry between the degrees of free-dom b , b and a temporally local gauge invariance, itcomplements the above explicit derivative expansion andsheds more light on the origin of the decoupling of Isingand Goldstone physics. The power counting applied here is based on the effective rela-tivistic Ising and Goldstone low energy actions with dynamicalexponent z = 1, and not the original nonrelativistic theory. For this purpose, let us first discuss the temporally lo-cal gauge invariance of the Bose-Hubbard Hamiltonian[31] in the presence of an infinite three-body repulsion,which is equivalent to the constrained model under con-sideration here. This adds a local term to the standardBose-Hubbard Hamiltonian, H = lim γ →∞ (cid:2) H BH + γ (cid:88) i ˆ n i (ˆ n i − n i − (cid:3) . (71)The temporally local gauge invariance results simplyfrom the fact that the Hamiltonian is not explicitly timedependent (while it is spatially non-local, such that aspatially local gauge invariance does not exist). Conse-quently, the constrained Bose-Hubbard action must takethe form S c = (cid:90) dτ (cid:88) i ( a † i ( ∂ τ − µ ) a i + H [ a † , a ]) (72)such that the temporally local gauge invariance is ex-pressed as an invariance under a i → exp i λ ( t ) a i , µ → µ + i ∂ τ λ ( t ) . (73)Since our construction must conserve this property, wealso require this invariance for the theory defined with(30). On the level of the effective action and in Fourierspace, this invariance translates into the Ward identityfor the effective action − ∂∂µ δ Γ δb † / ( q ) δb / ( q ) (cid:12)(cid:12)(cid:12) b / =0; q =0 (74)= ∂∂ (i ω ) δ Γ δb † / ( q ) δb / ( q ) (cid:12)(cid:12)(cid:12) b / =0; q =0 , i.e. the coefficient of the linear time derivative mustequal the derivative with respect to the chemical poten-tial. Therefore, in a derivative expansion of the effectiveaction, which is appropriate at low energies, we have:Γ= (cid:90) b † [ z ∂ τ + y ∂ τ + m + ... ] b + (cid:96) ( b † + b ) (75)+ b † [ z ∂ τ + y ∂ τ + m + .... ] b + h ( b † b + b b † ) + ... The presence of a condensate for b , θ (cid:54) = 0, generates off-diagonal terms in the b inverse propagator, i.e. (cid:96) (cid:54) = 0.Here we restrict to the spatially local part of the effec-tive action, since this is the sector where the couplingof Ising to Goldstone mode emerges. The Ward identity(74) implies z / = − ∂m / /∂µ =: g / .Furthermore, using solely the global gauge invariance,we can make the connection between g and g . Indeed,we have a phase locking in the K (21) K (10) † + h.c. term.As a consequence of these terms, the phases of b and b cannot transform independently, and we have b ,i → exp i λb ,i , b ,i → exp 2i λb ,i , (76)0leading to the additional Ward identity2 ∂∂µ δ Γ δb † ( q ) δb ( q ) (cid:12)(cid:12)(cid:12) b / =0; q =0 = ∂∂µ δ Γ δb † ( q ) δb ( q ) (cid:12)(cid:12)(cid:12) b / =0; q =0 , (77)or g = 2 g . In sum, we have the following relations: z = g = 2 z = 2 g . (78)Next we discuss properties of the “compressibility”coupling g ( n ) = − ∂m /∂µ | n , which fixes how stronglythe bound state excitation couples to the chemical po-tential. In the limits n = 0 , n = 0, we find g >
0, while at n = 2we obtain g <
0. These opposite signs can be expected,as at n = 0 the excitations are well-defined dimers, whileat n = 2 we face well defined di-holes. If we do not re-define the chemical potential, then adding a di-hole isenergetically equivalent to delete a dimer. Under themild assumption that the compressibility is a continuousmonotonic function of n (our description is tailored todescribe the DSF phase including the phase border, andtherein we do not expect additional phase transitions),then g ( n ) must have a unique zero crossing. We notethat we should use the above derivative prescription asan operational definition of g ; in principle, there couldbe a µ -independent constant adding to the full mass orgap term of b . For b , such a situation takes actuallyplace and we have an additional mass or gap term U .As a consequence of Eq. (78), a zero crossing of g also implies a zero crossing of the coefficients z , g , z .Thus, the leading frequency dependence is not linear, butquadratic, and the analogous statement is valid in thetime domain, where the leading behavior is a quadrati-cally appearing time derivative.With this result, we now discuss the possible form ofthe coupling of the Ising to the Goldstone mode. Asabove, we decompose linearly into massive and phasemode, and absorb the phase fluctuations into dressed b fields, b → ˜ b = b (1 − i ˜ κπ/ / ) , ˜ κ = h/(cid:96) . Indeed thelow energy effective action can only depend on derivativecouplings associated to the phase mode π , due to theglobal U (1) invariance under transformations π → π + λ .The transformaton cancels the cubic term in Eq. (75)associated to phase fluctuations, while the contributionassociated to the real part σ can be dropped at low en-ergies since the amplitude is massive. At the same time,the b part in the dressed frame now readsΓ = (cid:90) ˜ b † [ z ∂ τ + y ∂ τ + m + ... ]˜ b (79)+i z ˜ κ∂ τ π ˜ b † ˜ b + i y ˜ κ∂ τ π ˜ b † ˜ b + ... Thus, for g = 0, Eq. (78) also implies that the cu-bic derivative coupling z ˜ κ with canonical dimension(3 − d ) / quadratic time derivative. This coupling hascanonical dimension (1 − d ) /
2, and thus is irrelevant near a Gaussian fixed point for d >
1. Similarly, a potential U (1) symmetric coupling term g (cid:48) (cid:82) ( ∂ τ π ) φ has canon-ical dimension 1 − d . Both therefore do not lead to aColeman-Weinberg phenomenon. In consequence, Gold-stone and Ising physics effectively decouple at low ener-gies, giving rise to a second order Ising transition.We summarize our result. Based on the zero crossingof g , phase locking and temporally local gauge invari-ance we find:(i) At the zero crossing point, the nonrelativistic timederivative terms vanish. In the sense of a derivative ex-pansion, the next relevant term is ∂ τ , in which case thetheory acquires a relativistic space-time isotropy in a d +1dimensional space-time. This is physically sound, as thispoint has a special kind of (di-)particle-hole symmetry, inthat the hybrid excitation consists of a superposition of“dimers” and “di-holes” to equal parts. However, we notethe absence of a particle-hole symmetry in the conven-tional sense – such a situation only occurs in the pertur-bative limit J/ | U | →
0, as discussed in Sec. IV. Beyondthe leading order perturbation theory, this symmetry isbroken. One manifestation of the absence of this sym-metry is the asymmetry of the critical line in the phasediagram, cf. Fig. 2.(ii) The cubic coupling of Goldstone to Ising mode alsovanishes at this point. Only terms which are irrelevantin d > n = 1, there is a d + 1 dimensional Ising quan-tum critical point. Examples of physical realizations ofIsing quantum critical points in nature are actually rare.Several systems exhibit Ising type phase transitions withdiscrete symmetry breaking, like the ASF-DSF transitionin the continuum Feshbach model [6, 8] and or a transi-tion between superconductors with different pairing sym-metries [30], but in these cases in the long wavelengthlimit a Coleman-Weinberg phenomenon takes place. Acubic coupling of the Goldstone mode with linear timederivative to the Ising density is actually quite generic innonrelativisic systems, where the Ising mode emerges asan effective degree of freedom describing the transitionfrom one ordered phase to the other. Here we have iden-tified a mechanism that suppresses this coupling. Oneof the few other examples for Ising quantum criticalityis possibly provided by the model magnet LiHoF [32],though the issue is debatable due to the long range in-teractions in the material, preventing an exact mappingto the Ising model.The fact that qualitative aspects of the critical behav-ior are changed in the vicinity of the particle-hole sym-1metric point n = 1 bears some resemblance to the physicsat the tip of the Mott lobe in the repulsive Bose-Hubbardmodel. There, the behavior changes from the nonrela-tivistic O (2) (or XY) universality class with dynamicalexponent z = 2 to the relativistic O (2) model with z = 1[33]. C. Estimate of the Correlation Length
To get an impression of the perspective to observe Isingquantum criticality in this system experimentally, we es-timate the correlation length. This quantity is accessiblewith current experimental technology [34], and has beenmeasured in continuum Bose gases to characterize criticalbehavior.The Coleman-Weinberg phenomenon manifests itselfin the presence of “runaway” trajectories on the RGflow diagram. We therefore can estimate the correlationlength at the first order phase transition as a scale l ∗ ,at which the runaway trajectory with the correspondinginitial conditions hits the boundary of the stability re-gion of the system [17]. The instability is characterizedby the quartic Ising coupling λ turning negative, i.e. thecondition λ ( l ∗ ) = 0.The scaling properties of the action (70) are deter-mined by three parameters: V = ( ξ/ξ + ) (cid:112) Z ϕ /Z , U =4! λ/ξ (cid:112) Z ϕ , and K = κ ξ /ξ Z (cid:112) Z ϕ with the corre-sponding RG equations derived in Ref. [28]. The quan-tity V scales to zero, therefore we can put V = 0 fromthe very beginning. Then the RG equations for the re-maining constants K and U read1 K dKdl = ε − U − K, (80)1 U dUdl = ε − U − K − K U , (81)where ε = 3 − d . To solve these equations, we first in-troduce new functions k = K exp( − εl ) , u = U exp( − εl ),and a new variable x = exp( εl ). The equations then havethe following form, ε dkdx = − (cid:18) uk + 52 k (cid:19) , (82) ε dudx = − (cid:18) u + 6 uk + 24 k (cid:19) . (83)Writing u = kf ( k ) and, therefore du/dk = f + kf (cid:48) , weobtain k dfdk = du/dxdk/dx − f = f + 28 f + 1922 f + 20 = ( f + 12)( f + 16)2( f + 10) . (84)This equation can easily be solved with the result kk = (cid:18) f + 16 f + 16 (cid:19) (cid:18) f + 12 f + 12 (cid:19) , (85) where f = u /k is the initial value for the function f when k = k .It follows from Eq. (83) that ε dudx = − k (cid:2) sf + 48 f + 192 (cid:3) = ε ddx [ kf ( k )] = ε dfdx (cid:20) f dkdf + k (cid:21) and, after using Eq. (84), we obtain ε dfdx = − k f + 16)( f + 12)= − k ( f + 12) ( f + 16) ( f + 16) . The solution of this equation reads − x − ε ≡ − exp( εl ) − ε = 83 k f + 12 (cid:34) − (cid:18) f + 16 f + 16 (cid:19) (cid:35) (86)and, together with Eq. (85), provides a general solutionof the RG equations (82) and (83) and, therefore (80)and (81).The above solution allows us to find the scale l ∗ , atwhich the RG flow reaches the border of stability, U ( l ∗ ) =0. In 3 D we obtain (after taking the limit ε = 3 − d → − l ∗ = 83 k ( f + 12) (cid:34) − (cid:18) f (cid:19) (cid:35) with f = u /k . As a result, close to the Ising criticalpoint, k →
0, we get l ∗ ∼ k − ∼ (1 − n ) − .This result indicates a rather broad critical domain indensity around the true Ising critical point, in which thecorrelation length extends over the whole system. Suchextended quasi-critical behavior can be expected for afluctuation induced first order transition, which resultsexclusively from the competition of very long wavelengthdegrees of freedom, and therefore should be weak. Forexample, already at filling n = 1 / , − / / , / a priori decoupled atomic and dimer excita-tions.2 VI. CONCLUSION
In this paper, we have performed a detailed analyti-cal investigation of the phase diagram of the attractivelattice Bose gas with a 3-body hardcore constraint. Forthis purpose, we make use of a method presented in [11]which allows to exactly map the constrained model to atheory for two unconstrained bosonic degrees of freedomwith conventional polynomial interactions. Within thisframework, we particularly focus on effects tied to inter-actions, which cannot be addressed within a mean fieldplus spin wave approach. While our analysis confirms therough features of the phase diagram obtained from a sim-ple mean field approach – the presence of an Ising-typephase transition from an atomic to a dimer superfluid,numerous interaction driven effects are identified. Thesearise on various length scales, ranging from the fluctua-tion induced formation of the dimer (or di-hole) boundstate on top of the vacua at n = 0 and n = 2 on themicroscopic level over a an understanding of the beyondmean field effects causing nonuniversal shifts in the phase boundary and giving rise to the proximity of the systemto a bicritical point with enhanced SO (3) symmetry instrong coupling, down to the assessment of the true na-ture of the phase transition at very long wavelength. Thisunderpins the fact that short and long range correlationscan then be treated within a unified formalism. Acknowledgements – We thank E. Altman, A. Auer-bach, H. P. B¨uchler, M. Fleischhauer, M. Greiter, A. Mu-ramatsu, N. Lindner, J. M. Pawlowski, L. Radzihovsky,S. Sachdev, J. Taylor and C. Wetterich for interestingdiscussions. This work was supported by the AustrianScience Foundation (FWF) through SFB F40 FOQUS,and project I118 N16 (EuroQUAM DQS), by the Euro-pean union via the integrated project SCALA, by theAustrian Ministry of Science BMWF via the UniInfras-trukturprogramm of the Forschungsplattform ScientificComputing and of the Centre for Quantum Physics, bythe Russian Foundation for Basic Research, and by theArmy Research Office with funding from the DARPAOLE program. [1] N. Syassen et al. , Science , 1329 (2008).[2] J. J. Garcia-Ripoll et al. , New J. Phys. , 013053 (2009).[3] A. J. Daley, J. Taylor, S. Diehl, M. Baranov, P. Zoller,Phys. Rev. Lett. , 040402 (2009); Erratum ibid. ,179902 (2009) .[4] M. Roncaglia, M. Rizzi, J. I. Cirac, Phys. Rev. Lett. ,096803 (2010).[5] Y.-J. Han et al. , Phys. Rev. Lett. , 070404 (2009); A.Kantian et al. , Phys. Rev. Lett. , 240401 (2009).[6] L. Radzihovsky, J. I. Park, P. B. Weichman, Phys. Rev.Lett. , 160402 (2004).[7] L. Radzihovsky, P. B. Weichman, J. I. Park, Ann. Phys. , 2376 (2008).[8] M. Romans, H. Stoof, S. Sachdev, Phys. Rev. Lett. ,020405 (2004).[9] H. C. N¨agerl, private communication (2009).[10] S. Diehl, M. Baranov, A. Daley, P. Zoller, Phys. Rev.Lett. , 165301 (2010).[11] S. Diehl, M. Baranov, A. Daley, P. Zoller, Phys. Rev. B , 064509 (2010).[12] M. Kohno, M. Takahashi, Phys. Rev. B , 3212 (1997).[13] G. G. Batrouni, R. T. Scalettar, Phys. Rev. Lett. ,1599 (1999); F. Hebert, G. G. Batrouni, R. T. Scalettar,G. Schmid, M. Troyer, A. Dorneich, Phys. Rev. B ,014513 (2002).[14] S.-C. Zhang, Phys. Rev. Lett. , 120 (1990).[15] S. Coleman, E. Weinberg, Phys. Rev. D , 292 (1974); D. J. Bergman and B. I. Halperin,Phys. Rev. B , 2145 (1976).[16] E. Altman and A. Auerbach, Phys. Rev. Lett. , 250404(2002).[17] D. Amit and V. Martin-Mayor, Field Theory, the Renor-malization Group and Critical Phenomena: Graphs toComputers , World Scientific Publishing Company (2005).[18] S. D. Huber, E. Altman, H. P. B¨uchler, and G. Blatter, Phys. Rev. B , 085106 (2007); S. D. Huber, B. Theiler,E. Altman, and G. Blatter, Phys. Rev. Lett. , 050404(2008).[19] V. W. Scarola, S. Das Sarma, Phys. Rev. Lett. , 033003(2005); V. W. Scarola, E. Demler, S. Das Sarma, Phys.Rev. A , 051601 (2006); C. Trefzger, C. Menotti, andM. Lewenstein, Phys. Rev. Lett. , 035304 (2009); L.Radzihovsky and S. Choi, Phys. Rev. Lett. , 095302(2009).[20] M. E. Fisher, D. R. Nelson, Phys. Rev. Lett. , 1350(1974).[21] D. Petrosyan, B. Schmidt, J. R. Anglin, M. Fleischhauer,Phys. Rev. A , 033606 (2007); B. Schmidt, M. Bortz,S. Eggert, M. Fleischhauer, D. Petrosyan, Phys. Rev. A , 063634 (2009).[22] E. Altman, E. Demler, and M. D. Lukin, Phys. Rev. A , 013603 (2004).[23] P. T. Ernst, S. G¨otze, J. S. Krauser, K. Pyka,D.-S. L¨uhmann, D. Pfannkuche, K. Sengstock,arXiv:0908.4242 (2009).[24] A. A. Burkov, A. Paramekanti, Phys. Rev. Lett. ,255301 (2008); R. Ganesh, A. Paramekanti, A. A.Burkov, Phys. Rev. A , 043612 (2009).[25] G. Vidal, Phys. Rev. Lett. , 040502 (2004); F. Ver-straete, V. Murg, and J. I. Cirac, Adv. Phys. , 143(2008).[26] A. J. Daley, S. R. Clark, D. Jaksch, and P. Zoller, Phys.Rev. A. , 043618 (2005).[27] A. J. Daley, C. Kollath, U. Schollw¨ock, and G. Vidal, J.Stat. Mech.: Theor. Exp. P04005 (2004); S.R. White andA.E. Feiguin, Phys. Rev. Lett. , 076401 (2004).[28] E. Frey, L. Balents, Phys. Rev. B Quantum Phase Transitions , Cambridge Uni- versity Press, Cambridge (1999).[32] D. Bitko, T. F. Rosenbaum, G. Aeppli, Phys. Rev. Lett.
940 (1996). [33] M. P. Fisher, P. B. Weichman, G. Grinstein, and D. S.Fisher, Phys. Rev. B , 546 (1989).[34] T. Donner, S. Ritter, T. Bourdel, A. ¨Ottl, M. K¨ohl, andT. Esslinger Science315