Relativistic Coulomb excitation within Time Dependent Superfluid Local Density Approximation
Abstract
Within the framework of the unrestricted timedependent density functional theory, we present for the first time an analysis of the relativistic Coulomb excitation of the heavy deformed open shell nucleus U. The approach is based on Superfluid Local Density Approximation (SLDA) formulated on a spatial lattice that can take into account coupling to the continuum, enabling selfconsistent studies of superfluid dynamics of any nuclear shape. We have computed the energy deposited in the target nucleus as a function of the impact parameter, finding it to be significantly larger than the estimate using the GoldhaberTeller model. The isovector giant dipole resonance, the dipole pygmy resonance and giant quadrupole modes were excited during the process. The one body dissipation of collective dipole modes is shown to lead a damping width MeV and the number of preequilibrium neutrons emitted has been quantified.
pacs:
25.70.De, 21.60.Jz, 25.70.z, 24.30.CzCoulomb excitation represents an ideal method to probe the properties of large amplitude nuclear motion, because the excitation process is not obscured by uncertainties related to nuclear forces. The excitation probabilities are governed by the strength of the Coulomb field only and they can be fully expressed in terms of the electromagnetic multipole matrix elements Winther and Alder (1979); Baur and Bertulani (1986); Emling (1994); *abe; Bertulani and Ponomarev (1999); Lanza et al. (1999); Hussein et al. (2004). From the theoretical point of view, Coulomb excitation can be treated as a textbook example of a nuclear system being subjected to an external, timedependent perturbation. However, in order to be able to probe nuclear collective modes involving multiphonon states for example Boretzky et al. (2003); Ilievski et al. (2004), a large amount of energy has to be transferred to the nuclear system. Thus the interaction time should be relatively short and the velocity of the projectile has to be sufficiently large for an efficient excitation of nuclear modes of frequency , the collision time has to fulfill the condition that . Here is the impact parameter, is the projectile velocity, and is the Lorentz factor. One of the best known examples of collective nuclear motion is the isovector giant dipole resonance (IVGDR). A reasonably good estimate of the IVGDR vibrational frequency is for spherical nuclei. It implies that the excitation of a heavy nucleus to such energies requires a relativistic projectile.
We report on a new and powerful method to study relativistic Coulomb excitation and nuclear large amplitude collective motion in the framework of Time Dependent Superfluid Local Density Approximation (TDSLDA). This is a fully microscopic approach to the problem based on an extension of the Density Functional Theory (DFT) to superfluid nuclei and timedependent external probes, where all the nuclear degrees of freedom are taken into account on the same footing, without any restrictions and where all symmetries (translation, rotation, parity, local Galilean covariance, local gauge symmetry, isospin symmetry, minimal gauge coupling to electromagnetic (EM) fields) are correctly implemented Bulgac (2013); Stetcu et al. (2011). The interaction between the impinging U projectile and the U target is very strong (, where is the fine structure constant), which thus require a nonperturbative treatment, and the excitation process is highly nonadiabatic. We assume a completely classical projectile straightline motion since its de Broglie wavelength is of the order of fm for . In evaluating the EMfield created by the uranium projectile with a constant velocity along the zaxis, we neglect its deformation. The projectile produces an EMfield described by scalar and vector LienardWiechert potentials. These fields couple to a deformed U target nucleus residing on a spatial lattice, see Ref. sup . The interaction leads to a CM motion of the target as well as to its internal excitation and full 3D dynamical deformation of the target. In order to follow the internal motion for a long enough trajectory that allows the extraction of useful information, we perform a transformation to a system in which the lattice moves with the CM. The required transformation for each single particle wave function reads , with describing the CM motion and the momentum operator. The equation of motion (simplified form here) for becomes
(1) 
where is the CM velocity and the total current density.
The target nucleus is described within the SLDA and its time evolution is governed by the TD meanfieldlike equations (spin degrees of freedom are not shown):
(3)  
The singleparticle Hamiltonian and the pairing field are obtained selfconsistently from an energy functional that is in general a function of various normal, anomalous, and current densities. The external electromagnetic (EM) field has the minimal gauge coupling (and similarly for the timecomponent) in all terms with currents, as well as in the definition of the momentum operator in Eq. (1), details in sup . In the current calculation, the Skyrme SLy4 energy functional Chabanat et al. (1998a); *Sly4E was adopted, with nuclear pairing as introduced in Ref. Yu and Bulgac (2003), which provides a very decent description of the IVGDR in U Stetcu et al. (2011). The coupling between the spin and the magnetic field was neglected. The Coulomb selfinteraction between protons of the target nucleus is taken into account using the modification of the method described in Ref. Castro et al. (2003), so as not to include contributions from images in neighboring cells. For the description of the numerical methods see Refs. Bulgac and Roche (2008); Bulgac et al. (2011) and many other technical details can be found in sup .
The DFT approach to quantum dynamics has some peculiar characteristics. Unlike a regular quantum mechanics (QM) treatment one does not have access to wave functions, but instead to various onebody densities and currents. Within a DFT approach quantities trivial to evaluate in QM become basically impossible to calculate. For example, by solving the Schrödinger equation one can evaluate at any time the probability that a system remained in its initial state from , where is the solution of the Schrödinger equation (or some of its approximations). Within DFT one has access to the onebody (spin)density and onebody (spin)current with no means to compute the probability . One can calculate for example a quantity such as , but there is no obvious way to relate it to the probability . One might try to define instead through the overlap of the initial and current “Slater determinants” constructed through the fictitious singleparticle wave functions entering the DFT formalism, which is a rather arbitrary postulate. One can find quite often in literature various formulas used within DFT treatment of nuclei, which are simply “copied” from various mean field approaches, without any solid justification provided. Restrictions inherent to a DFT approach, prevent us from being able to calculate various quantities, which within a QM approach are easy to evaluate. Even though we evaluate accurate densities and currents well beyond the linear regime, within a DFT approach we cannot separate for example the emission of one and two photons from an excited nucleus, which however could be estimated relative easily within a perturbative linear response approach such a (Q)RPA. On the other hand a DFT approach has unquestionable advantages, allowing us to go far into the nonlinear regime and describe large amplitude collective motion.
The incoming projectile excites various modes in the target nucleus and the axial symmetry of the initial ground state is lost. Because U is highly deformed the energy of the first is 45 keV, corresponding to a very long rotational period, and thus during simulation time considered here ( sec.) it can be considered fixed. The identification of these modes requires certain care, since during the collision the system is beyond the linear regime and the analysis using the response function is not applicable in general. However, the information about the excited nuclear modes is carried in the subsequent EM radiation leading to nuclear deexcitation. Deexcitation to the ground state via photon emission requires times of about sec., which is four orders of magnitude longer than in the current calculations. However, it is possible to compute the spectrum of the preequilibrium neutrons and gamma radiation, which allows the identification of the excited nuclear modes. We can accurately treat the system as a classical source of electromagnetic radiation and the time dependence of the proton current governs the rate of emission, see Refs. sup ; Baran et al. (1996); Oberacker et al. (2012):
(4) 
with . Here , is the spherical Bessel function of order , and is the proton current. The emission rate is plotted in Fig. 1. The magnitude of this quantity indicates that the total amount of radiated energy during the evolution time (about fm/c) is rather small compared to the total absorbed energy and does not exceed MeV, which is about of the deposited energy reported in Table 1 below. This implies that the effect of damping of nuclear motion due to the emitted radiation can be neglected for such short time intervals. Consequently, the decreasing intensity of the radiation, see Fig. 1, is merely related to the rearrangements of the intrinsic structure of the excited nucleus caused by damping of collective modes due to the onebody dissipation mechanism. It has to be emphasized that within the framework of the presented approach one is able to extract only a small fraction of the spreading width , which is due to the onebody dissipation mechanism. The twobody effects require e.g. stochastic extension of TDSLDA which would allow for a dynamic hopping between various meanfields, and thus could account for collisional damping as well.
TDSLDA provides the EM power spectrum sup ; Baran et al. (1996); Oberacker et al. (2012), arising from the multipole expansion in Eq. (4). This quantity is different from what one would compute within a linear response approach or first order perturbation theory, see, e.g., Refs. Winther and Alder (1979); Baur and Bertulani (1986); Emling (1994); *abe; Bertulani and Ponomarev (1999); Lanza et al. (1999); Hussein et al. (2004), which provides the excitation probability only , where is the transition density and  the external field. is proportional to the excitation probability, here in the nonlinear regime, and the subsequent photon emission probability as well. A typical example of the emitted EM radiation for a given impact parameter is shown in Fig. 2, here due to the internal excitation of the system alone. The EM radiation due the CM motion has been calculated separately (see Table 1 and sup ).
In Fig. 2(a) the emitted radiation shows a well defined maximum at energy MeV which corresponds to the excitation of IVGDR. We have applied a smoothing of the original calculations using the halfwidth of MeV. Therefore, the original separate peaks split due to the deformation of U merge into a single wider peak. However, at larger frequencies another local maximum exists which we associate with the isovector giant quadrupole resonance (IVGQR). In order to rule out other possibilities we have repeated the calculation by retaining only the dipole component of the electromagnetic field produced by the projectile sup . The results are shown in Fig. 2(b). In this case, the highenergy structure above 20 MeV disappears, evidence that the high energy peak is related to the IVGQR. Noticeable contribution to the total radiation is coming from the quadrupole component of radiated field.
At low energies a change of slope occurs at about MeV, present at the same energy for all impact parameters and orientations, see Ref. sup . It indicates a considerable amount of strength at low energies, giving rise to an additional contribution to the EM radiation. We attribute this additional structure to the excitation of the pygmy dipole resonance (PDR). The inset of the figure 2 shows the spectrum of emitted radiation due to this mode. The contribution to the total radiated energy coming from the PDR is rather small and reads: , , and keV for impact parameters , , and fm, respectively. It corresponds to about %, %, %, % of the emitted radiation (due to internal motion) respectively. The relatively small amount of E1 strength obtained in our calculations, in the region where the PDR is expected, agrees with recent measurements Hammond et al. (2012).
The comparison between the average energy transferred to the internal motion of the target nucleus for three values of the impact parameter obtained within TDSLDA and within a simplified GoldhaberTeller model Goldhaber and Teller (1948) presented in Table 1 shows that significantly more energy is deposited by the projectile within the TDSLDA. The GoldhaberTeller model is equivalent to a linear response result, assuming that all isovector transition strength is concentrated in two sharp lines, corresponding to an axially deformed target. An exact QRPA approach would therefore severely underestimate the amount of internal energy deposited, one reason being the nonlinearity of the response, naturally incorporated in TDSLDA. Another reason is the fact that the present microscopic framework describing the target allows for many degrees of freedom to be excited, apart from pure dipole oscillations. At the same time, the CM target energy alone is approximately the same as obtained in a simplified point particles Coulomb recoil model of both the target and projectile.
12.2  39.29  0.668  0.911  0.960  17.58  24.68 

