Master equation vs. partition function: canonical statistics of ideal Bose–Einstein condensates


Abstract

Within the canonical ensemble, a partially condensed ideal Bose gas with arbitrary single-particle energies is equivalent to a system of uncoupled harmonic oscillators. We exploit this equivalence for deriving a formula which expresses all cumulants of the canonical distribution governing the number of condensate particles in terms of the poles of a generalized Zeta function provided by the single-particle spectrum. This formula lends itself to systematic asymptotic expansions which capture the non-Gaussian character of the condensate fluctuations with utmost precision even for relatively small, finite systems, as confirmed by comparison with exact numerical calculations. We use these results for assessing the accuracy of a recently developed master equation approach to the canonical condensate statistics; this approach turns out to be quite accurate even when the master equation is solved within a simple quasithermal approximation. As a further application of the cumulant formula we show that, and explain why, all cumulants of a homogeneous Bose–Einstein condensate “in a box” higher than the first retain a dependence on the boundary conditions in the thermodynamic limit.

PACS

  • 05.30.Jp;
  • 03.75.Fi;
  • 05.30.Ch

Keywords

  • Ideal Bose gas;
  • Canonical ensemble;
  • Statistics of occupation numbers;
  • Master equation;
  • Generalized Zeta functions

1. Introduction

It sometimes happens that a problem in physics can be solved by two quite different methods, each one requiring its own propositions and approximations. In such a case it can be rewarding to compare these two approaches in some detail: What is hard to derive within the first approach might be obvious within the second; both methods might lend themselves to different generalizations. In such a comparison it is not so much the result itself that matters, but rather the “best” way to obtain it; the hope is that a thorough understanding of the solution techniques will provide a key for attacking more demanding problems.

Such a situation arises when studying the ideal Bose gas in the canonical ensemble. While the grand canonical treatment of ideal Bose gases constitutes classic textbook material [1], [2] and [3], recent experiments with mesoscopic samples of dilute Bosonic atoms in thermally isolated traps (for a Review, see Ref. [4]) actually call for a microcanonical analysis. The canonical ensemble provides a convenient intermediate step: It still assumes that the trapped gas exchanges energy with a heat bath, as in a grand canonical setting, but it embodies the constraint that the number of gas particles, N, be fixed. Experiments with 4He in a porous medium [5] provide an interesting example of Bose–Einstein condensation in the canonical ensemble.

One approach to the canonical statistics of ideal Bose gases, presented in Ref. [6] and developed further in Ref. [7], consists in setting up a master equation for the condensate and finding its equilibrium solution. This approach has the merit of being technically lucid and physically illuminating, since it reveals important parallels to the quantum theory of the laser [8]. For deriving that master equation, in the approximation of detailed balance in the excited states, it was assumed that given an arbitrary number n0 of atoms in the condensate, the remaining Nn0 excited atoms are in an equilibrium state at the prescribed temperature T—or, in other words, that the occupation numbers of the excited-states subsystem thermalize significantly faster than the condensed ones, no matter how many particles the condensate contains. Whereas the existence of two separate time scales is evident in laser physics, where the atomic dynamics usually are fast compared to the escape of photons from the cavity, the physical motivation for such a separation of time scales seems less obvious in the present case. However, it was shown [9] that indeed the time scale of long-range coherent ordering is much greater than the collisional time scale responsible for the kinetic equilibration of above-condensate atoms. Thus, while this assumption might not be sufficient for describing the details of the process of equilibration, it does lead, in principle, to the correct canonical-ensemble equilibrium state [10]. Still, for evaluating that equilibrium state one has to make certain additional approximations, as detailed in Section 2, so that it is desirable to check conclusions drawn from the master equation approach against the results provided by independent techniques.

