Investigation of the Energy Correlations of Spallation Neutrons by Monte Carlo Simulations

Earlier works have suggested that the energy correlations in a spallation
source may influence the neutron noise measurements in an ADS. For the
calculation of this effect not only the generally known and used one-particle
spectrum is needed but also the so-called two particle spectrum, which describes
also the energy correlations. Since measured data are not available
for the energy correlation of the neutrons from a single spallation event, the
physical models of the MCNPX code have been used to investigate the effect.
The calculational model has been successfully validated with measurements
of the number distribution of spallation neutrons. The simulated one- and
two-particle energy distributions and spectra proved that the energy correlations
exist and have an important effect in low multiplicity spallation events
and in thin targets. On the other hand, for thick targets this effect appears
negligible and the factorization of the two-particle spectrum seems an acceptable
approximation. Further investigations are in hand to quantify the actual effect of the energy correlations on the neutron noise measurements.


Introduction
Spallation sources are considered as neutron sources for Accelerator-Driven Subcritical Systems (ADSs).In an ADS an accelerator is supposed to provide a high-energy (∼1 GeV) proton beam, which produces neutrons in a heavy metal (e.g., lead) target through the spallation process.Besides the high neutron energies (up to the energy of the proton) the spallation neutron source is also distinguished from other neutron sources by its high multiplicity: one proton can produce up to 40-50 neutrons.The high multiplicity increases the importance of the source in a subcritical system, especially in case of deeper subcriticality.Therefore, the precise description of the spallation source is inevitable for the modeling of an ADS.
A special case is the modeling of the neutron fluctuations in an ADS (e.g., neutron noise measurements for reactivity determination) as this requires also the description of the higher moments of the probability distributions.While one can easily find measured data about spallation sources for the average values (the first moment) in the literature, the higher moments are not available due to the smaller interest and the difficulties of such measurements.This paper makes an attempt to reproduce the higher moments and correlations needed for the accurate simulation of the neutron fluctuations in an ADS with the help of the physical models of the spallation process implemented in the MCNPX [1] highenergy Monte Carlo particle transport code.

Effect of the Spallation Source on the Neutron Fluctuations
2.1.Effect of the Multiplicity of the Source.The distribution of the number of neutrons q from a single source event can be described by the probability distribution p(q), which by definition has to fulfill the following criteria: Yamane and Pázsit et al. investigated the effect of the high multiplicity of the spallation source on different neutron noise measurement methods [4][5][6].Let us take the varianceto-mean ratio or Feynman-α measurement as an example.
According to the well-known formula, if a Poisson neutron source (e.g., Pu-Be) is placed into a subcritical core, the ratio of the variance and the mean value of the number of counts N in a detector during a certain measuring time Δt will be the following [6]: where is the detector efficiency, ρ is the reactivity of the subcritical system, β is the delayed neutron fraction, α is the prompt decay constant, and D νp is the so-called Divenfactor for the prompt neutrons, which can be calculated from the first and second factorial moments of the number distribution of prompt neutrons from fission (ν p ): As it is shown in [4,5] in case of a multiplicative source a (1 + δ) factor has to be introduced in (2): δ is defined as: where D q is the Diven factor for the number of the neutrons from a source event q and can be calculated analogously to D νp : A similar formula with the same (1 + δ) correction factor can be derived for the Rossi-α or autocorrelation measurement, too [6].From (5) one can arrive at the trivial conclusion that the higher the multiplicity of the source compared to the multiplicity from fission and the deeper the subcriticality, the more important the effect of the source on the neutron fluctuations.
The importance of this fact comes if one calculates the typical values of δ for a spallation source.The Diven-factor for spallation sources is somewhat higher than for fission.For example, for the number distribution in [3] measured for 35 cm thick Pb target (see in Section 5.1 in details) D q ≈ 1.4 can be obtained.Taking typical values for 235 U fission as in [4] one arrives at D q /D νp ≈ 1.75.The more important effect comes from the source multiplicity since in the same measurement as above q ≈ 20 was found for a spallation source, which gives q / ν p ≈ 8.3.Assuming these data δ ≈ −14.5 ρ, which means δ ≈ 0.76 in case of a subcritical system with k eff = 0.95.If one considers a deep subcriticality with k eff = 0.7 then δ ≈ 6.2 is obtained.These high values for δ emphasize the importance of the spallation source in the neutron fluctuations.The importance of the source multiplicity has been confirmed by experimental results from Kitamura et al. [2,6].In Figure 1 one can see the results of Feynman-α measurements in a subcritical system of k eff = 0.9874 with different neutron sources.The curves obtained with the Am-Be source (which is a Poisson source) and with the 252 Cf source cover each other due to the fact that the multiplicity of the 252 Cf source is about the same as the fission multiplicity, which results in a small δ value.On the other hand, when a randomly triggered pulsed D-T neutron generator (Random-PNG) is applied, despite the relatively small negative reactivity, the effect of the higher δ value can be observed due to the high number of neutrons from a pulse.

