Modeling overview

The Geant4 Binary Cascade is an intranuclear cascade propagating primary and secondary particles in a nucleus. Interactions are between a primary or secondary particle and an individual nucleon of the nucleus, leading to the name Binary Cascade. Cross section data are used to select collisions. Where available, experimental cross sections are used by the simulation. Propagating of particles is the nuclear field is done by numerically solving the equation of motion. The cascade terminates when the average and maximum energy of secondaries is below threshold. The remaining fragment is treated by precompound and de-excitation models documented in the corresponding chapters.

The transport algorithm

For the primary particle an impact parameter is chosen random in a disk outside the nucleus perpendicular to a vector passing through the center of the nucleus coordinate system an being parallel to the momentum direction. Using a straight line trajectory, the distance of closest approach $d_i^{min}$ to each target nucleon $i$ and the corresponding time-of-flight $t_i^d$ is calculated. In this calculation the momentum of the target nucleons is ignored, i.e. the target nucleons do not move. The interaction cross section ${\sigma}_i$ with target nucleons is calculated using total inclusive cross-sections described below. For calculation of the cross-section the momenta of the nucleons are taken into account. The primary particle may interact with those target nucleons where the distance of closest approach $d_i^{min}$ is smaller than $d_i^{min} < \sqrt{\frac{\sigma_i}{\pi}}$. These candidate interactions are called collisions, and these collisions are stored ordered by time-of-flight $t_i^d$. In the case no collision is found, a new impact parameter is chosen.

The primary particle is tracked the time-step given by the time to the first collision. As long a particle is outside the nucleus, that is a radius of the outermost nucleon plus $3fm$, particles travel along straight line trajectories. Particles entering the nucleus have their energy corrected for Coulomb effects. Inside the nucleus particles are propagated in the scalar nuclear field. The equation of motion in the field is solved for a given time-step using a Runge-Kutta integration method.

At the end of the step, the primary and the nucleon interact suing the scattering term. The resulting secondaries are checked for the Fermi exclusion principle. If any of the two particles has a momentum below Fermi momentum, the interaction is suppressed, and the original primary is tracked to the next collision. In case interaction is allowed, the secondaries are treated like the primary, that is, all possible collisions are calculated like above, with the addition that these new primary particles may be short-lived and may decay. A decay is treated like others collisions, the collision time being the time until the decay of the particle. All secondaries are tracked until they leave the nucleus, or the until the cascade stops.

The description of the target nucleus and fermi motion

The nucleus is constructed from $A$ nucleons and $Z$ protons with nucleon coordinates $\mathbf{r}_i$ and momenta $\mathbf{p}_i$, with $i = 1,2,...,A$. We use a common initialization Monte Carlo procedure, which is realized in the most of the high energy nuclear interaction models:

• Nucleon radii $r_i$ are selected randomly in the nucleus rest frame according to nucleon density $\rho(r_i)$. For heavy nuclei with $A > 16$ [2] nucleon density is
  

  $$\rho(r_i) = \frac{\rho_0}{1 + \exp{[(r_i - R)/a]}}$$ (25.1)

  
  
  where
  

  $$\rho_0 \approx \frac{3}{4\pi R^3}(1+\frac{a^2\pi^2}{R^2})^{-1}.$$ (25.2)

  
  
  Here $R=r_0 A^{1/3}$ fm and $r_0=1.16(1-1.16A^{-2/3})$ fm and $a \approx 0.545$ fm. For light nuclei with $A < 17$ nucleon density is given by a harmonic oscillator shell model [3], e. g.
  

  $$\rho(r_i) = (\pi R^2)^{-3/2}\exp{(-r_i^2/R^2)},$$ (25.3)

  
  
  where $R^2 = 2/3<r^2> = 0.8133 A^{2/3}$ fm$^2$. To take into account nucleon repulsive core it is assumed that internucleon distance $d > 0.8$ fm;

• The nucleus is assumed to be isotropic, i.e. we place each nucleon using a random direction and the previously determined radius $r_i$.

