Time-dependent Hartree-Fock calculation of 12C + 12C with a realistic potential

I

2.N

[

Nuclear Physics A270 (1976) 471-7488; (~) North-Holland Publishing Co., Amsterdam

I

Not to be reproduced by photoprint or microfilm without written permission from the publisher

TIME-DEPENDENT HARTREE-FOCK CALCULATION OF 12C+ t2C WITH A REALISTIC POTENTIAL * J. A. M A R U H N tt Physics Division, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37830 and

g. Y. CUSSON t,, Physics Department, Duke University, Durham, North Carolina 27706 Received 14 June 1976 Abstract: We have used a realistic single-particle K-matrix model to compute the head-on scattering of 12C+t2C at incident projectile lab energies o f 3.2, 6.4, 12.8, 19.2, 25.6, 32, 51.2 and 64 MeV/nucloon, above the Coulomb barrier, in the time-dependent Hartree-Fock approximation. Direct and exchange Coulomb forces as well as spin-orbit forces are included. A large deformed harmonic oscillator basis is used. Spatial density and current distributions at various times are shown. The outgoing energy is found to be Eo = 0.8Et~--28 (MeV), in the c.m. system. Fusion and fully relaxed scattering axe observed at low energy. Some compression is seen at higher energies but no shock waves can be detected. Consequences for heavy-ion reactions are discussed.

1. Introduction The dynamics of heavy-ion collisions is currently a topic of considerable interest to experimental and theoretical nuclear physicists alike. Experimental effects such as deep inelastic scattering i), fully relaxed collisions and characteristic plots of outgoing kinetic energy versus scattering angles 2) strongly suggest the occurrence of physical processes unlike anything else previously studied. It is to be expected that the understanding of these experimental effects will require parallel new developments in the theoretical methods for handling macroscopic but quantal dynamical effects in colliding nuclear fluid droplets. Since the static nuclear fluid droplet concept is not unlike the Thomas-Fermi model of the atom one is soon led to consider suggestions such as Dirac's 1930 proposal a) that one use the time-dependent Hartree-Fock (TDHF) approximation to study the dynamics of a moving Fermi gas, represented by a single time-dependent Slater determinant. The problem of computing the time evolution in three dimensions, using realistic nuclear interactions or reaction matrices remains however a formidable computational task still out of reach of present facilities. One must therefore develop the program in stages. The first step is then to consider one-dimensional motion using t Work supported in part by E R D A . *t Research sponsored by the US Energy Research and Development Administration under Contract with Union Carbide Corporation. ,tt Consultant at Oak Ridge National Laboratory. 471

472

J . A . M A R U H N A N D R. Y. C U S S O N

simplified nuclear interactions. This has been done recently 4) with the conclusion that the T D H F method does indeed predict, at least in one dimension, the occurrence of highly inelastic scattering, capture, and large amplitude collective motion. Using a simplified Skyrme interaction it has also recently been shown, in the case of 160 + 160, that one can observe exciting effects when solving the full three-dimensional scattering problem 5). The restriction to an oversimplified interaction is lifted in the present paper, but at the cost of considering head-on collisions only, which reduces the problem to a manageable two-dimensional one (axial symmetry and cylindrical coordinates are used). The interaction to be used is a recently developed realistic phenomenological single-particle K-matrix model 6) together with Coulomb direct and exchange potentials, spin-orbit and Brueckner rearrangement potentials. This K-matrix model has been used to study the statics of heavy-ion potentials 7) in the static Hartree-Fock approximation [ref. 7) will be referred to as CHK in this paper]. It was shown, in CHK, that the static potentials can simultaneously reproduce the asymptotic experimental potentials and the total energy versus quadrupole moment in the fused, compound nuclear region. Fig. 1 is taken from CHK and shows the static cluster potential energy for 12C + 12C versus the effective separation R defined as

R

= ?2it

'

It -- A l + A 2 ,

(1)

where A

Q°t = E <2 Z 2 - X 2 - y2>,,

(2)

i=l

is the total mass quadrupole moment of the system. We note that R is the separation between the centers whenever the fragments are spherical. The reader is refered to CI-IK for further explanations concerning this figure. Here we note that although we will present graphs for the development of R versus time, we will not be able to show the quantity which would correspond to the time dependence of Vet, due to damping of the collision. We will return to this matter later on. We show, in this paper, the T D H F predictions for the head-on scattering of two 12C ions at incident lab energy per projectile nucleon E~,b/Ap varying between 3.2 and 60 MeV/nucleon. The quantity Elab/Ao is preferred as a description of the collision because it is given by h2 Eiab/ap = ~m Ik'¢d2' (3) where k,c~ is relative momentum of a nucleon at rest in the projectile with respect to a nucleon at rest in the target, and is therefore independent of the masses of either target or projectile. Various compact graphical mefhods of displaying the otherwise unmanageable wealth of results are employed. Sect. 2 discusses the method of calculation, especially in those areas where it~differs most from the static calculation of

TDHF CALCULATION I

