THE UNIVERSITY OF MICHIGAN INDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING ONE-SPEED NEUTRON TRANSPORT PROBLEMS IN PLANE GEOMETRY Norman J. McCormick A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in the University of Michigan Department of Nuclear Engineering 1964 February, 1965 IP-693

ACKNOWLEDGEMENTS I gratefully thank Dr. I. Kuscer for his interest and encouragement and for his extensive help in completing Chapters III and VII of this thesis. I also appreciate the helpful discussions with Dr. M. R. Mendelson concerning the material of Chapter V. To Professors G. C. Summerfield and F. C. Shure, co-chairmen of my doctoral committee, I express my appreciation for their patience and assistance. I also wish to acknowledge the advice of Professor P. F. Zweifel, who first suggested a thesis topic in the field of neutron transport theory, and who served as a co-chairman of my committee while on the other members of my doctoral committee, Professors K, M. Case, R. K. Osborn, and A. Z. Akcasu. The financial support from an Atomic Energy Commission Special Fellowship in Nuclear Science and Engineering for three years from 1961 to 1964 is hereby acknowledged. I am also grateful for a fellowship from the University of Michigan Horace H. Rackham School of Graduate Studies for the present school year. The assistance of the Industry Program of the College of Engineering in the final preparation and reproduction of the manuscript is also acknowledged. The patience and understanding of my wife during these years is also greatly appreciated,

TABLE OF CONTENTS Page ACKNOWLEDGEMENTS o....................................... i.i LIST OF FIGURES....................................... iv LIST OF APPENDICES....................................o. o v CHAPTER I INTRODUCTION..................................... 1 II REVIEW OF CASE'S METHOD............................ 7 A. Isotropic Scattering.............................. 7 Bo Linearly Anisotropic Scattering................... 17 III ORTHOGONALITY RELATIONS........ 24 A. Isotropic Scattering........................... 24 B, Linearly Anisotropic Scattering........ O o0 00 35 IV HALF-SPACE APPLICATIONS WITH LINEARLY ANISOTROPIC SCATTERING........a................................. 46 V SLAB ALBEDO PROBLEM FOR ISOTROPIC SCATTERING,,... 66 VI QUADRATICALLY ANISOTROPIC SCATTERING IN A NONABSORBING MEDIUM......... o............................. 92 VII GENERAL ANISOTROPIC SCATTERING.......................o 110 VIII SUMMARY........................................ 135 APPENDICES........................................o...o.... 137 REFERENCES..........O.....*........................ 191 iii

LIST OF FIGURES Figure Page 2.1 Behavior of 0()................................... 12 5.1 Coordinate System for Slab Geometry................... 67 A. 1 The Path C for Equation (A.2)......................... 139 C.1 The Path C for Equation (C.5)........................ 148 C.2 The Path C Used to Obtain Equation (C.10)............ 150 Co3 The Path C Used to Obtain Equation (C.19)............. 154 E.1 The Path C for Equation (E.5)........................ 160 E.2 The Path C Used to Obtain Equation (E.10)............. 161 E.3 The Path C Used to Obtain Equation (E.18)............. 166 G.1 The Path C Used to Obtain Equation (G.4).............. 176 H,1 Coordinate Systems Used in the Chandrasekhar and Case Formulations.................................... 180 iv

LIST OF APPENDICES Appendix Page A DERIVATION OF EQUATIONS (2.30) AND (2.31)............. 138 B DETERMINATION OF THE GENERAL-RANGE, HALF-RANGE, AND FULL-RANGE WEIGHT FUNCTIONS....................... 142 C EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHAPTER IV.......................... 147 D MODIFIED ORTHOGONALITY RELATIONS FOR ISOTROPIC SCATTERING............a......... 156 E EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHAPTER V....................... 158 F GENERALIZATION OF SLAB BOUNDARY CONDITIONS............ 169 G EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHAPTER VI....................e.... e 175 H RELATIONSHIP BETWEEN CHANDRASEKHAR AND CASE QUANTITIES AND NEW INTEGRAL EXPRESSIONS FOR CHANDRASEKHAR'S XAND Y-FUNCTIONS....................................... 180 v

CHAPTER I INTRODUCTION The macroscopic behavior of neutrons in a system is governed by the linearized Boltzmann or transport equation | V0 a - -v) _ _) v t) -(V <(v ~ v r j t(v )_ + _(r t) ( where the notation is that of Reference 1. It will be assumed that all neutrons, on the average, are of one speed and that the medium is homogeneous. With these assumptions, the time-independent equation is'V + I] ( r<) = c f(f >~f) /(I j)d )g (1.2) where distances are measured in units of mean free paths and c is the mean number of secondaries per collision. The scattering function f(Q' - Q) is defined so that f(_' -Q _)dQ is the probability that a neutron scattered from direction ~' will have velocity direction in dO about Q after the collision. For an isotropic medium f(Q_ - _) = f(Q' Q~ ) and we write f(o i) = X (I fi) (1.3) truncating the expansion after N terms. Assuming a one-dimensional geometry and symmetry in the azimuthal direction about the x-axis and -1

-2applying the spherical harmonics addition theorem,(l) the transport equation becomes N I [,Ax + {] Y(x M) Y 2 X Q') /(x74)k (1.4) -=0 -I For various problems, the angular density at any position x and angle cos-1l, 4(x,k), will be obtained from this equation. Various exact methods have been used to solve the above integrodif ferential equation. Chandrasekhar, (2) with a method developed by Ambarzumian, (3) has used Equation (1.4) to study photon transfer ignoring polarization. The method is based upon recognizing physical principles of invariance which the angular density must obey. These invariance principles and the transport equation lead to solutions for the outgoing densities from slabs of finite or infinite width under appropriate boundary conditions. An extension of the approach used by Ambarzumian and Chandrasekhar to inhomogeneous regions was made by Bellman and Kalaba.(4) The method of solving the transport equation by use of invariance concepts is known as the method of invariant imbedding. Another exact approach used to solve Equation (1.4) is the Wiener-Hopf method of Fourier transforms.(5'6) The transport equation is first transformed into an integral equation and solutions are then obtained for various boundary conditions. Case(7) has obtained exact solutions of Equation (1.4) for isotropic scattering of neutrons by still another method. The technique

consists of expanding the angular density in terms of the eigenfunctions of the homogeneous transport equation. The expansion coefficients are determined by applying boundary conditions; the method of solution is thus analogous to the eigenfunction expansions used to solve other partial differential equations of physics. Case's treatment of the transport equation has an advantage over the method of invariant imbedding since the angular density is given for any arbitrary distance inside the plane-parallel medium and may also be used to give the surface quantities obtainable from the older method. The great virtue of the Case method is that it requires the difficult mathematics, for a given class of problems, to be done only once; the solutions for various physical applications can then be found immediately. Case used his method to solve full-space and half-space problems with isotropic scattering and briefly discussed the extension to problems in a two adjacent half-space geometry. Explicit solutions to problems in the latter geometry were obtained by Mendelson and Summerfield.(8) The method has also been used to attack problems which have no known closed form solutions. Mitsis(9) obtained a solution for the critical slab problem and results for special problems in other than slab geometry,(l) while Mendelson examined various reflected slabs. References 12 to 14 also deal with problems of this type, but explicit results were not obtained. The method of Case has also been used to examine energydependent problems(15 18) and time-dependent problems,(7l19)20) but we shall not be concerned with problems of this type here. Our main

-4concern will be with problems where the neutrons can be scattered in an anisotropic manner in the laboratory coordinate system. Mika(21) extended the completeness theorem of Case to consider problems with N'th order anisotropic scattering and Shure and Natelson(22) later showed how certain explicit results can be obtained. In addition to the above problems in neutron and photon transport theory, a similar eigenfunction expansion technique has been used to study plasma oscillations.(23'24) The method has also been used in gas dynamics to obtain solutions to the slip-flow problem.(25) The generality of the approach to problems in linear transport theory has been discussed in Reference 26. It is the purpose of this thesis to solve some of the simpler problems in one-speed transport theory by the method of Case. When some of these problems have been solved, it is felt that Case's method, with its inherent advantages, can then be used to solve more difficult problems not readily amenable to attack by the older methods. Chapter II is devoted to a review of Case's method in situations were the scattering is isotropic or linearly anisotropic. Here, results are first derived for the partial-range case where a < L < P. Although these relations have no known physical application at present, they conveniently tie together the full-range (a = 1, P = 1) and halfrange (a = 0, P = 1) cases. Chapter III contains the derivation of general orthogonality relations which may be used to solve problems in half-space and fullspace geometries when the scattering is isotropic or linearly anisotropic.

-5The orthogonality relations in the full-range case are simply those obtained previously by Case(7) and Mika(21) for isotropic and linearly anisotropic scattering, respectively, and will not be considered further. Since the isotropic scattering relations for the half-range case yield answers to half-space problems which were previously obtained,'7) no half-space applications for isotropic scattering will be explicitly exhibited. Chapter IV, however, contains the solutions to the half-space Milne,(22) albedo, constant isotropic source, and Green's function problems when the scattering is linearly anisotropic. Also included are simplified expressions for the emergent angular distributions and for the densities and net currents on the surface of the half-space. All equations in this chapter may be specialized to give the results for isotropic scattering. In Chapter V, orthogonality relations are used to attack the slab albedo problem. Because the solution obtained is not in closed form and is consequently more involved, we will restrict our attention to the case of isotropic scattering in this chapter. In an appendix, some of the results of this chapter are used to derive new integral expressions for the X- and Y-functions of Chandrasekhar,(2) which appear in the application of invariant imbedding techniques to slab problems. In another appendix, some results are used to consider the slab albedo problem for a slab surrounded on one side by a vacuum and on the other side by a specularly and diffusely reflecting face. In Chapter VI, we examine a very specific problem in anisotropic scattering, namely a medium which conservatively scatters

-6neutrons in a quadratically anisotropic manner (i.e., N = 2 and c = 1). Although this problem is not of much interest in neutron transport theory, it has been extensively investigated in the past with regard to photon transport theory. It is of some interest to us, however, since it provides additional insight into the nature of the half-space orthogonality relations in a situation which is no more complicated than that covered in Chapter IV where N = 1 and c < 1 The Milne(22) and albedo problems are then explicitly solved. Chapter VII contains a speculative approach for obtaining the orthogonality relations for general-order anisotropic scattering with c < 1.

CHAPTER II REVIEW OF CASE'S METHOD A. Isotropic Scattering If the medium is assumed to scatter neutrons isotropically in the laboratory system, the homogeneous transport equation is [r ~]yc(/) Cj' x,/')+l (2.1) Separation of variables may be achieved through the ansatz - x/7/ t ( XA) - e (2.2) with the allowed values of the constant v to be determined. Because the transport equation is homogeneous, the eigenfunctions may be normalized in an arbitrary manner and the normalization condition is selected such that* It~~~ ~~~~~p~,t~,O =1~ {(2.3) For -1 < v < 1, this leads to the eigenfunctions )=;s ) C12 P e t+ ) (2 (2.4) * This is permissible except for the trivial solution v(u) 0. -7

-8where X(v) is determined from (2.3) and (2.4) as AV)=I- ) = ( -zL -g(-l'- l (2.5) 6(v-4) is the Dirac 5-function. The symbol P is a reminder that we must take the Cauchy principal value of any integral over v or I. For c < 1, the argument principle(27) shows that there are also two real eigenvalues + V outside of the interval -1 < v < 1 with corresponding eigenfunctions O+(O) given by (2.6) The eigenvalues + vo are determined from the normalization condition and are roots of the equation Q (+ od) - O (2.7) where () = - (.8) The function A(z) is analytic in the complex plane excluding a cut from -1 to 1 on the real axis. For large z 4( Aei)= A(o)- 2 /,:2, _- / -. (2.9) i coDJl —

-9For c = 1, the two roots coelesce at vo = co and it is easy to verify that two linearly independent solutions of the transport equation are Vl(X,p) = and 72(x,>1) = 2(x-). In this case, lim z2A(z) For c > 1, there are again two discrete z - c>1 3h roots + vo, but in this case vo is imaginary. The Plemelj formulae* are of fundamental importance in the solution of the neutron transport equation by Case's method. Consider the function F(z) defined by F(z), H (2.10) F(z) is analytic in the entire complex plane, cut along the real axis from a to 8. For z approaching v where a < v < P,**'the boundary values of F(z) as z approaches v from above/below are denoted by Ft(7)=, F (-/ Le) (2.11) E -+0+ and may be calculated from F ~(P) =C4)c2 fz) - (2.12) These are the Plemelj formulas. They may also be written as *Reference 28, page 42. **The behavior of F(z) in the neighborhood of the endpoints ac is discussed in Reference 28, Chapter 4.

-10 - ~ -.i:: +F-()] -- /:C,64A Tr LF (v) + ( = P -- (2.13) i+ ) FIV)) = f(V) (2.14) These equations may be used to relate i(v) to the boundary values of A(z). Use of Equation (2.12) gives Ago = q fle) L T cU -j t agv)u(2.15) and addition of the two equations in (2.15) gives X(7')= ~[A+'v #A-i9] (2.16) Case(7) has shown that the functions (2.4) and (2.6) form a complete set of functions of i for Ct in the interval [-1,1] The solution to (2.1) is therefore given by the expansion -q e + f) Afr)+,, A)e ch/.(2.17) Through the application of boundary conditions of a specific problem to (2.17), one ultimately arrives as an equation of the type - /44i+ ~) + -- A ~ V)(2.18)

-11where -1 <a < a < 1 t, e [a]. Case has also proved this more general "partial-range" completeness property of the eigenfunctions Ov(p) by directly solving the singular integral Equation (2.18). A sufficient assumption in the proof of completeness is that V?(I,) satisfy the Holder condition(28) |I$(1) - t(i2)1 < M I1 - i L21p with < p < 1 for a < L < f and vary as MIiL - aIK or MIg - BIK -1 < K < 0, at the endpoints.* In the proof of the completeness of the expansion in (2.18), a general X-function, denoted here by Xg(z), is constructed such that: 1) it is analytic in the complex plane cut from a to P 2) it is non-vanishi.r in the finite cut plane 3) it obeys the ratio condition _____) A\ (s) C en <. < ej (2.19) X8 T() Am(v) 4) it varies as iz - aIK and I - IK as z - a and, respectively, where -1 < K < 0 o** One function, denoted by Xgo(z), which has properties 1-3 is X~ (z)= J [t) A- ]d ](2.20) * The proof is immediately extendible to the case where t(P) is any distribution (for example a s-function) which is the weak limit of a sequence of functions satisfying the Holder condition (Case,private communication). **The function is permitted to have a weak singularity at the endpoints, ca and I, of the cut. It may not have a zero at the endpoints, however, because of the manner in which the X-function is used in the constructive proof of completeness of the eigenfunction expansion.

-12This may also be written as @ 6)0 E (2.21) X1o (2) - p/ ff'where Q(p) - arg A+(L). From (2,5) and (2.15), rF c islet) =,<1 Z, (2.22) The function @(u) is an odd function of ju and has the approximate appearance shown in Figure 2.1. -I / r -1 1 _: Figure 2,1. Behavior of @(u) (approximate),. Another function which obviously has properties 1-3 is X'9o(g) 0, (2.23) e _ a

-13since 5 - z is non-vanishing in the cut plane and conditions 1 and 2 are obviously satisfied. This function also satisfies* property 4 since near the point z = 5, ) - Z) (2.24) Likewise, the function,Xo()' o 4, (2.25) has the properties 1-4. Any proper linear combination of (2.20), (2.23), and (2.25) also satisfies the four conditions. This linear combination, the general function Xg(z), thus is given by h(x) by b(X)={ (2.27) lO x O The requirement that sl, 0 and s2 / O in Equation (2.26) guarantees that Xg(z) will have the proper behavior at thate endpoints of the cut. *Reference 28, Chapter 4.

-14Define (v, ) = 2. ) ) C Z (2.28) where C_0(V) = -11)~~~ 2,9 A0~,(2.29) is the particular function corresponding to Xgo(z) of (2.20). By the use of the Cauchy integral formula, it is shown in Appendix A that X} ()-X (~~) - kJ' k (2.30) For z = v, c < v < A, this equation gives (Appendix A) 2 X(X) = P e v ". (2.31) Following Case,(7) for the special case of az = 0 and p = 1, arbitrarily define the half-range X-function, denoted by X(z), as the particular linear combination of Xg(z) with so = 0 and s2 = 1 such that x Jo=, ) (2.32) s2 =1

-15The corresponding y-function is given by tA(2/) = 9 )A(() _'il0) " / (2.33) -0 O s$ =I From (2.32), X(z) X -l/z for large z and, therefore, (2.30) is X(E) = J j )~ (2.34) 0 and (2.31) is (t) () c2 p Ov ).. (2.35) Also, since lim -zX(z) = 1, Z -e 00 _";;?,0 B y t I ~~, J~ Is~,~- /(2.36) Another identity of importance for the half-range case is(7) A_ 2) (2.37) X z X'(- 2) = (.o2_t 2) A l 0)

-16so that, from (2.33) CaV =(2'),~2 0 (2,38) (N.Z - v ) A tm) X (-V) The insertion of (2.38) into (2.34) results in the equation cX(. (2.39) X- 2 Mm\(-) J(N - ~ )2) x-)(/ + ) and this can be used for numerical evaluation of the X(-..) function.(l] If a = -1 and P = 1, then the full-range X-function, also denoted by X(z),* was chosen by Case(7) as the particular linear combination of Xg(z) with so = 0 s = s2 = so that X(Z)- Xs ( (l) (2,40) so = s = Corresponding to this particular function, )= C XI-/l_) _- i', -I =< 1.o /, (2-41) 2 AY7-) — 2 So= 0 The following relations (similar to the half-range relations) hold: X(l= )) (2.42) -I *No confusion will arise if the half-range and full-range X-functions are both denoted by the same symbol since the full-range function is used only as an intermediate step in the derivation of Chapter III and no final answers are expressed in terms of it,

