One-speed transport disadvantage factor calculations for general anisotropic scattering

Journal of Nuclear Energy. Vol. 24, pp. 23 to 34.

Pergamon Press 1970.

Printed in Northern Ireland

ONE-SPEED TRANSPORT DISADVANTAGE FACTOR CALCULATIONS FOR GENERAL ANISOTROPIC SCATTERING G. W. ECCLESTON,JR. and N. J. MCCORMICK Department

of Nuclear Engineering, University of Washington, Seattle, Wash. 98105, U.S.A. (Received 22 September 1969)

Abstract-The singular eigenfunction expansion method has been used to obtain the solution for the onespeed disadvantage factor in slab cells of fuel and moderator which may scatter neutrons in a highly anisotropic manner. New two-media cross-product integral relations, together with previously known one-medium orthogonality relations, have been utilized to obtain the solution. Sample calculations on slab cells of natural uranium and water show that, for cases in which the anisotropy of the scattering is quadratic or more, the disadvantage factor differs from the isotropic result by 7-27 per cent depending upon the order of the anisotropy, the effective mass of the scatterer, and the slab cell width. In all calculations, scattering in the fuel has been assumed to be isotropic. 1. INTRODUCTION

A MEASURE of the differences in the one-speed neutron fluxes in the fuel and moderator of an infinite slab lattice is provided by the disadvantage factor, 5. This factor is defined as the ratio of the space and angle-averaged flux in the moderator to that in the fuel, lb

r

-d I a dx s -1 i=

1 I

ra

yz(x,

,4 +

,

fl

dx a 0 -J J y,(x, -1

(1)

P) dp

where the subscripts 1 and 2 imply fuel and moderator, respectively, and “a” and “A” are their optical half-widths as shown in Fig. 1. FERZIGERand ROBINSON(1965) utilized the method of singular eigenfunction expansions developed by CUE (1960) to rigorously calculate 5 for isotropically scattering slab cells. BONDand SIEWERT(1969) later extended the method in order to calculate the effect upon 5 of a two-term scattering kernel of the moderator. The calculational procedure used here extends the technique of Bond and Siewert to include the effect upon 5 of a general anisotropic scattering kernel in both the fuel and moderator regions; consequently, a quantitative understanding of the neglect of higher order scattering terms may be obtained in sample calculations. The necessary theoretical development is based upon new two-media full-range cross-product integral relations developed in Appendix A. These relations are generalizations of those first derived by KUSZELL(1961). It should be pointed out that alternate procedures to the eigenfunction method might have been used to account for the effects of anisotropic scattering upon 5. For example, integral transport methods such as the collision probability method or the discrete method of CARLVIK(1967) can take into account the effects of anisotropic scattering. Also the invariant imbedding method (MINGLE, 1969) can be used. VANMAS~ENHOVE and GROSJEAN (1965) have also developed an exact analytical solution to the infinite-medium Green’s function problem for arbitrary anisotropic scattering. 23

G. W. ECCLESTON,JR. and N. J. MCCORMICK

...

Fuel, / medium 1 j itu+ I 0

Moderator, ) medium 2 ; +---h-t!

. ..

I-x b

a

FIN. 1.-Slab cell geometry.

The complexity of the method, however, appears to make it unattractive for multimedia problems. In Section 2, the basic formalism of the eigenfunction method is presented. Section 3 contains the generalization of the procedure of Bond and Siewert for obtaining the disadvantage factor by utilizing full-range, two-media cross-product integral relations to determine the unknown expansion coefficients of the eigenfunction expansion. Disadvantage factors obtained for typical cells of THEYS (1960) are discussed in Section 4. Calculations were performed on these cells involving increasing orders of anisotropic scattering in the moderator in order to assess the importance of accurately accounting for the scattering when computing the disadvantage factor. Scattering in the fuel was assumed to be isotropic. 2. BASIC

FORMALISM

The homogeneous transport equation for one-speed plane geometry problems with no azimuthal dependence and with general anisotropic scattering is (MIKA, 1961)

whereY,(x, p) is the angular density of medium a and where the order of anisotropy is denoted by N,. The oLccare the Legendre components of the scattering kernel; ,u is the cosine of the angle between the neutron direction and the x-axis, where x is measured in units of mean free paths. Using the completeness theorem of MIKA (1961), the solution to the angular density may be written in the form

