Microscopic predictions of ﬁssion yields based on the time dependent GCM formalism

. Accurate knowledge of ﬁssion fragment yields is an essential ingredient of numerous applications ranging from the formation of elements in the r-process to fuel cycle optimization in nuclear energy. The need for a predictive theory applicable where no data is available, together with the variety of potential applications, is an incentive to develop a fully microscopic approach to ﬁssion dynamics. One of the most promising theoretical frameworks is the time-dependent generator coordinate method (TDGCM) applied under the Gaussian overlap approximation (GOA). Previous studies reported promising results by numerically solving the TDGCM + GOA equation with a ﬁnite di ﬀ erence technique. However, the computational cost of this method makes it di ﬃ cult to properly control numerical errors. In addition, it prevents one from performing calculations with more than two collective variables. To overcome these limitations, we developed the new code FELIX-1.0 that solves the TDGCM + GOA equation based on the Galerkin ﬁnite element method. In this article, we brieﬂy illustrate the capabilities of the solver FELIX-1.0, in particular its validation for n + 239 Pu low energy induced ﬁssion. This work is the result of a collaboration between CEA,DAM,DIF and LLNL on nuclear ﬁssion theory.


Introduction
Data on fission fragment yields play an essential role both in the conception of new generations of reactors as well as in fuel cycle technology. Other growing activities such as neutron interrogation would also benefit from a better knowledge of differential fission distributions involving the mass and charge but also the spin of the fragments. Over the last decades, part of the theoretical effort on fission has been focused on predicting induced fission yields based on the time-dependent generator coordinate method (TDGCM) associated with the Gaussian overlap approximation (GOA) [1,2]. This approach provides a fully quantum mechanical description of the time evolution of the fissioning system. Recent applications provided the first fully-microscopic estimate of the mass and kinetic energy distributions of fission fragments for neutron-induced fission in the actinides region [3][4][5]. Note that these studies used the same solver based on the discretization of the TDGCM+GOA equations using finite differences in a regularly meshed hyper-cube. Because of the computational resources required by this numerical scheme, fission dynamics was only studied in 2 it is well-known from both semi-phenomenological and fully microscopic approaches that at least four or five collective variables play a role in the dynamics of fission [6][7][8][9]. Although possible in principle, extending the current scheme to N > 2 collective spaces would rapidly become prohibitive computationally. These reasons were the primary motivation to upgrade to a more flexible and more scalable discretization scheme. In this work, we give a quick overview of the new solver FELIX-1.0 [10] for the TDGCM+GOA equation. We briefly outline the numerical methods implemented before discussing some preliminary benchmark results on fission fragment yields obtained within this approach for the neutron-induced fission on 239 Pu.

TDGCM+GOA approach to compute fission yields
The time-dependent generator coordinate method (TDGCM) provides a fully quantum mechanical description of the time evolution of the fissioning system. It starts from the general observation that the time evolution of the many-body fissioning system is governed by the time-dependent Schrödinger equation. In our applications, the nuclear HamiltonianĤ involved contains an effective two-body potential such as, e.g., the Skyrme or the Gogny interaction. In the TDGCM method, one chooses the following anszatz for the many-body state |Ψ(t) associated with the fissioning nucleus [1,2], In fission, the states |Ψ(q) are solutions to the static Hartree-Fock-Bogoliubov (HFB) equations under a set of constraints q while the functions f (q, t) are to be determined. The constraints, also called collective variables, are in this work expectation values of the quadrupole and octupole moments. The quadrupole moment q 20 characterizes the elongation of the nucleus whereas the octupole moment q 30 is related to the asymmetry between the pre-fragments. Inserting Eq.(1) in a variational form of the time-dependent Schrödinger equation yields a timedependent integro-differential equation for the weight function f (q, t) known as the Hill-Wheeler equation. Numerically solving this equation has not yet been attempted in the case of fission as it would require a tremendous amount of computation. Instead, a widespread approach consists in assuming that the norm kernels Ψ(q)|Ψ(q ) can be approximated by a Gaussian form factor [11]. This Gaussian overlap approximation (GOA) yields a local, time-dependent Schrödinger-like equation for the dynamics of the fissioning system in the space of the collective variables q, where • The unknown complex function g(q, t) is related to the weight function f (q, t) appearing in Eq. (1). It contains all the information about the dynamics of the system. In the limit of sharply peaked norm kernels, the quantity |g(q, t)| 2 can be interpreted as the probability density for the system to have the collective characteristics q at time t; • The potential field V(q) and the inertia tensor B kl (q) are fully determined by the knowledge of the effective HamiltonianĤ and the generator states |Ψ(q) .
The equation (2) is referred to as the TDGCM+GOA equation.
Our microscopic approach to compute the fission yields is a three-step process. First, we compute the generator states |Ψ(q) for a range of values of the collective variables. This is the most timeconsuming part, since any realistic two-dimensional fission calculation involves of the order of 50,000 points of the collective space. From these generator states, we can compute both the potential energy surface and the inertia tensor. A second step consists in numerically solving the TDGCM+GOA equation (2) starting from an initial function g(q, t = 0). This initial condition is arbitrarily built as a wave packet localized at low deformation, in the first well of the potential energy surface. The system begins its evolution with an average initial energy of 1 MeV above the inner barrier to simulate an induced fission. In the third and final step, we extract fission product yields from the time-dependent solution of the dynamic equation by defining a frontier in the collective space, which separates scissioned configurations from configurations where the system is still a whole nucleus. Computing the probability for the system to cross through any point of this line provides a first estimate of the fission yields. This raw data is then convoluted with a Gaussian to take into account effects such as quantum fluctuations of particle number. Figure 1 gives a schematic picture of the whole process.