-17PJ-2 -', 1 (2.43) I / TV)+ = 0 (2.44) (2.45) Xci) = (avow EZ)At) (2.46) (VCZ - 7)A (cx) Tb" cv) (-Z/> ) 72 -c. (2.47) B. Linearly Anisotropic Scattering For this case, the homogeneous transport equation is [/ n + 1] iC (-] ( 2.48) Separation of variables is obtained by use of the ansatz of (2.2) with the result

(>- - ) [ A b/rJ']9,f')*t'(2.49) Multiplying (2.49) by dci and integrating from -1 to 1 gives,,^Xef ( -C)J 4P) (2.50) -I -I If the normalization condition is selected as 9i V ) I ) (2.51) -I the eigenfunctions for -1 < v < 1 are where 2_____ ( v d _)4& i t (253 ) -I -I and where d(vt) is defined by d() = I + b(, (- ( c)J/. (2e54)

-19(21) +o For c < 1, there are again two discrete roots,21 + V of the equation A(~72o) =0 (2.55) where = c;E = /-'it JicrC l (2.56) -I -I The corresponding discrete modes are _ = () Z d(9g/A) (2.57) if c f 1. For c = 1, two linearly independent eigenfunctions which satisfy the transport equation are jl(X,1l) = 1 and ~2(x,1 ) = 1 (x -I ) and the two discrete roots again coelesce at infinity. The value of A(z) as z -, c will be needed later and can be found from A Ac — ) i dI)fr (2.58) The result is = (I - c)(' CLb') (2.59)

-20Also, the value of X(v) may again be related to the boundary values of A(z). Use of (2.12) gives A (7V) = 2(. #) ~ a c JS),(2.60) Once again the solution of the transport equation may be found by expanding the angular density in terms of the eigenfunctions of the transport equation. Since there are only two discrete eigenfunctions, the angular density may be expressed as(21) -X/~o x/o r -x/V (7(x<)= )g+=e.+ k fv4)9)MEe 9JA (2.61) and the expansion coefficients will be found from an equation of the form of (2.18). For the general range of ca < t < f, let Xg(z) again denote a linear combination of functions analogous to (2.20), (2.23), and (2o25): X (o ) =X ()[t+ o Xte 0 (262) o For linearly anisotropic scattering, the constants ti, i = 0,1,2 may depend upon c and bI. This Xg(z) satisfies the following identities:

-21(2.63) Xt( (a) - () -= ot (2,64 For the hal~-range, 0 <,u < 1, the X- and 7-functions are now defined by X () (2.65) t2 -a and'b"(2 == 2 ___o O — &7t: _ &() 0i2t. (2.66) to =0 t2- Also, since there is only one pair of roots, + vo, Equations (2.37) and (2.38) remain true. Equations (2.63) and (2.64) are now written as X(X) = 0 c- ttdyA)J' (2.67) tb4(21) h(fJ'Sj p)i JfA>M) A Oc it I (2.68) - z ~-:'

-22/ =| / r)Jfdt)d (2.69) (22) As was pointed out by Shure and Natelson, a technique which minimizes the number of iterations when evaluating the X-function is to use the identity CZ-) ^AI- d__r_) (2.70) where -'(L ) [ Xto) ~ x.(2.71) This equation could also be used instead of (2,39) for the isotropic scattering case if bl were set equal to zero, Equations of the form of (2.40), (2.41), (2.46), and (2.47) are valid for full-range applicationso Additional full-range relations areo X (B) = TYW) (2-72) Th'Xzv) = 7/ p ) -_ I (2073) /" - 2/~~~~~~~~~~~~~~

-t LA cM cr cM II II To - - -a I — I

CHAPTER III ORTHOGONALITY RELATIONS As was pointed out in Chapter II, the angular density is given by (2.17) or (2.61); the solution of every problem of interest reduces to the determination of the expansion coefficients from the equation it q+t/) + + A(V)g)J;^, & (2.18) It is the purpose of this chapter to derive orthogonality relations between the eigenfunctions 0+([) and Xv(p). Use of these relations will facilitate the determination of a+ and A(v) v A, Isotropic Scattering Case(7) has shown that the following orthogonality relations are true in the full-range case: -I / +P t od)= 0 (3.2) (.) Equation (3,1) is an abbreviation for the statement rl.- )[f4 A I I] = A(>) )9A(')Avf)- (3.4)I -24

-25In the evaluation of Equation (3.1), there is the term 2 1 / ~ -I M p 7 A("YA (35) PZ/ The order of integration in this term may be interchanged by the use of the Poincare-Bertrand formula (28) pr I, [pf (;). ] = Pfdv f, P f )(,').. JLL L. LL + 1 T (,, ) t IT 2 ) (3.6) As was pointed out in Reference 299 a mnemonic form of (3.6) may be written as P P i.....' S Pt- F- -1) (3 -7) provided the integration over Ip is carried out first. The order of integration in -1[ A('V) 4)ig1r

-26may therefore be interchanged with the use of (3.7), We will use this approach in the subsequent determination of the orthogonality relations. The orthogonality conditions of interest are obtained by observing that (2.31) may be rewritten in the form'2rX2(X) =f -1), (3.8) From (2.4), 2 -T= (v-,)9r) (3.9) and use of (309) in (3.8) gives We now seek weight function 7g ) which (10) and We now seek a weight function yg(W) which will make Xv(p) and V (0), a < v < v and a < v' < 0, orthogonal over the range a < i < The orthogonality relation and the weight function for the general-range case, Yg(P), are derived in Appendix B or can be observed by inspection of (3.10). The proper weight function for the half-range and full-range cases of interest may then be obtained as in Appendix B where additional orthogonality properties are obtained. The weight function for the half-range case may also be directly obtained without first considering the general-range case (see Reference 1).

-27With our knowledge of the general-range weight function, we therefore write Equation (3.10) for v and v', multiply the first by Xv ()7Yg(L)di and the second by Ov(kl)7g([i)dll, and integrate over o Subtraction yields (-9,P /))9/fA 6P)+f=, rQg3=0)' A v (3.11) provided Xa D) O. (3.12) The normalization relation between,v(iL) and Xv' (i), ca < v < f and (c < v' < B, may now be found. Use of (3.7) gives tat multli(g fni an integratiLon + A3+.)A/ \,) oc - i) i,-3 )) so that multiplication of (3.13) by the weight function and integration over i yields = it + %(74A (v)AY(a) &(tv') ~(3.14) where (3.8) has been used. The weight function 7g(PL) for the general range is still somewhat arbitrary since it is a linear combination

-28of functions corresponding to those of (2~26). One can therefore select a particular one which meets certain additional requirements (see Appendix B)o For c - O and P = 1, the eigenfunction 0+(j) normally occurs for problems where the medium of interest is the right halfspace. From (2~26) it is seen that for so 1 and s2 = vo - 1, a permissible Xg(z) function which, satisfies (3.8) and (3.12) is X (Z.) _i- (7~Io= (3.15) which is simply (vo z) times the X-function defined for the half~ space by Equation (2,32)~* This would correspond to (vo0 ) times the y(~i) function of (2033). Use of (308) and. (3o14) and the halfspace weight function Yg(it) = (vo - [)y(i), determined from Appendix B, gives 924IV) ( = ( (3516) *One might wonder whether it would be simpler to forget about using the half-space X-function used by Case (with so = 0 and s2 = 1) and simply use the X 9-function defined in (3.15) to obtain solutions to half-space transport problems. The Xg-function varies like +1 as z- oo so Equation (2,30) would be X2 ( ) - I =%,. ) WThis equation yields no auxiliary equation of the form of (2,36) and consequently the analysis might be somewhat easier because there would be no simplifying identity to be used, A major disadvantage of the use of Xg(z) would be that the limit of Xg(z) as c-> 1 is infinite so numerical values would. be difficult to tabulate unless the function Xg(z)/vo was considered, The particular half -range function of Case will be retained here because of its general acceptance,

-29I ~I,) fa (7J/O /^~) = 0 (3.17) 0 where O < v < 1, O < v' < 1 With this choice of rg(GI), 0'( and Ov(0) are mutually orthogonal over the range 0 < 1i < 1 Other formulas of interest follow in a straight-forward manner and are shown below, for O < v < 1 and O < v' < 1:; 2, VA(t )+ = C; X(-PO) (3.18) 2+ir, ~~,(Vb~)^GLP)(iD= (~) Xf )(-Jo) (3.19) (4t'+ ),;-v( - >)~ =' ) X'(-z,) (3.20) r Z |?/,,P +(2-,)t,)l. - z+[f7'+ 0)] (3.23) 0~~~~~~~~~~Ct

-30NO~r)(io-,<)t>),2t = (3.24) 0 J ~/'t) -(~/O^s e - T. (3.25) Note that Ov(O) contains no singularity. We have defined 7(i) by I. (0)= tllif*, (3.26) where, for isotropic scattering, O= (3.27) from (2.36). In general, for any range c to A, the integrated moments of the yg-function will be defined as =iSX~~g Jr N? ) +r ( ~3(.28) The above half-range equations may also be used in problems where there is no absorption if they are modified by dividing by vo and letting c -> 1 (vo - co) o Equations (3.16) to (3.25) then reduce to the relations:

-31~'S4%) Yfri) T)d7 = (T)AC(v)A7 z) &(b-v') (3.29) 0 toe) o) + - 0 (3.30) I Jy~,'d/~'t>) A(3-31) I vk)9fvM)ThvAX/-) (3.32) 0 - elf >)+ 2 Xt(-) (3.33) along with (3.26) and (3.27). We therefore see that the half-space weight function when c = 1 is simply 7(i). Using the relation(5) i O2 (I ~ c) = C. )(3.34) Equation (2,38) yields the result ) = 3;0 -<, 14 (3.35)

-32Equations (3,29) to (3.33) are sufficient to determine the expansion coefficients completely since the discrete eigenfunctions are 2 and ~x ~) o It is interesting to consider the full-range case as a possible check on (3.1) to (3.3). From (2.26), with so = 1 and s= s2 = ivo2 - 1), a permissible Xg(z) which satisfies (3.8) is 2 _ z2 This is simply (vo 2) times the X-function defined in (2.40) for the full-range case. Therefore the full-range weight function (v2 - 2)7y(i), derived in Appendix B, satisfies (3.8) and gives )Z/)) = (o-i6'fr)A+(I) A(7) &(-_V)(3-37) i ~t) iP ) (9/0"- I) A o (3.38) 9+ k) f ) = Q. (3539) -I The other equation of interest is NO )= +c X(+s) (34) _ I

-33 - From (2.9) and (2.47) it is seen that (<~ /) /A) Z(I -c (3.41) so Equations (3.1) to (3.3) are immediately obtained and (3.40) reduces to (.42) - 2 i [joZ_! z4zj 2) This last result was also obtained by Case.(7) A few other relations for special problems are (/,4IP)+ = U/(j-C)(3.43),2 (5.44) -I yzt~ I= 720(V -C) (3.45) (3.46)

-34When c -1, the orthogonality relations of (3.1) to (3.3) and (3.43) to (3.46) simply reduce to the identities,1 I t+) gv l&)+ =;//\+(ZI)/-(2/) S(J-1/1) (3.47) _I /A/O (349) (,Jr = O (3-50) (31.o) where (3.34) has been used. The weight function in the full-range case thus appears the same for c < 1 and c = 1. Both the half-range and full-range cases are special cases of the "two-media" problem composed of two adjacent half-spaces of scattering media with different values for the mean number of secondaries per collision, These results are exhibited in Reference 29.

B. Linearly Anisotropic Scattering In order to find orthogonality relations which parallel those of the case of isotropic scattering, one observes that (2.64) may be written as z X(X = ) (),A) O (3.52) From Equation (2.52), 2V = (Zv / ) (3-53) Therefore, x'(0) (2 d) A 1')4'' AO o u (3.54) In the same manner as before, one uses (3.54) to get A-) B )A f= 0) (3.55) provided X (0) 0. (3.56)

Notice that this equation is of no use to us in its present form, however, since an eigenfunction expansion of the form of (2.61) does not contain functions of the type v(jl)e /V/d(vu) but of the type -X/V ~v(ll)e. If (3.56) is multiplied by the symmetric function d(v v') however, and the identity m(i;X) d~) -,(/,) dI(vM) + b, (I- ~)( [;)[') (3 57) is used, then (i~ v') Ire, V) fop r + 6 (v>) )= (3.58) Use of (3.53) gives Ok-C) 0) 0 (3(59) Z Z =. j The first term in (3.59) has the form we desire. The extra term suggests that Equation (3.59) may be a bi-orthogonality relation if properly rewritten. That is, we expect that we will be able to write (2-i')[tt>) B-4 8 =0) (5.60) 2 )'0)

-37with B a constant to be determined. B is obtained via ~~~~-B td, c (-c)T (3.61) Use of (2.64) and d(W/) = d ) + 6, - ) (3.62) gives tt Z2= O) +6) so Equation (3.61) becomes B (3.64) X O(D) + b, (I - c)(') Equations (3.60) and (3.64) combine to give a bi-orthogonality relation for the continuum eigenfunctions. Using this form of the bi-orthogonality relation, consider the normalization condition for the continuum eigenfunctions. Use of (3.7) gives t 1 8=-Y'LV 4(ZZM)9g P) - d ] + A,) A. ) &f-v) &f4-v9. (3.65)

-38Equation (3.65) can be used to give /i)[9, jA ) B 2 j55) = C.(1)4A1(V)&') (3.66) provided + 2. B,0t.t ~d, = O ~(3.67) The relation d (iy) = d (y') -, (I-c) 7("2'?-) (3.68) and Equations (3,61) and (3.63) reduce (3.67) to,[- 2 ( 6 ~ (0 z i(r-c);~) O3.69) Use of (3.53) reduces:3o69)to te identity 0 0. For the half-range case, Equation (2.62) yields X~(c)-Xs,()to...9..L] t,,. (370)

-39where to and t2 are dependent upon c and bl. Since (3.56) stipulates that Xg(oo) O, to $ 0 and (3.70) becomes X(Z) = X(z)[to0A t- (3.71) where (2.65) has been used. Also, using (2.66), T (v) = - u)[t o&t2 Di~ ] V (3.72) We wish to find to and t2 such that the continuum eigenfunction Ov(O) is also bi-orthogonal to 0+(), i.e., that I I~tp)[iP) B ]B tPJ = 0 (3.73) Using (3.71) and (3.72), (3.64) becomes 6i (, -c)[to+',t - t= oA 3- _ 4 bi(Vc)Etot{.a{0i7] ) = 1 (3.74) to, (t, )(- c) (to. tZ. where we have defined - = =..o. (3.75) 0 (0d and have used (2.69). Use of (3.72), (3.74) and the relation 1(o/,) d() - (d 7,,)cd j,,,) -, ( ( 4 )(v:)(3.76)

-40in (3.73) gives to +b, (I-c)(t0tz)7Jvok -A,,- n(l.)J ~)[ to. tz ] Jr, (3.77) Jo By inspection, we see that if to = 1 and t2 = vo - 1 and if (2.68) is used, (3.77) is the identity 0 0-. The half-space weight function 7g(L) = (Vo - )yr(iL) is therefore identical to that for the isotropic scattering case. Equation (3.74) can therefore be rewritten in the form 6B = ( (c)(v-; )0 =Ov=I * (3.78) For linearly anisotropic scattering, the relations of interest are, for 0 < v < 1 and 0< v < v' < 1: t tF[t )+ 2 ](io r)P>+= (N- v() A v)A V-V)(3.79 ) ltriry,+B~ 2 ](M~- )T>)+= Q (3.80)

-4'O+[h t) + B j(0;Y) 5) = 0 (3.81) JR P)[u?, +" ) t" = X(o2 d+ Z) X 7- (3.82) -P+ d (.,,g) t~CIt~)[Rt+)~B C 2-]( /o/) drw = IJ) d(7'Ao) d (-,) (3 83 (PPlo)[X>) to 2-](7JQ-J^)W>) = C 4 Aa)dl s9{ — 9) (3.84) j~+) [ 2,,4 d,i-2) (3.85) 2 d kIt) eta,)+ ~2JLOJ Y( a Ad)d = ) (A ( () d(,,o-v)d( v,) (3.86) 0~tt 2 ] / < Z [ d(iY)) B 2(',,,)~'~)up -'2, o),,) (s.86)

-42t/[ +t)~ + -dh/0))>lr= i i~,t) (3.89) Thus the solutions to half-range problems with linearly anisotropic scattering will be in terms of the calculated functions X(-[l) (or?t') ) and T(i) i = 0,1, and the variables c and bl Since dividing Equations (.79) to (.89) by and taking the limit as dividing Equations (3.79) to (3.89) by vo and taking the limit as c - gives t ) ~ t' )T = Th +()AV) A-V) &(s7-"') (3.91) atlt, t>) r = O (3.92) ~JI~, +) =- Z 3(3.93) rl'L~kR~Y~a03 t~,i~ln,-~ i / ~t) d V (3.94)

-43= zly" ) = X[(-') (3.95) ~) = I{, (3.96) Notice that the form of these equations is identical to that of Equations (3.29) to (3.33) which were obtained for the case of isotropic scattering with c = 1. The above results could also have been obtained by examining (2.49) and (2.50) and observing that, for c = 1 the transport equations for linearly anisotropic scattering and isotropic scattering are identical. For the full-range case, 7g(G) = (vo2 - >2)7(>) is the weight function and Xg(oo) = 1 from (3.36). The moments of the gamma function are now given by Y (i) = v 2y(i) _ 7(i+2) so (3.64) becomes 0 bj (-c)[ioZb"O'-;r) " ='I, ~ = X. (3.97) Use of (2.74) and (2.75) reduces (3.97) to (0) 8A Be -, =,( - =-I, = I (3.98) From Equation (2.47), Y7(1) is an odd function of AL so X (EVEN) o (3.99) "i~ (EE)

-44Therefore, B= O = I) e I) (3.oo00) which means that the orthogonality form of Mika(21) is obtained. The relations are: 7' ~~Yt') 9tliQ)cjt = viA4() A1) v-;") (3.101) r + 0JJ = 0 (5.102) (3.107) I J/^t5VO = +}JO J 7g/o)(l~c)(l~ C. ) (3-107)!'/'?+~)S = - ""r'~'~'~g)( — I ) (3o108)

