## 1 Introduction

A neutron flux is a function of time, space, direction and energy of the constituent neutrons. The behavior of neutrons in a medium is accurately expressed by the neutron *transport equation* [1]. However, in a situation involving an extensive amount of calculations as in most reactor-related cases, the neutron *diffusion equation* - which ignores the angular dependence of the neutron flux - is often used to approximate the neutron transport equation to reduce the computational load. However, it is important that, in applying the diffusion approximation, the transport boundary conditions should be adjusted accordingly so that the angular dependence aspect is removed. In this study we mathematically and physically analyze several boundary conditions - such as the zero net-current, the Marshak, the Mark, the zero scalar flux, and the Albedo condition - in diffusion theory and show how they approximate the reflective and vacuum conditions in transport theory [2]. In the process, we clarify these diffusion boundary conditions by adopting proper terminologies and point out the physical implication of each condition. Finally, the numerical results obtained by applying the transport and the diffusion equations are compared against one another with some explanation on their differences and significances. The safe management and disposition of radioactive waste is a fundamental aspect of the nuclear fuel cycle. The accurate calculation of neutron flux using appropriate boundary conditions assures higher computational modeling capability for the performance of storage, transportation, and disposal options for a range of fuel cycle alternatives [3].

## 2 Transport Equation and Transport Boundary Conditions

### 2.1 Transport Equation

The one-dimensional, time-independent Boltzmann neutron *transport equation* is given by

where *φ*(*x*,*μ*,*E*) is the neutron *angular flux*, *x* and *E* are the spatial and neutron’s energy variables, respectively, *μ* is the directional cosine (cosθ) between the neutron’s motion and the *x*-axis as shown in Fig. 1. Moreover, *σ* and *σ _{s}* are the macroscopic total cross section and scattering cross section, respectively, and

*Q*represents the neutron source. To simplify our discussion, we make the assumption of isotropic scattering and source within the given energy group. Under such assumption Eq. (1) becomes

where *ϕ* is the neutron *scalar flux* defined as(3)

in which we have normalized the angle by ${\int}_{-1}^{1}d\mu =2$

### 2.2 Transport Boundary Conditions

Here, we discuss the vacuum and reflective conditions as examples of transport boundary conditions. A similar analysis can be applied in analyzing other types of boundary conditions that are frequently used in expressing the neutron transport behavior.

#### 2.2.1 Reflective Condition

The *reflective condition* describes the neutron’s mirrorlike reflection at interfaces. If we apply the condition on the left boundary of Fig. 1 for example, we have

where *μ* is the incidence angle and *μ*' is the reflection angle of neutron. This condition is frequently used, as it imposes symmetry and thereby reduces the size of the medium in which calculations are carried out. Eq. (4) can also be written more generally

referred to as the transport *Albedo condition*. In particular when *α*=1, Eq. (5) reduces to Eq. (4) [4].

#### 2.2.2 Vacuum Condition

The *vacuum condition*, sometimes called the *free surface condition*, refers to a situation in which no neutrons flow into the medium from the outside as in the case of vacuum for example. Therefore, the vacuum condition applied to the right boundary in Fig. 1 can be expressed as

If the value of the incoming angular flux *f _{μ}* is known then Eq. (6) becomes(7)

and it is called the *known incoming condition* or *nonhomogeneous condition*.

## 3 Diffusion Equation and Diffusion Boundary Conditions

### 3.1 Diffusion Equation

Here, we briefly describe the procedure of deriving the diffusion equation from the transport equation. A similar technique will be used later to derive the diffusion boundary conditions from the transport boundary conditions. First, we expand the angular flux *φ*(*x*,*μ*) in terms of the directional cosine *μ* by the *P*_{1} -approximation

Here, *P*_{1} indicates the two-term expansion using Legendre polynomials *P _{n}*(

*μ*) as the expansion function. From the orthogonal property of Legendre polynomials(9)

the expansion coefficients are given by ${\varphi}_{n}\left(x\right)=\frac{1}{2}{\displaystyle {\int}_{-1}^{1}{P}_{n}\left(\mu \right)\phi \left(x,\mu \right)d\mu}$. As *P*_{0}(*μ*) = 1, *P*_{1}(*μ*) = *μ* we can easily obtain the first two expansion coefficients as(10)