A classic technique that almost seems to suggest itself for studying canonical statistics is the saddle-point method: Starting from the known grand canonical partition function one employs the saddle-point approximation for extracting its required canonical counterpart, which then yields all desired quantities by taking suitable derivatives. It turns out that this program requires some caution, since the customary form of the saddle-point approximation, as advocated by Schrödinger in his famous treatise [11], is not correct in the condensate regime; the close approach of the saddle-point to the ground-state singularity of the grand partition function forbids the usual Gaussian approximation. Fortunately, this difficulty can be overcome with the help of a strategy that has been suggested by Dingle [12] and worked out in detail in [13] and [14]; this approach furnishes canonical occupation numbers and their fluctuations, say, for reasonably large, experimentally relevant particle numbers N and all temperatures. However, accurately solving the equation for the saddle-point requires some numerical skills.

In this paper, therefore, we will follow still another avenue for checking the master equation approach. Instead of aiming for approximations valid at all temperatures—as provided by the proper saddle-point method—we will restrict ourselves to temperatures below the onset of Bose–Einstein condensation. In this case there exists a transparent approximation which allows us to bring the canonical partition function into a simple form—namely into that of a system of independent harmonic oscillators [15], [16], [17] and [18]. Using this, we derive a formula which expresses all the cumulants κk(β) of the canonical distribution pN(ex)(M;β)—i.e., the probability distribution for finding M of the N atoms in an excited state at inverse temperature β, and hence n0=NM atoms in the condensate—in terms of the poles of a generalized Zeta function determined by the system's single-particle energies. This result is of substantial interest, since systems that lend themselves to an explicit calculation of all higher moments beyond the usual Gaussian second order are quite rare in physics. Our approach enables us, first of all, to make contact with the master equation, and to disperse possible objections against it: The first and the second moment of pN(ex)(M;β), as obtained from the master equation within a “quasithermal” approximation, merge exactly into the corresponding results provided by the bona fide partition function; the higher moments are systematically approximated in a mean field-like fashion. The Zeta-function technique also allows us to clarify an important issue: In the case of the homogeneous Bose gas “in a box”, the “higher” statistical properties (i.e., the cumulants κk(β) with order k⩾2) do depend on the particular boundary conditions even in the thermodynamic limit. When it comes to a detailed discussion of the statistics of occupation numbers, simply taking the homogeneous Bose gas with the convenient periodic boundary conditions might, therefore, not be sufficient.

We restrict this paper to the discussion of orthodox equilibrium statistical mechanics and do not consider the method of canonical quasiparticles, which was introduced [17] and [18] on the basis of a particle-number conserving formalism established by Girardeau and Arnowitt [19] and [20]. That method is very efficient for the solution of various non-equilibrium problems in the canonical ensemble, such as the dynamics of condensate formation, vortex evolution, and the kinetics of correlations, both for the ideal and the interacting Bose gas.

The main part of this paper is organized as follows: In Section 2, we summarize the master equation approach [6] and [7], to the extent that it will be needed later. In Section 3, we then explain why the canonical partition function of the ideal Bose gas reduces in the condensate regime to a partition function of independent harmonic oscillators, and quantify the accuracy of the master equation by comparing the canonical cumulants κk(β) obtained from this approach with those following from the “independent oscillator” point of view, and with results of exact numerical calculations. In Section 4, we elaborate the connection between the cumulants and the poles of the generalized Zeta functions, while Section 5 contains the explicit evaluation of our general formula for an isotropic harmonic trapping potential, and for the “box” potential with periodic or Dirichlet boundary conditions, respectively. The detailed analysis of the higher cumulants, in particular, of those of third and fourth order (skewness and flatness), yields interesting insight into the non-Gaussian nature of the condensate fluctuations. A brief discussion concludes the paper in Section 6. Since the required techniques for handling the generalized Zeta functions might not be generally known, and in order to keep the paper self-contained, Appendix A offers the relevant mathematical details.

2. Essentials of the master equation approach

