Monte Carlo simulation for statistical Hauser-Feshbach theory

Monte Carlo simulations for particle andγ-ray emissions from a compound nucleus based on the Hauser-Feshbach statistical theory are performed. The Monte Carlo method is applied to the neutron induced nuclear reactions on 56 Fe, and the results are compared with a traditional deterministic method. The neutron and γ-ray emission correlation is examined by gating on an 847 keVγ-ray that is produced by an inelastic scattering process. The partial γ-ray energy spectra for differentγ-ray multiplicities is inferred using this Monte Carlo method.


Introduction
Nuclear reaction model codes like GNASH [1], have rather complicated algorithms to generate reaction-wise (or exclusive) energy spectra for individual particle emission process. In the GNASH case, all the decay probabilities at each intermediate stage are recorded, and the data are processed with a supplemental FORTRAN program, RECOIL [2], to calculate the exclusive spectra. A straightforward and more accurate approach to understand all the particle and γ-ray emission probabilities in the nuclear reaction, is to utilize a Monte Carlo technique. Instead of performing the complicated integration over all the probabilities at each intermediate stage in the course of a compound nucleus decay, we simulate the particle emission process by the Monte Carlo method, and record all the particles and γrays with their energies until the residual nucleus reaches its ground state.
The Monte Carlo technique also facilitates new insights into the nuclear reaction, as correlations between different processes can be studied by gating on a particular γ-ray so as to provide useful information on spin and parity distributions in the continuum state [3]. Coincidence experiments, i.e. a scattered neutron detection gated by a particular γ-ray transition, will provide an important benchmark of the modeling of nuclear properties. This kind of investigation, in both theory and experiments, will be more important in the near future for acquiring reliable nuclear data. The Monte Carlo method may help answer questions regarding nuclear reaction mechanisms that have been difficult to address by traditional techniques.
Based on this Monte Carlo study for the Hauser-Feshbach model [4], we can also envision a more rigorous simulation of neutron and γ-ray emissions from fission fragments. In the past, a Monte Carlo simulation for the prompt fission neutron emission has been performed [5][6][7][8], assuming that the neutron emission spectra from the fission fragments follow a Weisskopf probability distribution. Our present a e-mail: kawano@lanl.gov   To describe a formula to calculate the exclusive spectrum, we first need to define a decay probability from a given compound state (c n , k n ), where c n denotes the compound nucleus, k n is the index of excited state with the energy bin width of ∆E, and the suffix n describes the n-th compound nucleus. The index k runs from the highest excitation toward the ground state. We define a probability P(c n , k n , c m , k m ), which is a decay probability of the state (c n , k n ) by emitting a particle or a γ-ray to form a (c m , k m ) state. The probability P can be calculated with the Hauser-Feshbach theory, where T (c m k m → c n k n ) is the particle or γ-ray transmission coefficients from the final (c m , k m ) to initial (c n , k n ) states. Note that all the quantum numbers of T are already summed. An algorithm to calculate the exclusive spectrum may depend on how much complexity we introduce in evaluating the decay path in Fig 1. Increasing the γ-ray multiplicities results in more summations over the number of intermediate stages. For example, the path (a) has the multiplicity of one, while it becomes 3 for the cases of (b) and (c). Because the particle emission probabilities are much larger than the γ-ray emission probability, it is often good enough to consider the path (a) only. However, we allow the γ-ray competition at each compound nucleus in this study.

Monte Carlo technique
Once the probabilities P(c n , k n , c m , k m ) are constructed with the Hauser-Feshbach formula, a Monte Carlo simulation for the particles and γ-ray decay processes is straightforward. Starting with the initial state (c 0 , k 0 ), a following state can be determined by throwing a dice based on the probability P(c 0 , k 0 , c m , k m ), where we know that We repeat this random sampling procedure until the nucleus reaches a stable state, and record all the particles and γ-rays emitted during the sequence. In this procedure, since the particle and γ-ray multiplicities are unlimited, we can reconstruct very accurately the exclusive particle emission spectra by tallying the recorded Monte Carlo histories. The theoretical model calculations are always in the center-of-mass system (CMS), despite the fact that particle transport calculations require all the quantities in the laboratory system (LAB). The Monte Carlo technique allows us to calculate recoils of residual nuclei in the laboratory frame. Although we do not discuss the system conversion between LAB and CMS here, we would like to note that quantities obtained with the Monte Carlo simulation are directly connected to what we see in a real experiment. Hereafter all the particle and γ-ray energies are in the CMS.