-45In the limit as c -, 1, the above relations take the form of Equations (3.47) to (3.51). Using (2.56), A'(+ vo) in (3.104) is given by d/\s) += [d(3JTI -,) |_ c Cllv)z(-vO) (3.109) ds =+.o V2~e'oao) (V2-I)d(o7o) t-0 0 0 7 VO

CHAPTER IV HALF-SPACE APPLICATIONS WITH LINEARLY ANISOTROPIC SCATTERING In this chapter, we shall obtain explicit solutions for four common time-independent half-space problems. These problems are the Milne, albedo, constant isotropic source, and Green's function problems. The solution of the Milne problem gives the neutron distribution in a source-free half-space with no incident distribution at the surface and a source of neutrons at infinity. The problem was originally studied by astrophysicists who were attempting to find the emergent photon distribution from the surface of a star which was approximated by plane geometry and a large source at the star's center. The results for this problem have been reported earlier by Shure and Natelson,(22) who obtained their answers without the use of the bi-orthogonality relations of Chapter III. For the albedo problem, we desire the angular density at any position and direction due to an incoming distribution on the surface of the half-space. This incident distribution will be assumed to be azimuthally symmetric and expressible in terms of a Dirac delta function. There are no other sources of neutrons for the system. Thus, this problem is the half-space penetration problem. The third problem concerns a constant isotropic source throughout the half-space with the stipulation that there is no incoming distribution at the surface. This may be interpreted physically as a half-space homogeneous reactor, surrounded by a vacuum, which is -46 -

-47subcritical but is operating in a steady-state condition due to its internal source. In the half-space Green's function problem, there is no incident distribution on the surface and no source anywhere inside the medium except for a single plane source at x = xO emitting neutrons in an arbitrary direction. The solution of this problem can be used to obtain the solution for a half-space, with no incident distribution on the surface, containing an arbitrary source. None of these four problems is ever encountered in reality, but the solutions can be used as approximations to some physical situations. The solutions are particularly valuable, however, when used as a means of checking the accuracy of the approximate methods used to solve the transport equation. The angular density, which satisfies Equation (2.48), is -~x/V x-/xl' (x/) = X4M )e -+x/+ (2.61) If x > 0, the boundary conditions for the half-space a) Milne, b) albedo, c) constant isotropic source, and d) Green's function problems are given by

-48O,, ~?..O a) 6 9; o),,- b) I-?~ )C. I,7?0 c)* I-c O, 0 d) (4.d) and x/'4 ) X/VO a)** O b) x-CLO )0 c) o d). (4.2) *The transport equation for the constant isotropic source problem is Vi]4x- k A) - + I, L ax ],[ ] where (0,1) = 0, k > 0, and lim (x) - 2o The transformation w - X,-t o 1-c T(x,) =?(x,) + qo/(l-c) gives the homogeneous transport equation (2.48) for 4(x,[) and the boundary conditions?(0,k) = -o,kL > 0 1-c and lim rf(xk) = 0. Thus the constant isotropic source problem is X-o 00 reduced to the solution of the albedo problem with a constant incident distribution on the surface of magnitude -qo/(1-c). We will calculate the angular distribution'$(x,) and not the total angular distribution ** From Mika21 we know that the asymptotic distribution far away from a plane source varies as 0+(k)e-x/vo. This means that the angular density varies like 0 (k)ex/vO when approaching the source.

-49An additional boundary condition which must be satisfied is o a) 0 b) 0~/ T-'(XO/) = 0Ic O c) d) (4.3) which takes into account the discontinuity introduced by the source in the Green's function problem. Since boundary conditions (4.2) and (4.3) are valid for all 1i, they are less restrictive than (4.1) and will be considered first. Application of (4.2) and (4.3) to (2.61) enables us to write (2.61) as (/(x7) =- f(x) + ++ +NVe) fA(7),4)e Chid (4.4) where: 9~)e )ci# a) 0 b) 0 c) C( (Xo,/oAO )d), (4.5)

-50The function G(xo,llo - x,~) will be abbreviated by G(x,L) and is the solution to the infinite-medium Green's function problem for linearly anisotropic scattering.(21) Applying (4.1) to (4.4) gives the equation (for t > O) C ~ R }ca) 1r |6~- 7&>~A 0) b) q A(V)9M)?, c c) I-c -C (~2) ~ d) * (4.6) Use of the half-range bi-orthogonality relations from Chapter III enables us to directly obtain the expansion coefficients a+ and A(v) from (4.6). The results are: qt= R/[ (- X(_o) 7(Z I0 (4.7) A()v= S (iV - V)Y v) A(z /v) (4.8) (-') X() d(i)o ) c a) ( /A ID) ID)[2+Co) + & 2 ] R - ci t~,lAt)dTo) b) )) (Ofp) (i/Yt~,>a)[qy/A) ~ B 2'] d) (4.9)

-51-'"cvoX(-%) (7.)..(I"7co)d'(% - ); c. 1 I a) I O d(Vo;) d(- &o)b | c #(;I c) t' ~ )i~.-r r z ) T [ 2] dd). (4.10) Notice that Equations (4.9b) and (4.lOb) illustrate explicitly the form of the half-range bi-orthogonality relations used to obtain them since the integration is over the function 8(L-io). Thus, the correct bi-orthogonal form for half-space problems may be obtained by inspection if the discrete and continuum expansion coefficients for the half-space albedo problem are already known. This is obviously true for any type of isotropic or general-order anisotropic scattering. In fact, it was this observation by Kuscer which first led to his discovery of the existence of orthogonal forms for half-space problems. The solutions of the four problems are now complete since the angular density is known from Equations (4.4), (4.5), and (4.7) through (4.10). The neutron densities and net currents, defined by f(x) tx<) dr -I1

-52are easily obtained by integration of (4.4). Upon use of (2.51), one finds that' -x/ pA () eX/ -I -x/1 f' -X/V The existence of the factor eX/V in the integral term in The existence of the factor e-X/v in the integral term in Equation (4.4) appears to make it impossible to evaluate the integral analytically unless x = 0. The incoming distribution at the surface is specified by the boundary conditions and hence is completely known. It is desirable, however, to examine the distribution emerging from the surface since the results are very important and may be greatly simplified. Using (4.4), the emerging angular distribution may be written as (~-;)A f(OA ~ + ) + fAv)ptd,,O (4.13) where, with (//>M) = /)- (o 7)) (4.14) the expansion coefficients are:

-53( /a - Y) T(u) f t /) A -IV) a,= - (D _ Jo?') d~') r~')d' From Appendix C, however, we have the result Jo (,o2-V'(v) /Wv) A-r f) I,'~ (:,o+ )')X(,)ll7")''; d(-Jol^)dD" 8> (4.17) (7 +'1(AAr') dX2oJ) X;A)

-54Use of (4.17) in (4.16) yields, after some algebraic manipulation, (~O,,)-= f(o,-) + x' I V {j)' (Vo o (,) d (') / /A >, 0. (4.18) This general result can now be specialized, by use of (4.1), (4.5), and (4.14), to give the results c~vZ X(-vo) d (V/oo C (t) - a)..........(~4v0)C ) ~ c:.3A-. I a) (VXZ-7Z) X(7) dl(i7).~~ [ d (?-B >o)] b) (t'o+j) )(/o.A)X(/A) (I - ) ) (z+/o)X(yt o) f(O,4) I r Jd(-7-7) d.(v,) d -c:I c)'~ro ~) —~- ot/ d) (4.19) — ), q4 (co, 7*) - -o' ) d(~ub) d O~'aID~bS~)) c~.~_ )

-55Equation (4.19a) was obtained earlier by Shure and Natelson.(22) Equation (4.19b) agrees with the result of Chandrasekhar(2) when the conversion equations of Appendix H are used. The densities and net currents on the surfaces of the halfspace may be obtained from either (4.18) or (4.19) after multiplication by do and [dt and integration. The equations are y(0)= 9Y(o,) + atIA5-) o 0'-(0) f(o ( ) - y ) )(4.20) o 0 and (4.1) must be used. We will proceed by considering the general equation in (4.18) and then specializing the results to obtain specific answers to the four problems under consideration, Equations (4.18) and (4.20) combine to give (k10'Y,7r')ckM2)V [j)0 d 1 (4.21) d( J X (;,AA),J

-560 0 I rl d+ (,oo) f J/, (4.22) Use of partial fractions shows that I I I Jo X() _ +- J) X3t) 0J x(y) +_% o) The three integrals needed in (4.17) and (4,18) are evaluated in Appendix C and are /I d~ cb I ~

-57where, from (2.37) and (2.59), we know that X = 9/[(i-c)(i.-b) (4.25) Equations (4.23) and (4.24) reduce (4.21) and (4.22) to ) (0) -I + 2a ~/ (, Ia~'( k(') ~>^ 2J'(4.26) c XQ)J (Q) + (0. y+A ( )0,) ) (4.27) Equation (4.15) was also used in the derivation of (4.26) and (C.15) from Appendix C was used in the derivation of (4.27). The specific results for the four half-space problems reduce to r1/2 2 [(V-)(I 3+)] X )() c.AL a) j,(o) = ad(4, o2)..)(. —.).. o dv.W/o)S J,$,) b) _I-c I c [ I

-58f(O) -t F (Oj4d, q H(Vo4V) X/o) _ 2.'AA o,/ ~') d.44a d) (4.28) 010)=.,io)wanted to solve in this chapter The solution to the half-space Green)'s ji/ a (a (o function problem, however, is not in its simplest form since it still (7) problem. G(x,pu) - G(x0,,k0 -oxp) and is given by problem. G(x~p = G(XoP -. x) p) and i s given by

-59(-x-xo)/7 /, -(x -)/a (x -xo)/I ~ -(x- Xo)/V -4_ qp_ ~)e Af-)jl) ze -/ Xo.X (4.3o) where G(~xo/t) - 0 Lj(X0,)/ - c (x;4j( ) = ______ (4.31) Using the second boundary condition and the full-range orthogonality equations, (3.101) to (3.104), G(x,,u) may be written as T(-xo (-o)/7 4. Ixtn)= (,0 hpie 4)g #v>Q/~0) aff(iiA > e d7./ (4 32) G ~ ~ ~ ~z(x,,AV') =32)) X >XO) where N+ is given by (3.104) and (3.109). For the half-space Green's function problem with the source located in the right halfspace, we therefore have'XV/V- x o/zI (0o)= - ) -~_ 0 7~ -'() e- aJ o (4.330 N_ V, A+(i'')A-(9

-60Equation (4.33) may be used to reduce the Green's function results of (4.9), (4.10) (4.19), (4.28), and (4.29) to the following forms: R C4 d, v,-,, j Xv), X (-Vo) d-yo).(,,) )e 4 d, ((ov) t _ / ~_ r (,',)X-v'),' d/ —')e d4/ (4 34) A+ (v/') A-Mv') 5 =- ~J1o X(-"o)..w dr) (,,,o) (-,,z) ~- o)e L N-_ d(o2) d()-Vo), ~~,,b',,)('V X(-r ~v['N,')+ e. d,']1 Vcz- d(Vo 8o) (o.. )] (4.356)'4- X,/( I ]_~ J - - ~,,,, o ) x f - r') e ~,,, /' + ~'/o d(-v'~) dt3o4 o)',3,, _~~dVo odf- ) No"

-61- xo/4 ( ) 1 N_ l, to). (I-c.) v _/o) (Vo +,') X(-,Ve 1)e dC (4.38) tr(,,)dtJo, A'(V') A -(v').... For isotropic source, the above expressions are simplified by integrating over duo which amounts to replacing 0(0o) by 1 wherever it appears. The preceeding equations were derived with the aid of (3.82) through (3.85) plus the following identities valid for 5 = vo or v', 0 < v' < 1 J cf w.(49 Another quantity of interest in neutron transport theory with linearly anisotropic scattering is the Milne problem extrapolation distance, zo, defined as the distance from the surface of the halfspace at which the asymptotic neutron density vanishes. The asymptotic contribution to the angular density comes from the discrete modes since e d/vo>> e/v for 0 < v < 1. Inserting f(x,k) from (4.5) into (4.12) and using (2.51) gives the Milne problem asymptotic density as

-62x/%O -X/Vo f~(x)= e + q+ e (4.40) Therefore,? (-. o) = O = e + Q+ e (4.41) or If0- ~~~~~~ 2.~~~~~ X (~,(4.42) Using (4.7) and (4.9a) in (4.42) gives aVo L- X(-4o) d(-Vo)J (4.43) It is also interesting to obtain an integral form for zo which does not contain the X-function. From (2.21) and (2.65), the halfspace X-function is given by X(S) =T 29l - J (4.44) where, for linearly anisotropic scattering, -Ir..,.d..)(4.45)

using (4.45), one may show that doe#)= ie A_ [c + 60 (I_ C )9 2 df )?(4.46) Following the analysis of Mitsis (lo) write X (-) F0 I~ l= + TL' d (4.47) and integrate by parts. Since (0) = and @(1) L 2 ]o (4.48) Using (4.4-) and. 48 ") beUsing 4.43) and (4.48), the Milne problem extrapolation distance may be written as

-64=0 ~ - dP(V'Oi) + A~J, I + I + 4, (({L~di~rg I (4.49) The value of A ^()At(+) is given from (2.60) as AtJ,4)A3)= A,.)% +rrc ~[~ t)] (4.50) where X(Ci) is given by (2.53). In this chapter, we have yet to consider the numerical evaluation of the quantities appearing in the preceeding equations. The value of y([), 0 < t < 1, can be related to the value of X(->) by Equations (2.38) and (2.59). The value of A+(Cl)A-(1) is found from (4.50) and the equation S4)[i= J#X)[ 4i57%] i1(l-c)~/A (4.51) Thus, all results may be expressed in terms of vO, X(+ vo), X(-k) for 0 < t < 1,and 7(i) for i = 0,1. The first step is to evaluate vo from (2.55) and (2.56). Then the values of X(+ vo) and X(->) may be obtained by use of the integral equation (2.70). The constants 7(i) can then be found.

It is also possible to determine X(+ vo) without the use of the integral equation (2.70). From (2.38), A()X( —o)X() = ~ A(z-) I _ A^() (4 AIov) X( vO~ x~o) = MO z -2d + io (4.52) where (3.109) gives the value of A'(+ vo). From (4.43), = - e (4.53) X(7o) where the value of zo defined by =C _=7z + d(-vi) (4.54) D - O 2 dfroi) ) is given by the integral in (4.49). Combining (4.52) and (4.53) yields the answer I3 -o Vo I/,2 (I - Co1. 1(-io)= E |/2[ (VP>[ V) ] cb, ( - c) dS(-vo) 1 (4.55) ~-,,o - ~2,, V'2 5 2' - K( ) wh ere K is tabulated in Case et al.. (5) where K is tabulated in Case, et al. (5)

CHAPTER V SLAB ALBEDO PROBLEM FOR ISOTROPIC SCATTERING The orthogonality relations presented in Chapter III enable one to solve several one-medium half-space and full-space geometry problems in closed form. Four specific half-space problems were solved in Chapter IV. Unfortunately, this program can not be accomplished in the case of a slab and, as will be shown, the expansion coefficients in the solution for the angular density are in terms of a Fredholm integral equation. The problem of interest consists of a finite plane-parallel medium surrounded by vacuum on both sides with neutrons incident upon one face. The incident neutrons are assumed to be azimuthally symmetric about the x-axis. As in the half-space applications, assume that c the mean number of secondaries per collision, is < 1. The problem requires the solution of the following equations: f 1I] (xPA )'() = xY') d<', O (2.1) where tY(J/~) = CA,) J / >'Opt (5.1) L d/)/4 = 0,'O,(5.2) The coordinate system is shown in Figure 5.1. -66

-67Cos cos 0 O d Figure 5.1. Coordinate System for Slab Geometry. This problem has been studied by Halpern and Luneburg(5O) who obtained an expression for the asymptotic density in thick slabs vV (51) by means of Laplace transform techniques. Later, Kuscer also obtained the asymptotic angular density, density, and net current for a thick slab by appropriately combining the solutions of the halfspace Milne and albedo problems. More recently, Zelazny and Kuszell(l3 used Case's method to consider a bare slab with arbitrary source distributions on each surface. Because of the generality of the problem, however, explicit solutions were not obtained. In the present analysis, the angular density is expanded using (2.17) and the problem is reduced by use of the boundary conditions and orthogonality conditions to solving two sets of two coupled equations for appropriate linear combinations, b+ and B+(v), of the two discrete and two half-range continuum coefficients. One equation in each set is a nonhomogeneous Fredholm integral equation for

-68B+(v) which is solved approximately by Neumann interation. This leads to explicit solutions for the angular density, scalar density, and net current in zeroth- and first-order approximations. The zeroth-order results are equivalent to those obtained by Kuscer. Equation (2.17) may be written as -I/O X/HvO'cx,): + e,7)) + + 5_- _-,)3) 1 AcJ) ~)e JTJ 1 4(~Y) fzd d (5.3)'0 where A(+ v) are defined for 0 < v < 1 and are the half-range continuum mode expansion coefficients and a+ are the discrete coefficients. It should be noted that either Xv(>) or 0 v(>) is nonsingular, depending upon whether i is negative or positive. Using (5.1) and (5.2) gives I -i I 4-v A( -J,~_,,~)d-J - A(v],,~)d?/,O s (5.4) ~ o~0 and - Cta + ) e' d - O4 l, ), 0I -/ A(wVL) AAerJ- j 4-Jf#tR I (5-i/ c'5) TO Of

