Quantum molecular dynamics simulation formultifragmentation resulting from an expandingnuclear matter

doi:10.4236/ns.2011.37082

Paper Menu >>

Journal Menu >>

Vol.3, No.7, 594-599 (2011) Natural Science

http://dx.doi.org/10.4236/ns.2011.37082

Quantum molecular dynamics simulation for

multifragmentation resulting from an expanding

nuclear matter

Atef Abdel-Hafiez1*, Mohamed. Abdo Khalifa2, Albiomy Abd El-Daiem3

1Department of Experimental Nuclear Physics, Nuclear Research Center, Atomic Energy Authority, Cairo, Egypt;

*Corresponding Author: [email protected]

2Mathematics Department, Faculty of Science, Tanta University, Egypt.

3Physics Department, Faculty of Science, Sohag University, Sohag, Egypt.

Received 20 December, 2010; revised 20 April, 2011; accepted 3 May, 2011.

ABSTRACT

Quantum molecular dynamics (QMD) is used to

investigate multifragmentation resulting from an

expanding nuclear matter. Equation of state, the

structure of nuclear matter and symmetric nu-

clear matter is discussed. Also, the dependence

of the fragment mass distribution on the initial

temperature (Tinit) and the radial flow velocity (h)

is studied. When h is large, the distribution

shows exponential shape, whereas for small h,

it obeys the exponentially falling distribution

with mass number. The cluster formation in an

expanding system is found to be different from

the one in a thermally equilibrated system. The

used Hamiltonian has a classical kinetic energy

term and an effective potential term composed

of four parts.

Keywords: Quantum Molecular Dynamics;

Multifragmentation

1. INTRODUCTION

Multifragmentation has been a long-standing topic in

nuclear physics. In particular, a fragment mass distribu-

tion has been extensively investigated both theoretically

and experimentally as a phenomenon which might have

some connection with nuclear phase transition [1-6].

Among the theoretical approaches, the statistical model

[7] is a familiar one which predicts the fragment mass

distribution.

According to this model, the exponentially falling dis-

tribution with mass number is derived under the condi-

tion that each fragment can be regarded as a thermally

equilibrated droplet. This exponentially falling distribu-

tion with mass number behavior has attracted many au-

thors’ interest because the exponentially falling distribu-

tion with mass number is a signature of a second order

phase transition [8-13].

In order to investigate the relation between the cluster

formation and the fragment mass distribution, many

heavy ion collision experiments have been conducted

[3-6]. In these experiments, two types of fragment mass

distribution have been observed. The exponentially fal-

ling distribution with mass number is one of them and is

originating from a spectator region where thermal equi-

librium is often assumed [3,4]. On the other hand, in a

participant region, the source of fragments undergoes

strong compression and each nucleon has large collec-

tive momentum at its decompression stage. Thus, we

cannot regard this source of fragments as thermal equili-

brated object that the statistical model assumes. In this

region, the mass distribution does not obey the power

law but obeys an exponential distribution and the shape

of which strongly depends on the collision energy [5,6].

There are several factors to prevent the direct under-

standing of cluster formation in heavy ion collision ex-

periments. The finiteness of heavy ion collision is a

source of complication. For instance, the surface effect

during fragment formation might affect the distribution

significantly. Moreover, it is doubtful whether thermal

equilibrium is achieved even in a spectator region. To

check the basic mechanism of fragmentation, therefore,

it is desirable to adopt a model which does not assume a

thermal equilibrium and is free from the complication

arising from the finiteness of the system.

The analysis of nuclear matter is useful to avoid the

complication coming from finiteness. Instability of nu-

clear matter against multifragmentation has been studied

based on the Boltzmann-Langevin approach [14], clas-

sical molecular dynamics [11], kinetic equation [1], and

linear response theory [15]. These are simulations of an

infinite system studying under what condition the insta-

A. Abdel-Hafiez et al. / Natural Science 3 (2011) 594-599

595

bility against the multifragmentation occurs at high ex-

citation. These studies, however, have not taken into

account the collective expansion of the system explicitly.

In [16], the effect of collective expansion has been

treated in the framework of Boltzmann-Langevin ap-

proach and the fragmentation is discussed for a finite

system. The instability to multifragmentation is also

treated in the quantal RPA approach [17].

The aim of this work is to investigate multifragmenta-

tion by using an expanding nuclear matter. When the

matter expands slowly, it simulates a spectator region

where each nucleon moves slowly. On the other hand,

the rapid expansion simulates a participant region where

the system has a rapid collective motion. The rate of

expansion is controlled by only one parameter (h). By