Note that the two coefficients of the *P*_{1}-expansion above correspond precisely the *scalar flux**ϕ*(*x*) and the *current J*(*x*) by definition in diffusion theory. Hence from Eq. (8) the angular flux can be written as

If we substitute Eq. (12) into the transport equation [Eq. (2)] and integrate over the interval -1 ≤ *μ* ≤ 1:

$\begin{array}{l}\frac{1}{2}{\displaystyle {\int}_{-1}^{1}\mu \frac{d}{dx}[\varphi (x)+3\mu J(x)]d\mu +\frac{1}{2}{\displaystyle {\int}_{-1}^{1}\sigma (x)[\varphi (x)+3\mu J(x)]d\mu}}\\ =\frac{1}{2}{\displaystyle {\int}_{-1}^{1}{\sigma}_{s}(x)[\varphi (x)d\mu +\frac{1}{2}{\displaystyle {\int}_{-1}^{1}Q\left(x\right)d\mu}}\end{array}$

then we have

Here, *σ _{a}* =

*σ*-

*σ*represents the macroscopic absorption cross section. Next, if we multiply both sides of Eq. (2) by

_{s}*μ*and integrate over -1 ≤

*μ*≤ 1:

$\begin{array}{l}\frac{1}{2}{\displaystyle {\int}_{-1}^{1}{\mu}^{2}\frac{d}{dx}\varphi (x)+3\mu J(x)]d\mu}\\ +\frac{1}{2}{\displaystyle {\int}_{-1}^{1}\mu \sigma (x)[\varphi (x)+3\mu J(x)]d\mu}\\ =\frac{1}{2}{\displaystyle {\int}_{-1}^{1}\mu {\sigma}_{s}(x)\varphi (x)d\mu +\frac{1}{2}{\displaystyle {\int}_{-1}^{1}\mu Q\left(x\right)d\mu}}\end{array}$

we obtain

Solving Eq. (14) for *J* yields the famous *Fick*’*s law*:

In the above the diffusion coefficient is defined by *D* = 1/3*σ*. Substituting Eq. (15) into Eq. (13) yields(16)

The last equation is the well-known neutron *diffusion equation* for the scalar flux *ϕ*(*x*). This equation is an approximation to the neutron transport equation that does not employ the angular variable *μ* [5].

### 3.2 Diffusion Boundary Conditions

Since the transport boundary conditions (i.e., the reflective and vacuum conditions expressed by Eq. (4) and Eq. (6), respectively) contain the angular variable *μ*, they cannot be used in the current form for the solution of the diffusion equation which precludes the angular variable. The transort conditions should be also modified so that *μ* is eliminated. The procedures involved are as follows.

#### 3.2.1 Diffusion Approximation to Reflective Condition

The *P*_{1}-approximation can be applied similarly as in the case of the reflective condition [Eq. (4)] in transport theory. The result is

Since *μ*' = -*μ* in one-dimensional reflection, Eq. (17) implies *J*(*a*) = 0 which, based on Eq. (15), is equivalent to(18)

This condition is called the *zero net-current condition*, which is an approximation of the reflective condition in transport theory.

#### 3.2.2 Diffusion Approximation to Vacuum Condition

Many versions of the diffusion model can be derived from the transport vacuum condition depending on how the angular variable *μ* is approximated. Here, we show the derivation procedures for the Marshak and Mark conditions, two of the most frequently used conditions in diffusion calculations.

##### 3.2.2.1 The Marshak Condition

The transport vacuum condition applied to the right boundary in Fig. 1 implies that there is no angular flux flowing into the medium. Thus, in a diffusion model, the vacuum condition can be expressed by setting the neutron current heading into the medium from outside to be zero. Therefore, using the definition of current [Eq. (11)], we can write the vacuum condition as

where the superscript ‘-’ denotes the negative x-direction. Following the treatment similar to that presented earlier the *P*_{1}-approximation of φ(*b*,*μ*) in Eq. (19) leads to(20)

which completes the *Marshak condition*

##### 3.2.2.2 The Mark Condition

For the Marshak condition, *J*^{-} on the boundary was appoximated by the integral of the angular flux over the interval -1 ≤ *μ* ≤ 0. However, for the Mark condition, *J*^{-} is expressed directly by the specific angular quadrature used in the *discrete ordinates method*, the so-called *S _{N} method* [1,2] . The

*P*

_{1}-approximation to

