A Classical and Semi-analytical Analysis of Isotropic Scattering of Neutrons through a Protective Wall of Iron-via Two Energy Groups-for Two Examples of Thicknesses

Introduction

In this paper, we introduce an analytical method of calculating transmission, or conversely shielding coefficients, and contained therein an effective method of predicting local internal dose due to beam. Here in, the examples of an iron wall (i.e. of the main Fe-56) with thicknesses of 0.865cm and 2.88cm are examined for predictions of transmission ration and gross back-scatter ratio. In service to the field of health physics, the essential background topics in medical physics, and the reviews of radiological safety of reactors of nuclear reactor operations, various radiation related calculations have been conducted with reasonable caution for over 70 years. These calculations include: radiation shielding calculations, penetration assessments [1], and radiation dosimetry calculations [2]. Tying into the physics of interaction of radiation with matter, a beam which is comprised of either energetic neutrons or high-energy photons (i.e. gamma rays or hard X-rays) have conventionally been looked upon as candidates for interception via the standard X-section (i.e. cross-section) inspired models, for which one uses the typical expression

$$I(z) = I(0)e^{-\mu z}$$

where $z$ is the depth and $\mu$ is the attenuation parameter or 'coefficient'.

The "simplistic" attenuation formula is an appropriate name for this short Equation (1). Many shielding calculations done by medical physicists and health physicists over the years (since at least 1970) have been done using an extension of this "simplistic attenuation formula" coupled with a buildup factor [3]. The buildup factor is necessary if we are to use an attenuation formula as our principal tool of analysis to account for scattered neutrons or energetic photons, respectively [4]. Some of these “n’s” and “y’s” are almost (or totally) elastically scattered, but some of these are down-scattered in energy. In keeping with the aspiration for excellence from the era of the “Space Age”, some computer/electrical engineers and nuclear engineers have enhanced these efforts by conducting Monte Carlo simulations of the transport and penetration of neutrons or photons through walls and various barriers with various respected codes (packages with X-section libraries) such as MCNP, EGS4, EGSNrc, and the versatile but highly tedious GEANT4 [5, 6, 7]. In this short paper, we introduce and present a deterministic and semi-analytical method for doing transport analysis on neutrons and isotropically scattering ‘hard’ photons which are placed in two energy and, with future ambitions, into 3 energy groups. There are advantages for doing such 2-group and higher multi-group analysis of radiative particles (i.e., neutrons and photons). These advantages are that we can more directly keep track of what percentage of radiative particles are close to the original high energy and how many are at significantly lower energy. An inspection of the profile of any build up function shows that the function is slightly larger than 1.0 at entry, then it rises to perhaps 2 or 3 within roughly one mean free path of the fast primary particles, and finally approaches the asymptote of 1.0 as the penetration depth gets progressively larger.

Review of Buildup and Discussion of Methods of Analysis

The attenuation formula in the introduction expresses the particulate intensity, not the energetic intensity in our convention. This corresponds to the choice of analysis of particulate flux rather than energetic flux of radiation. For those with a non-nuclear background: “Flux” is used by the medical physics community and health physicists with a definition considerably different from that of the flux of electric fields. Our flux has units of “particles” per cm$^2$ per second. See Frank Attix’s text [8] if this is unclear. Particulate intensity is less than or equal to the scalar flux of the particles. Very often intensity is defined as the magnitude of net current of transported particles per cm$^2$ of surface per sec. Indeed, in a case where equally many particles approach and penetrate a wall bidirectionally the net current is zero. However, the scalar flux in such an example is much larger than zero. Admittedly, one can make some inferences on the approximate ratio of down-scattered particles at a given depth as a function of position by inspecting the build up function, should it be available in published tables for a given shape and material. Here the intensity with a buildup coefficient can be expressed as

$$I(z) = B(E,z)I(0)e^{(-\mu z)}$$

(2)