I

I

I

473 i

CLUSTER POTENTIAL FOR 12C~'12C

-150

I

vs R

COULOMB

-160 FROM EXR"',~ /¢-

-170

STRONG ABSORPTION

ff

-I 8 0 -190 -200 -210

0

I 2

t

I 6

4

I 8 R (fro)

I 10

r 12

14

Fig. 1. Static cluster potential energy o f 12C--~-llC versus the effective radius R, as obtained in ref. 7). The T D H F calculation is initialized with the two a=C ions at the strong absorption radius R = 8 flu, a n d with a velocity potential o f the f o r m z ( r ) = -- 1 / 2 k o V ' a - T ~ , a = 0.6 fm. The same K-matrix model is used for b o t h the static potential and the T D H F calculation.

CHK Sect. 3 presents our results in a compact graphical format whenever possible, while sect. 4 discusses both these results and the problems involved in extending the realistic K-matrix calculation to finite impact parameters in three dimensions. 2. Method of calculation

The TDHF method assumes that one can represent the time development of a heavy-ion reaction by a single Slater determinant (SSD) with occupied orbits ~k](r, t), 2 = 1 , . . . , ½Z or ½N. We have assumed that spin is saturated and that the Hamiltonian is axially symmetric. The symmetry remains axial throughout the scattering because the initial densities of the fragments are axially symmetric and the collisions are head-on. The orbits 2 are classified using the z-projection of angular momentum j~. Spin saturation implies that both orbits with +j= and -Jz are occupied. The label is the isospin label. We need to solve the TDHF single-particle (s.p.) Schr6dinger equation 1)

h'~Ua(r, t)

=

,.h

v

~

t),

(4)

where the mean field s.p. Hamiltonian h" is a function of the instantaneous densities D'(r, t) with the general form ha

h'(k,D(r)) = ~m Ik'12+v~(k'D(r))+As(k'D(r))+v~c°"~+v't'°"

(5)

474

J . A . M A R U H N A N D R. Y. C U S S O N

Here the first term is the usual free kinetic energy, v~ is the nuclear interaction potential, As is the Brueckner rearrangement potential, O~o~ is the direct plus exchange Coulomb potential and v~. ,, is the spin-orbit potential. The Hamiltonian (5) and the oscillator basis used to expand the eigenfunctions are identically the same in the present work as in CHK, including the parameters of the interaction. We will therefore not discuss the details o f h in this work. We only recall that h is chosen to saturate nuclear matter at an equilibrium fermi momentum kFc = 1.36 fm-1, with an equilibrium binding energy of 16.4 MeV/A, a nuclear compressibility of 150 MeV and a symmetry energy coefficient of 66 MeV. The Hamiltonian was studied in detail in ref. 6) and found to give good results for the binding energies, radii and single particle energies of spherical nuclei from 160 to 2°Spb, while retaining a realistic energy dependence, as determined from optical model considerations. The Hamiltonian was also shown, in CHK, to give good results for static heavy-ion potentials such as 12C--1-12C. Fig. 1 shows the static potential of C H K versus the effective separation. The potential v~(k, D(r)) depends on the density mainly linearly but in order to obtain proper saturation properties, higher order terms are also included. The density is time dependent and is obtained as D'(r, t) = Z n,~l~]( e, t)l z,

0 _-< na < 1,

(6)

where dn~. = 0, dt

E n~. = A t +A2. ~.,,

(7)

These equations represent a slight generalization of the T D H F method since a SSD total wave function implies that all the na are zeros or one. The generalization (6) and (7) is employed because it is not convenient to deal with t z c as a SSD and maintain its spherical symmetry unless one uses a pure spherical shell-model configuration of (s~) 4 (p~)8. As we will shortly see a more sophisticated initital wave function was chosen. The wave function ~kx(r, t) must in practice be represented as a finite number of degrees of freedom. This can be done by choosing ffx(r~, t), k = 1, 2 . . . . K, the values of ~ on a finite space grid as the independent quantities 1, s), or as we do here, by expanding 0k in a deformed harmonic oscillator basis ~k(r) as t) = y

(8)

k

F r o m here on the isospin label is suppressed although it is understood that separate equations for protons and neutrons are to be solved. The basis states qS,(r) have been discussed in detail in CHK. For the present application we only note that k stands for n,, n, and IJ=l, the oscillator quantum numbers in the z, and r = ~/x 2 +y2 direction, and the projection of total angular momentum along the z-axis j , = Ix +½trr We have assumed throughout that orbits are occupied in time reversed pairs with + J r Our

TDHFCALCULATION

475

basis is energy selected, with the frequency parameters ho~z = 7.77 MeV,

ho9r = 12.6 MeV.

(8)

All states up to a kinetic energy of 125 MeV were included which resulted in matrices having dimensions of 45 and 34 components for the angulaI momentum values JJzl = ½ and ~-, mainly relevant for this problem. The basis can accomodate two ~2C nuclei separated by a distance of 8 fm between centers with an accuracy of approximately 0.2 MeV on the total static energy of the system. The T D H F Schr6dinger eq. (4) can be rewritten in terms of the expansion coefficients C~(t) as