Effect of the Energy
Correlations of the Source.The aboveshown importance of the spallation source indicates that minor effects may be worthwhile for deeper investigation, as well, such is the energy correlation between the neutrons from a spallation source.Neutron fluctuations are often handled in an energy independent one-group theory, but the high-energy spectrum assumed in an ADS seeks for energy dependent multigroup treatment.Such treatment requires also the description of the energy distribution of the source neutrons.This can be described in full detail by giving the f q (E 1 , E 2 , . . ., E q ) energy distribution of the q-particle emission event for each q.Being conditional probability density functions, these functions also have to be normalized to unity: One can also define n-particle distributions (where n < q) by averaging over the remaining energy coordinates: In practice the f q (E) one-and the f q (E 1 , E 2 ) two-particle distributions have relevance, since they are sufficient for the determination of the first and second moments.If the energies of the neutrons originating from a single source event are totally independent from each other, the q-particle distribution can be factorized and substituted by the oneparticle distribution: which obviously also means that: In some cases the energy independence can be a good approximation but generally it cannot be assumed.The spectrum of the neutrons from the source χ(E) can be calculated by averaging the one-particle distribution f q (E) weighted by the expected number of neutrons from a qparticle emission event qp(q): This spectrum is known from experiments for both fission and spallation sources.However, Pázsit et al. showed [7] that in an energy-dependent approach the calculation of the second moment of the neutron fluctuations requires also the two-particle spectrum χ(E 1 , E 2 ), which involves the energy correlations between the source neutrons.This can be expressed with the help of the two-particle energy distributions: This function is usually unknown.In the case of fission it is generally replaced by the one-particle spectrum assuming that χ(E 1 , E 2 ) = χ(E 1 )χ(E 2 ).However, from (12) one can conclude that this assumes not only that ( 10) is valid but also the independence of f q (E) from q: This approximation is often used (e.g., in [8] where a general theory of neutron noise measurements is derived assuming different neutron sources), although the conditions in (10) and ( 13) do not apply for a spallation source.In order to decide on the applicability of this approximation the determination of the two-particle spectrum χ(E 1 , E 2 ) is needed.As measurement data are not available for this purpose, one possibility is to use the physical models describing spallation and calculate a Monte Carlo estimation of χ(E 1 , E 2 ).