*φ*(

*b*,

*μ*) in Eq. (6) yields

If we use the simplest *S*_{2} angular quadrature *μ* = -1/$\sqrt{3}$ (taking negative value since -1 ≤ *μ* ≤ 0), Eq. (22) becomes to(23)

which leads directly to the *Mark condition*

The Marshak condition of Eq. (21) and the Mark condition of Eq. (24) will not lead to significantly different results as 1 /$\sqrt{3}$ = 0.5774 ≈ 1/2 [6].

#### 3.2.3 Albedo Condition

*Albedo condition* is widely used in commercial diffusion codes because of its flexibility in being able to accommodate various types of boundary conditions. In reactor theory, the value of Albedo α is defined as the ratio of the incoming current to the outgoing current as

Substituting the relation ${J}^{\pm}=\frac{1}{4}\varphi \pm \frac{1}{2}J$ into Eq. (25) yields

If we set *α* = 0 in the above expression, Eq. (26) reduces to *J* = *ϕ*/2, which is equivalent to the Marshak condition. On the other hand, if we let *α* = 1, then Eq. (26) simplifies to *J* = 0, the zero net-current condition. Finally, when *α* = -0.072 (even though the negative value of *α* has no physical meaning by its definition), *J* = *ϕ*/${J}^{\pm}=\frac{1}{4}\varphi \pm \frac{1}{2}J$ is equivalent to the Mark condition. Often, the Albedo condition is also expressed as

In this form *β* = 0 corresponds to the zero-net current condition, *β* = 0.5 the Marshak condition, *β* = 1/$\sqrt{3}$ the Mark condition, and finally, *β* = ∞ (i.e., *ϕ* =0) the zero scalar flux condition.

## 4 Numerical Results

In this section, we compare the various diffusion boundary conditions discussed above by the numerical results for the sample problem shown in Fig. 2. The length of the medium is 15 cm, which consists of three different regions: a neutron source region, a scattering-dominant region, and an absorption-dominant region. The numerical values of the macroscopic cross sections and the neutron source rate for each region are tabulated in Table 1. The reflective and vacuum conditions are applied on the left and right boundaries, respectively. The reflective condition is approximated by the zero net-current condition and the vacuum condition is approximated by the Marshak, Mark, and zero scalar flux conditions. The finite difference method (FDM) is used for the solution of the diffusion equation. The transport equation is also solved to serve as a reference calculation. The classical diamond-difference method (DD) and the S_{16} angular quadrature set are used for this purpose. To reduce the spatial discretization effect, the medium is divided by 60 fine meshes each with the size of 0.25 cm.

Scalar fluxes computed with different models are shown in Fig. 3. All the diffusion model results show a good agreement with that of the transport model in regions of high neutron population but there are noticeable mismatches at region interfaces and near the vacuum boundary. This is expected, because the diffusion models have not incorporated the directional character of the neutron flux, i.e., the *μ*-dependence. On the left boundary, we can see that the transport reflective condition is well-approximated by the zero net-current condition. On the right boundary the results obtained from applying the Marshak and Mark conditions show similar trends as expected, whereas there are clear discrepancies for the zero scalar flux condition. As is known in reactor theory, the zero scalar flux condition is not a robust approximation to the vacuum condition. This is the reason for the need in the past to employ the socalled *extrapolated distance* in application of the zero scalar flux condition. Table 2 is a tabulation of the scalar fluxes integrated over the length of the medium, which is useful in computing the total reaction rate in the medium. Notice that the relative error of integrated fluxes between the transport and diffusion models are only about 0.78%. This apparent inconsistency can be explained by noting that, in Fig. 3, the flux level tends to be very low near the vacuum boundary compared with the rest of the medium. Note that the y-axis is drawn on a logarithmic scale. From these findings we may conclude that the diffusion models with application of well-formulated boundary conditions lead to results which do exhibit some discrepancies in the flux distribution but are accurate enough for the calculation of total reaction rate in the given medium.

## 5 Discussion

Various diffusion models were formulated for the transport boundary conditions using the *P*_{1} and the *S*_{2} approximations. While these approaches are not new in themselves, we were able to present a direct comparison of these boundary conditions based on common mathematical and physical grounds as well as verification of these conditions’ validity. Analysis of the multi-dimensional case or other boundary conditions such as periodic and white boundary conditions, though more complex in nature, may be considered for future work.