Z hkt(t)C~(t) = ih ~ C~(t), t Ot

(9)

where the hk~(t) are the (time-dependent)matrix elements . All the time dependence occurs through the time dependence of D on t. At any one time to the matrix elements hu(to) are obtained exactly as in C H K , e.g. using numerical integration techniques (given D(r, to) ). A solution of eq. (9) requires an initial set of C~(0), to be discussed shortly. Right now we outline our method of integrating eq. (9) in time. Several methods are available to step eq. (9) in time, provided the time dependence of the Hamiltonian matrix is not too violent. The most straightforward one consists in using a power series expansion technique to derive predictor-corrector equations s). A simpler but somewhat less accurate technique has been used here, in order to save computation labor. It goes as follows: let C~(to) be known and C~(t o +At) be desired. We suppose for now that the Hamiltonian matrix likt(t o +½At) can be obtained. The method consists in replacing the exact eq. (9) by the approximate equation

Z]tkt(to+½At)C~(t) = ibeX(t),

t o < t < to+At.

(i0)

1

This equation has a constant ~kt matrix and can be solved exactly once the eiyenvalues Ea of the instantaneous mean Hamiltonian/i k, are known. Let = e/,Iu,

p = 1. . . . N,

(11)

1

be the eigenvalue equation for gu, where N is the basis dimension. An equation of the type (11) is easily solved for each ILl value involved. It is easy to show that eq. (10) is then solved by Ckx(t) = E ~kt(t)C~(to), (12) l

where =

g, r~e_,,_ i~,,tt- ,o)/hrl _,,.

(13)

Ig

We note that the S-matrix is exactly unitary so that no instabilities can occur at any time. The use of an ~ evaluated half-way through the time step also contributes to

476

J.A. MARUHN AND R. Y. CUSSON

the stability, especially towards energy conservation. One could obtain h at the halfstep by first using eq. (12) to advance C~ to the half-step using h at the beginning of the step. A much faster way consists of using the continuity equation

OD(r, t) + V . J(r, t) = 0,

(14)

0t where the current J(r, t) is the usual

J(r, t) = ~ h ~a nx(~O*(r,t)V~(r, t)-(V~0*(r, t))~a(r, t)]. For a small

(15)

At the eq. (14) yields

O(r, t+At) ~ D(r, t)-AtV . J(r, t).

(16)

The eq. (16) is also used to pick the size of interval At over which one replaces the exact Schr~dinger equation (9) by the approximate one (10). Letting IIV'gl[ represent the largest value of [V. J] in the region of interest, at the beginning of a step, we pick At as At = IIADll IIV "~ll'

IIADII = 0.0014 nucleons/fm3,

(17)

which suffices to insure that the maximum change in D(r, t) during the step will be no more than 2 % of the equilibrium value of 0.07 nucleons/fm3. Having obtained D atthe half-step using eq. (16) one then obtains ~ at the half-step and diagonalizes it to obtain ~ and F~. These are used in eqs. (12) and (13) to obtain C~ at the end of the step, from which D is constructed at the end of the step. This method assures total energy conservation to about 0.01 MeV per step and an overall drift of the total energy of 1 MeV or less during the complete scattering. The total number of steps depends on the initial velocity of the collision and the collision time in such a way as to be nearly constant for most of the events studied here. About 200 steps were needed for scattering events and nearly 400 steps were taken to describe the capture event at 3.2, MeV/A. The initial wave functions ~Oa(r,0) were obtained as

0) =

(18)

where ~b°a(r) represents a solution of the static Hartree-Fock problem and x(r) is a velocity potential common to all the nucleons in the system. The static problem was solved using the code of CHK for the case of two 12C spherical ions held at a distance of 8 fm from each other. The initial orbits were mostly s~, p~ (Jz = ½) and p~ (J, = ½) orbits of very slightly deformed right and Ieft positioned 12C nuclei. Because of the residual pairing force, also included in the static calculation, about 2 % of the wave functions were contained in some excited Jz = ½ orbits. The velocity potential is

T D H F CALCULATION

477

chosen as

z ( r ) = --½kr~l x/a2 + z 2,

a = 0.6 fm,

(19)

where krcl is the relative momentum given in table 1. The velocity field 1

v(r) --

O(r)

J(r),

(20)

due to initial wave functions given by eq. (18) is

,(r)

= h_ vz(,)

hkr©I z 2m 4a-r-+

=

m

.

(21)

This shows that the value a = 0 would give the desired uniform velocity field v = __.Vo on each side of the z-axis. This was not done because it represents a singular velocity potential. It must be kept in mind that we did not obtain the initial wave functions for two shifted separate 12C ions but rather for one 24-nucleon compound system consisting of two nearly but n o t completely separated clusters. Thus the passage from ~b°(r) to ~Ja(r, 0) in eq. (18) must be made by a basis transformation of the form

