# Green’s Function Approach to the Bose-Hubbard Model

###### Abstract

We use a diagrammatic hopping expansion to calculate finite-temperature Green functions of the Bose-Hubbard model which describes bosons in an optical lattice. This technique allows for a summation of subsets of diagrams, so the divergence of the Green function leads to non-perturbative results for the boundary between the superfluid and the Mott phase for finite temperatures. Whereas the first-order calculation reproduces the seminal mean-field result, the second order goes beyond and shifts the phase boundary in the immediate vicinity of the critical parameters determined by high-precision Monte-Carlo simulations of the Bose-Hubbard model. In addition, our Green’s function approach allows for calculating the excitation spectrum both for zero and finite temperature and for determining the effective masses of particles and holes.

## I Introduction

Ultracold bosonic gases trapped in the periodic potential of optical lattices represent tunable
model
systems for studying the physics of quantum phase transitions fisher (); jaksch (); bloch2 (); sachdev (). They are described by
the Bose-Hubbard Hamiltonian which decomposes into two parts: . The first term
,
with the on-site energy and the chemical potential , describes the repulsion of more than one boson
residing on a lattice site. It is local and diagonalizable in the occupation number basis for any lattice site.
The second term with the hopping matrix element
if the lattice sites and are nearest neighbors and otherwise
describes the hopping between two sites due to the quantum-mechanical tunneling effect. The competition between the two
energy scales and determines the existence of two different phases. When the on-site energy is small compared to the
hopping amplitude, the ground state is superfluid (SF) as the bosons are delocalized in a phase coherent way over the whole
lattice. In the opposite case, where the on-site interaction dominates over the hopping term, the ground state
is a Mott insulator (MI) where each boson is trapped in one of the respective potential minima.

This characteristic quantum phase transition of the Bose-Hubbard model has been
studied extensively both with analytical monien (); monienpade (); stoof1 (); stoof2 (); santos (); Holthaus5 (); Holthaus6 () and numerical
monienDMRG (); QMC3D (); QMC2D2 () methods for zero temperature, while less literature exists on the finite-temperature
properties of this transition schmidtPRL (); buonsante (); pelster (). In this letter we work out an analytical Green’s
function approach in order to determine the MI-SF phase boundary and the excitation spectrum in the Mott phase
both for zero and finite temperature with a hopping expansion. Our findings compare well with the latest findings
of Quantum Monte Carlo simulations and allow to propose a thermometer for bosons in optical lattices.
In the following we restrict ourselves to a spatially homogeneous system and neglect the effects arising from the
additional harmonic confining potential, which is present in all experimental settings, like the formation of a shell
structure bloch1 (). However, these effects could be taken into account by applying the local density approximation
where the external potential is taken into account in form of a spatially dependent chemical potential.

## Ii Green’s Function Approach

As all single-particle properties of a quantum many-body system are contained in its Green function, we base our calculation on this quantity. Because we are interested in describing a system at non-zero temperature, we use the imaginary-time formalism abrikosov (); zinnjustin (). Therein, the single-particle Green function is defined as the thermal average of the time-ordered product of a creation and an annihilation operator in Heisenberg representation

(1) |

with , and we have introduced the partition function . Because it is not possible to obtain analytic expressions for the eigenstates and eigenenergies of the full Bose-Hubbard Hamiltonian, we cannot calculate the Green function exactly. Instead, we aim at a perturbative treatment and calculate this quantity as a power series in the hopping matrix element . As that parameter is small for the Mott phase, where the lattices are deep and the interaction between particles is strong, we refer to this treatment as a strong-coupling expansion. In order to employ this perturbative expansion for finite temperature we make use of the Dirac interaction picture and write the imaginary-time evolution operator in form of a Dyson series,

(2) |

where the time dependence of the Dirac-picture operators is determined by the local Hamiltonian . With its help we can write the Green function as

(3) |

