Nuclear Fission Dynamics: Past, Present, Needs, and Future

Significant progress in the understanding of the fission process within a microscopic framework has been recently reported. Even though the complete description of this important nuclear reaction remains a computationally demanding task, recent developments in theoretical modeling and computational power have brought current microscopic simulations to the point where they can provide guidance and constraints to phenomenological models, without making recourse to parameters. An accurate treatment compatible with our understanding of the inter-nucleon interactions should be able to describe the real-time dynamics of the fissioning system and could justify or rule out assumptions and approximations incompatible with the underlying universally accepted quantum-mechanical framework. Of particular importance are applications to observables that cannot be directly measured in experimental setups (such as the angular momentum distribution of the fission fragments, or the excitation energy sharing between the fission fragments, or fission of nuclei formed during the r-process), and their dependence of the excitation energy in the fissioning system. Even if accurate predictions are not within reach, being able to extract the trends with increasing excitation energy is important in various applications. The most advanced microscopic simulations of the fission process do not support the widely used assumption of adiabaticity of the large amplitude collective motion in fission, in particular for trajectories from the outer saddle toward the scission configuration. Hence, the collective potential energy surface and inertia tensor, which are the essential elements of many simplified microscopic theoretical approaches, become irrelevant. In reality, the dynamics of the fissioning system is slower than in the case of pure adiabatic motion by a factor of three to four times and is strongly overdamped. The fission fragment properties are defined only after the full separation, while in most of the current approaches no full separation can be achieved, which increases the uncertainties in describing fission-related observables in such methods.


THE PAST
In a matter of days after Hahn and Strassmann [1] communicated their yet unpublished results to Lise Meitner, she and her nephew Otto Frisch [2] understood that an unexpected and qualitatively new type of nuclear reaction has been put in evidence and they dubbed it nuclear fission, in analogy to cell divisions in biology. Until that moment in time nuclear fission was considered a totally unthinkable process [3,4], "as excluded by the small penetrability of the Coulomb barrier [5], indicated by the Gamov's theory of alpha-decay" [2]. Meitner and Frisch [2] also gave the correct physical interpretation of the nuclear fission mechanism. They understood that Bohr's compound nucleus [6] is formed after the absorption of a neutron, which eventually slowly evolves in shape, while the volume remains constant, and that the competition between the surface energy of a nucleus and its Coulomb energy leads to the eventual scission. Meitner and Frisch [2] also correctly estimated the total energy released in this process to be about 200 MeV. A few months later Bohr and Wheeler [7] filled in all the technical details and the long road to developing a microscopic theory of nuclear fission ensued. In the years since, a few more crucial theoretical results have been firmly established: (i) the defining role of quantum shell effects [8,9] and in particular the special role played by the pairing type of nucleon-nucleon interaction in shape evolution [10,11]; (ii) the fact that the subsequent emission of neutrons and gammas can be described quite accurately using statistical methods [12,13]; (iii) and that the non-relativistic Schrödinger equation should be adequate as well, as no genuine relativistic effects, such a retardation, are expected to play any noticeable role in fission dynamics.
Whether a fissioning nucleus undergoes either spontaneous fission or induced fission the time it takes the nucleus to evolve from its ground state shape until outside the barrier or past the outer saddle is very long in case of neutron induced fission ≈ O(10 −15 ) s. in comparison with the time the nucleus slides downhill until scission, which is estimated to be O(10 −20 ) s. Therefore, the saddle-to-scission stage of nuclear shape evolution is the fastest and arguably the most nonequilibrium stage of the nucleus dynamics from the moment a neutron has been captured and the stage which plays a crucial role in in determining the FFs properties. The many intricacies of the fission process, the multitude of aspects, which required a deep theoretical understanding, ultimately rooted in the quantum nature of this phenomenon defied the efforts of many generations of theorists, and a huge plethora of mostly phenomenological models have been put forward, often based on contradictory assumptions between models. The extensive range of assumptions, range from adiabatic evolution on top of which one adds (relatively weak) dissipation and fluctuations, to strongly overdamped motion, when the role of the collective inertia becomes irrelevant, to full statistical equilibrium near the scission configuration, and to mixing quantum and classical descriptions.
Microscopically inspired models are typically based on the (ill suited choice of words, as the our analysis shows) adiabatic approximation, which is often conflated with slow motion. The class of adiabatic transformations, during which only mechanical work is performed and no heat transfer or entropy production occurs, are only a subclass of slow motion or quasistatic processes. Theorists believed that the nuclear shape evolution until the moment of scission was so slow that individual nucleons had a sufficient time to adapt to avoided single-particle level crossings [14] and the entire nucleus would follow the lowest "molecular term, " using the FIGURE 1 | The schematic evolution of the single-particle nucleons levels (Upper panel) and of the total nuclear energy (Lower panel) as a function of deformation parameter q [10,16,17]. The thick line represents the Fermi level and the up/down arrows depict the Cooper pairs of nucleons on the Fermi level only, in time-reversed orbits (m, −m). This figure is reproduced from Bulgac et al. [17] under the terms of the American Physical Society copyright agreement.
Born-Oppenheimer chemical terminology [15] (see Figure 1). Following this assumption at first the generator coordinate method (GCM) has been introduced by Hill and Wheeler [14] and Griffin and Wheeler [18] and later on a related alternative approach, the adiabatic time-dependent Hartree-Fock (ATDHF) method [19][20][21][22]. GCM is still one of the most popular tools still in use in the microscopic theories of fission [22][23][24][25][26][27][28][29]. The GCM and ATDHF method have been shown to be basically equivalent [30], when GCM is defined with complex generator coordinates [31]. As however Goeke and Reinhard [30] succinctly state: "Usually the |q is obtained by an educated guess using the preconceived knowledge of the process." (Typically |q stands for a generalized Slater determinant, aka Hartree-Fock-Bogoliubov many-nucleon wave function.) Even though many efforts have been dedicated to find a better way to choose the collective or generator coordinates the methods proved to be too difficult to implement in practice and the quality of the decoupling between collective and intrinsic degrees of freedom either too not very good or difficult to assess [32]. An exact separation between collective and intrinsic degrees of freedom (DoF) it is equivalent to an adiabatic evolution of the set of collective DoF. Then the collective DoF would always follow the lowest "molecular orbital" and only work would be performed on the intrinsic DoF, and thus with no heat transfer, and the intrinsic system would remain "cold" during the entire evolution.