c~(o) = Z sI, c'L

(20)

1

where C~° represents the expansion coefficients of ~b° and S~ is the matrix

(21)

=f

The sum in I has to be truncated to the size of the basis used and the transformation (20) will remain accurate provided the matrix Skz does not have large far-off-diagonal elements. This in turn implies a smooth x(r). The value a = 0.6 represents a suitable compromise between the need for a uniform velocity field (small a) and a smooth X (large a). The accuracy of the transformation (20) can be verified by checking the unitarity of the resulting Ck(0) coefficients. At the highest energy studied the particle number was off by less than 0.1 particles and was exact to five digits at the lowest energy. Another check of the initial transformation from the static configuration to the moving one (eq. (20)) can be obtained by computing the total collective kinetic energy /'col given by

T~l = -[d3r IJ(r)lZ . d

(22)

2nO(r)

This should equal (A 1 +A2)(h2/2m)(½k,,i) 2 for a perfectly uniform velocity field. The actual initital values are given in table 1 and indicate a discrepancy of about 3 % in the initial total collective energy for the slow collision, increasing to about 11 ~o for the fastest collision. The increase comes partly from basis truncation effects at

478

J. A. MARUHN AND R. Y. CUSSON TABLE1 The 12C-{-12Chead-on collisions at various incident energies

Elab/Ap ") k~,l

Eceo~n'(in)c) E~n.(out) c) Rteff(min) d) t*(Rmta)

(MeV) (fin) -1 3.2 6.4 12.8 19.2 25.6 32.0 51.2 64.0

(MeV)

0.393 18.6 0.556 36.6 0.786 72.6 0 . 9 6 2 107.9 1.11 143.8 1.24 180.3 1.57 262.1 1.70 341.8

(MeV)

(fm)

(fm/c)

b) b) 35.5 60.6 84.8 107.2 187.7 248.1

4.16 3.64 3.32 2.99 2.68 2.53 2.48 2.45

56 45 40 33 28 25 22 19

to(Re) ") x*

xo

DR o

(z --- O, p = 0) (N/fm3) (fm/c) (fro) (fro) 140 80 70 57 56 47 40

4.63 5.25 6.59 6.67 6.49 6.51 7.25 6.78

16.2 13.2 14.1 13.3 14.6 15.5 14.3

0.153 0.156 0.163 0.176 0.177 0.175 0.172 0.184

*) Nominal value only; see third column for actual value used. b) System has effectivelycaptured; if it decays later the energy will be just the Coulomb energy at scission. ¢) Coulomb energy of about 6.5 MeV should be added to find asymptotic energy. ~) This radius is Rteff as discussed in the text. °) Time required to return to initial R,u. =

4Q°t/2p,

higher momenta. This 11% difference between ideal and actual kinetic energy does not mean that a similarly large error in the total energy is made once the system is propagating in time. Even at the higher velocities the changes in the total energy during a step remained of order 0.01 MeV as remarked earlier. This stability comes from working in an orthonormal basis of wave functions. The initial discrepancy does however, give some idea of the influence of the basis on the development of the dynamics.

3. Results The results of the present calculation are summarized in table I, whose first entry is the nominal lab energy per projectile nucleon. As explained earlier this is used, via eq. (3), to determine the relative motion m o m e n t u m given in col. 2 of the table. The actual total collective energy, eq. (22) is given in col. 3 and is the asymptotic initial c.m. energy minus the Coulomb energy at a separation of 8 fm ( ~ 6.5 MeV). The calculation is stopped in most cases when the effective radius, given by eq. (1), returns to the initial starting v a l u e R o = 8 fro, or when a steady state seems to be established. The resulting outgoing c.m. collective energy is given in col. 4, and is plotted in fig. 2. The figure shows that once a threshold of about 37 MeV above the Coulomb barrier is passed the outgoing energy is approximately a linear function of the incoming one. The formula T(out) = 0 . 8 T ( i n ) - 2 8 MeV represents the dashed line drawn through the points. The threshold corresponds to a lab bombarding energy of ~ 81 MeV for a 12C beam on a t2C target at rest. Unfortunately, head-on collision of identical fragments cannot be detected experimentally due to the Coulomb interaction and the

TDHF CALCULATION i

i

i

