Constraining the Woods-Saxon potential in fusion reactions based on the
  neural network
The accurate determination of the nuclear interaction potential is essential
for predicting the fusion cross sections and understanding the reaction
mechanism, which plays an important role in the synthesis of superheavy
elements. In this work, the neural network, which combines with the
calculations of the fusion cross sections via the Hill-Wheeler formula, is
developed to optimize the parameters of the Woods-Saxon potential by comparing
the experimental values. The correlations between the parameters of Woods-Saxon
potential and the reaction partners, which can be quantitatively fitted to a
sigmoid-like function with the mass numbers, have been displayed manifestly for
the first time. This study could promote the accurate estimation of
nucleus-nucleus interaction potential in low energy heavy-ion collisions.
Zepeng Gao, Siyu Liu, Peiwei Wen, Zehong Liao, Yu Yang, Jun Su, Yongjia Wang, Long Zhu
Constraining the Woods-Saxon potential in fusion reactions based on a physics-informed neural network

###### Abstract

The accurate determination of the nuclear interaction potential is essential for predicting the fusion cross sections and understanding the reaction mechanism, which plays an important role in the synthesis of superheavy elements. In this letter, a Physics-Informed Neural Network (PINN) method, which combines the neural networks with the calculations of the fusion cross sections via the Hill-Wheeler formula, is developed to optimize the parameters of the Woods-Saxon potential by comparing the experimental values. The correlations between the parameters of Woods-Saxon potential and the reaction partners, which can be quantitatively fitted to a sigmoid-like function with the mass numbers, have been displayed manifestly for the first time. This study could promote the accurate estimation of nucleus-nucleus interaction potential in low energy heavy-ion collisions.

+
Footnote †: journal:

## 1 Introduction

In recent decades, the study of heavy-ion fusion reactions has garnered growing interest in the field of nuclear physics due to their importance for extending the period table of elements, as well as the understanding of the interplay between nuclear structure and the reaction dynamics [1; 2; 3; 4; 5]. However, the complexity of these reactions presents significant challenges for both experimental and theoretical investigations. One of the key challenges is accurately determining the nuclear interaction potential, which is crucial for investigating the reaction mechanism, but remains an arduous task due to several factors, including the quantum many-body problem, the form of nuclear force, among other complex influences [6; 7; 8; 9].

As one of most successful phenomenological forms of nuclear potential, the Woods-Saxon potential has been widely used to describe the nuclear interaction in heavy-ion fusion reactions, especially for the synthesis of the superheavy nuclei. However, the accuracy is limited by the uncertainties in the parameters. Much effort has been made to constrain the parameters in recent years [10; 11; 12; 13]. Nevertheless, due to a complex function with multiple and correlated parameters as well as existing computational limitations, the majority of studies primarily apply constraints on the parameters for specific systems, such as focusing solely on \({}^{12}\)C or \({}^{16}\)O induced reactions [14; 15]. Therefore, it is imperative to accurately and efficiently constrain and optimize the Woods-Saxon potential parameters from a global perspective.

Machine learning, being adept at uncovering underlying patterns from vast amounts of data, as well as accurately fitting and predicting data, have garnered significant interest across various fields. In the realm of nuclear physics, they hold immense potential for addressing challenges in terms of nuclear theories and experiments [16; 17; 18; 19]. For example, machine learning methods have demonstrated successful applications in predicting nuclear masses [20; 21], charge radii [22; 23], half-lives [24; 25], nuclear energy density functional [26], fission product yields [27; 28], as well as facilitating studies on nuclear reactions [29; 30; 31; 32; 33; 34; 35; 36; 37; 38; 39]. In recent years, although great progress has been made in applying machine learning methods to address issues in nuclear physics, it is worth to note that most of those studies are data-driven approaches. Data in a real process are governed by physical laws, thus, it is imperative to integrate the fundamental principles of physics into machine learning methods, enabling them to be guided by the laws of physics rather than relying solely on data-driven approaches.