THE PRESENT
Unfortunately until recently this crucial aspect of large amplitude collective nuclear dynamics, whether the large amplitude collective motion in fission is indeed adiabatic was never tested and, as our results unequivocally show, the adiabatic assumption is strongly violated. Surprisingly at first glance, we have recently shown that the evolution of the nuclear shape is in reality significantly slower that the adiabatic assumption would predict. One would naively expect that if the motion is even slower then at an avoided level crossing (see Figure 1), the probability that the system would follow the lower "molecular term" is even greater and thus the adiabatic assumption would be even more likely be valid. An analogy with a classical system can help and demonstrate just the opposite. If a railroad car is released on top of a hill, it will convert basically all the gravitational potential energy difference when reaching the bottom of the hill into kinetic energy. Thus, only pure mechanical work on the intrinsic DoF of the railroad car will occur with essentially no heat transfer or entropy production. This is an adiabatic process. However, if one were instead to block the wheels of the railroad car, the friction will slow down the car and almost the entire gravitational potential energy difference will be converted into heat, the wheels will get red hot, thus increasing the "intrinsic energy" of the car, and the speed of the car at the bottom of the hill would be rather small.
The only practical theoretical framework to consider in the treatment of the dynamics of large nuclei is the (Time-Dependent) Density Function Theory [(TD)DFT], which has been formulated a long time ago [33][34][35][36][37][38]. One of the main difficulties consist in constructing the energy density functional, for which no rigorous recipes exist, DFT and TDDFT arrive at the mathematical conclusion that the stationary or timedependent solution of the many-body Schrödinger equation is in a one-to-one-to-one correspondence with the number density, an (arbitrary) applied one-body external potential, and that the number density can be obtained by solving the much simpler DFT or TDDFT equations. There is a continual debate in nuclear physics that DFT is not applicable to nuclei, which are self-bound isolated systems. At the same time however, no one would argue that DFT cannot describe neutrons and protons in the neutron star crust and deeper into the star. Neutrons are delocalized in the neutron star crust and below it. In the rod, slab, and tube phases in the neutron star crust and below both protons and neutrons are delocalized and in this respect they are similar to electrons in solids. One can imagine that in the future one might produce a nuclear trap using some kind of γ -lasers, similarly to what nowadays experiments are made with cold atoms. Until then one can mentally imagine that one can put an isolated nucleus in a spherical infinite square well-potential with a radius about 3. . . 5× the nuclear radius (or even a harmonic potential) and compare the results of such a DFT treatment of the nucleus with the widely accepted DFT alter ego, the NEDF approach, a tool of choice in theoretical nuclear calculations. The results of the these two approaches are numerically indistinguishable under these conditions, and therefore the debate alluded above is merely pedantic. One deficiency of a pure DFT approach is that the number density alone cannot disentangle between a normal and a superfluid system, and one needs an order parameter, as one does in the case of magnetization. The practical local density approximation (LDA) [34], which is the local formulation of DFT, has to be augmented with the anomalous density [39,40], and it was dubbed the superfluid LDA (SLDA). In Bulgac [39,40] one can find detailed reviews of the developments, verification, and validation of TDSLDA for a variety of physical systems, ranging from cold atoms, nuclei, and to neutron star crust.
A (TD)DFT framework for nuclear structure and dynamics should satisfy several requirements (in this order of importance): (i) the DFT and the Schrödinger description of observables should be identical, as both in ultimate instance rely on the same inter particle interactions; (ii) both DFT and Schrödinger equations should describe correctly Nature, thus we need accurate interactions between nucleons; (iii) the numerical implementation of the (TD)DFT should faithfully reproduce the theory. At present we definitely do not have acceptable answers to the requirements (i) and (ii) and rely instead to a significant amount of phenomenology. The numerical implementation of DFT without any drastic physical restrictions became possible only relatively recently, with the advent of supercomputers [39,40].
The TDSLDA is formulated in terms of Bogoliubov quasiparticle wave functions (qpwfs). The evolution of the qpwfs is governed by the equations: where we have suppressed the spatial r and time coordinate t, and k labels the qpwfs (including the isospin) [u kσ (r, t), v kσ (r, t)], with σ =↑, ↓ the z-projection of the nucleon spin. The single-particle (sp) Hamiltonian h σ σ ′ (r, t), and the pairing field (r, t) are functionals of various neutron and proton densities, which are computed from the qpwfs, see Jin et al. [41] for technical details. No proton-neutron pairing is assumed in the present study, and the pairing field is singlet in character. A TDSLDA extension to a more complex pairing mechanisms is straightforward. While a nuclear system evolves in time one can uniquely separate the energy into collective kinetic energy and intrinsic energy contributions [17] using the nuclear energy density functional (NEDF) E τ (r, t), n(r, t), ... , in a similar manner as in hydrodynamics: Above n(r, t) is the number density, τ (r, t) is the kinetic density energy, and p(r, t) = mn(r, t)v(r, t) are linear momentum and local collective/hydrodynamic velocity densities, and ellipses stand for various other densities. p(r, t)/n(r, t) is the position of the center of the local Fermi sphere in momentum space. The first term in Equation (2) is the collective/hydrodynamic energy flow E coll and the second term is the intrinsic energy E int in the local rest frame. For the sake of simplicity we have suppressed the spin and isospin DoF, even though they are included in all the actual calculations. The collective energy E coll (t) is not vanishing only in the presence of currents and vanishes exactly for stationary states. The inertia tensor in E coll (t) in the case of irrotational collective motion is fully equivalent to the Werner-Wheeler inertial tensor [23]. The intrinsic energy E int (t) is determined only by the fermionic matter distribution. The qualitative new result established in Bulgac et al. [17] and also illustrated here in Figure 2 is that in an fully unrestricted TDSLDA the collective flow energy is almost negligible until scission, in total discrepancy with what one would have naively expected if the adiabatic assumptions would be satisfied.
Our simulations point to an unexpectedly small E coll from saddle-to-scission, corresponding to a collective speed v coll /c ≈ 0.002 · · · 0.004, significantly smaller than the Fermi velocity v F /c ≈ 0.25 (see Figure 2). Since in TDSLDA one simulates the one-body dynamics exactly, it is natural to discuss adiabaticity at the mean-field level. The transition rate between sp states is suppressed if the time to cross an avoided level-crossing configuration satisfies the restriction t ≪h/ ǫ ≈ 400 fm/c, where ǫ = 1/ρ sp (ǫ F ) is the average sp energy level spacing at the Fermi level. Since on the way from saddle-toscission the time required is 1 . . . 3 × 10 3 fm/c and several dozen of avoided level crossings occur [16,43], this condition is clearly violated. Somewhat surprisingly, the adiabatic assumption is also violated even in the case of SLy4 NEDF (see Bulgac et al. [42] and Figure 2), when the saddle-to-scission time is O(10 4 ) fm/c as well. The collective motion is thus expected to be strongly overdamped. From saddle-to-scission the nucleus behaves as a very viscous fluid, the role of collective inertia is strongly suppressed, and the trajectories follow predominantly the direction of the steepest descent with the terminal velocity determined by the balance between the friction and the driving conservative forces (see Figure 2).
This result serves as the first microscopic justification for the assumption of the overdamped Brownian motion model [44][45][46][47][48][49] and partially to the scission-point model [50][51][52][53]. In both these phenomenological models it is assumed that the preformed FFs are in statistical equilibrium and that the collective energy flow is either vanishing or very small. The main difference is that in the scission-point model there is no mechanism to ensure that all equilibrium scission configurations could be reached dynamically, while the nucleus evolves from the saddleto-scission. Moreover, the relaxed FF properties are defined only after the FFs become sufficiently well-separated, see below. It is equally unexpected that in the case of enhanced pairing, when the pairing condensates retain their long-range order throughout the entire saddle-to-scission evolution, the collective dynamics remains strongly overdamped.
While evolving from saddle-to-scission a nucleus encounters a large number of avoided level crossings and instead of following the lowest potential energy surface, as would happen in an adiabatic evolution, many transitions to higher excited levels  occur. Similarly to what is known for decades in chemistry [54][55][56], one should consider not only the lowest potential energy surface, when initially the system could be found on the lowest with unit probability, but all "molecular terms, " which become populated during the evolution (see Figure 3). The separation between such potential energy surfaces reaches values of order O(10) keV or even less, and only if the system traverses a level crossing in a time much longer than ≈h/10 keV ≈ 20, 000 fm/c or longer the nucleus will remain on the lowest potential energy surface. On the other hand the saddle-to-scission time is O(10 3 ) fm/c and this is why the adiabatic assumption is in the final analysis strongly violated. Why do pairing correlations play an important role in fission dynamics? In a static nuclear configuration each single-particle level is doubly-degenerate, due to the Kramers degeneracy. At the level crossing both nucleons on the highest occupied level (homo-highest occupied molecular orbital) would have to transition simultaneously to the lowest down-sloping level (lumo-lowest unoccupied molecular orbital) to ensure that the local momentum distribution remains approximately spherical, as otherwise it would acquire an oblate shape [10,11,14], while the shape of the nucleus becomes more prolate. Nucleon-nucleon interactions at low momentum transfer can be modeled with a reasonable accuracy with a zero-range δ-interaction, which favors the transitions between pairs of time-reversed orbitals, exactly as the Kramers degenerate orbitals (see Figure 1). The up-sloping levels are characterized by larger projections of the angular momentum on the fission axis, |m| ≈ k F r 0 A 1/3 , and these levels should be depopulated, since in a FF the largest angular momenta are smaller, ≈ k F r 0 (A/2) 1/3 . While evolving from one level crossing to the next, the entire evolution is likely rather well-reproduced by a simple one-body dynamics, as each single-particle level occupation probability changes little. What one-body dynamics lacks is the contributions arising from the Boltzmann collision integral. However, at each level crossing the two-correlated nucleon pairs will undergo a collision, and at low energies transitions between pairs of time-reversed orbitals expected to dominate the collision rate. One should take with a grain of salt this simplistic picture of "collisions" and jumps between sp levels, as nothing happens instantaneously or at one point in space in quantum mechanics. In the presence of a Bose-Einstein condensate of nucleon Cooper pairs the nucleus has a superfluid component and pair transfers are enhanced due to the Bose enhancement factor. The dynamics of the nuclear systems then approaches the evolution of classical inviscid (no viscosity) or perfect fluid. An illustration of this behavior was exemplified in Figure 4 in Bulgac [40] and in Figure 4. When the magnitude of the pairing field was artificially increased from a realistic value to a value 3. . . 4× larger the evolution time from saddle-to-scission decreased by a factor of ≈ 10 and at the same time the long range coherence of the pairing field across the entire nucleus survived. For realistic values of pairing strengths during the descent from the saddle-to-scission both proton and neutron pairing fields fluctuate both in space and time, long range order basically vanishes, but quite often it is revived.
Another important aspect which emerged from our TDSLDA fission simulations [17,42], which is of significant importance in the implementation of various phenomenological models, concerns the shapes of FFs at scission and when they reach their relaxed shapes (see Figure 5). This behavior is apparently confirmed indirectly by experiments. In Langevin or Fokker-Planck [60][61][62][63][64][65], TDGCM [26,66], and scission-point [50][51][52][53] models the calculation of the FFs yields consider only a very limited range of nuclear shapes. In particular in such simulations one never introduces the octupole FF moments. Our results, as well recent analysis by Scamps and Simenel [67], clearly demonstrate that the FFs emerge at scission octupole deformed and also with a significantly larger quadrupole deformation than the relaxed values. Moreover, even after scission the FFs there is a significant Coulomb interaction between them, which leads to the excitation of both low energy and giant resonances in FF [68,69]. This interaction enables additional excitation energy exchange between the FFs after scission, and it also affects their total kinetic energy, a behavior also seen in our simulation, but yet not documented. In statistical scission-point models there is no dynamics, and only the competition between FFs configurations at the scission point are considered [50][51][52][53], a model to which our results lend partial support. However, the only shapes considered are quadrupole deformation of the relaxed FFs, which clearly is not what our dynamical simulation demonstrate.
Our simulations put in evidence another very important aspect, the mechanism of the excitation energy sharing between the FFs. As a rule the heavy FF emerges in the end cooler than the light FF, even though they have been in contact for quite a long time before scission. Moreover, when increasing the initial energy excitation of the fissioning nucleus we have established that only the heavy FF becomes hotter and that is reflected in the average neutron multiplicity number of emitted, and results which is in apparent agreement with experimental findings (see Figure 6). Experimentally it is extremely difficult to infer the excitation energies of the FFs, which are a crucial input in various statistical codes [70,72,73]. Bertsch et al. [74] argue that the FFs spin distribution, which determines the prompt gamma angular distribution, can be used to infer information of the excitation energy sharing between FFs. Randrup et al. [73] point to a pronounced anti-correlation betweenν(A) and mean total kinetic energyTKE, which can be used reduce uncertainties in data analysis. Schmidt and Jurado [75,76] suggested a phenomenological model, the "energy-sorting" mechanism based on the empirical constant temperature parametrization of the nuclear level densities due to Gilbert and Cameron [77]. In this "energy-sorting" model the FFs before scission have different temperatures, with a lower temperature of the heavy fragment, which would generate an energy flow from the light/hotter to the heavy/cooler fragment. Our simulations demonstrate however that the near scission the two FFs have properties quite different (shape, excitation energy, pairing correlations) from the properties of relaxed fragments. It is therefore problematic to relate the properties of excited isolated nuclei with the properties of FFs in contact before rupture.
As we have mentioned above, the number of collective DoF and their character is a long standing problem in microscopic inspired theoretical and phenomenological models. The choice FIGURE 6 | We compare here the average neutron multiplicityν(A) emitted by FFs using a CGMF simulation [70], which assumes an E n dependence for the energy sharing extracted using the excitation energy sharing between the FFs in our calculation with NEDF SeaLL1, as a function of the equivalent incident neutron energy in 235 U(n,f) reaction along with available experimental data [71]. Note that in this figure the parametrization was based on 240 Pu calculations, while 236 U calculations are in progress. of collective DoF and their character is guided by the authors' intuition, their computational and other abilities, educated guess, imagination, and/or available resources [30]. It was never proven or demonstrated that a GCM representation of the nuclear wave function is ever accurate. The wave function of a manyfermion system in TDGCM is constructed according to the following prescription where (x|q) are (generalized) Slater determinants depending on nucleon spatial coordinates, spin, and isospin x = (x 1 , . . . , x A ), x k = (r k , σ k , τ k ) and parameterized by the collective coordinates q = (q 1 , . . . , q n ), and where f (q, t) is the collective wave function. There is no criterion or small parameter which controls the accuracy of such a representation, that is We have presented rather simple arguments that most likely such an accurate representation does not exist in general, particularly in fission large amplitude collective motion [17]. A rather simple estimate of the possible number of different shapes, and therefore of independent terms in Equation (3), shows that it is basically infinitely small in comparison with the number of possible terms in an exact multi-configurational Slater determinant representation of exact . Can we check whether GCM is a good approximation for the saddle-to-scission evolution? Now we definitely can and we have the answer, and since where T is the temperature of the intrinsic system, we know that the collective motion is strongly overdamped, that there is an irreversible energy flow from the collective/shape DoF to the intrinsic DoF, and that during the saddle-to-scission evolution the temperature of the intrinsic system increases, as does its (entanglement) entropy as well. We should also remember, that the evolution of the fissioning nucleus from saddle-to-scission is a truly non-equilibrium one, and the notion of a "slowly evolving" nuclear temperature T might be questionable. This suggest that in phenomenological calculations a more physically motivated choice for the potential energy surface would be one in which with increasing deformation one would increase the "temperature T, " so as to keep the "collective potential" energy V(q, T) essentially equal to the initial excitation energy of the compound nucleus. Basically this is the prescription that Randrup et al. [48,49] have implemented lately. There were indirect indication in microscopic approaches that the number of collective DoF can vary along the fission path. In the overwhelming majority of fission studies the only deformation of the nucleus close to the ground state is the axial symmetric quadrupole moment Q 20 , in spite of the fact that the Bohr-Mottelson five-dimensional collective Hamiltonian is one tool of choice to describe the low energy structure of openshell nuclei [58,78,79]. Ryssens et al. [80] have shown in a beyond mean field calculation of the 1-dimensional fission path the nucleus 240 Pu is axially deformed around the ground state configuration, it becomes triaxially deformed when it reaches the fission isomer region, and before it reaches the outer fission barrier it breaks the axial symmetry as well. We have recently developed a framework into which one can include fluctuations and dissipation in a quantum approach [81], which allows for shape fluctuations of any kind. What we observed that by allowing a bending mode to become active (see Figure 7), the mass and TKE distributions acquire shapes in quite nice agreement with observations. It was discussed quite some time ago that such bending modes might be responsible for the angular distributions of the FFs [82] and this type of distributions can be also extracted from TDSLDA calculations in the near future [74,83]. Phenomenological or even TDGCM approaches do not consider so far such bending modes, which are definitely physically relevant for a large number of observables.