Modeling of Spallation
The spallation process is modeled in high-energy particle transport codes (e.g., MCNPX, LAHET [9], FLUKA [10], and HERMES [11]).Such codes use Monte Carlo methods both in the physical models to determine the outcome of nuclear interactions and for the transport of particles between interactions.The transport between interactions is done in the conventional way except that total cross-sections are calculated from physical models.Whenever a collision site is sampled, physical models are used to calculate the outcome of the interaction.The physical models are composed of two basic parts.The first part is the so-called intranuclear cascade (INC) model.This describes the nuclear interaction as series of nucleon-nucleon collisions inside the nucleus.During this process secondary particles are generated as "knocked-out" nucleons.The process is followed while a thermal equilibrium of the nucleons is achieved, and the nucleus is left in an excited state.The second part of the model describes the deexcitation of the nucleus, which leads to either evaporation (simultaneous emission of several particles) or fission.
It is important to note that the modeling of the spallation process is fully analogous.Therefore, one can expect that these models provide good estimation not only for the mean values of the distributions but also for the higher moments, which are important for the simulation of neutron fluctuations.

Calculational Methods
The MCNPX code has been chosen for the calculation.Geometry and parameters of a published measurement of spallation neutron yields have been used in the simulation in order to validate the physical models.In the experiments performed by Hilscher et al. [3] a lead target of variable thickness was bombarded by 1.22 GeV protons.The target was surrounded by a 4π neutron detector, the so-called Berlin Neutron Ball (BNB) (see Figure 2).This is a spherical shell filled with liquid scintillator and 0.4 weight% of Gd in order to make it sensitive to neutrons.24 photomultiplier tubes are attached to the outer surface of the shell to detect the scintillation events.
In the simulation only the target has been modeled as a natural lead (ρ = 11.34 g/cm 3 ) cylinder with a diameter of 15 cm and a height of 0.2 cm, 5 cm, and 35 cm, respectively.The proton beam arrives axially at the middle of the cover plane of the target.
The MCNPX code includes the Bertini [12], the Isabel [13], the INCL4 [14] and the CEM03 [15] INC models, and the ABLA and Dresner [16] deexcitation models.CEM03 consists of an intranuclear cascade model, followed by a preequilibrium model and an evaporation model.For Bertini and ISABEL and the INCL4 the Dresner evaporation model with Rutherford Appleton Laboratory (RAL) fission has been used.The transport of protons, neutrons, and pions was followed in the simulations.
During the simulations data have been collected about neutrons leaving the target with the help of the PTRAC event file of the MCNPX, in which (upon user request) data are recorded about certain events.A program has been developed to process this PTRAC file, extract the energy of the neutrons escaping the target, and reconstruct the requested distributions.It has to be noted again that in this way fully analogous estimators have been used, which do not bias the higher moments of the distributions.The number of simulated source events has been chosen large enough to arrive at acceptable statistics with a few percent of relative error.

Validation to the Measured Number Distributions.
In order to validate the applied models, first the p(q) number distributions have been calculated.For this purpose the number of source events producing q neutrons leaving the target (N q ) has been determined from the PTRAC event file and divided by the total number of source events (N): In [3] fitted functions are given for the differential crosssection of the source protons of the neutron multiplicity dσ/dn, where multiple reactions initiated by the same source proton are counted as one.This quantity needs to be converted to probability distribution per incident particle in order to be comparable with the simulated results: where P reac is the probability that a source proton enters a nuclear interaction in the target, while σ tot is the microscopic total cross-section of the source protons for entering a nuclear interaction in the lead target.p(0) needs to be increased by the survival probability of source protons 1 − P reac .All parameters have been taken from [3].The fitted function is the sum of an exponential and a Gaussian distribution: where the parameters are the ratio of the exponential and the Gaussian part (S E and S G , resp.), the mean value of the exponential and the Gaussian distribution (T n and M max n , resp.), and the spread s of the Gaussian term.
Comparison of the measured and simulated data using four different physical models can be seen in Figures 3, 4 and 5 for the target thickness of 0.2 cm, 5 cm, and 35 cm, respectively.The main conclusion of the comparison is that the simulated results give good approximation of the measured ones and follow the same behavior: the thicker the target, the lower the exponential part, the wider the Gaussian distribution, and more neutrons are produced.This effect is due to the increasing number of secondary reactions.
One can also observe that there is a systematic overestimation of the number of neutrons produced compared to the measured values.This is most probably due to the fact that during the evaluation of the measurements an average efficiency of 85% was used for the neutron detection efficiency, which had been measured by a 252 Cf source, although the spallation source emits also much higher energy neutrons for which the efficiency is lower [3].Concerning the different models one can observe that the very similar Bertini and Isabel provide approximately the same results.CEM03 also provides a similar distribution but with higher probabilities in the Gaussian part.The distribution obtained from INCL4 is shifted toward smaller multiplicities.The results prove that the physical models are suited to reproduce the distributions in all details and therefore they preserve not only the mean values but also the higher moments.The MCNPX default option Bertini model has been chosen for the detailed investigation of the energy correlations.

