Predictions of New Neutron-Rich Isotopes at N = 126 in the Multinucleon Transfer Reaction 136Xe + 194Ir

Aiming to produce new neutron-rich nuclei at N = 126, the multi-nucleon transfer reactions 136Xe+194Ir, 208Pb are investigated by using the GRAZING model and the three-dimensional time-dependent Hartree-Fock (TDHF) approach. Deexcitation processes of the primary fragments are taken into account in both models. Comparison with the experimental data of 136Xe+208Pb at Ec.m. = 450 MeV indicates that the isotopic production cross sections around the entrance channel can be well reproduced by both models. The results of GRAZING indicate that 136Xe+194Ir is a promising candidate for producing new neutron-rich isotones with N = 126. The limitation of using TDHF to investigate multi-nucleon transfer reactions is also discussed.


INTRODUCTION
Neutron-rich nuclei are of great importance for understanding the astrophysical r-process [1]. For instance, those around the closed neutron shell N = 126 can provide information on the solar rabundance distribution [2]. The neutron-rich radioactive isotopes can also be used as projectiles for synthesizing super-heavy nuclei which is one of the most interesting challenges in nuclear physics. However, new heavy neutron-rich nuclei out of limits of the present nuclear chart are hardly to be produced via traditional ways, such as fusion reactions with stable beams at low energies or fragmentation of heavy projectile at relativistic energies.
In recent years, theoretical predictions indicate that multi-nucleon transfer (MNT) reactions would be a possible route to produce heavy neutron-rich nuclei far away from the stability line [3,4]. The experimental results of 136 Xe + 198 Pt at the incident energy of 8 MeV/nucleon [2] show that the production cross sections of neutron-rich nuclei with N = 126 are orders of magnitude larger than those obtained in fragmentation reaction of 208 Pb (1 AGeV) + Be [5]. But in some other experiments of MNT reactions with stable beams such as 64 Ni + 207 Pb [6] and 136 Xe + 208 Pb [7,8], no new neutron-rich nuclei is detected. Based on the N/Z equilibrium concept, theoretical predictions show that using neutron-rich radioactive beams can remarkably improve the production cross sections of neutron-rich isotopes along N = 126 [9,10]. However, due to the fact that the experimental intensities of radioactive beams are orders of magnitude lower than stable beams, the advantage of using radioactive beams may be canceled. Opportunities will arise in the near future on the second generation radioactive beam facilities like the European EURISOL [11]. At the present time, to find an optimum stable projectile-target combination are highly appealed.
To describe MNT reactions, many theoretical models are developed. For example, the multidimensional Langevin model shows great success on predicting the isotopic production cross sections, even for those nuclei far away from projectile and target [3,4,12]. Semi-classical GRAZING model [13,14] and the complex Wentzel-Kramers-Brillouin (CWKB) model [15] can well describe transfer process in peripheral collisions. In GRAZING model, the reactants move on classical trajectories in the combined field of Coulomb repulsion and nuclear surfacesurface attraction. Surface modes of the colliding nuclei are taken into account. Independent single nucleon transfer between the projectile-like and target-like nuclei during the collision is governed by the quantum coupled equations [13,14]. Neutron evaporation is considered for the excited primary fragments. The mass, charge, energy, and angular momentum distributions of the products produced in grazing collisions can be obtained. For details of GRAZING model and its applications (see, e.g., [16][17][18][19], and references therein). The dinuclear system model can reproduce the experimental data related with quasi-fission or deep-inelastic collisions [10,[20][21][22][23]. The improved quantum molecular dynamics (ImQMD) model [24][25][26] is capable of describing collisions from central to very peripheral regions on a microscopic basis, the widths of isotopic distributions can be reproduced in ImQMD since stochastic two-body collisions are taken into account [27][28][29]. The time-dependent Hartree-Fock (TDHF) theory shows success in describing few-nucleon transfer process [30][31][32][33][34]. Recently, the stochastic mean-field (SMF) approach beyond TDHF [35] has been proposed to investigate the damped MNT reaction 136 Xe + 208 Pb, the experimental broad mass distribution can be reproduced.
In this work, the GRAZING model and TDHF theory incorporating with GEMINI++ 1 [36] are adopted to investigate the production of neutron-rich nuclei in MNT reactions. The TDHF theory [37] is based on the independent particle picture and is a good approximation to the nuclear many-body problem. It is capable of describing low-energy heavy-ion reactions and provides insight on the average behavior of the dynamics. The state-of-art TDHF calculations are performed in a threedimensional (3D) framework without any symmetry constraints due to the advances in computational power. It has been applied for investigations on various subjects, for instance, collective vibration [38,39], fusion reaction [40,41], fission dynamics [42][43][44], dissipation mechanism [45][46][47][48]. Recently, 3D TDHF is applied to MNT reactions [30][31][32][33]49]. Fluctuation and dissipation can not be properly described in TDHF since twobody collisions and internucleon correlations are not included. Thus, widths of the mass or isotopic distributions in MNT reactions at incident energies above the Coulomb barrier are underestimated in TDHF. To consider these effects in TDHF is beyond the scope of this paper and further studied will be carried out in the future. This paper is outlined as follows. In section 2, the TDHF approach and particle-number projection (PNP) method as well as the numerical details are introduced. In section 3, numerical results of the isotopic production cross sections in both GRAZING and TDHF approach are shown. Dynamical properties of the reactions are also discussed for TDHF. Finally, a summary is drawn in section 4.

