Density dependence of the symmetry energy and the nuclear equation of
  state: A Dynamical and Statistical model perspective
Density dependence of the symmetry energy and the nuclear equation of state: A
Dynamical and Statistical model perspective
D.V. Shetty, S.J. Yennello, and G.A. Souliotis
Cyclotron Institute, Texas A&M University, College Station, TX 77843, USA
(Dated: October 25, 2018)
The density dependence of the symmetry energy in the equation of state of isospin asymmetric
nuclear matter is of significant importance for studying the structure of systems as diverse as the
neutron-rich nuclei and the neutron stars. A number of reactions using the dynamical and the
statistical models of multifragmentation, and the experimental isoscaling observable, is studied to
extract information on the density dependence of the symmetry energy. It is observed that the
dynamical and the statistical model calculations give consistent results assuming the sequential
decay effect in dynamical model to be small. A comparison with several other independent studies
is also made to obtain important constraint on the form of the density dependence of the symmetry
energy. The comparison rules out an extremely “ stiff ” and “ soft ” form of the density dependence
of the symmetry energy with important implications for astrophysical and nuclear physics studies.
PACS numbers: 21.30.Fe, 25.70.-z, 25.70.Lm, 25.70.Mn, 25.70.Pq
I. INTRODUCTION
The fundamental goal of nuclear physics is to under-
stand the basic building blocks of nature - neutrons and
protons - and the nature of interaction that binds them
together into nuclear matter. Studying the nature of
matter and the strength of nuclear interaction is key to
understanding some of the fundamental problems such
as, How are elements formed? How do stars explode into
supernova? What kind of matter exists inside a neu-
tron star? How are neutrons compressed inside a neutron
star to density trillions of times greater than on earth ?
What determines the density-pressure relation, the so-
called equation of state?
The key ingredient for constructing the nuclear equa-
tion of state is the basic nucleon-nucleon interaction. Un-
til now our understanding of the nucleon-nucleon inter-
action has come from studying nuclear matter that is
symmetric in isospin (neutron-to-proton ratio, N/Z ≈ 1)
and matter found near normal nuclear density (ρo ≈ 0.16
fm−3). It is not known how far this understanding re-
mains valid as one goes away from the normal nuclear
density and symmetric nuclear matter. Various interac-
tions used in “ ab initio ” microscopic calculations pre-
dict different forms of the nuclear equation of state above
and below the normal nuclear matter density, and away
from the symmetric nuclear matter [1, 2, 3, 4, 5, 6]. As a
result, the symmetry energy, which is the difference in en-
ergy between the pure neutron matter and the symmetric
nuclear matter, shows very different behavior above and
below normal nuclear density [6] (see Fig. 1).
In general, two different forms of the density depen-
dence of the symmetry energy have been predicted. One,
where the symmetry energy increases monotonically with
increasing density (“ stiff ” dependence) and the other,
where the symmetry energy increases initially up to nor-
mal nuclear density and then decreases at higher densities
(“ soft ” dependence). Constraining the form of the den-
sity dependence of the symmetry energy is important not
0 0.5 1
ρ / ρ0
0 1 2 3
ρ / ρ0
SkLya
var AV
+δv+3-BF
DD-TW
DD-ρδ
FIG. 1: (Color online) Symmetry energy as a function of den-
sity predicted by microscopic “ ab initio ” calculations. The
left panel shows the low-density region, while the right panel
displays the high-density range. The figure is taken from Ref.
only for a better understanding of the nucleon-nucleon
interaction, and hence its extrapolation to the structure
of neutron-rich nuclei [7, 8, 9, 10], but also for deter-
mining the structure of compact stellar objects such as
neutron stars [11, 12, 13, 14, 15, 16, 17]. For example, a
“ stiff ” form of the density dependence of the symmetry
energy is predicted to lead to a large neutron skin thick-
ness compared to a “ soft ” dependence [8, 10, 18, 19, 20].
Similarly, a “ stiff ” dependence of the symmetry energy
can result in rapid cooling of a neutron star, and a larger
neutron star radius, compared to a “ soft ” density depen-
dence of the symmetry energy [20, 21, 22]. The nuclear
Equation Of State (EOS) is therefore a fundamental en-
tity that determines the properties of systems as diverse
as atomic nuclei and neutron stars, and the knowledge of
which is of significant importance [16, 17, 23].
Experimentally, the best possible means of studying
the nuclear equation of state at sub-normal nuclear den-
http://arxiv.org/abs/0704.0471v1
sity is through intermediate-energy heavy-ion reactions
[26, 27]. In this kind of reaction, an excited nucleus
(the composite of the projectile and the target nucleus)
expands to a sub-nuclear density and disintegrates into
various light and heavy fragments in a process called
multifragmentation. By studying the isotopic yield dis-
tribution of these fragments one can extract important
information about the symmetry energy and its density
dependence. Current studies on the nuclear equation of
state are limited to beams consisting of stable nuclei. It
is hoped that in the future radioactive beam facilities
such as, FAIR (GSI) [24], SPIRAL2 (GANIL) and FRIB
(USA) [25] will provide tremendous opportunities for ex-
ploring the nuclear EOS in regions never before studied
(i.e., extreme isospin and away from normal nuclear den-
sity).
In this work, we have made an attempt to study
the density dependence of the symmetry energy using
two different theoretical approaches for studying multi-
fragmentation, namely the dynamical and the statistical
model approaches of multifragmentation. In section II,
the isoscaling technique used to study the density de-
pendence of the symmetry energy, and their different in-
terpretations in terms of statistical and dynamical ap-
proach, are presented. In section III and IV, a brief
description of the experiment and the experimental re-
sults are presented. The dynamical and the statistical
approaches used to interpret the experimental results are
presented in section V. A comparsion between the two
approaches with other independent studies is presented
in section VI. Finally, a discussion and summary, and
conclusions are presented in section VII and VIII, re-
spectively.
II. SYMMETRY ENERGY AND THE ISOTOPIC
YIELD DISTRIBUTION
It has been shown from experimental measurements
that the ratio of the fragment yields, R21(N ,Z), taken
from two different multifragmentation reactions, 1 and
2, obeys an exponential dependence on the neutron num-
ber (N) and the proton number (Z) of the fragments; an
observation known as isoscaling [28, 29, 30]. The depen-
dence is characterized by the relation,
R21(N,Z) = Y2(N,Z)/Y1(N,Z) = Ce
(αN+βZ) (1)
Where, Y2 and Y1 are the fragment yields from the
neutron-rich and the neutron-deficient systems, respec-
tively. C is an overall normalization factor, and α and β
are the parameters characterizing the isoscaling behavior.
Isoscaling is also theoretically predicted by the dynam-
ical [31, 32, 33, 34, 35] and statistical [36, 37, 38, 39]
models of multifragmentation. In these models, the dif-
ference in the chemical potential of systems with different
neutron-to-proton ration (N/Z) is directly related to the
isoscaling parameter α. The isoscaling parameter α, is
related to the symmetry energy Csym, through the rela-
tion,
4Csym
where, Z1, A1 and Z2, A2 are the charge and the mass
numbers from the two systems and T is the temperature.
This relation provides a simple and straight-forward con-
nection between the symmetry energy and the fragment
isotopic yield distribution.
It must be mentioned that although the above equation
derived from the statistical and the dynamical models of
multifragmentation appears similar in form, the physical
meaning of the terms involved in this equation differ for
each model.
1) In statistical models, the Z/A in Eq. (2) corresponds
to the charge-to-mass ratio of the initial equilibrated frag-
menting system. Whereas, in dynamical models, it cor-
responds to the charge-to-mass ratio of the liquid phase
at a certain time (≈ 300 fm/c) during the dynamical evo-
lution of the colliding systems.
2) The interpretation of the symmetry energy Csym, in
dynamical and statistical models also differs significantly.
The dynamical models relate the symmetry energy in the
above equation to that of the fragmenting source. The
statistical models, on the other hand, relate Csym to that
of the fragments formed at freeze-out.
These conceptual differences between the statistical
and the dynamical models are due to the radically differ-
ent approaches taken in the interpretation of the multi-
fragmentation process. The different interpretation has
also lead to conflicting results from the use of Eq. 2, due
to the different sequential decay effects predicted for the
primary fragments by each model.
The isoscaling parameter α, in Eq. 2 corresponds to
the hot primary fragments which undergo sequential de-
cay into cold secondary fragments. These secondary frag-
ments are the ones that are eventually detected in exper-
iments. The experimentally determined isoscaling pa-
rameter must therefore be corrected for the sequential
decay effect before comparing it to the theoretical mod-
els. It has been observed that while statistical model
calculations show no significant change in the isoscaling
parameter after sequential decay [40], dynamical models
give contrasting results; with some showing no significant
changes [41], while others showing a change of as much
as 50% [42].
In this work, we adopt both theoretical approaches
with their respective interpretations to study the den-
sity dependence of the symmetry energy. In particular,
we use the Antisymmetrized Molecular Dynamics (AMD)
model [31, 43] and the Statistical Multifragmentation
Model (SMM) [36] for this study. A comparison between
the two can provide useful insight into the physical mean-
ing of the above equation in the two models.
III. EXPERIMENT
A. Experimental Setup
The experiments were carried out at the Cyclotron
Institute of Texas A&M University (TAMU) using the
K500 Superconducting Cyclotron and the National Su-
perconducting Cyclotron Laboratory (NSCL) at Michi-
gan State University (MSU). Targets of 58Fe (2.3
mg/cm2) and 58Ni (1.75 mg/cm2) were bombarded with
beams of 40Ar and 40Ca at 33 and 45 MeV/nucleon for
the TAMU measurements [44], and targets of 58Fe (∼ 5
mg/cm2) and 58Ni (∼ 5 mg/cm2) were bombarded with
beams of 40Ar and 40Ca at 25 and 53 MeV/nucleon for
the NSCL measurements [45]. The various combinations
of target and projectile nuclei allowed for a range of N/Z
(neutron-to-proton ratio) (1.04 − 1.23) of the system to
be studied, while keeping the total mass constant (A =
98). In a separate experiment at TAMU, beams of 58Ni
and 58Fe at 30, 40, and 47 MeV/nucleon were also bom-
barded on self-supporting 58Ni and 58Fe targets.
The beams in the TAMU measurements were fully
stripped by allowing them to pass through a thin alu-
minum foil before being hit at the center of the target
inside the TAMU 4π neutron ball [46]. Light charged
particles (Z ≤ 2) and intermediate mass fragments (Z >
2) were detected using six discrete telescopes placed in-
side the scattering chamber of the neutron ball at angles
of 10◦, 44◦, 72◦, 100◦, 128◦ and 148◦. Each telescope con-
sisted of a gas ionization chamber (IC) followed by a pair
of silicon detectors (Si-Si) and a CsI scintillator detector,
providing three distinct detector pairs (IC-Si, Si-Si, and
Si-CsI) for fragment identification. The ionization cham-
ber was of axial field design and was operated with CF4
gas at a pressure of 50 Torr. The gaseous medium was 6
cm thick with a typical threshold of ∼ 0.5 MeV/nucleon
for intermediate mass fragments. The silicon detectors
had an active area of 5 cm × 5 cm and were each sub-
divided into four quadrants. The first and second silicon
detectors in the stack were 0.14 mm and 1 mm thick,
respectively. The dynamical energy range of the silicon
pair was ∼ 16 - 50 MeV for 4He and ∼ 90 - 270 MeV
for 12C. The CsI scintillator crystals that followed the
silicon detector pair were 2.54 cm in thickness and were
read out by photodiodes. Good elemental (Z) identifi-
cation was achieved for fragments that punched through
the IC detector and stopped in the first silicon detector.
Fragments measured in the Si-Si detector pair also had
good isotopic separation. Fragments that stopped in CsI
detectors showed isotopic resolution up to Z = 7. The
trigger for the data acquisition was generated by requir-
ing a valid hit in one of the silicon detectors.
The calibration of the IC-Si detectors were carried out
using the standard alpha sources and by operating the
IC at various gas pressures. The Si-Si detectors were cal-
ibrated by measuring the energy deposited by the alpha
particles in the thin silicon and the punch-through ener-
gies of different isotopes in the thick silicon. The Si-CsI
detectors were calibrated by selecting points along the
different light charged isotopes and determining the en-
ergy deposited in the CsI crystal from the energy loss in
the calibrated Si detector.
The setup for the NSCL experiment consisted of 13
silicon detector telescopes placed inside the MSU 4π Ar-
ray. Four of which were placed at 14◦, each of which
consisted of a 100 µm thick and a 1 mm thick silicon
surface-barrier detector followed by a 20 cm thick plastic
scintillator. Five telescopes were placed at 40◦, in front
of the most forward detectors in the main ball of the 4π
Array. They each consisted of a 100 µm surface-barrier
detector followed by a 5 mm lithium drifted silicon de-
tector. More details can be found in Ref. [45]. Good
isotopic resolution was obtained as in TAMU measure-
ments.
B. Event Characterization
The event characterization of the NSCL data was
accomplished by detection of nearly all the coincident
charged particles by the MSU 4π Array. Data were ac-
quired using two different triggers; the bulk of which
was obtained with the requirement of a valid event in
one of the silicon telescopes. Additional data were taken
with a minimum bias 4π Array trigger for normalization
of the event characterization. The impact parameter of
the event was determined by the mid-rapidity charge de-
tected in the 4π Array as discussed in Ref. [47]. The
effectiveness of the centrality cuts was tested by compar-
ing the multiplicity of events from a minimum bias trigger
with the multiplicity distribution when a valid fragment
was detected at 40◦ [48]. The minimum bias trigger had
a peak multiplicity of charged particles of one, whereas
with the requirement of a fragment at 40◦, the peak of
the multiplicity distribution increased to five.
The event characterization for the TAMU data was
accomplished by using the 4π neutron ball that sur-
rounded the detector assembly. The neutron ball con-
sisted of eleven scintillator tanks segmented in its me-
dian plane and surrounding the vacuum chamber. The
upper and the lower tank were 1.5 m diameter hemi-
spheres. Nine wedge-shaped detectors were sandwiched
between the hemispheres. All the wedges subtended 40◦
in the horizontal plane. The neutron ball was filled with
a pseudocumene-based liquid scintillator mixed with 0.3
% (b.w.) of Gd salt (Gd 2-ethyl hexanoate). Scintilla-
tions from thermal neutrons captured by Gd were de-
tected by twenty 5-inch phototubes : five in each hemi-
sphere, one on each of the identical 40◦ wedges and two
on the forward edges. The efficiency with which the neu-
trons could be detected is about 83%, as measured with
a 252Cf source.
The detected neutrons were used to differentiate be-
tween the central and peripheral collisions. To un-
derstand the effectiveness of neutron multiplicity as
a centrality trigger, simulations were carried out us-
ing a hybrid BUU-GEMINI calculations at various im-
pact parameters for the 40Ca + 58Fe reaction at 33
MeV/nucleon. The simulated neutron multiplicity distri-
bution was compared with the experimentally measured
distribution. The multiplicity of neutron for the impact
parameter b = 0 collisions was found to be higher than
the b = 5 collision. By gating on the 10% highest neu-
tron multiplicity events, one could clearly discriminate
against the peripheral events.
To determine the contributions from noncentral im-
pact parameter collisions, neutrons emitted in coinci-
dence with fragments at 44◦ and 152◦ were calculated at
b = 0 fm and b = 5 fm. The number of events were ad-
justed for geometrical cross sectional differences. A ratio
was made between the number of events with a neutron
multiplicity of at least six, calculated at b = 0 fm, and the
number of events with the same neutron multiplicity at b
= 5 fm. The ratios were observed to be 19.0 and 11.1 at
44◦ and 1.3 and 2.2 at 152◦ for 33 and 45 MeV/nucleon
respectively. At intermediate angles, high neutron multi-
plicities were observed to be outside the region in which
b = 5 fm contributes significantly. At backward angles
the collisions at b = 5 fm made a larger contribution to
the neutron multiplicity.
In addition to the neutron multiplicity distribution, the
charge distribution of the fragments was also used to in-
vestigate the contributions from central and mid-impact
parameter collisions. The b = 5 collisions produced es-
sentially no fragments with charge greater than three in
the 44◦ telescope.
In an earlier work [44], some analysis of the fragment
kinetic energy and charge distributions were presented.
It was shown that at a laboratory angle of 44◦ the ki-
netic energy and the charge distributions are well repro-
duced by the statistical model calculation. Using a mov-
ing source analysis of the fragment energy spectra, it was
also shown that the fragments emitted at backward an-
gles originate from a target-like source, while those emit-
ted at 44◦ originate primarily from a composite source.
In this work, we will concentrate exclusively on data from
the laboratory angle of 44◦, which corresponds to the cen-
ter of mass angle ≈ 90◦, to study the symmetry energy
and the isoscaling properties of the fragments produced.
The choice of this angle enables one to select events which
are predominantly central and undergo bulk multifrag-
mentation. The contributions to the intermediate mass
fragments from the projectile-like and target-like sources
can therefore be assumed to be minimal.
IV. EXPERIMENTAL RESULTS
A. Fragment isotopic yield distribution
The experimentally measured relative isotopic yield
distributions for the Lithium (left) and Carbon (right)
elements, in 58Ni + 58Ni (star symbols), 58Ni + 58Fe
(square symbols), 58Fe + 58Ni (circle symbols) and 58Fe
FIG. 2: Relative yield distribution of the fragments for the
Lithium (left) and Carbon (right) isotopes in 58Ni + 58Ni
(stars and solid lines), 58Fe + 58Ni (circles and dashed lines),
58Ni + 58Fe (squares and dashed lines), and 58Fe + 58Fe (tri-
angles and dotted lines) reactions at various beam energies.
+ 58Fe (triangle symbols) reactions, are shown in Fig. 2
for beam energies of 30, 40 and 47 MeV/nucleon. Sim-
ilarly, the isotopic yield distributions for Lithium (left),
Berillium (center) and Carbon (right) elements, in 40Ca
+ 58Ni (star symbols), 40Ar + 58Ni (circle symbols) and
40Ar + 58Fe (square symbols) reactions, are shown in fig-
ure 3 for beam energies of 25, 33 and 45 MeV/nucleon.
The isotope distribution for each element in Fig. 3
shows higher fragment yield for the neutron rich iso-
topes in 40Ar + 58Fe reaction (squares), which has the
largest neutron-to-proton ratio (N/Z), in comparison to
the 40Ca + 58Ni reaction (stars), which has the small-
est neutron-to-proton ratio. The yields for the reaction,
40Ar + 58Fe (circles), which has an intermediate value of
the neutron-to-proton ratio, are in between those of the
other two reactions. A similar feature is also observed for
the 58Ni + 58Fe, 58Fe + 58Ni and 58Fe + 58Fe reactions
shown in Fig. 2. The fragment yield distributions there-
fore show the isospin dependence of the composite system
on the fragments produced in the multifragmentation re-
action. One also observes that the relative difference in
the yield distribution between the three reactions in each
figure decreases with increasing beam energy. This is due
to the secondary de-excitation of the primary fragments,
a process that becomes important for systems with in-
creasing neutron-to-proton ratio and excitation energy.
FIG. 3: Relative yield distribution of the fragments for
Lithium (left), Berillium (center) and Carbon (right) isotopes
in 40Ca + 58Ni (stars and solid lines), 40Ar + 58Ni (circles and
dashed lines), and 40Ar + 58Fe (squares and dotted lines) re-
actions at various beam energies.
B. Isotopic and Isotonic scaling
As discussed in section II, the ratio of isotope
yields in two different systems, 1 and 2, R21(N,Z) =
Y2(N,Z)/Y1(N,Z), follows an exponential dependence
on the neutron number (N) and the proton number (Z)
of the isotopes in relation known as isoscaling.
In Fig. 4, we show the isotopic yield ratio as a func-
tion of neutron number N , for Ar + Fe, Ar + Ni and
Ca + Ni systems at beam energies of 25, 33, 45 and
53 MeV/nucleon. The left column shows the ratio for
the 40Ar + 58Fe and 40Ca + 58Ni pair of reaction and
the right column shows the ratio for the 40Ar + 58Ni and
40Ca + 58Ni pair of reaction. One observes that the ratio
for each element shows linear behavior in the logarithmic
plot and aligns with the neighboring element quite well.
This feature is observed for all the beam energies and
both pairs of reactions studied. One also observes that
the alignment of the data points varies with beam ener-
gies as well as the pairs of reaction. To have a quantita-
tive estimate of this variation, the ratio for each element
(Z) was simultaneously fit using an exponential relation
(shown by the solid lines) to obtain the slope parame-
ter α. The values of the parameters are shown at the
top of each panel in the figure. The value of the slope
FIG. 4: Experimental isotopic yield ratios of the fragments as
a function of neutron number N , for various beam energies.
The left column correspond to 40Ar + 58Fe and 40Ca + 58Ni
pair of reactions. The right column correspond to 40Ar + 58Ni
and 40Ca + 58Ni pair of reactions. The different symbols
correspond to Z = 3 (circles), Z = 4 (open stars), Z = 5
(triangles), Z = 6 (squares) and Z = 7 (filled stars) elements.
The lines are the exponential fits to the data as explained in
the text.
parameter α is larger for the 40Ar + 58Fe and 40Ca +
58Ni reactions, which has a larger difference in the N/Z
of the systems in the pair, compared to the 40Ar + 58Ni
and 40Ca + 58Ni reactions, which has a smaller difference
in the corresponding N/Z. The α value furthermore de-
creases with increasing beam energy. A similar feature
is also observed in Fe + Fe, Fe + Ni, Ni + Fe and Ni
+ Ni systems. Fig. 5 shows the isotope yield ratios and
the isotone yield ratios for the Fe + Fe and Ni + Ni re-
actions for the 30 MeV/nucleon beam energy. A relative
comparison of how the isoscaling parameter α, evolves as
a function of beam energy and the isospin of the system
is shown in Fig. 6. The figure clearly shows that the α
value decreases with beam energy from 25 MeV/nucleon
to 53 MeV/nucleon. In addition, there is also a clear
drop in the α values with the decrease of the N/Z values
of the system.
FIG. 5: Experimental isotope yield ratios (top) and isotone
yield ratios (bottom) from 58Fe + 58Fe and 58Ni + 58Ni re-
actions as a function of N and Z for 30 MeV/nucleon beam
energy. The solid lines are fit to the data as discussed in the
text.
FIG. 6: Experimental isoscaling parameter α, as a function
of the beam energy. The solid circles are for the 40Ar + 58Fe
and 40Ca + 58Ni reactions. The open triangles are for 58Fe +
58Fe and 58Ni + 58Ni reactions. The solid stars are for 40Ar
+ 58Ni and 40Ca + 58Ni reactions. The open squares are for
58Fe + 58Ni and 58Ni + 58Ni reactions.
V. THEORETICAL MODEL COMPARISON
A. Dynamical AMD model
The Antisymmetrized Molecular Dynamics (AMD)
[31, 43] is a microscopic model that simulates the time
evolution of a nuclear collision. The colliding system in
this model is represented in terms of a fully antisym-
metrized product of Gaussian wave packets. During the
evolution, the wave packet centroids move according to
the deterministic equation of motion. The followed state
of the simulation branches stochastically and successively
into a huge number of reaction channels. The interactions
are parameterized in terms of an effective force acting be-
tween nucleons and the nucleon-nucleon collision cross-
sections. The advantage of using a dynamical model to
study the nuclear equation of state is that it allows one
to understand the functional form of the density depen-
dence of the symmetry energy at a very fundamental level
(i.e., from the basic nucleon-nucleon interaction).
Recently [31], the fragment yields from heavy ion col-
lisions simulated within the Antisymmetrized Molecular
Dynamics (AMD) calculation were reported to follow a
scaling behavior of the type shown in Eq. 1. A linear re-
lation between the isoscaling parameter α and the differ-
ence in the isospin asymmetry (Z/A)2 of the fragments
(as given in Eq. 2), with appreciably different slopes,
was predicted for two different forms of the density de-
pendence of the symmetry energy; a “ stiff ” dependence
(obtained from Gogny-AS interaction) and a “ soft ” de-
pendence (obtained from Gogny interaction).
In this section, we compare the experimentally de-
termined isoscaling parameter with the predictions of
the AMD model calculation. The isospin asymmetry
of the fragments for the present systems was estimated
at t = 300 fm/c of the dynamical evolution using the
AMD calculation. The values for the fragment asym-
metry (Z/A)2, were obtained by interpolating between
those calculated for the 40Ca + 40Ca, 48Ca + 48Ca and
60Ca + 60Ca systems by Ono et al. [31]. These systems
are symmetric and nearly similar in charge and mass as
studied in the present work. Fig. 7 shows the AMD cal-
culation of the fragment asymmetry, (Z/A)2 at t = 300
fm/c, as a function of initial asymmetry at time t = 0
fm/c, for two different choices of the nucleon-nucleon in-
teraction, Gogny and Gogny-AS. The asymmetry values
for the 40Ca + 40Ca, 48Ca + 48Ca and 60Ca + 60Ca sys-
tems of Ref. [31] are shown by solid and hollow square
symbols for the Gogny and Gogny-AS interaction, re-
spectively. The lines are the linear fits to the calcula-
tions. The interpolated values for the present systems
are shown by the solid circles and triangles, and hollow
circles and triangles for the Gogny and Gogny-AS inter-
action, respectively.
We note that the AMD calculations carried out in Ref.
[31] and shown in Fig. 7 correspond to the beam energy
of 35 MeV/nucleon. The interpolated values of the asym-
metries for the present systems obtained from Fig. 7 are
FIG. 7: AMD calculations of the fragment asymmetry
(Z/A)2, at t = 300 fm/c for the Gogny (solid line and solid
squares) and Gogny-AS (dotted line and hollow squares) in-
teractions at 35 MeV/nucleon. The calculations are taken
from Ref. [31] for the systems shown by the square symbols.
The lines are linear fit to the square symbols. The other sym-
bols are the interpolated values for the systems studied in this
work.
therefore for the beam energy of 35 MeV/nucleon. In or-
der to compare the experimentally determined isoscaling
parameter to that of the calculations, we therefore make
use of the experimental isoscaling parameter for the beam
energy of 35 MeV/nucleon using Fig. 6.
Fig. 8 shows a comparison between the experimentally
observed α and those from the AMD model calculations
plotted as a function of the difference in the fragment
asymmetry for the beam energy of 35 MeV/nucleon. The
solid and the dotted lines are the AMD model predictions
for the “ soft ” (Gogny) and the “ stiff ” (Gogny-AS)
form of the density dependence of the symmetry energy,
respectively. The solid and the hollow symbols (squares,
stars, triangles and circles) are the results of the present
study for the two different values of the fragment asym-
metry, assuming Gogny and Gogny-AS interactions, re-
spectively. Also shown in the figure are the scaling pa-
rameters (asterisks, crosses, diamond and inverted tri-
angle) taken from various other works in the literature
[36, 49]. It is observed that the experimentally deter-
mined α parameter increases linearly with increasing dif-
ference in the fragment asymmetry of the two systems
as predicted by the AMD calculation. Also, the data
points are in closer agreement with those predicted by
the Gogny-AS interaction (dotted line) than those from
the usual Gogny force (solid line).
In the above comparison between the data and the cal-
FIG. 8: Isoscaling parameter α, as a function of the differ-
ence in fragment asymmetry for 35 MeV/nucleon. The solid
and the dotted lines are the AMD calculations for the Gogny
and Gogny-AS interactions, respectively [31]. The solid and
the hollow squares, stars, triangles and circles are from the
present work as described in the text. The other symbols cor-
responds to data taken from [49] (asterisks) and [36] (crosses,
diamonds, inverted triangles).
culation, the corrections for the isoscaling parameter α
due to the sequential de-excitation of the fragments are
not taken into account. The slightly lower values of the
isoscaling parameters (symbols) from the present mea-
surements with respect to the Gogny-AS values (dotted
line) could be due to the small secondary de-excitation
effect of the fragments not accounted for in this compar-
ison. Recently, it has been reported by Ono et al. [42],
that the sequential decay effect in the dynamical calcu-
lations can affect the α value by as much as 50 %, and
the ability to distinguish between the “ stiff ” and the “
soft ” form of the density dependence of the symmetry
energy diminishes significantly. The calculations by Ono
et al., were carried out for the above studied systems
using the AMD model. However, dynamical calculation
carried out by Tian et al. [41], using Isospin Quantum
Molecular Dynamic (IQMD) model shows no significant
difference between the primary and the secondary α. The
sequential decay effect from the IQMD calculation was
also carried out for the same systems and beam energy
as studied by Ono et al. [42] using the AMD model. The
contrasting results between the two dynamical calcula-
tions for the same systems and energy currently present
significant amount of uncertainty in reliably estimating
the effect of sequential decay using dynamical models.
One reason for this could be the large discrepancy that
exists in the determination of the primary fragment exci-
tation energy from the current dynamical models. It has
been shown using another dynamical model (stochastic
mean field calculation) (see Liu et al. [33]), that it re-
quires a significantly lower value of the primary fragment
excitation energy (by as much as 50%), to be able to
reproduce the experimentally observed fragment isotope
distribution.
In the above comparison between the data and the cal-
culation, we have therefore assumed the effect of the se-
quential decay to be negligible. A correction of about 10
- 15 %, as determined and well established from various
statistical model studies [28], results in a slight increase
in the α values bringing them even closer to the dotted
line. The observed agreement of the experimental data
with the Gogny-AS type of interaction therefore appears
to suggest a “ stiff ” form of the density dependence of
the symmetry energy.
Figure 9 shows various forms of the density dependence
of the symmetry energy in isospin asymmetric nuclear
matter used by Chen et al. [50], and those used in the
present dynamical model analysis. The dot-dashed, dot-
ted and the dashed curve corresponds to the momentum
dependent Gogny interactions used by Chen et al., to ex-
plain the NSCL-MSU isospin diffusion data. Assuming
that the density dependence of the symmetry energy can
be parametrized as,
Csym(ρ) = C
(MeV ) (3)
where Cosym, is the value of the symmetry energy at sat-
uration density and γ is the parameter that characterizes
the stiffness of the symmetry energy, the above depen-
dences used by Chen et al. can be written as, Esym ≈
31.6 (ρ/ρ◦)
γ , where, γ = 1.6, 1.05 and 0.69, respectively.
The solid curve and the solid point in Fig. 9 correspond
to those from the Gogny and Gogny-AS interactions used
to study the isoscaling data in the present work.
By parameterizing the density dependence of the sym-
metry energy that explains the present isoscaling data,
one obtains, Csym(ρ) ≈ 31.6 (ρ/ρ◦)
γ , where γ = 0.69,
from the dynamical model analysis.
B. Statistical Multifragmentation Model
The Statistical Multifragmentation Model (SMM) [51,
52] is based on the assumption of statistical equilibrium
at a low density freeze-out stage. All breakup channels
composed of nucleons and excited fragments are taken
into account and considered as partitions. During each
partition the conservation of mass, charge, energy, mo-
mentum and angular momentum is taken into account,
and the partitions are sampled uniformly in the phase
space according to their statistical weights using Monte
Carlo sampling. The Coulomb interaction between the
fragments is treated in the Wigner-Seitz approximation.
Light fragments with mass number A ≤ 4 are considered
FIG. 9: Different forms of the density dependence of the nu-
clear symmetry energy used in the dynamical analysis of the
present measurements on isoscaling data and the isospin dif-
fusion measurements of NSCL-MSU [50]. The curves are as
described in the text.
as elementary particles with only translational degrees
of freedom (“nuclear gas”). Fragments with A > 4 are
treated as heated nuclear liquid drops, and their individ-
ual free energies FA,Z are parametrized as a sum of the
volume, surface, Coulomb and symmetry energy.
For the present study we make use of the SMM ver-
sion adopted by Botvina et al. [36]. In this version,
the secondary de-excitation of large fragments with A
> 16 is described by Weisskopf-type evaporation and
Bohr-Wheeler-type fission models [51, 53]. The decay
of smaller fragments is treated with the Fermi-breakup
model. All ground and nucleon-stable excited states of
light fragments are taken into account and the popula-
tion probabilities of these states are calculated according
to the available phase space [53]. The sequential decay
effect on the isoscaling parameter in this version of SMM
has been established to be small and in good agreement
with other versions of the statistical models.
Unlike dynamical calculations, the form of the density
dependence of the symmetry energy is not known a priori,
but has to be deduced from the systematic correlations
between the isoscaling parameter, temperature, symme-
try energy and the density of the multifragmenting sys-
tem. To build this correlation, we make use of the frag-
ment yield distributions measured in 58Ni, 58Fe + 58Ni,
58Fe reactions at 30, 40 and 47 MeV/nucleon to study
the isoscaling parameter α, as a function of the excita-
tion energy of the fragmenting source. The parameter α
was obtained from the ratio’s of the isotopic yields for two
different pairs of reactions, 58Fe + 58Ni and 58Ni + 58Ni,
and 58Fe + 58Fe and 58Ni + 58Ni as discussed in section
4.2. The excitation energy of the source for each beam
energy was determined by simulating the initial stage of
the collision dynamics using the Boltzmann-Nordheim-
Vlasov (BNV) model calculation [54]. The results were
obtained at a time around 40 - 50 fm/c after the projec-
tile had fused with the target nuclei and the quadrupole
moment of the nucleon coordinates (used for identifica-
tion of the deformation of the system) approached zero.
These excitation energies were also compared with those
obtained from the systematic calorimetric measurements
(see Ref. [55]) for systems with mass (A ∼ 100), and sim-
ilar to those studied in the present work, and are in good
agreement. Fig. 10 (a) shows the experimental isoscaling
parameter α, as a function of the excitation energy for
Fe + Fe and Ni + Ni, and Fe + Ni and Ni + Ni pairs
of reactions. A systematic decrease in the absolute val-
ues of the isoscaling parameter with increasing excitation
energy is observed for both pairs. The α parameters for
the 58Fe + 58Fe and 58Ni + 58Ni are about twice as large
compared to those for the 58Fe + 58Ni and 58Ni + 58Ni
pair of reactions.
The experimental isoscaling parameter was compared
with the predictions of the Statistical Multifragmentation
Model (SMM) [51, 56] calculations to study their depen-
dence on the excitation energy and the isospin content.
The initial parameters such as, the mass, charge and ex-
citation energy of the fragmenting source for the calcula-
tion was obtained from the BNV calculations as discussed
above. The possible uncertainties in the source param-
eters due to the loss of nucleons during pre-equilibrium
emission was accounted for by carrying out the calcula-
tions for smaller source sizes. The break-up density in
the calculation was taken to be multiplicity-dependent
and was varied from approximately 1/2 to 1/3 the sat-
uration density. This was achieved by varying the free
volume with the excitation energy as shown in Ref. [51].
The form of the dependence was adopted from the work
of Bondorf et al. [57, 58], (and shown by the solid curve in
Fig. 10 (d)). It is known that the multiplicity-dependent
break-up density, which corresponds to a fixed interfrag-
ment spacing and constant pressure at break-up, leads
to a pronounced plateau in the caloric curve [57, 58]. A
constant break-up density would lead to a steeper tem-
perature versus excitation energy dependence.
The symmetry energy in the calculation was varied un-
til a reasonable agreement between the calculated and the
measured α was obtained. Fig. 10 (a) shows the com-
parison between the SMM calculated and the measured
α for both pairs of systems. The dashed curves corre-
spond to the calculation for the primary fragments and
the solid curves to the secondary fragments. The width
in the curve is the measure of the uncertainty in the in-
puts to the SMM calculation. It is observed that, within
the given uncertainties, the decrease in the α values with
increasing excitation energy and decreasing isospin dif-
ference ∆(Z/A)2, of the systems is well reproduced by
4 6 8 10
0 2 4 6 8 10
E* (MeV/nucleon)
4 6 8 10
0 2 4 6 8 10
FIG. 10: (Color online) Isoscaling parameter α, temperature,
symmetry energy and density as a function of excitation en-
ergy for the Fe + Fe and Ni + Ni (inverted triangles), and
Fe + Ni and Ni + Ni (solid circles) reactions at 30, 40 and
47 MeV/nucleon. a) Experimental isoscaling parameter as
a function of excitation energy. The solid and the dashed
curves are the SMM calculations as discussed in the text. b)
Temperature as a function of excitation energy. The solid
stars correspond to the measured values and are taken from
Ref. [55]. The solid and the dashed curve corresponds to the
Fermi-gas relation. The dotted curve corresponds to the one
obtained from Eq. 4. c) Symmetry energy as a function of
excitation energy. d) Density as a function of excitation en-
ergy. The solid stars correspond to those from Ref. [69]. The
open triangles are those from Ref. [70]. The solid curve is
from Ref. [57].
the SMM calculation. One also notes that the effect of
sequential decay on the isoscaling parameter is small as
observed in several other studies [40, 59] using statistical
models.
We show in Fig. 10 (b), the temperature as a func-
tion of excitation energy (caloric curve) obtained from
the above SMM calculation that uses the excitation en-
ergy dependence of the break-up density to explain the
observed isoscaling parameters. These are shown by the
solid and inverted triangle symbols. Also shown in the
figure are the experimentally measured caloric curve data
compiled by Natowitz et al. [55], from various measure-
ments for this mass range. The data from these mea-
surements are shown collectively by solid star symbols
and no distinction is made among them. The Fermi-gas
model predictions with inverse level density parameter
Ko = 10 (solid and dashed curve), is also shown. It is
evident from the figure that the temperatures obtained
from the SMM calculations are in good agreement with
the overall trend of the caloric curve. Somewhat lower
value for the temperature is observed when the break-up
density of the system is kept constant at 1/3 the nor-
mal nuclear density. By allowing the break-up density
to evolve with the excitation energy, a near plateau that
agrees with the experimentally measured caloric curves
is obtained. This assures that the input parameters used
in the SMM calculation for comparing with the data are
reasonable.
The symmetry energies obtained from the statistical
model comparison of the experimental isoscaling param-
eter α, are as shown in Fig. 10 (c). A steady decrease in
the symmetry energy with increasing excitation energy is
observed for both pairs of systems. Such a decrease has
also been observed in several other studies [60, 61, 62, 63].
We have also estimated the effect of the symmetry en-
ergy evolving during the sequential de-excitation of the
primary fragments [60, 64]. These are reflected in the
large error bars shown in Fig. 10 (c).
The phase diagram of the multifragmenting system is
two dimensional and hence the excitation energy depen-
dence of the temperature (the caloric curve) must take
into account the density dependence too. Often this de-
pendence is neglected while studying the caloric curve.
In the following, we attempt to extract the density of the
fragmenting system as a function of excitation energy. It
has been shown by Sobotka et al. [65], that the plateau
in the caloric curve could be a consequence of the ther-
mal expansion of the system at higher excitation energy
and decreasing density. By assuming that the decrease
in the breakup density, as taken in the present statisti-
cal multifragmentation calculation, can be approximated
by the expanding Fermi gas model, and furthermore the
temperature in Eq. 2 and the temperature in the Fermi-
gas relation are related, one can extract the density as a
function of excitation energy using the relation
Ko(ρ/ρo)2/3E∗ (4)
0 0.5 1 1.5
FIG. 11: (Color online) Symmetry energy as a function of
density for the Fe + Fe and Ni + Ni pair of reaction (inverted
triangles), and Fe + Ni and Ni + Ni pair of reactions (solid
circles) for the 30, 40 and 47 MeV/nucleon. The solid curve is
the dependence obtained form the dynamical model analysis
as explained in the text.
In the above expression, the momentum and the fre-
quency dependent factors in the effective mass ratio are
taken to be one as expected at high excitation energies
and low densities studied in this work [66, 67, 68].
The resulting densities for the two pairs of systems are
shown in Fig. 10 (d) by the solid circles and inverted
triangles. For comparison, we also show the break-up
densities obtained from the analysis of the apparent level
density parameters required to fit the measured caloric
curve by Natowitz et al. [69], and those obtained by Viola
et al. [70] from the Coulomb barrier systematics that are
required to fit the measured intermediate mass fragment
kinetic energy spectra. One observes that the present re-
sults obtained by requiring to fit the measured isoscaling
parameters and the caloric curve are in good agreement
with those obtained by Natowitz et al. The figure also
shows the fixed freeze-out density of 1/3 (dashed line)
and 1/6 (dotted line) of the saturation density assumed in
various statistical model comparisons. The caloric curve
obtained using the above densities and excitation ener-
gies (shown by solid stars, circles and the triangles) with
Ko = 10 in Eq. 4, is shown by the dotted curve in Fig.
10 (b). The small discrepancy between the dotted curve
and the data (solid stars) below 4 MeV/nucleon is due
to the approximate nature of Eq. 4 being used.
It is therefore evident from figure 10 (a), (b), (c) and
(d) that the decrease in the experimental isoscaling pa-
rameter α, symmetry energy, break-up density, and the
flattening of the temperature with increasing excitation
energy are all correlated. One can thus conclude that the
expansion of the system during the multifragmentation
process leads to a decrease in the isoscaling parameter,
decrease in the symmetry energy and density, and the
flattening of the caloric curve.
TABLE I: Parameterized form of the density dependence of the symmetry energy obtained from various independent studies.
Reference Parametrization Studies
Fuchs et al. [85] 32.9(ρ/ρo)
0.59 Relativistic Dirac-Brueckner calculation
Heiselberg et al. [82] 32.0(ρ/ρo)
0.60 Variational calculation
Danielewicz et al. [81] 31(33)(ρ/ρo)
0.55(0.79) BE, skin, isospin analog states
Tsang et al. [79] 12.125(ρ/ρo)
2 Isospin diffusion
Chen et al. [50] 31.6(ρ/ρo)
1.05 Isospin diffusion
Li et al. [80] 31.6(ρ/ρo)
0.69 Isospin diffusion
Piekarewicz et al. [77, 78] 32.7(ρ/ρo)
0.64 Giant resonances
Shetty et al. [73, 74, 75] 31.6(ρ/ρo)
0.69 Isotopic distribution
Famiano et al. [87] 32.0(ρ/ρo)
0.55 neutron-proton emission ratio
Tsang et al. [28] 23.4(ρ/ρo)
0.6 Isotopic distribution
From the above correlation between the symmetry en-
ergy as a function of excitation energy, and the density
as a function of excitation energy, we obtain the sym-
metry energy as a function of density. This is shown by
the inverted triangles and solid circles in Fig. 11 for the
Fe + Fe and Ni + Ni, and the Fe + Ni and Ni + Ni
pair of reactions. The temperature in the present work
remains nearly constant for the range of excitation ener-
gies studied, the observed decrease in the symmetry en-
ergy with increasing excitation energy is therefore a con-
sequence of decreasing density. This is also supported
by microscopic calculations which shows an extremely
slow evolution of the symmetry energy with temperature
[71, 72]. The evolution is practically negligible for the
temperature range studied in this work. The solid curve
in Fig. 11 corresponds to the dependence Csym(ρ) =
31.6 (ρ/ρ◦)
0.69 MeV, obtained from the dynamical Anti-
symmetrized Molecular Dynamic (AMD) calculation, as
discussed in the previous section. It is thus observed that
the dynamical and statistical models lead to similar den-
sity dependence of the symmetry energy.
VI. COMPARISON WITH OTHER
INDEPENDENT STUDIES
In the following, we compare the form of the den-
sity dependence of the symmetry energy obtained from
the present experimentally measured isoscaling param-
eter using the statistical and the dynamical multifrag-
mentation models with several other recent independent
studies. Fig. 12 shows this comparison. The green solid
curve corresponds to the one obtained from Gogny-AS
interaction in dynamical AMD model that explains the
present results [73, 74], assuming the sequential decay
effect to be small. The inverted triangle and the circle
symbols also correspond to the present measurements ob-
tained by comparing with the statistical multifragmenta-
tion model [75]. The red dashed curve corresponds to the
one obtained recently from an accurately calibrated rela-
tivistic mean field interaction, used for describing the Gi-
ant Monopole Resonance (GMR) in 90Zr and 208Pb, and
the IVGDR in 208Pb by Piekarewicz et al. [76, 77, 78].
The pink dot-dashed curve correspond to the one used
to explain the isospin diffusion results of NSCL-MSU us-
ing the isospin dependent Boltzmann-Uehling-Uhlenbeck
(IBUU) model by Tsang et al. [79]. The blue dot-dashed
curve also corresponds to the one used for explaining the
isospin diffusion data of NSCL-MSU by Chen et al. [50],
but with the momentum dependence of the interaction
included in the IBUU calculation. This dependence has
been further modified to include the isospin dependence
of the in-medium nucleon-nucleon cross-section by Li et
al. [80], and is in good agreement with the present study.
The shaded region in the figure corresponds to those ob-
tained by constraining the binding energy, neutron skin
thickness and isospin analogue state in finite nuclei us-
ing the mass formula of Danielewicz [81]. The yellow
solid curve correspond to the parameterization adopted
by Heiselberg et al. [82] in their studies on neutron stars.
By fitting earlier predictions of the variational calcula-
tions by Akmal et al. [83, 84], where the many-body and
special relativistic corrections are progressively incorpo-
rated, Heiselberg and Hjorth-Jensen obtained a value of
Cosym = 32 MeV and γ = 0.6, similar to those obtained
from the present measurements. A similar result is also
obtained from the relativistic Dirac-Brueckner calcula-
tion, with Cosym = 32.9 MeV and γ = 0.59 [85]. The
Dirac-Brueckner is an “ab-initio” calculation based on
nucleon-nucleon interaction with Bonn A type potential
instead of the AV18 potential used in the variational cal-
culation of Ref. [84]. The density dependence of the sym-
metry energy has also been studied in the framework of
expanding emitting source (EES) model by Tsang et al.
[28], where a power law dependence of the form Csym(ρ)
= 23.4(ρ/ρo)
γ , with γ = 0.6 was obtained. This depen-
dence (shown by the black dotted curve) is significantly
softer than other dependences shown in the figure. The
solid square point in the figure correspond to the value
of symmetry energy obtained by fitting the experimental
differential cross-section data in a charge exchange reac-
tion using the isospin dependent CDM3Y6 interaction of
FIG. 12: (Color online) Comparison between the results on the density dependence of the symmetry energy obtained from
various different studies. The various curves and the symbols are described in the text.
the optical potential by Khoa et al. [86].
An alternate observable, the double neutron/proton
ratio of nucleons taken from two reaction systems using
four isotopes of the same element, has recently been pro-
posed as a probe to study the density dependence of the
symmetry energy [87]. This observable is expected to be
more robust than the isoscaling observable. It was shown
recently [87] that the experimentally determined double-
ratio for the 124Sn + 124Sn reaction to that for the 112Sn
+ 112Sn reaction, results in a dependence with γ = 0.5
(shown by black dashed curve in Fig. 12), when com-
pared to the predictions of the IBUU transport model
calculations. This observation is in close agreement with
other studies discussed above. However, this dependence
has been obtained by using the momentum independent
calculation of Ref. [88]. A more recent calculation [89]
using a BUU transport model that includes momentum
dependent interaction show significantly lower values for
the double neutron/proton ratio of free nucleons com-
pared to the one reported by Famiano et al.
The parameterized forms of the density dependence of
the symmetry energy obtained from all the above men-
tioned studies are as given in Table I. The close agree-
ment between various independent studies show that a
constraint on the density dependence of the symmetry
energy, given as Csym(ρ) = C
sym(ρ/ρo)
γ , where Cosym ≈
31 - 33 MeV and γ ≈ 0.55 - 0.69 can be obtained.
VII. DISCUSSION
We make the following observations from the above
comparison between the statistical and the dynamical
model analysis :
1) Assuming a negligibly small sequential decay effect,
the form of the density dependence of the symmetry en-
ergy obtained from the dynamical model analysis is in
good agreement with the one obtained from the statistical
model analysis: As mentioned earlier, the sequential de-
cay effect among various dynamical model calculations is
still a subject of debate. The statistical models however
consistently show small sequential decay effect. If the se-
quential decay in both the dynamical and the statistical
model is determined by the excitation energy, charge (Z)
and mass (A) of the fragments, and not by the process
that leads to these fragments, the de-excitation of the
fragments must lead to a same amount of change in the
isoscaling parameter (either a large change or no change
at all). It is therefore unrealistic to assume that the se-
quential decay effect is different in the dynamical and the
statistical model calculations. One comparison by Hudan
et al. [90], show good agreement between the experimen-
tally determined primary fragment excitation energy and
those calculated using the statistical multifragmentation
model (see table II of Ref. [90]). Furthermore, if dy-
namical and statistical models are merely two different
ways of interpreting the same multifragmentation process
(i.e., one simulating the entire process from the forma-
tion to the breakup stage, and the other simulating only
the later breakup stage), the isoscaling parameter from
both interpretation must lead to consistent results. It
is well known, and as discussed in section II, that both
interpretations predict isoscaling in multifragmentation.
As discussed in section V A, the apparent disagreement
between the sequential decay effect in statistical and dy-
namical models, could be due to the large discrepancy
that exists in the determination of the primary fragment
excitation energy from current dynamical models.
It has been argued [91] that the effect of sequential
decay on the isoscaling parameter α, in statistical multi-
fragmentation model depends not only on the excitation
energy but also on the value of the symmetry energy. The
fragments in their primary stage are usually hot and the
properties of hot nuclei (i.e., their binding energy and
mass) differ from those of cold nuclei. If hot fragments
in the freeze-out configuration have smaller symmetry
energy, their mass at the beginning of the sequential de-
excitation will be different and this effect should be taken
into account. At smaller values of the symmetry energy
the sequential decay effect can be large. In order to study
the effect of symmetry energy evolution on isoscaling pa-
rameter during sequential decay, we have adopted in this
work a phenomenological approach of Buyukcizmeci et
al. [64]. Fig. 13 and 14 shows the primary and the sec-
ondary isoscaling parameter as a function of symmetry
energy calculated from the statistical multifragmentation
model (SMM) for the Ar + Ni and Ca + Ni, and Ar + Fe
and Ca + Ni pair of reactions. The various panels from
top to bottom correspond to different system excitation
energies. Fig. 13 shows the result of the calculations
where the symmetry energy is kept fixed, and Fig. 14
shows the result for the calculations where the symme-
try energy is varied during the de-excitation process. The
dashed lines in the figure correspond to the primary frag-
ments (Eq. 2) and the solid lines to the secondary frag-
ments. It is observed that there is no significant change
in the primary and the secondary alpha.
2) The result of the statistical model analysis is in good
agreement with other independent studies : A comparison
between the density dependence of the symmetry energy
obtained from the statistical model analysis (for which
the sequential decay effect is known to be small) and
other independent studies shows excellent agreement.
3) The isoscaling parameter probes the property of in-
finite nuclear matter : The symmetry energy obtained
from dynamical model analysis (shown by the solid curve
in Fig. 11) relates to the volume part of the symmetry en-
ergy as in infinite nuclear matter, whereas, the symmetry
energy obtained from the statistical model analysis (solid
circles and inverted triangles in Fig. 11) relates to the
FIG. 13: SMM calculated isoscaling parameter α as a func-
tion of symmetry energy for various excitation energies. The
open circles joined by dotted lines correspond to the primary
fragments and the open stars joined by solid lines to the sec-
ondary fragments. The left column shows the calculation for
40Ar + 58Ni and 40Ca + 58Ni pair, and the right column for
the 40Ar + 58Fe and 40Ca + 58Ni pair.
fragment that is finite and has surface contribution. The
similarity between the two can probably be understood
in terms of the weakening of the surface symmetry free
energy when the fragments are being formed. During
the density fluctuation in uniform low density matter,
the fragments are not completely isolated and continue
to interact with each other, resulting in a decrease in
the surface contribution as predicted by Ono et al. [43].
It must be mentioned that by fitting the nuclear masses
with mass formula, a volume contribution to the symme-
try energy of about 31 - 33 MeV and surface contribution
of about 11 - 13 MeV was obtained by Danielewicz [92]
for nuclei at normal density. Using the constraint ob-
tained for the volume part of the symmetry energy from
the present study, and following the expression for the
symmetry energy of finite nuclei by Danielewicz, we write
the general expression for the density dependence of the
symmetry energy as,
SA(ρ) =
α(ρ/ρ◦)
1 + [α(ρ/ρ◦)γ/βA1/3]
where, α ≡ Cosym = 31 - 33 MeV, γ = 0.55 - 0.69 and
α/β = 2.6 - 3.0. The quantities α and β are the volume
and the surface symmetry energy at normal nuclear den-
sity. The above equation reduces to Eq. 3 for infinite
FIG. 14: Same as in fig. 13, but with the modified secondary
de-excitation with evolving symmetry energy.
nuclear matter in the limit of A → ∞, and to the sym-
metry energy of finite nuclei for ρ = ρ◦. The ratio of the
volume symmetry energy to the surface symmetry energy
(α/β) is closely related to the neutron skin thickness. De-
pending upon how the nuclear surface and the Coulomb
contribution is treated, two different correlations between
the volume and the surface symmetry energy have been
predicted [17] from fits to nuclear masses. Experimental
masses and neutron skin thickness measurements for nu-
clei with N/Z > 1 should provide tighter constraint on
the surface-volume correlation.
4) The density dependence of the symmetry energy ob-
tained using the statistical model approach is consistent
with other experimentally determined observables : In the
past, attempts have been made to study the density de-
pendence of the symmetry energy by looking at specific
observables and comparing them with the predictions of
the dynamical models. Such an approach attempts to ex-
plain the observable of interest without trying to simul-
taneously explain other properties, such as, the temper-
ature, density and excitation of the fragmenting system.
This has lead to a variety of different dependences with-
out any clue to what density is being probed. It might not
be straightforward to distinguish different realistic EOS
interactions using dynamical models, due to the large
uncertainties that currently exist in the sequential decay
effects for these models. But the idea of extracting infor-
mation on the symmetry energy from the point of view
of the basic nucleon-nucleon interaction is a very pow-
erful approach. On the other hand, the determination
of the density dependence of the symmetry energy from
statistical model analysis by simultaneously explaining
the isoscaling parameter, caloric curve and the density
as a function of excitation energy is a reverse approach.
This approach attempts to explain the experimental ob-
servables without any prior knowledge of the governing
interaction and arrives at a dependence which can then
be compared with those predicted from the basic inter-
actions.
5) Symmetry energy determined from the multifrag-
mentation study is lower than that of normal nuclei
: Theoretical many-body calculations [1, 93, 94, 95]
and those from the empirical liquid-drop mass formula
[96, 97] predict symmetry energy near normal nuclear
density and temperature to be around 30 - 32 MeV. As-
suming a negligible evolution of the symmetry energy as
a function of temperature, as shown in Ref. [71, 72], the
present statistical model analysis yields symmetry energy
of the order of 18 - 20 MeV at half the normal nuclear
density.
6) The above constraint on the density dependence of
the symmetry energy has important implications for as-
trophysical and nuclear physics studies :
a) Neutron skin thickness: It has been shown [98] that
an empirical fit to a large number of mean-field calcu-
lations yield neutron skin thickness for 208Pb nucleus,
Rn - Rp ≃ (0.22γ + 0.06) fm, where γ is the exponent
that determines the stiffness of the density dependence
of the symmetry energy. From the above comparison, as-
suming only those density dependences of the symmetry
energy which have symmetry energy at normal nuclear
density 31 - 33 MeV, one obtains a neutron skin thick-
ness of 0.18 - 0.21 fm. An accurate determination of
the neutron skin thickness from the parity-violating elec-
tron scattering measurement [18] will provide a unique
observational constraint on the thickness of the neutron
skin of a heavy nucleus. The above constraint also leads
to symmetric matter nuclear compressibility of K ∼ 230
b) Neutron star mass and radius : The constraint also
predicts a limiting neutron star mass of Mmax = 1.72
solar mass and a radius, R = 11 - 13 km for the “ canon-
ical ” neutron star. Recent observations of pulsar-white
dwarf binaries at the Arecibo observatory suggest a pul-
sar mass for PSRJ0751+1807 of M = 2.1+0.4
−0.5 solar mass
at a 95% confidence level [99].
c) Neutron star cooling : Furthermore, it predicts a
direct URCA cooling for neutron stars above 1.4 times
the solar mass. If such is the case, then the enhanced
cooling of a M = 1.4 solar mass neutron star may provide
strong evidence in favor of exotic matter in the core of a
neutron star.
These results have important implications for nuclear
astrophysics and future experiments probing the prop-
erties of nuclei using beams of neutron-rich nuclei. The
above constraint was obtained by studying the low den-
sity behavior of nuclear matter. Measurements at den-
sities above normal nuclear matter should further con-
straint the form of the symmetry energy. Such measure-
ments should yield consistent results when extrapolated
to low densities.
VIII. SUMMARY AND CONCLUSIONS
In summary, a number of reactions were studied to in-
vestigate the density dependence of the symmetry energy
in the equation of state of isospin asymmetric nuclear
matter. The results were analyzed within the framework
of the dynamical and the statistical models of multifrag-
mentation. It is observed that a dependence of the form
Csym(ρ) = 31.6 (ρ/ρo)
0.69 MeV, agrees reasonably with
the experimental data indicating that a “ stiff ” form of
the symmetry energy provides a better description of the
nuclear matter EOS at sub-nuclear densities. A compar-
ison with several other independent studies shows that
a constraint on the density dependence of the symmetry
energy given as Csym(ρ) = C
sym(ρ/ρo)
γ , where Cosym ≈
31 - 33 MeV and γ ≈ 0.55 - 0.69, can thereby be obtained.
The present observation has important implications for
astrophysics, as well as, nuclear physics studies to be car-
ried out at future radioactive beam facilities worldwide.
IX. ACKNOWLEDGMENTS
This work was supported in part by the Robert A.
Welch Foundation through grant No. A-1266, and the
Department of Energy through grant No. DE-FG03-
93ER40773. We also thank A. Botvina for the SMM
code, and the Catania group for the BNV code.
[1] A.E.L. Dieperink, Y. Dewulf, D. Van Neck, M. Waro-
quier, and V. Rodin, Phys. Rev. C 68, 064307 (2003).
[2] R.B. Wiringa, V. Fiks and A. Fabrocini, Phys. Rev. C
38, 1010 (1988).
[3] C.H. Lee, T.T.S. Kuo, G.Q. Li and G.E. Brown, Phys.
Rev. C 57, 3488 (1998).
[4] B. Liu, V. Greco, V. Baran, M. Colonna, and M. DiToro,
Phys. Rev. C 65, 045201 (2002).
[5] N. Kaiser, S. Fritsch and W. Weise, Nucl. Phys. A 697,
255 (2002).
[6] C. Fuchs and H.H. Wolter, Eur. Phys. J. A 30, 5 (2006).
[7] B.A. Brown, Phys. Rev. Lett 85, 5296 (2000).
[8] C.J. Horowitz and J. Piekarewicz, Phys. Rev. Lett. 86,
5647 (2001).
[9] R.J. Furnstahl, Nucl. Phys. A 706, 85 (2002).
[10] K. Oyamatsu, I. Tanihata, Y. Sugahara, K. Sumiyoshi,
and H. Toki, Nucl. Phys. A 634, 3 (1998).
[11] J. Lattimer, C. Pethick, M. Prakash and P. Hansel, Phys.
Rev. Lett. 66, 2701 (1991).
[12] C. Lee, Phys. Rep. 275, 255 (1996).
[13] C.J. Pethick and D.G. Ravenhall, Annu. Rev. Nucl. Part.
Sci. 45, 429 (1995).
[14] J.M. Lattimer and M. Prakash, Phys. Rep. 333, 121
(2000).
[15] W.R. Hix, O.E.B. Messer, A. Mezzacappa, M. Lieben-
dorfer, J. Sampaio, K. Langanke, D.J. Dean, and G.
Martinez-Pinedo, Phys. Rev. Lett. 91, 201102 (2003).
[16] J.M. Lattimer and M. Prakash, Science, 304, 536 (2004).
[17] A.W. Steiner, M. Prakash, J.M. Lattimer and P.J. Ellis,
Phys. Rep. 411, 325 (2005).
[18] C.J. Horowitz, S.J. Pollock, P.A. Souder and R. Michaels,
Phys. Rev. C 63, 025501 (2001).
[19] C.J. Horowitz and J. Piekarewicz, Phys. Rev. C 66,
55803 (2002).
[20] J.R. Stone, J.C. Miller, R. Koncewicz, P.D. Stevenson
and M.R. Strayer, Phys. Rev. C 68, 034324 (2003).
[21] J. Lattimer et al., Astrophys. J. 425, 802 (1994).
[22] P. Slane, D.J. Helfand and S.S. Murray, Astrophys. J. L
571, 45 (2002).
[23] P. Danielewicz, R. Lacey and W.G. Lynch, Science, 298,
1592 (2002).
[24] GSI Conceptual Design Report,
http://www.gsi.de/GSI-Future.
[25] RIA homepage, http://www.orau.org/ria.
[26] B.A. Li, C.M. Ko and W. Bauer, Int. J. Mod. Phys. E 7,
147 (1998).
[27] V. Baran, M. Colonna, V. Greco and M. DiToro, Phys.
Rep. 410, 335 (2005).
[28] M.B. Tsang, W.A. Friedman, C.K. Gelbke, W.G. Lynch,
G. Verde, and H.S. Xu, Phys. Rev. Lett. 86, 5023 (2001).
[29] H.S. Xu et al., Phys. Rev. Lett. 85, 716 (2000).
[30] D.V. Shetty, S.J. Yennello, E. Martin, A. Keksis, G.A.
Souliotis, Phys. Rev. C 68, 021602 (2003).
[31] A. Ono, P. Danielewicz, W.A. Friedman, W.G. Lynch
and M.B. Tsang, Phys. Rev. C 68, 051601 (2003).
[32] Q. Li, Z. Li, H. Stocker, Phys. Rev. C 73, 051601 (2006).
[33] T.X. Liu et al., Phys. Rev. C 69, 014603 (2004).
[34] M. Colonna and F. Matera, Phys. Rev. C 71, 064605
(2005).
[35] C.O. Dorso, C.R. Escudero, M. Ison, and J.A. Lopez,
Phys. Rev. C 73, 044601 (2006).
[36] A.S. Botvina, O.V. Lozhkin and W. Trautmann, Phys.
Rev. C 65, 044610 (2002).
[37] Al.H. Raduta, Ad.R. Raduta, Phys. Rev. C 65, 054610
(2002).
[38] W.A. Friedman, Phys. Rev. C 42, 667 (1990).
[39] D.H.E. Gross, Phys. Rep. 279, 119 (1997).
[40] M.B. Tsang et al., Phys. Rev. C 64, 054615 (2001).
[41] W.D. Tian et al., arXiv: nucl-th/0601079 (2006).
[42] A. Ono, P. Danielewicz, W.A. Friedman, W.G. Lynch
and M.B. Tsang, arXiv: nucl-ex/0507018 (2005).
[43] A. Ono, P. Danielewicz, W.A. Friedman et al., Phys. Rev.
C 70, 041604 (2004).
[44] H. Johnston et al., Phys. Rev. C 56, 1972 (1997).
[45] S.J. Yennello et al., Phys. Lett. B 321, 15 (1994).
[46] R.P. Schmitt et al., Nucl. Instrum. Methods Phys. Res.
A 354, 487 (1995).
[47] C.A. Ogilvie, D.A. Cebra, J. Clayton, S. Howden, J.
Karn, A. Vander Molen, G.D. Westfall, W.K. Wilson,
and J.S. Winfield, Phys. Rev. C 40, 654 (1989).
[48] S.J. Yennello et al., Proc. of the 10th Winter Workshop
on Nuclear Dynamics, Snowbird, UT, edited byW. Bauer
http://www.gsi.de/GSI-Future
http://www.orau.org/ria
http://arxiv.org/abs/nucl-th/0601079
http://arxiv.org/abs/nucl-ex/0507018
(1994).
[49] E. Geraci et al., Nucl. Phys. A 732, 173 (2004).
[50] L.W. Chen, C.M. Ko and B.A. Li, Phys. Rev. Lett. 94,
032701 (2005).
[51] J.P. Bondorf et al., Phys. Rep. 257, 133 (1995).
[52] A.S. Botvina et al., Nucl. Phys. A 584, 737 (1995).
[53] A.S. Botvina, A.S. Iljinov, I.N. Mishustin, J.P. Bondorf,
R. Donangelo, and K. Sneppen, Nucl. Phys. A 475, 663
(1987).
[54] V. Baran, M. Colonna, M. Di Toro, V. Greco, M.
Zielinska-Pfabe and H.H. Wolter, Nucl. Phys. A 703, 603
(2002).
[55] J.B. Natowitz, R. Wada, K. Hagel, T. Keutgen, M. Mur-
ray, A. Makeev, L. Qin, P. Smith, and C. Hamilton, Phys.
Rev. C 65, 034618 (2002); and references therein.
[56] A.S. Botvina and I.N. Mishustin, Phys. Rev. C 63,
061601 (2001).
[57] J.P. Bondorf, R. Donangelo, I.N. Mishustin and H.
Schultz, Nucl. Phys. A 444, 460 (1985).
[58] J.P. Bondorf, A.S. Botvina and I.N. Mishustin, Phys.
Rev. C 58, 27 (1998).
[59] W.P. Tan et al., Phys. Rev. C 64, 051901 (2001).
[60] J. Iglio, D.V. Shetty, S.J. Yennello, G.A. Souliotis, M.
Jandel, A. Keksis, S. Soisson, B. Stein and S. Wuenschel,
Phys. Rev. C 74, 024605 (2006).
[61] D.V. Shetty, A.S. Botvina, S.J. Yennello, A. Keksis, E.
Martin and G.A. Souliotis, Phys. Rev. C 71, 024602
(2005).
[62] A. LeFevre et al., Phys. Rev. Lett. 94, 162701 (2005).
[63] D. Henzlova, A.S. Botvina, K.H. Schmidt, V. Henzl, P.
Napolitani, and M.V. Ricciardi, arXiv: nucl-ex/0507003
(2005).
[64] N. Buyukcizmeci, R. Ogul, and A.S. Botvina, Eur. Phys.
J. 25, (2005) 57.
[65] L.G. Sobotka, R.J. Charity, J. Toke and W.U. Schroder,
Phys. Rev. Lett. 93, 132702 (2004).
[66] R. Hasse and P. Schuck, Phys. Lett. B 179, 313 (1986).
[67] S. Shlomo and J.B. Natowitz, Phys. Lett. B 252, 187
(1990).
[68] S. Shlomo and J.B. Natowitz, Phys. Rev. C 44, 2878
(1991).
[69] J.B. Natowitz, K. Hagel, Y. Ma, M. Murray, L. Qin, S.
Shlomo, R. Wada, and J. Wang, Phys. Rev. C 66, 031601
(2002).
[70] V.E. Viola, K. Kwiatkowski, J.B. Natowitz, and S.J. Yen-
nello, Phys. Rev. Lett. 93, 132701 (2004).
[71] Isospin Physics in Heavy Ion Collisions at Intermediate
Energies, edited by B.-A. Li and W. Schroder (Nova Sci-
ence, New York, 2001).
[72] B.A. Li and L.W. Chen, Phys. Rev. C 74, 034610 (2006).
[73] D.V. Shetty, S.J. Yennello, A.S. Botvina, G.A. Souliotis,
M. Jandel, E. Bell, A. Keksis, S. Soisson, B. Stein, and
J. Iglio , Phys. Rev. C 70, 011601 (2004).
[74] D.V. Shetty, S.J. Yennello, and G.A. Souliotis, Phys.
Rev. C 75, 034602 (2007); arXiv: nucl-ex/0505011
(2005).
[75] D.V. Shetty, S.J. Yennello, G.A. Souliotis, A.L. Keksis,
S.N. Soisson, B.C. Stein, and S. Wuenschel, Phys. Rev.
C (2007) (Submitted); arXiV: nucl-ex/0606032 (2006).
[76] B.G. Todd-Rutel and J. Piekarewicz, Phys. Rev. Lett 95,
122501 (2005).
[77] J. Piekarewicz (Private Communication).
[78] J. Piekarewicz, Proc. of the International Conference on
Current Problems in Nuclear Physics and Atomic En-
ergy, Kyiv, Ukraine, (May 29 - June 3, 2006).
[79] M.B. Tsang et al., Phys. Rev. Lett. 92, 062701 (2004).
[80] B.A. Li and L.W. Chen, Phys. Rev. C 72, 064611 (2005).
[81] P. Danielewicz, arXiv : nucl-th/0411115 (2004).
[82] H. Heiselberg and M. Hjorth-Jensen, Phys. Rep. 328, 237
(2000).
[83] A. Akmal and V.R. Pandharipande, Phys. Rev. C 56,
2261 (1997).
[84] A. Akmal, V.R. Pandharipande and D.G. Ravenhall,
Phys. Rev. C 58, 1804 (1998).
[85] C. Fuchs (Private Communication); E.N.E. van Dalen, C.
Fuchs, and A. Faessler, Nucl. Phys. A 744, 227 (2004).
[86] D.T. Khoa and H.S. Than, Phys. Rev. C 71, 044601
(2005).
[87] M.A. Famiano et al., Phys. Rev. Lett. 97, 052701 (2006).
[88] B.A. Li, C. Ko, and Z. Ren, Phys. Rev. Lett. 78, 1644
(1997).
[89] B.A. Li, L.W. Chen, G.C. Yong, and W. Zuo, Phys. Lett.
B 634, 378 (2006).
[90] S. Hudan et al., Phys. Rev. C 67, 064613 (2003).
[91] M. Colonna and M.B. Tsang, Eur. Phys. J. A 30, 165
(2006).
[92] P. Danielewicz, Nucl. Phys. A727, 233 (2003).
[93] W. Zuo, I. Bombaci, and U. Lombardo, Phys. Rev. C 60,
024605 (1999).
[94] M. Brack, C. Guet, and H.B. Hakansson, Phys. Rep. 123,
276 (1985).
[95] J.M. Pearson and R.C. Nayak, Nucl. Phys. A 668, 163
(2000).
[96] W.D. Myers and W.J. Swiatecki, Nucl. Phys. 81, 1
(1966).
[97] K. Pomorski and J. Dudek, Phys. Rev. C 67, 044316
(2003).
[98] C.J. Horowitz, Eur. Phys. J. A 30, 303 (2006).
[99] D.J. Nice et al., Astrophys. J. 634, 1242 (2005).
http://arxiv.org/abs/nucl-ex/0507003
http://arxiv.org/abs/nucl-ex/0505011
http://arxiv.org/abs/nucl-ex/0606032
http://arxiv.org/abs/nucl-th/0411115