We consider an ideal Bose gas that consists of N particles and is stored in some arbitrary trap. The trapping potential determines the discrete single-particle energies εν; the index ν=0,1,2,… labels the individual single-particle eigenstates. It is further assumed that this system is kept in thermal contact with a harmonic-oscillator heat bath of temperature T; at sufficiently low T a fraction of the particles undergoes Bose–Einstein condensation and occupies the ground state ν=0. We require that the spectral density of the bath be flat, so that for each transition of a gas particle, involving an energy difference ενεσ, there is always the same number of frequency-matched oscillators which can provide or accept this quantum. The average occupation number of such a heat-bath oscillator is written as

equation1
Full-size image (<1 K)
where β=1/kBT, with kB denoting Boltzmann's constant.

Starting from the equation of motion for the total density matrix, tracing out the bath degrees of freedom in the usual manner, and making the Markov approximation of very fast relaxation in the bath compared to the dynamics of the system of excited and condensed Bose particles in the trap, one obtains a master equation for the Bose gas. This equation can be further reduced to a master equation for the condensate alone, if one assumes that the state of the excited particles is determined by detailed balance [6] and [7]. The condensate master equation describes the evolution of the probability pn0 that the condensate contains n0 out of the N particles, and includes information about the excited-state occupations only through the “cooling coefficients”

equation2
Full-size image (<1 K)
and the “heating coefficients”
equation3
Full-size image (<1 K)
where 〈nνn0 is the canonical expectation value for the number of particles occupying the νth excited state, subject to the condition that there be n0 condensate atoms. In terms of these coefficients, the master equation reads
equation4
Full-size image (<1 K)
The constant κ obviously carries the dimension of an inverse time. It embodies the spectral density of the bath and the coupling strength of the bath oscillators to the gas particles [7], and determines the rate of condensate evolution since there is no direct interaction between the particles of an ideal Bose gas. Here, as in Ref. [7], we study the equilibrium solution of Eq. (4), so that the question how fast the equilibration proceeds is not important.

The content of this master equation (4) can be grasped in an intuitive manner, even without a formal derivation. The cooling coefficient Kn0 describes all those processes that add another particle to a condensate consisting so far of n0 atoms: The new particle originates from one of the excited states ν⩾1, each of them endowed with the occupation number 〈nνn0. The particle's transition to the ground state is accompanied by the emission of a phonon of energy ενε0 into the heat bath, as expressed by the additional factor (ην0+1) in Eq. (2). Likewise, the heating coefficient (3) summarizes all those processes through which an n0-particle condensate loses one of its constituents: By absorbing a phonon with the appropriate energy from the heat bath (—multiplicity factor ην0—), the particle can make a transition to any of the excited states. Since, prior to the arrival of the new particle, the occupation number of such a state was 〈nνn0, the factor (〈nνn0+1) properly accounts for the “emission” of the new particle into this mode. As indicated in Fig. 1, there are four contributions that affect the rate of change of pn0: an (n0−1)-particle condensate acquiring an additional particle through cooling, an (n0+1)-particle condensate losing one particle through heating, and an n0-particle condensate either picking up or losing one particle; these four contributions make up the right-hand side of Eq. (4).

Visualization of the condensate master equation (4). There are four ...
Fig. 1. 

Visualization of the condensate master equation (4). There are four contributions that affect the rate of change of pn0, the probability of finding n0 particles in the condensate: An (n0−1)-particle condensate (upper horizontal line) gaining one particle through cooling (upper downward arrow), an (n0+1)-particle condensate (lower horizontal line) losing one particle through heating (lower upward arrow), and an n0-particle condensate gaining or losing one constituent. These four contributions add up to give the right-hand side of Eq. (4).