However, the buildup coefficient (i.e. B(E,z)) does not make clear just what percentage of scattered $\gamma$-ray or X-ray radiation is scattered so as to retain most of its energy and how much has been “demoted” to photons with 50% less or more of energy per radiative particle. This reality (of radiative particles often undergoing elastic or nearly elastic scattering) holds for ‘free’ neutrons which scatter off nuclei with an Atomic Number greater than 4. For example, when a ‘n’ with a kinetic energy of 1.0MeV collides with an Fe-56 nucleus, it has a 99% chance of undergoing elastic or quasi-elastic scattering (retaining most of its 1MeV). This same ‘n’ has a probability of less than 0.8% of down-scattering to a “low” energy neutron with less than 0.201MeV of kinetic energy. The chance of capture at 1MeV is less than one in a 1000. Fast neutrons are not easily captured, they usually are just scattered. This information is given based on inspecting generally available public claims of two group data and, more importantly, by our careful conducting of MCNP [5] simulations in which we reproduce the conditions of scattering and so called “buildup” of scattered neutrons in rectangular walls of Iron (Fe-56). We are authorized to use and has extensive experience with the very versatile MCNP Monte Carlo code. Regarding neutrons, the buildup coefficient data either is not widely published for nuclear engineers or not readily available. Thus, the case of neutron shielding analysis offers a major service for performing multi-group energetic neutron flux and dose calculations. In this paper we stick with 2-group n’s (i.e., neutrons). The following formula is an overly simplistic & inadequate expression often used for metals with atomic number less than 84, away from the “uranic” family

$$\Phi(z) = \Phi(0)e^{(-\sum_{tot,at,1MeV}z)}$$

(3)

However, Equation (3) gives a very incomplete story of the local neutron flux. Just as Equation (1) of I(z), which equals I(0)exp(-\mu z), gives an incomplete story of local current/intensity of X-rays and $\gamma$-rays, when one fails to include the Buildup factor included in Equation (2) where B(E,z) in Equation (2) is the buildup factor for the calculation of local ‘flux’ of $\gamma$-rays or photons at various depths (e.g. z) of penetration. It is more insightful to replace I(z) from Equation (1) with $\Phi(z)$ in Equation (3) where Flux is more appropriate than ‘n’ or photon current in the geometry of a wall or box since many of the neutrons or photons no longer travel straight forward along the z-axis after one or two collisions of scatter occur.

If the nuclei/atoms of a medium which is entered by the beam of n’s or photons is a pure absorber, then Equation (1) and Equation (3) are acceptable solutions for penetration and the differential equation which explains the transport of the particle(s) is given by:

and where $n_{\text{rad}}$ is the number density from radiation. If the medium is Boron-10 and the neutrons are at low energy, then Equation (4) would be a realistic equation for modeling transport of the neutrons, because B-10 is nearly a pure absorber. However, most materials are not pure absorbers.

If we presume that there are two energy levels for neutrons and photons (i.e., fast n's and slow n's), then it is appropriate to write a double energy-group Maxwell Boltzmann Transport Equation (MBTE) in order to express what is going on for transmission and for energy demotions (i.e., down scatters). This pair of equations is written below for the transport of neutrons with two possible energy levels [9, 10, 11]. Note in Equation (5) that $\Psi_{[2]}$ is proportional to $1/(4\pi)\Phi_{[2]}$ for fast neutrons (of grp-2). Here $\Phi_{[2]}$ is the scalar flux which includes neutrons in the energy grouping of 0.201MeV up to 10MeV (at least when we model & analyze iron shielding). $\Phi_{[2]}$ is the scalar flux which includes all neutrons of energy 0.20MeV and lower. Here is the two-group MBTE:

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = \int \Sigma_{[r,2]} \cdot \Psi_{[2]}(\bar{r}, \theta) \cdot \left( \Sigma_{[s,2]} + \Sigma_{[j,2]} \right) d\theta$$

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = 0.$$

Suppose that $\Sigma_{[s,2]}$ and $\Sigma_{[s,2]}$ are equal to zero, as we might imagine for a medium of "super-Boron". Then, Equation (5a) can be rewritten as:

$$\bar{n} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = 0.$$