-69Changing variables in (5.5) from L to -ii gives - d/1Vo O2/o - a ~)e~t - a_ (5M)e Al~(8)2> d =feA(-v) 4)e di)>O (5.6) J o We are now faced with the problem of trying to find the expansion coefficients from two simultaneous singular integral equations for:1 > 0, Equations (5.4) and (5.6). The latter equation is not even in the standard form as solved by Case because of the exponential term ed/v in the integral containing the singular eigenfunction V(). We therefore resort to adding and subtracting, respectively, (5.4) and (5.6) to obtain 0 0 where we have defined our new expansion coefficients as /Jo b+= + e B+jj)= A A (-m)e (5.8) Note that the right hand side of (5.7) is, formally at any rate, a halfrange expansion of the (unknown) left hand side. Indeed, the solution will be obtained in terms of the half-range X-function.

-70Equation (5.7) may also be obtained by using (5.3) to consider the angular densities?(x,0) and 4(d-x, -). Adding and subtracting the two angular densities corresponds to finding twice the symmetric density in the slab, denoted by I+(x,$), and the non-symmetric density, denoted by t_(x,>~), and gives -x/o (x-:)/o = Us+Ir)e b A) e ~,fa. b(vr)~ ~J)e d() rt it(x-d)/1 +J +(t) 2 )e 1c (5.9) D The boundary conditions (5.1) and (5.2) combine to give Y,, (oX) -; o. > 0, (5.10) Applying (5.10) to (5.9) gives (5.7). Thus, b+ and B+(v) may be interpreted as the expansion coefficients for twice the symmetric angular density and b_ and B.(v) may be interpreted as the expansion coefficients for the non-symmetric angular density. The orthogonality relations of Chapter III may be used to obtain b+ and B+(v) from (5.7). Multiplying (5.7) by h+(g)(vo - pi)?y(7)dL, integrating over ~ from 0 to 1, and then using (3.17), (3.19) and (3.20) gives

-71b[+X (-.o) + e -X(-..] a a) Use of the relationship(7'9) X(~ O) gives b=-C T o) 2 B 2)e XI-V) (l6, (_1,, - (5..3) -' " ) ~(-V) e + e eJ/~o] Here, the Milne problem extrapolation distance, zo, is defined by ~ = cio r' II~ -]tZ]4,/ ~ (5.14) + -, and is tabulated.(5) Multiplication of (5.7) by fv(i)(vo-,)7(~)du and integration over i from O to 1 gives, upon application of (3.16), (3.17), (3.18), and (3.21):

-72B+ ()(cJo-v)7 A+f-) A = (ii-Jy)"~fo) ~(,4)A -dlv' 7 br e c. X(-M) ( J'BZ('. l(v)(. +') X- v'9S (5.15) Use of partial fractions, (2.38), and (5.11) reduces (5.15) to B y~) = X()(v,;'-l) (,-') z,o),,A.) -b t~P)Xlvo),) Xt-,, Xt v Tfs+(t)ih gytt)~X/-vJS2] (5O16) Equation (5,16) ma~y a1o b@ obtained dir~eoyr without th us@ of partial fractirona, by use of a epecia eet of itotropio rlationg which do not have 4(k) and,v(4) mutually ortho$onale Th' n-e:adary relations are derived and applied to (5.7) in Appendix D. The equations in (5,16) are nonhorhOgeneosti IredholM equations for B+(O) And can be Solved by Xeuimahn itePatioh, p:Ovidtd tht retuitAnt ge~ieo Dwnvt~e es Defining Kz (~,A) = A+ ) uA') (V;+)

-73and -A) X(7A+V,)ol#)(V~c[ 2 Y/A) -)qT) 8, t+H~n)nDAt ) L'M 0) ) -a, ni t t>)+.- 4/ ) X( 7o) (5.18) Equation (5.16) becomes = ) c(l-) _ ) C (5.19) (32) Following Tricomi2 since K+(v,[p) is square integrable, the following relation must hold for the Neumann series solution to converge:. l- c I I I/2 ~c(1- 6)[tt;1<~~2( z~~i2] 1<(5.20) -d (5.21) a conservative estimate of the above inequality is

-74Here we have used the fact that ni~~~,n~~r J(r. (5.23) A*()A ) = A - ( —), (Z I )T ( where g(c,o) is a tabulated function.(5) For c = 0, the inequality is satisfied, From Reference 5, IhA AX Tr Z (5.24) for c very small and non-zero so that the Neumann series for B+() converges for c> Z e 0 113 e (5.25) Tr2 Because g(c,p) is very sharply peaked near l = 1 for small c it is felt that the convergence is much better than the conservative estimate shown above. For other values of c, convergence is guaranteed if < < 3.61 e. (5.26) Although we have no proof, it is reasonable to surmise that the Neumann series actually converges for all c < 1 and all d, Convergence is rapid only for reasonably large d, say d 2 3 ~

-75From (5.13) and (5.19) the b+ and B+(O) may be found to arbitrary orders of approximation. The zeroth-order solutions, b+(0) and B+ () (), valid for large d, are obtained by noting that the integral terms with ed/V vanish more rapidly than terms proportional to e-d/vo. For the n'th-order approximation, use the B+ (n-1)() to get the b+( and then the b+( and the first n+l terms of the series solution of (5.19) to determine B+(n)(P). Equation (5.8) may be used to find a+(n) and A(n)(+ A). Notice that this scheme gives zeroth-order results which contain a continuum component. Diffusion theory results may always be obtained by ignoring the integral terms as was done by Mitsis(9) for his lowest order approximation to the critical slab problem. For the zeroth-order (thick slab) approximation, we use (5.8) and (5.13) to find that'io+zd)lo Jiuo (0) ---— e (5.27) -C X7 VO)=l [z?+. )/ o] Equation (5.19) reduces to B+(O) = B+ (l), which gives Al f) X(-a)6(<D w)(I-(l-C) zFlo lo tt~) ) (D) C.) L' CO l(o) (o) - 9?4+)X(v0) - c )

-76(0) _O"f)~,) ef'o A (~) ) _ c)( ) I T e ()7)lX(J) A+M)A-A ) (0) -JlVo t+ e ~ _ )X(-vo), (5.28) The first-order discrete coefficients are (I) {D) c I[ ] Io) = I* o Z Z o i + CAo ZZO + eTzo,Io J] (5.29) where Int and nm are defined as I _ Ie(,-G r, e XZ-(,) bI7 = 0,, 2,3 Imn _D (3)n A+f) A-(i) 0efr -odrcniumcef eI, -J/ z() hr- A(I-c)..... Ah.... (I.3o) The first-order continuum coefficients are

-77All = ^(oD*) C__(__)X__C)_(MO_ _ )_ (X + e +,!'_ e(v[ cI,(,) -)( 2) (x)r - 4 C+ e(1 l)[ -,, ])OI (5.31) +(O Xl-v,)[~_,, -(I, -c) where the A(O *)(+ ) are given by (5.28) with the a+(0) coefficients replaced by a+). The zeroth- and first-order solutions of the slab albedo problem are now complete since all of the expansion coefficients in (5.3) have been determined. The density, p(x) = f 4f(xi')d', and 1 the net current, j(x) = Jl I'(xW,t')dp' are easily obtained from (5.3) upon use of the normalization condition (2.3). The equations are - X/Ho X/7o,p(X),+e _ e -x/- X/ +f A( d),+ A. v)e d S X (x)= (Ic) e a

-78The expressions for the reflected and transmitted angular densities on the slab surfaces may be simplified by the use of Cauchy's integral formula. The results from Appendix E are (for l > 0): X(7) L o V'/ + O,)(-)O4 ) XrJ AI) - X(^ 4[+t (+)X&.) - _e. 9 )X(V) t+A MV)Q-)XX (5-33) In Appendix H. these results are used to correlate the X(-i) function of Case with the X(i) and Y(4) functions of Chandrasekhar.(2) Equation (5.33) was derived in a different manner in Reference 33. When d is large, asymptotic expressions for the above angular densities can easily be obtained since e-d/vo > e-d/v and the integral terms may be ignored. Use of (5.12) and (5.27) in (5.33) gives the results - ~io [ I _ o 4(a[(Zo+d)/vD] -J, A __ _ _ ___________o) (5.35-)

-79The reflected angular density for the half-space albedo problem is Jay (2 t) IDgtl) X(-w) - (5.35) (co) X(-> )(ub +o) ( -C) Xo(z ),A) and this agrees with the result for isotropic scattering obtained from (4.19b). Equation (5.33) can be used to find the scalar densities and net currents on the slab surfaces. From Appendix E: Ote X O+ [ X(-to) + A/ X(-v)] a (0) Z 144j' i 2X(-') -aho d/ldo - a'0g (~)'~ D~C~a + =+ -JZ X1*4O) + a_ 2 XT-a/);.('d) - o X( - v. at.+e,Z0X(-) Q_ A (ux) +(A (v) w XI- v) Ji (5.36) - ~1,'o0

-80Since X(O) = [vo2l0-c)] - 2 use of (2.38), (5.12), and (5.27) gives the asymptotic forms of the above equations (where the integral terms are ignored): 0o? D coI[(Z~ +d)/ s p d t /- C (U~-no X ) x(;"o) 7(/0Z 0) X(>O) These results can be specialized to obtain the following half-space albedo problem results: n f(d) - o d-2 -A;Gcr~ a-2)= o,(5.38) 1440

-81The results in (5.37) agree with those obtained by Kuscer. (31) Kuscer observed that the neutron density in a thick slab could be written as a combination of the asymptotic densities from the halfspace albedo and Milne problems, denoted by pa(x) and Pm(x), in the forms P(x)l = Ng (d-x) x-d (5.39) g where M and N are constants which must be determined. The values of M and N are obtained by observing that the two equations are equal somewhere in the interior of the slab. Using the half-space results from Chapter IV with isotropic scattering (where d(pq) = 1 and y(0) = 1 ) in conjunction with (5.12), we have..... - o o/ - X/ [ -2.o/4o -x _ -Nzl(- 0) - 1( x)/ (5.40) o) e 3e~~~~~~~~~~~~~~

-82Comparing the coefficients of e-X/Vo and eX/Vo yields the result (2 o+J)-,,O r A0) czo X(-VO)AA4([(Z20+4)/1] -/ (5.41) Z o xV.Xo) (HZ tA )(1 - C) Xto4o)At;4[( Zt0+ d)//] One can immediately obtain the results of (5.37) for p(O)!large d and p(d)llarge d using the isotropic forms of (5.28). Kuscer also noted that analogous relations are valid for 4(x,p) and j(x). In fact, the results of (5.34) and (5.37) are obtained by using (4.19), (4.29), and the equations,'o,:t' Yo: -.. - M ~,(:,~ -r) t'tJ,,)|I = N' (o, 7) 0) - (~) - _ M a, (~) i(e)|, - -N / (0 (5.42) (The minus sign in the last equation is necessary to take the reversal p -, -p. into account.) The zeroth-order results are obtained from the "large d" results by the addition of the integral terms. Use of (5.28) and (5.30) to evaluate the integral terms yields the answers

-83f (0) -k?o _ 0o o I [(Z+d)/,o] 9b d v1C X(; o) (o l l l);o Y o -~ oI0o ~I~)(d)Of(d)| + V i7(>4I ziX(,Mo)'(0) (0)- ~ ~ (dd(. 2..d)/.] (:0),,"o) o In order to estimate the correction terms in (5.43) due to the integral terms, we need to asymptotically evaluate In o and JnO Use of the identity,ta % = 2 I,"+y and (5.23) enables us to use (5.30) to write = c(I-c) f/. n e X2(v) Ji = 1)2, 3, (5.44) /-/ Z.

-84For c = O or 1, In_. For any c, the major contribution to the integral occurs for v 1. Following closely the work of Case et al,,(5 set v = [l+x]- and observe that I J. d~l-JX~( )') el ci.~ -~ X2(_ _ _ _ ____ex _ (5.45) This may be rewritten as where x 2e2/. We therefore have -d 0 #0) (5.47) where f)e 9 Y (5.48) and X = ln xod = ln 2d - 2/c. The function f(\) is tabulated(5) and its ultimate asymptotic behavior for all c is(5) _ _ rc 254 (59

-85Equations (5.47) and (5.49) may be combined to yield the asymptotic result I X)I,,I —)X(_i) e — 2,54 )C ~ ~ d)Z nI O,2 (5 50) The evaluation of JnO may be achieved by using (2.4) and (5.30) to obtain ho = ( I-C 6.)tm0 e X +-/#oM)(Iko- 0) )) S (5.51) where h+i -d/V, S c (I-c) J i.) e xZ( V) 802) ad (C) (5.52) For c = 0 or 1, S: O0 If 1o is small, the asymptotic evaluation of (5.52) can be carried out in the same manner as before since the integrand is not singular over the region where the major contribution to the integral occurs. After setting v = (l+x)-1, we have ( C) ______________ _ -dx C.) )( ) [e-ol +) I~ — (oA (5.53)

-86If ~ << 1, then 1-[o _O__X oX/(1Ao) I __J I_-/o e I-/Mo(I+x) I>-/o I-/~o so -' (I X'Z (VoJo - 0. e-Dx dX 1TP,. (5 54) where D=- ~ (5.55) I -iwo If D >> 1, the integral in (5.54) equals D f(A) where A = in 2(d - o. One therefore has ____",. D -, ce o, (5.56) If) is such that the above restrictions are not satisfieday If pLo is such that the above restrictions are not satisfied, it may be written in a form convenient for numerical computation as I = C/>o)E dJ Wv(57) /.o + -W o) 5

-87where W(v) - Xc)l-)) X..,+l -d( (5.58) e, XZ( —L/) (5-.8) In addition to the previous results already derived in this chapter, a few specialized answers may also be obtained from the analysis. We will now proceed to the determination of: 1. The expansion coefficients for the half-space albedo problem for isotropic scattering. 2. The transmission and reflection coefficients for a thick slab. 3. The current through a slab when c = 1 4. The expansion coefficients for a bare slab with a constant surface source (plane source) on one face. 1. In the limit as the slab thickness becomes infinite, all approximations reduce to the zeroth-order approximation and, from (5.27) and (5.28), the expansion coefficients are: (0) A ) O c74_ -X,,) )

-88These results may also be obtained from Equations (4.7) to (4.10). 2. The transmission coefficient for a slab is defined as _________ ((d)( 6 T = =(5.60) and the reflection coefficient (commonly called the "albedo" in halfspace problems) is given by f -i0 J()O) (5.61) The values of T and R may be found from (5.43) for a thick slab. 3. When c = 1, the current is constant. From (5.27), (5.28), and (5.32), we find that T (ero) = s,lt (5.62) This zeroth-order result for the current agrees with the answer obtained by Kuscer(31) when Equation (H.5) of Appendix H is used. Using (5.29) and (5.31), one obtains Y, (X) _, (5.63) CJl (Zo. d) X/- o)

-89where the correction term to the zeroth-order approximation is given by I -d/9 K= *.e (5.64) C->l [ 3 (2 D ~J) C (I- C) 4. The expansion coefficients for the bare slab with an isotropic surface source on one face may be obtained from the previous results by integration. Using the Equation (D.7), we observe that I' tk~o ~Ad) to t Cot At od = ~ (5.65) Also, from (2.36), 0)" j ) = (5.66) J0 Integrating (5.27), (5.28), (5.29) and (5.31) over k~, one therefore obtains _(0o) (2 ao+ 2,)/ eo (o) =e (0) __ At =-e ci0X( 7zo)X>[(2to+ei/9, o X(- Avo ) I [R z) / [+;)X4)) to) + a_ p_)xi-v S/Vs~~

-90A ( p)(r) _-e / X (/\ (2o7A+(,m )A (,4 ) [iX o) e - (0) ] 9+~+% e cvO) = I +,o I[ 1#~ a-(') ) A(03~ Aopz-) ( ) (xa-) (xi) e X(-vo) - I A (-r) = A (J^) + e,o ((-+) xe )Tho~-~~ t) A [qX()[A: + ^,) + c_ X(-) /oI)[, Vl-:I1 (5.67)

-91The analysis of the slab problem in this chapter can be used for different boundary conditions. In Appendix F, other slab problems are considered where there is an incident distribution on one face and a specularly or diffusely reflecting boundary at the other face.

CHAPTER VI QUADRATICALLY ANISOTROPIC SCATTERING IN A NONABSORBING MEDIUM As was mentioned in Chapter I, problems with this type of scattering are not of much interest in neutron transport theory, but are examined in photon transport theory because the scattering function is a generalization of the Rayleigh scattering function.(2) Because of our assumption that the transport equation and boundary conditions possess azimuthal symmetry, however, we shall derive only the azimuthallyindependent part of the angular distribution, i.e. the part which is obtained after averaging the complete solution over the azimuthal angle. Some discussion of the transport equation subject to azimuthallydependent boundary conditions is found in Appendix F. Using (1.4) with N = 2 and c = 1, Once again we use the ansatz V(Xc) = -,r)e (6.2) to separate the variables in (6.1) and obtain c(i7);) )= t +;[/AP + bZ lo 4 P. - (6.3) -92-

-93 - Multiplying (6.3) by dp and ALd~L and integrating shows that (6.4) and so fPu)9(e)4,d Lft>'d)J (6.5) Therefore, (5.3) reduces to ();z ) 14>,, =, 2 >+' (6.6) where f(i) is defined by the equation We again impose the normalization condition Jt+)JRq.,,o =t ( dr =((6.8) which enables us to write (6.6) as At2P) = p;/r+ +;,)i,) /). (6.9) The value of \(v) may be obtained from the two preceeding equations and is _R'AA)d, (6.10) Z~~'-A

