Specific heat and bimodality in canonical and grand canonical
versions of the thermodynamic model
G. Chaudhuri ∗† and S. Das Gupta‡
Physics Department, McGill University, Montréal, Canada H3A 2T8
(Dated: October 22, 2018)
Abstract
We address two issues in the thermodynamic model for nuclear disassembly. Surprisingly large
differences in results for specific heat were seen in predictions from the canonical and grand canon-
ical ensembles when the nuclear system passes from liquid-gas co-existence to the pure gas phase.
We are able to pinpoint and understand the reasons for such and other discrepancies when they
appear. There is a subtle but important difference in the physics addressed in the two models.
In particular if we reformulate the parameters in the canonical model to better approximate the
physics addressed in the grand canonical model, calculations for observables converge. Next we
turn to the issue of bimodality in the probability distribution of the largest fragment in both
canonical and grand canonical ensembles. We demonstrate that this distribution is very closely
related to average multiplicities. The relationship of the bimodal distribution to phase transition
is discussed.
PACS numbers: 25.70Mn, 25.70Pq
∗ On leave from Variable Energy Cyclotron Center, 1/AF Bidhan Nagar, Kolkata 700064, India
†Electronic address: gargi@physics.mcgill.ca
‡Electronic address: dasgupta@physics.mcgill.ca
http://arxiv.org/abs/0704.0288v3
mailto:gargi@physics.mcgill.ca
mailto:dasgupta@physics.mcgill.ca
I. INTRODUCTION
In models of statistical disassembly of a nuclear system formed by the collision of two
heavy ions at intermediate energy one assumes that because of multiple nucleon-nucleon
collisions a statistical equilibrium is reached. The temperature rises. The system expands
from normal density and composites are formed on the way to disassembly. As the system
reaches between three to six times the normal volume, the interactions between composites
become unimportant (except for the long range Coulomb interaction) and one can do a
statistical equilibrium calculation to obtain the yields of composites at a volume called the
freeze-out volume. The partitioning into available channels can be solved in the canonical
ensemble where the number of particles in the nuclear system is finite (as it would be in
experiments). In some experiments, the number of particles can fluctuate around a mean
value. In such a case a sum of several canonical calculations could be appropriate. Even when
the number of particles is fixed one can hope to replace a canonical model calculation by
a grand canonical model calculation where the particle number fluctuates but the average
number can be constrained to a given value. The case we will look at corresponds to
this situation. Usually the grand canonical model is more easily solved. Hence it is more
commonly used although in the case of nuclear physics (where particle numbers are typically
≈ 200 or less) the use of the canonical ensemble would be more appropriate.
Apart from ease of calculation, there is another reason why the grand canonical model is
very useful. Known properties of nuclear interactions predict that if nuclear systems were
arbitrarily large (consider a fictitious system where the Coulomb interaction is switched off)
disassembly of nuclear systems would show features of liquid-gas phase transition [1]. Since
the grand canonical ensemble is expected to become accurate for large systems, this would
seem to be a suitable framework to describe bulk properties. The canonical model, built with
a constant particle number in mind, can not be pushed to arbitrarily large particle number
although it can be implemented for fairly big systems containing thousands of particles. One
can then extraprolate from finite particle number systems to the infinite particle number case.
When this was done, a huge difference showed up between grand canonical and canonical
results for cv (specific heat per particle at constant volume). In the grand canonical model
cv was merely discontinuous at phase transition [2] but in the canonical model cv would go
to infinity as the system particle number approached infinity [3].
This discrepancy was examined in detail for a system with 2000 particles [4]. The analysis
showed that in the co-existence region, even when the average number of particles is 2000,
fluctuations in the number of particles is huge if one uses the grand canonical ensemble.
The cv with a fixed number of particles is sharply peaked at a temperature T where T is
a function of the number of particles. Because in the grand canonical ensmble the particle
number fluctuations are significant when the average number is 2000, the resulting cv is
smeared out. What this analysis did not answer are two significant questions: (a) under
what conditions are the grand canonical results valid and (b) under what conditions can the
canonical and grand canonical results agree? We answer both these questions in this work.
The apparent paradox is explained.
We also address the issue of bimodality in the probability distribution of the largest
fragment as a function of the mass number of the largest fragment for a finite system when
the thermodynamic model is used. Clearly the canonical model is appropriate here but
we also study bimodality in the grand canonical ensemble, where the average value of the
number of particles is constrained to a typical value expected in heavy ion collisions. In
the thermodynamic model it is easy to devise systems which may (as in the nuclear matter
case) or may not (i.e., by switching off surface tension term in binding energy formula) have
a phase transition (some other models may not have this versatility). Thus in this model
the connection of bimodality to phase transition can be directly established. It is interesting
to note (as described in detail later) that both canonical model and the approximation of a
finite system with the grand canonical model show features of bimodality but quantitatively
the results are quite different.
The plan of the paper is as follows. In section II we set up the the formulae for the grand
canonical ensemble. In section III the methodology of the canonical ensemble is presented.
We point out that, although not obvious, the two ensembles actually addressed different
physics causing the difference in results for specific heat. In section IV we reformulate the
canonical model to better approximate the grand canonical model. We show that results
for cv progressively converge. We then turn to the question of bimodality in the probability
distribution of the largest fragment as a function of the fragment mass number. Formulae
for the probability distribution are given in section V. In section VI we present results
for this distribution and discuss the bimodality in both the canonical and grand canonical
ensembles. The connection between the probability distribution of the largest fragment and
average multiplicity is established in section VII. Summary is presented in section VIII.
As in [1, 2, 4] we use one kind of particle and no Coulomb interaction. This is adequate
for the purpose of this study and offers considerable numerical simplifications. Numerous
applications of the canonical [5] and the grand canonical models [6] with two kinds of particles
exist where fits experimental data are the main issues.
II. FORMULAE IN THE GRAND CANONICAL MODEL
If we have na particles of type a, nb particles of type b, nc particles of type c etc. all
enclosed in a volume V and interactions between particles can be neglected, the grand
partition function for this case can be written as
Zgr =
i=a,b,c...
(1 + eβµiωi + e
.......) =
i=a,b,c,..
exp(eβµiωi) (1)
Here the µi is the chemical potential and ωi the canonical partition function of one particle
of type i. The average number of particles of type i is given by ∂(lnZgr)/∂(βµi) :
ni = e
βµiωi (2)
It is possible that one of the species can be built from two other species. In reverse, a heavier
species can also break up into two lighter species. If α number of particles of type a can
combine with β number of particles of type b to produce γ number of particles of type c,
then chemical equilibrium implies [7] that the chemical potentials of a, b and c are related
by αµa + βµb = γµc.
In our model we have N nucleons in a volume V (which is significantly larger than the
normal nuclear volume) but these nucleons can be singles or form bound dimers, trimers etc.
Chemical equilibrium implies that a composite with k bound nucleons has a chemical poten-
tial kµ where µ is the chemical potential of the monomer (nucleon). Thus our ensemble has
monomers, dimers, trimers etc. upto some species with kmax bound nucleons. In the actual
world of nuclear physics kmax terminates around 250 because of Coulomb interaction but in
the model pursued here we may terminate it arbitrarily at 1 (monomers only), 2(monomers
and dimers), 3 or any large kmax. It was demonstrated in [8] that liquid-gas type phase
transition occurs for large kmax > 2000.
The total number of nucleons will be denoted by N . Of course, the grand canonical
ensemble works best when N is very large, ideally infinite.
We now look into ωi, the partition function of one composite of i nucleons.This factors into
two parts, a traditional translation energy part and an intrinsic part: ωi = zi(tran)zi(int)
where
zi(tran) =
exp(−βp2/2mi)d
(2πmiT )
3/2 (3)
The intrinsic part zi(int) of course contains the key to phase transition. If we regard each
composite to exist only in a ground state with energy e
i , then zi(int) = exp(−βe
i ). We
use e
i = −iW +σi
2/3 where nuclear physics sets W=16 MeV and σ = 18 MeV. This simple
model itself will lead to the main results of this paper. Because of the surface term, energy
per particle drops as i grows. Let us denote by F the free energy of the N nucleons where N
is the total number of nucleons; E be the energy and S, the entropy: F = E−TS. At finite
temperature F will go to its minimum value. The key issue is how the system of N nucleons
breaks up into clusters of different sizes as the temperature changes. At low temperature
E and hence F minimises by forming very large clusters (liquid). But as the temperature
increases S will increase by forming larger number of clusters thus breaking up the big
clusters. Gaseous phase will appear. How exactly this will happen requires calculation and
these show that the system goes through a first order liquid-gas phase transition [2, 8]. As
is the common practice, we used here a slightly more sophisticated model for zi(int). We
make the surface tension temperature dependent in conformity with usual parametrisation
[9]; σ(T ) = σ0[(T
c − T
2)/(T 2c + T
2)]5/4. Here σ0 =18 MeV and Tc=18 MeV. At T = Tc
surface tension vanishes and we have a fluid only. For us this is unimportant as our focus
will be the temperature range 3 to 8 MeV. Also in zi we include not only the ground state
but also the excited states of the composite in the Fermi-gas approximation [1, 9]. The
expression for zi(int) is now complete and easily tractable.
Let us now summarise the relevant equations. For k = 1 (the nucleon which has no
excited states)
(2πmT )3/2 exp(µ/T ) (4)
and for 1 < k ≤ kmax
(2πmT )3/2k3/2 exp[(µk +Wk + kT 2/ǫ0 − σ(T )k
2/3)/T ] (5)
Here nk is the average number of composites with k nucleons. In the rest of the paper, for
brevity, we will omit the qualifier “average”.
A useful quantity is the multiplicity defined as
kmax∑
nk (6)
The number of nucleons bound in a composite with k nucleons is knk and obviously N =
∑kmax
k=1 knk. The pressure is given by
kmax∑
T (7)
This follows from the identity pV = T lnZgr.
Quantities like N, V, nk are all extensive variables. These equations can all be cast in
terms of intensive variables like N/V = ρ, nk/N etc so that we can assume both N and V
approach very large values and fluctuations in the number of particles can be ignored. Thus
for a given temperature and density we solve for µ using
(2πmT )3/2
(exp(µ/T ) +
kmax∑
k5/2 exp[(µk +Wk + kT 2/ǫ0 − σ(T )k
2/3)/T ]) (8)
The sum rule N =
∑kmax
k=1 knk changes to 1 =
knk/N . The energy per particle is given by
kmax∑
Ek (9)
where Ek =
T for k=1 and for 1 < k ≤ kmax
T + k(−W +
) + σ(T )k2/3 − T [∂σ(T )/∂T ]k2/3 (10)
The term T [∂σ(T )/∂T ]k2/3 arises from the temperature dependence of the surface tension.
The effect of this term is small. Eqs. (9) and (10) follow from the identity E = µN− ∂
lnZgr.
From what we have described so far it would appear that V in eqs.(3) to (9) is the
freeze-out volume V , the volume to which the system has expanded. Actually if the freeze-
out volume is V then in these equations we use Ṽ which is close to V but less. The
reason for this is the following. To a good approximation a composite of k nucleons is an
incompressible sphere with volume k/ρ0 where the value of ρ0 is ≃ 0.16 fm
−3. The volume
available for translational motion (eq.(3)) is then Ṽ = V − Vexcluded where we approximate
Vexcluded ≃ N/ρ0 = V0 the normal volume of a nucleus with N nucleons. Similar corrections
are implicit in Van der Waals equation of state. This is meant to take care of hard sphere
interactions between different particles. This answer is approximate. The correct answer is
multiplicity dependent. The approximation of non-interacting composites in a volume gets
to be worse as the volume decreases. We restrict our calculation to volumes V greater than
2V0. This is how the calculations reported proceed. We choose a value of V0/V = ρ/ρ0 from
which V0/Ṽ = ρ̃/ρ0 = ρ/(ρ0 − ρ) is deduced. This value of ρ̃ is used in eq.(8) to calculate µ
and all other quantities. We plot results as function of ρ/ρ0. If we plotted them as function
of ρ̃/ρ0 the plot would shift to the right.
III. THE CANONICAL MODEL SOLUTION
The statistical equilibrium model as described above can be solved for a given fixed
number of particles when the number of particles N is finite. No spread in the number of
particles, which is inherent in the grand canonical ensemble, needs to be made. Extensive
use of the canonical model has been made to fit experimental data [5] so just an outline will
be presented for completeness. Among other applications, the canonical model can be used
to study finite particle number effects on phase transition characteristics.
Consider again N identical particles in an enclosure V and temperature T . These N
nucleons will combine into monomers, dimers, trimers etc. The partition function of the
system in the canonical ensemble can be written as
Here ωi is the one particle partition function of a composite which has i nucleons. We already
encountered ωi in section II: ωi = zi(tran)zi(int) with zi(tran) and zi(int) given in detail.
Other forms for ωi can be used in the method outlined here. The summation in eq.(11) is
over all partitions which satisfy N =
ini. The summation is non-trivial as the number of
partitions which satisfy the sum is enormous. We can define a given allowed partition to be
a channel. The probablity of the occurrence of a given channel P (~n) ≡ P (n1, n2, n3....) is
P (~n) =
∏ (ωi)ni
. (12)
The average number of composites of i nucleons is easily seen from the above equation to
ni = ωi
Since
ini = N , one readily arrives at a recursion relation [10]
kωkQN−k (14)
For one kind of particle, QN above is easily evaluated on a computer for N as large as
3000 in matter of seconds. It is this recursion relation that makes the computation so easy
in the model. Of course, once one has the partition function all relevant thermodynamic
quantities can be computed. For example, eq. (7) still gives the expression for pressure
although one could correct for the center of mass motion by reducing the multiplicity by 1:
p = T (M − 1)/Ṽ . The chemical potential can be calculated from µ = F (N) − F (N − 1)
where the free energy is F (N) = −T ln QN which is readily available from the calculation.
IV. GENERATING GRAND CANONICAL RESULTS FROM THE CANONICAL
ENSEMBLE
We first consider pressure (eq.(7)) in the grand canonical ensemble. The V in eq.(7)
cancels out the V in eqs. (4)and (5) and thus pressure is given in terms of intensive variables
directly. We may assume that this is truly the pressure in infinite systems (V and nk
arbitrarily large in which case fluctuations in the grand canonical ensemble can be ignored).
However the grand canonical answer does depend upon the value of kmax. In Fig.1 we have
used kmax=2000, a value large enough so that liquid-gas transition type features emerge
(the flatness of pressure against ρ). For kmax significantly lower, the flatness disappears
(see details in [8]). In Fig.1 we also show several canonical model results all with the same
kmax=2000 but different values of N . For N=2000, the canonical results are quite different
from the grand canonical results except for very low densities. In particular a region of
mechanical instability is seen which can give rise to a region of negative cp, the specific heat
per particle at constant pressure (see [5] for detailed discussion). In the same figure, we
have also shown prssures in the canonical model when N=100,000 and 500,000. We see that
the pressure approaches the grand canonical value as N increases (the periodicity obvious
in the curve for N= 100,000 arises from the fact the largest composite has k =2000 but we
will not get into a detailed analysis here). The conclusion here is that the grand canonical
value of pressure in Fig.1 refers to a system which has N = ∞ where the largest cluster has
k = 2000. This is quite different from the usual canonical model result which would have
N = kmax = 2000. To address the physics of the grand canonical model we keep kmax still
at 2000 but need to keep on increasing the value of N . Then the canonical results converge
towards the grand canonical values.
For this given problem we have approached the grand canonical result from a canonical
ensemble. One can consider the reverse problem: getting the canonical result starting from
the grand canonical model. It is of course obvious that the correspondence would be exact
provided one uses the appropriate V in the grand canonical ensemble and then projects
from it the part which has an exact N . This is because the grand canonical ensemble is a
particular weighted sum of canonical ensembles with different N ’s.
Let us now turn to Fig.2 which deals with cv, the specific heat per particle at constant
volume. We again keep kmax = 2000. One finds that the cv in a canonical calculation for
N = kmax = 2000 produces a very sharp peak. The grand canonical expression for energy
per particle (eq.(9)) is an intensive quanity and so is its derivative cv. We expect this grand
canonical result for cv is valid for N very large. Comparison shows that the grand canonical
result differs drastically from the canonical N=2000 result in a very narrow window when
the system passes from the co-existence to a pure gas phase. Can we make the results
converge by successively increasing the value of N in the canonical model? The answer is
“yes” as Fig.2 demonstrates. We see that the canonical result with kmax=2000 approaches
the grand canonical result with kmax=2000 as the number of particles N in the canonical
calculations is progressively increased beyond N=2000.
To summarise: the grand canonical model is applicable when N >> kmax. The limit
N >> kmax; kmax → ∞ is robust (as shown in [8]) and produces a first order phase transition.
This model is distinct from the canonical model N = kmax; kmax → ∞. There is no scaling
in this latter model:N = kmax, both N and kmax very large is not equivalent within a factor
of scaling to a system with 2N = 2kmax. We are unable to provide a robust limit for the
canonical model of N = kmax;N → ∞. Fig. 3 shows the progression of the E/N and
pressure curves as N = kmax increases from 2000 to 50,000.
V. BIMODALITY: THE BASIC FORMULAE IN GRAND CANONICAL AND
CANONICAL ENSEMBLE
In event by event analysis in experiments, one can in principle ascertain the largest mass
(or the largest charge) emerging in each event from multifragmentation. The probability
distribution of this largest mass can be plotted as a function of the value of mass of this
largest fragment. It is shown that a bimodality in this distribution at a certain temperature
is a signature of a first order phase transition: that is, if the system were infinitely large
it would have a first order phase transition [12, 13]. Thus from a finite system one can
have a signal for phase transition. We will now see how the probability distribution of the
largest fragment as a function of the mass of the largest fragment can be computed in the
thermodynamic model in the two ensembles. First the grand canonical ensemble.
The grand canonical ensemble works best for a large system and we have already seen in
the previous section that application of this model to multifragmentation of finite nuclei can
lead to serious errors in some temperature (equivalently energy) window. Nonetheless, let
us proceed to see how results can be derived. We fix a value for ρ/ρ0 (in figs. (4) and (5) we
have kept this at 0.25) and choose the appropriate value of the volume so that the average
number N =
∑kmax
1 knk is 150, the system whose results we show. The heaviest composite
allowed in the model kmax is also 150.
From eqs. (1) and (2) one can derive that the probability that a particular composite
with k nucleons does not occur at all is
and the probability that it occurs at least once or more is
enk−1
Note that our nk here is the same as nk’s of eqs. (4) and (5), the average values in the
grand canonical ensemble. The probability that k is the highest mass fragment in an event
is then given by (k < kmax)
Pm(k) =
enk − 1
e−(nk+1+nk+2+.....nkmax) (15)
From the above eq. one readily derives
Pm(k + 1)
Pm(k)
enk+1 − 1
enk − 1
enk =
enk+1 − 1
1− e−nk
If nk+1 > nk then Pm(k + 1) > Pm(k). If further both nk+1 and nk are small compared to 1
Pm(k + 1)
Pm(k)
Let us turn to the calculation of the probability disribution of the largest fragment as a
function of the mass of the largest fragment in the canonical model. A detailed formulation
when two kinds of particles are present was given in a recent paper [11] but for competeness,
we review the development.There is an enormous number of channels in Eq.(11). Different
channels will have different values for the largest fragment. For example there is a term
in the sum of Eq.(11). In this channel all the fragments and hence also the largest fragment
has mass 1. The probability of this channel occurring is (from Eq.(12)) Pm(1) =
The full partition function can be written as QN = Q(ω1, ω2, ω3, .......ωkmax). If we construct
a QN where we set all ω’s except ω1 to be zero then this QN(ω1, 0, 0, 0..............) =
and this has the largest mass 1. Consider now constructing a QN with only two ω’s :
QN(ω1, ω2, 0, 0, 0, ......). This will have the largest mass sometimes 1 (as
is still there)
and sometimes 2 (as, for example, in the term
(N−6)!
). It then follows that
Pm(k) =
QN (ω1, ω2, ...ωk, 0, 0, 0..)−QN (ω1, ω2, ...ωk−1, 0, 0, 0...)
In the above the first term in the numerator takes care of the occurrence of all partitions
where the largest fragment is between 1 and k and the second term takes care of all the
partitions where the largest fragment is between 1 and k−1. The difference, divided by QN
is the desired answer.
Since one has the general formula for the probability Pm(k), one can compute the average
value of the mass of the largest fragment as well as the root mean square deviation. In fact,
these have been measured in some experiments [14] and have recently been calculated [11].
But we will not need this for this paper.
VI. REPRESENTATIVE RESULTS
The probability distribution Pm(k) of the largest fragment as a function of k where k
is the largest fragment in an event is shown in Fig.(4) where the freeze-out density ρ/ρ0 is
0.25 and the dissociating system has N=150 (for the grand canonical ensemble the average
value is 150). The canonical and grand canonical results are quite different but both display
bimodality (there are two maxima with similar heights), the grand canonical at temperature
≈ 5.9 MeV and the canonical at temperature ≈ 6.2 MeV. In Fig.(5) we have compared the
nk’s of the two models. Near the end value 150 the differences are very substantial at all
temparatures. At lower values of k they agree very well at T = 6.8 MeV, quite well at T=6.2
MeV but gets worse at lower temperatures becoming quite different at T=5.0 MeV. These
differences have been noted and discussed before [15].
VII. CONNECTION BETWEEN PROBABILITY DISTRIBUTION OF THE
LARGEST FRAGMENT AND AVERAGE MULTIPLICITY
The very first experiments in heavy ion collisions measured nk, the average multiplicity
against k. One of the earliest postulates was the following. At low energy nk first falls with
k but after reaching a minimum rises again. This is the so-called “U” shape. This shape
at lower temperature is an indication that the system will undergo a liquid-gas type phase
transition. As the energy of collision increases, the height of the maximum on the heavier
side will decrease, will then disappear (this marks the phase transition temperature). At
higher energy, nk decreases monotonically with k. This is discussed in many places including
[1, 3]. Basically then one looks at the behaviour of nk as a function of k and energy as one
signature of phase transition. Since bimodality in the probability distribution is also a
signature of phase transition, we hope to get a connection between Pm(k) and nk.
For bimodality one requires that after the minimum following the first maximum, Pm(k)
will rise again witk k. Similarly in conjectures involving the multiplicity, nk, after reaching
a minimum must rise again with k. These two features are intimately related. In the grand
canonical model this is very simple to prove. Equation (16) shows that if nk+1 > nk then
Pm(k+1) > Pm(k) and bimodality can happen. The reverse is not true; nk+1 < nk does not
imply that Pm(k + 1) is less than Pm(k).
There is similar connection in the canonical model. Here it can be proven that on the
heavier side k > N/2, a rise of nk with k guarantees that Pm(k) will rise with k. In fact it
is even more direct than that. For k > N/2, we have an equality; Pm(k) = nk. This can
be proven from eq. (18) but there is a an easier proof. We can rewrite Pm(k) as a sum of
terms:
Pm(k) = P
m(k) + P
m(k) + P
m(k) + ........ (19)
where in each of the terms in the right hand side k is the highest mass that occurs but in
P 1m(k) the composite k occurs only once, in P
m(k) it occurs twice, in P
m(k) it occurs three
times and so on. Specifically,
P 1m(k) =
; N − k =
i=1 ini
P 2m(k) =
; N − 2k =
i ini
It is clear what the structures for higher terms in the series will be. It is then also obvious
nk = P
m(k) + 2P
m(k) + 3P
m(k) + .......+O.C. (20)
In the above, O.C. stands for other channels where mass k occurs but it is not the highest
mass in the channel. If k > N/2 then only P 1m(k) exists. Thus for k > N/2 we have
Pm(k) = nk. If nk rises with k in this region then so does Pm(k). For k ≤ N/2, the
relationship is nk ≥ Pm(k), with nk usually significantly larger than Pm(k). For bimodality
to exist we need to have nk rising with k in some region k > N/2.
The relation Pm(k) = nk for k > N/2 is not limited to the thermodynamic model only.
It is true in any number conserving model.
VIII. SUMMARY
This paper had two goals. One, to resolve and understand the difference between grand
canonical and canonical values of specific heat and pressure in thermodynamic models as
applied to heavy ion collisions. This issue we believe is resolved. Second, to understand
the link between bimodality in the distribution of the heaviest fragment and the average
multiplicity of fragments (which has also been linked with aspects of phase transition).
We think we have gained an understanding. In a later publication we expect to show
more results for bimodal distributions for realistic cases with two kinds of particles and the
Coulomb interactions included. Calculations for a particular case have already appeared
[11].
IX. ACKNOWLEDGEMENT
This work is supported by the Natural Sciences and Engineering Research Council of
Canada.
[1] S. Das Gupta, A. Z. Mekjian and M. B. Tsang, Advances in Nuclear Physics, vol. 26, 91,
(2001).
[2] K. A. Bugaev, M. I. Gorenstein, I. N. Mishustin and W. Greiner, Phys. Rev. C62, 044320
(2000).
[3] S. Das Gupta and A. Z. Mekjian, Phys. Rev. C57, 1361 (1998).
[4] C. B. Das, S. Das Gupta and A. Z. Mekjian, Phys. Rev. C68 , 031601(R) (2003)
[5] C. B. Das, S. Das Gupta, W. G. Lynch, A. Z. Mekjian and M. B. Tsang, Phys. Rep. 406, 1
(2005)
[6] S. Das Gupta and A. Z. Mekjian, Phys. Rep. 72, 131 (1981)
[7] F. Reif, Fundamentals of statistical and thermal physics,(McGraw Hill, New York, 1965) Chap-
ter 8.
[8] G. Chaudhuri, S. Das Gupta and M. Sutton, Phys. Rev. B74, 174106 (2006).
[9] J. P. Bondorf, A. S. Botvina, A. S. Iljinov, I. N. Mishustin and K. Sneppen, Phys. Rep. 257,
133 (1995).
[10] K. C. Chase, and A. Z. Mekjian, Phys. Rev C52, R2339, (1995)
[11] G Chaudhuri and S. Das Gupta, Phys. Rev. C75, 034603 (2007)
[12] F. Gulminelli and Ph. Chomaz, Phys. Rev. C 71, 054607 (2005).
[13] M. Pichon et al., Nucl. Phys. A779, 267 (2006).
[14] J. B. Elliott et al., Phys Rev C67,024609 (2003)
[15] C. B. Das, S. Das Gupta, and B. K. Jennings, Phys. Rev. C70,044611(2004)
0 0.1 0.2 0.3 0.4 0.5
 ρ /ρ0