In fact, a physics-informed neural network (PINN) model was proposed by M. Raissi _et al._, which has been utilized for solving partial differential equations (PDEs) [40]. PINN is an innovative neural network model that combines physical principles with traditional data-driven approaches. By incorporating physical information as computational constraints, PINNs enhance the accuracy and robustness of predictions made by purely data-driven neural networks. It has also effectively tackled several crucial problems in nuclear physics. For nuclear masses, M. R. Mumpower _et al._ integrated the constraint of the Garvey-Kelson relations existing in the nuclear binding energy into the machine learning model, resulting in remarkably accurate predictions without reference to any underlying theoretical model [41]. Furthermore, F. P. Li _et al_ used PINN to successfully reproduce the equation of state of a strongly-interacting quantum chromodynamics [42]. In Ref. [43], theoretical models of bottomonium productions in high-energy nuclear collisions was also characterized by using PINN. In general, PINN
are capable of achieving objectives that are beyond the reach of conventional neural networks.

In this letter, we will explore the capability of the PINN to constrain the parameters of the Woods-Saxon potential model by comparing the calculations of fusion cross sections using the Hill-Wheeler formula with the corresponding experimental values. Subsequently, we will explicitly demonstrate the correlations between the parameters and colliding partners. The outcomes of this investigation will be highly significant in enhancing our comprehension of heavy ion fusion reactions and will offer a fresh perspective for exploring nuclear physics co-nundrums.

## 2 Theoretical method

### Fusion cross section and Woods-Saxon potential

At a given center-of-mass energy, the fusion cross section \(\sigma_{\text{fusion}}\left(E_{\text{c.m.}}\right)\) can be expressed as the sum of the cross sections at each partial wave \(J\):

\[\sigma_{\text{fusion}}\left(E_{\text{c.m.}}\right)=\frac{\pi\hbar^{2}}{2\mu E_{ \text{c.m.}}}\sum\limits_{J=0}^{J_{\text{thresh}}}(2J+1)T\left(E_{\text{c.m.}} \right.,J)\,, \tag{1}\]

where \(\mu\) is the reduced mass, \(T\) is the penetration probability and \(J\) is the incident angular momentum. Note that we have specifically selected relatively light combinations of target and projectile, where the occurrence of quasi-fission is relatively minimal. The Hill-Wheeler formula [44] is a well-known analytical expression for the penetration probability:

\[\begin{split}& T_{\text{HW}}\left(E_{\text{c.m.}},J\right)\\ &=\left\{1+\exp\left[\frac{2\pi}{\hbar\omega(J)}\left(\frac{ \hbar^{2}J(J+1)}{2\mu R_{\text{B}}^{2}(J)}+B-E_{\text{c.m.}}\right)\right] \right\}^{-1},\end{split} \tag{2}\]

where \(R_{\text{B}}\) and \(\hbar\omega(J)\) correspond to the position and curvature of the barrier under the \(J^{\text{th}}\) partial wave, respectively. The barrier curvature \(\hbar\omega(J)\) can be calculated using the following formula:

\[\hbar\omega(J)=\left.\sqrt{-\frac{\hbar^{2}}{\mu}\frac{\partial^{2}}{\partial R ^{2}}V(R,J)}\right|_{R=R_{\text{B}}(J)}. \tag{3}\]

The interaction potential \(V(R,J)\) consists of a long-range Coulomb potential, a short-range nuclear potential, and a centrifugal potential, which can be expressed as follows [45]:

\[V(R,J)=V_{\text{C}}(R)+V_{\text{N}}(R)+V_{\text{R}}(R,J). \tag{4}\]

The fusion cross sections at the sub-barrier, especially the deep sub-barrier region strongly depends on the coupling of the nuclei structures. We would like to emphasize that in order to minimize the uncertainties and focus on the parameters of Woods-Saxon potential, the fusion reactions occurring above the barrier are only investigated. The Coulomb potential and centrifugal potential can be expressed by the forms:

\[\begin{split}& V_{\text{C}}\left(R\right)=\frac{Z_{1}Z_{2}e^{2}}{R}, \\ & V_{\text{R}}\left(R,J\right)=\frac{\hbar^{2}J(J+1)}{2\mu R^{2}}. \end{split} \tag{5}\]

The Woods-Saxon nuclear potential can be expressed as [46]:

\[V_{\text{N}}\left(R\right)=\frac{-V_{0}}{1+\exp\left[\left(R-R_{\text{P}}-R_{ \text{T}}\right)/a\right]}, \tag{6}\]

with

