## Solutions of the Diffusion Equation – Non-multiplying Systems

**diffusion theory**is widely used in core design of the current Pressurized Water Reactors (PWRs) or Boiling Water Reactors (BWRs). This section is not about such calculations, but provides an

**illustrative insights**, how can be the neutron flux distributed in any diffusion medium. In this section we will solve diffusion equations:

in various geometries that satisfy the **boundary conditions** discussed in the previous section.

We will start with simple systems and increase complexity gradually. The most important assumption is that **all neutrons** are lumped into a** single energy group**, they are emitted and diffuse at **thermal energy (0.025 eV)**.

In the first section, we will deal with neutron diffusion in **non-multiplying system**, i.e., in system where fissile isotopes are missing and therefore the **fission cross-section is zero**. The neutrons are emitted by external neutron source. We will assume that the system is uniform outside the source, i.e. **D** and **Σ _{a }**are constants.

**neutron source**(with strength

**S**) as an

_{0}**infinite plane source**(in y-z plane geometry). In this geometry the flux varies so slowly in y and z allowing us to eliminate the y and z derivatives from ∇

^{2}. The flux is then a

**function of x only**, and therefore the

**Laplacian**and

**diffusion equation**(outside the source) can be written as (everywhere except x = 0):

For x > 0, this diffusion equation has **two possible solutions** **exp(x/L)** and **exp(-x/L)**, which give a general solution:

**Φ(x) = Aexp(x/L) + Cexp(-x/L)**

Note that, B is not usually used as a constant, because B is reserved for a **buckling parameter**. To determine the coefficients A and C we must apply boundary conditions.

From **finite flux condition** (0≤ Φ(x) < ∞), that required only reasonable values for the flux, it can be derived, that **A must be equal to zero**. The term exp(x/L) goes to ∞ as x ➝∞ and therefore cannot be part of a physically acceptable solution for x > 0. The physically acceptable solution for x > 0 must then be:

**Φ(x) = Ce ^{-x/L}**

where C is a constant that can be determined from source condition at x ➝0.

If S_{0} is the source strength per unit area of the plane, then the number of neutrons crossing outwards per unit area in the positive x-direction **must tend to S _{0} /2 as x ➝0**.

**neutron source**(with strength S

_{0}) as an

**isotropic point source**situated in spherical geometry. This point source is placed at the origin of coordinates. In order to solve the

**diffusion equation**, we have to replace the

**Laplacian**by its

**spherical form**:

We can replace the **3D Laplacian** by its **one-dimensional spherical form**, because there is no dependence on angle (whether polar or azimuthal). The source is assumed to be an **isotropic source** (there is the spherical symmetry). The flux is then a function of **radius – r** only, and therefore the diffusion equation (outside the source) can be written as (everywhere except r = 0):

If we make the **substitution** **Φ(r) = 1/r ψ(r)**, the equation simplifies to:

For r > 0, this differential equation has** two possible solutions** **exp(r/L)** and **exp(-r/L)**, which give a general solution:

Note that, B is not usually used as a constant, because B is reserved for a **buckling parameter**. To determine the coefficients A and C we must apply **boundary conditions**.

To find constants A and C we use the identical procedure as in the case of infinite planar source. From** finite flux condition** (0≤ Φ(r) < ∞), that required only reasonable values for the flux, it can be derived, that **A must be equal to zero**. The term exp(r/L)/r goes to ∞ as r ➝∞ and therefore cannot be part of a physically acceptable solution for r > 0. The physically acceptable solution for r > 0 must then be:

Φ(r) = Ce^{-r/L}/r

where C is a constant that can be determined from **source condition** at x ➝0.

If S_{0} is the source strength, then the number of neutrons crossing a sphere outwards in the positive r-direction **must tend to S _{0} as r ➝0**.

So that the solution may be written:

**neutron source**(with strength S

_{0}) as an isotropic line source situated in an infinite

**cylindrical geometry**. This line source is placed at r = 0. In order to solve the

**diffusion equation**, we have to replace the Laplacian by its

**cylindrical form**:

Since there is no dependence on angle Θ and z-coordinate, we can replace the 3D Laplacian by its **one-dimensional form**, and we can solve the problem only in **radial direction**. The source is assumed to be an isotropic source. Since the flux is a function of radius – r only, the diffusion equation (outside the source) can be written as (everywhere except r = 0):

This differential equation is called the **Bessel’s equation** and it is well known to mathematicians. In this case, the solutions to the Bessel’s equation are called the **modified Bessel functions** (or occasionally the **hyperbolic Bessel functions**) **of the first and second kind, **I_{α}(x) and K_{α}(x) respectively.