Calculation of the One-Particle Energy Distribution and
Spectrum.For the calculation of the one-particle energy distributions an energy group structure (E i ) has been set up and the number of neutrons from source events producing q neutrons has been determined in each energy bin (M i q ).The one-particle energy distribution f q (E) can be estimated as: where M q = N q q is the sum of the neutrons from source events producing q neutrons.The average one-particle neutron spectrum χ(E) has also been determined: where M = q M q is the total number of neutrons produced.One can see in Figures 6, 7, and 8. the results of the simulations for the 0.2 cm, the 5 cm, and the 35 cm thick targets, respectively.It is obvious from the figures that the condition (13) does not stand for a spallation source and the energy distribution of the neutrons clearly depends on the number of neutrons produced.The explanation of this effect is related to the energy conservation: the fewer the number of neutrons produced, the higher the average excitation energy per one neutron is.This is the reason why at low neutron numbers (q < 10) the probability of a high energy neutron is much higher than at higher neutron numbers.
The comparison of the results obtained with different target thicknesses shows that in case of the thicker targets the difference between the one-particle energy distributions for low q and for high q decreases.In the case of the 35 cm thick target (see Figure 8) for neutron numbers higher than 25, there is no difference between the distributions, which are very close to the average spectrum χ(E).This can also be seen in Figure 9 where the portion of the low-energy (<100 MeV) neutrons can be seen compared to the highenergy (>100 MeV) neutrons for the different targets.Up to about 25 neutrons/proton this ratio goes together in the three cases, and the spectrum becomes softer as the number of neutrons increases.But above this value it steeply increases for the thin (d = 0.2 cm) target while saturates for the thick ones.This is due to the fact that in the thin target all the neutrons come from one single nuclear interaction while in the thick targets the higher multiplicity is the result of secondary reactions and an average spectrum is produced.This averaging effect of the secondary reactions reduces the differences between the energy distributions for high neutron numbers in thick targets.

Calculation of the Two-Particle Energy Distribution and
Spectrum.It has been shown above that the condition formulated in (13) does not stand for a spallation source, and therefore the factorization of the two-particle spectrum χ(E 1 , E 2 ) is not possible.Now it is also important to investigate whether the energies of the neutrons from a source event are independent from each other and approximations ( 9) and ( 10) can be used.For this purpose all the q(q − 1)/2 possible pairs have been created from the q neutrons escaping the target after a single source event and collected into energy bins according to the two energies E i and E j to determine the number of neutrons in each bin (M i, j q ).In this way the twoparticle energy distribution for q outcoming neutrons can be obtained as: Figure 6: One-particle energy distributions f q (E) for different number of produced neutrons q and average spectrum χ(E) calculated for the 0.2 cm thick target.Figure 7: One-particle energy distributions f q (E) for different number of produced neutrons q and average spectrum χ(E) calculated for the 5 cm thick target.while the two-particle spectrum is given as: In order to express the deviation of the above functions from the case when conditions ( 9), (10), and ( 13) are valid, the following covariance functions have also been calculated: In Figure 10 one can see the two-energy distributions f q (E i , E j ) for different neutron numbers q and the corresponding covariance functions C fq (E 1 , E 2 ) for the thin target.It can be observed that for low neutron numbers the condition (10) is not valid as the distributions are clearly antisymmetric for the E 1 = E 2 axis.This represents an "anticorrelation" of the neutron pairs: if a neutron is emitted with higher energy, then the other one tends to have lower energy.This effect diminishes at higher neutron numbers, due to the much higher degree of freedom the neutron energies become independent from each other.The same effect can be observed for the thick targets, but due to the averaging effect of the secondary particles the "anti-correlation" diminishes even faster.Concerning the two-particle spectra in Figure 11 one can see that the correlation observed in Figure 10 influences also the spectrum.In the case of the thin target the two-particle spectrum χ(E 1 , E 2 ) is also asymmetric.This is demonstrated by the spectral covariance function C χ (E 1 , E 2 ), as well.For the thick targets this asymmetry is not obvious, but the covariance functions show that the effect exists, although in a much smaller extent than for the thin target.This is again due to the many more secondary particles produced in thicker targets.