14.6  19.2  0.608  0.567  0.963  10.32  14.51 
16.2  12.87  0.547  0.411  0.963  7.41  10.43 
20.2  5.41  0.444  0.199  0.961  3.43  4.84 
12.2  25.11  0.588  0.5  0.941  12.94  18.17 
14.6  13.16  0.498  0.306  0.942  7.22  10.16 
16.2  8.97  0.470  0.217  0.939  5.02  7.07 
20.2  3.8  0.367  0.106  0.934  2.16  3.05 
12.2  24.21  0.591  0.407  0.930  12.36  17.33 
14.6  12.58  0.513  0.245  0.929  6.65  9.34 
16.2  8.5  0.464  0.175  0.926  4.49  6.32 
20.2  3.5  0.353  0.085  0.919  1.78  2.51 
The average energy radiated due to the internal excitation represents only part of the total radiated energy. (One should remember that a straightforward DFT approach provides no measure for the variance.) Also, because of the spreading of the strength due to onebody dissipation only a fraction of the energy (where is the EMwidth alone and the total width of the IVGDR) is emitted as a pulse, as shown in Fig. 1. A subsequent pulse of reduced amplitude is to be expected after a delay fm/c. Since our simulations times are much shorter we are not able to see emission of the second photon, as reported in experiment Boretzky et al. (2003); Ilievski et al. (2004), where two photons were measured in coincidence. In our calculations we have followed the nuclear evolution during approximately fm/c after collision. The other component of the EM radiation arises from the CM acceleration as a result of collision (Bremsstrahlung), during the relatively short time interval . The radiation emitted from the internal motion has a much longer time scale.
We can estimate the cross section for the emission of radiation by assuming that the asymptotic transition probability for a given impact parameter is given by . Here , and are the total energies transferred to the target nucleus during the collision at the impact parameter and for the three independent orientations. Our simulation yields mb. A detailed comparison of intensities of radiation for various impact parameters and orientations is shown in Table 1. It is evident that although the intensity of radiation decreases with increasing impact parameter, the ratio between the intensities due to the internal modes with that of the CM motion remains fairly constant and depends slightly on orientation For the orientation perpendicular to the beam and parallel to the reaction plane the target nucleus is the most efficiently excited which results in a larger ratio.
The average energy deposited in the target nucleus is of the order of
the neutron separation energy. In Fig. 3
we plot the total number of neutrons inside a (smoothed) sphere of
radius fm which is slightly larger than the nuclear diameter (see
Ref. sup for details). For all these impact parameters neutrons
can leak from the excited system. Since more energy is deposited in
the nucleus with perpendicular orientation with respect to the beam,
the rate of emitted neutrons is larger in that case. However, the simulation box is large enough (about 40 times bigger than the nucleus) so that during the evolution the calculations are not affected by the emitted neutrons.
We thank G.F. Bertsch for a number of discussions and reading the manuscript. We acknowledge support under U.S. DOE Grants DEFC0207ER41457, DE FG0208ER41533, NSF grant PHY1415656, Polish National Science Centre (NCN) Grant, decision no. DEC2013/08/A/ST3/00708, and in part by the ERANETNuPNET grant SARFEN of the Polish National Centre for Research and Development (NCBiR). Part of this work was performed under the auspices of the National Nuclear Security Administration of the US Department of Energy at Los Alamos National Laboratory under contract No. DEAC5206NA25396. This research used resources of the National Energy Research Scientific Computing Center, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DEAC0205CH11231 and of the Oak Ridge Leadership Computing Facility at the Oak Ridge National Laboratory, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DEAC0500OR22725. Some of the calculations reported here have been performed at the University of Washington Hyak cluster funded by the NSF MRI Grant No. PHY0922770. KJR was supported by the DOE Office of Science, Advanced Scientific Computing Research, under award number 58202 ”Software Effectiveness Metrics” (Lucille T. Nowell).
References
 Winther and Alder (1979) A. Winther and K. Alder, Nucl. Phys. A 319, 518 (1979).
 Baur and Bertulani (1986) G. Baur and C. Bertulani, Phys. Lett. B 174, 23 (1986).
 Emling (1994) H. Emling, Prog. Part. Nucl. Phys. 33, 729 (1994).
 Aumann et al. (1998) T. Aumann, P. F. Bortignon, and H. Emling, Annu. Rev. Nucl. and Part. Sci. 48, 351 (1998).
 Bertulani and Ponomarev (1999) C. Bertulani and V. Ponomarev, Phys. Rep. 321, 139 (1999).
 Lanza et al. (1999) E. Lanza, M. Andrés, F. Catara, P. Chomaz, and C. Volpe, Nucl. Phys. A 654, 792c (1999).
 Hussein et al. (2004) M. Hussein, B. Carlson, and L. Canto, Nucl. Phys. A731, 163 (2004).
 Boretzky et al. (2003) K. Boretzky, A. Grünschloß, S. Ilievski, P. Adrich, T. Aumann, C. A. Bertulani, J. Cub, W. Dostal, B. Eberlein, T. W. Elze, H. Emling, M. Fallot, J. Holeczek, R. Holzmann, C. Kozhuharov, J. V. Kratz, R. Kulessa, Y. Leifels, A. Leistenschneider, E. Lubkiewicz, S. Mordechai, T. Ohtsuki, P. Reiter, H. Simon, K. Stelzer, J. Stroth, K. Sümmerer, A. Surowiec, E. Wajda, and W. Walus ((LAND Collaboration)), Phys. Rev. C 68, 024317 (2003).
 Ilievski et al. (2004) S. Ilievski, T. Aumann, K. Boretzky, T. W. Elze, H. Emling, A. Grünschloß, J. Holeczek, R. Holzmann, C. Kozhuharov, J. V. Kratz, R. Kulessa, A. Leistenschneider, E. Lubkiewicz, T. Ohtsuki, P. Reiter, H. Simon, K. Stelzer, J. Stroth, K. Sümmerer, E. Wajda, and W. Waluś (LAND Collaboration), Phys. Rev. Lett. 92, 112502 (2004).
 Bulgac (2013) A. Bulgac, Annu. Rev. Nucl. and Part. Sci. 63, 97 (2013).
 Stetcu et al. (2011) I. Stetcu, A. Bulgac, P. Magierski, and K. J. Roche, Phys. Rev. C 84, 051309(R) (2011).
 (12) “Supplementary online material,” .
 Chabanat et al. (1998a) E. Chabanat, P. Bonche, P. Haensel, J. Meyer, and R. Schaeffer, Nucl. Phys. A635, 231 (1998a).
 Chabanat et al. (1998b) E. Chabanat, P. Bonche, P. Haensel, J. Meyer, and R. Schaeffer, Nucl. Phys. A643, 441 (1998b).
 Yu and Bulgac (2003) Y. Yu and A. Bulgac, Phys. Rev. Lett. 90, 222501 (2003).
 Castro et al. (2003) A. Castro, A. Rubio, and M. J. Stott, Canadian Journal of Physics 81, 1151 (2003), arXiv:physics/0012024 .
 Bulgac and Roche (2008) A. Bulgac and K. J. Roche, J. Phys.: Conf. Series 125, 012064 (2008).
 Bulgac et al. (2011) A. Bulgac, Y.L. Luo, P. Magierski, K. J. Roche, and Y. Yu, Science 332, 1288 (2011), arXiv:1011.5999 .
 Baran et al. (1996) V. Baran, M. Colonna, M. D. Toro, A. Guarnera, and A. Smerzi, Nuclear Physics A 600, 111 (1996).
 Oberacker et al. (2012) V. E. Oberacker, A. S. Umar, J. A. Maruhn, and P.G. Reinhard, Phys. Rev. C 85, 034609 (2012).
 Hammond et al. (2012) S. L. Hammond, A. S. Adekola, C. T. Angell, H. J. Karwowski, E. Kwan, G. Rusev, A. P. Tonchev, W. Tornow, C. R. Howell, and J. H. Kelley, Phys. Rev. C 85, 044302 (2012).
 Goldhaber and Teller (1948) M. Goldhaber and E. Teller, Phys. Rev. 74, 1046 (1948).