WHAT ARE THE NEEDS?
Nuclear fission is a complex process, in which a heavy nucleus evolves from a compact shape to a configuration in which two or more fragments are produced, and is accompanied by emission of prompt neutrons and gamma rays (and, eventually, electrons and antineutrinos). Most of the energy is released in form of kinetic energy of the fragments, while prompt neutrons emitted before beta decays play the main role in applications like energy production.
The dynamics of the nuclear system from the formation of compound nucleus until after the acceleration of FFs and prompt particle emission is too rapid to be experimentally resolved. On the other hand, the time scale of the weak interaction, which governs the decay of FFs in experiments toward stability ranges from seconds to minutes, and thus the dynamics involving beta decay toward stability, including delayed neutron and gamma emissions, can be decoupled from the initial more rapid part that is only governed by the strong interactions. Consequently, the fragment properties directly influence the properties of prompt fission and gamma rays, which have been the subject of comprehensive, albeit not exhaustive, experimental investigations over the years. Thus, at Los Alamos National Laboratory, experimental campaigns have investigated the prompt fission neutron spectra for neutroninduced fission of major actinides ( 235,238 U, 239 Pu) for a large range of incident neutron energies using the ChiNu experimental setup [84,85]. On the other hand, experiments that measure the average neutron multiplicity as a function of pre-fission neutron mass [86][87][88][89][90][91], which can be used to guide energy sharing in fission fragments, are more scarce, and for a limited number of reactions (usually spontaneous fission or neutroninduced fission with thermal neutrons). The surrogate reaction technique is used by Lawrence Livermore National Lab and Texas A&M to measure the 239 Pu and 241 Pu prompt fission neutron multiplicity (average and distribution) as a function of equivalent incident neutron energy [92,93]. Complementary, significant resources have been devoted toward measuring the prompt fission gamma rays produced in the decay of fission fragments, after the neutron emission, at Los Alamos, using the DANCE calorimeter [94][95][96][97], and worldwide employing highresolution detectors [98][99][100][101][102]. Such measurements complement existing measurements of prompt-gamma rays by Verbinski et al. [103], and many others [104][105][106]. Furthermore, data on total gamma production can be useful in evaluating the prompt fission gamma properties, if no other data exist, fission dominates, and the other gamma-producing channels can be modeled with reasonable accuracy [107].
In addition to detecting and measuring the products of decay of FFs, a concerted effort has been directed toward the direct measurement of the FF distributions. While simpler, the "2E" measurements [108][109][110][111] have the shortcoming that additional information is required to identify the mass of the fragment, hence the poor 4-6 amu resolution. Reducing such large uncertainties was the main argument for the construction of the SPIDER "2E-2v" experimental setup [112,113] at the Los Alamos Neutron Science Center, which should start taking data in the near future. In this setup the two fragment kinetic energies and velocities are measured simultaneously and can achieve one mass unit resolution, as also demonstrated by similar setups, COSI FAN TUTTE at ILL, France and VERDI at JRC-Geel, Belgium. The Sofia experiment [114] at GSI has produced an extensive amount of experimental data with accuracies in mass and charge <1 mass unit, but suffer from the fact that a large range of excitation energies that cannot be disentangled is produced in Coulomb excitations. LLNL and LANL activation measurements at TUNL have provided invaluable data on the incident energy dependence of cumulative fission product yields [93,115], and ongoing efforts are directed toward testing the Bohr hypothesis of independence between entrance and outgoing channels in a compound reaction by comparing fission yields in the 239 Pu(n,f) and 240 Pu(γ ,f) reactions [116] that lead to the same fissioning system, at least in mass and charge. At CARIBOU, two projects are underway: (i) measurement of fission product properties (isomeric yield ratios, gammaray decay branching ratios, and β-delayed neutron emission properties) and (ii) improvements of the antineutrino spectrum simulations by performing measurements of the β-decay data using GAMASPHERE [117]. Other experiments related to the properties of FFs have concentrated on measuring the average total kinetic energy (TKE) of the fission fragments at the LANSCE WNR facility for energy from a few hundred keVs to 200 MeV [111,113]. The High Rigidity Spectrometer [118] at FRIB [119] will allow to explore the fissioning of very neutron-rich nuclei, study their shell structure, the presence of superheavy, the limits of stability and the equation of state of neutron-rich nuclear matter, the neutron skins, and shed light on the nature of the r-process, so far an unknown territory.
Treating the evolution of a heavy nucleus from a compact configuration until the start of beta decays in the FFs, including prompt neutron and gamma emissions, is a complicated task, computationally unfeasible within a unified microscopic approach. Within this reality, one has to consider a mixture of approaches, in which the initial part of the fission process is treated within a microscopic framework that can inform more phenomenological treatments that model the emission of prompt neutrons and gamma rays, and whose results can be directly compared against a large set of experimental data. And the input from microscopic models does not necessarily need to be restricted to the fission process. The systematics of several physical quantities used in the phenomenological models of neutron and gamma emission is based on data for stable nuclei. Since the FFs are nuclei far from stability, it is important to investigate within more microscopic models the validity of various systematics far from stability, as prompt neutron and gamma observables are sensitive to ingredients like level densities, optical models, or gamma strength functions.
Several models have been proposed to describe the shape dynamics, and to some extent many of them are able to reproduce experimental quantities like the pre-neutron emission mass distributions, irrespective of the approximations involved. But, as noted above, the direct measurement of the FFs before neutron emission is not possible. Therefore, even before one goes into details regarding the validity of the approximations involved in theoretical models, one needs to consider the corrections involved in the analysis of the experimental data. The post-neutron emission FF mass distributions are obtained using information on average prompt neutron multiplicity as a function of the fragment mass,ν(A), either from measurements, where available, or from phenomenological models. Given the scarcity of the data, and the fact that these phenomenological models have been built in a systematics based on little data, the uncertainties arising fromν(A) could be considerable. Another assumption is that no neutrons are emitted during scission or at neck rupture. It is conceivable that the number of neutrons emitted during the scission dynamics is not necessarily large, but how large is large? If the fraction of scission neutrons is significant, as some phenomenological models predict, citing the experimental prompt fission spectrum as evidence to support the claim [120][121][122][123][124][125], what is the impact on the experimental analysis? While these models are not universally accepted, a microscopic approach should be able to answer such a question, and the assumption of small numbers of neutrons emitted during scission should not be implicitly included in the model. Moreover, it would be rather impossible to asses the fraction of scission neutrons if the model does not follow the dynamics of the process until full separation.
All the codes that model the prompt neutron and gamma emissions have been built in the assumption that the neutron emission proceeds after the full acceleration of the fission fragments [70,[126][127][128][129][130]. This assumption can have significant consequences, as at full acceleration the neutrons are maximally boosted. If the emission occurs instead during the acceleration, then the prompt neutron spectrum in the lab frame can noticeably be altered.
In addition to validating assumptions in phenomenological codes, a reliable theoretical model should provide information on observables that cannot be measured, but are essential in modeling the neutron and gamma emission. Thus, while the total excitation energy available in the fragments can be determined from the Q-value and TKE, no additional constrains on how this energy is shared between the FFs are available. As the most efficient way to lower the energy in FFs is via neutron emission,ν(A) can be used to parameterize the excitation energy sharing. However, such data are scarcely available, usually for a limited number of spontaneous and neutron-induced fission with thermal neutrons, and thus guidance from microscopic models on how the excitation energy is shared with increasing incident neutron energy is required not only for minor actinides, but even for induced fission reaction of major actinides 235,238 U and 239 Pu isotopes.
Gamma emission competes with the neutron emission when the excitation energy in the nucleus is around the neutron separation energy. The spin of the fragments, which cannot be directly measured experimentally, determines the exact strength on the neutron-gamma competition. Only indirect information can be extracted regarding the spin distribution of the FFs from properties of prompt gamma rays, like average multiplicity, prompt gamma fission spectra [131], or isomeric ratios [132]. The measurements of gamma rays are also scarce with increasing the incident neutron energy. Experimentally, one observes an increase of the total prompt fission gamma energy released in fission [133], which has been interpreted as an increase of the average spin with the incident energy. However, only little experimental data exist and no microscopic model has been able to reproduce such trends.
The main message of this section is that experimentally one cannot isolate and measure properties of post-scission FFs, and any full model will have to include simulations of prompt neutron and gamma emission, in order to compare and validate against experimental data. Given that many quantities used in the modeling of prompt neutron and gamma emission are taken from data systematics, which in general is available only for nuclei close to stability, the reliability of the theoretical model near scission and beyond, until full separation, is very important. Such models should allow for a full separation of the fragments, and employ approximations that are validated and under control, as to reduce any uncertainties regarding the FF properties.

WHAT LESSONS HAVE WE LEARNED SO FAR AND WHAT IS THE MOST LIKELY PATH TO THE FUTURE?
While pairing is not the engine driving the fission dynamics, pairing provides the essential lubricant, without which the evolution may arrive rather quickly to a screeching halt [134][135][136]. So far we have not considered whether proton-neutron pairing might have a role in fission. It is very unlikely that a condensate of proton-neutron pairs exists in heavy nuclei, but as we have learned from our simulations, proton-neutron transitions with L = 0 between single-particle orbitals could be important as the neutron-neutron and proton-proton transitions with (S,L,T) = (0,0,1) at low excitation energies, even in the absence of a condensate of such pairs, see the discussion concerning Figures 1, 4.
TDSLDA framework for fission dynamics, while it does not incorporate fluctuations, has provided a lot of insight into the real quantum dynamics, and it revealed extremely valuable information into nuclear processes and quantities, which are either not easy or impossible to obtain in laboratory or observations: FFs excitation energies and angular momenta distributions prior to neutron and gamma emission, element formation in astrophysical environments, as well as other nuclear reactions in a parameter free approach. In particular, the excitation energy sharing mechanism between FFs and its evolution with the initial excitation energy of the compound nucleus was not accessible until now within a dynamic approach. Fluctuations, which are essential in order to reproduce mass and charge yields for example, can be now incorporated into a pure quantum framework [81].
The quality of the agreement with experimental data is surprisingly good, especially taking into account the fact that no attempt was made to reproduce any fission data. Basically all phenomenological NEDFs satisfy the most important requirements to describe the gross properties of nuclear fission: saturation, realistic surface tension and symmetry energy, Coulomb energy, realistic pairing and shell corrections energy. Nevertheless, the quality of existing NEDFs needs improving. One can make a strong argument that we have now a clear path from more phenomenology and adjusted parameters to more fundamental theory and increased predictive power [137].
Perhaps the most important aspect we have observed in all our simulations is the strong violation of the adiabatic assumption in fission large amplitude collective motion. Basically since the 1950's the adiabatic assumption was the main simplification included in all microscopic frameworks, GCM, ATDHF and a large majority of phenomenological models as well, such as the Langevin and the Fokker-Planck equations, where one needs to introduce a potential energy surface and an inertia tensor in the space of the collective variables. If the collective motion is overdamped, the inertia tensor becomes irrelevant, and moreover, considering only the lowest potential energy surface is physically unacceptable (see Figures 3, 2) and the corresponding discussion in section 2. As we have established in Bulgac et al. [81] the fluctuations, or equivalently the role of two-body collisions, does not affect this conclusion.
A somewhat unexpected result was the character of the energy sharing mechanism between the fission fragments, the fact that the heavy fragment is cooler than the light fragment and it has less excitation energy [17,42]. And this conclusion is not a result of the fact that the heavy FF is closer to the double magic 132 Sn nucleus. In the original experiment of Hahn and Strassmann [1] the heavy fragment had a charge closer to Z ≈ 52 − 56, a fact recognized also by Meitner and Frisch [2,138], and explained by Scamps and Simenel [67], this is due to a stabilization of the octupole deformation in FFs, also observed in our simulations [17,42]. With increasing excitation energy of the fissioning nucleus the heavy FF appears to be the only one who absorb the increase and emits more neutrons (see Figure 6) and the accompanying discussion. The character of the excitation energy sharing mechanism has major consequences on the predicted spectrum of emitted neutrons and gammas.
Another important outcome was the clear indications that many more collective DoF appear to be relevant in fission than have ever been considered in either phenomenological or microscopically inspired models, such as GCM and ATDHF frameworks. At scission, both FFs are octupolly deformed. Moreover, the FFs attain their relaxed shapes only after the separation between them is ≈ 5...6 fm. No phenomenological or microscopically inspired approach on the market follow the FFs to such large separations. As we stressed in section 3 the spectrum of emitted neutrons is affected by the number of neutrons emitted before the full FFs acceleration. Another DoF, which appears to be relevant as well is the bending mode (see Figure 7 and Bulgac et al. [81]), the inclusion of which likely is going to influence the angular momentum distributions of the FFs [82].
We have pointed to several directions into which phenomenological models and theoretical models, such as GCM and ATDHF would have to be altered, in order to describe nuclear fission in a manner more consistent with theoretical expectations inferred from unrestricted quantum mechanical simulations (see section 3). We did not cover or mention all phenomenological models on the market and not all microscopically inspired theoretical frameworks, as this is not a review of such approaches, for which we recommend [129,139].

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.

FUNDING
The work of AB and SJ was supported by U.S. Department of Energy, Office of Science, Grant No. DE-FG02-97ER41014 and in part by NNSA cooperative agreement DE-NA0003841. The work of IS was supported by the US Department of Energy through the Los Alamos National Laboratory. Los Alamos National Laboratory is operated by Triad National Security, LLC, for the National Nuclear Security Administration of U.S. Department of Energy (Contract No. 89233218CNA000001).

ACKNOWLEDGMENTS
We thank many people with whom we had discussion over the years and for their input: I. Abdurrahman, G. F. Bertsch, P. Magierski, K. J. Roche, and N. Schunck. We thank G. F. Bertsch and J. Randrup who made a number of suggestions on our manuscript, which has been released as a preprint [140]. The