Fourier transform of single particle wave functions in extremely deformed nuclei: Towards the momentum distribution of scission neutrons

. The Fourier transform of single particle wave functions in cylindrical coordinates is applied to the study of neutrons released during scission. We propagate the neutron wave packets in time through the bi-dimensional time dependent Schrödinger equation with time dependent potential. We separate the parts of these wave packets that are in the continuum and calculate their Fourier transforms at di ﬀ erent times: immediately after scission (T = 1 × 10 − 22 s) and at several intervals afterwards (until T = 50 × 10 − 22 s). The momentum distributions corresponding to these Fourier transforms are then estimated. The evolution of these distributions in time provides an insight into the separation of the neutron from the ﬁssioning system and asymptotically gives the kinetic energy spectrum of that particular neutron.


Introduction
The Fourier transform (FT) appears in many nuclear physics applications. For example, FT is required in theoretical approaches of direct nuclear reactions of stripping, pick-up and knock-out types. In particular, the FTs of the single-particle wave-function (WF) represent the momentum transfers that occur during the reaction and give the probability that such an event takes place. Another application of the FT of a WF is in the field of nuclear fission. The momentum distribution of the scission neutrons (those that are emitted during the neck rupture) is given by the FT of the tail of their wave packets that is in the continuum. From the momentum distribution one can deduce an important observable: the kinetic energy distribution. This may allow us to distinguish the scission neutrons from the neutrons evaporated from fully accelerated fragments. Although the nuclei involved in the direct reactions are often deformed, the calculations have been performed in spherical coordinates neglecting the, undoubtedly important, role of the deformation. This is mainly because a deformed WF is commonly expressed in cylindrical coordinates and the procedure to calculate its FT is quite complicated. Recently, we have developed adequate numerical procedures to calculate the FT of a WF in cylindrical (ρ, z) coordinates [1]. After a short description of the numerical methods used, we apply them to scission neutrons and present some physical results in the frame of the Dynamical Scission Model [2,3].

Continuous case
The standard form of the Fourier transform of a function Ψ( r) in three dimensions is: It gives the probability in the momentum space i.e., the probability that the nucleon has its momentum p in the volume element d 3 p.
In our applications we consider nuclear forms with axial symmetry so that we use the cylindrical coordinates (ρ, z) and omit the dependence on the angle θ. It can be shown that the Fourier transform in this case has the form: where J 0 is the zero-order Bessel function of the first kind. The transform with respect to the variables ρ, R is called the zero-order Hankel transform. Thus, the Fourier transform in cylindrical coordinates implies a combination of Hankel and one-dimensional Fourier transforms. The form used in quantum mechanics is based to the relations R = k ρ /(2π) and Z = k z /(2π). The variables ρ, z belong to the position space, while k ρ , k z belong to the momentum space. The factor 1/2π is introduced as a phase convention.

Discrete grid
Suppose now that the function f (ρ, z) is known only on the nodes of a discrete grid. Then, the Fourier transform will be also a discrete function, approximation of the continuous transform. Let ρ j = ρ 0 + jΔρ, j = 0, 1, . . . , M − 1 and z k = z 0 + kΔz, k = 0, 1, . . . , N − 1 be the points (uniformly spaced) defining the grid for which f (ρ j , z k ) are given. The function f (ρ, z) is supposed to be zero or negligible for ρ > ρ M−1 , z < z 0 and z > z N−1 . The Fourier transform will be also calculated on a grid with the mesh points R m , Z n . For any R m from the discrete grid, let us consider the integral: The integral on the interval [ρ 0 , ρ M−1 ] is represented as a sum of integrals on subintervals of the size 2Δρ. To approximate these partial integrals we have deduced a quadrature formula (see [1]) whose coefficients are expressed in terms of Bessel and Struve functions [4] of order 0 and 1.
Let us now consider the integral: The variable Z takes a set of N values defined as The value of F(R m , Z n ) is approximated by using the Fast Fourier Transform -FFT [5]. In conclusion, the values F(R m , Z n ), representing the Discrete Fourier Transform in cylindrical coordinates, are obtained by the discrete Hankel (ρ-axis) and Fourier transforms (z-axis).

Application to the scission neutrons
We apply the Fourier transform of single particle wave functions in cylindrical coordinates to the study of neutrons released during scission. The neutron wave packets are propagated through the bi-dimensional time dependent Schrödinger equation (TDSE) with time dependent potential (TDP): Δ is the Laplacean, V is the potential, C is a constant and the operators S a , . . . , S d represent the spinorbit coupling. We use a realistic mean field potential of Woods-Saxon type adapted to nuclear shapes described by generalized Cassini ovals [6]. The wave functions have two components corresponding to spin "up" and spin "down". For a given projection Ω of the total angular momentum on the symmetry axis they assume the form: Due to the axial symmetry, only the coordinates ρ and z remain.

