Collision probability calculations of multiplication factors for lattices

Journal of The Franklin Institute DEVOTED

TO SCIENCE AND THE MECHANIC

ARTS

June

Volume 291, Number 6

1971

Collision Probability Calculations of Multiplication Factorsfor Lattices by

DOUGLAS

C. HUNT

Dow Chemical Company Rocky Flats Division, Golden, Colorado A mu~tigroup

ABSTRACT: mine

the multiplication

Jissile components. long cylinders) arranged

0~’

the

slabs).

ma,y, however,

in a linear, planar

Lattices

computation

u.niform neutron-source density

energy

transport The lattices similar

method and,

neutron

(cubes or

compon.ents

component

collisions

in the lattice,

in each lattice component.

in a 16-group

is speci$ed predicts

satisfactorily

thermal expansion

in$nitely

components

The lattice

the permissible

60 deter-

(i.e. spheres,

one-dimensional treated.

U&d

may

be

arrangements

starting

from

which is obtained from

format

a

The initial sourceS,

solutions.

to this degree

lattices.

iS

lattices of identical

geometry.

dens&y distribution

distribution

equation

generatiOn

must be one-dimensional of nearly

or cuboidal manner;

traces successive

SUCCSSSiVt?

and unreJlected periodic

be approximately

being limited by the lattice component The

Of

??%%Od

of unmoderated

The lattice components

or infinite

long cylinders

Ve%iOlt

factors

the critical

of conjklence,

The method

pitch

of

is used to compute

two experimentally

is also used to fina! the reactivity

for several just-critical

for large arrays of small elements

arrays.

This

measured

the critical parameters temperature

coeffkient

which have a high interaction

is found

of other

coeficient

of

to be positive

factor.

I. Introduction

The critical parameters of finite arrays and lattices of fissile material are generally determined using either sophisticated computer codes (such as the 05R Code) (1) or simple scaling techniques, such as the density analogue method (2). In this paper, an approximate calculational technique is outlined which, subject to certain restrictions, offers a simple and accurate method for determining the multiplication factor of finite periodic arrays. The restrictions are that the average path length for neutron travel in the arrays before escape should be less than about two transport mean-free paths and that the probability of scattering collisions is much less than unity. Ways of modifying the above restrictions are discussed in Ref. (5).

417

Douglas C. Hunt The calculations are done using the collision probability method (3-6). Briefly, in this method, successive neutron collisions are traced in the lattice starting from a uniform neutron-source density distribution in each lattice component and assuming each lattice element has the same neutron content. This initial spatial source distribution is assumed to persist in all successive neutron collisions. The initial source-density energy distribution is specified in a 16-group format obtained from S, transport equation solutions. Critical points are found in the usual manner by determining the multiplication factor as a function of the system property of interest and graphically determining where a multiplication factor of unity occurs. The method is tested by determining the critical spacing of 2 by 2 by 2 and 3 by 3 by 3 lattices measured by Thomas (7). Having thus established a degree of confidence for the method, other critical spacing calculations are made varying the number of array elements and the size of the array elements. As a further application of the method, the expansion reactivity temperature coefficient of several jzcst critical arrays are computed. II.

Calculation&

Multiplication

Procedures

Factor Expression

The multiplication by Hunt (3). Coll&on

factor expression and the associated notation are given

Probability Expressions

The array collision probabilities used in the calculations are obtained from the following expression [derived by Rothenstein (S)]

(1) In Eq. (l), the 1 is the mean chord length of an array element or four times its volume-to-surface ratio, Q is a constant whose value depends on the array element geometry and on its blackness, and C is the array interaction factor Dancoff and Ginsburg (9) suitably summed and averaged over all possible visible neighbor interactions in the array. Rothenstein gives values of Q for arrays of one-dimensional elements (spheres, infinitely long cylinders and two dimensionally infinite slabs) in the limits as ZZTRj1-+O and YlTR,i E+ co. The values of Q for the white limit (CTRj Z+ 0) are used in all computations. The white limit values slightly overestimate collision probabilities, but generally give good agreement with exact values for the entire range of array element blackness. This result is expected since uniform source density-collision probability calculations are only valid in the white limit. Equation (1) assumes an infinite periodic lattice of fuel elements immersed in an infinite moderator bath. It is, however, valid for finite lattices when C is interpreted as the finite lattice interaction factor (10). Equation (1) is derived by replacing averages of products by products of averages in the exact expression for array collision probabilities and in

418

Journal of The Franklin

Institute

Collision Probability Calculations of Multiplication

Factors for Lattices

addition

assumes that the Dancoff factor is small compared to unity. An expression for array collision probabilities has been obtained in a similar manner by the author by replacing averages of products by products of averages, and assuming a uniform and isotropic interaction flux between array elements. This expression is

alternate

where P$ is the collision probability for neutrons in the array element in which they originated and P$ P$, . . . are collision probabilities for neutrons in unshielded neighbors, singly shielded neighbors, etc.* Thus, e,,i and P;j are given by the following equations, with similar expressions for multiply shielded collision probabilities : P$ = c, l&J

I( 1 - P,oj)

and P$ = c,qQqj

l(1 -Pz,j,

(1 -P&J.

(3)

(4)

In Eqs. (3) and (4), C, and C,, are Dancoff factors for unshielded and singly shielded neighbors, respectively. The advantage of Eq. (2) is that the exact values of the self-collision probabilities (PE,i) may be used and also that no restriction on the magnitude of the Dancoff correction factors is present. The advantage of Eq. (1) is the simple manner in which the array element dimensions and cross-sections appear. Calculation of Interaction Factors The array interaction factor C in Eq. (l), which is identically equal to the C, of Eq. (3), is computed in terms of C,, the interaction factor for pairs of lattice elements. The net interaction factor is then obtained from C =

z8 f,CS.

(5)

The s is an array element symmetry index and f, is the fraction of the array elements with symmetry s. The values of C, are computed according to

The C,,$ is the interaction factor at center-to-center element spacing i, and nS,irefers to the number of elements at the ith spacing from an array element having symmetry s. Equations (5) and (6) assume array elements are spaced sufficiently far apart so that partial shielding of one element by another is not important. Expansion

Reactivity Temperature Coeficient for an Array

On the assumption that reactivity changes only through density changes of the array elements, the multiplication factor expression may (see Ref. * Unshielded, singly shielded, those which have zero, one, two, the neighbor element.

Vol.291,No. 6,June 1971

doubly shielded, . . . neighbor array elements are . . . array elements between the source element and

419

Douglas C. Hunt 3) be differentiated

with respect to temperature

to give

(7) where

Using Eq. (1) for P,j [ = Pc(&R,3x)], dP.JdT is given as dP.

-g =

and denoting

C,R,il

by pi, then

1

2Q (QP~+~)~ '

CL&(1- Pc,J2

(9)

In Eq. (9), 01is the linear temperature expansion coefficient for the array elements and x is a characteristic array element dimension. If the array elements are spheres, the variation of C, with the ratio of sphere radius (rs) to center-to-center sphere spacing (do may be quite accurately approximated by the cubic equation (11):

Cc,,{=

A(z)++j2+D($,

A = 0.0017, Equation

B = 0.207s,

D = 0.2032.

(9) may then be written for arrays of spherical elements as

(11) In the above Eq. (1 l), g = ~fJw

gs = ~%,A%i

and

go,<= B(z)“+2D(-$.

(12)

The terms in Eq. (12) are defined as in Eqs. (5), (6) and (10). Extrapolation

Procedures

The limit indicated in the multiplication factor according to the procedure described by Hunt (3).

expression

is obtained

Scaling Techniques The experimental data (7), which are compared with the predictions of the present method, pertain to arrays with cylindrical elements. Since the present method is directly applicable only to arrays of one-dimensional lattice elements, the Thomas arrays must be suitably scaled before a comparison is possible. According to Eqs. (2), (3) and (4), the following equivalence relations must be satisfied to scale one-lattice element geometry to

420

Journa.l of

The Franklin

Institute

Collision Probability Calculations of Afultiplication

Factors for Lattices

another and maintain the same &, value: P$ = constant, C,l = constant, C,, 1 = constant.

(13)

Satisfying the first condition in Eq. (13) is the most important, since the terms involving the C, and C, in Eq. (2) (the P$ and P$) are significantly smaller than the P$. Thus, the scaling method is to maintain a constant element spacing and to require that P$ is constant for all energy groups. The Thomas arrays are accordingly scaled to equivalent spherical element arrays, using the collision probability values of Carlvik (12) for cylinders and those of Case et al. (15) for spheres. It develops that this scaling technique also approximately satisfies the conditions for the constancy of CVl and C,,l required by Eq. (13). ZZZ. Results

Comparison of Theory and Experiment The scaling procedure discussed in the previous section is applied to the 2 by 2 by 2 and 3 by 3 by 3 Thomas arrays. The dimensions of the Thomas array cylindrical elements equal a height of 10.77 cm and a diameter of 11-519 cm. Each element is composed of enriched uranium metal (235U = 93 per cent and YJ = 7 per cent) having a density of 18.76 g/cm3. The radius of the equivalent spherical lattice element is 6.061 cm. The comparison of the pitches* of the computed equivalent lattices and of the corresponding experimental lattices is shown in Table I. The method appears to somewhat TABLE

I

Calculated Results an,d Comparison with Experimental

Geometry*

2X2X2 array? 3X3X3 array$.

Computed dimension

Minimum center-tocenter element spacing

and Theoretical

Values

Computed critical dimension

Measured critical dimension

(cm)

(cm)

Computed critical dimension (Monte Carlo)

13.88

713.84-f

13.97 (14)

18.42

‘18.07t

18.54 (14)

* The lattice elements are uranium metal and have a density of 18.76 g/cm3. The uranium contains 93 per cent 235U and 7 per cent zaeU. t The array elements are 6.061 cm radius spheres. $ Density and array element geometry scalings are applied to the measurements.

overpredict the observed pitches. This is due to the white-limit assumption used to compute the array collision probabilities according to Eq. (1). * The pitch of a periodic lattice is the minimum spacing.

Vol.291,No.6,June1971

center-to-center

lattice element

421

Douglas C. Hunt Lattice Pitch vs. Number and Size of Lattice Elements The results of the calculations of critical lattice pitch are given in Figs. l-4. In Fig. 1, the nature of the convergence of c&1 l~~&,,~~ to keff is shown for several typical calculations. In general, the value of keff is about 15-20 per cent greater than the ;rfz1 kn,eif value for the maximum N considered in a computation. In Fig. 2, the graphical solution of Eq. (1) for the critical (ke.* = 1) lattice pitch is shown. It is to be noted that keff varies nearly linearly with lattice pitch. Figures 3-4 give the variation of the computed I

I

I

I

Lattice ,pitcc,

(5by5

by

La&e (2 by 2 by 2) Sphere radius

v

0 FIG.

1.

I

I

I

2

I

I

3 4 Number of collisions

I

I

I

I

5 , N

6

7

8

Convergence of the multiplication factor for N collisions to the system multiplication factor.

60

0

0.86

0.92

0.96 System

I.00

multiplication

FIU. 2. Graphical determination

I.04

1’08

I-12

factor (k,tt)

of critical lattice pitch.

values of critical lattice pitch with the number of lattice elements in a cubic latticeandwith theradiusof alattice element.InFig. 3, thedataarerepresented by plotting lattice pitch vs. number of lattice elements parametric in lattice element radius. In Fig. 4, the lattice pitch is plotted vs. lattice element radius, parametric in the number of lattice elements. For this representation, the limit of the lattice element radius as the lattice pitch tends to infinity is the bare sphere critical radius. Also, shown in Fig. 4 is the line representing

Journal of The Franklin

Institute

Collision Probability Calculations of Multiplication

Factors for Lattices

lattices of elements in contact. Along this line, critical points are computed using the density analog method. As would be expected, for lattices with a large number of elements, the density analog points agree well with the collision probability predictions.

IO -

(4.01”

-

0

I 20

,

X=Experimenlal

,

I

I

40

,

60 Number-of

I

points ,

80 lattice

I

loo

I

I,

120

140

elements

FIG. 3. Critical pitch vs. number of lattice elements for cubic lattices. Ice 90 I30 g

70

; u -5

50

60

g ‘5

40

4

30

(5bv5 bv5)-

20 IO

Radius

of

lattice

element

, cm

FIG. 4. Lattice pitch vs. lattice element radius for critical cubic lattices.

Temperature CoefJicient Calculations In Figs. 5, 6 and 7, the results of the temperature coeflicient calculations are presented. The calculations are all made for lattices which are just critical. In Fig. 5, the trend toward the convergence of ~~=‘=,dlc,,&dT to dk,,,/dT [see Eq. (7)] is shown for several of the calculations. The curves show less curvature than those of C~=‘=,/C~,,~~ vs. N; the dlc,,,/dT values being generally about 50 per cent higher than the ~~=‘=,dlc,,,JdT for the maximum N calculated. In Fig. 6, the values of dk,,,ldT for 5 x 5 x 5 lattices are shown plotted vs. lattice element radius and vs. the ratio of lattice element radius to lattice pitch. The temperature coefficient of reactivity changes from negative to positive as the lattice interaction factor increases (the curve of lattice interaction

Vol. 291, No. 0, June 1971

423

Douglas C. Hunt

12 -

All Curves (5 by 5‘ by 5 lattice)

Number

of cdlisions,

N

FIG. 5. Convergence of reactivity temperature coefficient of expansion for N collisions to system reactivity temperature coefficient of expansion. Radius

of lattice

cm

Arrows indicate axis from which curves should be read

15v;‘ 0 -z

,

element

-

u

IO-

-1-5 -

Lattice

element

radius/lattice

pitch

reactivity temperature coefficient of expansion vs. element radius and the ratio of element radius to lattice pitch for critical lattices (5 by 5 by 5) ; also lattice element radius vs. lattice interaction factor for critical lattices (5 by 5 by 5). FIG.

6. System

_,J/o 0

FIU.

424

7. System reactivity

20

40 Number

60

of

60 lattice

100

I20

elements

temperature coefficient of expansion vs. number of lattice elements for critical cubic lattices.

Journal of The Franklin Institute

Collision Probability Calculations of Multiplication

Factors for Lattices

factor vs. lattice element radius also appears on the figure). This effect, which is important from a nuclear safety standpoint, can be easily explained. As the temperature rises, the lattice element density decreases. The density decrease in each lattice element causes more neutron leakage from the lattice, thus decreasing the lattice reactivity. At the same time, the temperature increase causes the lattice elements to increase in size, thus increasing the interaction probability and making the lattice more reactive. For critical lattices in which the interaction factor is relatively large, the second effect predominates causing a positive reactivity temperature coefficient. Finally, in Fig. 7, d&.JdT is plotted vs. the number of lattice elements for 6.061 cm radius lattice elements. IV.

Discussion

The collision probability method as described has three primary sources of uncertainty. First, an ampproximate expression is used to compute array collision probabilities, and, second, the exact collision probability expression assumes that the same uniform and isotropic neutron density is present in all array elements. Lastly, the multiplication factor values are obtained by an extrapolation technique which must entail some error. The amount of error due to the approximate collision probability expression may be estimated by comparison with the results of Rothenstein (8). Rothenstein computes escape probabilities for hexagonal lattices of long uranium-238 rods, immersed in water, using a Monte Carlo code and assuming an initial uniform neutron-source density distribution. His results are given in Table II. The predictions of Eq. (1) are seen to be acceptable over a wide range of element blackness and for two diverse (and relatively large) interaction factor values. TABLE Comparison

0.5 o-5 3.0 3.0 6.0 8.0

1)

Carlo Collision

Probability Calculations

(Eq. 5, & = +)

Escape probability (Monte Carlo)

Difference (per cent)

Interaction factor

0.559 0.672 0.164 O-233 0.087 0.096

0.572 O-691 o-173 0.250 0.088 0.100

2-2 2.7 5.2 5.2 1.1 4.0

0.46 0.21 O-46 0.21 O-46 0.21

Escape probability mm,,

II

of Results from Eq. 5 and Monte

The uncertainty due to the uniform flux approximation is discussed by Harper (13). He concludes that the uniform flux approximation is valid if the neutrons, on the average, traverse less than two mean-free paths before escape from a system and if the probability of a scattering collision is small.

Douglas C. Hunt The extrapolation technique uncertainty is estimated by making sucoessive three-point fits to the multiplication factor expression at N = 1, 2, 3 ; 2, 3, 4 ; and 3, 4 and 5. It is estimated that the inclusion of more collisions in the calculations would reduce the predicted critical spacings by 1 per cent. It is concluded that the approximations involved in the method ‘should not seriously affect the computations. Hence, the results of the calculations, and particularly the interesting positive temperature coefficient prediction, are considered valid.

(1) D. C. Irving, R. M. Freestone, Jr. and R. B. K. Kern, )“OBR, A general purpose Monte Carlo neutron transport code’!, ORNL-3622, Union Carbide Corp., Oak ’ Ridge National Lab., Feb. -1965. (2) C. L. Schuske, “Method for calculating uranium-235 metal storage arreys”, RFP-108, The Dow Chemical Co., Rocky Flats Division, May 7, 1958. (3) D. C. Hunt, “Collision probability criticality calculations”, J. Franklin Inst., Vol. 289, p. 175, 1970. (4) D. C. Hunt, “Collision probability multiplication factor calculations”, Trans. Am. Nut. Sot., Vol. 12, p. 212, June 1969. (5) D. C. Hunt, “Collision probability calculations of multiplication factors for lattices”, RFP-1186, The Dow Chemical Co., Rocky Flats Division, Jan. 1969. (6) T. Lefvert, “A multiple scattering collision probability method for neutron transport problems”, submitted to Nucl. Sci. Engng, 1970. (7) J. T. Thomas, “Neutron physics division annual progress report”, ORNL-3499, Vol. I, Oak Ridge National Lab., Dec. 10, 1963. (8) W. Rothenstein, “Collision probabilities and resonance integrals for lattices”, Nucl. Sci. Engng, Vol. 7, pp. 162-171, 1960. (9) S. M. Dsncoff and M. Ginsburg, “Surface resonance absorption in a close-packed lattice”, CP-2157, Univ. Metallurgical Lab., Oct. 11, 1944. (10) W. Rothenstein, private communication, June 1969. (11)H. K. Clark, “Comparison of a simple theoretical treatment of critical arrays of fissionable units with experiment”, DP-898, Savannah River Lab., Jan. 1964. (12) I. Carlvik, “Collision probabilities for finite cylinders and cuboids”, Nucl. Sci. Engng, Vol. 30, pp. 156-152, 1967. (13) R. C. Harper, “The collision probability code-MINOS”, J. N&. Energy, Vol. 21, pp. 767-785, 1967. (14) J. T. Michalczo, “Multiplication factor of uranium mete1 by one-velocity Monte Carlo calculations”, NW& Sci. Engng, Vol. 27, pp. 567-663, 1967. (15) K. Case, F. de Hoffmann end G. Placzek, “Introduction to the theory of neutron diffusion”, Los Alamos Scientific Lab., June 1953.

426

.Journal

of The Franklin

Institute