+

s

’ &(+$Jy, -1

ru)e+

dv.

(3)

The B,(+J and &(Y) are expansion coefficients to be determined in terms of the following functions which were defined by Mika:

Disadvantage factor calculations

v+ *kl+l,a(17) - bm,(r)

25

+ k&l,&) = 0, 1 2 2

(64

MO,

(6b)

go, = 1,

g1cz=

h,, = 21 + 1 - OLa

(7)

(8) Here q denotes either the continuum eigenvalue, - 1 < v -=L1, or the discrete eigenvalues fvj;.or,j = 1 to M,, where AI, I IV, + 1. The value of vjolis calculated from the dispersion relation A,(fVj,) = 0 (9) where A,(Z) = 1 - ; j;rge

dp.

(10)

Numerical use of equation (10) cannot accurately determine the largest discrete eigenvalue for cases in which the mean number of secondaries, UJ,,,approaches unity and h,, approaches zero. To determine the eigenvalues for these cases, a reduced form of the expansion* of HOLTE (1948) should be used, -

1

= hohl 1

VI2

1

576h12 ) +. . . . h22h32h4

(11)

Finally, a useful relation between A, and 1, is obtained by considering the boundary values of the function A,(z) along its cut in the complex plane, (-1, 1). These boundary values are

n,*(v) = A,(v) f 3. DISADVANTAGE From

equation

Fg,(v,

v).

(12)

FACTOR

(3), the solutions for the angular densities are taken to be

Y,(x, P) = MvII)+r(rI1,

1u)e-5/yl~+ ~I(-vll)+l(-vrI, i-

Iu)e3c’yll ’ AI(v)&(v, ,u)e-“1”dv s -1

(13a)

and Y2(x, P) = A2(v12)+2(v12,cl)e(aS-A-z)lvl~+ A2(-~IY12)+2(_v12,~)e-(@+A-z)/vlz +

’ &(v)$~(v, p)e(a+A-z)iv dv + 1 s -1

(13b)

where the last term of equation (13b) is a particular solution of the inhomogeneous form of equation (2) containing an assumed uniform isotropic source resulting from neutrons slowing down in the moderator. The source is normalized to an arbitrary * Holte’s expansion was for a more general problem and has been restricted to the present form by KuS~R

(1956).

G. W. ECCLFSTON, JR. and N. J. MCCORMICK

26

magnitude. 1, are

The boundary conditions onYI(x, ,u) andY,(x,

,u), valid for -1

(14a)

%(O, P) = YI(O, -P) Yz(a + A, P) =Yz(a

I ,U I

+ A, ---EC)

(14b)

Yl(Q, P) = Ys(G P).

(14c)

Equations (13a) and (13b) have been restricted to cases with only one pair of discrete eigenvalues, +vltc, for each medium. This is because multiple pairs of discrete eigenvalues did not happen to arise in any disadvantage calculations performed here (although more than one pair of discrete modes can occur for model scattering kernels such as exp [lo@,, - l)], as was considered by FRANCIS et al. (1963), where p. is the scattering angle). In problems in which more than one pair of discrete eigenvalues exist, equations (13) would contain more than one pair of discrete expansion coefficients and could not be solved without altering the solution procedure of Bond and Siewert used here. Boundary conditions (14a) and (14b) require symmetry of the expansion coefficients in equation (13), i.e. A,(yr,) = U--y,,)

and

dI = 1,2.

A=(Y)= A,(+,

(15)

Applying boundary condition (14~) to equations (13) and using equation (15) gives ‘_A1(v)~l(v, ,u)ewa’” dv =

1 + A&12)[#2(v12, ,4eA”‘1~ + +a(-%,

P)e-A”le] + s’ A,(v)#,(v, -1

,u)eA"'dv.

(16)

Equation (16) is used to obtain a form of the disadvantage factor expression amenable for numerical calculations and to determine the expansion coefficients. An expression for the disadvantage factor which is useful for numerical calculations was obtained by FERZIGERand ROBINSON(1965) by taking the zeroth and first moments of equation (16) with respect to ,u and expressing A1(vlJ and A2(v12) in terms of AI(v) and A,(v). Their result may be written as

(s9&4:) +K-- 11+$4;)