changing h we can investigate both spectator and par-

ticipant regions. In the case of h = 0, the system is just in

a static thermal equilibrium state with temperature (T).

Such a system is similar to what is assumed in the statis-

tical model. Our main purpose is to know how the frag-

ment mass distribution depends on the expansion veloc-

ity and the temperature.

The used quantum molecular dynamics (QMD) and

the nuclear matter at the saturation density



0 and a

temperature Tinit.

On the basis of QMD, a series of simulations is car-

ried out. There are four steps in the calculation.

1) The nuclear matter at the saturation density



0 and

a temperature Tinit is prepared by Metropolis sampling

method.

2) To this thermal matter with Tinit , a radial flow ve-

locity (h) is imposed.

3) The time evolution of the system is calculated from



0 to 0.001



0 by QMD.

4) A fragment mass distribution is calculated at

0.001



The Hamiltonian has a classical kinetic energy term

and an effective potential term composed of four parts:

nuclsurface Pauli coulmb

HVVVV





 (1)



2(1 )

nucl ii

ij ijij

ij iij iij

ij iij

Ccc





























 

























(2)

 

5/3

surface

ij i

Vrrr









 (3)



,i,j

.00

exp

pauli

Aij ij

ij i

qp c

















 









pp pp(4)

 

e11

222

coulombi j

ij i







 



 

 









 rrr r

(5)

where



i is the nucleon spin and



i is the isospin [



i =1/2

for protons and –1/2 for neutrons], 



i is an overlap of

one body density



i for ith nucleon with other nucleons

defined as

 

iij

ji ji

 







rr (6)

and the explicit form of



i is given by

 



exp 2

2π











rrR

(7)

where L is the squared width of the wave packet, L =

1.95 fm2. This wave packet with fixed width is charac-

teristic of the QMD model. It has the advantage of deal-

ing with the many body correlation directly but imposes

considerable restrictions on the instability of the nuclear

matter.

We approximate the fermionic nature of the system by

introducing the Pauli potential between two identical

nucleons [18,19], instead of antisymmetrization. The

nuclear potential Vnucl has several components: The

first two terms are Skyrme type nuclear interactions. The

term with Cs0 is a symmetry potential where Ci is 1 for

protons and 0 for neutrons. The last two terms with Cex(1)

and Cex(2) are momentum dependent potentials originat-

ing from the exchange term of the Yukawa interaction.

The term Vsurface stands for a surface interaction to de-

scribe the fragment property. The values of the parame-

ters in the interaction are fitted to give good agreements

with the nuclear matter saturation properties, the real

part of optical potential, the rms radii of nuclei, and their

binding energies.

The incompressibility constant K at the saturation

point is simulated to the following parabolic form











Equation of state and the structure of nuclear matter.

A. Abdel-Hafiez et al. / Natural Science 3 (2011) 594-599

596

In this section we study properties of nuclear matter at

several conditions. It is desirable to use a cell large

enough to include several periods of structure and to

avoid the spurious effects of the boundary condition on

the structure of matter. Though our calculation with typi-

cally 1024 particles in a cell is not fully satisfactory in

this respect, we consider it is enough for semi qualitative

discussions at the beginning of this study. Actually, the

global quantities, e.g., the ground state energy of the

system, are well saturated at this number of particles in a

cell.

In Figure 1, we show the energy per nucleon of the

infinite system as a function of the particle number (NA)

in a cell for four average densities from



= 0.4



0 to 2.2



0. For all densities, the energy of the system has al-

ready approached the asymptotic value above 256 parti-

cles within 100 KeV/ nucleon. In this study we simulate

an infinite system by the periodic boundary condition

mainly with 1024 or 2048 particles in a cell, and inves-

tigate the ground state properties of nuclear matter.

Figure 2 shows the energy per nucleon as a function

of the average density. The solid squares indicate the

energy of ‘‘uniform’’ nuclear matter while open squares

are the results of energy-minimum configurations.

The ‘‘uniform’’ matter energy is calculated as follows:

First we distribute nucleons randomly and cool the sys-

tem only with the Pauli potential. The Pauli potential is

repulsive and does not spoil the uniformity of the system.

Then we impose the other effective interactions and cool

only in the momentum space, fixing the positions of

particles. The system turns out to be approximately uni-

form with this procedure. Note that simulated ‘‘uniform’’

matter is not exactly the same as ideal nuclear matter

since the latter is continuous and completely uniform.

Both cases of uniform and energy-minimum configura-

tions have almost the same energy per nucleon for the

higher densities as is seen in this figure. Below the satu-