If our wall is very broad and if the distribution to approximation depends only on coordinate $z$, then this simple 1st order equation is completely equivalent to Equation (4) on attenuation on the 2nd page. In this paper, we presume that scattering is isotropic, which is often a good approximation for the scattering of neutrons. Therefore, we establish that $\Sigma_{[s,2]}$ and $\Sigma_{[a,2]}$ are completely independent of $\theta$, $\theta'$, and any angle. This 2-group version of the MBTE is an example of a pair of integrodifferential equations. Generally, it is easier to solve a purely integral equation, such as a Fredholm Int. Equation [12]. Holding on to the presumption of isotropic scattering, Equations (5a) and (5b) can be subjected to a special integral transformation via Green's functions in order to re-express them as the following integral equations.

$$\frac{d\left(I_{[2]}(z)\right)}{dz} + (\mu) \cdot I_{[2]}(z) = 0, \quad \text{for } \mu = \sigma_{\text{absorb}} n_{\text{rad}}$$

and where $n_{\text{rad}}$ is the number density from radiation. If the medium is Boron-10 and the neutrons are at low energy, then Equation (4) would be a realistic equation for modeling transport of the neutrons, because B-10 is nearly a pure absorber. However, most materials are not pure absorbers.

If we presume that there are two energy levels for neutrons and photons (i.e., fast n's and slow n's), then it is appropriate to write a double energy-group Maxwell Boltzmann Transport Equation (MBTE) in order to express what is going on for transmission and for energy demotions (i.e., down scatters). This pair of equations is written below for the transport of neutrons with two possible energy levels [9, 10, 11]. Note in Equation (5) that $\Psi_{[2]}$ is proportional to $1/(4\pi)\Phi_{[2]}$ for fast neutrons (of grp-2). Here $\Phi_{[2]}$ is the scalar flux which includes neutrons in the energy grouping of 0.201MeV up to 10MeV (at least when we model & analyze iron shielding). $\Phi_{[2]}$ is the scalar flux which includes all neutrons of energy 0.20MeV and lower. Here is the two-group MBTE:

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = \int \Sigma_{[r,2]} \cdot \Psi_{[2]}(\bar{r}, \theta) \cdot \left( \Sigma_{[s,2]} + \Sigma_{[j,2]} \right) d\theta$$

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = 0.$$

Suppose that $\Sigma_{[s,2]}$ and $\Sigma_{[s,2]}$ are equal to zero, as we might imagine for a medium of "super-Boron". Then, Equation (5a) can be rewritten as:

$$\bar{n} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = 0.$$

If our wall is very broad and if the distribution to approximation depends only on coordinate $z$, then this simple 1st order equation is completely equivalent to Equation (4) on attenuation on the 2nd page. In this paper, we presume that scattering is isotropic, which is often a good approximation for the scattering of neutrons. Therefore, we establish that $\Sigma_{[s,2]}$ and $\Sigma_{[a,2]}$ are completely independent of $\theta$, $\theta'$, and any angle. This 2-group version of the MBTE is an example of a pair of integrodifferential equations. Generally, it is easier to solve a purely integral equation, such as a Fredholm Int. Equation [12]. Holding on to the presumption of isotropic scattering, Equations (5a) and (5b) can be subjected to a special integral transformation via Green's functions in order to re-express them as the following integral equations.

$$\frac{d\left(I_{[2]}(z)\right)}{dz} + (\mu) \cdot I_{[2]}(z) = 0, \quad \text{for } \mu = \sigma_{\text{absorb}} n_{\text{rad}}$$

and where $n_{\text{rad}}$ is the number density from radiation. If the medium is Boron-10 and the neutrons are at low energy, then Equation (4) would be a realistic equation for modeling transport of the neutrons, because B-10 is nearly a pure absorber. However, most materials are not pure absorbers.