\[R_{i}=r_{0}A_{i}^{1/3},i=\text{P},\text{T}, \tag{7}\]

where the depth \(V_{0}\), the radius parameter \(r_{0}\), and the diffuseness parameter \(a\) are the major parameters of the Woods-Saxon potential. By changing the parameters of the Woods-Saxon potential, we can modify the position, height as well as curvature of the potential barrier, which ultimately affects the prediction of the fusion cross sections. Therefore, we can evaluate the optimization results of the parameters by comparing the predicted cross sections via PINN with experimental values.

### The physics-informed neural network

PINN primarily involves grafting a physical constraint onto the output of a neural network. In this letter, we perform calculations of fusion cross sections and train the neural network by comparing it with experimental data, as illustrated in Fig.1. Initially, \(Z\) (proton number) and \(A\) (mass number) of the projectile and target nuclei as the input features are fed into a multi-output neural network, which then generates three parameters corresponding to the Woods-Saxon potential. Moreover, the nuclear interaction potential at various angular momenta and the associated input parameters required for the Hill-Wheeler formula are obtained. By integrating the capture probability over all angular momenta, the corresponding fusion cross sections can be derived. By comparing the fusion cross-section calculated using Hill-Wheeler formula with parameters \(V_{0}\),\(r_{0}\) and \(a\) obtained from PINN to the experimental values, we can establish a loss function, as expressed below:

\[\mathcal{L}(i)=\left\{\sum_{e}\left[\lg\left(\sigma_{\text{NN}}^{(i)}(e) \right)-\lg\left(\sigma_{\text{exp}}^{(i)}(e)\right)\right]\right\}^{2}, \tag{8}\]

where \(e\) represents the \(E_{\text{c.m.}}\) associated with the experimental value in the \(i\)th projectile-target combination. It is important to note that for each set of parameters of Woods-Saxon potential, a fusion excitation function is calculated and compared with all experimental values associated with the projectile-target combination. Consequently, the resulting loss is not solely determined by each individual energy, but rather by considering the combined outcomes of all the \(E_{\text{c.m.}}\) with the specific system. Obtaining the loss allows for the training of the multi-output neural network, which will be discussed in the subsequent context.

Foremost, it is crucial to clarify that \(\mathcal{L}(i)\) is not attributed to the parameters of Woods-Saxon potential, but rather to the fusion cross sections. Moreover, there is no explicit functional association between the fusion cross sections and the parameters of Woods-Saxon potential, so that it is not feasible to directly obtain the loss of three parameters using Backpropagation algorithm. However, the qualitative association can be analyzed and the existence of a positive correlation between them has been established in advance. Therefore,
loss of the parameters of Woods-Saxon potential to the loss of fusion cross sections and assign suitable weight coefficients in a straightforward manner, as expressed below:

\[\mathcal{L}(i)^{(j)}=\alpha^{(j)}\mathcal{L}(i),j=V_{0},r_{0},a. \tag{9}\]

The weight coefficients, by taking into account the empirical ranges, were assigned as 1, 0.001, and 0.01 for \(V_{0}\), \(r_{0}\) and \(a\), respectively. It is worth to emphasize that the weight coefficients weakly influence the behavior of correlations found in this work. After obtaining the loss of the parameters, the Backpropagation algorithm can be employed for training the multi-output neural network.

### Parameter setting of the multi-output neural network

The multi-output neural network consists of three neural networks, wherein each network comprises three hidden layers with 32, 64, and 128 neurons, respectively. "Tanh" and "Adam" were selected as the activation function and optimizer. The learning rate \(\alpha\) started at \(10^{-4}\) and gradually decreases to \(10^{-6}\) during the training process. Furthermore, a total of 343 reactions from the online dataset [47] were chosen for this study, where \(Z_{1}Z_{2}\leq 1600\) and there were at least three experimental data points above the Coulomb barrier. To select the data, the Coulomb barrier can be estimated by the latest fitting formula [48] as:

\[V_{\text{B}}=\frac{Z_{\text{P}}Z_{\text{T}}e^{2}}{0.9782\left(A_{\text{P}}^{1/ 3}+A_{\text{T}}^{1/3}\right)+4.2833}. \tag{10}\]