Figure 1. Particle number (NA) dependence of the en-

ergy per nucleon (E/A) of the infinite system for four

average densities from



= 0.4



0 to 2.2



Figure 2. The energy per nucleon of symmetric nuclear

matter [E/A = 0.5] at zero temperature as a function of the

average densities.

ration density



0, the energy per nucleon of the energy-

minimum configuration is lower than the uniform case.

The deviation amounts to about 5 MeV. As we see in the

following, this is due to the structure change of matter

from uniform to nonuniform such as a clusterized one.

From the left, the open squares are the results with

soft [K = 210 MeV] and hard EOS [K = 380 MeV] ob-

tained by the damping equation of motion searching the

energy-minimum configuration in the full phase space.

The solid squares indicate results obtained from the spa-

tially uniform distribution. The kinetic energy of the

electron is not included. We use 1024 particles in a cell

for all cases.

If a system is finite and has high temperature and/or

pressure, it expands spontaneously. However, an infinite

matter cannot expand without an initial collective mo-

tion even if it has high temperature and pressure. There-

fore, we give a collective momentum to each constituent

nucleon so that the matter could expand homogeneously.

This collective momentum is expressed in terms of one

parameter h as

coll i

ph p





R (8)

where i

R is the position of i th particle [the center of

wave packet] measured from the center of the primitive

cell and pF is the Fermi momentum at



0. The distance is

normalized by using



–1/3. The nearest neighbor pair of

nucleons, on average, have a relative momentum hpF

under this condition. Time evolution from



0 to 0.001



by QMD.

The time evolution of the nuclear matter with a given

temperature Tinit and radial flow velocity with [h,Tinit]

by QMD has been followed by Shinpei Chikazumi, et al.

[20]. Figure 3(a) shows a snapshot of an initial state at



0. Figure 3(b) shows a part of the whole system at



= 0.05



0 during the expansion. The expansion is stopped

when the average density reaches



= 0.001



0 shown in

Figure 3(c). In this procedure, the normal periodic

A. Abdel-Hafiez et al. / Natural Science 3 (2011) 594-599

597

Table 1. Effective interaction parameter set.

K = 210 MeV[soft] K = 380 MeV[hard]



[MeV] –121.9 –21.21



[MeV] 197.3 97.93

 1.33333 1.66667

Cs0 [MeV] 25.0 25.0

Cex(1) [MeV] –258.5 –258.54

Cex(2) [MeV] 375.6 375.6



1 [MeV] 2.35 2.35



2 [MeV] 0.4 0.4

VSF [MeV] 20.68 20.68

CP [MeV] 115.0 115.0

p0 [MeV] 120.0 120.0

q0 [fm] 2.5 2.5

L [fm2] 1.95 2.05

(a)



= 1.000



0 (b)



= 0.05



0 (c)



= 0.001



Figure 3. Snapshots of expanding nuclear matter at different

densities. (a) An initial state with



= 1.0



0, (b) An intermedi-

ate state at



= 0.05



0 during the expansion. (c).

boundary condition is not applicable to the expanding

matter. Therefore, we introduce an extended periodic

boundary condition.

2. FRAGMENT MASS DISTRIBUTION

The matter leave to expand until the average density

becomes small enough so that each fragment can be

identified. All the fragments are isolated when the aver-

age distance between nearest neighbor pair reaches

about 5 fm at 0.05



0. However, we noticed that the in-

trinsic expansion of fragments is not yet stopped at



0.05



0 especially when h is large. Therefore, we set the

final density equal to 0.001



0. At this density, all the

fragments are stabilized so that we can identify them

only by their positions.

We show our main result of the fragment mass dis-

tribution. The distribution strongly depends on the radial

flow velocity h whereas the initial temperature has little

effect on the distribution. First, we concentrate on the

distribution for Tinit = 30 MeV and investigate how the

radial flow affects the distribution.

Figure 4 shows the distributions calculated for h = 0.1,

0.5, 1.0, 2.0 with the initial temperature Tinit = 30 MeV.

In this figure, all the distributions follow straight lines

except for h = 0.1. These straight lines in the semi loga-

rithmic scale mean that the number of large fragments

decreases exponentially as fragment size increases. This

exponential shape can be understood as the manifesta-

tion of random distribution[21]. When h is large enough

to suppress the interaction, the final distribution is sim-

ply an enlarged copy of the initial configuration. Under

this limit, the particles are randomly distributed irrespec-

tive of the interaction because the initial configuration is

prepared at saturation density and the particles are ran-

domly distributed under any temperature. The value of h

directly corresponds to the radial flow velocity in the