The above discussion also highlights an important difference between the canonical ensemble studied here and a microcanonical ensemble. In a canonical setting, the energy of each transition from or to the condensate is provided or accepted by the heat bath, implying that such transitions do not lead to correlations among the other gas particles. This is no longer the case under microcanonical conditions; here each amount of energy associated with a particle's transition to or from the condensate has to be compensated by the other particles themselves. Assuming the presence of an external heat bath enormously simplifies the analysis, but it also affects the results: For instance, microcanonical fluctuations of the number of condensate particles are smaller than their canonical counterparts [16] and [21].

For evaluating the master equation (4), one now needs an approximation to the cooling and heating coefficients  and . In Ref. [7], a quasithermal approximation for the conditional occupation numbers of the excited states was suggested:

equation5
Full-size image (<1 K)
where
equation6
Full-size image (<1 K)
Summing Eq. (5) over ν⩾1, one finds
equation7
Full-size image (<1 K)
so that this ansatz incorporates the particle number constraint in an exact manner; moreover, 〈nνn0 is taken to be proportional to the phonon occupation number ην0. Within this approximation, the cooling coefficient (2) is given by
equation8
Full-size image (<1 K)
where the “cross-excitation parameter” [7]
equation9
Full-size image (<1 K)
has been introduced. Likewise, Eq. (3) now yields
equation10
Full-size image (<1 K)
With these approximations  and , the steady-state solution
equation11
Full-size image (<1 K)
to the condensate master equation (4) takes the form
equation12
Full-size image (<1 K)
with ZN=1/pN fixed through the normalization condition Full-size image (<1 K), giving the formal representation
equation13
Full-size image (<1 K)
of the canonical partition function. Interestingly, the distribution (12) closely resembles a negative binomial distribution; the only difference lies in the fact that here Nn0 is restricted to integers in the range from 0 to N, whereas all non-negative integers figure in the true negative binomial case. All moments of this distribution (12) can be calculated analytically. In particular, one finds
equation14
Full-size image (<1 K)
as the canonical expectation value of the number of condensate particles; its variance
equation15
Full-size image (<1 K)
the third centered moment
equation16
Full-size image (<1 K)
and the fourth centered moment
equation17
Full-size image (<1 K)
Thus, assuming the validity of the quasithermal approximation (5), the master equation (4) provides a complete description of the canonical statistics of ideal Bose–Einstein condensates.

For later reference, let us specialize the moments , ,  and  to the condensate regime, that is, to temperatures so low that p0, the probability for finding no particle at all in the ground state, is practically zero. Recall first that within the grand canonical ensemble the occupation numbers take the form 〈nνgc=[z−1exp(βεν)−1]−1, with the fugacity z being tied to the ground state energy, z≈exp(βε0). Since for large N these grand canonical occupation numbers equal their canonical counterparts 〈nν〉 (where the absence of the index n0 distinguishes these expectation values from the conditional occupation numbers entering the heating and cooling coefficients  and ), one has

equation18
Full-size image (<1 K)
so that the parameters Full-size image (<1 K) and η, defined in  and , can be expressed through these occupation numbers 〈nν〉 as
equation19
Full-size image (<1 K)
Accordingly, taking the limit p0→0, we obtain for the centered moments
equation20
Full-size image (<1 K)
equation21
Full-size image (<1 K)
equation22
Full-size image (<1 K)
equation23
Full-size image (<1 K)
Eq. (20) is self-evident, but Eq. (21) contains an important insight: Since 〈nν〉(〈nν〉+1)=Δnν2 is the fluctuation of the occupation number nν of the νth excited state (ν⩾1), Eq. (21) expresses the fact that in the condensate regime these occupation numbers are uncorrelated stochastic variables; one simply has to add up their variances to obtain the variance Δn02 of the ground-state occupation number Full-size image (<1 K). This finding is a first hint that within the canonical ensemble there exist independent degrees of freedom in the condensate regime; this will become more clear from the “independent oscillator” point of view outlined in the following section.

3. Cumulants for the ideal Bose–Einstein condensate