[

I2Co+IzC,

i"u

i

i

i

u

i

479 i

i

[

i

i /

OUTGOING vs INCOMING

TOT. c.m. COLLECTIVE KIN. ENERGY

50C Tc'm'(out) (MeV)

T(out) =0.ST(in

)-28

(MeV)

//o t "

200

/,

100

)," /

/

/

/

/-

./

O_O/ n I

I I I 100

L I £ I I 200

I

, 1 i 300

l

Tc'm'(in){MeV) Fig, 2. Total outgoing c.m. collective kinetic energy of the system versus the incoming energy. Less than 20 MeV o f the outgoing collective energy is internal collective energy of the fragment; the rest is relative motion kinetic energy. The difference between incoming and outgoing energy is internal fluctuation (thermal) energy of the fragments. The points correspond to the events studied here and are summarized in table 1. The dashed line obeys the e q u a t i o n / ' ( o u t ) = 0.ST(in)--28 MeV.

Pauli requirement that the cross section be symmetric about 90 ° in the c.m. system. +All head-on collision effects are masked by the zero-degree Rutherford singularity. However, one could look at angles near zero or 180 ° for an inelastic group with the energy predicted here, since the angular dependence of the energy loss is a smooth function 2, s). The energy dependence of the minimum effective total radius is given in col. 5 of table 1. Upon referring to fig. 1, which shows that the ground state of 2+Mg has R ~ 3.2 fm, we see that one must go to very high bombarding energy to penetrate past the ground-state radius. This is due to the incompressible nature of the nuclear fluid and to the substantial collective energy loss during the scattering. The time t* required to reach this minimum radius is given in col. 6, in units 9) of fm/c, (1 fm/ c = 0.33 x 10 -23 sec). To obtain a better idea of the scaling of the reaction as one changes the bombarding energy, we show, in col. 8, the quantity x*

X* = __h krclt, '

(23)

m

which is a characteristic distance traveled by a hypothetical nucleon, until the turning point of the collision. Since the distance between the two centers was initially 8 fm, one would expect x* to be less than or equal to 8 fm. The observed value begins near 5 fm and slowly increases towards the 8 fm limit. This provides some indication that the finite size of the ions becomes less important as we go to higher energies. Col. 7 shows the time, to, required to return to the initial radius. Apart from the first two

480

J . A . M A R U H N A N D R. Y. CUSSON

energies, where one has effectively captured the corresponding distance Xo = (h/ m)k,,ffo, given in col 9 is nearly constant. One might interpret this as saying that the finite size effects which bring Xo down from its 16 fm nominal value, are nearly unaffected by relative motion energies up to 64 MeV/A. This may be a feature of the present K-matrix model, whose momentum dependence appears to implement this behavior. The last col of table 1 shows the value of the central density at the minimum radius. Again one sees that it increases only very slowly above the equilibrium value of 0.16 nucleons/fm3. As we will shortly see the maximum density does not actually occur at the minimum radius but somewhat earlier. More information on the collision process can be obtained from figs. 3 and 4, which show the time dependence of the radius R, the total central density Dt(z = O, p = 0), and the total collective energy ~ . l for the lowest and highest energies considered here. At the lowest energy the capture effect is seen through the oscillations in R. We note that we are unable to exclude the possibility that the system may actually come apart at some later time. In this case we would actually have a resonance instead of permanent capture. Such resonances have been observed in one dimension 4). Fig. 3 shows that the collective energy fluctuates with several periods, some as short as 40 fm[c. As we wi.ll see shortly several collective excitation modes are excited and energy is continuously being transferred from kinetic to potential energy of the various modes. The oscillations of the central density in fig. 3 again represent the net effect of many modes propagating through the compound 81{

~ZC +~C E~./Ao= 3.2 MeV/A • RADIUS;TOT. COL. EN., TOT. CENTRAL DENSI]Y vs TIME (fro/c) • R(fm)

6 5 centro

, /

;

o6

', "J

,

,:

~,- ~ 1 " ,' TJoI/4 ( MeV ) I I lO0 200 t (fm/c)

~., t 500

Fig. 3. Effective radius, central density and total c.m. collective kinetic energy versus time, for a bombarding energy o f 3.2 MoV/)I, above the Coulomb barrier. The system appears to have captured a n d is oscillating between a highly necked configuration with R m 7-7.5 fm and an ellipsoidal configuration near R = 3--4 fro. Se~ figs. 10--12 for the details o f the density and current distribution at various times. The small fluctuations in the central density are due to the non-uniform central density o f the clusters and do not appear to have any special physical meaning.

TDHF CALCULATION

48.1

i

12C*12C,EL/Ap=64 MeV/A RADIUS,TOT.COL.EN.,TOT.CENTRAL /R(fm)

._.~. DENSITYvS TIME(fro/c)/ :\, /.--. '~,,.

5j-,

t-

I

O/ 0

I

/

",..

"

/, ....... "~.~"

/

I I0

-

/

\ .

equil, central

.

"'a

.,