If we presume that there are two energy levels for neutrons and photons (i.e., fast n's and slow n's), then it is appropriate to write a double energy-group Maxwell Boltzmann Transport Equation (MBTE) in order to express what is going on for transmission and for energy demotions (i.e., down scatters). This pair of equations is written below for the transport of neutrons with two possible energy levels [9, 10, 11]. Note in Equation (5) that $\Psi_{[2]}$ is proportional to $1/(4\pi)\Phi_{[2]}$ for fast neutrons (of grp-2). Here $\Phi_{[2]}$ is the scalar flux which includes neutrons in the energy grouping of 0.201MeV up to 10MeV (at least when we model & analyze iron shielding). $\Phi_{[2]}$ is the scalar flux which includes all neutrons of energy 0.20MeV and lower. Here is the two-group MBTE:

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = \int \Sigma_{[r,2]} \cdot \Psi_{[2]}(\bar{r}, \theta) \cdot \left( \Sigma_{[s,2]} + \Sigma_{[j,2]} \right) d\theta$$

$$\bar{n}_{\theta} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta) + \left( \Sigma_{[a,2]} + \Sigma_{[s,2]} \right) \cdot \Psi_{[2]}(\bar{r}, \theta) = 0.$$

Suppose that $\Sigma_{[s,2]}$ and $\Sigma_{[s,2]}$ are equal to zero, as we might imagine for a medium of "super-Boron". Then, Equation (5a) can be rewritten as:

$$\bar{n} \cdot \nabla \Psi_{[2]}(\bar{r}, \theta

thicker samples MCNP simulations have been verified with converging persuasiveness that the escape ratio is in the range of ½ through 0.65, depending on thickness of shielding. $I_{[battered]}$ refers to the particulate current density of n's which have undergone at least one collision before escaping from either the left-hand side (LHS) or right-hand side (RHS) of the slab. Many books casually refer to the 'I' in $I_{[battered]}$ and $I_{[entering, beam]}$ as intensity. In this paper, 'I' shall mean current density of particles (mostly neutrons here) per cm$^2$ (or m$^2$) per second. If the beam and scattered particle are all monoenergetic, then, of course, energetic I, or energetic intensity, is given by the product of 'I' (of 'n') and Energy $[of n']$.

Computations and Results

Following our earlier work Steinfelds EV, et al. [13, 14], we can iteratively solve integral equations (6a) and (6b) in the case where we have one very broad rectangular slab. If there is a monoenergetic beam of fast neutrons which enter the designated slab of shielding, then before the first iteration we find that $\Phi_{[2][0]}(z)$ goes as $I(0)e^{(-\sum_{i=1}^{n}M_i z)^2}$, as in Equation (3) and where $I(0)$ and $\Phi[0]$ were interchanged. $\Phi_{[2][0]}(z)$ is that part of the flux of neutrons at z that have never undergone scattering. The first index, which holds [2], shall denote that this is flux/current of the fast group of neutrons. In the first iteration (upon the event of scatter) we find the formula of $\Phi_{[n][1]}(z)$ for the n-group's neutrons where n equals either 1(slow) or 2(fast). Upon integration this flux can be found analytically for the n-group of neutrons. Consider Equation (7)

$$\bar{\Phi}_{[2]}(\bar{r}) = I(0)e^{-(\sum_{i=1}^{n}M_i z) / \cos(\theta)}$$