The canonical condensate distribution (12), and the moments , ,  and  derived from it, hinge on the validity of the quasithermal approximation (5). This ansatz appears quite reasonable, and it manifestly takes the constraint of fixed particle number into account, but its level of accuracy still needs to be quantified. In order to assess the accuracy of this quasithermal approximation, we will now tackle the problem of canonical condensate statistics from a different angle, staying strictly within the framework of orthodox equilibrium statistical mechanics, and trying to avoid any uncontrolled approximation. The price to pay for the mathematical rigor that will be attained in this and the following sections is that we have to remain restricted to temperatures below the onset of Bose–Einstein condensation, whereas the distribution (12) applies at all temperatures.

We start from the familiar representation [1], [2] and [3]

equation24
Full-size image (<1 K)
of the canonical N-particle partition function, where the prime indicates that the summation runs only over those sets of occupation numbers {nν} that comply with the constraint ∑νnν=N. Singling out the ground-state energy 0, we write the total energy of each such configuration as
equation25
Full-size image (<1 K)
so that E denotes the true excitation energy of the respective configuration. Grouping together configurations with identical excitation energies, we have
equation26
Full-size image (<1 K)
with Ω(E,N) denoting the number of microstates accessible to an N-particle system with excitation energy E, that is, the number of microstates where E is distributed over N or less particles: Ω(E,N) counts all the configurations where E is concentrated on one particle only, or shared among two particles, or three, or any other number M up to N.

Detailed understanding of canonical (or microcanonical) statistics, however, requires knowledge of the number of microstates with exactlyM excited particles, for all M up to N; these numbers are provided by the differences

equation27
Φ(E,M)≡Ω(E,M)−Ω(E,M−1),M=0,1,2,…,N.
By definition, Ω(E,−1)=0. Given Φ(E,M), the canonical probability distribution for finding M excited particles (and, hence, NM particles still residing in the ground state) at inverse temperature β is determined by
equation28
Full-size image (<1 K)
This distribution is merely the mirror image of the condensate distribution pn0 considered in the previous section, i.e., we have pN(ex)(M;β)=pNM. Since the definition (27) obviously implies
equation29
Full-size image (<1 K)
it is properly normalized,
equation30
Full-size image (<1 K)
In the following, we will derive a general formula which gives all cumulants of the canonical distribution (28), provided the temperature is so low that a significant fraction of the particles occupies the ground state—that is, provided there is a condensate.

To this end, we recall that the set of canonical M-particle partition functions ZM(β) is generated by the grand canonical partition function, which has a simple product form [1], [2] and [3]:

equation31
Full-size image (<1 K)
here, z is a complex variable. However, every single M-particle partition function enters into this expression with its own, M-particle ground-state energy 0. In order to remove these unwanted contributions, we define a slightly different function Ξ(β,z) by multiplying each ZM(β) by (zeβε0)M, instead of zM, and then summing over M, obtaining
equation32
Full-size image (<1 K)
On the other hand, in view of Eq. (26) this function also has the representation
equation33
Full-size image (<1 K)
Therefore, multiplying Ξ(β,z) by 1−z and suitably shifting the summation index M, we arrive at a generating function for the desired differences Φ(E,M):
equation34
Full-size image (<1 K)
Going back to the representation (32) now reveals that multiplying Ξ(β,z) by 1−z means amputating the ground-state factor ν=0 from this product and retaining only its “excited” part; we therefore denote the result as Ξex(β,z):
equation35
Full-size image (<1 K)
Moreover, from Eq. (34) it is obvious that the canonical moments of the unrestricted set Φ(E,M) (that is, the moments pertaining to allΦ(E,M) with M⩾0) are obtained by repeatedly differentiating Ξex(β,z) with respect to z, and then setting z=1:
equation36
Full-size image (<1 K)
So far, all rearrangements have been exact.