-94Equation (6.10) may also be written as \(i)= A ()+ {i) A(6.11) where A_______ 4;1,4 ~~(6.12) and where A(z) is given by A(z) = - - (6.13) The equation A(z) = 0 is satisfied only in the limit of large z since Ath)= - Z + +32. sztLI [- bS ],(6.14) This means that there are two roots at v = X in a manner similar to the case of isotropic and linearly anisotropic scattering with c = 1 The two linearly independent discrete eigenfunctions of the transport equation are found directly from the Boltzmann equation as Lf.jIXV) I 4[ _ w-X (6.15)

-95Since all the eigenfunctions of the homogeneous transport equation are now determined, we can use the completeness theorenm) to write the angular density of (6.1) as t8X) 2 9r " (Xv,,.JAtV)X ) e. (6.16) There is an additional difficulty encountered in the special case of c = 1 which we have not previously considered. Substituting (6.15) into (6.8) shows that the normalization condition is satisfied by *l(X,.) but is not satisfied by?2(x,IL). Therefore, we are still lacking a normalization condition for V2(x,u). From Shure and Natelson(22) we realize that when c = 1, the steady-state homogeneous transport equation satisfies the equation 1 &=J/(X }, = CONSTANT J (6. 17) Thus, this is the normalization condition we seek. Using (6.16) and (6. 17) gives I I 2 =1 3 - b',,r 9 0J1(X- 30*(6.18) By use of (6.4), this reduces to By use of (6.4), this reduces to 5 _ _. (6.19) 3- b,

-96Now re-examine Equation (6.1). It may be rewritten as 4 6,J? ~ (6.20) It is necessary to make (6.20) a homogeneous equation, This may be done by the substitution of P/%x)- =_ it) ~ 6,7x. (6.21) which reduces (6,20) to AX l]t(X)= Lj [14 + Zt6 ~ ~ I x/)6/ (6.22) Therefore, the transport equation for quadratically anisotropic scattering with c = 1 and b1 not necessarily zero may be reduced to the same equation with b! = 0 by the use of (6.21)o The two problems we will solve in this chapter are the halfspace a) Milne and b) albedo problems for quadratically anisotropic scattering with c = 1 and b1 not necessarily equal to zero. The boundary conditions for the angular density, 4(x,,u), are ( 0 )^,0 a) |oi~),/A,' ~ b) (6.23)

-97|2 i- 3- )AS X -o a)* FINITE AS X —-o b) (6.24) where 4(x,11) is given by (6.16). Applying (6.24) to (6.16) shows that V(x,p) may be written as t(x,~- x(~+ + fJ A(i)9s')e. h) (6.25) 0 where xt~~ (X ~ )a) o b). (6.26) The net currents are therefore given by 3- 6, a) j=- 3L b)** (6.27) *The slowest growing mode for large x is 4r2(x,), arbitrarily normalized by setting a2 = 1. **Perhaps the logical order of (6.24b) and (6.27b) should be reversed. Indeed there is no a priori reason to demand that * be finite as x -e A. On the other hand, since every neutron entering the medium eventually leaves at the interface, there is indeed no net current. The discrete mode ~2 is the only one which can contribute to the current so the boundary condition on * as x -> co is selected to make the net current vanish.

-98If we consider the solution to the transport equation where b = 0, Equation(6.22), the angular density is now given by' (Xiy)= X(/) + 3 6 ) dJA7/) J(6.28) The V'(x,0) is related to 4(x,9t) (for the case where b1 is not necessarily zero) by bx a) qd(OX,/) )a( (X )) b~ 4) (6.29) where (6.21) and (6.27) have been used. Now that the two problems we wish to solve have been formulated, we will proceed to the various identities satisfied by the X-function and obtain the orthogonality relations which will enable us to determine a1 and A(v) from the equations 4( ( b)t'~> a) $6O & o), /,,,, b).(6o30) When comparing the work of this chapter with that concerning linearly anisotropic scattering in Chapters II and III, one observes that the polynomial f(p) replaces the polynomials cd(v~) and cd(^t) in the prior work. Therefore, in analogy to (2.67) to (2.69),

-99 - the half-space X-function of Case (which has the properties attributed to all X-functions in Chapter II and vanishes like z-1 as z -, o ) satisfies ~~x~~I i'(ei*) = it(6.31) (-I) h(j)= vP p ~ tf~ )d', 4~ o, ~,~ I (6.32) (6.32), = ()f~ ) (6.33) where _-ci) V XYV) O-{>J 84'~ (6.34) Equation (6.9) may be used to rewrite (6.32) as JO )t>)9< = 0. (6.35) Also by analogy, Equation (3.65) becomes () tI )V = 1/! [0 j VID) J Z ) 94)] +A nt+r)/\-p) sl*n-8 zi~- A-1)m ~(6.36) Since f(p) has no zeros in [0,1], Equation (6.36) may be multiplied by y(t)dl/f(G) and integrated over ~ to get

-100( 1 QEI) A*()a-,) A (V-V-) (6.37) where (6.35) has been used. This equation is the normalization relation, with weight function y(t)/f(k) between the continuum modes. This is no guarantee, however, that the continuum eigenfunction Xv(p) will be orthogonal to the eigenfunction V1(0,IJ) = with weight function y(0)/f([). In fact, by use of the identity f/ () __, -) +. 3 (6.b8) (63 8) it can be shown that ~,2 it,) = ) (6.39) a0 f) 8(v) In (6.39), we have again used the definitions r; j6 =t~ r t)s, = 0, 1 (6,40) and -o (6.41) a (O)~

-101Equation (6.39) may be rearranged to give!+ a8 f(,;) J # 0 (6.42) Use of (6.35) and (6.37) shows that?,(v) A4v,) A -(i) 6(,- v - (6.43) In order to completely determine the continuum mode coefficient by use of the above relations, one also needs the equation -o z 8,[,, 3b'.,.r,,,) ] Ma, = [r[o, f() )] ( 4) 4] (6.44) (03"'0' fly) 4-'~'~4Following the approach of Shure and Natelson, (22) the right hand side of (6.44) may be simplified. The X-function under consideration obeys the equation (z) X(-) = A() X 2(0) (6.45)

-102where _ -) _+zb -i ZA(3 I rr _ A(6.46) is obtained using (6.14). The equation f(G) = 0 has roots I= + [I+T ]+ - + 4- - ~bo From (6.31), X(+-) - - I+z1 r) ( ~ o)d = 3 l (0) (6.47) 1Using (6.45) through (6.47), we find that 3A( () 2l) 21 )-ab.(r ) (6.48) Setting z = a in (6.13) gives a. _)) _ _b~r enz6 (6.49) so that (6.48) becomes i 6i. [ ](i* ~ ) - 3b4 [4(aX)] [#41 ] f( i) (6.50)

Thus, Equation (6.44) becomes z~~~~~J ~(6.51) 4- 7 (~) ( Ir)((- 62/5) and the four equations necessary to determine the expansion coefficients in the half-space Milne and albedo problems are (6.35), (6.42), (6.43), and (6.51). These expansion coefficients will be found in terms of b2 X(->) or Y(~), and 7(i) for i = 0,1. The function y(G) can be related to X(-k) by use of (6.34), (6.45), and (6.46); the result is: God)_, 3/^ 0 1 O f <(6.52) We now apply the orthogonality relations to Equation (6.30) to obtain al and A(v). Using (6.35), (6.30) reduces to 3-b 3 - b, 2 To0) b). (6.53) Equations (6.42), (6.43), and (6.51) reduce (6.30) to Equations (6.42), (6.43), and (6.51) reduce (6.30) to

-104Z (3- b, T () A(v)-/-(/) A to) (+ e o)" )f(.....b) (6.54) T(V) A+() A -(v) so the angular density f'(xp) for quadratically anisotropic scattering with c = 1 and b1 = 0 is now completely specified. Equations (6.28), (6.29), (6.53), and (6,54) give the solutions to the halfspace a) Milne and b) albedo problems for quadratically anisotropic scattering with c = 1. As was the case in Chapter IV, the integrals for the emergent angular distribution from the half-space may be done explicitly, From (6.28) and (6,29), (0}) V=?( (ot)= f(Q), AA) + +' (6.55) where, for he) - p (o,nr -ffI) n(6-56) the expansion coefficients are

-105AfC() -V)')_ 3b(6 57) ^i, _ A/(,,) A-( [ ( Using (6.57) in (6.55) gives (,~o,) - = J(Q,) + Io )?rv)....... x{Yyk,)+ t-3-lib28f />>~' Jd',*>0(6.58) From Appendix G,

1 063 V b' z,,') (v +,-OM70 /Jo VTg~)J(v)[ 9a(ty) A 3Jbf(/)b' +Y) I [ T 3[ -/1- 4 (6.59) -Le') (,~, + x'-~), and so (6,59) simplifies to 7(o, y) -/(o,-7) X, ^)O' [ 4 f c]) Since a) z (3 -bI) kb)1(6.61) the results of (6.60) can finally be reduced to

-1073; 0a) 2(3- b,) 0)-, X/,')(l- b2 /5 ) b[1o- 3bz 3/bo)..t..0 b)(6.62) 2 o) X b)(I- 1/5) X(;o) where (6.50) was used in the derivation of (6.62a) and (6.52) was used to obtain (6.62b). Equations (6.53a), (6.54a) and (6.62a) were obtained in a different manner by Shure and Natelson.* Equation (6.62b) agrees with the result of Chandrasekhar(2) when the conversion equations of Appendix H are used. To find the surface densities and net currents, we can integrate over (6.60) or (6.62). The net currents are already known by (6.27), however. We will find the surface densities by integration of (6.62). From Appendix G, the integrals needed are *Shure and Natelson(22) normalized to J = -1/3 instead of J = -1/(3-bl).

-o08he =* 2(I-b2/5) s 3b- b5)Z~lba/5)[, X( e 0) + 1'2 X'(' [I+ 4- + o - 1, _ 2a (6,15) X(6o) a x2tto) 3 140 (I + bZ/4),r)J - 4 (66 63 I-o(/ + 2)Z/4) The surface densities are given by f70 =JT(}M1 J i y2(o ym)+ (6,64) and, using (6.23) and (6o62), are found to be /I - b3/ - 4. )b2 (3 - b,) -'() (I + bz/4) x o) ( / (6.65)

-109For x ~ 0, the density may be found by use of Equations (6.8), (6.28), and (6.29). The result is Ix +I - x/iz b, x J x a) 3- b) o ~ b). (6.66)

GENERAL ANISOTROPIC SCATTERING The transport problem for a medium which scatters neutrons in an arbitrary manner was considered by Mika, (21) who proved the completeness theorem and derived full-range orthogonality relations. For this case where the scattering is N'th-order, Equation (1.4) is to be solved. In this chapter, we shall proceed in a manner similar to that of Chapter III and attempt to determine a bi-orthogonal form which makes an eigenfunction of the homogeneous 3oltzmann operator orthogonal to an adjoint function with respect to an appropriate weight function. No conclusive results, however, are obtained. The calculations are presented mainly to suggest a reasonable approach and to illustrate the manner in which difficulties arise. Using the ansatz of (2~2) in (1L4) and defining 5A (v) -/ Pi ge)+) (7-1) one obtains ('A-k) -P) = 2-E, QL)j b (7,2) AkO Multiplying both sides of (7L2) by (2+lj) P,(i.), integrating over and using the orthogonality and recurrence relations for Legendre polynomials, we find that the Ek(v) are themselves polynomials with the re currance relation -110

-111where h Ze + { cb (7X4) Again we choose the normalization So(Z) =. (7-5) One then obtains (i)= h = (ho _ ) 4SL~) -- i(oh h+. - {C o~h;,h. ~. h34. ~ j, A V1,1(7.6) A-3~~~~~~/ A /,

-112where the general term has been given by Kuscer. (3 T) The general term may also be written as [R/2] ~_, i!-262~' =,- I. (7-7) where I1 = 1 and where [l/2] is defined as the greatest integer in i/2. The function Q(v) is a polynomial of order i and is a generalization of the Legendre polynomial P([i) to which it degenerates if c = 0. Using (7.6), one has gRy(8) =X>$ jQ vtt~29 (7.8) where Abo - I (7.9) and where, for i > 1 "o ho h- +'+_-_o'and-4 h, _z (7.10)o o and so on.

-113 - For the general-order anisotropic scattering, define D(,p) =- (7.11),0 0 The eigenfunctions for -1 < v < 1 are At+=AVz p Z. r + aim) > )(7.12) where x(E)-= I- Cj Pf1 DI) J (7.13) The discrete eigenfunctions are 7mp) i ii ii i (7.14) In Equation (7.14), the eigenvalues + vj for j = 0,1,...M are the roots of A(~t 7 )- 0 (7.15) where Alz) =- I -' (7.16) z. -. ~~i,,' ~i~. ~~

-114and where M < N (8) Equation (716) is proved by showing that the integral I1 Sn~t(- SRIX) _ iA) XKR 4N vanishes. Since e(z) is a polynomial of order < Q in z, the integral reduces to a sum of integrals of the form which indeed vanish. The values of A(z) as z approaches the cut from above/below are At (v) - (rJ,) ~+: 1' Dtyv) -I I (y,,) Al (V + L-I eve ~ (7-17) The behavior of A(z) at infinity can be found by considering QN I nc) = R Ll,P- 2 J'I # ZF ZZ+ ]4- (7c18) The value of A(oo) is finite and is given by A(o)= h o -C h_, CbTE _ h, hL.,.(N Tr RM5) (70 19) 5 33 3 5,7

-115By inductions it: may be proved that(35) A/() = ] [a+1u] (7.20) L o For the case of general anisotropic scattering where N is finite, the solution to the Boltzmann equation becomes, by virtue of the completeness theorem of Mika: (21) M X/ M X/ + A (v) i)e dg (7i21) JI The function Xg(z) is again given by a linear combination of functions analogous to the forms in Equations (2.26) and (2.62): PZ) =- (a -a )W' ] ) 4(7.22) where Xgo(z) has the form of (2o20) and where the constants wi i = 0,1,2, depend upon c and bi, i = 0 to N. The integer constants Pi and P2 are found from the relations

-116o A+(.) -AL, AVa) I o06 +Pi, 4i06- Ale {) -.A () p < (7,23) since [,. A A+(d) —,on A1A).9(e~o~t)'~. ($- ~:) as ~ -^ —}, Xgo,(i) (e-t3 )AS Z e t'(7.24) These value of Pl and P2, for the special cases of interest, are given by: (X B P1 P2 -1 1 M+1 M+1 0 1 0 M+l -1 0 M+1 0 Two equations of interest for the N'th-order anisotropic scattering are @'

-117() - Xa(X_ ='__P_ a _ D_,_ _ (7.26) The particular relations in the half-range case, where wo 0 O and w2 = 1, are: X(z) - (7.27),vVrt) hlz/)?P 1)4) D>,r) " D74/ (7.28) J Ai) ) dr =. (7.29) The half-range X-function satisfies the identity =(-A() (7y30) M which means that cJ T(7J)= -, o o __ V (7.31) X(-v) Alo) (5,Z-vZ) i=0 For C = -1 and f = 1 the full-range X-function may be awritten as X(~) M A (7.32) AI)R ( V.2 ) 0':D

-118so that T () 2= - M, -I ~ y I (7033) A(O) ] (7.D_)V2) a=3 1 Other full-range relations corresponding to w 0 and wl w2 = 0o are: X(~) ( ) X { )'D (7354) PfT(/,)LY/4sW4)dt11 D(7-35) fi )DI x = =0 (7.36) /J) VDA/) = I (7037) In the present form, the preliminary analysis indeed appears formidable. It can. be greatly simplified, however, by assuming that there is only one pair of roots to the dispersion relation (7o15). The (2) (6) vv (36) Chandrasekhar, Davison, and Kuscer. Although it is tempting to argue that there are only two discrete roots for general anisotropic vv (36) scattering, Kuscer has shown that this is not necessarily true. However, for scattering thich is not too highly anisotropic, the assumption of only one pair of roots should be permissible and we will therefore make this simplification. This assumption

-119is reasonable since, according to Davison,* it is true in "all practical problems hitherto considered," In the full-range case, one verifies (21) directly from (7~2) that the orthogonality relations (3 1) to (3.3) for isotropic scattering remain valid. We proceed directly to the problem of finding orthogonality relations for the half-space problems. From our experience in Chapter III, we shall assume that the weight function for the half-range case is still 7g (4) = (Vo-p) y(4) o This assumption, however, is only a guess. Following the approach used in the case of linearly anisotropic scattering, one can prove that (9- v9) DlY') O~)(i-,)%)fr=O, 0 (7.38) provided Xg (m) f 0 0 The rearrangement of this equation into a generalization of the useful bi-orthogonality relation developed before is not as straight-forward since D(v,4) is not a symmetric function except for the two cases previously considered (N = 1 or N = 2 with c = 1). Thus it is much more difficult to:..remove the polynomials D(v,4) and D(v'1,) from the denominator of (7.38). Instead, we will approach the problem in a slightly different manner and examine a "recipe" for finding a bi-orthogonal form for half-space problems. From (7.12) or (7o14), (I->~ =,,'- D (7 39) * Reference 6, p. 242

-120arid writing the equation again for ~v, and operating on the two equations and subtracting gives (YVI D+A),) DI9,k")](ATR* (7.40) In order to obtain the desired bi-orthogonal form which will make the right hand side of (7.40) vanish, add the term 0 to both sides. Here, the m in B(mn;v'lt) denotes the mode of v in (7o40)o m = c for 0 < v < 1 and m = d for v = v. The mode of v' is denoted by n and n = c for 0 < v' < 1 and n = d for v' = vo0 Note that B(cn;v Il,) can not contain the variable of integration in the eigenfunction expansion (7.21) However, B(mn; v?94) can be an explicit function of v~ and of p, as is denoted by the notation selected, One nlow has ( provi de m w'p ) i constrj) + ]u )ted s = (7.41) provided B(mn;v',4) is constructed such that