0.003
0.004
0.005
0.006
0.007
0.008
0.009
0.011
T= 6 MeV k max=2000
FIG. 1: Pressures calculated in the canonical model compared with pressure calculated in the grand
canonical model. For all of these the largest composite allowed has 2000 nucleons (kmax = 2000)
and the temperature is 6 MeV. The grand canonical calculation is the solid curve. The canonical
calculations are done with N=2000 (dash-dot), N=100,000 (dots) and N=500,000(dash). As N
increases, agreement with the grand canonical result becomes better and better.
6 6.5 7 7.5 8
Temperature (MeV)
6 6.5 7 7.5 8
Temperature (MeV)
 ρ/ ρ0 = 0.25
 max 
=2000
 ρ/ρ0=0.25 k max=2000
FIG. 2: Specific heat at constant volume in grand canonical and canonical models. Two different
scales are needed to highlight differences in values. Again the canonical calculations are done
with N=2000 (dash-dot), N=100,000 (dots) and N=500,000 (dash). The difference between the
grand canonical result (solid) and the N=2,000 calculation is huge around 7 MeV (upper panel)
but for N=100,000 and N=500,000 the results are so close to the grand canonical values that
they are nearly indistinguishable in the scale of the upper panel. In the lower panel results for
N=100,000 and N=500,000 are compared with grand canonical values. Even in this expanded
scale the N=500,000 canonical results are indistinguishable in the curve from the grand canonical
results.
0 0.1 0.2 0.3 0.4 0.5
0.005
0.0075
6.5 7 7.5 8
Temperature (MeV)
T = 6 MeV
 ρ/ρ0 =0.25
