1. Introduction
Silicon-based materials, notably silicon carbide (SiC), are extensively utilized in aerospace applications owing to their exceptional radiation hardness and high-temperature resilience [1–4]. Nevertheless, the operational reliability of electronic and structural components remains significantly challenged by radiation-induced degradation and single-event effects [5, 6] within harsh radiative environments. These phenomena are primarily triggered by nuclear interactions between the constituent nuclei and cosmic-ray particles. In these environments, protons dominate both galactic (about 87%) [7] and solar (over 89%) [8] cosmic-ray fluxes. Given that $ ^{28}$Si constitutes approximately 92.2% [9] of natural silicon, it serves as the predominant target for proton-induced reactions in aerospace instrumentation. Such interactions generate primary knock-on atoms and secondary particles, which can subsequently trigger single-event upsets. Consequently, high-precision nuclear reaction data are indispensable for the predictive modeling of radiation-induced damage.
Energy spectra and double-differential cross sections (DDCSs) are essential for understanding reaction mechanisms, as they encode detailed energy-angular correlationsv of outgoing particles. However, the EXFOR database [10] provides only limited DDCSs of light charged particles (LCPs) for the p+$ ^{28}$Si reaction. This scarcity of data is insufficient to impose effective constraints on the model parameters to enhance the predictive power of the theoretical models. Consequently, measurements using natural Si targets offer a viable alternative, owing to the high isotopic abundance of $ ^{28}$Si.
The nuclear reaction codes, such as GNASH [11] and TALYS [12], which consider direct nuclear reactions as well as pre-equilibrium and compound nucleus emission, can reproduce the general features of the experimental data for the p+$ ^{28}$Si reaction. More detailed comparisons between the calculated results with the experimental data indicate that there is still room for improvement in these codes [13, 14]. In addition, a consistent analysis of the collective band structure and nucleon scattering data for n & p+$ ^{28}$Si reaction has been performed using the soft rotator model and coupled-channels method. However, the DDCSs for outgoing neutrons and light charged particles are not provided [15]. In this work, on the basis of spherical optical model, direct reactions in discrete levels and intranuclear cascade nucleons emission, pre-equilibrium statistical theory, evaporation model and Hauser–Feshbach theory with width fluctuation correction, which are all encapsulated within the MEND code [16], we self-consistently calculate the reaction cross sections, angular distributions, energy spectra and DDCSs of the neutron, proton, and the light charged particles (including deuteron, triton, $ ^3$He, and alpha). The main objective is to use the improved Iwamoto–Harada model to provide a comprehensive description of the energy spectra and DDCSs of the light charged particle (d, t, $ ^3$He and $\alpha$) emissions in the intermediate-energy range.
This paper is organized as follows. In section 2 , the theoretical framework is briefly described. In section 3 , the analysis and comparison of calculated energy spectra and DDCSs with experimental data are performed. A summary is succinctly given in section 4 .
2. Theoretical model
The optical model is one of the most fundamental theoretical models in nuclear reaction theory [17]. It is applied to describe nucleon-induced reaction cross sections and elastic scattering angular distributions, and to calculate the transmission coefficients and absorption cross sections for compound nucleus decay and the pre-equilibrium emission process. Phenomenological global optical model potentials (OMPs), such as Becchetti–Greenlees (BG)[18], Koning–Delaroche (KD)[19], and Chapel-Hill 89 [20] have been widely employed. In recent years, increased emphasis has been placed on the reliability and rigorous uncertainty quantification of these models [21]. Furthermore, promising progress has been made in developing microscopic OMPs rooted in chiral effective field theory (chiral-EFT), which provide a more fundamental description and facilitate robust theoretical uncertainty quantification [22].
In this paper, the phenomenological BG potential is adopted as the foundational framework [18]. To enhance its predictive accuracy for the p+$ ^{28}$Si system, the potential parameters have been specifically optimized. Following the formalisms established in recent evaluations [23, 24], the expression of the optical potential is defined as:
$\begin{align} V\left(r\right) = V_{\mathrm{R}}\left(r\right) + V_{\mathrm{C}}\left(r\right) + U_{\mathrm{SO}}\left(r\right) + \mathrm{i}\left[ W_{\mathrm{D}}\left(r\right) + W_{\mathrm{S}}\left(r\right) \right],\end{align}$
where the real part consists of the real central potential $V_{\mathrm{R}}(r)$ and the Coulomb potential $V_{\mathrm{C}}(r)$. The imaginary part includes the volume absorption potential $W_{\mathrm{D}}(r)$ and the surface absorption potential $W_{\mathrm{S}}(r)$. The term $U_{\mathrm{SO}}(r)$ represents the spin-orbit interactions.The radial dependence of the real central potential is described by the standard Woods–Saxon form factor:
$\begin{align} V_{\mathrm{R}}\left(r\right) = -\frac{V_{\mathrm{R}}\left(E\right)}{1 + \exp\left(\frac{r - R_{\mathrm{R}}}{a_{\mathrm{R}}}\right)}.\end{align}$
To account for the surface-peaked reaction mechanisms, the imaginary surface absorption potential is constructed using the derivative Woods–Saxon form:
$\begin{align} W_{\mathrm{S}}\left(r\right) = -4W_{\mathrm{S}}\left(E\right) \frac{\exp\left(\frac{r - R_{\mathrm{S}}}{a_{\mathrm{S}}}\right)}{\left[1 + \exp\left(\frac{r - R_{\mathrm{S}}}{a_{\mathrm{S}}}\right)\right]^2}.\end{align}$
Conversely, the imaginary volume absorption potential utilizes the standard Woods–Saxon shape:
$\begin{align} W_{\mathrm{D}}\left(r\right) = -\frac{W_{\mathrm{D}}\left(E\right)}{1 + \exp\left(\frac{r - R_{\mathrm{D}}}{a_{\mathrm{D}}}\right)}.\end{align}$
The complex spin–orbit coupling potential is formulated with a Thomas-like radial derivative, and its explicit dependence on the angular momentum quantum numbers is given by:
$\begin{align} U_{\mathrm{SO}}\left(r\right) & = -2\left(V_{\mathrm{SO}} + \mathrm{i}W_{\mathrm{SO}}\right) \frac{1}{a_{\mathrm{SO}} r} \frac{\exp\left(\frac{r - R_{\mathrm{SO}}}{a_{\mathrm{SO}}}\right)}{\left[1 + \exp\left(\frac{r - R_{\mathrm{SO}}}{a_{\mathrm{SO}}}\right)\right]^2} \nonumber\\ &\quad \times \left[ j\left(j+1\right) - l\left(l+1\right) - s\left(s+1\right) \right],\end{align}$
where $l$, $s$, and $j$ represent the orbital, spin, and total angular momentum quantum numbers of the interacting system, respectively.The Coulomb potential $V_{\mathrm{C}}(r)$ assumes a uniformly charged spherical distribution for the interaction between the incident projectile and the target nucleus:
$\begin{align} V_{\text{C}}\left(r\right) = \begin{cases} \dfrac{z Z \mathrm{e}^2}{2 R_{\mathrm{C}}} \left( 3 - \dfrac{r^2}{R_{\mathrm{C}}^2} \right), & \text{if } ~~r \lt R_{\mathrm{C}}, \\ \dfrac{z Z \mathrm{e}^2}{r}, & \text{if} ~~ r \unicode{x2A7E} R_{\mathrm{C}}, \end{cases}\end{align}$
where $z$ and $Z$ are the charge numbers of the projectile and the target nucleus, respectively. $\mathit{e}$ is the elementary charge, and $R_{\mathrm{C}}$ is the Coulomb radius.Furthermore, to accurately reproduce the scattering observables across different incident energy regions, the energy dependencies of the potential well depths are explicitly parameterized as follows:
$\begin{align}V_{\mathrm{R}}\left(E\right) & = V_0 + V_1 E + V_2 E^2 + V_3 \frac{N-Z}{A} + V_4 Z A^{-1/3},\end{align}$
$\begin{align}W_{\mathrm{S}}\left(E\right) & = \max\left\{0, W_0 + W_1 E + W_2 \frac{N-Z}{A} \right\},\end{align}$
$\begin{align}W_{\mathrm{D}}\left(E\right) & = \max\left\{0, U_0 + U_1 E + U_2 E^2\right\},\end{align}$
where $E$ denotes the incident laboratory energy. The variables $N$, $Z$, and $A$ represent the neutron, charge, and mass numbers of the target nucleus, respectively.The radii are given by 2 )–(5 ) characterize the diffuseness of the real, surface absorption, volume absorption, and spin-orbit couple of the optical potential, respectively.
$\begin{align} R_i = r_i A^{1/3}, \quad i = \mathrm{R, S, D, SO, C},\end{align}$
where $r_{\mathrm{R}}$, $r_{\mathrm{S}}$, $r_{\mathrm{D}}$, $r_{\mathrm{SO}}$ and $r_{\mathrm{C}}$ are the radii of the real part, the surface absorption, the volume absorption, the spin-orbit couple and the Coulomb potential, respectively. The parameters $a_{\mathrm{R}}$, $a_{\mathrm{S}}$, $a_{\mathrm{D}}$, and $a_{\mathrm{SO}}$ appearing in equations (In this work, all potential depths ($V_{\mathrm{R}}, W_{\mathrm{S}}, W_{\mathrm{D}}$, and $V_{\mathrm{SO}}$) as well as the incident energy $E$ are expressed in units of MeV. The corresponding geometric parameters, including the reduced radii $r_i$ and diffuseness $a_i$ are specified in fm.
Based on the formalism presented above, the OMP comprises 22 specific parameters: $V_0$, $V_1$, $V_2$, $V_3$, $V_4$, $W_0$, $W_1$, $W_2$, $U_0$, $U_1$, $U_2$, $V_{\mathrm{SO}}$, $W_{\mathrm{SO}}$, $r_{\mathrm{R}}$, $r_{\mathrm{S}}$, $r_{\mathrm{V}}$, $r_{\mathrm{SO}}$, $r_{\mathrm{C}}$, $a_{\mathrm{R}}$, $a_{\mathrm{S}}$, $a_{\mathrm{V}}$, and $a_{\mathrm{SO}}$. Among these, the values for $V_3$, $W_2$, $V_{\mathrm{SO}}$, $r_{\mathrm{SO}}$, and $a_{\mathrm{SO}}$ are fixed to the universal global parameters originally proposed by Becchetti and Greenlees [18]. Since the angular distributions and excitation functions exhibit non-uniform sensitivities to various optical potential parameters across different spatial angles and incident energy regions, this characteristic allows us to utilize the scarce experimental measurements to break parameter degeneracy. Consequently, it provides robust physical constraints to accurately evaluate and optimize the selection of specific phenomenological parameters. Specifically, the corresponding parameters in their sensitive regions are systematically optimized by reproducing the experimental reaction cross sections, elastic and inelastic scattering angular distributions.
In order to obtain a reliable set of proton OMP parameters for the p+$ ^{28}\text{Si}$ reaction, the APMN code [25] is employed in this work. To address the parameter optimization, the APMN code utilizes an improved steepest descent method to iteratively minimize the $\chi^2$ deviation between the theoretical calculations and the experimental measurements. Through this numerical algorithm, the optimal proton OMP parameters can be automatically searched to best fit the relevant experimental data, including reaction cross sections and elastic scattering angular distributions. The final optimized proton OMP parameters are summarized in table 1. The neutron optical potential parameters are derived from [20], while the optical potentials of charged particles, such as deuteron, triton, $ ^3$He and alpha, are partly adopted from our previous works [23, 26–28], respectively. By solving the Schrödinger equation with the optimized spherical optical potential, the cross sections, including the shape elastic scattering cross section ($\sigma_{\mathrm{se}}$), the reaction cross section ($\sigma_{\mathrm{R}}$), and elastic angular distributions are directly obtained.
Table 1. The optimal proton optical model potential parameters for the p+$ ^{28}$Si reaction. |
| Pameters | Value | Parameters | Value |
|---|---|---|---|
| $V_{0}$ | 55.72 049 | $V_{\mathrm{SO}}$ | 6.20 000 |
| $V_{1}$ | $-$0.32 252 | $W_{\mathrm{SO}}$ | $-$0.17 291 |
| $V_{2}$ | 0.00 007 | $r_{\mathrm{R}}$ | 1.10 715 |
| $V_{3}$ | 24.00 000 | $r_{\mathrm{S}}$ | 1.11 140 |
| $V_{4}$ | 0.52 157 | $r_{\mathrm{D}}$ | 1.14 639 |
| $W_{0}$ | 15.49 838 | $r_{\mathrm{SO}}$ | 1.01 000 |
| $W_{1}$ | $-$0.17 291 | $r_{\mathrm{C}}$ | 1.25 000 |
| $W_{2}$ | 12.00 000 | $a_{\mathrm{R}}$ | 0.77 581 |
| $U_{0}$ | $-$3.44 338 | $a_{\mathrm{S}}$ | 0.73 122 |
| $U_{1}$ | 0.20 210 | $a_{\mathrm{D}}$ | 0.45 981 |
| $U_{2}$ | $-$0.00 044 | $a_{\mathrm{SO}}$ | 0.75 000 |
The theoretical calculations for the direct inelastic scattering cross sections were performed within the distorted wave Born approximation (DWBA) theory [29, 30]. Initially, the optimized OMP parameters are utilized to generate the distorted wave functions, $\chi_{i,f}(\mathbf{r})$, for both the entrance and exit channels. By incorporating the deformation parameter $\beta$ of the excited levels which contribute to the inelastic scattering cross sections, the transition matrix element $T_\mathrm{if}$ between the initial and final states is evaluated. Consequently, the direct inelastic scattering differential cross section is rigorously expressed as:
$\begin{align} \frac{\mathrm{d}\sigma}{\mathrm{d}\Omega} = \frac{\mu_\mathrm{i} \mu_\mathrm{f}}{\left(2\pi \hbar^2\right)^2} \frac{k_\mathrm{f}}{k_\mathrm{i}} |T_\mathrm{if}|^2,\end{align}$
where $\mu_\mathrm{i(f)}$ and $k_\mathrm{i(f)}$ denote the reduced mass and the wave number in the initial (final) channel, respectively. The angle-integrated direct inelastic cross section $\sigma_{\mathrm{dir}}$, is subsequently obtained by integrating the differential cross section over the entire solid angle. The direct scattering contributions from 16 low-lying discrete levels are explicitly considered, with their corresponding physical properties summarized in table 2. Their excitation energies ($E^*$), spins ($J$), and parities ($\pi$) are strictly adopted from the Evaluated Nuclear Structure Data File (ENSDF) [31]. Notably, the deformation parameters ($\beta$) utilized in our calculations represent the nuclear potential deformations, which essentially differ from the pure electromagnetic deformations. Taking the recommended values from the RIPL-3 database as initial references, these empirical $\beta$ parameters were specifically adjusted and optimized to yield the best reproduction of the experimental inelastic scattering angular distributions. Thus, these $\beta$ values are model-dependent and should be used with caution. It is also worth noting that slightly different $\beta$ values can be derived from experimental data [32, 33]. When the excitation energy is higher, the contribution of particle emission is considered in the form of a continuous energy level.Table 2. The excited energies ($E^*$), spins ($J$), parities ($\pi$), and deformation parameters ($\beta$) of the 16 discrete levels that predominantly contribute to direct inelastic scattering. The data is adopted from ENSDF [31], and the optimal deformation parameters $\beta$ are recommended in this paper. |
| $E^*$ (MeV) | $J$ | $\pi$ | $\beta$ |
|---|---|---|---|
| 1.77 903 | 2.0 | + | 0.459 |
| 4.61 786 | 4.0 | + | 0.387 |
| 6.87 879 | 3.0 | — | 0.303 |
| 6.88 765 | 4.0 | + | 0.280 |
| 7.38 059 | 2.0 | + | 0.180 |
| 7.41 626 | 2.0 | + | 0.180 |
| 7.93 345 | 2.0 | + | 0.280 |
| 8.25 874 | 2.0 | + | 0.280 |
| 8.54 356 | 6.0 | + | 0.480 |
| 8.90 480 | 1.0 | — | 0.280 |
| 9.16 468 | 4.0 | + | 0.380 |
| 9.38 155 | 2.0 | + | 0.280 |
| 9.47 949 | 2.0 | + | 0.280 |
| 9.76 452 | 3.0 | — | 0.280 |
| 9.79 595 | 2.0 | + | 0.280 |
| 10.1816 | 3.0 | — | 0.280 |
The unified Hauser–Feshbach and exciton model [34] are used to describe the nuclear reaction equilibrium and pre-equilibrium decay processes. The main objective of this paper is to use the improved Iwamoto–Harada model to provide a comprehensive and improved description of the energy spectra and DDCSs of the light composite particle (including d, t, $ ^3$He and $\alpha$) emissions in the intermediate-energy range [35]. The energy spectrum of emitted composite particle b can be expressed as [36–41]
$\begin{align} & \frac{\mathrm{d}\sigma}{\mathrm{d}\varepsilon_{\mathrm{b}}} \nonumber\\ &\quad = \sum_{J\pi} \sigma_{\mathrm{a}}^{J\pi} \left\{\sum_{ {n = 3}}^{{ {n_{\max}}}} P^{J\pi}\left(n\right) \frac{W_{\mathrm{b}}^{J\pi}\left(n, E^*, \varepsilon_{\mathrm{b}}\right)}{W_T^{J\pi}\left(n, E^*\right)} + Q^{J\pi} \frac{W_b^{J\pi}\left(E^*, \varepsilon_b\right)}{W_T^{J\pi}\left(E^*\right)} \right\},\end{align}$
where $W_{\mathrm{b}}^{J\pi}(n, E^*, \varepsilon_{\mathrm{b}})$ is the emission rate of particle b at the $n$th exciton state with outgoing energy $\varepsilon_{\mathrm{b}}$, which is obtained by $\begin{align} W_{\mathrm{b}}^{J\pi}\left(n, E^*, \varepsilon_{\mathrm{b}}\right) & = \frac{2s_{\mathrm{b}}+1}{\pi^2 \hbar^3} \mu_{\mathrm{b}} \sigma_{\mathrm{b}}^{J\pi}\left(\varepsilon_{\mathrm{b}}\right) \varepsilon_{\mathrm{b}} F_{lm}^{\mathrm{b}}\left(\varepsilon_{\mathrm{b}}\right) Q_{\mathrm{b}}\left(n\right) \nonumber\\ & \quad \times \frac{\omega\left(p-1, h, E^* - \varepsilon_{\mathrm{b}} - B_{\mathrm{b}}\right)}{\omega\left(p, h, E^*\right)},\end{align}$
where $s_{\mathrm{b}}$, $\mu_{\mathrm{b}}$, $\sigma_{\mathrm{b}}^{J\pi}$ stand for the spin, reduced mass and the inverse cross section of b particle, respectively. $F_{lm}^{\mathrm{b}}(\varepsilon_{\mathrm{b}})$ is the formation probability of b particle with the configuration $[lm]$ from the improved Iwamoto–Harada model, while $F_{lm}^{\mathrm{b}}(\varepsilon_{\mathrm{b}}) = 1$ for nucleon. $Q_{\mathrm{b}}(n)$ stands for the combination factor to account the nucleon type composed the cluster for memory by excitation system. $\omega(n, E^*)$ stands for the exciton state density. The explicit expressions of $F_{lm}^{\mathrm{b}}(\varepsilon_{\mathrm{b}})$, $Q_{\mathrm{b}}(n)$ and $\omega(n, E^*)$ can be found in [35]. The upper limit $n_{\max}$ represents the dynamically determined equilibrium exciton number at which the system reaches statistical equilibrium.The total emission rate can be expressed as summing over all light charged particles
$\begin{equation} W_T^{J\pi}\left(n, E^*\right) = \sum \int W_{\mathrm{b}}^{J\pi}\left(n, E^*, \varepsilon_{\mathrm{b}}\right) \mathrm{d}\varepsilon_{\mathrm{b}}.\end{equation}$
The occupation probability $P^{J\pi}(n)$ of the $n$th exciton state in the $J\pi$ channel is given by solving the $j$-dependent exciton master equation to conserve the angular momentum in the pre-equilibrium reaction processes.$W_{\mathrm{b}}^{J\pi}(E^*,\varepsilon_{\mathrm{b}})$ is the emission rate of the emitted particle b at the equilibrium state with outgoing kinetic energy $\varepsilon_{\mathrm{b}}$, of which the contributions are from the discrete states (table 2) and continuum states of the excited nuclei. The Ignatyuk nuclear level densities [42] are used, which include the washing out of shell effects with increasing excitation energy, collective excitations as well as single particle ones, depart from more traditional ones, and are matched continuously onto low-lying experimental discrete levels (table 2). The Ignatyuk model for describing the statistical level density properties of excited nuclei is particularly appropriate for the relatively high energies. For energy levels above the highest excited states, the level density theoretical model is employed. In the low-energy region (below 20 MeV), the Gilbert–Cameron formula [43] is adopted. $E^*$ is the excited energy of compound nucleus, and $Q^{J\pi} = 1-\sum_{n = 3}^{n_{\max}}P^{J\pi}(n)$ is the occupation probability of the equilibrium state.
In order to simplify the calculations, the angular dependent part of the double differential cross sections for emission neutron, proton, deuteron, triton, $ ^3$He and alpha are obtained from Kalbach phenomenological approach [44, 45]. It is based on a systematical study of a wide variety of experimental data. The parameterization established for the DDCSs as a function of the total energy spectra is given by the equation:
$\begin{align} \dfrac{\mathrm{d^2}\sigma}{\mathrm{d}\Omega\mathrm{d}\varepsilon_{\text{b}}} = \dfrac{1}{4\pi}\frac{\mathrm{d}\sigma}{\mathrm{d}\varepsilon_{\mathrm{b}}}\dfrac{a}{\sinh \left(a\right)} \left[ \cosh \left(a \cos \theta\right) + f_\mathrm{PE} \sinh \left(a \cos \theta\right) \right].\end{align}$
In this expression, $\theta$ is the emission angle in the center of mass frame, and the term $a$ is the slope parameter depending on the incident particle type and energy, the target nucleus and the exit channel. The $f_\mathrm{PE}$ parameter is the fraction of particle emission away from equilibrium. It will be calculated using the formula: $\begin{equation} f_\mathrm{PE} = \frac{\left(\mathrm{d}\sigma/\mathrm{d}\varepsilon_b\right)_\mathrm{PE}}{\left(\mathrm{d}\sigma/\mathrm{d}\varepsilon_{\mathrm{b}}\right)} = \frac{\left(\mathrm{d}\sigma/\mathrm{d}\varepsilon_b\right)_\mathrm{PE}}{\left(\mathrm{d}\sigma/\mathrm{d}\varepsilon_b\right)_\mathrm{PE} + \left(\mathrm{d}\sigma/\mathrm{d}\varepsilon_b\right)_{EQ}},\end{equation}$
where the PE and EQ symbols refer respectively to pre-equilibrium and equilibrium emissions. Using the above calculated energy spectra equation (3. Theoretical results and analysis
Owing to the experimental data of the energy spectra for the p+$ ^{28}$Si reaction currently unavailable, the calculated results of this paper are compared with the experimental data obtained from the natural Si target [14]. Figure 1 shows the comparisons of calculated results with the experimental energy spectra of the proton (a), deuteron (b), triton (c), $ ^3$He (d), and $\alpha$ (e) at incident energies 26.5 MeV, 46.5 MeV, and 62.9 MeV, respectively. The calculated results are in rather good accordance with the experimental data. The energy spectra of deuteron, triton, and $ ^3$He emissions are primarily attributed to pre-equilibrium reactions. In contrast, the energy spectra of the proton and $\alpha$ are mainly dominated by equilibrium emission in the lower incident energy region. In figure 1(a), the structural characteristics observed in the proton energy spectra within the high emission energy region arise from the direct inelastic scattering processes. Similarly, for deuterons, as depicted in figure 1(b), the observed energy spectra within the high emission energy region may stem from the contributions of the knockout reaction mechanism which are not enough considered in this paper. Thus, several structural peaks of the outgoing deuteron in the present calculations do not match the experimental data very well. Additionally, the energy spectra derived from TALYS-2.0 are also shown in figure 1. One can see that our results are superior to those calculated by TALYS-2.0 except several structural peaks of the outgoing deuterons.
Figure 1. The energy spectra of the outgoing proton (a), deuteron (b), triton (c), $ ^3$He (d), and $\alpha$ (e) for the p+$ ^{28}$Si reaction compared with the experimental data at incident energy of 26.5 MeV, 46.5 MeV, and 62.9 MeV, respectively. The scattered symbols denote the experimental data taken from [14]. The red solid and blue dash-dot curves denote the calculated results of this work and the TALYS-2.0, respectively. The data is shifted downward by factors of $10^{0}$, $10^{-1}$, $10^{-2}$ respectively. |
The DDCSs of outgoing light charged particles for natural Si were measured [14]. The comparisons of the calculated DDCSs of the outgoing protons with the incident energies of 26.5 MeV, 46.5 MeV, and 62.9 MeV are shown in figure 2, respectively. The vertical axes of figure 2 employ a logarithmic scale to display the broad dynamic range of the DDCSs. Although the experimental uncertainties are symmetric in linear space, the logarithmic transformation results in asymmetric error bars due to nonlinear scaling. To ensure clarity while preserving scientific accuracy, only the upper error bars are displayed. The calculated DDCSs are in good agreement with the experimental data, except for the high-energy region of the outgoing particles, due to the structural effects in this region are very pronounced. Additionally, as illustrated in figure 2(a), the calculated results for proton emission energy in the range of 4–8 MeV are lower than the experimental data when the incident energy is 26.5 MeV. This is primarily due to the insufficient consideration of the proton emission from equilibrium state in this paper, especially when compared to other higher incident energy scenarios. The calculations also reveal that the contribution of pre-equilibrium reactions decreases with increasing outgoing angle. Furthermore, from the data of incident protons with higher energies as shown in figures 2(b) and (c), we can observe that the significance of the pre-equilibrium emission becomes increasingly prominent as the incident energy rises.
Figure 2. The double-differential cross sections of outgoing proton for p+$ ^{28}$Si reaction with different outgoing angle at incident energy of 26.5 MeV (a), 46.5 MeV (b), and 62.9 MeV (c), respectively. The scattered symbols denote the experimental data taken from [14]. The red solid curves denote the calculated results of this work. The data is shifted downward by factors of $10^{0}$, $10^{-1}$, $10^{-2}$, and so on. |
Figure 3 shows the DDCSs of outgoing deuterons with a different outgoing angle for the p+$ ^{28}$Si reaction. The experimental data from the natural Si target [14] are used to compare. Similarly to figure 2, figure 3 employs a logarithmic scale and displays solely upper error bars to maintain clarity in depicting asymmetric uncertainties near the axis origin. In figure 3, the structural characteristics are observed in the higher emission energy regions. In these regions, the calculated results cannot reproduce the structural peaks, as well as the energy spectra of the outgoing deuterons. This deviation likely stems from the exclusion of direct reaction in the (p, d) channel within the current framework. Because the direct effects of this reaction channel cannot be readily incorporated into the theoretical framework of this paper through the DWBA theory alone.
Figure 3. Same as figure 2, but for the outgoing deuteron. |
The calculated DDCSs of the outgoing tritons are compared with the experimental data [14] at incident proton energies of 26.5 MeV, 46.5 MeV, and 62.9 MeV, as shown in figure 4. Lower error bars are omitted in figure 4 for identical reasons as well as in figures 2 and 3. The magnitude and shape of the calculated results are in good agreement with the experimental data.
Figure 4. Same as figure 2, but for the outgoing triton. |
The calculated DDCSs of the outgoing alphas are compared with the experimental data [14] at incident proton energies of 26.5 MeV, 46.5 MeV, and 62.9 MeV, as shown in figure 5. Lower error bars are also omitted in this figure for identical reasons as well as figures 2–4. The magnitude and shape of the calculated results are reasonable, although the agreement is not as good as those of the proton, deuteron, and triton in the lower outgoing energy region. This could be attributed to the fact that the contribution of alpha particles primarily stems from the equilibrium emission process, while the contributions of other lighter charged particles mainly originate from pre-equilibrium processes.
Figure 5. Same as figure 2, but for the outgoing alpha. |
We also conducted model calculations for the p+$ ^{28}$Si reaction using the open-source TALYS-2.0 code [12]. Overall, the calculated results provide a reasonable reproduction of the experimental DDCSs of the light charged particles [14]. Nevertheless, the DDCSs of tritons were not published in previous papers, and there are no experimental DDCSs (despite the experimental energy spectrum data) of the outgoing $ ^3$He. To gain a deeper insight into the differences between the TALYS-2.0 code and this work, we further compared the DDCSs of the emitted light charged particles. Figure 6 shows the results of outgoing light charged particles at an incident energy of 62.9 MeV with an angle of 60$ ^\circ$ calculated by TALYS-2.0 and this work, respectively. The black points represent the experimental data derived from [14]. The red solid lines and the blue dash-dot lines denote the results from this work and TALYS-2.0, respectively.
Figure 6. The comparisons of the double-differential cross sections between the results of this work and those of the TALYS-2.0 for the outgoing proton (a), deuteron (b), triton (c), and $\alpha$ (d) for the p+$ ^{28}$Si reaction at incident energy of 62.9 MeV with an angle of 60$ ^{\circ}$. The scattered symbols denote the experimental data taken from [14]. The red solid and blue dash-dot lines denote the result of this work and TALYS-2.0, respectively. |
It should be acknowledged that global nuclear reaction codes, such as TALYS-2.0, exhibit remarkable universal applicability and generally provide reliable predictions across the entire nuclear chart. However, relying solely on global phenomenological parameters can inevitably lead to unpredictable uncertainties or limited accuracy for specific individual reactions. In this context, the present work demonstrates a highly effective and targeted evaluation methodology. By optimizing the specific OMP parameters based on limited but critical experimental measurements, and coupling them with robust physical models (especially the explicit treatment of light-charged particle pre-equilibrium emissions), significantly more accurate nuclear data for the p+$ ^{28}\text{Si}$ reaction are obtained. This tailored optimization approach is highly essential and practically valuable to satisfy the stringent precision requirements in advanced nuclear engineering applications.
4. Conclusions
This study establishes a comprehensive theoretical framework for the p+$ ^{28}$Si reaction up to 200 MeV using the MEND code, demonstrating robust predictive capability through systematic validation against experimental data as our previous works [46–50]. The direct, pre-equilibrium, and equilibrium reaction contributions for the different outgoing particles are analyzed in detail. The model successfully reproduces the energy spectra and the DDCSs of the light charged particles at incident energies of 26.5, 46.5 and 62.9 MeV except several structural peaks of the outgoing deuterons. The results show that the adequacy of this procedure is verified. Furthermore, the generated data on the basis of this work can provide essential nuclear inputs for spacecraft radiation shielding design and SiC electronics reliability simulation.