Dynamical model: UrQMD

The UrQMD model [39-42] is a typical transport model used for microscopic many-body nonequilibrium dynamics. In this model, each nucleon is represented by a Gaussian wave packet of a certain width L. Empirically, L = 2 fm² are chosen to simulate Au+Au collisions. The coordinates and momentum of each nucleon are propagated using Hamilton's equations of motion. The total Hamiltonian of the system $\langle H\rangle$ consists of the kinetic energy T and the effective interaction potential energy U, which includes the Skyrme potential energy Uρ, Coulomb energy U_Coul, and momentum-dependent potential energy U_md: $\langle H\rangle =T+U_{\rho }+U_{\text{Coul}}+U_{\text{md}}\mathrm{.}$ (1) To study the HICs at intermediate energies, the Skyrme energy density functional was introduced in the same manner as that in the improved quantum molecular dynamics (ImQMD) model [43, 44]. The local and momentum-dependent potential energies can be written as $U_{\rho \mathrm{,}\text{md}}={\displaystyle \int }{u}_{\rho \mathrm{,}\text{md}}\text{d}r$, where $\begin{array}{c}{u}_{\rho }=\frac{\alpha }{2}\frac{{\rho }^2}{{\rho }_0}+\frac{\beta }{\gamma +1}\frac{{\rho }^{\gamma +1}}{{\rho }_0^{\gamma }}+\frac{g_{\text{sur}}}{2{\rho }_0}{\left(\nabla \rho \right)}^2\\ +\frac{g_{\text{sur,iso}}}{2{\rho }_0}{\left[\nabla \left({\rho }_n-{\rho }_p\right)\right]}^2\\ +\left[A_{\text{sym}}\left(\frac{\rho }{{\rho }_0}\right)+B_{\text{sym}}{\left(\frac{\rho }{{\rho }_0}\right)}^{\eta }+C_{\text{sym}}{\left(\frac{\rho }{{\rho }_0}\right)}^{5/3}\right]{\delta }^2\rho \mathrm{,}\end{array}$ (2) and $u_{\text{md}}=t_{\text{md}}{{\displaystyle \ln }}^2\left[1+a_{\text{md}}{\left(p_i-p_j\right)}^2\right]\frac{\rho }{{\rho }_0}\mathrm{.}$ (3) Here, δ = $\left({\rho }_{\text{n}}-{\rho }_{\text{p}}\right)/\left({\rho }_{\text{n}}+{\rho }_{\text{p}}\right)$ denotes the isospin asymmetry defined by the neutron (ρ_n) and proton (ρ_p) densities. In the present work, the normal nuclear matter density ρ₀ = 0.16 fm^-3, α = -393 MeV, β = 320 MeV, γ = 1.14, g_sur = 19.5 MeV fm², g_sur,iso = -11.3 MeV fm², A_sym = 20.4 MeV, B_sym = 10.8 MeV, C_sym = -9.3 MeV, η = 1.3, t_md = 1.57 MeV, and a_md = 500 (GeV/c)^-2 are chosen, corresponding to a soft and momentum-dependent (SM) equation of state with the incompressibility K₀ = 200 MeV and the symmetry-energy slope parameter L = 80.95 MeV, which are in their commonly accepted regions [42, 45-49].