Hauser-Feshbach model calculation
The CoH Hauser-Feshbach code [9] has been upgraded recently. The version 3.0 code, CoH 3 , calculates all the nuclear reaction cross sections above the resonance range, including the total, shape elastic scattering, direct inelastic scattering, direct/semidirect capture cross sections, preequilibrium emission, and particle and γ-ray emission in the compound decay process. All the particle transmission coefficients are calculated internally, so that it does not require an external solver of optical model potentials. The code is adopted to generate the probabilities in Eq. (1).

Comparison of deterministic and Monte Carlo methods
As an example, we consider a 20-MeV neutron induced reaction on 56 Fe. To obtain a satisfactory statistics in the Monte Carlo simulation, we tallied 10 8 histories of particle and γ-ray decays for the 20-MeV neutron induced reactions on 56 Fe case. Figure 2 shows a comparison of the calculated exclusive spectra by both deterministic and Monte Carlo methods for (n,n'), (n,2n), (n,np) and (n,nα) reactions. Each spectrum is normalized per incident neutron per outgoing energy. Therefore (n,n') (n,2n) (n,np) (n,nα) Fig. 3. Exclusive γ-ray spectra for the (n,n'), (n,2n), (n,np) and (n,nα) reactions calculated with the Monte Carlo method. The spectra are normalized per outgoing energy per incident neutron.
hold, where σ n ′ and σ 2n are the inelastic and (n,2n) reaction cross sections, and σ R is the total reaction cross section calculated with the optical model. These two reactions are clearly separated at the neutron separation energy, and the two methods give almost identical results, except for tiny differences at low energies of (n.n') spectrum. The agreement for the (n,np) and (n,nα) reaction cases is also excellent. It is not so surprising that these two methods agree very well, as we started with the same probability distributions, together with the fact that the γ-ray transition probabilities are usually much smaller than the particle emission probabilities.

γ-ray emission spectra and multiplicities
The γ-ray emission spectra and the γ-ray multiplicities corresponding to a specific reaction can be obtained very easily with the Monte Carlo technique, by summing all the γ-rays emitted and dividing by the number of simulations performed, which is shown in Fig. 3. The calculated γ-ray multiplicities are 2.57 for the (n,n') and 1.67 for (n,2n) reactions.

Particle and γ-ray correlation
One of the notable benefits of the Monte Carlo method is that we can sample any kind of reactions by looking into the history record. As an example, we calculate 20-MeV neutron induced reactions on 56 Fe again, and observe an inelastic scattering process that produces an 847 keV γ-ray, which is a transition from the first 2 + state to the 0 + ground state in 56 Fe. In this simulation, we are gating on the inelastic scattering process by the partial 847 keV γ-ray, so  that only one neutron is allowed to be emitted, while multiple γ-rays can be emitted in one simulation. The energy spectra of inelastically scattered neutrons that produce the 847 keV γ-ray are shown in Fig. 4, for various m values. The γ-ray spectra, each corresponding to the neutron spectrum case in Fig 4, are shown in Fig. 5. The calculated average multiplicity is m=3.04. This implies that a two-step γ-ray cascading plus E2 transition from the 2 + state after neutron emission is the most probable process. Assuming E1 is the dominant γ-ray transition, which changes a parity at each decay step, the initial compound nucleus state has an even parity, and the intermediate state (after the first γ-ray emission) has a negative parity. This γ-ray cascading might be sensitive to a parity distribution of the compound nucleus.