BRIEF INTRODUCTION TO TDHF APPROACH
In this section, we briefly introduce TDHF formalism, PNP method and the coupling to GEMINI++. Details of the TDHF theory can be found in, e.g., [50,51], and references therein. In TDHF approach, time evolution for the single-particle states are described by a set of coupled non-linear equations where ψ α, ν (r, t) is the single-particle state and N is the total number of states.ĥ is the self-consistent mean-field Hamiltonian of single-particle motion and it is always related to Skyrme energy density functional (EDF) which depends on local densities [52]. We underline that there are no adjustable parameters in the TDHF approach. The uncertainty in TDHF calculations may arise from the uncertainty of fundamental nuclear properties, such as the distribution of shape deformation of the reactants in their ground states. This can be obtained by preparing reactants with deformation constraints. The above TDHF equations can be derived from the variation of the action [53] is the correlated many-body wave function of the system and is given by a single Slater determinate.
In the present work, the 3D unrestricted TDHF code Sky3D [51] with Skyrme SLy5 parametrization [54] is adopted for both the static and dynamic calculations. The static HF calculations are performed on 32 × 32 × 32 Cartesian grids with 1.0 fm grid spacing in all three directions. In dynamical calculations, the meshes are extended to 70 × 32 × 70 with the same grid spacing in static HF. The projectile and target are initially placed at a separation distance of 24 fm and then boosted with the associated center-of-mass energy E c.m. and the impact parameter b. The time step △t is set to be 0.2 fm/c and six-order Taylor series expansion is employed. The TDHF simulations are stopped when the separation distance of the primary fragments after collision reaches 30 fm. Since single-particle wave functions are partially exchanged between the projectile-like and target-like nuclei in the collision process, the outgoing states are not the eigenstates of the particle-number operators (Ẑ for protons andN for neutrons) but superpositions of them. One can only get the expectation (mean) values of the charge and mass numbers for each fragment after collision. If one wants to get the distributions of proton or neutron numbers in one of the primary fragments, one should project the many-body states on good particle numbers by introducing the PNP operator [50,55] where q = n, p labels the nucleon species, and N n means N neutrons while N p means Z protons. The subscript V denotes the region of coordinate space encompassing one of the primary fragments,N q V = N q α=1 V (r) and V (r) = 1 if r ∈ subspace V and 0 elsewhere. The integral is performed with an M−point uniform discretion. For convergence, M is set to be 300.
The probability to find N q particles in subspace V is then obtained accordingly. The cross section of primary fragment with neutron number N and charge number Z at a certain incident energy is where P(N, Z; b) = P(N; b)P(Z; b) represents the probability to find N neutrons and Z protons at the impact parameter b; b max is a cutoff impact parameter which depends on the incident energy and should be large enough to guarantee that most of the transfer cross sections are included. But it is not necessary to set b max to be too large, because the transfer probability is extremely small and elastic collision dominates in very peripheral collisions. In this work, we set it to be 10 fm for 136 Xe+ 208 Pb at E c.m. = 450 MeV and 13 fm for 136 Xe+ 194 Ir at E c.m. = 720 MeV. b ranges from 0 to b max with the interval b = 1 fm. Deexcitation of the primary fragments are considered by using the statistical code GEMINI++ with default parameters [36,56]. All possible sequential binary decay modes, from emission of nucleons and light particles through asymmetric to symmetric fission as well as the γ -emission are included in GEMINI++. The code needs information of a primary fragment including charge and mass number as well as the excitation energies and the angular momentum as the inputs. Detailed calculations of these quantities can be found in Jiang and Wang [34]. The deexcitation of a certain projectile-like fragment (PLF) or target-like fragment (TLF) should be repeated M trial times due to the statistical nature of GEMINI++. Here M trial = 1, 000 is used. The number of events in which final fragment with (N final , Z final ) is counted and denoted as M(N final , Z final ). Then the final production cross section is given as Owing to the intrinsic stochastic nature of the Monte Carlo method employed in GEMINI++, Type A standard uncertainties for the isotopic production cross sections are calculated in the simplest case. The deexcitation process of a certain fragment is performed ten times repeatedly with M trial = 1, 000 for each time. Average values of σ (N final , Z final ) and the uncertainties can be obtained straightforwardly. We find the uncertainties are very small. For the sake of simplicity, the deexcitation process is performed only once with M trial = 1, 000 in this work. 136 Xe + 208 Pb is a candidate reaction for the production of neutron-rich nuclei at N = 126. The experiment at E c.m. = 450 MeV was performed at Argonne in 2015 [8]. In Figure 1 we plot the calculated isotopic production cross sections of the PLFs in this reaction. The results of GRAZING are shown as blue solid lines while those of TDHF+GEMINI are presented as black dashed lines. The experimental data are shown for comparison. The production cross sections of the TLFs are also calculated but have already been published in another paper [34]. One can find that for Z = 51 − 53 and Z = 55 − 58, the magnitude of the peak values can be reproduced by TDHF+GEMINI. Whilst those for Z = 52, 53, 55, 56 can be reproduced by GRAZING. For Z = 54, the peak values are overestimated in both models. This is because the results of (quasi)elastic channels are not excluded in our calculations. One can also find that better predictions are obtained for proton pickup channels than stripping channels in both models. As the number of transferred nucleons increases, discrepancies between model predictions and the experimental data get larger. These isotopes far away from the entrance channel may be produced in strongly damped collisions. Such processes can not be well estimated by the two models: two-body collisions are not considered in TDHF [57,58] while the GRAZING model only takes grazing collisions into account. However, no new neutron-rich nuclei were detected in 136 Xe + 208 Pb. The experiment results of 136 Xe+ 198 Pt at E c.m. = 645 MeV [2] indicate that this reaction is a better candidate to produce neutron-rich nuclei with N = 126. Calculations of this reaction are performed by using both GRAZING and TDHF+GEMINI. Unfortunately, no neutron-rich nuclei with N = 126 is obtained in TDHF+GEMINI. The production cross sections of those nuclei predicted by GRAZING are given in Figure 2. Detailed discussions about the predictions of TDHF on this reaction will be reported elsewhere. At the end of this section, we will show some results of 136 Xe + 194 Ir given by TDHF+GEMINI and discuss the limits of TDHF approach on investigating MNT reactions when the incident energy is much above the Coulomb barrier. In the following, we concentrate on the predictions of GRAZING on neutron-rich nuclei with N = 126.

