Two-stage description of charge distributions in 56Fe+p spallation reactions at 0.3-1.5 GeV/A

Spallation reactions are studied in the framework of a two-stage model. The intra-nuclear cascade (INC) stage of the reaction is simulated with the code ISABEL and a Constrained Molecular Dynamics (CoMD) code. The de-excitation of the highly excited pre-fragments is described with the multi-sequential binary decay code MECO, based on a generalized Weisskopf-Ewing formalism. Emission of nucleons, gamma rays and IMFs in their ground, excited bound and unbound states is considered. Calculated cross sections are compared with the experimental charge distributions of Fe + p spallation reaction products studied at GSI with the fragment separator FRS in the energy-range 0.3 − 1.5 GeV/A. At all bombarding energies, a good description of the experimental data is obtained with ISABEL-MECO using a global set of parameters. A preliminary calculation with CoMD-MECO at 1GeV/A is discussed.


Introduction
In recent years, proton-induced spallation reactions have received a great deal of attention in basic and applied Nuclear Physics research. In basic research, they provide a framework for the development and testing of high-energy nuclear reaction models [1]. In Applied Nuclear Research they provide useful information for the development of accelerator-driven systems (ADS), the transmutation of nuclear waste, spallation neutron sources and the production of radioactive beams. In the case of targets with mass A~50-60, spallation reactions are used for the production of isotopes of medical interest, studies concerning the interaction of cosmicrays with interstellar bodies (radiation damage) and radiation protection in space [2].
A proton-induced spallation reaction (SR) proceeds in two stages [3]. In the first stage, the incident particle induces a sequence of collisions transferring energy and momentum to the target nucleons. As a result, we have the formation of an intranuclear cascade (INC) of highenergy (greater than 20 MeV) protons, neutrons and * Corresponding author, email: nnicolis@uoi.gr pions within the nucleus. This is a fast process lasting approximately 10 . During this process some of the energetic hadrons may escape from the target. Others deposit their kinetic energy in the nucleus leaving it in an excited state.
In the second stage, the produced hot pre-fragments equilibrate and deexcite. Sequential evaporation is assumed to occur with a typical time scale of 10 . Emission of nucleons, protons, alpha-particles and gamma-rays is dominant. Emission of heavier nucleon clusters in their ground or excited states is also possible. If the target is heavy enough, high-energy fission may compete with sequential evaporation. The deexcitation products of target-like and/or fission fragments may be radioactive.
Recently, detailed experimental and theoretical studies have been devoted on the spallation of with protons. In a high precision experiment, spallation residue cross sections were measured in + bombardments at 0.3 -1.5 GeV/A. The experiment employed reverse-kinematics with a liquid hydrogen target and the fragment separator FRS at GSI (Darmstadt) [4]. High resolution measurements of fully identified fragments and coincidence measurements of residues and light particles at 1.0 GeV/A have also been reported [5,6]. Efforts are made for the understanding of the odd-even structure observed in the charge distributions [7,8] and the possibility of a multifragmentation decay mechanism at the highest bombarding energies [6,9].
In the present work, we employ two different models for the INC stage coupled with the same evaporation code for the description of the experimental charge distributions of + spallation reaction data at 0.3-1.5 GeV/A [4,5]. In the next Section we describe our reaction models. There follows a comparison of our calculations with the experimental data, a brief discussion, and a summary.

Reaction models
We describe the INC stage with the Monte-Carlo code ISABEL [10,11]. ISABEL describes the target nucleus as a continuous medium with a diffuse surface. Nucleon-nucleon collisions may be elastic or inelastic. Linear trajectories are assumed between collisions, which occur with a criterion based on the mean-freepath. Free N-N cross sections are assumed. A full Pauli blocking mechanism forbids nucleons from falling below the Fermi surface.
We have also employed the Constrained Molecular Dynamics (CoMD) model described in Refs. [12][13][14][15]. CoMD is a semi-classical model originally designed for reactions near and below the Fermi energy. Nucleons are considered as Gaussian wave packets interacting with a phenomenological N-N interaction (Skyrme). Pauli principle is taken into account through an appropriate restriction in phase-space. Angular momentum conservation is obeyed. Recognition of fragment formation is made at each step of the time evolution.
The deexcitation of the pre-fragments is calculated in the context of equilibrium compound nucleus decay. Possible pre-equilibrium emissions between the two stages are not considered. For the equilibrium decay process, we employ the Multi-sequential Evaporation COde MECO [16]. MECO is a Monte-Carlo code sufficient to describe low-spin, high-excitation energy reaction systems. It describes the equilibrium decay of excited nuclei as a sequence of binary processes involving emission of fragments in their ground, excited bound and unbound states. Any number of user-defined channels may be considered. Emission of nucleons and heavy fragments up to symmetric mass divisions is described in the framework of a generalized Weisskopf-Ewing evaporation formalism. The applicability of this formalism is limited to systems with an effective fissility below the Businaro-Gallone point [17], where the saddle point becomes unstable toward asymmetric divisions.
In the statistical decay calculations, each fragment emitted in an excited-bound state is considered as a separate decay channel. We performed calculations with a total of 180 decay channels consisting of 42 ground states of fragments with Z=0 (neutrons) up to Z=12, 100 excited states for fragments with ≥ 5. The heaviest 28 fragments were allowed to decay in the continuum. We employ phenomenological level densities described with the composite Gilbert and Cameron formula [18]. Level-density parameters were considered as energyindependent according to Gilbert and Cameron. Adjustments were made in the energy shifts of the Fermi-gas formula, thus improving the reproduction of the charge-distribution yields. Inverse cross sections were calculated with the optical model for Z<4. For heavier fragments with ≥ 5, barrier penetration calculations were performed with the Christensen and Winther nuclear potential [19].
At each bombarding energy of the present study, 20000 events were generated from the INC models following a triangular distribution in the impact parameter range = 0 − 6.0 . These events were fed into MECO and their sequential decay was calculated. There followed a second stage of equilibrium decay calculations for the unbound binary decay fragments. Events from both equilibrium decay calculations were sorted to produce the isotopic distributions of spallation residues to be compared with the data.  Before applying our models to the high-energy data, we compared our procedure to the predictions of the reaction code TALYS [20], known to provide an accurate description of low-energy proton-induced reactions at proton beam energies up to 200MeV. A comparison of the mass and charge distributions predicted by TALYS was made with calculations based on ISABEL/MECO and CoMD/MECO for protons incident on with an energy of 50, 100 and 150MeV. Consistency between the three procedures was found for the cumulative A-and Z-distributions.

