A boundary element method is developed for the 1-D nodal diffusion equation in cylindrical geometry. This method retains the matrix qualities of the nodal formulation while providing an accurate computational benchmark for evaluating reactor analysis codes

Ganapol, B. D.

Ougouag, Abderrafi M.

UNT Digital Library

This is a preprint of a paper intended for publication in a journal or proceedings. Since changes may be made before publication, this preprint should not be cited or reproduced without permission of the author. This document was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States Government nor any agency thereof, or any of their employees, makes any warranty, expressed or implied, or assumes any legal liability or responsibility for any third party’s use, or the results of such use, of any information, apparatus, product or process disclosed in this report, or represents that its use by such third party would not infringe privately owned rights. The views expressed in this paper are not necessarily those of the United States Government or the sponsoring agency. INL/CON-06-11443PREPRINTThe Boundary Element Formulation of the One Group 1D Nodal EquationsAmerican Nuclear Society Annual MeetingB.D. Ganapol A.M. Ougouag June 2006 The Boundary Element Formulation of the One Group 1D Nodal Equations B.D Ganapol+Department of Aerospace and Mechanical Engineering, University of Arizona, Tucson, Arizona, 85721 A.M. Ougouag Reactor and Nuclear Physics Idaho National Laboratory, PO Box 1625, Idaho Falls, Idaho 83415-3885 The nodal method is a powerful numerical tool with which to treat neutron diffusion in heterogeneous media.1,2  The popularity of nodal methods stems from the common features they share with finite element and finite difference methods.  In particular, like finite element methods, in nodal methods, the flux is approximated over a coarse grid, possibly by piecewise continuous functions to provide an accurate representation of the flux and/or current map.  Like finite difference approximations, the coupling matrix that results from the application of the nodal approach is relatively sparse and can usually be treated by routine numerical inversion.  In addition, nodal methods are compatible with modern homogenization techniques that characterize fuel pin geometry in an efficient and accurate manner.  Finally, nodal methods (especially those of the analytical variety) ensure neutron balance and therefore have a degree of robustness not found in finite difference techniques. There are many derivations of nodal methods that can be found in the literature1 that relate the flux and/or current within a cell to values on the cell boundaries.  Here, we will again consider the 1D steady state neutron diffusion equation through the application of the boundary element method applied to a single cell.  The advantage of this approach is its succinctness and completeness, simultaneously yielding all possible analytical nodal equations in a 1D heterogeneous medium setting.  Also the approach is extendable to multigroup formulations and multidimensions which is to be the subject of future investigations.  While the solutions obtained are well known, what is new is the procedure from which they are derived providing a variety of solutions.  This effort represents a first step in the development of a highly accurate multidimensional benchmarking tool that will be used to assess multidimensional + Work performed under INL contract no. 00044817Pebble Bed Reactor (PBR) production codes.3To this end, analytical nodal formulation will be compared to a proposed 1D finite difference (FD) scheme. Reformulation of the Diffusion Equation We consider the one group, 1D diffusion equation in the jth region (thickness j' ) of a heterogeneous reactor, where the jth region means the interval (xj-1, xj):               222 j jjjjdB xdxS xq xD )     ª º« »¬ ¼       (1)             2 ajfjjjBDQ6  6{ª º« »¬ ¼with interfacial continuity of flux and current             1 1 1 1j j jj j j j jxx x   ) { )) { )  )           1 1 1 1 .j j jj j j j jJ J xJ J x J x   {{  A zero net current is assumed at the reactor centerline as is a zero flux at the outer surface. The derivation begins with the extension of the flux in region j to all space through                   j j jx x xU< { )                 1j j jx x x x xU { 4  4 where jU  is the characteristic function that ensures the flux outside the specific slab of interest is zero, since it does not concern us as long as the boundary fluxes are properly set.  Ĭis recognized as the Heavyside step function.  When a second spatial derivative is taken of j< ,and Eq.(1) is taken into account, there results Transport Methods: General 521              222112.j jj jjj jjj jdB xdxx x x x xdx x x x xdxx q xG GG GU <  c c   )     ) ª º« »¬ ¼ª º¬ ¼ª º¬ ¼ (2) Thus, artificial boundary sources have been introduced consisting of yet-to-be- determined boundary fluxes with their presence ensured by Dirac-delta functions and their derivatives.  In this way, particle balance is rigorously maintained.  It is this equation that is new. The Green’s Function and Solution The Green’s function, satisfying          222,j jdB G x x x xdxGc c   ª º« »¬ ¼,can be shown to be simply (using Fourier transforms)           1, sin2j jjG x x B x xBc c   . (3) When the Green’s function is integrated against Eq.(2) over all space, we find the solution                   1 11 1 1sin1 12 sincos sgn12 cos sgnjj j j jj j j j j jj j j j jj j j j jjxJ x B x xD B J x B x xx B x x x xx B x x x xQ x       )   )  <ª º« »¬ ¼ª º« »¬ ¼(4) where      1,x jj j jx jQ x dx q x G x xc c c{ ³ .Useful relationships between the boundary flux and current are obtained by introducing 1jx x H   and jx x H   into Eq.(4) to obtain            111sincos 2j j jjj jj j j j jJ BD BB Q x ' ) ' )   (5a)     111sincos 2 .j j j jj jj j j j jJ BD BB Q x  ' ) ' )              (5b) By restricting x to the interval 1,j jx xª º¬ ¼ , Eqs. (4) and (5) lead to the following relationships for --the interior flux to boundary fluxes             > @   > @    1 1j j j j jj j j j jx Q x xQ x x Q x   &   & ) ))        (6a) --the interior current to the boundary currents                   1 1j j j j cj j jj cj j j cjx D B Q x xJ Q x x Q xJ J    &   & ª º¬ ¼ª º¬ ¼   (6b) with               1,sgnxjjcj jxjdG x xQ x dxq x x xdxcc c c{ ³ .and                 1.sinsinsinj jjjj jjjj j jB x xxSB x xxSS B& {& {{ 'ª º« »¬ ¼ª º« »¬ ¼These are just two of the possible eight relationships that can be obtained from Eq.(4).  The solution is finalized, say for Eq.(6a), when the boundary fluxes are found recursively by requiring current continuity           1 2j j j jj jb fJ    ) ) )   (7) where the coefficients are known in terms of the nuclear parameters of regions j and j-1. The numerical solution to Eq.(7) is most effectively obtained using a tridiagonal solver.  Although a solution to Eq.(2) has previously been obtained through approaches based on the repeated use of finite domain Green’s functions for each homogeneous sub-domain, it is noteworthy that the present work obtains the relations directly from a single Green’s function for the entire heterogeneous domain.  Once the collected ensemble of all types of Eq.(6) for all cases are solved, the remainder of the method presented here allows the determination of fluxes (and Transport Methods: General522currents if desired) throughout the domain [per Eqs.(5a) and (5b)]. To emphasize the novelty of the development, it should be noted that the boundary element methods prescribes the most efficient form of the homogeneous solution when given two boundary conditions on either the flux or current.  In addition, a direct extension to multigroup is possible. A Benchmarking Application In this section, an application of the above analytical result to assess the performance of a proposed FD approximation is presented. Consider a grid of Ln intervals (slabs) representing a particular subdivision of a heterogeneous medium where the different material boundaries are maintained; however, the contiguous material regions can be further subdivided into subregions. Also let each material thickness be removed from the diffusion equation by the transformation / jx xo ' .  If Eq.(1) is integrated over the jth interval assuming a uniform source, then              12j j jjj j jh jJ x J xD B dx x hS   )  ³  (8) with 1 1 /j j Nh x x   {  .N is the number of discrete nodes per region (either fuel or water).  If a trapezoidal rule integration is assumed and if   1j jj j jJ x D h) )ª º « »¬ ¼=  ;consequently, since               1 21 1 1 1 j jj j j j jJ x J x D h     ) )ª º« »¬ ¼=  ,we find a three-point (finite-difference) recurrence for Eq.(8) in terms of cell edge fluxes  1 22 212 21 / / 21 / 2j j j j j jj j jjjx b fD D B hbB hJ c c c    c { ) ) )12 2/1 / 2j jjjD DB hJ c {22 2.1 / 2jjjh qfB hc { For zero current at the centerline, the recurrence for j = 2 can be shown to be as above except with 2 2 2 1 22 2 1 2.    0 b b bf f fJ JJc c c co  oc c c co As an example, consider a heterogeneous half symmetry subcritical MTR reactor of slabs of alternating water and BWR fuel of thicknesses 1.1158 and 3.231cm respectively. A uniform source exists in each water channel to simulate the thermalized fast flux.  The flux variation from Eq.(6a) is shown in Fig. 1 for Ln = 30, 60, 120 slabs.  The stability and fidelity of the analytical solver is evident as seen from the repetitive nature of the flux approaching an infinite heterogeneous MTR reactor. Figure 2 shows the relative error of the FD scheme in comparison to the analytical flux x0 50 100 150 200 250 300Flux0.00.51.01.52.02.53.0Ln = 30Ln = 60Ln = 120           Fig. 1. MTR reactor with increasing                       number of fuel cells. for increasing discretization.  For the largest number of discrete nodes, the FD scheme can deliver accuracy better than 10-3 which is an indication of a robust FD scheme.  The motivation for this comparison is the 2D application.  It is anticipated to eventually be able to use a 2D finite difference formulation with convergence acceleration as a legitimate benchmark—but you will have to wait for this one since convergence acceleration will first be applied to the above 1D scheme to assess its accuracy.Transport Methods: General 523X0 10 20 30 40 50 60 70Relative Error1e-71e-61e-51e-41e-31e-21e-1N= 2N= 8N= 32N= 50Ln = 30      Fig. 2. Relative error comparison. References [1] R. D. Lawrence, Progress in Nuclear Energy, Vol. 17, No. 3, 271-301, (1986). [2] R. J. J. Stammler and M. J. Abate, Methods of Steady-State Reactor Physics in Nuclear Design, (Academic Press), Chapter XI, “Nodal Equations for 3D Reactor Calculations”, (1983).  [3] Modular Pebble Bed reactor Project University Research Consortium Annual Report, INEEL/EXT-2000-10034, July 2000. Transport Methods: General524

The Boundary Element Formulation of the 1-Group, 1-D Nodal Equations

https://digital.library.unt.edu/ark:/67531/metadc889290/m2/1/high_res_d/911668.pdf

The Boundary Element Formulation of the 1-Group, 1-D Nodal Equations

Abstract

Similar works

Full text

Available Versions

UNT Digital Library