Supplemental material to:
Relativistic Coulomb excitation within Time Dependent Superfluid Local Density Approximation
I. Stetcu, C. Bertulani, A. Bulgac, P. Magierski and K.J. Roche
August 17, 2019
Density functional in TDSLDA and coupling to electromagnetic field
Here we present various definitions and conventions which we have used in the manuscript. The density functional is constructed from the following local densities and currents:

density:

spin density:

current:

spin current (2nd rank tensor):

kinetic energy density:

spin kinetic energy density:

anomalous density:
where
(5) 
The coupling of the nuclear system to the electromagnetic field:
(6)  
(7) 
requires the following transformation of proton densities and currents:

density:

spin density:

current:

spin current (2nd rank tensor):

spin current (vector):

kinetic energy density:

spin kinetic energy density:
As a result the proton single particle hamiltonian has the form:
(9) 
and
(10)  
(11)  
(12)  
(13) 
Numerical Implementation
We build a spatial threedimensional Cartesian grid in coordinate space with periodic boundary conditions, and derivatives evaluated in momentum (Fouriertransformed) space. This method represents a flexible tool to describe large amplitude nuclear motion as it contains the coupling to the continuum and between singleparticle and collective degrees of freedom. For the present problem, we have considered a box size of with the lattice constant 1 fm. The time step has been set to with a total time interval of about . The projectile is initially placed at such a distance from the target nucleus that the collision occurs after . Even though inially the projectle is far enough from the target and hence the EM fields are weak, spurious excitations produced by a sudden switch of the EM interaction at are possible. They were avoided by multiplying the the EM potentials in Eq. (LABEL:LW) by the smoothing function , where fm, fm. This ensures that the EM field varies smoothly within the distance , but stay approximately equal to its physical value within the distance .
Coulomb potential on the lattice
Here we describe the method used to describe the Coulomb selfinteraction of the target nucleus.
Consider the charge distribution :
(14)  
(15) 
Note that above we have defined as (note in the formula). After the Fourier transform one gets:
(16) 
The above prescription generates however the spurious interaction between neighbouring cells.
Therefore we define the modified potential ( denote number of equidistant lattice points in each direction, , , is lattice constant):
(17) 
Clearly the Fourier transform is:
(18) 
and moreover
(19) 
where in the last term is the Fourier transformed density on the lattice . In practice it means that one has to perform forward and backward Fourier transforms on the lattice three times bigger in each direction.
This may however be avoided if one realizes that the Fourier transform of the density in the larger lattice can be expressed through the Fourier transforms in the smaller lattices:
(20) 
and we need to perform 27 FFTs on the smaller lattice for of the following quantities:
Subsequently we obtain the potential through the relation:
(21)  
Dipole component of the electromagnetic field produced by the projectile
Coordinates (see fig. 7):
(22) 
Interaction potential:
(23) 
where
(24) 
and . is the charge density of the nucleus at location and are proton wavefunctions . The vector potential is given by
(25) 
In order to extract the dipole component we used the interaction Hamiltonian:
(26) 
where one subtracts the second term which is responsible for the c.m. scattering (i.e. monopole field).
Consequently the dipole term reads:
(27) 
where r is the coordinate of one of the protons in the target. A sum over has to be performed, i.e.
(28) 
where is the part of the interaction which does not involve the intrinsic structure of the nucleus:
(29) 
Electromagnetic radiation from a nucleus described within TDSLDA
Let us consider the proton density and current (we use Gauss units):
(30)  
(31) 
where
(32)  
(33) 
Maxwell equations:
(34)  
(35)  
(36)  
(37) 
and spatial Fourier transforms:
(38)  
(39)  
(40)  
(41) 
Hence clearly:
(42)  
(43) 
where
(44)  
(45)  
(46) 
and
(47)  
(48) 
Clearly
(49)  
(50)  
(51) 
Therefore one has a freedom to choose (gauge transformation) whereas is the gauge invariant part of the vector potential.
We choose the Coulomb gauge:
(52) 
Hence
(53)  
(54)  
(55) 
and only perpendicular components of electric and magnetic fields are responsible for emission of radiation. The important equation in this case is the fourth Maxwell equation:
(56) 
Since the lhs represents the vector of type therefore:
(57) 
and
(58) 
Substituting the potential :
(59)  
(60) 
where
(61) 
Therefore in the Coulomb gauge
(62) 
and in the far zone