where $\Phi_{[2][n]}(z)$ is generated iteratively from the input of $\Phi_{[2][n-1]}(z)$ into the RHS of Equation (6). The writing out of the summation of $\Phi_{[2][n]}(z)$ with respect to iteration index n is akin to the early methods of perturbative quantum electrodynamics done by Hans Bethe circa 1949 [15], but only without the need to subtract infinities out. Renormalization methods are not required to solve Equations (6), analytically or computationally. The first iteration, $\Phi_{[2][1]}(z)$, contains logarithms of z, EI functions of z, and terms of exp($-z$) z for several factors. The "2" refers to the 1-MeV neutrons and can be efficiently approximated by a much more tractable sum of a polynomials and log expressions. $\Phi_{[2][1]}(z)$ is the flux of fast neutrons which have been scattered only once. We allow batter-1 to refer to radiative particles which have undergone collisions only once. Following through with the iterations, $\Phi_{[2][2]}(z)$ is the flux of fast neutrons which have been battered twice. Likewise, $\Phi_{[2][3]}(z)$ is the flux of fast neutrons which have been battered three times. In general, $\Phi_{[m][n]}(z)$ is the flux of neutrons of group 'm' which have been battered n times without being either demoted from nor promoted from its energy-group "m". 'Fast' albedo is the ratio of the sum of fast scattering cross sections divided by the total fast cross section of neutrons. So long as the albedo is less than or equal to 1, the series of $\Phi_{[n][n]}(z)$ is guaranteed to converge. Five to eight iterations have proven to be required for sufficiently convergent solutions. We use the approximation that the current of escaping neutrons which have penetrated the wall is given by the sum of the 'unbattered' flux plus half (i.e. 0.5) of the sum of fluxes comprised of scattered neutrons, where 0.5 is the analytically and geometrically guaranteed minimum. Judiciously, we occasionally replace 0.5 with 0.52 or even 0.55, where this is the dimensionless factor of escape which relates surface $\Phi$ value to current density of scattered neutrons which depart from the surface.

For the sake of space and the context of this paper, we will focus on predictions for scattering deterministically of the transmission of neutrons through sample slabs of iron at various thicknesses and on the prediction for the portion of neutrons which are returned backwards (due to back-scatter) from the slab. Our deterministic method distinguishes between the population of fast scattered neutrons and of slow neutrons (where the choice of '1' for 'm' designates the down-scattered (i.e., slow) neutrons). We also conducted Monte Carlo simulations of neutrons from a beam which approach the same wall of iron material. The two Monte Carlo (MC) codes used are: MCNP (which was developed and updated by LANL) and a custom Monte Carlo code which was developed early in 2014 [13, 14] and had proven to be valid for the modeling of isotropic scattering. In Table 1 the name assigned to our custom MC code is SMUSKE, where SMUSKE is designed to simulate isotropic scattering or radiative particles in arbitrarily chosen rectangular geometries, as done in the past.

Two examples are given below for a broad rectangular slab of homogeneous cast-iron where a beam of neutrons approaches the barrier at normal incidence. Our deterministic predictions, the predictions of MCNP, and the predictions of the MC code, SMUSKE, are included in Table (1). In the common style, MFP shall denote Mean-Free-Path of neutrons in response to the total cross section for the fast (0.9 MeV to 1.0 MeV) neutrons. Our deterministic algorithm is named "Integrgiterator", which stands for Iterative Integral Equation Solver. "Deterministic Iterator" is a slightly lengthy alias for our abbreviation of Integrgiterator. As a reminder, Equations (6a) and (6b) are Fredholm integral equations in terms of mathematical structure. "Mfpm" in column one means the value of MFP/(Wall Thickness). "Grp#' in the 2nd column of the table refers to the energy group number of the neutrons. Grp# 2 includes all neutrons which are in the energy range of 0.9 MeV through 1.0 MeV. These are the designated fast neutrons. Grp# 1 denotes the slow neutrons, which include the ‘slow’ Group and the ‘fast’ Group of neutrons, therefore, we do not include such an intermediate group of negligible population.

Table 1: Summary of Forward Scattering and Back-scattering Predictions based on the Monte Carlo Codes of SMUSKE and MCNP and on the IntegIterator Deterministic Algorithm. In the first half of Table (1) the thickness of the rectangular slab of iron is 1 MFP long, which is 2.884 cm. According to SMUSKE, out of 10,000 incoming fast neutrons, 3132.8 fast n’s travel out or backward from the wall, and slow n’s escape back toward the source of the beam. Accordingly, SMUSKE predicts that 6,691.7 fast n’s escape forwards out through the iron wall, and 84 slow n’s escape forward. Comparing the Deterministic Iterator and SMUSKE’s predictions starting with the thick wall at 2.884 cm first, we find that for the collective back-scattered n’s, which escape out of wall backward, the percent differences between the Deterministic Iterator and SMUSKE are: 45.6% for slow n’s (i.e., of grp.1) and 14.6% for fast n’s (i.e., of grp.2). For the rate of forward escape, or transmission, the percent differences between IntegIterator and SMUSKE are: 37% for slow n’s (i.e., of grp.1) and -4.85% for fast n’s (i.e., of grp.2). We also consider the wall which is at 3 tenths of an Mfpm at 0.8652 cm, giving a rate of forward escape, or transmission, the percentage of disagreement between our IntegIterator and SMUSKE as: 103% for slow n’s and -2.03% for fast n’s (i.e., of grp.2).