where the time-ordering operator acts also on the time variables which are resulting from the expansion of the Dirac imaginary-time evolution operator in (2). When we now expand perturbatively in powers of the tunneling matrix element , the th order contribution in (3) turns out to depend on the -particle Green function of the unperturbed system as . In order to determine we cannot use standard Wick’s theorem because the unperturbed Hamiltonian is not quadratic in the Bose operators. Although the lack of Wick’s theorem in the present situation prevents us from using the powerful perturbative technique based on standard Feynman diagrams, we can nevertheless simplify by decomposing it into cumulants. To this end, we follow an approach reviewed by Metzner metzner () in the context of the Hubbard model which describes electrons in a conductor. This decomposition is based on the important observation that the on-site Hamiltonian is local. Consequently, the unperturbed Green functions are also local and can be decomposed into time-dependent cumulants . For instance, we have with the cumulant

(4) |

where is the unperturbed partition function for a single-site system with the on-site energy eigenvalues . With this decomposition, we can represent diagrammatically: We denote an -particle cumulant at a lattice site by a vertex with entering and leaving lines with imaginary-time variables, so we have, for instance, for the first two cumulants

(5) | |||||

(6) |

Furthermore, the hopping matrix element is symbolized by a line connecting two vertices:

(7) |

With all this, we can set up the diagrammatic rules for calculating the th order contribution of the Green function in .
First: Draw all possible combinations of vertices with total internal and one entering and one leaving line.
Second: Connect them in all possibles ways and assign time variables and hopping matrix elements to the lines.
Third: Sum all site indices over all internal lattice sites and integrate all internal time variables from to .

We also note here that we have to sum all site indices over the whole lattice, no matter whether two sites
in a diagram coincide or not. We make use of the translational invariance in
imaginary time and transform all expressions to Matsubara space. In the second diagrammatic rule the integrals
over the time variables have to be replaced by sums over all bosonic Matsubara frequencies with
integer where the sum of the incoming frequencies must equal the sum of the outgoing ones.

The formalism developed so far allows for calculating the Green function to any given order in
the tunneling matrix element . But because one
of our main goals is to describe the phase transition between the Mott insulator and the superfluid phase, and it is
well known in the theory of critical phenomena
that such a transition is characterized by diverging long-range correlations zinnjustin (); verena (), we must
employ a non-perturbative method which is archieved by summing an infinite subset of diagrams. In
order to perform this task, we introduce the sum of all one-particle irreducible diagrams including their respective symmetry
factors and multiplicities as

(8) |

The self-energy which describes the movement of a single particle in a many-body enviroment is defined as
abrikosov (). For our case, it can be written in the
form ,
where is the -dimensional lattice dispersion.
For instance, the first order in the self-energy which corresponds to a summation of all simple chain diagrams,
yields with
.
We note that performing the limit
in the Green function allows to obtain the quasi-momentum distribution needed to explain
experimental time-of-flight pictures blochvisibility (); alexpaper (); bloch3 ().

## Iii Phase Boundary

The boundary of the phase transition is characterized by diverging long-range correlations zinnjustin (); verena (). Thus, we set and solve for the value of where the Green function diverges which is only possible for . This yields up to first order in

(9) |

which coincides with the finite-temperature mean-field result pelster (); buonsante (); gerbier (). The zero-temperature
limit of (9), i.e. , agrees with the seminal mean-field result of
Ref. fisher (). The reason for this agreement is that each approximation becomes exact in the limit of infinite
spatial dimension. In order to see that one must suitably scale the hopping parameter vollhardtferm (); dmft1 (); vollhardtneu ().
When we define , the contribution of the th order chain diagram is proportional to
because there exist possibilities in a chain diagram
for every hopping line to connect to neighboring sites. The lowest-order term neglected by that summation is the one-loop
diagram in (8)
which has two internal lines but only one free index and is, therefore, proportional to . Thus, it
vanishes in the limit of .

Taking now the one-loop diagram into account yields

(10) |

