{'fraction_non_alphanumeric': 0.03532644682823467, 'fraction_numerical': 0.008839888756456098, 'mean_word_length': 4.335629747394453, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Behavior of heat capacity of an attractive Bose-Einstein Condensate approaching collapse

We report calculation of heat capacity of an attractive Bose-Einstein condensate, with the number N of bosons increasing and eventually approaching the critical number Ncr for collapse, using the correlated potential harmonics (CPH) method. Boson pairs interact via the realistic van der Waals potential. It is found that the transition temperature Tc increases initially slowly, then rapidly as N becomes closer to Ncr . The peak value of heat capacity for a fixed N increases slowly with N, for N far away from Ncr . But after reaching a maximum, it starts decreasing when N approaches Ncr . The effective potential calculated by CPH method provides an insight into this strange behavior.


I. INTRODUCTION
Bose-Einstein condensation (BEC) is the transition process, in which a macroscopic fraction of bosons goes into the lowest energy state, as the temperature is lowered below a certain critical temperature T c [1]. It was predicted by Einstein in 1925, based on Bose's explanation of black body radiation. A great deal of activity, both theoretical and experimental, has been seen in this field, since the experimental realization of BEC in 1995. Although a number of static, dynamic and thermodynamic properties have been studied [2,3], not much attention has been paid to the heat capacity of attractive condensates. The main motivation of this work is to fill this gap.
In laboratory experiments, the condensate is trapped by a confining potential, usually a harmonic oscillator potential. An attractive condensate (e.g. 7 Li condensate) has a negative value of the s-wave scattering length a s and collapses, when the number of particles N in the condensate exceeds a critical number N cr . On the other hand a repulsive condensate (e.g. 87 Rb condensate) corresponds to a s > 0 and is stable for any N , since repulsively interacting bosons are contained in the externally applied trap. The situation is quite different for an attractive condensate: attractive bosons tend to come to the center of the trap, which is balanced only by the kinetic pressure, resulting in a metastable condensate. The total attraction increases as the number of pairs N (N − 1)/2, while the kinetic pressure increases as N . Thus for N larger than a critical value N cr , the net attraction dominates and a collapse occurs.
In this communication, we report the calculation of heat capacity of an attractive condensate containing a fixed number of 7 Li atoms, interacting via the realistic van der Waals potential, appropriate for the experimental scattering length. The features are markedly different from those of repulsive condensates, only which have so far been investigated. In a repulsive condensate, the heat capacity C N (T ) for a fixed number N of bosons in the trap, as also the critical temperature T c , smoothly approach a constant value as N → ∞ [4,5]. As a function of T , the heat capacity for a given N increases to a maximum (C N ) max , then falls rapidly to a saturation value 3N k B as T increases (k B is the Boltzmann constant). These features are qualitatively similar to those of a trapped non-interacting condensate [3,5]. However for an attractive condensate, there are important changes in the nature. This is due to the fact that the number of available energy levels of the system is limited, especially when N → N cr , while for any N , there are infinitely many energy levels for the repulsive or non-interacting condensate. In the limit of high T , both the repulsive and the non-interacting condensate behave as the corresponding trapped Bose-gas, resulting in a saturation in C N (T ). On the other hand, an attractive condensate also shows similar behavior, only if it is allowed to absorb energy internally through rotational motion involving large orbital angular momenta.
We provide an understanding of this peculiar nature based on the many-body picture. For the theoretical calculation, we adopt the correlated potential harmonic (CPH) method [6,7] to solve the many-body problem approximately. This technique is based on the potential harmonics (PH) expansion method [8]. The laboratory BEC must be very dilute to preclude three-body collisions, which lead to molecule formation and consequent depletion. Hence only two-body correlations are relevant. The PH is a subset [8] of the full hyperspherical harmonics (HH) basis [9], that involves only two-body correlations. Hence the PH basis is a good approximation for expanding the condensate wave function. It reduces the bulk of the numerical procedure immensely, while retaining the most important basic features of the condensate. However, the leading members of the PH basis do not have the correct short separation behavior of the interacting Faddeev component. This causes a very slow rate of convergence of the PH expansion basis. To correct for this, we include a short-range correlation function in the expansion basis. This correlation function is obtained as the zero-energy solution of the two-body Schrödinger equation [7]. It is a correct representation of the short separation behavior and also incorporates the s-wave scattering length, a s , through its asymptotic behavior [3]. The technique has been shown to reproduce known results, both experimental and theoretical [10]. These include the following: ground state properties (energy, wave function, condensate size, one-body density, pair-correlation, etc.) as also multipolar moments of both repulsive and attractive condensates, correct prediction of the critical number and collapse scenario of attractive condensates, thermodynamic properties of repulsive condensates, properties of condensates in finite traps, etc.
We can understand the behavior of heat capacity of attractive condensates in terms of the energy levels of the system produced by the CPH method. This method generates an effective potential in which the condensate moves. For an attractive condensate, the effective potential has a metastable region (MSR), separated from a deep well on the inner side by an intermediate finite barrier. A finite number of energy levels are supported by the MSR. As the number N of atoms increases, the MSR shrinks and the number of energy levels reduce drastically. As temperature increases, particles are distributed in higher energy levels, according to Bose distribution formula. Thus at low temperatures the internal energy and C N (T ) increase with temperature. At higher temperatures, the bosons have fewer levels to occupy, causing C N (T ) to differ from the repulsive case. There is also a dominant effect of the drastically reducing number of energy levels as N → N cr .
The paper is organized as follows. For easy readability and to introduce our notations, we briefly review the correlated potential harmonic method in Section II. Section III provides our numerical procedure. Results and discussion are presented in Section IV. Finally we draw our conclusions in Section V.

II. CORRELATED POTENTIAL HARMONICS (CPH) METHOD
We adopt the correlated potential harmonics method [6,7] to solve the many-body problem of the BEC. We briefly recapitulate the technique in the following. Interested readers can find details in the cited references.
For the relative motion of a system of N identical spin-less bosons, we introduce (N − 1) Jacobi vectors where x i is the position vector of the i-th particle. The Schrödinger equation governing the relative motion of the system trapped in a harmonic well, is where N = N − 1 and the trapping potential V trap and the interatomic interaction V are expressed in terms of the Jacobi vectors. The energy of the relative motion is E R . Next, we introduce hyperspherical variables corresponding to the set of N Jacobi vectors. First, a hyperradius is defined as The remaining set of (3N − 1) 'hyperangles' consists of 2N polar angles of N Jacobi vectors and (N − 1) angles defining their relative lengths [9]. In the hyperspherical harmonics expansion method (HHEM) ψ is expanded in the complete set of hyperspherical harmonics (HH), which are the eigenfunctions of the grand orbital operator [hyperangular part of the N dimensional Laplace operator, given by the sum in the first term of Eq.
(2) and projection on a particular HH result in a set of coupled differential equations. Imposition of symmetry of the wave function and calculation of the matrix elements become increasingly difficult and tedious as N increases. In addition, the degeneracy of the HH basis increases very rapidly [9] with the increase in the grand orbital quantum number K. Hence a convergent calculation using HHEM with the full HH basis is extremely computer intensive and unmanageable for N > 3. This is the price one pays for keeping all many-body correlations in ψ.
However all these complications can be avoided and a much simpler computational procedure can be formulated for the laboratory BEC, which is designed to be extremely dilute (typical number density is ∼ 10 15 cm −3 ) in order to avoid recombination through three-body collisions. Thus three-body correlations and three-body forces are totally negligible. We can then express ψ as a sum of two-body Faddeev component ψ ij for the (ij)interacting pair [8] Note that the assuption of two-body correlations only makes ψ ij a function of the pair separation vector and the hyperradius only. One can then expand ψ ij in a subset of HH, called the potential harmonics (PH) subset, which is sufficient for the expansion of the interaction potential V ( r ij ) as a function in the hyperangular space for the (ij)-partition. Since the labeling of the particles is arbitrary, we can choose r ij = ζ N . Then the corresponding PH, P lm 2K+l (Ω ij N ) (the argument is the full set of hyperangles for the (ij)-partition) is independent of { ζ 1 , . . . , ζ N −1 } and a simple analytic expression is possible [8]. Expansion of the Faddeev component in the PH basis reads where P lm 2K+l (Ω ij N ) is a potential harmonic [8]. The rdependent factor in front is included to remove the first derivative with respect to r. Substitution of this expansion in the Faddeev equation for the (ij)-partition and projection on the PH corresponding to the (ij)-partition give a set of coupled differential equations in r. Note that any realistic two-body potential, V ( r ij ) can be used. A realistic interatomic potential has a very strong repulsion (arising from the nucleus-nucleus repulsion) at very short separations. Consequently, corresponding ψ ij must be vanishingly small for small values of r ij . But the leading PH (corresponding to K = 0) in the expansion in Eq. (5) is a constant and does not have this behavior. Hence convergence of the expansion in Eq. (5) will be very slow. To improve the rate of convergence, we include a shortrange correlation function η( r ij ) in the expansion basis, so that Eq. (5) is replaced by The short-range correlation function is chosen to have the expected behavior of ψ ij ( r ij , r) for small r ij in the following manner. The small r ij behavior of ψ ij will be that of a zero-energy pair interacting via V ( r ij ), since the energy of the interacting pair is practically zero. We obtain η( r ij ) by solving the zero-energy two-body Schrödinger equation Inclusion of the short-range correlation function η( r ij ) enhances the rate of convergence greatly, which has been checked in our numerical calculation.
The laboratory BEC is very dilute; hence the average separation of the atoms is very large compared with the range of interatomic interactions. Moreover, the atoms scatter with almost zero energy. Hence the effective two-body interaction is represented by the s-wave scattering length a s . In our calculation, we take V ( r ij ) to be the van der Waals potential with a hard core: V ( r ij ) = − C6 rij 6 for r ij ≥ r c and = ∞ for r ij < r c . The correlation function obtained by solving Eq. (8) quickly attains its asymptotic form C(1 − as rij ) for large r ij . The asymptotic normalization is chosen to make the wavefunction positive at large r ij . The hard core radius r c is adjusted so that the calculated a s is the actual experimental value of the scattering length [3]. This procedure assures that the realistic two-body interaction appropriate for the condensate has been incorporated.
Substitution of the expansion, Eq. (7) in Eq. (6) and projection on the PH corresponding to the (ij)-partition result in where 2 , l being the orbital angular momentum contributed by the interacting pair. f 2 Kl is a constant representing the overlap of the PH for interacting partition with the full set of all partitions, which can be found in Ref. [8]. The correlated potential matrix element V KK ′ (r) is given by [7] Here h αβ K and W l (z) are respectively the norm and weight function [11] of the Jacobi polynomial P αβ K (z). Note that the inclusion of the short-range correlation function, η(r ij ) makes the PH basis non-orthogonal. Numerical solution of Eq. (8) shows that η(r ij ) differs from a constant value only in a small interval of small r ij values. Hence the dependence of the overlap N )η(r kl ) > on the hyperradius r is quite small. Disregarding derivatives of this overlap with respect to the hyperradius, we approximately get Eq. (9), with V KK ′ (r) given by Eq. (10). The effect of the overlap being different from unity is represented by the asymptotic constant C of η(r ij ). The emerging physical picture is: the effective interaction between pairs of atoms at very low energy becomes V (r ij )η(r ij ). This is justified, since at very low kinetic energy, the atoms have a very large de Broglie wave length and do not approach each other close enough to "see" the actual interatomic interaction. In the limit of zero energy, the scattering cross section becomes 4π|a s | 2 and the effective interaction is governed by the s-wave scattering length a s , through the asymptotic form of η(r ij ).
Introduction of the PH basis and inclusion of the short-range correlation function, referred to as the correlated potential harmonic (CPH) method, simplifies the many-body problem dramatically. A fairly fast computer code can solve Eq. (9) with upto 15000 particles in the condensate. This technique has been tested against known results, both experimental ones and theoretical ones calculated by other authors, for repulsive as well as attractive condensates [6,7,10].

A. Solution of coupled equations
Although Eq. (9) can be solved by an exact numerical technique using the Numerov method, we adopt the hyperspherical adiabatic approximation (HAA) [12], which apart from simplifying the computations greatly, provides an effective potential in the hyperradial space, in which the condensate moves. This effective condensate potential provides a physical picture for the internal mechanism of the condensate.
In the HAA, one assumes that the hyperangular motion is much faster than the hyperradial motion, since the latter corresponds to the breathing mode. Therefore, one can solve the former adiabatically for a fixed value of r and obtain the solution as an effective potential for the hyperradial motion, as in Born-Oppenheimer approximation. The hyperangular motion is solved by diagonalizing the potential matrix V KK ′ (r) together with the hyper-centrifugal potential [second term of Eq. (9)]. The lowest eigenvalue ω 0 (r) [corresponding eigen column vector being χ K0 (r)], is the effective potential for the hyperradial motion [12]: The third term is an overbinding correction. Eq. (11) is solved by the Runga-Kutta method, subject to appropriate boundary conditions to get E R and the hyperradial wave function ζ 0 (r). The many-body wave function can be constructed in terms of ζ 0 (r) and χ K0 (r) [12]. Total energy is obtained by adding the center of mass energy (1.5 o.u.) to E R . Energy levels, E nl , are characterized by the quantum numbers (n, l), where n represents the excitation quantum number for a given orbital angular momentum l. The HAA has been tested for nuclear, atomic and molecular systems and shown to give better than 1% accuracy, even for the long-range Coulomb potential [13]. In our case, the van der Waals potential has a shorter range and HAA is expected to be better. Moreover, in a BEC, the dominant confining harmonic oscillator potential is smooth and the corresponding hyperradial equation is completely decoupled. Hence in a BEC, the HAA is expected to be far better. We tested this by solving the CDE, Eq. (9), with the interatomic potential for the ground state by the renormalized Numerov method [14,15], which is an exact numerical algorithm for solving a set of coupled differential equations.
where β = 1/k B T and µ is the chemical potential. The latter is determined from the constraint that the total number of particles is N . Clearly, µ has a temperature dependence. The total number of bosons in the trap is fixed and at any temperature it can be written as At a particular temperature T , µ is determined from the constraint Eq. (13). The total energy of the system at T is given by The specific heat at fixed particle number N is calculated using the relation Using (12), (14), (15) one can obtain the heat capacity as For an ideal non-interacting bosonic gas containing N bosons in a three-dimensional isotropic harmonic well, the critical temperature T 0 c is well defined [2]. µ remains equal to the energy of the single particle ground state for T < T 0 c and start decreasing rapidly for T > T 0 c . In the standard text book treatment [1] the sums in eqs. (13)(14) are replaced in the semi-classical approximation by integrals over energy, assuming a continuous energy spectrum. In a harmonic trap the energy spectrum is discrete and this assumption is not valid, particularly for small N at low energies. A correct treatment [5] shows that µ decreases slowly from its maximum value (equal to the ground state energy) as T increases from zero, the rate of decrease becoming suddenly rapid at some temperature close to the reference temperature T 0 c corresponding to the same value of N [4]. Thus, in this case the critical temperature is not well defined. In the correct treatment, the heat capacity also becomes a smooth function of T , attaining a maximum at a temperature, at which µ suddenly becomes a rapidly decreasing function of T . In the limit of large N , the behaviors of µ(T ) and C N (T ) curves approach those of the text book treatment. The transition temperature T c for a finite interacting system is defined as the temperature at which C N (T ) is a maximum [5] ∂C N (T ) ∂T For the numerical calculation for a chosen particle number N , the CPH equations are solved for a large number of energy levels -typically n running from 0 to 300 and l running from 0 to 200, subject to an upper energy cutoff value, E UL (see below), so that E nl < E UL . Using previously calculated values of E nl , Eq. (13) is solved for µ, at a chosen temperature T , by a modified bisection method. Next E UL is increased and the process repeated, until convergence in µ is achieved. Using this upper energy cutoff, C N (T ) is determined using Eq. (16).

IV. RESULTS AND DISCUSSIONS
We consider the attractive condensate of 7 Li atoms in the trap used in the experiment at Rice University [16].
The harmonic trap used was axially symmetric with ν x = ν y = 163 Hz and ν z = 117 Hz. For simplicity, we consider an isotropic trap with ν = (ν x ν y ν z ) 1 3 . The experimental value of a s is −27.3 o.u. We use oscillator unit (o.u.) of length ( mω ) and energy ( ω). As mentioned earlier, we choose van der Waals (vdW) potential for the interatomic interaction, with known value [3] of C 6 = 1.715 × 10 −12 o.u. The value of r c is obtained by the procedure discussed following Eq. (8), so that calculated a s has the experimental value [3]. Its numerical value is 5.338 × 10 −4 o.u. The calculated effective potential, ω 0 (r) is plotted as a function of r in Fig. 1 for N = 1300. In the r → 0 limit, ω 0 (r) becomes N cr , for N = 500, 1300, 1400, 1410, 1420, and 1426 respectively. From Fig. 2, one notices that both the depth and width of MSR decrease as N increases towards N cr . Hence the number of bound energy levels supported by the MSR decreases rapidly with N (see also Fig. 7).
However, the effective potentials shown in Fig. 1 and Fig. 2 are obtained for l = 0. For higher l, the effective potential has a higher IB, arising from the l-dependent terms of the hyper-centrifugal repulsion [see Eq. (9)]. Thus the position of the MSR rises higher in energy as l increases, as can be seen in Fig. 3 for l = 0, 1, 2, 3, 4 for a condensate containing 1420 atoms. Hence particles with l > 0 can attain higher energy levels. Inclusion of these levels will have a profound effect on the heat capacity. If such energy states were ignored (i.e., only the energy levels supported by the l = 0 MSR considered), the heat capacity would reduce drastically and T c would increase indefinitely as N → N cr , since all the atoms would be forced into the few remaining energy levels available for internal excitation. In our calculation, we have retained all energy levels supported by a given l. A question arises as to whether the metastable condensate can have large l values. Intuition indicates that, with increase of temperature, the system can absorb energy only by increasing its rotational kinetic energy, thereby increasing the stability of the metastable system with enhanced centrifugal repulsion. Increase of kinetic energy due to faster linear motion alone would cause the system to fall in the DNW near the center of the condensate. Compared with the non-interacting or repulsive condensates, the attractive condensate has a clear distinction, viz. while the number of hyperradial excitations for a given l in the former is not limited, it is drastically limited in the latter. Thus for an attractive condensate, there are fewer energy states, in which the system can reside. This causes (C N ) max to increase initially for N ≪ N cr (when energy levels are not greatly restricted), but as N → N cr , it starts decreasing, after attaining a maximum. Fig. 4 shows how C N (T ) depends on T , for selected values of N . It is seen that the transition temperature T c increases gradually with N , but (C N ) max increases up to N = 1300 and for larger N , it starts decreasing. When N ≪ N cr , the nature is similar to that of a repulsive condensate [4], since in this case, the number of available energy levels are still large enough (the top most energy level -including l = 0in the MSR has an energy much greater than k B T ), so that the top most levels are still practically unoccupied and there is scope for further internal excitation as T increases. Consequently T c increases gradually with N , as in the repulsive case. As N approaches N cr , the number of energy levels supported by the MSR decreases rapidly and there is less scope for absorbing energy internally as T increases. Hence (C N ) max decreases and T c increases faster, as N increases towards N cr . At higher temperatures, higher l states are excited, which push atoms further outwards, increasing the average interatomic separation. Consequently, the system behaves ultimately as a non-interacting Bose gas. Thus the asymptotic value of C N (T ) becomes 3N k B . In Fig. 4, we plot the dimensionless quantity C N (T )/(N k B ) against T (in All the curves appear to converge to the Bose gas limit. But a closer scrutiny shows that the curves for N = 1350 and 1400 show a slight downward trend. This is due to a limitation in the higher energy cut-off used in our calculation. In levels as N approaches N cr . We demonstrate this in Fig.  7 for the l = 0 energy levels as N increases from 1200 to N cr . The decrease in the number of available energy levels forces a larger fraction of the bosons to be in the ground state as T increases [see Eq. (13)]. This causes a decrease of (C N ) max and an increase of T c , as N → N cr .
A possible scenario of the attractive condensate as its temperature is gradually increased is the following. The standard definition of critical number N cr referred to in the literature, corresponds to the l = 0 condensate at zero temperature. As T is gradually raised, the system absorbs energy by occupying higher available energy levels upto the top of the MSR. However, there is a finite life time of atoms in higher energy levels due to tunneling through the IB into the DNW. Thus there will be a decrease in the number of atoms in the MSR. The rate of loss of atoms will increase with the energy of the level, as also with N (increase of N will lower the IB). The usual definition of heat capacity, C N (T ), is the rate of change of internal energy with respect to T [Eq. (15)], for a fixed number N of atoms in the condensate. This definition is unambiguous for a repulsive or non-interacting condensate, since in these cases there is no loss. However, since the loss is appreciable for an attractive condensate for N close to N cr , or if T is such that k B T is comparable with the highest excitation energy allowed by the IB, this definition demands that atoms be pumped into the condensate at the same rate as the loss rate from the condensate. Such an experimental procedure has not been adopted yet. However, for N ≪ N cr and T ≪ T c , the highest appreciably occupied levels will have negligible tunneling probability into the DNW. Under these conditions, the metastable attractive condensate is fairly long lived and the standard definition of C N (T ) is acceptable. Hence standard experimental techniques can be adopted. Thus our results presented in Figs. 4 -6 are experimentally verifiable in the small N , small T limit. We have presented, for theoretical completeness, results for N close to N cr and for T beyond T c . The question of how the system can absorb energy internally for such values of N and T was already discussed above.

V. CONCLUSIONS
In this work, we report a detailed calculation of the heat capacity C N (T ) of an attractive Bose-Einstein condensate containing N atoms of 7 Li. The correlated potential harmonics method, which is appropriate for the dilute BEC, has been used. The effective potential, in general, supports a large number of energy levels. At T = 0, the lowest energy level accommodates all the bosons. As temperature increases, particles are distributed in higher energy levels, according to Bose distribution formula. Thus the internal energy of the system increases. Heat capacity C N for a fixed number N of particles in the condensate is defined as the temperature derivative of the total internal energy. For a repulsive condensate trapped by an ideal harmonic oscillator, the effective potential has no upper cut off. Hence the energy levels are not limited in energy. Consequently, total internal energy and C N increase as temperature increases upto T c . For T > T c , the ground state occupation becomes suddenly microscopic (negligible) and the system behaves like a harmonically trapped Bose gas. Hence C N decreases rapidly above T c , reaching its asymptotic value 3N k B . Thus C N first increases, reaches a maximum value (C N ) max and then decreases rapidly to its asymptotic value, as T increases from zero.
For N < N cr bosons with mutual attraction, a metastable condensate is formed in the metastable region (MSR) of the effective potential. On the left of the MSR, an intermediate barrier (IB), followed by a deep narrow well and finally a strongly repulsive core appear, as one approaches the center of the condensate. For N ≪ N cr , the IB is very high and the minimum of the MSR is very low, so that the metastable well is sufficiently deep compared with thermal excitation energy k B T c at T c , and a large number of energy levels are supported. Hence for T ≤ T c , even the most thermally excited particles do not feel the effect of the IB and C N increases gradually, as in the repulsive case.
With increase of temperature, the system with a fixed N absorbs energy internally by increasing the occupation probability of higher energy levels supported by the metastable region of the effective potential. Atoms in energy levels close to the top of the intermediate barrier have appreciable probability to tunnel into the deep narrow well, causing the condensate to loose atoms. But such levels are not occupied with any appreciable probability if N ≪ N cr and T ≪ T c . Hence such a condensate is quasi-stable and C N (T ) calculated for a fixed N is appropriate. When the rate of atom loss from the condensate is appreciable, standard definition of heat capacity at constant N requires feeding the attractive condensate with additional atoms at a rate such as to compensate for the loss rate. Although this is not the usual experimental technique, we investigate such cases for a complete theoretical study. In such a situation, there are drastic changes. The peak value of C N (T ) (the temperature at which this occurs is the transition temperature T c ) initially increases gradually with N , then after reaching a maximum, decreases fairly rapidly near N cr . On the other hand, for small N , T c increases almost linearly up to N ∼ 1300. For larger N , the transition temperature increases rapidly with N . We provide an explanation of this behavior, based on the microscopic mechanism of absorption of internal energy, as T increases. As N increases towards N cr , the depth and width of the metastable well decrease rapidly. As a result, the number of energy levels supported by the metastable well decreases rapidly. This tends to increase T c , since fewer energy levels are available for absorption of internal energy, and bosons are forced to be in lower energy levels as T increases. Rapid reduction of available energy levels as N → N cr , causes quicker saturation of internal energy of the condensate. Consequently, the maximum of C N (T ) decreases rapidly as N → N cr .

VI. ACKNOWLEDGEMENT
We would like to thank Dr. Parongama Sen for useful discussions. SG acknowledges CSIR (India) for a Ju-