participant region. Actually, this kind of heavy ion colli-

sion is investigated by a simulation of a finite system

with a radial flow velocity [22]. Figure 5 shows the

same quantity as Figure 4 in double-logarithmic scale.

In our method, we can discuss most clearly the effect of

radial flow because our system is an infinite expanding

nuclear matter which is free from any finite size effect of

the parent nucleus.

Figure 4. The fragment mass distribution obtained in our

QMD simulation in case of semi logarithmic scale.

Figure 5. The fragment mass distribution obtained in our

QMD simulation in case of double logarithmic scale.

A. Abdel-Hafiez et al. / Natural Science 3 (2011) 594-599

598

When h is small, the interaction must play an impor-

tant role in fragmentation. In this case, h = 0.1 follows a

straight line whereas the others do not. The straight line

in double-logarithmic scale refers to the so-called expo-

nentially falling distribution with mass number. The ex-

ponentially falling distribution with mass number has

been discussed in many works in connection with critical

phenomena. Usually, the exponentially falling distribu-

tion with mass number is thought to be a manifestation

that the system undergoes a second order phase transi-

tion. In the model assuming thermal equilibrium, the

power law appears only when the second order phase

transition happens. However, we must remember that the

expanding system is completely in nonequilibrium.

Even when h = 0.1, the system does not have enough

time to reach thermal equilibrium where the temperature

can be determined. Therefore this exponentially falling

distribution with mass number we obtained must be ex-

plained in a different way.

3. EXPANDING FRAGMENTATION VS.

ISOTHERMAL FRAGMENTATION

It is significant to compare two types of the power law,

i.e., dynamical and isothermal ones. The isothermal frag-

mentation can be investigated by Metropolis sampling

method. Like creating initial configuration, we prepare

1000 samples with a fixed temperature T at



= 0.05



instead of



0. From investigation of isothermal pres-

sure, we have already known that the critical tempera-

ture is around 8 MeV. Figure 6 shows how fragment

mass distribution depends on the temperature T = 5, 8,

18 MeV at



= 0.05



0. For T = 5 MeV, the distribution

shows the U-shape where large fragments and small

fragments exist but there is no middle size fragments. T

= 8 MeV is just a critical temperature in which the

power law appears. The exponent



, i.e., the slope of

distribution, is around 2.5 which is consistent with the

Fisher’s droplet model. As the temperature increases

Figure 6. The comparison between isothermal fragmentation and expanding fragmentation.

The three figures [[i], [ii], and [iii]] are isothermal distributions [t = 5, 8 and 18 MeV] at



0.05



0. The two figures [[iv] and [v]] are the distributions obtained by the expansion [h =

0.10].

A. Abdel-Hafiez et al. / Natural Science 3 (2011) 594-599

599

beyond the critical temperature Tc = 8 MeV, the distribu-

tion deviates from the exponentially falling distribution

with mass number and becomes an exponential shape. In

the bottom figure [T = 18 MeV], thermal motion com-

pletely overcomes attractive interaction. Therefore, a

random distribution appears like the rapid expansion

case.

4. CONCLUTIONS

In this paper we have presented the QMD approach

for studying multifragmentation resulting from an ex-

panding nuclear matter. We can reproduce well the finite

nuclear properties for various mass ranges by inclusion

of the Pauli and surface potential. We have investigated

the EOS of nuclear matter by simulating an infinite sys-

tem with the used QMD. The fragmentation during ex-

pansion can be classified according to the speed of ex-

pansion h. Also, the symmetric nuclear matter is dis-

cussed.

Our QMD model contains a further possibility for the

simulation of the dynamical evolution of infinite nuclear

matter such as supernova explosions, the glitch of the

neutron star, and the initial stage of the universe. An

intensive and systematic study of nuclear matter with the

present model will be important since it contains fewer

assumptions than the foregoing models as to the struc-

ture of matter.

REFERENCES

[1] Schmelzer, J., Ro¨pke, G. and Ludwig, F.-P. (1997) Nu-

clear multifragmentation processes and nucleation theory.

Physical Review C, 55, 1917-1927.

doi:10.1103/PhysRevC.55.1917

[2] Strachan, A. and Dorso, C.O. (1999) Temperature and

energy partition in fragmentation. Physical Review C, 59,

285-294. doi:10.1103/PhysRevC.59.285

[3] Gilkes, M.L. et al., (1994) Determination of critical

exponents from the multifragmentation of gold nuclei.

Physical Review Letters, 73, 1590-1593.

doi:10.1103/PhysRevLett.73.1590