-121-, V NP.(A) -(r)(-) 6No Z +CV (;28' citB~niTV) 0H)%<=.(7.42) - = 0 Equation (7.11) has been used in obtaining the above equation. To simplify (7.42), define ~jP)X/I<)( H-in ^) ~iS^ N R;^ ) i7N (7.43) where Ki,N(V) is a polynomial in v of order < N (v may be either a continuum or discrete mode). Equations(7.42) and (7.43) combine to give Z. E [(VI);, ), ( IN (X] izo +(Z~ % 2* 5~iJ2>j^),2vt)(~aIO/)r) 0, (.4

-122After dividing by v-v, this may be written in the form t 8 1vn n j >//^) go) K( N-) 1< (7 2 wfhere LN (,) = L=o.L (7,46) The function LN(v v') is a symmetric polynomial of order < (N-i) in v and v1 so therefore B(mn;v pi) must be a polynomial of order < (N-l) in v' if (7.45) is to be valid. In order to simplify the solution of (7.45), we will. try to expand B(mn;vl,,) in the form N-I 13(muj hZ>2^);>. (n7) (7.47) where Z = 0 and where Qi Ln; ) =I eidlnn)vl,(7.48) Use of (7.47) reduces (7.45) to the equation N-i [ —O

-123after (7.43) is used. From (7.47) and (7.48), the final form of B(mn;v',4) (such that (7.41) is valid) may be written as N-I N-I (m3 (r)i;n ) =Z PL~)j jijn) (vm (7 50) where the values of qij(mn) may be determined from the equation N-1 N-I I KIV)j l(r) ) La (z') (7.51) provided they exist. Notice that no where have we proved the existence of the qij(mn) which may depend only upon m and n, c, v and 7) and bi, i = 0 to N. We will merely assume that they may be found from the above speculative analysis. In the solution for qij(mn), the left and right hand sides of (7.51) may be grouped in terms multiplied i j by v1 and (v:), i, j = 0 to N-1. Equation (7.51) then dissolves itself into N2 equations for qij(mn)o Because of the symmetry on the right hand side, however, we may judiciously group the equations into N systems of N equations each. Thus the values of qij(mn) may be found after evaluating N'th-order determinants. The need for evaluating N'th-oarder determinants when considering the transport equation with N th-order anisotropic scattering and arbitrary c is not unique to this particular method of solution; it was also encountered by Kuscer in his work on the Milne problem for anisotropic scatteringo(3)

-124Using the bi-orthogonality relation formulated above, one can now calculate the normalization condition for the continuum eigenfunctions Where 0 < v < 1 and 0 < wv < 1 An equation analogous to (3~65) is I [ Zl, dOrln) -Z'\^+A )\ ^-),'-v"),e-,').(752) Using (7043) and (7,52), rC (- -',) T~J)/ lh/v),(-'6, (7, 53) provided (7o50) is valid for m = c and n = c As can be seen from the above analysis, a first step in the solution for qij (mn) is to determine a general expression for KiN(V)o Using (7,11), (7o12). and (7.43)9 Z: K C I-)- P: (v),o-,) (,J,),v) N tIP' ) P () )) L (7"54) Qn O

-125In order to use (7.28) to eliminate the first term on the right hand side of the above equation, we utilize (7.8) and also the expansion(37) Px 5 = V (2') -z (7055) J~=o to obtain _____ - t.([Q/2J R-2 — -s( )?(L/ -2Az)! E-2A1I-r r (7.56) (2' -!) ( -k)! (-ZA)! where Z = 0. Use of the equation and (7.54) and (7.56) gives the general form for Ki N(v), 0 < V <1

.N [,z] R -2,-i + (.v, | - (i b IV 5 l=O t=O - -5 N~ L/ -2 - V ~ ~~ (7.57) I: tk=o r —oZ 0 4V I Notice that (7~57) is also valid for v = vo since the result (7.58) is also obtained by use of (7.14) and (7.43). This is a consequence of the selection of the half-range weight functior as (vO-P) 7(U) Specific values of KiN(v) may be obtained from (7.57) and are:

127K ( 7v) ID) K0, 2z) -= o,, (L) _ 42 [ _- h.J-h 1 I (V) (= K - ( -)~ [''h).Oi+ (0)] bz, 3 T K ) 2 [ 3 1Z) )] 3b, JbD 8 ([ )b0A 3 1 ( 3 z) ) I ]96 ( T' ( + ~] (2- Y)[3o Xl0) 3l)?d() +?7/h'' bI,) IO) e (+ [; ~ },-,] [~,~ v' (,) )

-128Because of the complexity of the Ki N(v), it is desirable to determine a recursion relation to check the polynomials obtained by use of (7 57). Using the recurrence relation for Legendre polynomials, multiplying by 6v(l)(v -4)y(7G), and integrating from 0 to 1 gives (n 2l)Z K4-1 N nC( (Z 1 ) oth)(uE07^)t>) TV I An4 -l ()= 0) I0N-I (7.60) Using partial fractions, the middle term may be reduced to so the desired recurrence relation so the desired recurrence relation, is (i +I),- N(V) l (z ) v k'+,(-) + k' I+,,, d) = l (7.61) Equation (7,61) is satisfied by the K iN(~) of (7o59). Up to now, we have only shown how to construct the qij(mn)o We have not yet said anything about which orthogonality relations must be constructed. From (7041), we see that for 0 < v < 1 and 0 < vl < 1)

-129then (9v)J 9 [ t ( vI /) i +] (7/o i )* 0~ (7 62) and the normalization condition is given by (7,53). In (7.41), if v = v 0 and 0 < v < 1, then' c](oT/)I l (7.63) We hope that -(cc. jv'<) = B' (Jc; ) (7.64) so that both the continuum and discrete eigenfunctions are orthogonal to [EvG('i) + B(mc; v',4)c'] m c = d. If the weight function (which 2 we assumed was given by (Vo-l)y7(l) ) is correct, then Equation (7~64) will be valid. This would make the determination of the continuum mode expansion coefficient for half-space problems trivial after the other necessary relations were found, Even if (7.64) is not true, however, the present analysis could be used to find B(cc; vt, 4) so that the determination of the continuum mode coefficient using (7,53) might still be possible after the integral rl -)

-130were evaluated. It must be noted that orthogonal forms are sometimes useful even if the discrete eigenfunction is not orthogonal to an eigenfunction of the form OV ( +) + B(cc; v', k) 2 For an example of this, see Appendix D. Another useful bi-orthogonality equation could be obtained from (7.41) by setting 0 < v < 1 and v' = v to obtain the result The co[nstrution oi B (cVd vo ~ ]) w = t, (7.65) The construction of B(cd; vo, p) would permit the discrete expansion coefficient in half-space problems to be determined after the other necessary relations were derived. In general, since Ki N(v) for v = vo is an easier expression to work with than if 0 < v < 1, B(cd; vo, ) will be more easily found than B(cc; v, A). The asymptotic distribution in any half-space transport problem with general anisotropic scattering, which is obtained from the discrete mode expansion coefficients, is therefore more accessible to evaluation than the non-asymptotic distribution, One advantage of the present method of construction of the bi-orthogonality relations (which could have been used in Part B of Chapter III) over the construction which would have to be followed using (7.38) and proceeding in a manner analogous to that of Chapter III is the following result:

-131Theorem. If B(dc; v, ~) is constructed such that t tt)[t B(dc.i/) ( Z ](-dA)a C)lt 0. (7.63) 0 then it Rl)[p~~, t)c~ + BLd% ]=a~-,)t~)J =~. o(7.65) The converse is also true. The proof is trivial and consists of showing that 2-.'llt B(62i Sov~) (idles)@>)+ (7.66) Use of (7.50) and (7.43) reduces (7.66) to N-I N-I Z K, ~) t- ( ) N-1 N-I 7> (7/)E T )./ (7.67) =o 0= Equation (7.51) can be used in (7.67) to give the equation Lo,{ (2O)V) = LN (1~) /(7.68)

-132But from its definition in (7.46), LN(v,v') is a symmetric function so the theorem is proved. An example of the theorem is seen in Equations (3.80) and (3.81). We have now formulated a "recipe" for finding bi-orthogonalit relations of the form of Equation (7.41). In the analysis, however, we have assumed that the qij(mn), i,j = 0 to N - 1, exist and can be found. Rather than proceeding to the existence proof, let us examin the results for the special cases N = 0,1,2. For N = 0, (7.50) becomes B(mn; v', ) =. For N = 1, the following relations are true: B (cc j Y, = t oo (cc) to B B(d;,) - $oo (dc) c-B B(c-dj vo)= (c d) =e f(7.69) where B is given in Equation (3.78). Note that it is only coincidental that qoo(dc) = qoo(cd) for linearly anisotropic scattering since the theorem in no way guarantees that B(dc; v', t) = B(cd; vo), For N = 2, the analysis becomes quite involved and explicit answers have not been obtained. Thus, even if the constants qij(mn) do exist for general N, the computation required in order to find them is indeed formidable. It should also be remembered that only if (7.64) is valid did we assume the proper weight function such that the eigenfunctions (4), < v < 1, and 0+(p) are necessarily both orthogonal to

-133-,(~L) + B(cc; v',. It thus appears pointless to attempt to prove the existence of the qij(mn) as presently defined if we want both the continuum and discrete eigenfunctions orthogonal to the same function. With the ideas of the previous paragraph in mind, however, it is easy to use the present analysis to formulate the proper funcCV' tion v',(() + B(cc; v', I ) 2 such that J0 (Iv) A t()-) A (V- V') (7-70) I,)+[ O (7.71)'()C.') + A(, M,) ] 0A. (7.72) One simply proceeds from Equation (7.39) writing everything in terms of the unknown weight function yg(O). Since yg(k) must reduce to (Vo - l)Y(~) when N = 0,1, Equation (7.39) can be used to write the general form of 7g( ) as ) =JTh)[w= o(u 7) + ] * D -', (7.73)

-134The existence of all the constants qij(cc) and qij(dc) for ij = O to N - 1 must be verified and the values must be obtained. Then the values of wo and w2 are determined by requiring that B (ccj i,) Bc(.cj >) (7.74) The existence of B(cd; vo, l) necessary to make (7.72) valid is guaranteed by the theorem. To summarize, we have attempted to determine a bi-orthogonal form which makes an eigenfunction of the homogeneous Boltzmann operator orthogonal to an adjoint function with respect to an appropriate weight function. We have observed that by proceeding in a manner similar to that of Chapter III, the determination of the bi-orthogonal form (if it exists) is quite involved. Perhaps this problem could better be attacked from a more abstract point of view using the properties of the adjoint Boltzmann operator.

CHAPTER VIII SUMMARY The method developed by Case was used to consider several time-independent problems of neutron transport theory for a homogeneous medium in plane geometry. This method utilizes a direct expansion of the neutron angular density in terms of the eigenfunctions of the homogeneous transport equation. The solutions to these problems were greatly facilitated by the use of half-range orthogonality relations (for O0 < p.< 1 ) between the eigenfunctions of the transport equation. Both the half-range orthogonality equations and the previously known full-range equations (for -1 < p < 1 ) were shown to be special cases of a more general partial-range orthogonality equation. These orthogonality relations were derived and were used to obtain solutions for the following problems: 1. The half-space Milne, albedo, constant isotropic source, and Green's function problems for linearly anisotropic scattering (i.e., for a half-space medium which scatters neutrons in such a manner that the scattering function is linear in the cosine of the scattering angle). 2. The half-space Milne and albedo problems for quadratically anisotropic scattering when there is no absorption. 3. The albedo problem for isotropic scattering in a slab of finite thickness. In each case, simplified expressions were derived for the emerging angular density and the density and net current on the surface(s). -135

-136A speculative approach for finding the half-range orthogonality relations for general-order anisotropic scattering was also discussed. Some of the results from the slab albedo problem were used to determine new integral expressions for the X- and Y-functions of Chandrasekhar. These yielded convenient asymptotic values valid for large slab thicknesses. Also included was a generalization of the slab albedo boundary conditions to the case of a slab bounded on one side by a vacuum and on the other side by a specularly and diffusely reflecting face.

APPENDICES -137

APPENDIX A DERIVATION OF EQUATIONS (2.30) AND (23.51) Extensive use will be made of the Cauchy integral formula,(27) which states that if f(z) is analytic on an open region S in the complex plane and if C is a closed curve such that C and its interior region lie within S, then for every point z inside of C, where g denotes the positive orientation of C Derivation of Equation (2.30) The derivation of (2.30) for the function Xg(z) is similar to the derivation for the half-space X(z) (7) The function Xg(z) - Xg(m) is analytic everywhere in the entire complex plane cut from a to 5 on the real axis and vanishes at infinity. One can therefore use (A.1) to write X (E) X (A.2) where C = C1 + C2 is defined as the contour shown in Figure A.1: -138

-139C2 Figure A.1. The Path C for Equation (A.2). If we deform the contour C2 to a very large circle and observe that the integrand vanishes faster than l/z, we find that there is no contribution from the integral over C2. (This reasoning will be tacitly used throughout the remaining appendices.) If we also allow C1 to shrink and enclose the cut, we then have x -+ l,"~~ I~ g Xt) X >(A.3) since there is no special contribution from the contour in the neighborhood of the endpoints (in view of the fact that Xg(G) possesses at

-140worst a weak singularity as p. -c,~ ). From (2.28), and this can be rewritten, using (2.15) and (2.19), as Y _ X+) - X- ) ~ < < A (A24) 2rrL Ao) Equations (A,3) and (A.4) combine to yield X () _-Xt() =f_,., (2.30) Derivation of Equation (2o31) For z v, Cx < v < 3, the use of (2.13) reduces (2.30) to — x o) - x +- X Pj; ar t (A,5) Multiplication by cv gives +,. %O -I C-).. 7 d = _P_,_ _ (A.6) — =Z-:~. —.. 1 (.,,,)+ -~~~ " x.',:, —"'-/,

-141Use of (2.15) changes the first term on the left hand side of (A.6) to C x 7 () +X /,)] = A g ( v) X ((A, and use of (2.19) gives C X +(t) +X -(v] xy )X- (>)X$(_ (A.8) This may be rewritten as... T X (V,)+ 4 [ + i )] [ L ][ L Tr Equations (2.15) and (A.4) reduce (A.9) to 4[C a' ) X(v)] = (v)) (A.10) and (A.6) and (A.10) yield the result X()d]) (2 OD)=Z p t,,) 0 e,(,<M (2.31) z..~~~~~~~~~'oe X-~ =

APPENDIX B DETERMINATION OF THE GENERAL-RANGE, HALF-RANGE, AND FULL-RANGE WEIGHT FUNCTIONS Determination of the General-Range Weight Function In the constructive proof of the general-range completeness theorem,(7) the function 4(i) was expanded as S^/f ) | A (i) 0Vf ) 1 7J (B.1) where the value of the expansion coefficient, A(v), was obtained by using the set of equations(7) -z-A() -I N- N,) N"'?rrG X()f T - (B.2) Using the Plemelj formulas (2.13) and (2o14) and (Bo2) gives T ()[ XL(v) + (v ] () Use of (2,19) shows that the function (X(z)-1 obeys the following relations along the cut a < v < f in the complex planeo -142

-143cI -T [X-) = A -(2 ) A j I r _ Z X/~-(~c) A(7v) 2( B (B. 4) ~ xiv) X7 (v) = T (' These equations may be verified by using (2.19) to obtain X*{~)AAt) t A~/ )7 [ Nf(v) + ] + = A"H) __ _ _ _ _ (B.5) L X tlt x -(Vt) x Att -tV)AV) and then using (2.15) and (2.28). Substituting the results of (B.4) into (B.3) yields the result A( v) () A ( v)Af-) ) _ 7' j. r' L - z v ar, S (B.6) and this can be rewritten as A(v1) (B.9A4')A h') ~-n/) 9 ) PY) S (B.7) Use of (B.1) in (B.7) gives

-144Atv)% ( )A*()A-(') = itr B~)?~ )' t')?~,).,d d~.(B.8 ) Equation (B08) can also be written in the abbreviated form f o t=)9%e s t) a4l = ( 2)A r)A,'c(-Y(y,9V- 7')(B.9) to parallel the form of (351). This is exactly the equation we are looking for since it shows that Yg(i) is the general-range weight function which makes av(>) and Ov,('L) orthogonal over the range a < Cc < e Determination of the Half-Range Weight Function When we are considering a = 0 and B = 1, we know that 53t)E{) HZ ~P~p _ O ~ i2g ), (2.35) This may be rewritten as a n 0 / (B.10) and (B.10) may be written as

-145II J i)+ A) g;:i, F) ( HoTV) as = ~ (B. 11) The eigenfunction X(4) is already orthogonal to Xv(>) by (B.9). To make 0+(0) also orthogonal to Gv(W), O < v < 1, the weight function for the half-range case must be Yg(~G) = (Vo - )()Y().* (This weight function is the desired one for problems were the region under consideration is in the right half-space because the discrete coefficient a+ is unknown. The weight function should be (vo + >)7(p) if the region is located in the left half-space.) Determination of the Full-Range Weight Function For a = -1 and f = 1, we have the equation (V)A(7/) =Pi d- (2.43) which may be rewritten as A 0000 ~~~~~A ~~(B.12) _ I or as f'[ [ (E =,// (B.13) *Of course 7g(l-) = K(vo-l)Y() where K is any constant is also a satisfactory weight function.