where is the second-order self-energy with the value of the one-loop diagram . The analytic formula for the phase boundary resulting from this formula is shown in Fig. 1. This result coincides in the limit of with the effective potential formalism employed in Ref. santos (). We want to emphasize that, unlike the first order (9), the one-loop corrected result depends on the system dimension in an non-trivial way and that the corrections are larger in two than in three dimension which is consistent with the fact that the first-order results becomes exact in the limit . Note that even higher-order corrections for the effective potential formalism have been obtained in Refs. Holthaus5 (); Holthaus6 (), which turned out to be indistinguishable from the Quantum Monte-Carlo data in Ref. QMC3D ().

For finite temperature, the phase boundary is shifted towards larger values of
as thermal fluctuations suppress quantum correlations which are responsible for the formation of the superfluid. This effect,
which occurs both in first and in second order, is most important between the Mott lobes as fluctuations are strongest when
the average particle number is not near an integer value. The correction to the mean-field result arising from the one-loop
diagram, visible in the difference between first and second order curve, plays an important role only near the tip of the
Mott lobe. This feature stems from the fact that quantum fluctuation are notably increased when the system
approaches the quantum critical point at the tip of the lobe sachdev (). Thus, we can say the quantum phase diagram
consists of a thermally dominated region where the influence of thermal fluctuations is large and a quantum dominated region
where the quantum corrections from the one-loop diagram are most important.

## Iv Excitation Spectrum

For a system in the Mott phase three different types of excitations exist. 1.
The addition of a particle from the environment (particle excitation). 2. The removal of a particle (hole excitation). 3.
The creation of a particle-hole pair (pair excitation). The last one is most important from a physical point of view
because it is also possible in an isolated system as realized in the current experiments. In the superfluid phase,
there are additional excitations corresponding to fluctuations of the phase of the superfluid order parameter. However,
these features can not be investigated within the present formalism but with the related effective action method barry ().

In order to obtain the spectrum of the quasi-particles, which is given by the poles of the real-time Green function
for finite temperature we must analytically continue our imaginary-time result. This is achieved by performing the
replacement with abrikosov (). The first and second order, respectively,
yield excitation spectra of which the former one agrees in the limit with the mean-field result from
Ref. stoof1 (). In the following we restrict ourselves to the most interesting case . Both finite
temperature and one-loop corrections are most effective for small wave numbers, the former because the effect of temperature
is the suppression of quantum correlations which mainly exist for long wavelengths, the latter because the dominant
fluctuations near a quantum critical point are the ones with vanishing wave number as shown in Fig. 2

All spectra show a characteristic gap which vanishes at the critical point, i.e. at the value of at the tip of a
Mott lobe. In Fig. 3 its temperature dependence is shown. As this gap is a quantity which is experimentally accessible
bloch1 (), it could serve as a method to determine the temperature of bosons in an optical lattice.
The quasi-particle can be ascribed an effective mass which is shown in Fig. 4 for the particle- and for the
hole-excitations. They both become massless at the critical point which is a result of the U(1) symmetry breaking at the
second-order phase transition.

## V Conclusion and Outlook

We have presented a powerful formalism to calculate the Green function for the
Bose-Hubbard model in the Mott phase. It allowed us to improve the mean-field phase boundary in an analytic way both for
finite and zero temperature where the former result deviates from recent Quantum Monte-Carlo studies by only 3% for
in Ref. QMC3D (). For finite temperature first analytic results beyond mean-field theory have been
presented and the importance of both thermal and quantum fluctuations in the different regions of the phase diagram has been
clarified. In addition, we have calculated the excitation spectrum of the quasi-particles and derived its effective masses
and compared them to numerical findings. We have also investigated the characteristic energy gap and determined its
temperature dependence.
More applications of this approach have recently allowed in Ref. Grass () to obtain more insight
into the superfluid phase by combining the present technique with the effective action approach from Ref. barry ().

## Acknowledgement

We cordially thank B. Bradlyn, H. Enoksen, A. Hoffmann, H. Kleinert, F. Nogueira, R. Graham, and F.E.A. dos Santos for stimulating and helpful discussions.