[4] Mastinu, P.F. et al., (1996) Circumstantial evidence for

critical behavior in peripheral Au+Au collisions at 35

MeV/nucleon. Physical Review Letters, 76, 2646-2649.

doi:/10.1103/PhysRevLett.76.2646

[5] Petrovici, M. et al., (1995) Cluster formation during

expansion of hot and compressed nuclear matter pro-

duced in central collisions of Au on Au at 250 A MeV.

Physical Review Letter ,74, 5001-5004.

http://dx.doi.org/10.1103/PhysRevLett.76.2646

[6] Reisdorf, W. et al., (1997) Central collisions of Au on Au

at 150-A/Mev, 250-A/Mev and 400-A/Mev, Nuclear

Physics A, 612, 493-556.

doi:10.1016/S0375-9474(96)00388-0

[7] Das Gupta, S., Mekjian, A.Z. and Tsang, B. (2001)

Symmetry and surface symmetry energies in finite nuclei.

Advances in Nuclear Physics, 26, 89.

[8] Mishustin, I.N. (1998) Fluctuations near the deconfine-

ment phase transition boundary. Nuclear Physics A, 630,

111-125. doi:10.1016/S0375-9474(97)00748-3

[9] Belkacem, M., Latora, V. and Bonasera, A. (1995) Criti-

cal evolution of a finite system. Physical Review C, 52,

271-285. doi:10.1103/PhysRevC.52.271

[10] Colonna, M. and Chomaz, Ph. (1998) Quantal effects on

growth of instabilities in nuclear matter. Physical Review

B, 436, 1-9. doi:10.1016/S0370-2693(98)00882-X

[11] Finocchiaro, P., Belkacem, M., Kubo, T., Latora, V. and

Bonasera, A. (1996) The nuclear liquid-gas phase transi-

tion within fermionic molecular dynamics. Nuclear

Physics A, 600, 236-250.

doi:10.1016/0375-9474(96)00040-1

[12] Ohnishi, A. and Randrup, J. (1997) Quantum fluctuation

effects on nuclear fragment formation. Physical Review B,

394, 260-268. doi:10.1016/S0370-2693(97)00032-4

[13] Sugawa, Y. and Horiuchi, H. (1999) Nucleon-induced

fragment formation with antisymmetrized molecular dy-

namics. Physical Review C, 60, 064613-054621.

doi:10.1103/PhysRevC.60.064613

[14] Colonna, M. andChomaz, Ph. (1994) Unstable infinite

nuclear matter in stochastic mean field approach. Physi-

cal Review C, 49, 1908.

[15] Baran, V., Colonna, M., Di Toro, M. and Larionov, A.B.

(1998) Fragment isotope distributions and the isospin

dependent equation. Nuclear Physics A, 632, 287-303.

doi:10.1016/S0375-9474(98)00006-2

[16] Guarnera, A., Colonna, M. and Chomaz, Ph. (1996) Iso-

spin effects in nuclear fragmentation. Physics Letters B,

373, 267-274. doi:10.1016/0370-2693(96)00152-9

[17] Ayik, S. and Colonna, M. (1995) Quantal effects on

growth of instabilities in nuclear matter. Physical Review

B, 353, 417-421. doi:10.1016/0370-2693(95)00596-D

[18] Dorso, C.O., Duarte, S. and Randrup, J. (1987) Classical

simulation of the fermi gas. Physics Letters B, 188, 287-294.

doi:10.1016/0370-2693(87)91382-7

[19] Maruyama, T. Niita, K., Oyamatsu, K. Maruyama, T.,

Chiba, S. and Iwamoto, A. (1998) Quantum molecular

dynamics approach to the nuclear matter below the satu-

ration density. Physical Review C, 57, 655-665.

doi:10.1103/PhysRevC.57.655

[20] Chikazumi, S., Maruyama, T., Chiba, S., Niita, K. and

Iwamoto, A. (2001) Quantum molecular dynamics simu-

lation of expanding nuclear matter and nuclear multi-

fragmentation. Physical Review C, 63, 024602.

doi:10.1103/PhysRevC.63.024602

[21] Holian B.L. and Grady, D.E. (1988) Fragmentation by

molecular dynamics: The microscopic ‘‘big bang”. Phy-

sical Review Letters, 60, 1355-1358.

doi:10.1103/PhysRevLett.60.1355

[22] Reisdorf, W.et al., (1997) Central collisions of Au on Au

at 150-A/Mev, 250-A/Mev and 400-A/Mev. Nuclear

Physics A, 612, 493-556.

doi:10.1016/S0375-9474(96)00388-0