In this work, the in-medium nucleon-nucleon (NN) elastic cross-section is considered as a density- and momentum-dependent medium correction factor $\mathcal{F}(\rho \mathrm{,}p)$ times the free nucleon-nucleon cross-section ${\sigma }_{\text{NN}}^{\text{free}}$, which is available from the experimental data. ${\sigma }_{\text{NN}}^{\text{in}-\text{medium}}=\mathcal{F}\left(\rho \mathrm{,}p\right)*{\sigma }_{\text{NN}}^{\text{free}}\mathrm{,}$ (4) with $\mathcal{F}\left(\rho \mathrm{,}p\right)=\{\begin{array}{l}{f}_0\mathrm{,}{p}_{\text{NN}}>1\text{G}eV/c\mathrm{,}\hfill \\ \frac{\lambda +(1-\lambda )e^{-\frac{\rho }{\zeta {\rho }_0}}-f_0}{1+{(p_{\text{NN}}/p_0)}^{\kappa }}+f_0\mathrm{,}{p}_{\text{NN}}\le 1\text{GeV}/c\mathrm{.}\hfill \end{array}$ (5) The parameters f₀ = 1, λ = 1/6, ζ = 1/3, p₀ = 0.3 GeV/c, and κ = 8 were adopted, with p_NN being the momentum in the two-nucleon center-of-mass frame. Further details regarding the correction factor $\mathcal{F}(\rho \mathrm{,}p)$ can be found in Refs. [41, 42, 50]. In addition to the in-medium NN cross-section, Pauli blocking plays an important role in determining the rate of collisions. The Pauli blocking treatment is the same as that described in Refs. [41, 19]; it involves two steps. First, the phase space densities fi and fj of the two outgoing particles for each NN collision are calculated. $\begin{array}{c}{f}_i=\frac{1}{2}{\displaystyle \sum _k\frac{1}{{(\pi \hslash )}^3}}\exp \left[-\frac{{(r_i-r_k)}^2}{(2{\sigma }_r^2)}\right]\\ \exp \left[-\frac{{\left(p_i-p_k\right)}^2\cdot 2{\sigma }_r^2}{{\hslash }^2}\right]\mathrm{.}\end{array}$ (6) Here, ${\sigma }_r^2$ denotes the width parameter of the wave packet. ${\sigma }_r^2$ = 2 fm² is used in this study. k represents the nucleons of the same type around the outgoing nucleon i (or j). The following two criteria were simultaneously considered: $\frac{4\pi }{3}{r}_{ik}^3\frac{4\pi }{3}{p}_{ik}^3\ge 2{\left(\frac{h}{2}\right)}^3;$ (7) $P_{\text{block}}=1-\left(1-f_i\right)\left(1-f_j\right)<\xi \mathrm{.}$ (8) Here, rik and pik denote the relative distance and momentum between the nucleons i and k. The above conditions were also considered for nucleon j. The symbol ξ denotes a random number between 0 and 1.

In addition, the isospin-dependent minimum spanning tree (iso-MST) algorithm was used to construct the clusters [36-38]. Nucleon pairs with relative distances smaller than R₀ and relative momenta smaller than P₀ = 0.25 GeV/c are considered to be bound in a fragment. Here, $R_0^{\text{pp}}$ = 2.8 fm and $R_0^{\text{np}}$ = $R_0^{\text{nn}}$ = 3.8 fm were used for proton-proton and neutron-neutron (neutron-proton) pairs, respectively. A fairly good agreement between the recently published experimental data (including collective flows, nuclear stopping power, and Hanbury–Brown–Twiss interferometry) and UrQMD model calculations was achieved with appropriate choices of the above parameters [45, 35, 42, 51, 52].

Statistical model: GEMINI++

The statistical GEMINI++ decay model [30, 53] is a Monte Carlo method used to simulate the decay of excited fragments, which includes light-particle evaporation, symmetric and asymmetric fission, and all possible binary decay modes. It has been widely used to treat the de-excitation of the primary fragments in studies of nuclear reactions at low and intermediate [23, 25, 54-60] and high energies [61-63]. The parameters of the shell-smoothed level density in GEMINI++ were set to their default values, which were k₀ = 7.3 MeV and $k_{\infty }$ = 12 MeV. A more detailed discussion of the GEMINI++ parameters can be found in Refs. [60, 53, 64].

There are four inputs for the GEMINI++ code (with default parameter settings): the mass number A, charge number Z, excitation energy E^*, and angular moment $\overrightarrow{L}$ of the primary fragment. After the UrQMD transport and iso-MST coalescence processes, the excitation energy of a fragment is calculated as follows: $E^{*}=E_{\text{bind}}^{\text{excited}}-E_{\text{bind}}^{\text{ground}}\mathrm{,}$ (9) where $E_{\text{bind}}^{\text{excited}}$ is the binding energy of the primary fragment, that is, the sum of the potential and kinetic energies of all nucleons that belong to a fragment. $E_{\text{bind}}^{\text{ground}}$ is the binding energy of the ground state obtained from nuclear data table AME2020 [65, 66]. When the excitation energy of the fragment was less than or equal to zero, it was regarded as a bound state and did not decay further. The angular momentum of the primary fragments was calculated using classical mechanics [56]. $\overrightarrow{L}={\displaystyle \sum _i{\overrightarrow{r}}_i\times {\overrightarrow{p}}_i}\mathrm{,}$ (10) where ${\overrightarrow{r}}_i$ and ${\overrightarrow{p}}_i$ are the coordinate and momentum vectors, respectively, of the i-th nucleon in the primary fragments in the c.m. frame of the fragment. The total angular momentum is the sum of all the nucleons in the primary fragments. We verified that the results remained almost unchanged if the contribution of the angular momentum in GEMINI++ was not included.

Observables

In this study, we focus on collective flows and nuclear stopping power, which are commonly used observables, and propose investigating the properties of dense nuclear matter. The directed (v₁) and elliptic (v₂) flows are defined as [67] $v_1=\langle \frac{p_x}{\sqrt{p_x^2+p_y^2}}\rangle \mathrm{,}{v}_2=\langle \frac{p_x^2-p_y^2}{\sqrt{p_x^2+p_y^2}}\rangle \mathrm{.}$ (11) Here, px and py are the two components of the transverse momentum p_t = $\sqrt{p_x^2+p_y^2}$. Angle brackets indicate the average of all considered particles from all events. The nuclear stopping power R_E, which characterizes the transparency of the colliding nuclei, can be defined as [68] $R_{\text{E}}=\frac{{\displaystyle \sum }{E}_{\perp }}{2{\displaystyle \sum }{E}_{\parallel }}\mathrm{,}$ (12) where $E_{\perp }$ ($E_{\parallel }$) is the c.m. transverse (parallel) energy, and the sum runs over all considered particles. Another quantity vartl, which was proposed by the FOPI Collaboration [69], is widely used to measure stopping power. It is defined as the ratio of the variances of the transverse to those of the longitudinal rapidity distribution and is expressed as $vartl=\frac{\langle y_x^2\rangle }{\langle y_z^2\rangle }\mathrm{,}$ (13) where $\langle y_x^2\rangle$ and $\langle y_z^2\rangle$ are the variances in the rapidity distributions of the particles in the x and z directions, respectively.

Rapidity distribution

The rapidity distributions of hydrogen (Z = 1) and helium isotopes (³He and ⁴He) in ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon and impact parameter b = 5.5–7.5 fm are shown in Fig. 1. The left, middle, and right panels show the calculated results when the switching time of the UrQMD model was set to 200, 250, and 300 fm/c, respectively. When considering the decay of the excited primary fragments, the yields of Z = 1 particles are slightly increased, whereas those of ³He and ⁴He clusters are clearly increased. The yields of ³He and ⁴He are approximately 1.5 and 6 times greater than those obtained from the simulations without considering sequential decay. It is known that the yields of light clusters such as ³He and ⁴He obtained using QMD-type codes are significantly smaller than experimentally measured yields. As discussed in our previous work [44], the yield of ³He calculated using the UrQMD model was approximately three times lower than that of the experimental data. The inclusion of GEMINI++ improved the description of light cluster production to some extent; however, the yields of light clusters were still underestimated. Additionally, the inclusion of GEMINI++ may have led to an overestimation of the yield of free nucleons. To solve these problems further, additional issues should be considered, such as the spin degree of freedom and the dynamical production process [7].

Fig. 1Full-screen View

(Color online) Rapidity distributions of Z = 1 [panels (a1–a3)], ³He [panels (b1–b3)], and ⁴He [panels (c1–c3)] particles as functions of the reduced rapidities (y₀ = y_z/y_pro) from peripheral (b = 5.5–7.5 fm) ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon at different reaction times. The solid and dashed lines represent the results obtained from UrQMD simulations with and without considering the decay of the primary fragments using GEMINI++, respectively.

Collective flows

The collective flows of Z = 1 particles as functions of y₀ in ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon are shown in Fig. 2 (directed flow v₁ [left panel] and elliptic flow v₂ [right panel]). The simulations were terminated at 200 fm/c. Clear differences are observed in the distributions of the collective flows between simulations with and without considering the decay of the primary fragments, especially in the case of elliptic flow v₂ from peripheral (b = 5.5–7.5 fm) collisions. From Fig. 1, the difference in the yields of Z = 1 particles between simulations with and without considering sequential decay is due to the particles produced from the decay of excited primary fragments, which have larger mass. These heavier excited primary fragments are usually produced at the target/projectile rapidity; thus, the directed flow at the target/projectile rapidity is significantly affected by sequential decay. Increasing the impact parameter results in more large-mass primary fragments being produced and the effects of sequential decay on the collective flows becoming more evident.

Fig. 2Full-screen View

(Color online) Reduced rapidity y₀ distributions of the directed (elliptic) flow v₁ (v₂) from semi-central (peripheral) ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon for Z = 1 particles. The blue circles denote the calculated results without primary fragment decay, whereas the olive down-triangles represent the calculated results with primary fragments de-excited.

To quantitatively evaluate the influence of sequential decay on the collective flows, as shown in Fig. 3, the v₁ slope (top panels) and v₂ (bottom panels) at mid-rapidity for Z = 1 particles (left panels) and free protons (right panels) are calculated at different stopping times and compared with the experimental data. The slope of directed flow v₁₁ and elliptic flow v₂₀ at mid-rapidity are extracted by assuming $v_1(y_0)=v_{10}+v_{11}\cdot y_0+v_{13}\cdot y_0^3$ and $v_2(y_0)=v_{20}+v_{22}\cdot y_0^2+v_{24}\cdot y_0^4$ in the range of $\left|y_0\right|<0.4$, in the same manner as that in the experimental report [69-71]. The left panels show the results from the simulations at E_lab = 50 MeV/nucleon, whereas the right panels show the results from the simulations at E_lab = 400 MeV/nucleon. There is an obvious difference in v₁₁ and v₂₀ between the simulations with and without considering the decay of primary fragments at E_lab = 50 MeV/nucleon. By considering the decay of the primary fragments, the value of v₁₁ decreases, while that of v₂₀ increases, and both are closer to the experimental data. This can be understood from the fact that both directed and elliptic flows are stronger for large-mass fragments [70, 69, 33]. The decay of fragments is isotropic in the statistical model. Thus, the light nuclei produced from the decay of primary fragments retain information about the primary fragments. Furthermore, at E_lab = 50 MeV/nucleon, v₁₁ first decreases with time and then saturates (increases) after 200 fm/c in the case without (with) the sequential decay effect. This can be explained by the competition between sequential decay and dynamical effects. At approximately t = 150 fm/c, the nucleons are located in a low-density environment, and the final state interactions are attractive; thus, a decrease v₁₁ is observed. When considering GEMINI++ at t = 300 fm/c, there are only a few large-mass fragments that can emit hydrogen isotopes, and the effect of sequential decay on v₁₁ gradually weakens. A comparison of the sequential decay effects at E_lab = 50 MeV/nucleon with those at E_lab = 400 MeV/nucleon shows that the sequential decay effects are more pronounced at lower beam energies. This is understandable because more large-mass fragments, which contribute to the production of light clusters, are formed at lower beam energies.

Fig. 3Full-screen View

(Color online) Directed flow slope v₁₁ and elliptic flow v₂₀ from ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 50 MeV/nucleon (for hydrogen isotopes (Z = 1)) and 400 MeV/nucleon (for free protons with u_t0 gt;0.8) with different stopping times. The INDRA and FOPI experimental data taken from Refs. [69-71] are indicated by shaded bands.

Figure 4 shows the beam energy dependence of v₁₁ (top panel) and v₂₀ (bottom panel) for semi-central ¹⁹⁷Au+¹⁹⁷Au collisions. The value of v₁₁ from the simulations without considering sequential decay is higher than that of the experimental data at lower energies. When sequential decay is included, v₁₁ is driven down, and v₂₀ is pulled up at low energies, which is consistent with the results shown in Fig. 3. Again, both v₁₁ and v₂₀ are considerably closer to the experimental data for all investigated energies. Moreover, the effects of sequential decay on the collective flows decrease with increasing beam energy. With increasing beam energy, the reaction becomes more violent, and the multiplicities of large-mass primary fragments decrease.

Fig. 4Full-screen View

(Color online) Beam energy dependence of directed flow slope v₁₁ (panel a) and elliptic flow v₂₀ (panel b) from semi-central ¹⁹⁷Au+¹⁹⁷Au collisions. As shown, the solid stars represent the data from INDRA Collaboration [70] for Z = 1 particles, and the open stars are the results from FOPI Collaboration [69, 71] for free protons with u_t0 gt;0.8 cut. The lines with different symbols are the results of calculations with and without considering the sequential decay. Where error bars are not shown, they fall within the symbols.

For a systematic investigation of the effects of sequential decay on the collective flows, the transverse 4-velocities u_t0 dependence of the v₁ of light charged particles in semi-central ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 120, 150, and 250 MeV/nucleon are calculated with and without considering the decay of primary fragments and compared with the experimental data [69], as shown in Fig. 5. The effects of sequential decay on v₁ (u_t0) of particles with mass number A lt;4 are relatively weak; it only affects v₁ at low u_t0. By considering the sequential decay effect, v₁ of α particles is obviously influenced and can reproduce the experimental data reasonably well at higher u_t0. This is because after the decay of excited primary fragments is considered, the multiplicity of α particles is enhanced, as shown in Fig. 1. The α particles produced by the decay of the excited primary fragments have a large flow effect because they inherit the flow information of the excited primary fragments. The remaining discrepancies in u_t0-dependent v₁ of light clusters may be due to simplifications in the initial wave function of particles (nucleons and possible clusters) and quantum effects in two-body collisions, as well as the lack of a dynamical production process, which requires further study [7].

Fig. 5Full-screen View

(Color online) Directed flow v₁ of protons (squares), deuterons (circles), A = 3 clusters (triangles), and α-particles (rhombuses) versus u_t0 in ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 120, 150, and 250 MeV/nucleon with b = 3.4–6.0 fm. Results from calculations with and without sequential decay are represented by solid and open symbols. FOPI experimental data (red symbols) are taken from Ref. [69]. The shaded areas represent the corresponding error bars.

Nuclear stopping power

The nuclear stopping power characterizes the transparency of the colliding nuclei and can provide insight into the rate of equilibration of the colliding system. Thus, it is meaningful to explore the effects of sequential decay on the nuclear stopping power. The yield distributions of deuterons and ⁴He clusters as functions of the reduced longitudinal ($y_z/y_{\text{pro}}$) and transverse ($y_x/y_{\text{pro}}$) rapidities for central ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon with and without considering sequential decay are shown in Fig. 6. The values of the calculated R_E are shown in each panel, and they slightly decrease; that is, the nuclear stopping power is depressed when considering the sequential decay. Notably, the values of R_E are slightly influenced by sequential decay; however, the yield spectrum of each particle is obviously influenced, especially for ⁴He clusters. In conjunction with the results presented in Fig. 1, it can be concluded that sequential decay has a strong effect on the particle yield and yield distribution; however, the nuclear stopping power is only weakly influenced because the de-excitation of the primary fragments in the GEMINI++ code is isotropic.

Fig. 6Full-screen View

(Color online) Yield distributions of ²H and ⁴He particles as functions of the reduced longitudinal and transverse rapidities from central (b = 0–2 fm) ¹⁹⁷Au+¹⁹⁷Au collisions at E_lab = 40 MeV/nucleon. Results from calculations with and without the sequential decay are shown in the left and right panels, respectively. The corresponding values for nuclear stopping power R_E are also indicated in each panel.

The nuclear stopping power R_E (E_lab = 40 MeV/nucleon, left panel) and vartl (E_lab = 400 MeV/nucleon, right panel) for free protons from central (b = 0–2 fm) ¹⁹⁷Au+¹⁹⁷Au collisions calculated at different stopping times are shown in Fig. 7. Similar to the collective flows shown in Fig. 3, the effects of sequential decay on the nuclear stopping power at lower energies are more obvious than those at higher energies. In all cases, the results were saturated and were close to or covered by the experimental data above 200 fm/c.

Fig. 7Full-screen View

(Color online) R_E and vartl for free protons from central (b = 0–2 fm) ¹⁹⁷Au+¹⁹⁷Au collisions with E_lab = 40 and 400 MeV/nucleon at different stopping times. The symbols sets are identical to those in Fig. 1. The shaded bands are the corresponding INDRA and FOPI experimental data taken from Refs. [72, 73].

Figure 8 shows the degree of nuclear stopping (R_E or vartl) in the central ¹⁹⁷Au+¹⁹⁷Au collisions as a function of the beam energy. It is clear that both R_E and vartl decrease and are much closer to the experimental data when the GEMINI++ code is applied, and the difference in the values of the nuclear stopping power gradually disappears with increasing beam energy. With the same impact parameter (centrality), more violent reactions occurred at higher beam energies, and less heavily excited primary fragments were produced. The difference with free protons is that they are produced from the decay of heavy excited primary fragments, which usually have a weak stopping power (smaller values of R_E and vartl) [74, 73]. As a result, the nuclear stopping power was slightly reduced by including the GEMINI++ code.

Fig. 8Full-screen View

(Color online) Beam-energy dependence of R_E (INDRA) and vartl (FOPI) for free protons from central (b = 0–2 fm) ¹⁹⁷Au+¹⁹⁷Au collisions. The symbols sets are identical to those in Fig. 5. Calculations with and without sequential decay are compared with the INDRA and FOPI experimental data, which are taken from Refs. [72, 73], respectively.