RESULTS AND DISCUSSIONS
In order to find another optimum projectile-target combination for producing new exotic neutron-rich nuclei with stable reactants, we carry out a systematic study on 15 MNT reactions with projectile around 136 Xe and target around 198 Pt at various incident energies by using GRAZING. The results indicate that 136 Xe + 194 Ir is a surrogate reaction to produce neutron-rich nuclei around N = 126. The production cross sections of isotones with N = 126 in 136 Xe + 194 Ir are compared with those of 136 Xe + 208 Pb and 136 Xe + 198 Pt. The center-ofmass energy for all the three systems is set to be E c.m. = 645 MeV (this energy is the same as the experiment of 136 Xe + 198 Pt [2], and we will show later that it is also an optimum energy for 136 Xe + 194 Ir). The simulation results are plotted in Figure 2 as open symbols. The experimental data of 136 Xe + 198 Pt taken from [2] are shown as black solid squares for comparison. One can first see that more N = 126 isotones with charge number Z 78 can be produced in both 136 Xe + 194 Ir and 136 Xe + 198 Pt rather than 136 Xe + 208 Pb. Particularly, the system 136 Xe + 194 Ir has huge advantages for producing N = 126 isotones with Z 74. Those nuclei are out of the limits of the present nuclear landscape and are of great interest for nuclear and astro-nuclear physics. One can also find that the experimental data are underestimated by GRAZING. This can be understood since only peripheral collisions are treated in GRAZING and those nuclei with large number of nucleon transferred are produced in damped collisions at small impact parameters.
To find an optimum incident energy to produce more neutron-rich nuclei with N = 126 in 136 Xe + 194 Ir, the production cross sections of N = 126 isotones with Z = 72 − 77 at various incident energies from E c.m. = 1.07V B to 1.85V B (V B is the Bass barrier [59] and it is around 410 MeV for this reaction) are calculated by using GRAZING. The results are presented in Figure 3. It can be seen that very neutron-rich isotones with Z 74 can not be produced if E lab /A 5.5 MeV. The cross sections of all these isotones increase with the increasing incident energy when E lab /A < 6.5 MeV. However, when E lab /A is in the range of 7-9 MeV, plateau-like structures are observed for all the curves. This phenomenon indicates that the production cross sections of those nuclei are insensitive to the incident energy if 7 E lab /A 9 MeV.
In Figure 4 we show the isotopic production cross sections of the TLFs with Z = 72 − 77 in 136 Xe + 198 Pt and 136 Xe + 194 Ir at E c.m. = 645 MeV (E lab /A ≈ 8 MeV) as black dashed and red solid lines, respectively. It can be found that for nuclei on the proton-rich side, larger production cross sections are obtained in 136 Xe + 194 Ir. For isotopes with Z = 76 and 77, more neutron-rich ones are produced in 136 Xe + 198 Pt. Whilst for Z = 75 the two systems give comparable results for isotopes on the neutron-rich side. As Z decreases, the advantage to produce more neutron-rich nuclei arises in 136 Xe + 194 Ir. Note that the results of Z = 72 in 136 Xe + 198 Pt are not given by GRAZING. The above results of the two reactions are difficult to be understood through the N/Z equilibrium process because N/Z is around 1.52 for 136 Xe and 194 Ir while it is about 1.54 for 198 Pt and 208 Pb, respectively. It might be related with the Q gg value effect, where Q gg is the ground state-to-ground state Q value. Other quantum effect like shell effect might also play a role. Further studies should be carried out to check these discussions.
Finally, in Figure 5 we show the density contour plots of 136 Xe + 194 Ir at some special reaction stages in TDHF for two different initial configurations at E c.m. = 720 MeV (E lab /A ≈ 9) and b = 6 fm. Isodensities at half the saturation density (ρ 0 /2 = 0.08 fm −3 ) is plotted as black solid lines. We should mention  136 Xe + 198 Pt are taken from Watanabe et al. [2] and shown as black solid squares. The lines are drawn to guide the eye.  136 Xe + 194 Ir as a function of the incident energy. The lines are drawn to guide the eye. that in TDHF the projectile and target are both deformed in their ground states with β 2 = 0.064 and 0.154, respectively. Since the deformation of 136 Xe is very small, we only take into account the orientation effect of the deformed 194 Ir at the beginning of the dynamical calculations. The two initial configurations are named "tip collision" (the symmetry axis of 194 Ir is set parallel to the bombarding direction: z-axis) and "side collision" (the symmetry  136 Xe + 198 Pt and 136 Xe + 194 Ir at E c.m. = 645 MeV, which are plotted as black dotted and red solid lines, respectively. axis of 194 Ir is set parallel to x-axis). It can be found that in both configurations, the composite system is strongly elongated before rupture of the neck. The primary fragments produced at the end of the dynamical calculations also have large deformations. Similar results are obtained for TDHF+GEMINI in 136 Xe + 198 Pt at E c.m. = 645 MeV. Such large deformation in the exit channel makes the primary fragments have very large excitation energies. This leads to strong evaporation of nucleons in the deexcitation process. So no neutron-rich isotope with N = 126 is observed in TDHF+GEMINI for 136 Xe + 194 Ir and 136 Xe + 198 Pt reactions when the incident energy is much larger than the Coulomb barrier. It is well-known that mean-field model such as TDHF is suitable for low-energy reactions, however, the lack of two-body collisions limits the predictive power of TDHF on estimating the yields of MNT reactions when the incident energy is much above the barrier.

SUMMARY
Employing the semi-classic GRAZING model and the microscopic TDHF approach, we have investigated the MNT reactions of 136 Xe + 194 Ir, 208 Pb. Neutron evaporation is considered in GRAZING while GEMINI++ is coupled with TDHF to deal with the deexcitation process. The calculated production cross sections of the PLFs for 136 Xe + 208 Pb at E c.m. = 450 MeV are compared with the experimental data. The model predictions can well describe the yields of nuclei near the projectile. The predictions of GRAZING on 136 Xe + 194 Ir show that this reaction has the advantage to produce more neutron-rich nuclei with N = 126 compared with 136 Xe + 198 Pt. The latter one has been carried out at GANIL [2] and the results are inspiring. No new neutron-rich isotope is observed for these two reactions at energies much above the barrier in TDHF+GMINI which might be interpreted as the lack of two-body collisions in the mean-field theory. Further investigations on 136 Xe + 194 Ir by using other theoretical models with two-body collisions included, such as ImQMD, are in progress.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/supplementary material.

AUTHOR CONTRIBUTIONS
XJ performed numerical calculations. XJ and NW wrote the manuscript.