-146In order to make both the 0+(l) and 0_(l) orthogonal to,v(p)' -1 < v < 1. we must use (2.44). Multiplying (2.43) by v and using partial fractions and (2.44) gives t-1Sk~a~ 0. ) ~(B.14) Therefore we also have the equation Muti; - ] pli n g4. +w) T =O. (B.15) Multiplying (B.13) by vo and adding/subtracting (B.15) gives il _+ i 7)(i_ Z)~) w_ 0 = (B.16) 2 2 Therefore the full-range weight function is (vo 2-(1) This weight function may be reduced to that obtained by Caseo(7) Multiplication of (B.16) by 2(1-c)/c and use of (2.9) and (2.47) gives the result =/ t it),) S O ~~(B.17) so the function t is also a weight function for the full-range case — in agreement with Equations (3,1) to (353),

APPENDIX C EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHAPTER IV Derivation of Equation (4.17,) The integral which must be evaluated in order to simplify the emergent angular density is rI — iT ft_ )[piJ +B ]d (C.1),/A'70 Use of the identity d-dv )= d(2zt') d(8) - 7A/ (-) )( )(C.2) and (2.52) gives spr (+Gi) cI(vu')d9 (Z N - b,(l-)] dv (C )A Jo (v- Vo) 7(v) A t(v) A-(m) -147

-148For linearly anisotropic scattering, the function (X(z))1 has the following behavior on the cut 0 < v < 1 in the complex plane: I 1'-1 () 2"t L' X(2) X -X T(V) A+(7) A{V) f + - -A _V ( ) 0 J < <l (C.4) 2 L X+(X)X-(>) J ((c) +(V ) (i) The derivation of (C.4) is analogous to that used to obtain (B.4). Using the results of (C.4) as a guide, consider the function [(z-KI')(Vo-z)X(z)]-1 which is analytic in the complex plane except for a simple pole at z = vo and a simple pole at z = p' imbedded in the cut from 0 to 1. The function vanishes as z - X oo Letting C3 go to infinity in the contour shown, C3 C7 \I, Figure C The Path C for uation (C Figure C.l. The Path C for Equation (Co5).

-149the use of the Cauchy integral formula gives (c 5) (- t')(i/o-)X(~) ZrL (- (,0)(vo-')X( )('-.5) Equations (2.12) and (C.5) can be used to obtain the result z, "[' + 1' ] + I.... ______ _'_ _ r)(_[X(C.)6) + L -. (c.6) (V y >,) (vo - ) Xlao) Here, the sign of the second term on the right hand side is negative because the contour is clockwise around the point z ='. The third term comes from the residue at z = vo o Use of (Co4) in (Co6) shows that, for z = -,.p I i ifi _+ A. (I.7) a.vb.) (VO, t. s.., I V-oATV..IV,.V..

150 Use of the identity r) _ d() b,(-.)') (C.8) and Equation (C.7) reduces (Co,3) to I = _,/,t,') +. (.. ~O, (.,,,, ) (w) A 7A(. A,-/, v,) v+W ().) AC..) N, o (on - t) f(i) A+v-)( ) Now consider the function [(vo-z)X(z)]- i It is analytic in the complex plane except for a simple pole at z vo and the cut from 0 to 1 and varies like +l as z -, oo The contour for the Cauchy integraL formula is shown in Figure Co2. I C3 Fi-r Co Tt I E Figure C02. The Path C Used to Obtain Equation (ColO)o

-151Again, the curve C3 is allowed to become large. The result obtained when z = -kL is I -+ (SO+f) X(-') (+D ) X(4) 0 (7/0-.) (v~j)T iv)A +(v)A-) The last function to be considered is b, (I-c) (Vo-z + dVO 7) d(Zz) X(a) It is analytic in the complex plane except for the simple poles at + 4-/7[bl(7-c)]- + C * and the cut from 0 to 1 and it behaves like -(vo - v) as z - oo. The Cauchy formula then yields the following result if z = v: b, _ I - c) ()Jo -z)/o + d I(g0) d(J~O,) X-Vo) i (2)7[ b(1 ( v ) + d(Io2)] ( ('-J- ) T(v) A+(' ) A-()') _+ F.....[b.(.-c)(7O-)- - d.i)] _ [b -4(7{-.-J )] (C 11) 2 b,(I-) (I- o x(6) 1(-1) )X(- -) J.l --,, *No pole is imbedded in the cut since bl(l-c) is always > -1. If bl were less than -1, the scattering function would be negative for some angles and this is not physically possible.

-152The first term on the left hand side of (C.11) reduces to -I[X(vo )-1 To simplify the right hand side, one needs to know the values of X(+ a). Following the approach of Shure and Natelson,() one may write (2o67) as Y(~t)- b8l- )f )d (i-c)+W (i'~<) (C.12) Using (2.37) and (C.12), it is evident that X(-')X( —, (I- ['] A W Cd) (C.-l) and using the definition of + a gives the identity [T(0)] 2do - Ai(to) /2() 2 (C.14) From (2.56) it may be shown that A(a) 1-c so we have the identity CbO [2()] d(vi7) d('o) _ (___ _ (C, 15) where (2.59) has been used. After considerable algebraic manipulation, one finds that the contribution in (C.11) from the residues at + a and - a cancel. Thus, Equation (C.11) reduces to I (~13'[b,(l-,10- AdiSI (inv ),(c.l6) O ( -~ to) OT (V) A\+(vY) A Xl O)

-153Putting (C.10) and (C.16) into (C.9) eventually leads to the result I1= (-vo Y/4)[) - BA' + (-) y) (c 17) (vo~) (,o -r') d (a, o-) X(tr) Derivation of Equation (4.24) In order to derive the equations in (4.24), we need the following generalizations of Equations (Ao4) and (A.10) to the case of linearly anisotropic scattering: d() (e - X)~ [x~++ X-(v] - 7.(v) m~)v 0 v4 1 (C.18) We consider the function

-154where E is a small negative number. The function vanishes as z c- oo and is analytic in the entire complex plane except for simple poles at z = e and z = + a. It also has a cut from 0 to 1 because of the X-function. Using the Cauchy integral formula for the region of Figure C.3 shows that Figure C.3. The Path C Used to Obtain Equation (C.19). d( Z - e) (E 2V- 2T) ( - ) (voE'- E 2) XC() ( 1.) d(go) (/ -z)

Taking the limit as e -e0 and using (C.12) and (C.18) reduces (C.19) to Ioz Xo 0) _ t~)d,o) d(-. ) (C.20) Z d ( ) Letting z = -A' and z = + vo gives the first two equations of (4.24). To evaluate the third integral in (4.24), consider the function (>/a - z~) X/Z) and follow through a similar analysis to obtain rd Ir -(l-C)(-C3b) i: Joc I X((- 2 2)X(c) s (ttin dz(-) += 0 ]i (C.21) Setting z = 0 in (C.21) gives the last equation of (4.24).

APPENDIX D MODIFIED ORTHOGONALITY RELATIONS FOR ISOTROPIC SCATTE.RING In the derivation of relations (33016) to (3,25), the weight function (Vo.)y(Lt) was used in order to make ~j(~) and v&(B) 9 O < v <. i orthogonal. If the weight function (V o4-)Y(P) were used9 then 0(-) and ~V(') wouLd have been orthogonal over the half-range of u With this in mind) the following generalities of (3o16) through (3o21) may be obtained for 0 < v < 1 and 0 < v <K 1 -.P7,7~A,~:-,.I1-'tJ,): -,-,~, (D o72).tt' -jv),,2(D.4) 20 Jo-4 Equations (D.2), (D,4), and (D,5) simplify to J41~~~~~ — +~~~~ ~(Do7) -156

-157= 2 X_) (D.9) Adding (D.2) and (D.3) gives C -1) + (D.10) Adding both equations of (Do 1) results in the answer 7V)?it= T h'JA) A+(i)A /4) A (2-i') (D.ll) and adding both equations of (D.6) yields t = X(-7) -v/ ) (D.12) From (D.10) it is seen that neither 0+(k) nor 0_(j) is orthogonal to Xv(4), 0 < v < 1, with respect to the weight function y( ). This is of no concern, however, in problems which do not have solutions in closed form such as the slab problem. In fact, Equations (5. 11) and (5016) are most easily obtained using the above relations. Multiplying (5,7) by 7(kp)dkt and integrating over the half-range of it gives (5.11) after Equations (D.7) through (D.9) are used. Use of (Do10) through (D.12) gives (5.16) directly —without the partial fraction decomposition which is needed if (3.16), (3.17), (3.18), and (3.21) are used on Equation (5.7).

APPENDIX E EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHA.PTER V Derivation of Equation (5331). The reflected and. transmitted angular densities at the surfaces may be written, for > > O0 as W(~, -r)= qt 2nt) + Ct M?) +J MV)) es 1 0 + Be,a___ (Eol) d../o d/ -/ ~(dg) = f,g4)e +a f)e A1 v) )e 4d' ji [+.l'Vv) 13_,f (E~2) where (5o3) and (538) have been used. Use of (5o 16) in the above equations gives -+ Q- [ -X(~~)- ]-,..... +J A-v2[eJAi_ Xvy v() + VA(v') A-, J) -158

-1599'1 &MO) d/ () MA +1v] A, -() l]/ (E.3) and 9'(d)l4) = ~h) - 9?~vt4) c/o I c/v3 [2P) ofV T()\(v, / ) ] /A )_ ) XZ )(( — T')] (E. 4 ) Consider first the function [(z-vo)X(z)]- which is analytic everywhere in the complex plane except for a simple pole at vo and the cut from O to 1 o It behaves like -1 as z -, o. Considering the contour in Figure E.1 and using the Cauchy integral formula yields I a C:~ + |=,. i d (E.5) c,+ ( —.

-a60C3 / ~ ~I oz /I C2,I Figure Eo 1 The Path C for Equation (Eo5), after letting C3 go to infinity. Use of (Bo4) with z - 9 0 < a < 1, gives _+.. Z. (E,6) (WO x) X3(7) (VotXit) (NO06v)tj/) A similar equation may be obtained by letting vo - -v0 Multiplying the equations by + ~ X(+ vo) gives the result 2 -, H ) A+ (V)A(v)). CA x( o) (+ - X(.{.) (E.7)

-161The second function to be considered is [(v+z)X(z)]-. It is analytic in the complex plane except for a simple pole at -v and the cut from 0 to 1 and it behaves like -1 as z -.o Proceeding in the same manner as above, one obtains I I, - + - = (E.8) XY X(-Y) X(-v) x(lJ )o(oo_ v) Multiplying by 2 X(-v) and then adding the term.() to each side allows us to write (E.8) as 9>4)?-X(,') Z -., A v+k' A-y'J.-s -m-. "-), M (E.9) The last function to be considered is [(z-io)X(z)]- It is analytic everywhere except for a pole at po which is imbedded in the cut from 0 to 1. It also behaves like -1 as z -co. Figure Eo2. The Path C Used to Obtain Equation (E.10),

-162The Cauchy integral formula and Equation (2.12) give.+ I _ p, (A+? )X()+ 2iio (Vf'o)_)_ /2' [ { + ~i / l (E.lo) after letting the outer contour go to infinity and setting z = - Use of (B.4) and multiplication by -Y7(o) reduces (E.10) to the form?A...^>..... =..... _ o) (E.11) Jo i} A A+(AA ^-~) (o) X -) Equations (EI7), (E.9), and (E.ll) reduce (E,o3) to ~ [)= -( I Mo) ( /, 94)X i') +cta(,)_-x ) +nA - EZ) Xa-vd ( ) i n X- El) gv e and (E,7) and (Eg) in (E.4) gives( 2)

-163j-J0vo d/vo - 4A(~,) e t7,>)X(-72)1 q X/-vol,o V X'(' - A/'v)e (E. 13) - &. + Now examine Equations (5.4) and (5.6). Multiplication of both equations by 7(k)dt and integration over ~1 from 0 to 1 gives A (-2,) X(-,)dJ - ____bc~/oD ca-'o (') - e 2 X(-vo) + e X(Vo) A(v) 2 Y(-V)dw O (E.14) where (D.7), (D.8), and (D.9) have been used. The reflected and transmitted angular densities, for t > 0,are therefore given by rX ) -, 1 tfl 4 )

-164 - X(7))% +, (5G33) The same equations can be derived in another manner(33) by use of results obtained in the constructive derivation of (5.11) and (5. 16)o Derivation of Equation (5.536) Using (5o1), (5.2), and (5o33)9 the surface densities and net currents are?o- >Xq-,,),o) +0jX(r J0 -.+) U0T X(A,) J +t4(vhe xt-w~tPdS (E. 5) +, &/Q kek~es~~p"~c- o~ A ^-,(-(,) xf~) j

-165and t qn)) MOXF C,f J0 )XLi0+ ) ] + _ -2 -,O X [-~ " J AM e X(-tS t'm/o I) (OE. 16 i' x(/-) )o) dJ (d) - ++. z I __ _I (El) 2j>(A ( /7) Xt L -t ~['

-166In order to evaluate the above integrals, consider the function (vo _z2)X(z)(zE)l where E is a small negative number, The function behaves like +1 as z - co and is analytic in the cut plane except for a simple pole at z = e. Applying the Cauchy integral formula to the analytic region enclosed by the contour shown in Figure E.3 gives N(yi- ) _J' _) L —ZX T[x ) A-x- ] _ 70 X() (E. 18) [ 2 Figure E,3, The Path C Used to Obtain Equation (E8) >(I- [ X) N =,_A I (E.1 ) TV~ ~ I C)

-167Using this and taking the limit as e - 0 in (E.18) yields -"_ zi-)(,:...) X(.) _I X'O) (E.20) X~ d- = 2(1-C) (JoZ- z ) XIZ) I_ As z = -c,i + Vo, and +t, one obtains - (- c -I + ] C dr =- 4MF>)d, = 2(1-c) C _ coXlo) (E.21) J- X(-,r) -'c Z M z Equations (E.15), (E.16), (E.17), and (E.21) combine to give _ 2(lC)2 )n~ + C X() -0{ICA X1v)-f )+(-V)X(-V)+.d].......[+C X(.JO +._eS A(),)e X(-+/JJ+/ /(o)- xto)L o',I l Z, Xv)-( D)

-168- i__ c)[M2 I_ MD) + a X4V) +4 X +A c XID) f) X) Zxro -O z (-D) (- -v)-Z xr) (d) - V(cqO X(-V) + a e X()+J A()e d X-) [' X( )(X)'"cJ/Y~ c Ul~ dl"~~ c... V XrC] (x) - q+e ) X(-+ pJ (E.22) and these results reduce to (5.36) when the equations in (E.14) are used.

APPENDIX F GENERALIZATION OF SLAB BOUNDARY CONDITIONS Using the solution of the slab albedo problem obtained in Chapter V, it is possible to write the solutions to related slab problems. For example, suppose we wish the solution to the albedo problem for the slab which specularly reflects neutrons (or photons) at the surface located at x = d. The boundary conditions on the angular density, 7r(x,k), are: (Irr jdj t) = Y/rd), ~-I a I (F.1) Using an approach which is essentially the same as the method of images in electrostatic theory, the result may be written down immediately. The angular density is given by (tr(X/)- Yk(X<)+ it X/) 0O X. d (F.2) after the slab thickness in the expansion coefficients is doubled to account for the image slab. Here, *a(XA,) is the solution of Chapter V and aa' (X,) =?a(d-x,- ). As another example, suppose that one desires to solve the slab albedo problem under the condition that the face at x = d reflects neutrons (photons) in a diffuse manner. The boundary conditions on the angular density, *d(x,)), are: -169

O lJ,) - &rA, o) A)>O Vd (ed, 7) - sfo)0'(d z); o,>,(F,3) The function f(o0) is a function which prescribes the angular distribution of the diffusely reflected neutrons. The constant S is to be determined by the requirement that the number of neutrons leaving the slab at x = d be equal to the number entering, i.e., 50 -//, c't (~)+,,die or J5 = oif Ad (9<4A~~)d.(F.4) In (F.4), we have chosen the normalization j uof(o)do = 1. The diffuse angular density is given by (x (X) + S% (F,5) and (F,4) and (F.5) combine to determine the value of S as I fJo c Atj (F.6) {.J A. J AA r i i i 1 ~~gf'~fo

-171Equation (F.6) may also be written as 5 A- T (F.7) I rjto R fIo) *4 fo where T and R are given by (5.43), (5.47), (5.60), and (5.61). As a guide to the study of photon transport in the atmosphere, the solutions to the slabs with specularly and diffusely reflecting faces may be used. Assume that one desires to know the one-speed steady-state photon distribution (ignoring polarization effects) in a homogeneous atmosphere illuminated by an incident distribution from the sun. Scattering in the atmosphere is assumed to be isotropic. Photons striking the earth's surface, however, are permitted to be either specularly or diffusely reflected back into the atmosphere or absorbed. The angular density in the atmosphere, T(x,i), is then given by PE(xt) I'F(X<)i +N tdI (x)) (F.8) where M and N are the coefficients of specular and diffuse reflection of the earth's surface, respectively. Clearly, M + N < 1. If the diffuse reflection obeys Lambert's law, i.e., if the diffuse reflection is isotropic, then f(10o) = 2 and the problem is solved once M and N are known. A problem arises when one attempts to apply the above theory to a model of the earth's atmosphere, however. When the scattering in

-172the atmosphere is most nearly isotropic, the atmosphere is clear and the mean free path thickness (dependent upon the photon energy) is of the order of unity, For a slab thickness this small, a high order of approximation is required for the solution of the Fredholm integral equations containing the continuum expansion coefficients. This becomes quite unwieldly since even the first-order approximation involves the numerical evaluation of the complicated integrals of (5.30)~ In attempting to apply the zeroth- or first-order approximation to the model of the earth's atmosphere for thicker slabs, one must be concerned with the transfer of monoenergetic radiation through a cloudy or foggy homogeneous atmosphere. This has the advantage that inside the atmosphere, polarization effects may be ignored which make the model more realistic. The primary failure of the above theory in this case, however, is that a foggy medium scatters radiation in a highly anisotropic manner. Another important drawback to the present formulation of the problem of photon scattering in an earth's atmosphere is that azimuthal symmetry was assumed. If the angle of the sun's monodirectional beam of incident radiation is such that uo i 1, the above model is satisfactory, For other values of 4o, however, the transport equation without azimuthal symmetry would have to be considered since the boundary condition for the sun's radiation is of the form O ( J+ = s A,4A- ) A 0 (F.9)