For making contact with the actual N-particle system under consideration, we now have to restrict the summation index M: If the sum over M did not range over all particle numbers from zero to infinity, but rather were restricted to integers not exceeding the actual particle number N, then Eq. (36), together with the representation (35), would yield precisely the non-normalized kth moments of the canonical distribution (28). As it stands, however, exact equality is spoiled by the unrestricted summation. At this point, there is one crucial observation to be made: In the condensate regime the difference between the exact kth moment, given by a restricted sum, and the right-hand side of Eq. (36) must be exceedingly small [16]. Namely, if there is a condensate, then Φ(E,N)/Ω(E,N) is negligible, since the statistical weight of those microstates with the energy E spread over all N particles must be insignificant—if it were not, so that there was a substantial probability for all N particles being excited, there would be no condensate! Consequently, we have Φ(E,M)/Ω(E,N)≈0 also for all M larger than N: In the condensate regime it does not matter whether the upper limit of the sum over M in Eq. (36) is the actual particle number N, or infinity;

equation37
Full-size image (<1 K)
Within this approximation—which will remain the only approximation in the entire argument!—the amputated function Ξex(β,z) provides, by means of Eq. (36), the moments of the excited-states distribution (28), as long as one stays in the condensate regime.

The rationale behind this reasoning can be interpreted in a twofold manner. Intuitively, one may divide a partially condensed Bose gas into the excited-states subsystem, and a supply of condensate particles. The approximation (37) then means replacing the actual condensate, consisting of a finite number of particles, by an infinite particle reservoir [21]; this point of view goes back to Fierz [22]. Of course, the added condensate particles do not take part in the dynamics, so that the statistical properties of the excited subsystem remain unchanged. For our purposes, another interpretation is more telling: For k=0, and utilizing the representation (35) of the moment-generating function Ξex(β,z), the approximation (37) brings Eq. (36) into the form

equation38
Full-size image (<1 K)
—stating that in the condensate regime, where (37) holds, the canonical partition function of the Bose gas on the right-hand side equals that of a system of harmonic oscillators, with frequencies ενε0 (ν⩾1), on the left-hand side of Eq. (38). While the original particles are indistinguishable, and subject to Bose statistics, the substituting oscillators obey Boltzmann statistics; the partition function of the oscillator system being the simple product of the geometric series representing the partition functions of the individual oscillators. It should be noted that this (almost-) isomorphism of a partially condensed, ideal Bose gas and a Boltzmannian harmonic oscillator system holds for any form of the single-particle spectrum, not only for harmonic traps.

Since, according to Eq. (36), the function Ξex(β,z) generates the moments of the canonical distribution (28) in the condensate regime, its logarithm generates the cumulants κk(β):

equation39
Full-size image (<1 K)
As is well known from elementary statistics, the kth order cumulant of the sum of independent stochastic variables is the sum of their individual kth order cumulants; moreover, all cumulants higher than the second vanish exactly for a Gaussian stochastic variable [23]. In the present case, taking the required derivatives of Full-size image (<1 K) and using expression (18) for the canonical occupation numbers 〈nν〉 in the condensate regime, we find for k=1,…,4:
equation40
Full-size image (<1 K)
equation41
Full-size image (<1 K)
equation42
Full-size image (<1 K)
equation43
Full-size image (<1 K)
The fact that κ3(β) and κ4(β) are non-zero signals the non-Gaussian nature of the fluctuations; the fact that all cumulants are given as simple sums over the oscillator index ν reflects the independence of the individual Boltzmannian oscillators.

Because the excited-states distribution (28) is related to the ground-state distribution pn0 through pn0=pN(ex)(Nn0;β), one also has

equation44
Comparison of , ,  and  derived here from orthodox statistical mechanics with the previous master-equation results , ,  and  thus shows that the master equation, solved within the quasithermal approximation (5), actually had yielded the first two cumulants exactly; the effect of the quasithermal approximation only consists in breaking higher-order moments into suitable combinations of first- and second-order ones:
equation45
In this respect, the quasithermal approximation is a kind of mean field approximation.