Next, we compare the escapee numbers of SMUSKE to those of MCNP for n’s: the rate of backward escape, returning to the beam source, when the wall is one MFP thick, the percentage of disagreement between the predictions of MCNP and SMUSKE are: -168% for slow n’s and 77% for fast n’s. For the rate of forward escape, when the wall is one MFP thick, the percentage of disagreement between the predictions of MCNP and SMUSKE are: -164% for slow n’s and 0.448% for fast n’s. For the rate of backward escape, returning to the external beam source, when the wall is 3/10 of a MFP thick, the percentage of disagreement between the predictions of MCNP and SMUSKE are: -172.4% for slow n’s and 75.1% for fast n’s. For the rate of forward escape, when the wall is 3/10 of a MFP thick, the percentage of disagreement between the predictions of MCNP and SMUSKE are: -172.6% for slow n’s and 0.349% for fast n’s.

On the other hand, on a more impressive note, there is only a -5.30% disagreement between the respective predictions of our Deterministic Iterator (i.e., IntegIterator) and those of MCNP for the number of forward transmitted fast neutrons through the wall which has 1 MFP (2.884 cm) of thickness. In addition, there is a -2.387% disagreement between the respective predictions of IntegIterator and MCNP for the number of forward transmitted neutrons through the wall with thickness of 3/10 of an MFP (0.8652 cm).

One can see that our deterministic IntegIterator makes predictions which are either as close as or almost as close in their respective predictions to those of SMUSKE to the predictions of forward escape and backward escape generated by using MCNP alone. IntegIterator and SMUSKE have comparable operation speeds if one is content with 1.5 percent statistical fluctuations of Monte Carlo from using SMUSKE. However, SMUSKE does not easily lend itself to decisive mapping of internal flux at high resolution. IntegIterator does, by default, offer the feature of internal mapping. We observe the effectiveness of SMUSKE and IntegIterator in spite of the disadvantage of using isotropic cross sections used in these codes. MCNP does consider angular probabilities in great detail for ‘n’-scattering for all of the well-known isotopes on the complete table of nuclides.

Indeed, in a nuclear historical context, the compilation of the ‘scatter-kernel’ data of MCNP was a project which spanned more than a decade. Thus, it would be very difficult to summarize such a large amount of angular data with approximations of the zeroth, first, and second order Legendre polynomials of the cosine of scatter-angle with any tractable and manageable database which could function without incurring “data strangulation” of an analytical iterator (e.g., our IntegIterator) written in Maplesoft, high- level Python, or similar language/package. Some ask, “why not just give the M.C. ‘jobs’ to GEANT4 to execute?”. With all due acknowledgement of the formidable abilities of GEANT4, such as readily offering options of angular data streaming, the three smaller codes SMUSKE, MCNP, and IntegIterator are all faster and easier to work with than GEANT4 for designated rectangular walls of metal bombarded by n’s.

Conclusions