• The initial momenta of the nucleons $p_i$ are randomly choosen between $0$ and $p^{max}_F(r)$, where the maximal momenta of nucleons (in the local Thomas-Fermi approximation [4]) depends from the proton or neutron density $\rho$ according to
  

  $$p^{max}_F(r) = \hbar c(3\pi^2 \rho(r))^{1/3}$$ (25.4)

  
  

• To obtain momentum components, it is assumed that nucleons are distributed isotropic in momentum space; i.e. the momentum direction is chosen at random.

• The nucleus must be centered in momentum space around $\mathbf{0}$, i. e. the nucleus must be at rest, i. e. $\sum_i {\bf p_i} = \bf0$; To achieve this, we choose one nucleon to compensate the sum the remaining nucleon momenta $p_rest=\sum_{i=1}^{i=A-1}$. If this sum is larger than maximum momentum $p^{max}_F(r)$, we change the direction of the momentum of a few nucleons. If this does not lead to a possible momentum value, than we repeat the procedure with a different nucleon having a larger maximum momentum $p^{max}_F(r)$. In the rare case this fails as well, we choose new momenta for all nucleons.

  This procedure gives special for hydrogen $^1$H, where the proton has momentum $p=0$, and for deuterium $^2$H, where the momenta of proton and neutron are equal, and in opposite direction.

• We compute energy per nucleon $e = E/A = m_{N} + B(A,Z)/A$, where $m_N$ is nucleon mass and the nucleus binding energy $B(A,Z)$ is given by the tabulation of [5]: and find the effective mass of each nucleon $m^{eff}_i = \sqrt{(E/A)^2 - p^{2\prime}_i}$.

Optical and phenomenological potentials

The effect of collective nuclear elastic interaction upon primary and secondary particles is approximated by a nuclear potential.

For projectile protons and neutrons this scalar potential is given by the local Fermi momentum $p_{F}(r)$

$$V(r) = \frac{p_{F}^2(r)}{2 m}$$ (25.5)

where $m$ is the mass of the neutron $m_n$ or the mass of proton $m_p$.

For pions the potential is given by the lowest order optical potential [6]

$$V(r) = \frac{-2 \pi (\hbar c)^2 A }{ \overline{m}_\pi } ( 1 + \frac{m_\pi}{M} ) b_0 \rho(r)$$ (25.6)

where $A$ is the nuclear mass number, $m_\pi$, $M$ are the pion and nucleon mass, $\overline{m}_\pi$ is the reduced pion mass $\overline{m}_\pi = (m_\pi m_N) / (m_\pi + m_N)$, with $m_N$ is the mass of the nucleus, and $\rho(r)$ is the nucleon density distribution. The parameter $b_0$ is the effective $s-$wave scattering length and is obtained from analysis to pion atomic data to be about $-0.042 fm$.

Pauli blocking simulation

The cross sections used in this model are cross sections for free particles. In the nucleus these cross sections are reduced to effective cross sections by Pauli-blocking due to Fermi statistics.

For nucleons created by a collision, ie. an inelastic scattering or from decay, we check that all secondary nucleons occupy a state allowed by Fermi statistics. We assume that the nucleus in its ground state and all states below Fermi energy are occupied. All secondary nucleons therefore must have a momentum $p_i$ above local Fermi momentum $p_F(r)$, i.e.

$$p_i > p_F^{max}(r).$$ (25.7)

If any of the nucleons of the collision has a momentum below the local Fermi momentum, then the collision is Pauli blocked. The reaction products are discarded, and the original particles continue the cascade.

The scattering term

The basis of the description of the reactive part of the scattering amplitude are two particle binary collisions (hence binary cascade), resonance production, and decay. Based on the cross-section described later in this paper, collisions will occur when the transverse distance $d_t$ of any projectile target pair becomes smaller than the black disk radium corresponding to the total cross-section $\sigma_t$

$${\sigma_t\over \pi} > d_t^2$$

In case of a collision, all particles will be propagated to the estimated time of the collision, i.e. the time of closest approach, and the collision final state is produced.

Total inclusive cross-sections

Experimental data are used in the calculation of the total, inelastic and elastic cross-section wherever available.