30 projectile-target combinations of the dataset were randomly selected as the testing set and validation set, respectively. Thus, the remaining 283 combinations were utilized as the training set and batch size was set to 70 during training.

### Training the neural network

The loss curves of training set, validation set and testing set, by training on 283 combinations through multiple optimizations and structural adjustments, are displayed in Fig. 2. The loss function can be defined as the average value of the individual loss function \(\mathcal{L}(i)\) computed for each projectile-target combination:

\[\text{Loss}=\frac{1}{N}\sum_{i}^{N}\mathcal{L}(i) \tag{11}\]

The loss for the validation and testing sets initially decreases and then plateaus or slightly increases, whereas the loss for the training set continuously decreases with the increase of epoch. Hence, a truncation, to prevent overfitting and to obtain the extrapolation ability in the neural network, was employed at epoch = 500. Furthermore, one can find that the loss for the training set is higher than that for the validation set, it does not necessarily imply underfitting, since the loss on the validation set has already plateaued. This discrepancy is mainly due to

Figure 1: The architecture of PINN used in this work. It combines neural networks and process of the calculation of fusion cross-sections via Hill-Wheeler formula.

Figure 2: The loss curve for training (solid line), validation (dash line) and testing set (dotted line) as functions of the epoch (training time).
the impact of data sampling on the training and validation sets. However, this effect is not significant for subsequent correlation analysis.

The extrapolation ability of the PINN can be demonstrated by its predictions on fusion cross section in the testing set, which has not been through the training process. The predicted cross section from PINN have been compared with the experimental data, as shown in Fig. 3. One can see that the extrapolation capability of the model is reliable, as the predicted cross section values are in good agreement with the experimental data for most of the reaction systems. However, in light projectile-target combination, there are some discrepancies between theoretical predictions and experimental data. This is primarily due to the strong structural effects that exist in such systems, which can lead to the possibility of cluster formation.

## 3 Results and discussion

Our aim is to explore the underlying correlations of the parameters of Woods-Saxon potential with the reaction systems. We show the parameters of Woods-Saxon potential constrained by the PINN in Fig. 4. The scatter plots for the parameters of mass number, proton number as well as Coulomb parameter (\(\chi=Z_{1}Z_{2}/(A_{1}^{1/3}+A_{2}^{1/3})\)) are shown. The strong correlations indicate that the Woods-Saxon potential parameters could be constrained and characterized by these quantities. It can be seen that all of these parameters fall within a reasonable range. The radius parameter \(r_{0}\) displays a rapid increase from 1.14 fm to approximately 1.18 fm, followed by a slower convergence to around 1.19 fm as the mass number \(A\) and proton number \(Z\) increase. Additionally, the depth parameter \(V_{0}\) exhibits a discernable trend of increasing with increasing \(A\), \(Z\), and \(\chi\), with values falling within the range of roughly 77-83 MeV. The surface diffuseness parameter \(a\) shows an even more pronounced correlation, with a rapid increase followed by a slower rise as \(A\), \(Z\), and \(\chi\) increase, within a range of approximately 0.4-0.8 fm. Note that the distributions of \(V_{0}\) and \(a\) exhibit a certain degree of diffuseness, primarily because these two parameters are shared by the entire combination and are not solely determined by either the projectile or the target nucleus. The Coulomb parameter is a measure of the direct correlation between the projectile-target combination, and there exists a strong linear correlation with the Coulomb barrier [48; 49], which directly affects the fusion cross sections and makes it highly relevant to the constrained \(V_{0}\) and \(a\).

The above trends and ranges between the three parameters of the Woods-Saxon potential and \(A\), \(Z\), and \(\chi\) are considered reasonably in physics [10; 15; 50]. The Woods-Saxon potential describes the nuclear force potential in the nucleus-nucleus interaction, and the nuclear force is a short-range strong interaction attraction. Among them, \(V_{0}\) controls the strength of the nuclear force, which increases with \(A\) in the region of light partner. However, the correlation between \(V_{0}\) and heavy collision partner is weak. This is because for large \(A\), the value of potential is approximately flat in the center and the variation is weak. Currently, the available experimental data are mostly in the reactions with stable projectiles and targets. Due to the curvature of the \(\beta\) stable line, the spread of neutron skin is shown with increasing isospin. Consequently, the trend that the parameters