Our deterministic iterative algorithm, IntegIterator, agrees reasonably well regarding the prediction of the energy distribution and the direction distribution to the corresponding distributions of energy and overall forward direction which are generated via our MCNP simulations, in spite of the deficiency of IntegIterator not being able to process anisotropic scattering in terms of the design of a ‘scatter’ kernel or helpful Green’s Function. This can be seen from the results posted in Table (1) and in the summary of Table (1) above. On the other hand, neither our deterministic code nor our isotropically designed SMUSKE agree extremely well with the predictions of back-scatter results from MCNP. All three methods of calculation agree well on the prediction of percent neutrons (sum of quasi-elastically scattered and down-scattered neutrons) which escape from the target slab of our ferrous metal informatively showing the trend of almost succeeding to sustain the population of free neutrons (albeit with records of past collisions and occasional energy reductions). In reference to “surviving” n’s, ‘free neutron’ means neutron not captured by any atom. It is evident from the slowing down that every one of the iron slabs in Table 1 incidentally functions as a neutron moderator. However, the directionality and ratio of down-scattering are subject to disagreement. As explained above, MCNP has extremely detailed cross section libraries. Iron-56 turns out to be one of these extremely anisotropically scattering isotopes. With the exception of the work of Chandrasekhar [16] and his use of H-functions for tracing intensities of scattered photons, there is little analytical work recorded on predicting flux densities theoretically which are solutions to the MBTE in which the scattering cross sections are anisotropic. The discrete ordinate method used to solve the MBTE is impressive in its flexibility and is deterministic, but that method is not analytical - thus being completely numerical at every step, unlike our IntegIterator which offers some “term- by-term” perspective for the theorist. It is ‘somewhat’ easy to do analysis of solutions of the MBTE when the neutron scatterers (i.e., the nuclei) are isotropic. Many experts of reactor physics and shielding analysis approximate the multi-group MBTE as a multi-group coupling of two, 3, or even 6 simultaneous diffusion equations of radiative particles. By its intrinsic nature, it is virtually impossible to do angle dependent “ray track tracing” of neutrons if one does modeling with a diffusion equation (or equivalently a pair of diffusion equations) rather than the MBTE.

It would be convenient from a clinical radiological treatment planners’ point of view to carefully look up the data for B(E,z), where B(E,z) is the buildup factor included in Equation (2). However, in regard to neutrons, buildup coefficient data either is not widely published or is not available to the broad national/international communities of health physicists or engineers. Moreover, a significant benefit of our deterministic method (and algorithm) of IntegIterator is the superior speed which it offers by its retention of the definitions of chosen geometries per wall to be bombarded and in its output extraction compared to the time and duties required at the conclusion of a corresponding run of an MCNP simulation. Similarly using GEANT4 is often even more involved and tedious than MCNP. Thus our 2-group IntegIterator algorithm and formulation is much faster than MCNP when predicting penetration ratios as well as the distribution of energy of penetrating neutrons [17].

In regard to the CPU (one or multi) time of a single ‘run’ of modeling transmissions and down-scattering of neutrons through rectangular walls, MCNP is sufficiently fast. However, the processing of the output data from the output files generated by the MCNP code requires considerable data processing which is done best either in a UNIX console environment or a DOS console environment. In order to conduct our local Flux calculations, we selected Maplesoft for the code and kernel of IntegIterator. Our ‘Maple’ version of IntegIterator can operate within the environments of Windows, XWindows of Mac OS, and Linux - as valid versions of Maplesoft can be placed in these OS’s. Much of our Maplesoft code can be translated into Matlab code, for the accommodation of the preferred software environment of many electrical engineers for Matlab. Another, great benefit of our deterministic code is that one can procure a polynomial approximation of the local dose of neutron flux at any depth within the metal. The structure of the spatial internal solution for flux is a combination of log and polynomial terms. With MCNP such a feat of mapping flux as a function of the depth in a wall would require updated writing of an Input file which is more than ten-fold more elaborate than the input file for the IntegIterator code for the same slab of metallic material. Internal depth profiling with SMUSKE also is a challenge, but less so than it is with standard MCNP input file declarations.

It is reasonable to anticipate a future effort of analysis of 3-group neutron flux distributions with respect to energy by radiation transport theorists. However, for now, we focus on constructing 2-group databases of neutron cross sections of various important materials besides just iron and boron and subsequently using the algorithm and code(s) of IntegIterator to predict collective forward escape and backward escape of neutrons which initially enter slabs of the respective materials of interest.