(=Q

A

I-K

i

(17)

!

where the following definitions have been introduced: K=A,O-A,O+Gcoth

!!fcoth_!!

Al

(1v12

v12

(18)

1

Ai

AjO= s0

cash

;

dv

a,j=

()

1

i 5

(194

xj =

1

A,j = 2.

A,1 = s0

W’b)

For cases when K = 0 in equation (17), the asymptotic reactor theory expression for

Disadvantage factor calculations

27

the disadvantage factor is obtained. Thus, by computing K we are adding the continuum correction to the disadvantage factor. To determine K requires finding solutions for A,(v) and A,(v). Since explicit expressions for these expansion coefficients could not be obtained, we used the formalism of Bond and Siewert to generate Fredholm integral equations in terms of the expansion coefficients. By iterating, AI(v) and A,(v) (and hence the integrals in equations (19)) could be obtained to a desired degree of accuracy. Following a procedure identical to that used by Bond and Siewert and utilizing equations (A-I)-(A-3) of Appendix A, one can show that A,(v) and AZ(v) are obtained from the solution of the following four coupled Fredholm equations: A”~(~~)Y~Q(Y~>T(Y~)

= S(ym

ool,A>

+ &WR(~,s

yll,A>

CW

(2Ob) A;W~Q12WW

=

S(v,

001,

A)

+

A”22~12)W52,

v, A) 1 +

A”2(v’)R(~‘,

v, A)

ii’Jv’)R(v’,

v, a)

dv’

I --A12W’)yQ12W(v)

=

S(v,

~02,

a)

+

A”,h)Rh,

(2Oc)

v, :) 1 +

dv’.

s0