hadron-nucleon scattering

For the case of proton-proton(pp) and proton-neutron(pn) collisions, as well as $\pi^=$ and $\pi^-$ nucleon collisions, experimental data are readily available as collected by the Particle Data Group (PDG) for both elastic and inelastic collisions. We use a tabulation based on a sub-set of these data for $\sqrt{S}$ below 3 GeV. For higher energies, parametrizations from the CERN-HERA collection are included.

Channel cross-sections

A large fraction of the cross-section in individual channels involving meson nucleon scattering can be modeled as resonance excitation in the s-channel. This kind of interactions show a resonance structure in the energy dependency of the cross-section, and can be modeled using the Breit-Wigner function

$$\sigma_{res}(\sqrt{s}) = \sum_{FS}{}{2J+1\over (2S_1+1)(2S_2+... ...i\over k^2}{\Gamma_IS\Gamma_FS\over (\sqrt{s}-M_R)^2+\Gamma/4},$$

Where $S1$ and $S2$ are the spins of the two fusing particles, $J$ is the spin of the resonance, $\sqrt(s)$ the energy in the center of mass system, $k$ the momentum of the fusing particles in the center of mass system, $\Gamma_IS$ and $\Gamma)FS$ the partial width of the resonance for the initial and final state respectively. $M_R$ is the nominal mass of the resonance.

The initial states included in the model are pion and kaon nucleon scattering. The product resonances taken into account are the Delta resonances with masses 1232, 1600, 1620, 1700, 1900, 1905, 1910, 1920, 1930, and 1950 MeV, the excited nucleons with masses of 1440, 1520, 1535, 1650, 1675, 1680, 1700, 1710, 1720, 1900, 1990, 2090, 2190, 2220, and 2250 MeV, the Lambda, and its excited states at 1520, 1600, 1670, 1690, 1800, 1810, 1820, 1830, 1890, 2100, and 2110 MeV, and the Sigma and its excited states at 1660, 1670, 1750, 1775, 1915, 1940, and 2030 MeV.

Mass dependent resonance width and partial width

During the cascading, the resonances produced are assigned reall masses, with values distributed according to the production cross-section described above. The concrete (rather than nominal) masses of these resonances may be small compared to the PDG value, and this implies that some channels may not be open for decay. In general it means, that the partial and total width will depend on the concrete mass of the resonance. We are using the UrQMD[13][14] approach for calculating these actual width,

$$\Gamma_{R\rightarrow 12}(M) = (1+r){\Gamma_{R\rightarrow 12... ...l+1)}}{M_R\over M} {p(M)^{(2l+1)}\over 1+r(p(M)/p(M_R))^{2l}}.$$ (25.8)

Here $M_R$ is the nominal mass of the resonance, $M$ the actual mass, $p$ is the momentum in the center of mass system of the particles, $L$ the angular momentum of the final state, and r=0.2.

Resonance production cross-section in the t-channel

In resonance production in the t-channel, single and double resonance excitation in nucleon-nucleon collisions are taken into account. The resonance production cross-sections are as much as possible based on parametrizations of experimental data[15] for proton proton scattering. The basic formula used is motivated from the form of the exclusive production cross-section of the $\Delta_{1232}$ in proton proton collisions:

$$\sigma_{AB} = 2\alpha_{AB}\beta_{AB}{\sqrt{s}-\sqrt{s_0} \ove... ...t ( {\sqrt{s0}+\beta_{AB}\over \sqrt{s}}\right )^{\gamma_{AB}}$$

The parameters of the description for the various channels are given in table25.1. For all other channels, the parametrizations were derived from these by adjusting the threshold behavior.

Table 25.1: Values of the parameters of the cross-section formula for the individual channels.

Reaction	$\alpha$	$\beta$	$\gamma$
$\rm pp\rightarrow p\Delta_{1232}$	25 mbarn	0.4 GeV	3
$\rm pp\rightarrow \Delta_{1232}\Delta_{1232}$	1.5 mbarn	1 GeV	1
$\rm pp\rightarrow pp^{*}$	0.55 mbarn	1 GeV	1
$\rm pp\rightarrow p\Delta_{*}$	0.4 mbarn	1 GeV	1
$\rm pp\rightarrow \Delta_{1232}\Delta^{*}$	0.35 mbarn	1 GeV	1
$\rm pp\rightarrow \Delta_{1232}N^{*}$	0.55 mbarn	1 GeV	1

The reminder of the cross-section are derived from these, applying detailed balance. Iso-spin invariance is assumed. The formalism used to apply detailed balance is

$$\sigma(cd\rightarrow ab) = \sum_{J,M}{}{\left < j_cm_cj_dm_... ...\right> \over \left< p_{cd}^2\right> }\sigma(ab\rightarrow cd)$$ (25.9)

Nucleon Nucleon elastic collisions

Angular distributions for elastic scattering of nucleons are taken as closely as possible from experimental data, i.e. from the result of phase-shift analysis. They are derived from differential cross sections obtained from the SAID database, R. Arndt, 1998.

Final states are derived by sampling from tables of the cumulative distribution function of the centre-of-mass scattering angle, tabulated for a discrete set of lab kinetic energies from 10 MeV to 1200 MeV. The CDF's are tabulated at 1 degree intervals and sampling is done using bi-linear interpolation in energy and CDF values. Coulomb effects are taken into consideration for pp scattering.

Generation of transverse momentum

Angular distributions for final states other than nucleon elastic scattering are calculated analytically, derived from the collision term of the in-medium relativistic Boltzmann-Uehling-Uhlenbeck equation, absed on the nucleon nucleon elastic scattering cross-sections:

$$\sigma_{NN\rightarrow NN}(s,t) = {1\over(2\pi)^2s}\left(D(s,t)+E(s,t) + (inverted t,u)\right)$$

Here $s$, $t$, $u$ are the Mandelstamm variables, $D(s,t)$ is the direct term, and $E(s,t)$ is the exchange term, with

$$D(s,t) = {(g_{NN}^\sigma)^4(t-4m^{*}2)^2 \over2(t-m_\sigma^2 )^2} + {(g_{NN}^\omega)^4(2s^2+2st+t^2-8m^{*2}s+8m^{*4}) \over (t-m_\omega^2)^2} +$$

$${24(g_{NN}^\pi)^4m^{*2}t^2\over (t-m_\pi^2)^2} - {4(g_{NN}^\sigma g_{NN}^\omega)^2(2s+t-4m^{*2})m^{*2}\over (t-m_\sigma^2)(t-m_\omega^2)},$$

and

$$E(s,t) = {(g_{NN}^\sigma)^4\left( t(t+s)+4m^{*2}(s-t)\right)\over 8(t-m_\sigm... ...(g_{NN}^\omega)^4(s-2m^{*2})(s-6m^{*2}))\over 2(t-m_\omega^2)(u-m_\omega^2) } -$$

$${6(g_{NN}^\pi)^4(4m^{*2}-s-t)m^{*4}t\over (t-m_\pi^2)(u=m_pi^2) }+ {... ...{NN}^\pi)^2 m^{*2} (4m^{*2}-s-t)(4m^{*2}-t) \over (t-m_\sigma^2)(u-m_\pi^2) } +$$

$${3(g_{NN}^\sigma g_{NN}^\pi)^2 t(t+s)m^{*2}\over 2(t-m_\pi^2)(u-m_\s... ...\omega)^2 t^2-4m^{*2}s-10m^{*2}t+24m^{*4}\over4(t-m_\sigma^2)(u-m_\omega^2) } +$$

$${(g_{NN}^\sigma g_{NN}^\omega)^2 (t+s)^2-2m^{*2}s+2m^{*2}t\over 4(t-... ...mega g_{NN}^\pi)^2 (t+s-4m^{*2})(t+s-2m^{*2})\over (t-m_\omega^2)(u-m_\pi^2)} +$$

$${3(g_{NN}^\omega g_{NN}^\pi)^2 m^{*2} (t^2-2m^{*2}t)\over (t-m_\pi^2)(u-m_\omega^2)}.$$ (25.10)