Dynamical scission model
According to the Dynamical Scission Model [2,3], we separate the calculation in two stages: 1) The scission process itself, i.e., the neck rupture and its absorption by the fragments. The corresponding nuclear configurations are defined by a set of deformations {α i } (just before scission) and {α f } (immediately after scission). The duration of this stage is relatively short (ΔT ≈ 10 −22 s) and therefore the potential in which the neutrons move changes rapidly.
2) The detachment from the fragments of the fraction of the neutrons that are left unbound at the end of the previous stage. Since, at this stage, their motion is much faster than that of the just separated fragments we maintain the fragments at the configuration {α f } and therefore keep the potential constant during the subsequent propagation in time.
In the 1 st stage we follow the time evolution of the wave functions corresponding to each neutron state that is occupied at {α i }. At each time step the nuclear shape and the potential is changing as the deformation parameter set varies uniformly from {α i } to {α f }. The initial wavefunctions (at t = 0) are taken as eigensolutions of the eigenvalue problem: HΨ = EΨ. To solve it, the infinite physical domain is replaced by a finite one, which is discretized on a grid with the mesh points ρ j and z k . At each point the partial derivatives in H are approximated by finite difference formulas [7]. The eigenvalue problem is transformed into an algebraic eigenvalue problem with a large sparse matrix, which is solved by an appropiate method (Arnoldi) [8]. The propagation in time is performed by a numerical scheme of Crank-Nicolson type [9]. At the end of this stage we calculate quantities like: • The probability that a neutron occupying the state |Ψ i before scission populates a state |Ψ f : • The probability that the neutron is emitted, given by: where v 2 i are the occupation probabilities of the initial eigenstates. • The emitted part of the wave packet: In the 2nd stage we follow these wave packets (that describe the unbound neutrons) in time (by TDSE with the potential fixed for the set {α f }) and calculate their FT at certain time intervals T .

The evolution of the emitted wave functions and of the momentum distributions
We have investigated the reaction 235 U(n th , f ) with the light fragment mass A L = 96 (which is the most probable experimental mass asymmetry). Two parameters of the Cassini representation have been considered: α (elongation) and α 1 (mass asymmetry). We have chosen: α = 0.985 (before scission) and α = 1.001 (after scission). The parameter α 1 results from the condition of fixed ratio of the fragment volumes that corresponds to A L =96. We have performed the Fourier transform for several wavefunctions corresponding to Ω = 1/2. At each point in the (k ρ , k z ) plane are associated an absolute value K = k 2 ρ + k 2 z and a probability P = |F(k ρ , k z )| 2 k ρ Δk ρ Δk z . To represent the Kdistributions as histograms, we divide the domain of K-values in equal intervals and group the points according to the interval to which they belong. Summing up the probabilities of the points in each group one obtains the probability P(K) that a neutron has its K-value in the respective interval. In the Figs.1 and 2 are shown the emitted wavefunctions with the indices 30 and 33 (as sum of square moduli of the two components) aside with the momentum histograms at different times T . The timedependent wavefunctions are represented relative to the function at T = 20, while the values on the ordinates of the histograms are the P(K) probabilities multiplied by 10 2 .
As one can see, at increasing time intervals the wavefunction values are diminishing, showing that the neutron is leaving the nucleus. At the same time, the K-distributions are shrinked to lower values, reflecting the fact that the neutrons are less and less present in the potential well. At large times they tend to a gaussian shape.
From the momentum distribution one can deduce the kinetic energy distribution (according to E = 2 2μ K 2 ). To obtain the whole kinetic energy spectrum one has to repeat such calculations for all the states at each Ω value. Here we present results only for Ω=1/2 states and from WF 11 to WF 41 . The states lying at the bottom of the well (WF 1 to WF 10 ) have negligible P em values. Note however that Ω = 1/2 gives the dominant contribution (65%). Figure 3 shows kinetic energy spectra for two individual states (30 and 33) around the Fermi level and the sum for all Ω=1/2 states with the index above 10 (having energies greater than −23 MeV). The total spectrum presents a maximum around 0.7 MeV and an exponentially decreasing tail till 5 MeV in qualitative agreement with the measured spectrum [10] of all prompt fission neutrons (PFN). For a quantitative comparison with experimental data, states with higher Ω values (3/2 to 9/2) has to be included in the calculations.

Summary
In this work we have presented numerical procedures and some applications of the Fourier transform on wave functions given in cylindrical coordinates. As initial functions for our procedure we have chosen single-particle states in extremely deformed nuclei associated with the wave packets of the neutrons released during scission. We propagate the neutron wave packets in time through the bi-dimensional TDSE with TDP. We select Ω=1/2 states and calculate their Fourier transforms at different times T and the momentum distributions corresponding to these transforms. From such distributions one can follow the separation of the neutron from the fissioning system and obtain the kinetic energy spectrum of each emitted neutron. The shape of the energy spectrum obtained summing over all emitted neutrons resembles the measured PFN spectrum. For a more accurate comparison the effect of reabsorption of the unbound neutrons by the imaginary potential has to be taken into account.