(204

In equations (20) the following definitions have been utilized: Ac,(v,,) = 4(4

(214

cash (44

L&(V)= A,(v) cash (u/v) A”2h2) A2(v)

=

A2(v12)

=

A,(v)

cash

cash

@lb) WC)

(Ah,)

(214

(A/v)

(214 Qtv,J =

Fg&,

Q[F]

cw Z=YjE

2XY

N&Y, 4 = -x2 - y2

0

y) tanh f

S(x, y, z) = x(1 - y) tanh (z/x) T(x) = tanh (u/x) + tanh (A/x)

Y

- y&,(x, x) tanh

” (X

)I

t2w CW (21i) (21j)

G.

28

W. ECCLE~~ON, JR.

and N.

J. MCCORMICK

In equation (21j), N,,, equals the larger of N,, NB. The quantity Q(v~~)of equation (21f) is inversely proportional to the “effective source strength” given in series form by I~ijNii and USSELI(1968). 4. NUMERICAL

ANALYSIS

Numerical analysis to obtain the disadvantage factor for general anisotropic scattering was convenient to perform. This was because only the Q- and &-values obtained earlier by Bond and Siewert had to be modified according to the definitions of equations (21e), (21f), and (21j). In cases involving either a very weakly absorbing or a slightly multiplying medium, it is necessary to determine the discrete eigenvalues from equation (11). In this situation the functions gra(vra) of equation (6a) and Q(Y~~) of equation (21f) are calculated using the procedure of Appendix B. All calculations of the disadvantage factor were done in double precision using the 21-point Gauss quadrature of KRONRAD(1965). Confidence of the accuracy in some calculations was obtained using 81-point Gauss quadrature since differences between the two results were of order 1O-6. A measure of the importance of anisotropy in physical applications was obtained by analyzing two slab cells consisting of natural uranium and water. Two cells of THEYS(1960) were used. Scattering in the strongly absorbing fuel was assumed to be isotropic for all calculations, with o,,r = 055370. Anisotropic scattering in the moderator was assumed to be entirely due to elastic scattering (isotropic in the center-of-mass system) with an effective mass, m, of either 1.00 or 2.06. The scattering coefficients were obtained from the equation (DAVISONand SYKES,1958) i [(m2 - 1 + s2)1/2+ s]2 2m(m2 _ 1 + $2)‘7”;-Pi(s) ds.

(22)

Evaluation of the integral in equation (22) is especially easy for m = 1, as is discussed in Appendix C. The use of m = 1.00 permits a determination of the maximum effect that anisotropy can have upon disadvantage factor calculations in the slab cells, while the use of m = 2.06 from KRIEGERand NELKIN(1957) provides a more realistic indication of the anisotropy effects of water. Table 1 contains the discrete eigenvalues as calculated from equation (11) for the various orders of anisotropic scattering in the moderator. To obtain a converged value of the disadvantage factor, the expansion coefficients and integrals in equation (20) were obtained to within a value of lo+ between successive iterations. In all cases the disadvantage factor converged to a value within

TABLE I.-DISCRETE EIGEIWALUES FOR ANISOTROPIC SCATTERING IN WATER Discrete eigenvalue, y1

Effective mass. In

0

Order of anisotropy of scattering, N 2 1 3

1

6.33217

10.87660

10.88966

10.88966

2.06

6.33217

7,684tl

7.68542

7.68542

4

10.88966 7.68542

Disadvantage factor calculations 10” between successive iterations. condition (14c), given by

s s

29

In addition, the first two moments of the interface

1

Yda>IU)P~ dp

-1

1

-1

-1

, k = 0,l

%(a, P)P” do

were computed and were found to be consistently less than 2 x 10e4 for cases involving one and two moderator expansion terms. The disadvantage factor for the two cells is given in Table 2. Also listed in the table is the percentage difference between the calculated value of 5 and the result for isotropic scattering. The percentage differences in the computed values of the disadvantage factors of Table 2 for cell A are seen to be approximately half those of cell B. In cells of type A and B, the percentage of the fuel and moderator in a unit cell is the same and hence the cells differ only in that system B is more lumped. Since the scattering in the fuel is isotropic, any effect of scattering upon the disadvantage factor arises in the relative widths of the moderators. The moderator width in cell A is 1.63 1 mean free paths versus a width of 3.262 mean free paths in cell B, and hence multiple anisotropic scatterings in the moderator are more predominant in cell B than in cell A. This causes a larger percentage difference from the isotropic result for 5 of cell B than cell A when anisotropic scattering terms are taken into account. It is seen from Table 2 that the effects of the linear and quadratic anisotropy terms upon the computed value of the disadvantage factor are quite important when the moderator is calculated using m = 1. By using 4th-order anisotropic scattering in cell B, the percentage difference of 5 from the isotropic result was found to be 27.8 per cent. This difference corresponded to an error in the isotropic result of 38.6 per cent relative to the 4th-order result. However, the quadratically anisotropic expansion result differed from the 4th-order result by only 1.4 per cent. Thus, for cases in which the severity of scattering anisotropy is less than that for m = 1, all scattering terms higher than the quadratic will have only a very small effect upon the calculated disadvantage factor. (See Appendix C for a numerical verification of this fact.) Hence at most only quadratically anisotropic scattering was considered for the m = 2.06 calculations shown in Table 2. The percentage differences between the computed values of 5 for anisotropic scattering and for isotropic scattering are less when m = 2.06 than when m = 1.00. By using 2nd-order anisotropic scattering in cell B, the percentage difference of 5 from the isotropic result was found to be 13.7 per cent. This difference corresponded to an error in the isotropic result relative to the 2nd-order result of 15.9 per cent; use of a lst-order anisotropic expansion gave an error relative to the 2nd-order result of 6.7 per cent, and hence the inclusion of the 2nd-order term of the scattering kernel appears warranted. Finally, comparisons of the results for anisotropic scattering obtained by the singular eigenfunction method used here with those results obtainable by theAmouya1 -Benoist-Horowitz method employed by THEYS(1960), the integral transport theory approach of CARLVIK(1967), or the S, method, are expected to show agreement much the same as the comparisons for isotropic scattering exhibited by BOND and SIEWERT(1969).

0.96275 0.96275

0.99163

1.98326

O-99163

@99163

1.98326

1.98326

0.99163

0.99163

0

%a

0.99163

aI2

AND

WATER

0

0,24224

CELLS

0

0

-0.37186

0

0

0

WP2

SLAB

ANISOTROPIC

1.1405

1.1986

1.0554

1.0664

1.1634

1.2317

Cell A*

Disadvantage

WITH

to thicknesses of 0.2 and O-7 cm, respectively. to thicknesses of 0.4 and l-4 cm, respectively.

0

0.0

0

0

0

W3S

0

1.23954

1.23954

0

0

War

Legendre scattering coefficients of moderator

URANIUM

* Cell dimensions are a = 0.1434 and A = 1.631 corresponding t Cell dimensions are a = 0.2868 and A = 3.262 corresponding

2.06

1

1,2.06

Moderator effective mass, m

TABLE 2.-DISADVANTAGE FACTORIN

14049

1.4988

1.1752

1.1924

1.3599

1.6284

Cell Bt

factor, 5

SCATTERING

7.4

2.7

14.3

13.4

5.5

0

Cell A*

13.7

8.0

27.8

26.8

16.5

0

Cell Bj’

‘A Difference from isotropic result

E A

P

E

; .+

6

?

8 “2

s

m 8

3

F

Disadvantage

factor calculations

31

Acknowledgments-This work was supported in part by National Science Foundation Grant GK-1559 and is based upon a portion of the dissertation submitted by the first author to theuniversity of Washington in partial fulfillment of the requirements for the Master of Science degree. The authors appreciate helpful discussions with Dr. A. G. Gibbs, Professor C. E. Siewert, and Mr. H. A. Larson and a numerical program provided by Mr. G. R. Bond. REFERENCES BOND G. R. and SIEWERTC. E. (1969) Nucl. Sci. Engng 35,277. CARLVIKI. (1967) Rep. No. AE-279, Aktiebolaget Atomenergie, Stockholm. CASE K. M. (1960) Ann. Phys. (N. Y.) 9, 1. CHANDRASEKHAR S. (1950) Radiative Transfer, Oxford University Press, London. DAVISONB. and SYKESJ. B. (1958) Neutron Transport lkeory, Oxford University Press, London, p. 234. FERZIGERJ. H. and ROBINSONA. H. (1965) Nucl. Sci. Eng. 21, 382. FRANCISN. C., BROOKSE. J. and WATSONP. A. (1963) Trans. Am. Nucl. Sot. 6,283. HOLTE G. (1948) Arkiv. f mat., astr. o. fys. 35A, No. 36. 1~0~0 E. and USSELIA. I. (1968) Nucl. Sci. Engng 34,39. KAONRAD A. S. (1965) Nodes and Weights of Quadrature Formulas, Consultants Bureau, Inc., New York. KRIEGERT. 5. and NELKIN M. S. (1957) Phys. Rev. 106,290. KUSCERI. (1956) J. Math. & Phys. 34, 256. KUSZELLA. (1961) Acta Physica Polonica 20, 567. LATHROPK. D. and LEONARDA. (1965) Nucl. Sci. Engng 22,115. MCCORMICKN. 5. (1969) Nucl. Sci. Engng. 37, 243. MCCORMICKN. J. and KUSCERI. (1966) J. Math Phys. 7,2036. MENDELSONM. R. (1966) J. Math. Phys. 7, 345. MIKA J. R. (1961) Nucl. Sci. Engng 11,415. MINGLE J. 0. (1969) Trans. Am. NucI. Sot. 12, 634. MUSKHELISHVILI N. I. (1953) Singulur Integral Equations, Noordhoff, Groningen, Holland. THEYSM. H. (1960) Nucl. Sci. Engng 7, 58. VANMASSENHOVE F. R. and GROSJEANC. C. (1965) Proc. Second ZnterdiscipIinary Conf of Electromagnetic Scattering, Gordon and Breach, New York, pp. 721-763.

APPENDIX

A

Two-media cross-product integrals” The new two-media cross-product integral relations derived here, which involve the eigenfunctions +,(r], y) for ~7= Y, -1 < v < 1, and for *vja,j = 1 to M,, are as followst:

s 1

_l

hdv, ~Mjdv’, AP dp = 4g(v, 4

where F,p(q, 7’) and Q,-&v) are defined in equations and (A-2), together with the equation

sv +

~~Qa,e(v)QJ - v’)

(21j) and (21e), respectively. Equations (A-l)

1 _-l Mrtqz,

~Mkq,, PIP dtL = fvjaQ@jd (A-3) s of MCCORMICKand KUSEER (1966), where Q(vja) is defined in equation (21f), are sufficient to determine the expansion coefficients in two-media problems. The form of equation (A-2) suggests that its straight-forward use in the determination of expansion coefficients from equations such as equation (16) may lead to singular integrals, in which case the resultant equation for the continuum expansion coefficient is a singular integral equation. (For the disadvantage factor problem of Section 3, however, BOND and SIEWERT(1969) were able to avoid these

* Preliminary results reported at the Conference on Transport Theory, Blacksburg, Virginia, January 20-24, 1969. t The one-medium orthogonality conditions of MIKA (1961) are obtained by restricting equations (A-l) and (A-2) to the case of one medium: s

;I +ix(~, ru)+i*(~‘, CL)~dp = 0,

3 St rl’

G. W. ECCLESTON,JR. and N. J. MCCORMICK

32

singular integrals by the judicious procedure used to derive equations (20).) Hence it appears that equations (A-l) and (A-2) are normally most convenient only in infinite lattice problems, such as the one considered here, or for cases not applicable to the “wide-region approximation” as used by MENDELSON(1966). For multimedia problems in which the thicknesses all exceed (say) 4 mean free paths so that the “wide-region approximation ” is applicable, principal value integrals may be avoided by calculating the H-function of CHANDRASEKHAR (1950) and using the integrals of MCCORMICK (1969) rather than using equations (A-l) to (A-3). The derivation of equations (A-l) and (A-2) follows the procedure used by MCCORMICK and KUSEER (1966) and begins by defining a general function Jc&,

