“CHAPTER 5” DIFFUSION THEORY

In many nuclear reactor physics texts (excluding Elmer Lewis’s recent text however), “Chapter 5” is dedicated to diffusion theory; hence, the title of this submission. Here, we will investigate analytical solutions to the most basic form of the monoenergetic 1D stationary diffusion equation. The intuitive approach taken radically departs from the usual method of solving the diffusion equation. In particular, we consider a general setting such that the method accommodates all solutions to the monoenergetic diffusion equations in 1D plane and curvilinear geometries. This is not your father’s diffusion theory and, for this reason, we anticipate it will eventually become the classroom standard.


With titles of The Diffusion of Neutrons Diffusion Theory: The Homogeneous One-velocity Reactor The Diffusion of Neutrons The One-Speed Diffusion Theory Model,
Chapters 5 in such well-known reactor physics texts of Edlund [1], Meghreblian and Holmes [2], Larmarsh [3] and Duderstadt and Hamilton [4] have been the primary source of our analytical solutions to the neutron diffusion equation since the 1950s. Little has changed since then, until recently. Our aim is to demonstrate a new paradigm where the monoenergetic neutron diffusion equations in 1D plane homogeneous and heterogeneous media as well as in curvilinear geometries and the multigroup approximation have a common origin enabling straightforwardly intuitive analytical representation. Here, we solidify the mathematics required to resolve the monoenergetic diffusion equation applied to classical 1D nuclear reactors. We derive a consistent theory for monoenergetic plane geometry for both multiplying and non-multiplying media alike. Multigroup solutions in heterogeneous media and in curvilinear geometries then immediately follow the pattern of 1D plane monoenergetic theory, but will not be considered here. The concept is simply to enable the independent solutions of the homogeneous diffusion equation to include boundary conditions to form the general solution. Once one knows the independent solutions then any particular solution follows through the method of Variation of Parameters (VoP). From matching of current at interfaces and specifying the flux and current at free surfaces, a three term recurrence results for the interfacial fluxes. Thus, we find the analytical solution as simply as the finite difference or finite element numerical solution but avoiding discretization error altogether. Our method provides nuclear engineering students with reinforcement of their study of second order ODEs including Wronskians, independent solutions to the homogeneous diffusion equation, particular solutions and solving 3-term recurrence. Arguably, our approach is the most efficient and consistent way of solving the neutron diffusion equation in 1D and, we therefore hope, will eventually become the educational standard adopted in all "Chapters 5" of future reactor physics texts.

SOLUTION TO THE 1-D-MONENERGETIC DIFFUSION EQUATION IN A HOMOGENEOUS MEDIUM
In the following, we consider the one-group characterization of a homogeneous reactor of width ∆ where a spatially dependent source is included for generality. With buckling (1) We leave the boundary conditions unspecified for the moment other than to assume they exist.
The general solution to Eq(1), written as the sum of the solutions to the homogeneous equation and a particular solution, is One uniformly observes, essentially in all reactor physics texts as well as in ODE math classes, the choice of two linearly independent exponentials for the solution to the homogeneous equation Introducing the known boundary conditions at x0 and x1 gives the following set of equations for the unknown coefficients a,b: to yield by Cramer's rule The simple re-arrangement, gives where the new independent solution is ( ) Thus, the solution automatically satisfies the boundary conditions, at least for the homogeneous ODE.
Once we have the solutions to the homogeneous ODE, VoP then provides the particular solution as where the Wronskian is and is non-zero for all buckling. Since the general solution becomes We could also write the solution as which is an alternative of many possible particular solutions.
Several principles of solving ODEs are in play here. Specifically, when solving an ODE, one has many practical choices for the two linearly independent solutions with some more convenient than others. As demonstrated, it is convenient to choose functions that automatically satisfy the boundary conditions as identities. Moreover, we choose the particular solution in the same way according to Eqs(4a) and (6b).
The advantage of Eq(5b), called the customized solution, is clear for the homogeneous medium and will become clearer for the heterogeneous medium. To solve any 1D problem in plane geometry, one need only assume Eq(5b) rather than the exponential solution, thus avoid solving for the coefficients a and b in the complementary solution, Eq(3a). No matrix algebra is necessary for the single medium with the solution entirely transparent. In addition, note that nothing was said about whether the medium was fissioning or not. The procedure as presented thus far should be equally applicable. If fissioning or not, the sign of B 2 can be either negative or positive. The case presented tacitly assumed B 2 was positive. If however, B 2 is negative, the sine switches branches and becomes the hyperbolic sine (sinh) and we write

SOLUTION TO THE 1-D-ONE-GROUP DIFFUSION EQUATION IN A HETEROGENEOUS MEDIUM
We now consider the one-group characterization of the heterogeneous reactor configuration of n slabs as shown here x-1 x0 x1 x2 xj-1 xj xn xn+1 Fig.2. Heterogeneous diffusing medium.

Now the diffusion equation becomes
We also require flux and current continuity at the interior slab interfaces The analysis begins with the development of the region wise solution to Eqs(8) in slab j where the interfacial fluxes are where for the two branches of the sine function are As above, the first two terms of Eq(9a) are the homogenous solutions and the third is the particular solution.
The solution automatically satisfies flux continuity at slab interfaces since Again, when we apply VoP, the following particular solution emerges: The advantage of the customized form of the solution is apparent when we satisfy the interfacial conditions of Eqs(8c,d) to give the following three-term recurrence relation for the coupling coefficients between two adjacent regions:

Including Boundary Conditions
Assume vacuum on the boundaries, which is defined as no incoming partial current, or [3] ( ) ( ) ( ) ( ) Assuming no sources in the first and last slabs using Eq(9a) with j=1 and n gives [ ] Thus, the extrapolation distance, to be determined iteratively, can be different at each surface. The physical boundaries now become x0 -d0 and xn + dn, where the flux vanishes.

AN EXAMPLE: ENTERING CURRENT AT x = 0, xn
Partial currents J+ and J-enter a heterogeneous slab, devoid of sources, at the left and right boundaries. The incoming currents at the boundaries are from Eqs(13a,b) As above, using Eq(9a) with j = 1 and n at x0 and xn respectively gives ( ) n n n n n n n n n n n n n n to be solved along with Eq(12a) using several possible methods--including transformation of the boundary value recurrence into two initial value recurrences, tridiagonal matrix inversion by the Thompson algorithm, matrix diagonalization. As a demonstration, a 101 slab heterogeneous medium, with J+, J-= 1, of alternating It is interesting to observe that both the analytical solution and FD approximation derive from the same recurrence relation and therefore represent essentially identical numerical effort. The question is therefore--Why would one use finite difference when the analytical is as easily obtained?

CONCLUSIONS
We have presented a new approach to analytically solving the 1D neutron diffusion equation. The power of the new approach is the two-slab connection between multiple slabs in a heterogeneous medium to give a readily solvable tridiagonal matrix for interfacial fluxes.
Unfortunately, space limitation precludes discussion of straightforward extensions to curvilinear and multigroup diffusion theory other than to write the following: + For curvilinear geometries: From the plane/point transformation, spherical geometry follows the identical form presented in plane geometry. Cylindrical geometry follows using Bessel functions as the independent solutions, but since Bessel functions do not combine like sines and cosines of different argument, an additional step forcing the linearly independent solutions to give the appropriate boundary conditions must be included. On a final note, for a homogeneous medium assuming medium width h is small, the coefficients of the recurrence relation can be expanded in a Taylor series to give for Eq(12a) ( )