FIG. 3: Pressure and energy per particle for the canonical model of N = kmax for N=2000 (solid),
N=10,000 (dash) and N=50,000 (dot).
0 50 100 150
0 50 100 150
0 50 100 150
0 50 100 150
T= 5.0 MeV
T=6.2 MeV T=6.8 MeV
T= 5.9 MeV
FIG. 4: Probability that the largest cluster has k nucleons plotted as a function of k in the grand
canonical(solid) and the canonical model(dot). Here N=150 in the canonical model and in the
grand canonical model the average value is set at N=150. The value of kmax is also 150. The
density is fixed at ρ/ρ0 = 0.25. In the grand canonical model bimodality is seen at about 5.9 MeV
and in the canonical model this appears at about 6.2 MeV.
1 10 100
1 10 100
1 10 100
1 10 100
T= 5.0 MeV
T=6.2 MeV T=6.8 MeV
T= 5.9 MeV
FIG. 5: For the same cases as above, the average multiplicity of each composite plotted as a
function of mass number k. Note that the grand canonical results (solid) approximate the canonical
results (dots) quite well at the highest temperature (except at very high mass numbers) but the
agreement worsens as the temperature is lowered.
	Introduction
	Formulae in the grand canonical model
	The canonical model solution
	Generating grand canonical results from the canonical ensemble
	Bimodality: the basic formulae in grand canonical and canonical ensemble
	Representative Results
	Connection between probability distribution of the largest fragment and average multiplicity
	Summary
	Acknowledgement
	References

  We address two issues in the thermodynamic model for nuclear disassembly.
Surprisingly large differences in results for specific heat were seen in
predictions from the canonical and grand canonical ensembles when the nuclear
system passes from liquid-gas co-existence to the pure gas phase. We are able
to pinpoint and understand the reasons for such and other discrepancies when
they appear. There is a subtle but important difference in the physics
addressed in the two models. In particular if we reformulate the parameters in
the canonical model to better approximate the physics addressed in the grand
canonical model, calculations for observables converge. Next we turn to the
issue of bimodality in the probability distribution of the largest fragment in
both canonical and grand canonical ensembles. We demonstrate that this
distribution is very closely related to average multiplicities. The
relationship of the bimodal distribution to phase transition is discussed.

Introduction
	Formulae in the grand canonical model
	The canonical model solution
	Generating grand canonical results from the canonical ensemble
	Bimodality: the basic formulae in grand canonical and canonical ensemble
	Representative Results
	Connection between probability distribution of the largest fragment and average multiplicity
	Summary
	Acknowledgement
	References