Figure 3: The experimental and predicted fusion excitation functions for testing set. The solid lines denote the predicted funsion cross sections from PINN, while the dots denote experimental data. The parameters constrained by PINN are also displayed.
\(a\) and \(r_{0}\) increase with increasing mass and charge number for both light and heavy partners can be seen.

Based on the qualitative constraints mentioned above, we expected to obtain a more specific, intuitive, and analytical quantitative expression. As a result, the strongly correlated parameters \(a\) and \(r_{0}\) were fitted with the mass numbers of projectile-target combinations using a sigmoid-like function, as shown in Fig. 5. The good performance of fitting using a sigmoid-like function also can be easily observed in upper panel of Fig. 5 where the predicted \(a\) is plotted versus the fitting results. The fitting results are better in regions where \(a\) is larger, which is consistent with the broadening of \(a\) observed in Fig. 4. The fitting accuracy for \(r_{0}\) is higher in lower panel of Fig. 5, mainly due to \(r_{0}\) depends only on the individual projectile or target nucleus instead of the combination. The figure displays the fitted parameters and coefficient of determination (\(R^{2}\)), which provide an indication of the goodness-of-fit. The coefficient of determination for both parameters is greater than 0.96, which suggests that the fitting is highly effective. The final quantitative constraints on \(a\) and \(r_{0}\) can be expressed as follows:

\[a=\frac{\mathrm{c_{1}}}{\mathrm{c_{2}+e^{\mathrm{c_{3}A_{1}}}+e^{\mathrm{c_{4 }A_{2}}}}}, \tag{12}\]

where, \(\mathrm{c_{1}}=1.0259\) fm, \(\mathrm{c_{2}}=1.3603\), \(\mathrm{c_{3}}=-0.0877\), and \(\mathrm{c_{4}}=-0.0352\).

\[r_{0}=\frac{\mathrm{d_{1}}}{\mathrm{d_{2}+e^{\mathrm{d_{3}A}}}}, \tag{13}\]

where, \(\mathrm{d_{1}}=19.7344\) fm, \(\mathrm{d_{2}}=16.5655\), and \(\mathrm{d_{3}}=-0.0343\). \(A_{1}\), \(A_{2}\), and \(A\) denote the mass numbers of light partner, heavy partner and total the whole system, respectively.

## 4 Summary

We have successfully incorporated the calculations of the fusion cross sections via the Hill-Wheeler formula into a neural network and constructed a physics-informed neural network framework to optimize the parameters of the Woods-Saxon potential in a more accurate and efficient way. Upon comparison with experimental fusion cross sections, the PINN demonstrated impressive predictive performance on the testing set, and effectively constrained the three parameters of the Woods-Saxon potential. Furthermore, the constrained values of \(a\) and \(r_{0}\) exhibit a strong correlation with the projectile-target combinations, which are accurately fitted by a sigmoid-like function for the first time. This function can be conveniently used to model heavy-ion fusion reactions. Moreover, PINN methods are expected to become the mainstream direction in nuclear physics. Compared with traditional data-driven methods, PINN is more efficient, more accurate, and can introduce prior knowledge into physical models to improve the interpretability of the model.

**Declaration of competing interest**

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

**Acknowledgments**

This work was supported by the National Natural Science Foundation of China under Grant No. 12075327; Fundamental Research Funds for the Central Universities, Sun Yat-sen University under Grant No. 23lgbj003; The Open Project of Guangxi Key Laboratory of Nuclear Physics and Nuclear

Figure 4: Scatter plots display the correlation between the horizontal and vertical labeled quantities. Each scatter point denotes one projectile-target combination.

Figure 5: Upper panel: the predicted \(a\) from PINN versus the fitted \(a\). Each scatter point denotes one projectile-target combination. Lower panel: the predicted \(r_{0}\) from PINN (dots) and the fitted \(r_{0}\) (solid line) as the functions of mass number. Each scatter point denotes a projectile or target nucleus. The fitting formula and corresponding coefficients, as well as the coefficient of determination, are also displayed.
Technology under Grant No. NLK2022-01; Guangdong Major Project of Basic and Applied Basic Research under Grant No. 2021B0301030006; Central Government Guidance Funds for Local Scientific and Technological Development, China (No. Guike ZY22096024).