For r > 0, this differential equation has two possible solutions **I _{0}(r/L)** and

**K**, the modified Bessel functions of order zero, which give a general solution:

_{0}(r/L)To find constants A and C we use the identical procedure as in the case of infinite planar source. From** finite flux condition** (0≤ Φ(r) < ∞), that required only reasonable values for the flux, it can be derived, that A must be equal to zero. The term I_{0}(r/L) goes to ∞ as r ➝∞ and therefore cannot be part of a physically acceptable solution for r > 0. The physically acceptable solution for r > 0 must then be:

**Φ(r) = C.K _{0}(r/L)**

where C is a constant that can be determined from** source condition** at r ➝0.

If S_{0} is the source strength, then the number of neutrons crossing a cylinder outwards in the positive r-direction must tend to S_{0} as r ➝0.

So that the solution may be written:

**various zones**of different composition. The consequence of this is that the

**diffusion coefficient**,

**absorption macroscopic cross-section**, and therefore, the neutron flux distribution, will vary per zone. For the determination of the flux distribution in various zones, the

**diffusion equations**in zone 1 and zone 2 need to be solved:

where **a** is the real width of zone 1 and** b** the outer dimension of the diffusion environment including the **extrapolated distance**. With problems involving two different diffusion media, the interface boundary conditions play crucial role and must be satisfied:

*At interfaces between two different diffusion media (such as between the reactor core and the neutron reflector), on physical grounds the neutron flux and the normal component of the neutron current must be continuous. In other words, φ and J are not allowed to show a jump.*

1., 2. Interface Conditions

It must be added, as **J** must be continuous, the flux gradient will show a jump if the **diffusion coefficients** in both media differ from each other. Since the solution of these two diffusion equations requires four boundary conditions, we have to use two boundary conditions more.

3. Finite Flux Condition

The solution must be finite in those regions where the equation is valid, except perhaps at artificial singular points of a source distribution. This boundary condition can be written mathematically as:

The presence of the **neutron source** can be used as a boundary condition, because it is necessary that all neutrons flowing through bounding area of the source must come from the neutron source. This boundary condition depends on the source geometry and for planar sourve can be written mathematically as:

For x > 0, these diffusion equations have the following appropriate solutions:

**Φ _{1}(x) = A_{1}exp(x/L_{1}) + C_{1}exp(-x/L_{1})**

and

**Φ _{2}(x) = A_{2}exp(x/L_{2}) + C_{2}exp(-x/L_{2})**

where the four constants must be determined with use of the four boundary conditions. The typical neutron flux distribution in a simple two-region diffusion problem is shown at the picture below.

**Nuclear and Reactor Physics:**

- J. R. Lamarsh, Introduction to Nuclear Reactor Theory, 2nd ed., Addison-Wesley, Reading, MA (1983).
- J. R. Lamarsh, A. J. Baratta, Introduction to Nuclear Engineering, 3d ed., Prentice-Hall, 2001, ISBN: 0-201-82498-1.
- W. M. Stacey, Nuclear Reactor Physics, John Wiley & Sons, 2001, ISBN: 0- 471-39127-1.
- Glasstone, Sesonske. Nuclear Reactor Engineering: Reactor Systems Engineering, Springer; 4th edition, 1994, ISBN: 978-0412985317
- W.S.C. Williams. Nuclear and Particle Physics. Clarendon Press; 1 edition, 1991, ISBN: 978-0198520467
- G.R.Keepin. Physics of Nuclear Kinetics. Addison-Wesley Pub. Co; 1st edition, 1965
- Robert Reed Burn, Introduction to Nuclear Reactor Operation, 1988.
- U.S. Department of Energy, Nuclear Physics and Reactor Theory. DOE Fundamentals Handbook, Volume 1 and 2. January 1993.

**Advanced Reactor Physics:**

- K. O. Ott, W. A. Bezella, Introductory Nuclear Reactor Statics, American Nuclear Society, Revised edition (1989), 1989, ISBN: 0-894-48033-2.
- K. O. Ott, R. J. Neuhold, Introductory Nuclear Reactor Dynamics, American Nuclear Society, 1985, ISBN: 0-894-48029-4.
- D. L. Hetrick, Dynamics of Nuclear Reactors, American Nuclear Society, 1993, ISBN: 0-894-48453-2.
- E. E. Lewis, W. F. Miller, Computational Methods of Neutron Transport, American Nuclear Society, 1993, ISBN: 0-894-48452-4.