1 g&G p) gj&‘, p) z’) = s -1 --.Fpdp z--p

(A-4)

where z and z’ are two nonidentical complex variables. Later, after performing some reductions in the expression for &p(z, z’), the values of z and z’ will be allowed to become the eigenvalues 77and Q’, respectively. The function g,&z, p) is the analytic continuation of the polynomial gor(q, cc) into the complex plane. Using partial fractions and equation (5), we arrange (A-4) to obtain z + z’

Jqdz, 2’) =

(A-5)

where use has been made of the definition &(r)

= /TIgz

Pr(/+

do.

(A-6)

The expression for &,(z) may be written as (MCCORMICKand Ku&?ER, 1966) Rl,(Z) = 2hragra(z) - 2h (z)Pr(z) 21+1 a *

(A-7)

Using equation (A-7) in equation (A-5) shows that J@(Z, r’) =

I z+zyy;

bzp - walg&)g&‘) - kz x Ig&‘,hL(4 - gc&,z’)A&‘N, z # z’

(A-8)

where we have used equation (7) to show that &

[eJz&, -

wzahpl= OZLl - ~Za.

(A-9)

Equating the definition of JQ(z, z’) from (A-4) with the result in equation (A-8) and multiplying by zz’/4 gives 1 ;g&, -. f -1

Z-P

p) fga(z~* p) 2’ - ,U

p dlu = Qdz> 2’). y&Z

22’12

- z’ - z

x rgsw, dhdz)

- g&,

z’)AgW,

z # z’

(A-10)

where equation (Zlj) has been used. Equation (A-10) represents the final form to be used to obtain the full-range cross-product integral relations. Equation (A-l) follows directly from equation (A-10) with the use of equations (8) and (9) whenever z and/or z’ is a discrete eigenvalue and with the use of equations (4a) and (12) whenever z or z’ is a continuum eigenvalue. (To perform the derivation for the case where (say) z tends to the eigenvalue Y, - 1 < Y < 1, we substitute z = v -+ ie into equation (A-10) and average the results from above and below the interval as E approaches zero in the limit and then use the Plemelj formulae.) Equation (A-2) is obtained from equation (A-10) by letting z and z’ become the continuum eigenvalues v and v’, respectively. When this is accomplished we have

Disadvantage

factor calculations

33

where use has been made of the Plemelj formulae and equation (12). Equation (A-11) involves the difference of two distributions and may be expressed as a product of two distributions using the Poincare-Bertrand formula (MUSKHELISHVILI, 1953).Using this formula and equations (4a) and (21e) gives equation (A-2) after rearranging.

APPENDIX

B