Here, in this first release, the in-medium mass was set to the free mass, and the nucleon nucleon coupling constants used were 1.434 for the $\pi$, 7.54 for the $\omega$, and 6.9 for the $\sigma$. This formula was used for elementary hadron-nucleon differential cross-sections by scaling teh center of mass energy squared accordingly.

Finite size effects were taken into account at the meson nucleon vertex, using a phenomenological form factor (cut-off) at each vertex.

Decay

In the simulation of decay of strong resonances, we use the nominal decay branching ratios from the particle data book. The stochastic mass of a individual resonance created is sampled at creation time from the Breit-Wigner form, under the mass constraints posed by center of mass energy of the scattering, and the mass in the lightest decay channel. The decay width from the particle data book are then adjusted according to equation 25.8, to take the stochastic mass value into account.

All decay channels with nominal branching ratios greater than 1% are simulated.

The escaping particle and coherent effects

When a nucleon other than the incident particle leaves the nucleus, the ground state of the nucleus changes. The energy of the outgoing particle cannot be such that the total mass of the new nucleus would be below its ground state mass. To avoid this, we reduce the energy of an outgoing nucleons by the mass-difference of old and new nucleus.

Furthermore, the momentum of the final exited nucleus derived from energy momentum balance may be such that its mass is below its ground state mass. In this case, we arbitrarily scale the momenta of all outgoing particles by a factor derived from the mass of the nucleus and the mass of the system of outgoing particles.

Light ion reactions

In simulating light ion reactions, the initial state of the cascade is prepared in the form of two nuclei, as described in the above section on the nuclear model.

The lighter of the collision partners is selected to be the projectile. The nucleons in the projectile are then entered, with position and momenta, into the initial state of the cascade. Note that before the first scattering of an individual nucleon, a projectile nucleon's Fermi-momentum is not taken into account in the tracking inside the target nucleus. The nucleon distribution inside the projectile nucleus is taken to be a representative distribution of its nucleons in configuration space, rather than an initial state in the sense of QMD. The Fermi momentum and the local field are taken into account in the calculation of the collision probabilities and final states of the binary collisions.

Transition to pre-compound modeling

Eventually, the cascade assumptions will break down at low energies, and the state of affairs has to be treated by means of evaporation and pre-equilibrium decay. This transition is not at present studied in depth, and an interesting approach which uses the tracking time, as in the Liege cascade code, remains to be studied in our context.

For this first release, the following algorithm is used to determine when cascading is stopped, and pre-equilibrium decay is called: As long as there are still particles above the kinetic energy threshold (75 MeV), cascading will continue. Otherwise, when the mean kinetic energy of the participants has dropped below a second threshold (15 MeV), the cascading is stopped.

The residual participants, and the nucleus in its current state are then used to define the initial state, i.e. excitation energy, number of excitons, number of holes, and momentum of the exciton system, for pre-equilibrium decay.

In the case of light ion reactions, the projectile excitation is determined from the binary collision participants ($P$) using the statistical approach towards excitation energy calculation in an adiabatic abrasion process, as described in [12]:

$$E_{ex} = \sum_{P} (E_{fermi}^P-E^P)$$

Given this excitation energy, the projectile fragment is then treated by the evaporation models described previously.

Calculation of excitation energies and residuals

At the end of the cascade, we form a fragment for further treatment in precompound and nuclear de-excitation models ([16]).

These models need information about the nuclear fragment created by the cascade. The fragment is characterized by the number of nucleons in the fragment, the charge of the fragment, the number of holes, the number of all excitons, and the number of charged excitons, and the four momentum of the fragment.

The number of holes is given by the difference of the number of nucleons in the original nucleus and the number of non-excited nucleons left in the fragment. An exciton is a nucleon captured in the fragment at the end of the cascade.

The momentum of the fragment calculated by the difference between the momentum of the primary and the outgoing secondary particles must be split in two components. The first is the momentum acquired by coherent elastic effects, and the second is the momentum of the excitons in the nucleus rest frame. Only the later part is passed to the de-excitation models. Secondaries arising from de-excitation models, including the final nucleus, are transformed back the frame of the moving fragment.