Comparison with experimental data
Symbols in Figure 2 show the experimental charge distributions at 0.3, 0.5, 0.75, 1.0 and 1.5 GeV/A from Ref. [4]. At 1.0 GeV/A, the IMF cross sections of Ref. [5] are also shown. The solid curves show the results of the ISABEL-MECO calculations. We obtain a good overall description of the Zdistributions except for the region with ≤ 13 at the 0.3 and 1.5 GeV/A. At all energies, the odd-even structure observed in the experimental data is reproduced. For the A-distributions we found the same kind of agreement between experimental and calculated cross sections. These calculations reproduce the centroids of the experimental isotopic distributions. The widths of the isotopic distributions are also reproduced for Ζ>20 and overestimated the lower Z's by less than one mass unit.
The symbols in Figure 3 show the experimental Zdistribution of spallation residues and IMF's at 1.0 GeV/A. The distribution marked with "1" represents the excited pre-fragments predicted by ISABEL. This distribution is structureless. If we allow these events to deexcite by neutron, proton, alpha and gamma emission, we end up with the distribution marked with "2". As Z decreases, the calculation underestimates the experimental cross sections by a lot. If we consider emission of IMFs in their ground states, we obtain the curve marked with "3". This calculation with MECO considers 42 fragments in their ground states with Z up to 12 and accounts for most cross sections. Including emission from excited-bound states and allowing for binary decays in the continuum brings the calculation with close agreement with the experiment (curve marked with "4"). From this Figure, we realize that the odd-even structure of cross sections appears because of the equilibrium decay. Evaporation of just n, p and alpha is sufficient to produce such an effect.  [4] and emitted IMFs [5] in + reactions at 1.0 GeV/A. Solid curves show the calculated distribution with ISABEL (curve 1) and ISABEL-MECO (curves 2, 3 and 4) using various assumptions for the number of equilibrium decay channels described in the text.
We have determined the contributions to the A-and Z-distributions from fragments emitted in their ground or excited states. The thick solid line on Figure 4 shows the total Z-distribution at 1.0 GeV/A calculated with ISABEL-MECO using 180 channels. Experimental data are shown with symbols. The shaded marked "1" shows the contribution from fragments emitted in their ground states. Region marked "2" shows the population from fragments emitted in excited-bound states. The twobumps marked "3" and "4" correspond to binary decays. Experimental data [4,5] are shown with symbols.
A preliminary calculation with CoMD for the Zdistribution at 1.0 GeV/A was performed. This is a stringent test for this model since the bombarding energy is well above its range of validity. It was realized that the excitation energy distribution of the prefragments, does not extend to energy as high as the ISABEL calculation does (cf. with Figure 1). Applying a correction factor of 1.3 to CoMD produces an excitation energy distribution close to ISABEL. The resulting calculation with CoMD-MECO is shown in Figure 5. We realize that there is a tendency to underestimate the distribution as Z decreases. Furthermore, there is excess of events predicted in the primary stage with zero excitation energy. The highenergy behaviour of the CoMD code is currently under investigation.  [4] and emitted IMFs [5] in + reactions at 1.0 GeV/A. The solid curve shows the calculated distribution with CoMD-MECO described in the text.

Discussion and Summary
We have obtained a good description of spallation residue cross sections with the combination of the intranuclear cascade model ISABEL coupled with the sequential binary decay code MECO. Assuming that treatment of the primary stage with ISABEL is realistic and neglecting the effect of possible preequilibrium emissions between the two stages of the reaction (not accounted for by ISABEL, if any), we may conclude that the treatment of the equilibrium decay stage is successful. The employed statistical model parameters are global and energy independent. In Figure 2 we see some underestimations in the mid-Z range of the charge distributions at 0.3 and 1.5 GeV/A. They may be related to specific nuclear parameters not accounted for by our average statistical model parameter set. In addition, the underestimation at 1.5 GeV/A may be related to the onset of multifragmentation. Inspection of the ISABEL events shows a small fraction of pre-fragments produced at 1.0 and 1.5 GeV/A with temperatures 3-4 MeV. Calculations with the statistical multifragmentation model SMM [21,22] are under way, keeping in mind that such an interpretation may not be straightforward [6,9].
Summarizing, spallation products of + reactions at 0.3-1.5 GeV/A were described as a twostage process with the ISABEL INC followed by the MECO sequential binary decay code. A good agreement with the experimental Z-distributions was obtained with an effective set of excitation-energy independent parameters for the equilibrium de-excitation. Contributions from multifragmentation may be considered at the highest bombarding energies (1.0-1.5GeV) in calculations based on more detailed statistical decay parameters and a proper consideration of multifragmentation decays. Improvements in the calculation of the INC stage with the CoMD model are under way.