## References

- (1) M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. Fisher, Phys. Rev. B 40, 546 (1989).
- (2) D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, Phys. Rev. Lett. 81, 3108 (1998).
- (3) M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch, and I. Bloch, Nature 415, 39 (2002).
- (4) S. Sachdev, Quantum Phase Transitions, Cambridge University Press, 1999.
- (5) J. K. Freericks and H. Monien, Phys. Rev. B 53, 2691 (1996).
- (6) N. Elstner and H. Monien, Phys. Rev. B 59, 12184 (1999).
- (7) D. van Oosten, P. van der Straten, and H. T. Stoof, Phys. Rev. A 63, 053601 (2001).
- (8) D. B. M. Dickerscheid, D. van Oosten, P. J. H. Denteneer, and H. T. C. Stoof, Phys. Rev. A 68, 043623 (2003).
- (9) F. E. A. dos Santos and A. Pelsterr, Phys. Rev. A 79, 013614 (2009).
- (10) N. Teichmann, D. Hinrichs, M. Holthaus, and A. Eckardt, Phys. Rev. B 79, 100503 (2009).
- (11) N. Teichmann, D. Hinrichs, M. Holthaus, and A. Eckardt, Phys. Rev. B 79, 224515 (2009).
- (12) T. D. Kühner and H. Monien, Phys. Rev. B 58, R14741 (1998).
- (13) B. Capogrosso-Sansone, N. V. Prokof’ev, and B. V. Svistunov, Phys. Rev. B 75, 134302 (2007).
- (14) B. Capogrosso-Sansone, Şebnem Güneş Söyler, N. Prokof’ev, and B. Svistunov, Phys. Rev. A 77, 015602 (2008).
- (15) H. Kleinert, S. Schmidt, and A. Pelster, Phy. Rev. Lett. 93, 160402 (2004).
- (16) P. Buonsante and A. Vezzani, Phys. Rev. A 70, 033608 (2004).
- (17) K. V. Krutitsky, A. Pelster, and R. Graham, New J. Phys. 8, 187 (2006).
- (18) S. Fölling, A. Widera, T. Müller, F. Gerbier, and I. Bloch, Phys. Rev. Lett. 97, 060403 (2006).
- (19) A. A. Abrikosov, L. P. Gorkov, and I. E. Dzyaloshinski, Methods of Quantum Field Theory in Statistical Physics, Dover Publications, 1963.
- (20) J. Zinn-Justin, Quantum Field Theory and Critical Phenomena, Clarendon Press, fourth edition, 2002.
- (21) W. Metzner, Phys. Rev. B 43, 8549 (1993).
- (22) H. Kleinert and V. Schulte-Frohlinde, Critical Properties of -Theories, World Scientific, 2001.
- (23) F. Gerbier et al., Phys. Rev. Lett. 95, 050404 (2005).
- (24) A. Hoffmann and A. Pelster, Phys. Rev. A 79, 053623 (2009).
- (25) F. Gerbier, S. Trotzky, S. Foelling, U. Schnorrberger, J. D. Thompson, A. Widera, I. Bloch, L. Pollet, M. Troyer, B. Capogrosso-Sansone, N. V. Prokof’ev, and B. V. Svistunov, Phys. Rev. Lett. 101, 155303 (2008).
- (26) F. Gerbier, Phys. Rev. Lett. 99, 120405 (2007).
- (27) W. Metzner and D. Vollhardt, Phys. Rev. Lett. 62, 324 (1989).
- (28) L. Amico and V. Penna, Phys. Rev. Lett. 80, 2189 (1998).
- (29) K. Byczuk and D. Vollhardt, Phys. Rev. B 77, 235106 (2008).
- (30) T.D. Grass, F.E.A. dos Santos, and A. Pelster, Phys. Rev. A 84, 013613 (2011).
- (31) B. Bradlyn, F. E. A. dos Santos, and A. Pelster, Phys. Rev. A 79, 013615 (2009).