FELIX-1.0, a new solver for the TDGCM+GOA equation
As already mentioned, previous studies of induced fission within the TDGCM+GOA framework relied on the same code tdgcm2d developed and used at CEA and LLNL over the last decades. This code uses a finite difference scheme to discretize the function g(q, t) in the collective space. Because of the high computational cost of this numerical method, fission dynamics could only be studied with two collective variables. However, it is well-known from both semi-phenomenological and fully microscopic approaches that at least four or five collective variables play a role in the dynamics of fission [6][7][8][9]. Although possible in principle, extending the current scheme to N > 2 collective spaces would rapidly become prohibitive computationally: Each additional dimension multiplies both the pre-calculation time of the potential energy surface and the evolution time of the TDGCM equations by approximately a factor 100. In order to fully exploit this microscopic technique, we have thus developed a new solver for the TDGCM+GOA equation called FELIX-1.0, with the following capabilities: • handling of more than two collective variables for a better description of the fissioning system, • full control on the numerical precision, • discretization of the collective space on irregular meshes to focus the computational effort on physically-relevant regions.
FELIX-1.0 relies on the Galerkin finite element method for the q discretization of the TDGCM equation. This well-established method allows the use of irregularly-spaced meshes in the collective space, which results in turn in calculations scaling much more efficiently with the number of collective variables. In addition, standard p-refinement techniques, i.e., the use of higher-degree polynomial bases in each finite element, give a better control on the numerical precision of the calculations. A Cranck-Nicholson unitary scheme is used for the time integration of the wave packet. This implicit method requires solving a sparse linear system at each time step of the evolution. The numerous matrix inversion operations involved are performed with a quasi-minimal residual (QMR) algorithm. Finally, we note that there are virtually no restrictions in the number of collective variables used. The code FELIX-1.0 has been tested up to N = 3 but should also work for higher dimensions. All this work is summarized in [10]. The first version of the FELIX-1.0 solver was released under an open source GPL-2 license and will be available in the Computer Physics Communications program library.

Preliminary results on n+ 239 Pu
We performed a set of preliminary calculations for the pre-emission mass yields of the low energy induced fission on 239 Pu. Our first objective was to check the convergence of the calculations. Four parameters control the numerical precision of the results: the mesh spacing parameter h, the time step δt, the tolerance for matrix inversions, and the convergence criteria for the initial state calculation.
Thanks to the efficiency of FELIX-1.0, we can perform multiple series of calculations to check the convergence associated with each one of these parameters. Figure 2 shows the calculated mass yields of the heavy fragment in n th + 239 Pu as a function of the mesh parameter h for the SkM* parametrization of the Skyrme force. Within a reasonable amount of computational time, the precision obtained reaches 2% in the peak region (A H ∈[127;150]) and is maintained below 6% in the wings of the distribution (A H < 127 and A H > 150). With this robust numerical scheme, we have solved the TDGCM+GOA equation for several potential energy surfaces based on both Skyrme and Gogny interactions. In a few cases, the computed mass yields show an encouraging 30% agreement with experiments in the mass range 127-150. On the other hand, the dispersion of our results with changes of potential energy surface may reach 100% at the peak. Although the interpretation of such discrepancy is still in progress, most of it seems to be related to the well-established presence of discontinuities [8] in the two dimensional collective space used. One possible way to overcome this drawback would be to increase the number of collective coordinates. Systematic, large-scale investigations of fission product yields in n th + 239 Pu for a variety of microscopic inputs are currently on-going and will be presented elsewhere.

Conclusion
In this work we use a fully quantum mechanical approach to compute fission fragment yields for neutron-induced fission with thermal neutrons. We presented the new code FELIX-1.0 to solve the dynamics of fission within the TDGCM+GOA framework. FELIX-1.0 gives full control on the numerical precision of fission product yields in neutron-induced fission, and it scalability also enables series of dynamical calculations on several potential energy surfaces. Preliminary results suggest an important sensitivity of our two-dimensional approach to the input potential energy surface.