\--~L,~Z5(N/'r.3

I 50

20

density

| 40

t (fro/c) Fig; 4. EffectiVe radius, central density and total c.m. collective kinetic energy versus time for a bombarding energy o f 64 MeV/A, the highost yaluestudied here, The density and current distribu-

tions for this case are shown in figs. 5 and 6. system; By Contrast,: the fig. 4 is much smoother. Of course, the time scale is now so short that most collective oscillations of the system get established only after the collision is over! One notes, in fig. 4, the location of the density peak beforethe occurrence of the minimum radius. This is the behavior one would expect from a shock phenomenon but as we will see. no other indications of shock behavior were obtained. In addition, the events at 6,4 and 51.2 MeV/A have been described elsewhere i o). They and all the other events entered in table 1 are similar to the first and last one described above. At this point it is clear that an exhaustive quantitative account of the complex behavior of this colliding 24-body system would be extremely time-consuming and might still not render the flavour of the overall macroscopic behavior o f the system. F o r this reason we now present some computer generated perspective drawings of the density in the z, p plane as a function of time D(z, p, t), Figs. 5 and 6 show the density for the 64 MeV/A event, at twelve different times from beginning to end. The first frame of fig. 5 shows the initial density along the vertical axis, the z-axis is along the horizontal axis and the p-axis is given in perspective at 45 °. As can be seen from this frame and fig. 4, the density overlap is initially about 10 % of the equilibrium value. The nuclei were started as close as possible to each other in order to minimize the basis size, and therefore the computing time. The frame at t --- 11 fm/c corresponds to the highest observed value of the density at 0.23 N/fro a, to be compared with the initial value of 0.16 N/fm a. The frame at t = 18 fm/c is t a k e n near the minimum effective radius. It is seen that the instantaneous appearance is very much that of a slightly deformed A = 24 system. The collective energy is almost zero (see fig. 4). The frame at t = 25 fm/c illustrates the density at the second maximum of fig. 4. The sequence of frames in fig. 6 show how the system proceeds to come apart. One

482

J . A . M A R U H N A N D R. Y. C U S S O N

notices how much more ragged the density contours are and how little correlation there is between the small secondary oscillations at t = 40 and 48 fm/c. This time is close to the nuclear relaxation time below which little can happen. The secondary oscillations in the t = 48 frame have a characteristic wavelength of about 2.5 fin. The calculation was not pursued for later times because the clusters had separated to such an extent by then that they were beginning to touch the edge of the basis. We presume that at later times, the central peak remaining in the last frame of fig. 6, would die away. A comparison of the t = 0 and the t -- 48 frame shows that the outgoing clusters are much larger and have lower density than the initial ones. It has been pointed out 1o) that this is the behavior expected from a system which has converted a substantial fraction of its kinetic energy into fluctuation (thermal) energy. We note that the initial collective energy of relative motion has not just been lost to internal collective energy of the higher multipoles. This is clear since Tt~o=includes collective energy from all sources, relative motion as well as internal collective motion. In that context we note that the outgoing collective energy of col. 4, table l, is not just relative motion collective energy as in the beginning. Up to 15-20 MeV of it is actually internal collective motion. Fig. 7 shows similar density drawings for the collision at 6A MeV/A. This collision is indicated as a capture in table 1. In fact, fig. 7 shows that we rather have an instance

t2C +t2C,EtlAp= 64 MeV/A D (z,p) ot vorious times

0,15

D

t=4fm/c

t=O

ODE

'

"9

~k

0

~

Z(fm)

*¢t = 18 fm/c

Fig. 5. Perspective three-dimensional plots o f the total density D(z, p) versus t h e coordinates z a n d p, at E/A = 64 MeV/A. T h e f o u r t h frame, at t = 11 fm/c s h o w s t h e highest density observed in this work, D ~ 0.23 nucleons/fro 3. T h e sketches o f figs, 5 a n d 6 are best interpreted in c o n j u n c t i o n with fig. 3.

T D H F CALCULATION

483

12G+lZC, E,JAp=64 MeV/A D (z,p) of various limes

0.15

D

I = 32 f m / c

(Nucl,~ 0.10 \ fm3/ 0.05

.

-9

I=36fm/c

0 Z~m)

"l

-,9

! =48 fm/c

t=4Ofm/c

! Fig. 6. Same as fig. 5 but at later time. One observes in tha last frame the fragments emerge in a highly

excited internal state. The region displayed on the graph is smaller than t h e o n e that can be represented accurately by the basis. Thus the cutoff in the last frames is for graphical convenienceonly.

of a so-called "fully relaxed" collision 1). The early part of the collision is quite similar to the previous one and is not shown. The first frame at t = 44 fmfc corresponds to the minimum radius configuration and has essentially zero collective energy ( ~ 3 MeV). The frame at t = 62 has a maximum in the collective energy lo) which then proceeds to decrease to less than 1 MeV at t = 133 fm/c. F r o m then on the system drifts apart very slowly and begins to pick up some Coulomb energy. The last two frames show this behavior if we note the quite large time interval between them. We can obtain further information on the events by simultaneously displaying the density and the current distributions. Figs. 9 and 10 show equidensity contour lines and current arrows at various times for the 25.6 MeV/A collision. The arrow marked "calibration" shows essentially how big a current corresponding to the largest current in the initial frame would be in the subsequent frames of that same energy. Thus when the calibration arrow is large we know that the distribution being displayed is actually very slow moving. This renormalization was done in order to display the intricacies of the currents near the turning points where the velocities change direction. The fifth frame of fig. 8 shows how the system comes to rest by ejecting mass in the radial direction under the pressure of the moving matter behind it. Frame no. 6 shows that the outer edges come to rest and reverse themselves before the central portion. This frame corresponds to nearly minimum radius (see table 1). The first