1735 Equation (1o3) would. then appear as aO= = a o iZE IX A)!Q pi'); [q(6 ~ ]]) bY ~ d, (F.10) which can be written more concisely as [I X 4wl]1 -t>>rJ (. < A Tkr)! 9t ]) d~at [ l/(P e~) 0] gr< _) (F.11) With the method of invariant imbedding, Chandrasekhar(2) was able to solve (F.11) with the boundary condition (F,9). After separating out the uncollided beam from the angular density, he transformed (F.11) into an inhomogeneous equation for the angular density of the radiation field which had suffered at least one scattering interaction. He then expanded the angular density of the transport equation in a finite Fourier cosine series in terms of the angle ($ -,)

-174Because he had separated out the uncollided beam of particles, the solution eventually reduced to the solution of N azimuthally-independent transport equations. An identical approach could be followed using Case's method. The solution to each azimuthally-independent, non-homogeneous transport equation is found by solving the Green's function problem.

APPENDIX G EVALUATION OF INTEGRALS NEEDED TO FIND THE SURFACE DISTRIBUTIONS OF CHAPTER VI Derivation of Equation (6.59) The integral which must be evaluated is ~3,Vbn {_,) f")) 1d. + (G.1),/~0 Use of (6.9) gives Ib L 2 (' ) f' ),' f -) 1' 5) 7 C (;/C 4TA))') Xt) d)/,' ~3 b~:o' /J dz( v 4 z$-......... I ((G.2) Jo TA(v) 4+ r)A-(. (IV /) The behavior of (X(z))-i on the cut O < v < 1 is given by -175

-176Fr |_ _' 1 ( 7 )Z Tre L X-(V) X-(") J _ 2 ()V (G-) i LX' X (v)J ^-/) ) ( ) Now consider the function [(z-I')X(z)] 1 It is analytic in the complex plane except for a pole at p' imbedded in the cut from 0 to 1. The function behaves like -1 as z -- co The contour to be used when applying the Cauchy integral formula is that of Figure G.l. Upon letting the contour C1 become very large, - +1- | = ds (G.4) Cl Figure G.l. The Path C Used to Obtain Equation (G.4).

-177Using Equations (2.12), (G.3), and (G.4) gives, for z = -y, I r' (X)2'f(v)dJ ~,) X(:) JO, (v ~ (V (O A+(v) A^ 2 (G.5) T'34A) A+4')t J A-k)m We now consider the function (z+v)/[f(z)X(z)]. It behaves like 4/(3b2) as z- oo and is analytic in the complex plane except b2 3b2 for simple poles at z = + 4(1 + )/ =+ a and the cut from O to 1. Applying the Cauchy integral formula gives the equation _(_)X(_ ) 34- _ O 1Xi) X(v) 4- 4 (+Fo ) 4 (- + ) (G.6) 3ba (2 6) (4t ) X(-T) After setting z = -F and using (6.47), (6.50), and (G.3), one finds that ( (J) f(I)(V )'i ) );

-178Equations (G.5) and (G.7) reduce (G.2) to ~t 3b'' 3]- } (6.59) since f(-Ol) = f(). Derivation of Equation (6.63) The behavior of the X-function on the cut is given by xwu + (A - X - V ( - 6 )X04V) Equation (G.8) and the Cauchy integral formula can be used to consider the function X(z)/f(z). This function is analytic everywhere except at the poles z = + a and the cut from 0 to 1. It vanishes for large z. Using the Cauchy formula, one therefore obtains t F X~s) + ~r-)' (G.9) Letting z = 0 and using (6.47) reduces (Go9) to dn Z[H43 5m bt' > 3] 3t? (6.63) J0X(M);" ~[l 3, -r ~

-179To evaluate the second integral of (6.63), we consider the function _) -ioI- 9I q(,)_ =L+ 4T +..... and use the Cauchy formula to find that += 3 J [I + ~4 + 4. / - 4. z (I- b,./5) X(-,)- 6) — -) _4C 4-l~l 6Z ( 3 4- ) (3-b)( zLetting e go to zero and setting z = -o0 gives z (I- b/s-) _ X(-,)_ 4- X(o) 3 +(I +bZ/4)To) ___+____9L- _ 4____ o]_ _, (6.63) (lIb./4)/o after Equations (6.47) and (6.50) are used.

APPENDIX H RELATIONSHIP BETWEEN CHANDRASEKHAR AND CASE QUANTITIES AND NEW INTEGRAL EXPRESSIONS FOR CHANDRASEKHAR S X- AND Y-FUNCTIONS Some of the emergent angular distributions obtained in Chapters IV and VI were also obtained by Chandrasekhar(2) using the method of invariant imbedding. This appendix contains the relations needed to check the agreement with these results, These half-space relations correlate the Chandrasekhar H-function with the Case Xfunction and also the moments of the H-function with the moments of the Case y-function, Also included in this appendix is a derivation of new integral expressions for ChandrasekharVs X- and Y-functions for slab problems using the results of the slab albedo problem of Chapter V. Correlation Between Chandrasekhar and Case Formulations One of the primary notational differences between results obtained by using Case's method and those of Chandrasekhar arises because of the differences in the coordinate systems used' -l -l cos co 1 0 o T1 O d Chandrasekhar Case Figure H.1. Coordinate Systems Used in the Chandrasekhar and Case Formulations. -l80O

-181Another important difference is the fact that Chandrasekhar's work evaluates only the emergent angular distributions I(O,0) and I(d,-C), of particles which have been scattered at least once whereas the approach of Case evaluates the entire angular distribution, r(x),), regardless of the past history of any particle. This means that Chandrasekhar's law of diffuse reflection and transmission and the slab albedo problem of our Chapter V, to which it corresponds, are related by the equations 2 -Y /I(, = _M 0 (H.1) 2IF~t,,? = tat ~ eL K(,ot) (H.2) provided the scattering is isotropic,* If the scattering is anisotropic, (H.1) and (H.2) relate the Case angular density to the azimuthallysymmetric part of the Chandrasekhar intensity, This distinction arises since the boundary conditionsof Chapters IV, V, and, VI assumed azimuthal symmetry whereas Chandrasekhar solved the problem wit'nout making this assumption, In the limit as the slab thickness becomes infinite, (Hol) can be used to correlate the emergent angular density for the half-space albedo problem with the azimuthally symmetric part of Chandrasekhar's solution to the problem of diffuse reflection. *The quantity 2/F is a normalization factor.

-182When attempting to correlate the results from Chandrasekhar's constant net flux problem, which is the half-space Milne problem for c = 1,* Equation (H.1) again is valid. In fact, (Hol) is valid for a m-dium which scatters particles in an anisotropic manner since the constant net flux problem has azimuthal symmetry. The solutions to Chandrasekhar's half-space problems are found in terms of the H-function, which is related to the Case Xfunction by(ll) I - (r~ x~r 2(r) (Hex) in the case of isotropic scattering. By the use of (2037) with z = O, this may be written as Ott4) = (14Mb') X>M) D(H.xo) This latterform is valid for arbitrary scattering law and conveniently illustrates the behavior as c -. 1 (vo -.) ) I _ X( -~.,C. (H.5) HkZ) Xl o) - *Note that the Milne problem answers of Chapter IV can not be used if c = 1 because the normalization was selected such that and this normalization is not correct for the discrete eigenfunction carrying the current, V2(x,L). In Chapter VI the normalization for the Milne problem is given by Equation (6.27a$. The normalization must be selected for the solution by Caseks method in order for it to agree exactly with the results of Chandrasekhar.

-183From (2.37) and (3.34), we note that X(O) = 31/2 for isotropic scattering and c = 1. For the linearly anisotropic scattering problems of Chapter IV, we have the result IX _ I __ (H.6) where we have used (2.37) and (2.59). For the quadratically anisotropic scattering problems of Chapter VI where c = 1, Equation (6.46) shows that 1/D X(O) [ b/53 (H.7) Thus, the relation between the H-function and the X-functions of Chapters IV and VI are 0(-00 ( CL> X(- M) IV) M)"' [(8 VI), (H.8) The solutions of Chandrasekhar are in terms of the H-functions and the moments of the H-function, defined by ^. = 8 yZ X~b) Jr L = OJ TV v T' @ (H.9)

-184Thus in order to complete the correspondence between the notations, we must relate ai to the y(i) of Case's formalism. For the y-function of Chapter IV, Equations (2.38), (2.59), and (H.8) give t 2+ I_ t( A be)Ax i" (H.lo) so 0 1:- =2(- c)(l- [- ] - (H. 11) The function 7 is obtained by setting z = 0 in (2o67)o In particular, the specific relation needed to correlate the emergent angular distributions is _Wo_ _ __- (H. 12) In Chapter VI, Equations (6o52), (Ho8), and (Ho9) can be combined to give the relation i,'d=,[_,_ _ )1__ -= (H 13) The only other relation needed is given by 3 6a'~z /6, - + /6 (H. 14) in (4 +ba),Z( b)

-185and this is proved using (6.50). It should again be noted that the source strength in the solution to the Milne problem is normalized by Chandrasekhar to a net current of -1 instead of -1/(3-b1) as was done in (6.27a). New Integral Expressions for Chandrasekhar's X- and Y-Functions The X- and Y-functions are defined by the coupled equations /~)- I {4 [/)X(,')- -n>...... A Y( =) / tj [>XlA) -By,)7t$)]+ (H.15) for isotropic scattering. From Section 62.1 of Chandrasekhar's book, (2) it is seen that his functions are related to the reflected and transmitted angular intensities at the surfaces of a slab illuminated by incident radiation from all directions io > 0 with magnitude proportional to 1/o. Using (H.1) and (H.2), the relations for X(A) and (YG) are y-ok~~~")~~ eo- e&ags-)](H,16)

-186where V(0,-I) and I(d,p) are given by (5.33). Adding and subtracting the above equations and using (5.8) and (5.33) gives Xk)~ 4)=~I + X2LfI I~ [' "o) r -J/< ~+ P.)x~voX +_ e L+ )Xr-) + B)e (H.17) From (2.34), X(o), ok (H. 18) Using partial fractions, (H.18), and (2.43), one can show that 0 _ 1 T1X) dXo) I ji >~~ob,) dP f )'1 - xr -- (H.19) /~o zvox 1VT:2O- x[:~) Using the above equations and (5. 12), Equation (H.17) becomes X(o) _____ /Vbd/v - < X e~e X(->J Qi(2)JJ (H.20) XCMJ %r )e (-)(Vd'

-187where, using (5.13) and (5.16), fi cI/V' f b+ _d =X(O) JA/ o " +) )e XF')4 (H.21) and 7+ /) e (J)X/-s/'voI')d d (H 22) QI) =1tt d,, -e The following identity may be proved: 2iMe' -dt o, () + xd'-v (H.25) e _e Use of (H.21) and (H.23) to eliminate q+ in Equations (H.20) and (H.22) gives two equations of the form

-188Xp(;+)iO)~{ c( fo r ~ + to) _O310 x(-~,(H.24) and x(-'(C-) 2 I- c + ( o F-~w (,oz,)2v,+,t ~v,/ e X(Jj j t wLb'~, c ( Er,~+J )] L (H.25)

-187where, using (5.13) and (5.16), C Xl(o) iJZ / 4+(,)e X( (H.2) = (~H.21) ~o Xli-o)[e 2_o/i + Cd/io] and Q(,),)( t X(-)[Y)e -Rrt~e~d/X~]t te f4/v (S)Xt w8JQ (8,92s,] (H.22) The following identity may be proved: Use 20 + -- ( to U ( tO ) (H.23) Use of (H.21) and (H.23) to eliminate qa in Equations (H.20) and (H.22) gives two equations of the form

-188i4 (H o -(24) and (~ } /t a (H. 24) X(-%) {t t- (, iI() + C / "' Z V )] ] Z o-,, 9v 2QR+,:dv (H.25)'U J

-189Equations (H.24) and (H.25) correlate the X-function of Case with the X- and Y-functions of Chandrasekhar. For large slab thicknesses, integral terms with e-d/v vanish more rapidly than terms proportional to csch (2zo + d)/vo and can be ignored in the lowestorder approximation. Ignoring the integral terms corresponds to obtaining the effect from the discrete modes, i.e., the asymptotic variation. Equation (H.24) then gives za +J + z>+J) ty) +!k ) = ) z ) / ( (H.26) so, for large d, _ ~T Xt (74o SE) =Ho -? LA (a), —1 O -t Y9 "Z")" \ ( / (H.27) The equations of (H.27) can be obtained more directly by using the asymptotic forms of (H.20) and (H.21). It will be noticed that the forms of (H.27) are familiar — they are just the asymptotic emergent surface densities, p(O) and p(d), obtained from (5.37) with Ho replaced by pt. This corresponds to the alternative interpretation of the X- and Y-functions given by Chandrasekhar,

-190The results in (H.27), although only approximate, are useful since the numerical work on the Chandrasekhar functions(38-40) considers d to be small. The integral equations of (H,16) are best for smaller thicknesses and require more iterations as the thickness increases. In the limit of infinite thickness, the equations in (Ho27) yield XfA) PY9A Y m (H. 28) where, in the first equation, (H.3) was used. This behavior was first noted by Chandrasekhar.

REFERENCES 1o K, M, Case and P. Fo Zweifel, Neutron Transport Theory, (to be published)o 2, S. Chandrasekhar, Radiative Transfer, Dover, 1960. 30 V. A. Ambarzumian Compt. Rend. Acad. Sci. U.R.S.S., 38, (1943) 229. 4. R. Bellman and R. Kalaba, Proc. Nat'l. Acad. Sci., 42, (1956) 629. 5. K. M Case, F. de Hoffmann, and Go Placzek, Introduction to the Theory of Neutron Diffusion, Vol. I, U.S. Gov't. Printing Office, 1953. 6. B. Davison, Neutron Transport Theory, Oxford, 1958. 7. K. M. Case, Ann. Phys., 9 (1960) 1; also, Recent Developments in Neutron Transport Theory, Michigan MemorialPhoenix Project Report, The University of Michigan (1961). 8. M. R. Mendelson and G. C. Summerfield, J. Math. Phys., 5, (1964) 668. 9. G. JO Mitsis, Nucl. Sci. Eng,, (1963) 55 or ANL-6459. 10o G. Mitsis, Transport Solutions to the Monoenergetic Critical Problems, Ph.D. Thesis, The University of Michigan (1963) or ANL-6787. 11. M. R. Mendelson, Solutions of the One-Speed Neutron Transport Equation in Plane Geometry, Ph.D. Thesis, The University of Michigan (1964). 12. R. Zelazny, J. Math, Phys. 2, (1961) 538, 13. R. Zelazny and A. Kuszell, Physica, 27, (1961) 797. 14. A. Kuszell, Acta Phys. Polon., 20, (1961) 567. 15. R. Zelazny and A. Kuszell, Ann. Phys. 16, (1961) 81. 16, J. H. Ferziger and A. Leonard, Ann. Phys., 22, (1963) 192. -191

-19217. R. J. Bednarz and J. R. Mika, J. Math. Phys., 4, (1963) 1285. 18. W. T. Fehlberg, Thermal Neutron Distributions Near Material Discontinuities Ph.D. Thesis, California Institute of Technology (1964). 19. R. L. Bowden, Time-Dependent Solution of the Neutron Transport Equation for a Finite Slab, Ph.D. Thesis, Virginia Polytechnic Institute (1963) or TID-18884. 20. I. KuscVer and P. F. Zweifel (submitted for publication to J. Math. Phys ) 21. Jo Ro Mika, Nucl. Sci. Engo, 11, (1961) 415. 22. F. Shure and M. Natelson, Anno Physo, 26, (1964) 274. 23. N. G. Van Kampen, Physica, 1, (1955) 949. 24. K. M. Case, Ann. Phys., 7, (1959) 349. 25. Co Cercignani, Ann, Phys., 2 (1963) 219o 26. R. E. Aamodt and K. M. Case, Ann. Phys., 21, (1963) 284. 27. See, for example, Chapter 3 of WO Kaplan, Operational Methods for Linear Systems, Addison-Wesley, 1962. 28. N, I. Muskhelishvili, Singular Integral Equations, Noordhoff, Holland, 1953 5 vv 29. Io Kuscer, N, Jo McCormick, and G. CO Summerfield, Ann, Phys. (to be published). 30. 0. Halpern and Ro. K. Luneburg, Phys. Revo, y6, (1949) 1811. vv 31, I, Kuscer, Can. JO Phys., 31, (1953) 11870 32. F. G. Tricomi, Integral Equations, Interscience, 1957o 33. N. JO McCormick and M, Ro Mendelson, Nucl. Sci. Eng,, 20, (1964) 462, 34. I. Kuscer, J. Math.& Phys,34, (1956) 256. 35. A, Leonard and T. W. Mullikin, J. Math. Phys., 2, (1964) 399.

-19336. I. Kuscer, Bulletin Scientifique, Yougoslavie, 2, No. 1, (1954) 15. 37. E. D. Rainville, Special Functions, MacMillan, 1960. 38. S. Chandrasekhar, D. Elbert, and A. Franklin, Astro;hys. J., 115, (1952) 244. 39. D. F. Mayers, Monthly Notices of the Royal Astronomical Society,!23, (1962) 471. 40. Y. Sobouti, Astrophys. J. Supplement No. 72, Vol. VII, (1963) 411.

UNIVERSITY OF MICHIGAN 3 9015 03465 848711111 3 9015 03465 8487