Special equations for 0~~ near unity The principal difficulty encountered when wccr is near unity arisesin the calculation of the quantities gla(vra), I = 1 to N. (Note that henceforth the subscript a denoting the weakly absorbing medium is omitted.) It is known (KuSEER, 1956) that the value ofg,(v,), I 2 1, actually tends to zero as v1 tends to infinity. However, straight-forward use of the recursion relation of equation (6a) to calculateg,(v,), especially for larger values of 1, involves dealing with terms of order vlt when v1 may be very large or tending to infinity. To circumvent this difficulty, we rearrange the recursion relation to yield a continued fraction expression which is more amenable to calculations when v1 is large. Equation (6a) is first rewritten in the form (1 i-

1) gt+,(h) -’ 1

-.glo VA

@- 1)

and this equation is applied recursively to obtain the continued fraction result

gzh) = gz-I(%)

$yz (C.F.),,

/>I

(B-2)

where (C.F.), is defined as the continued fraction

(B-3)

Since, by equation

(W, g&z) = 1, equation (B-2) may be used to obtain all g&5), 12 1. Indeed,

g*w = v1 ‘h

hl! fI (C.F.),, h

1 2.. * c3=1

12 1.

Two pertinent observations may be made concerning the above results: 1. Arbitrarily setting (C.F.), = l,j = 1 to I, is a good first-order numerical approximation when v1 is large. When (CF.), = 1, this corresponds to using the recursion relation of equation (B-2) to calculate g,(vl) instead of the recursion relation of equation (B-l), i.e. it corresponds to setting g,+,(v,) identically to zero. 2. Equation (B-3) may be used to show the recursion relation for (C.F.),, i.e. (C.F.), = [I - ;zl

(C.F.),+;I -I*

(B-5)

Thus, for numerical calculations, one can go out to a large I and assume that (C.F.),,, = 0; one then obtains (CF.), = 1 and (C.F.), 2 l,j < 1. The accuracy of (C.F.), tends to increase as one gets to smaller and smaller values of k. A second difficulty encountered when o ,, is near unity arises in the calculation of the discrete normalization constant, Q(vl), given by equation (21f). This term may be treated in two stages. First, since it is known (MCCORMICKand KIJ%ER, 1966) that

34

G. W. ECCLESTON,JR. and N. J. MCCORMICK

one can immediately show that the series (B-7 ) yields good accuracy for large Q. Second, it is also necessary to obtain an asymptotic series for g(~r, or) as defined in equation (5). Using equation (B-4) and the series representation for the Legendre polynomials, one finds

(B-8) where [l/2] is the maximum integer in l/2. The first few terms of this series are as follows:

g(v,, VI) = wg +

1

(CF.), f

1

-t s*

w2 &Q3

(5 - 3%~2) fi

k=l

v,-2)(c.F.),(c.F.),

(B-9)

(C.F.),

fi (C.F.), + . . . k=l Equations (B-7) and (B-8) may then be combined together to give Q(vr) values of vt.

of equation (21f) for large

APPENDIX C A more accurate calculation Although a 4th-order anisotropy expansion has been used in Section 4 and previously to treat elastic scattering when m = 1 [see, for example, LATHROP and LEONARD (1965)], the scattering frequency is negative over a portion of the angular range. This is because a low-order polynomial approximation is unable to adequately represent the condition of no backscattering. The range over which the scattering function is negative can be minimized by approximating the function by a higherorder polynomial; to investigate this effect, calculations for 30th-order anisotropic scattering of the moderator were carried out for cell A and the disadvantage factor was found to be 1.0513. This is to be compared to a disadvantage factor of 1.0554 for 4th-order anisotropic scattering. Thus, the loworder approximations do indeed give quite good results. The expansion coefficients, wla, used in this 30th-order calculation for m = l*OOwere obtained from equation (22) 2

= 2(21 t- 1) j; #r(s) ds.

(C-1)

It is possible to evaluate the integral in this expression to show, for integer k, that %?k+l -

00

=

26,,,

k>O

(C-2)

W2 - = 514, % k 2 2.

(C-4)

When generating the expansion coefficients, however, it is easiest to use the equation %k -=%k--2

(4k + 1)(2k - 3) 2(4k - 3)(/c + 1) ’

k 2 2.

(C-5)