484

J . A . M A R U H N A N D R. Y. C U S S O N

O.15it =44fm/c ,NuDcI,O.lo'

t2 C +l 2C, E~,/Ap=6.4MeV/A D (z,p) at various times

"

: 81 fm/c

~ : ~

~0

1=297fm/c

t =209 fm/c

-t=133fm/c

Fig. 7. Similar to fig. 5, but past the point o f minimum radius (see table 1) and for an energy o f 6.4 MeV/A. This event is a fully ielaxed Collision in which th e system is'essentially at rest at t ----133 fm/c and slowly drifts apart while acquiring the Coulomb energy (t = 209 and 297 fro/c).

oL

CAU BRATION 12C * 12C, E.e/Ap =25.6 MeV/~

4

J(z,p)& D(z,p) ot various times

CURRENT

p(fm) 2

. :'- i . :

C

.

.

.

.

. : .....

.

.

]"

t = 17 f m / c

.'- . . . .

".

.

..:.

.." ....

:

Z(fm)

t=lfm/c DENSITY

41

.

.:-'_

......

t=7fm/c

'

t=llfm/c

i

t = 22 fm/c

-"

t = 29 fm/c

Fig. 8. Density and current distributions versus Z and # at a bombarding energy o f 25.6 MeV/A. The calibration arrow in each frame is essentially the size o f arrow which would represent a current equal to the largest current plotted in the first frame o f this energy. The point o f minimum radius occurs near t ----27 fm/c (see table 1).

T D H F CALCULATION

485

--CALIBRATION 12C, +12C F,t/AP= 3.2 MeV/A ~t

CURRENT

J'(z.p)BD(z,p)at various times

P(fm) 5

t =8fm/c DENSITY

$]

Z (fm)

1=22fm/c

t=33fm/c

,,, M

:;~

t = 50fm/c

.-, ,-. - ~ . ~ :

t = 56fm/c

= 82 f m / c

Fig. II. Same as fig. 10, but at later time. This events appears to be capture event where the system oscillates between the stretched configuration at t -----123 fm/c (second frame) and the ellipsoidal configuration at t ---- 56 frrdc (fifth frame, fig. 10).

~CALIBRAT ON I I ~ " ' ! 2C + 2C,Ez/Ap=25.6MeV/A

6 4 2

'-~-4-:~'6 6

2 t = ,,'32fm/c DENSITY

4

2F

0

l

!J(z,p)&D(z,p)at

CURRENT ::: :::

:_:.: 4 ~ Z(fm)

various times !" -~

t= 3 6 f m / c

t =40fm/c

--

I

Fig. 12. Same as fig. 11, except at still later times. From t = 218 fm/c on, an asymmetry between the two fragments can be seen. The asymmetry first develops because o f minor numerical errors but its growth may be physical and indicate a tendency to break up in non-symmetrical fragments, which would contribute to the mass asymmetry o f the reaction.

486

J . A . M A R U H N AND R. Y. CUSSON •

~CALIBRATION '

6

12C -~t2C, E l / A p = 3 . 2 MeV/A

4

CURRENT ,-:.

::

P(fm) 6

-6 -~

-2

0

2

t=ll2fm/c DENSITY

4

4

J(z,p)&D(z,p)at

~rious times

i'-i

i

.

.

.

.

.

.

.

.

.

.

.

.

.

.

6 Z(fm) •

l= 1 4 7 f m / c

t =123fm/c

t =168fm/c

t =130fm/c

t =183fm/c

Fig. 9. Same as fig. 8, except at later times. The frame t = 58 fm/c has the same radius as the initial configuration. The last two frames show that the outgoing fragments arc bigger and have lower density than the incoming z=C ions.

CALlBRATI0 N ~2C + ~2C, E~/Ap = 3.2 MeV/A J(z,p)&D(z,p)at various times

CURRENT P(fm) -6-4-2

0

2

t=199fm/c DENSITY

t =2 5 5 f m / c

4 Z (fm)

t :218 fm/c

t = 2 6 4 fm/c

Fig. I0. Similar to fig. 8 but at the lowest bombarding energy studied here, Refer to fig. 3 to find the radius at the times shown.

t =252 fm/c

t =2 8 0 f m / c

El.b/Ap =

3.2 MeV/A.

TDHF CALCULATION

487