Method to Investigate the Effect on Noise Measurements
As it was shown above, energy correlations between spallation neutrons exist, although the effect is very small and seems diminishing as the target thickness increases.Therefore, it is important to quantify the effect of these energy correlations on actual noise measurements in order to decide whether this effect needs to be considered or can be neglected by the assumption in (13).In the followings a simulation method is proposed for this investigation.
Noise measurements can be simulated by Monte Carlo calculations as described in, for example, [17].In such calculations the data of the produced neutrons recorded during the above simulations can be used as the neutron source.In order to distinguish the effect of the different approximations, the calculations need to be performed by modified source data also, where the neutrons are redistributed between the source events: (i) by preserving the number distribution p(q) but sampling the neutron energies from the one-energy distributions f q (E) to investigate the effect of the assumption in (9), that is, the independence of the energy distribution of the neutrons from the same source event; (ii) by preserving the number distribution p(q) but sampling the neutron energies from the one-energy spectrum χ(E) to investigate the effect of the assumption in (13), that is, the independence of the energy spectrum of the produced neutrons from the number of neutrons produced in a source event; (iii) by preserving the number distribution p(q) but sampling the neutron energies from the two-energy spectrum χ(E 1 , E 2 ) to investigate whether the second moments are fully preserved in this way.Calculations are in hand with the above described method for a typical lead cooled ADS design to quantify the effect of the energy correlation on neutron noise measurements.

Conclusions
In the present paper the physical models implemented in the MCNPX code have been used for the investigation of the energy correlation in a spallation neutron source.After a successful validation of the model with a neutron multiplicity measurement, the one-and two-particle energy distributions and spectra of the spallation neutrons have been determined.It has been shown that the one-particle energy distributions highly depend on the number of neutrons produced in a source event.Fewer neutrons result in a higher probability of high-energy neutrons as the excitation energy per neutron is higher.It has also been observed that the averaging effect of the large number of secondary reactions in thick targets reduces this effect.The energy correlation between the neutrons from low multiplicity spallation events has been demonstrated also by the two-energy distributions.Finally, one can conclude that the energy correlation is an important effect in a thin spallation target and for the complete description of the neutron fluctuations the two-energy spectrum χ(E 1 , E 2 ) has to be used.On the other hand this effect reduces in thick targets and the factorization of the two-energy spectrum χ(E 1 , E 2 ) = χ(E 1 )χ(E 2 ) appears as an acceptable approximation.Further investigations are in hand to quantify the actual effect of the observed energy correlations on the higher moments of the neutron fluctuations in a subcritical system.

Figure 3 :Figure 4 :
Figure 3: Comparison of measured and simulated distribution of the number of neutrons produced by a source proton for the 0.2 cm thick target.

Figure 5 :
Figure 5: Comparison of measured and simulated distribution of the number of neutrons produced by a source proton for the 35 cm thick target.

Figure 8 :Figure 9 :
Figure8: One-particle energy distributions f q (E) for different number of produced neutrons q and average spectrum χ(E) calculated for the 35 cm thick target.