three frames of fig. 9 show how the radial motion is gradually converted to motion along the x-axis. The last three frames of fig. 9 illustrate the effects of the residual oscillations on the current distribution. The fifth frame at t = 58 fm/c is taken at the instant where the effective radius is back to the initial value. A comparison of this frame with the first frame of fig. 8 reveals the extent to which the clusters are now enlarged, presumably due to their heating up. One also notes in the last frame of fig. 9, the onset of octupole vibrations of the fragments 8). The collision at 25.6 MeV is, everything considered, relatively smooth and uneventful. Our last illustration of the detailing power of the T D H F method will be given by considering the current-density display for the capture event at 3.2 MeV/A. The fourth and fifth frames of fig. 10 show the approach to the minimum radius to be similar to the higher energy 25.6 MeV/A case. Complete reversal has set in the last frame at t = 82 fm[c and from fig. 3 we see that this is a time of maximum collective energy. The first three frames of fig. 11 then show the approach to the second stationary point at t = 130 fm/c. Again we note that the outer edges turn over first and generate a compression region inside the two clusters. The approaches to the second radial stationary point at t = 215 fm/c is shown in the subsequent frames, continued in fig. 12. Starting at t = 218 fm/c in fig. 12, we see that an oetupole vibration of the compound system is developing. The initial feeding of this mode is due to computer roundoff but its growth may indicate some slight instability towards mass asymmetry in the decay. The last frame at t -- 280 fm/c may again be compared with the first frame of fig. 12. The density at the center of the clusters is somewhat lower than the initial one and the neck region is still quite large. But the size of the cluster is not enlarged as in the last frame of fig. 9 (E/A = 25.6 MeV/A). There has been considerably less dissipation and the shape at t = 280 fm[c is similar to the one at t = 123 fm/c (frame 2, fig. 1I), the time of the second stationary point. A recurrent shape oscillation between these two configuration (t = 123 and 218) appears to be taking place, and is reminiscent of liquid drop motion. The secondary oscillations seen in the central density (fig. 3) appear to be related to the dip in the central density of the 12C clusters and do not show up here as interesting distinct modes. 4. Discussion

We have seen that a realistic K-matrix model can indeed be used to discuss heavyion collisions in two dimensions. Unfortunately, the complexity of the interaction has so far precluded its use in a three-dimensional calculation. It would be interesting to develop valid two-dimensional approximations, by going to a rotating frame 12) and devising a suitable self-consistent moment of inertia, to treat the low energy collisions. The K-matrix model used here would then be applicable and considerable computing time could be saved by using the recently developed s) fast Fourier transform scheme instead of the present harmonic oscillator basis. Clearly the present calculation is more of an exploratory nature, due mainly to the restriction to head-on collisions.

488

J . A . MARUHN AND R. Y. CUSSON

The neutrons and protons have been treated separately in this work. For such a light system the Coulomb force does not play a dominantrole and we have observed that the proton quantities follow the neutron ones to within a few percent, due tothe large attractive n-p force. The ability to treat these two Separately will, of course, become a mandatory feature in the study of heavier reactions. An interesting conclusion which can be drawn from the present results concerns the time scaling of the reaction. Over a fairly broad range of incoming energies the sequence of density distributions can be obtained to first approximation by scaling the time from some standard reaction say at 25.6 MeV/A. Thus if there are any shock waves to be observed, it must be at either higher energy or higher mass than considered here. We have not discussed here the matter of computing a potential energy versus time. Such a potential cannot simply be obtained by subtracting the collective energy from the total (conserved) energy due to the single particle dissipation effect observed here. An estimate of the fluctuation energy is needed along the path. Wong 13) has suggested a method based orr a sequence of corrections to the Thomas-Fermi energy. Preliminary results with this method are encouraging and may lead to useful potential energy estimates l*). Given such a potential one could then draw potential curves similar to those of fig. 1, but for different bombarding energies, and thereby obtain a dynamical prediction of the energy dependence of the optical potential. Such a project is all the more interesting that it does not require One to consider non head-on collisions, provided the angular-momentum dependence of the optical potential can be pinned down by some other method. We acknowledge stimulating discussions on the TDHF method with Drs. S. Koonin and J. W. Negele. References 1) L. G. Moretto and J. S. Svemek, Symp. on macroscopic features of heavy,ion collisions, Argonne, Illinois (April 1976) LBL report 5006 (1976) 2) J. W. Wilczynski, Phys. Lett. 4713 (1973) 484 3) P. A. M. Dirac, Prec. Cambridge Phil. Soc. 26 (1930) 376 4) P. Bonche, S. Koonin and J. W. Negele, Phys. Roy. C13 (1976) 1226 5) R. Y. Cusson, R. K. Smith and J. Marulm, Phys. Rev. Lett. 36 (1976) 1166 6) R, Y. Cusson, H. P. Trivedi, H. W. Meldner and M. S. Weiss, Phys. Rev., submitted 7) R. Y. Cusson, R. Hilko and D. Kolb, Nucl. Phys. A270 (1976) 437 8) S. Koonin, Phys. Lett. 61B (1976) 227 9) C. Y. Wong, unpublished report (ORNL) 10) R. Y. Cusson and J. Maruhn, Phys. Lett., in press 11) C. Y. Wong, Saclay Lectures (1975) (unpublished) 12) S. E~ Koonin, K. T. R. Davies, H. Feldmeier, S. J. Krieger, V. Marulm-Rezwani and J. W. Negele, to be published 13) C. Y. Wong, to be published 14) R. Y. Cusson, J. Maruhn and R. K. Smith, to be published