THE UNIVERSITY OF MICHIGAN CNDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING TIME DEPENDENT NEUTRON ENERGY SPECTRA IN PULSED MEDIA Surendra N. Purohit A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in The University of Michigan 1960 February, 1960 IP-422

ACKNOWLEDGEMENTS The author wishes to express his thanks to Professor P, F, Zweifel for his continuous interest and help, He also thanks Professor H. J. Gomberg and Professor R, K, Osborn for their encouragement during the early stage of this study, Thanks are also due to Dr. A, M, Weinberg and Mr. E, P. Blizard, who made the facilities of the Oak Ridge National Laboratory available to the author during the last stage of this study, Drs, L, Dresner: E, Inonu~ and Wo Hafele are thanked for their stimulating discussions during the authors's stay at Oak Ridge National Laboratory, The author expresses his gratitude also to his wife, Henriette, for her cooperation and typing of the manuscript, ii

TABLE OF CONTENTS Page ACKNOWLEDGMENTS w * * -:: ii LIST OF TABLES. * -.'6Wi a a' a. vi LIST OF FIGURES. 4. 4 -.*: a. *. *' a i Vii CHAPTER I INTRODUCTION - * - * i *.* -. * a:.. V.. I I Experimental Studies4,- -,..,.. *:,..,.. -.,.. 2 II. Theoretical Studies...*.,....,.* * *. 4.. *., 4 *. **. *.,.4 III, Treatment of the Problem........._.: -.,...**'....,**..8 II, THE MAT:HEMATICAL FORMULATION 0OF THE PROBLEMl.*...*****,*......*.. 12.I* Boltzmann Equation...2..,..,.,. -.*...-:..,,: i... 12 II Angular Distribution. *-*. * -.. -,,.,,,,:.*.:* 14 III. Space Distribution.. *....,:.,...,*..*....,,.,......,. - 17 IV.~ Time Dependent Proble *.19 V. Energy Independent Solut1ion:,...,....*,-.....-..,. 21 VI. The Scattering Integral*..*...**.. *,***.. -..:* -** * 2 III, FAST AND THER1AL NEUTRON SCATT:ERING,......:... r....., 24 I. The Sloving Down a4 A. The Heavy Element Scattering Integral (A >> 1)**- -'i*w -*:*-:: a.26 1. Determination of S(EE-) by Fermi Age Model 29 2. Determination of S(Eits) using the Ga-ussian Distribution. x,,x,.. * -, % *.~ * b.-.-., 3.0 I3B Discussion of the Space and Time Dependent Slowing down Solution for Hydroge4n:..,*.:,*.,*.*... - 32 II, Thermal N-eutron. Scattering... 4.,.*,....,..,..,... 36 iii

TABLE OF CONTENTS (CONT'D) CHAPTER Page IV., EIGEN VALUE SOLUTION OF THE TIME DEPENDENT PROBLEM*..*... 40 I General Formulation..,.,. *.....,* 40 II. Determination of Eigen Values xo and 1.............. 46 III. Relation Between the Diffusion Cooling Coefficient'and 1 48 V. a a - a a - a a i a a * a i * a *' *:r t. 0 s O. * * * * 4 IV. Determination of.,....,.,.............. 52 V. Determination of M2,.. *-*t.*.. * * 54 A. Diffusion Cooling Method.,...,....... 54 B. Infinite Medium 1/v Absorption Method,.......... 55 V. NEUTRON THERMALIZATION IN THE HEAVY GAS MEDIUM...*......* 56 I. Transformation of Integral Equation into Differential Equation.. aQ*..*. **- W 4 **..0 O..............w** 56 II* Determination of the Eigen Values...*..*..........* 59 III. Determination of the Eigen Functions....4.,..........* 61 VI* THE TIME DEPENDENT ENERGY SPECTRUM IN THE INFINITE MEDITUM 66 I, Details of the Numerical Method...,,....,, *.......... 67 II. The Infinite Medium, Zero Absorption Case..,*.....**,. 71 A. Asymptotic Energy Spectrum......,............,... 74 Bo The Thermalization Time......,.............-..-. 77 C. The Exponential Fit of the Distribution Function cp(E,t)o.,.............. O......*.,... 78 D. Comparison Between Analytical and Numerical Results.,,.. Z,^.....,-*........ *.4*9 ***9 *.. 84 III, The Infinite Medium l/v Absorberr.................... 85 VII. THE TIME DEPENDENT ENERGY SPECTRUM IN FINITE MEDIUM...... 88 I. Constant Diffusion Coefficient (Graphite)....*....... 88 A. Equilibrium Spectrum.... *.....*....,....,. 89 B* Average Speed.*.*...............~.....~......94 C. Average Eiergy.....**.................*.,.. 95 D,. Neutron Temperature...,,...,.... * *........ ~ 96 E.- Decay Constant......,....,.O...,....*.., I..V#.. 97 iv

TABLE OF CONTENTS (CONT' D) Page II. Energy Dependent Diffusion Coefficient (Beryllium)*.., 99 Ak Equilibrium Spectrum4..* *-* 0 * 99 BA Average Speed........,.... *o....... *.. ~ 103 C* Neutron Temnperature.....- *, -.........-........ *.... 103 D. Decay Constantt.*.. *,...,,,:,,........ o.*... 104 III, Diffusion Cooling Coeffic ient-,.. * * *,..............., * *o-. 104 A, Graphite (Constant Diffusion Coefficient.)...,...,. 105 B4 Beryllium (,Energy Dependent Diffusion Coeffici;ent) 106 VII:I. THE TIME DEPEDENDENT DISTRIBUTION IN THE PULSED MULTIPLYING.MEDIUM.......,..,. #~ *. *v. 4 *..o 4,, o..*, *.... *4'484* 108 I* General Formulation.***. * o.,.......,.... o..* 108 II* The Monoenergetic Case.,.,.*...,...4..*....*....49-4* 113 A. Heavy Element Case...a.. *,*i 4o..*.c, *, *...:*4.,4:. 113 B. Hydrogen Moderation...*..... 8,,.., o.,, *. *. *., 115 III* The Multienergy Case119 IV. Analysis of Experimeit1 s............. o...... ~.... 122 IX. DISCUSSION OF RESULTS. o a....o4,:.o.oo —.o. —...-,-* —.,-. 125 I. First: Eigen Value*..o,... 4,o.... -.;...,,... o0. 125 II. Determination of M2.4.o0.44.*4 o. a 4. 0a #a 4 o. -:, 126 III, Thermalization in a Heavy Gas Medium4. o..,,,,,,,,, 127 IV. -The Time Dependent Energy Spectrum in the Infinite Medium.**.*+. o -... o 0 a, 4.,44 4 o4 o: 4,...0., is 4 128 V. The Time Dependent Energy Spectrum in the Finite Medium,,...,...0* 4 O4 - a.,a.:,.. o Q.: o 0:o6 o:,: ~.~...* 129 VI. Conclusion*1.04 4a4 404:-8u a 444...o..4 *,a 131 REEERENCES *0, * - 44 48 * r 0 a a 133 APPENDICES APPENDIX I. PROPERTIES OF ASSOCIATED LAGUERRE POLYNOMIALS.... 136 APPENDIX II. TIME DEPENDENT ENERGY SPECTRA FOR DIFFERN T B2. 3 39 APPENDIX III, HIGtER SPATIAL MDDE.CONTRIBUTION i4.. - 4.... 144 APPENDIX IV. COMPARISON OF THE STEADY AND TIME-DEPENDENT ENERGY SPECTRUM IN THE FINITE MEDIUM.,*,,,,..*4*,.* 147 APPENDPIX V. DIFFUSION PARAMETERS OF BERYLLIUM AND GRAPHITE e,, 151

LIST OF TABLES Table Page 3.1 The Source Distribution S(E,ts)........................ 32 4.1 The Thermalization Parameter M2...................,... 50 4.2 The Thermalization Time Constant tth,............ 51 4.3 Comparison of o0 and X1 in Beryllium....*.,,,..**..**.... 52 5,1 Eigen Values in Heavy Gas Model........****. **.......*,* 61 5,2 Coefficients of the First Eigen Function.*.*....,...... 63 5*3 Eigen Functions in Heavy Gas Model..,......*...*,..... 64 6.1 Numerical Values for co and c1 in a Heavy Gas Medium... 80 6.2 The Values of kl/X2 for Different Absorptions......... 85 6.3 The Decay Constant for Different Absorptions in a Infinite Medium.n,.,...,.,,....,*....,..*-.*... 86 7*1 Equilibrium Energy Distribution (Graphite).**...b..*.. 92 7.2 Average Speed for Various B2 Values (Graphite).*,.*-.... 94 7.3 Average Energy for Various B2 Values (Graphite)..,...* 96 7.4 Neutron Temperature for Various B2 Values (Graphite)... 97 7.5 Decay Constants for Various B2 Values (Graphite)#..... 98 7-6 Average Speed (Beryllium)....,,.,,.,............,. 103 7.7 Decay Constants (Beryllium).....o..ooo............ 104 8.1 The Decay Constant for Different B2 Values..,,..,.*.... 122 8.2 Nonescape Probabilities..,,*..,, os,, e., *.,,,, 123 vi

I

LIST OF FIGURES Figure Page 3.1 Comparison of Time Dependent Distribution Function cp(x) for Mass A = 12.,..................... 27 3.2 Source Distribution S(E,1s') for Graphite (A = 12).... 31 4.1 l Variation of Xo/Xl with B2 for Beryllium.............. 53 6.1 Time Behavior of Neutron Energy Spectrum in Beryllium for Infinite Medium and Zero Absorption.... 72 6.2 Time Behavior of Neutron Energy Spectrum in Graphite for Infinite Medium and Zero Absorption. 73 6.3 Time Behavior of Energy Groups in Beryllium for Infinite Medium and Zero Absorption.................. 75 6.4 Time Behavior of Energy Groups in Graphite for Infinite Medium and Zero Absorption............... 76 6.5 Exponential Fit of the Infinite Medium Data for Beryllium................. 81 6.6 Comparison of Numerical and Analytical Eigenfunction cp1(E)....................................... 82 6.7 Comparison of Numerical and Analytical po(E)........ 83 6.8 Neutron Energy Spectrum at 300 4 sect in Infinite Medium of Graphite for Different Amounts of Absorptions.......................................... 87 7.1 Time Behavior of NeUtron Energy Spectrum in Finite Medium of Geometrical Buckling B2 = 18.5 x 10-3 cm-2 in Graphite........................................ 90 7.2 Time Behavior of Each Energy Group in the Finite Medium of Graphite of Geometrical Buckling B2 = 18.5 x 10-3 cm2-............................... 91 7.3 Equilibrium Neutron Energy Spectra in Graphite for Geometrical Bucklings B2 = 0 and B2 = 18.5 x 10-3 cm-2............................................. 93 7.4 Time Behavior of Neutron Energy Spectrum fin Finite Medium of Beryllium of Geometrical Buckling B2 = 7.18 x 10-2 cm-2..... 100 vii

LIST OF F'IGURES (CONT'D) Figure Page 7-5 Time Behavior of Each Energy Group in the Finite Medium of Beryllium of Geomettical Buckling B2 = 7,x 108-xl2cm2 o.,-* -.oo o0 b o o -2 oo 7.~6 Equilibrium Neutron Energy Spectra in Beryllium for Geometrical Buckling B2 0 and B2 = 7.18 x 102cm'2*4 102 APPENDIX 1.1 Time Behavior of Neutron Energy Spectrum in Finite Medium of Beryllium of Geometrical Buckling B2 = 2o96 x 10 -2 2o 140 11.2 Time Behavior of Neutron Energy S ectrum in Finite Medium of Geometrical Bucklings B = 3~96 x 10 -2cm-2. 141 II3 Time Behavior of Neutron Energy Spectrum in Finite Medium of Beryllium of Geometrical Buckling B2 = 5436 x 10-2cm-2 in Beryllium.,a.o..o**.o... *4***0 142 Ii.4 Neutron Energy Spectrum at 300 ~1 sec in Graphite for Various Geometrical Bucklings, B2,Q,.... 143 viii

CHAPITER I. INTRODUCTION The study of the behavior of neutrons from a pulsed source is of great physical interest. The wide use of pulsed neutrons to measure reactor physics parameters has increased the importance of theoretical studies. The parameters obtained. from these experiments are averaged over the neutron distribution existing inside the medium. Complete knowledge of the neutron distribution as a function of space, time and energy variables is required. to analyse and to interpret -the results of these experiments. The measurement of the neutron~ distribution as a function of space, time and energy variables in.all ranges is.a formidable task. The rigorous determination.of neutron distribution can be obtained analytically as well as numerically under certain physical approximations. A burst of high energy neutrons, when incident on a medium, undergoes three physical processes in succession: moderation, thermalization and diffusion. In the first stage, neutrons lose energy by elastic and inelastic collisions, This stage is called moderation of neutrons. In the second stage, neutrons gain as well as lose energy during the collisions and undergo thermalization. When the neutrons are completely thermalized., the exchange:of energy stops, and neutrons acquire an equilibrium energy disc tribution. In the last stage, neutrons.diffuse in the medium. The d.iffusion process is governed by the absorption.and the leakage properties of the medium. During diffusion, the shape of the equilibrium energy distribution is maintained. The process of diffusion is terminated by absorption or leakage of neutrons. Neutrons spend less time during slowing down than -1

during thermalization and diffusion, On the other hand, the distance.travelled is larger while slowing down. It becomes evident, therefore, that the study of the spatial distribution of neutrons is perferred for understanding the slowing down process. In order to understand thermalization and diffusion, however, time.dependent studies are desirable, Time-dependent studies of slowing down of neutrons above one electron-volt have been made by Ornstein and Uhlenbeck(36), Marshak(3), (50) (23) Waller (5), Kazarnovsky and others, Reference to their work shall be made in the third chapter, No references for such studies below one.electron-volt exist in the literature. This thesis.deals with' the study of time dependent energy spectrum below 1/4 -electron-volt. Analytical and numerical methods have been employed to study these problems. The characteristics of thermalization and diffusion processes, namely, the thermalization time constant, the thermalization time, the asymptotic energy spectrum, the decay constant, the average speed and the.diffusion cooling coefficient, have been determined. We briefly describe pulsed neutron studies undertaken by other authors, both experimentally and theoretically. I. Experimental Studies (4. 8) In 1953, von Dardel introduced pulsed neutron techniques to measure tho diffusion parameters: the absorption cross section and the diffusion constant, A beam.of pulsed neutrons is allowed to fall upon a finite medium. The integrated neutron density is measured as a function of time, and the decay constant is determined as a function of geometrical

-3buckling, B. The decay constant X is fitted, as a function of B as follows. A = U A tD- J-1-Be) (1.1) where: Za = absorption cross section D = diffusion coefficient v = speed of neutrons c = diffusion cooling coefficient. Zsv and Dv are averaged over the Maxwellian distribution. It was suggested by von Dardel that the energy distribution is not Maxwellian in the case of a finite medium, and that the deviation of the spectrum is determined by the diffusion cooling coefficient c. In the finite medium, there is a preferential leakage of neutrons of high energy, and, therefore, the energy spectrum shifts toward the low energy side, This decreases the mean energy associated with the spectrum, causing a phenomenon termed "diffusion cooling." Using pulsed neutron sources, the dififusion parameters of beryllium, beryllium oxide, graphite and water have been studied. Review articles by Amaldi(l) and also by von Dardel and Sjostrand, summarize these studies. Formula (1.1) is based upon the existence of a fundamental spatial mode in a medium characterized by geometrical buckling. The experiments of Campbell and Stelson(ll) suggest the existence of a single decaying spatial mode at long times. However, the need for an experiment which measures the time dependence as well as the spatial distributions as a function of energy has yet to be f-Plfilled.

In addition to the non-multipyinlg medium experiments, measurements have been made in multiplying medium by Campbell and Stelson(), (26) Kloverstrom and Komoto and others. II, Theoretical Studies Theoretical studies reported in the literature have been limited to the determination of the decay constant X and the diffusion cooling coefficient c. In these studies all authors have made the assumption that at sufficiently long times the decay of'the asymptotic energy spectrum in the pulsed medium can be represented by a single decay constant. Mathematically, this means: At t (-Etr = tE) C (1J2) H-rwitz and Nelkin(21), using perturbation theory, determined the diffusion cooling coefficient c for a heavy gas model. The value of c for graphite and for beryllium was calculated and found to be lower than the experimental value. (32) Nelkin repeated the calculations for c using a variational technique. He employed an energy distribution function c(E) given by the following modified Maxwellian distribution: IEME) Ma(E-e (1e3) In Equation (1.3), M(E) is the Maxwellian energy distribution. The energy is measured in units of KTm. The moderator temperature is given by Tm and the neutron temperature by Tn. The variational parameter P in the

-5trial function is defined as follows: p - Tv- TrT (1.4) -,h Nelkin obtained the following result for c: 2 -. C - (oi-Y<) \[7T) (U( FN.) (1.5) where: D = - M(E) is the Maxwellian.distribution vo is the speed corresponding to the most probable energy M2 is the thermalization parameter: 2 M2 = J f zs (E - E')M(E)(E'-.E) dEdE' 00 In deriving Equation (1.5), the transport mean free path was taken.as proportional to Ea. Singwi and Kothari (44) also derived a relation similar to 1.5, using an arbitrary energy dependence for the transport mean free path. According to them, c is given as follows: CI= C 0 (Ao ~ < ) ( i)(1.6) where: Ai(E) = axAtr(E E SXp(-E)dE Xtr(E) = transport mean free path. The above integral has to be determined from a knowledge of ~tr(E). For "the case of a constant diffusion coefficient, a is equal to zero, and A1/Ao is equal to two. Under these conditions, Equations (1.5) and (1.6) become identical.

-6The concept of neutron temperature was used to obtain both expressions for c given above. This concept of assigning temperature to neutrons has been criticized on physical grounds because "temperature" is defined thermodynamically only for equilibrium.distributions. Considerable advance has been made in neutron thermalization and diffusion theory by the application of Laguerre polynomials. This is dle to thie fact that the eigen fuhnctions of the heavy gas model differential operator are associated Laguerre polynomials of the first order, A large class of neutron thermalization and diffusion problems in heavy media were solved by Kazariovsky et al. (24) H'afele and Dresner(19) determined the diffusion cooling coefficient in all orders, for the heavy gas model. Their value of c is 33 per cent higher than the value obtained using Equation (1.5), for the case of a constant diffusion coefficient, Singwi( using the Laguerre polynomial expansion derived the same expression for the diffusion cooling coefficient as given by the variational method. Singwi's work represents a considers able advance in the theory of diffusion cooling, as it is not based upon the concept of neutron temperature, All the above results were obtained using the diffusion theory approximation, A transport correction for the diffusion cooling coefficient has been considered by Sjo'strand( 45) Singwi(43) and Nelkin(35). This correction was determined using single velocity neutrons; and it has been found to have the opposite sign of the multi-velocity correction mentioned above. The correction is given as follows: whetr = 15trnpr meanfeat2 where: Err.= transport mean free path,

The diffusion cooling coefficient has been related to the thermalization time constant (tth). We define the thermalization time constant as the time constant with which neutrons reach thermal equilibrium with the atoms of the moderator. On the basis of elementary thermodynamic considerations, Beckurts(4) gave the following expression relating the diffusion-cooling coefficient and the thermalization time constant: Ltt = lT Tnj= Thi (i8) aD/8Un represents the variation of the diffusion coefficient with the neutron temperature. Von Dardel(48) gave the following result for tth for the case of a monoatomic gas of mass M. where m = the ass of the neutron where: m = the mass off the neutron M = the mass of the moderator &b = the mean free path of the neutrons for rigidly bound nuclei. Nelkin(32), using the variational method mentioned above, derived the following result for the thermalization time constant, I =t T) DSo,io (1,10) In Chapter IV, we have derived a relation for the thermalization time constant which is identical with Nelkin's result, without using the concept of neutron temperature.

(2) Antonov et al( gave a two-energy group theory for the analysis of their pulsed neutron experiments in non-multiplying media, A theory of pulsed neutron experiments in multiplying media has been given by Krieger and Zweifel(29) In their treatment thermal neutrons have been assumed to be of a single groupo The problem of slowing down neutrons was treated by using the kernel methodo A two-gro-p analysis of the transient behavior of thermal fl1ux in a subcritical assembly has been made by Fultz(18) and (6) (25) applied to the experiments of Bengston et alo, Kirschbaum and (40) Reynolds o ITII Treatment of the Problem The objective of this thesis is the determination of the timedependent low-energy spectrum of a pulse of neutrons, and the study of the characteristics of neutron thermalization and diffusiono Analytical and. n-merical methods have been employed to carry out this study, Tfhe analytical method involves an eigen-value technique. In the numerical method, the time behavior of the neutron energy spectrum is generated for all times greater than the slowing down time. We describe here, in brief, both methods and the problems to which they are applied, A time-dependent problem can be reduced to a stationary problem by taking the Laplace transform, The transformed problem can then be solved by the eigen-value method, with the transform variable as a parameter, The complete determination of the neutron distribution involves the knowledge.of all the eigen-functions and eigen-values. In addition, the initial distribution at a specified time is needed to determine the amplitudes of the

-9eigen functions. Each eigen function is expanded as a sum of an infinite number of associated Laguerre polynomials of first order. (The choice of the Laguerre polynomial, however, is arbitrary.) By using this method, we shall determine the first two eigen values for any physical model. Determination of -the higher eigen values involvesa knowledge of energy transfer moments greater than two.. These are not readily available due to a lack of information about energy-exchange scattering cross sections. However, for the heavy gas model it is possible to determine all eigen values, since all of.the energy transfer moments are known. It.is not.always possible. to construct the.complete neutron distribution by using:two or three of the lower eigen values. These lower eigen values, however, help us in the study of the.characteristics of the last stage of neutron thermalization and diffusion. The first eigen value determines the rate of thermalization and the zeroth eigen value determines the final decay of a pulse -of neutrons in the finite medium. It.is of great physical interest to study analytically the first eigen value and its associated eigen function, Practical considerations limit the use of'the:eigen value method, although it is rigorous. As far as the analytical method is concerned, we shall restrict ourselves to the determination of lower eigen values. A multigroup method has been developed to obtain the time dependent low energy neVtron spectrum in a heavy gas-medium. A burst of neutrons of high energy is followed in time during the thermalization and diffusion periods. The initial source of neutrons at -the slowing down time is determined from Fermi age slowing down theory. Using this.source,

-10 - we generate the neutron energy distribution at different times. From these energy spectra, it is possible to study the characteristics of therzalization and diffusiono To obtain these spectra, the energy interval between 10 KT and zero has been divided into groups. The problem is solved for the heavy gas model, since this requires only the solution of the differential equation which is easy to handle, he electric analog computer was employed to generate these distributions, with tine as the continuous variable. Using -the nmerical methoad, we shall solve the infinite medium problem. The -thermalization time constant and the thermnalization time are determined fromn-the energy spectr-um curves, The results obtained analytically are compared with those obtained by the numerical method, The same numerical met1hod has been applied to the study of the finite medium problem. Two cases have been dealt with: the constant diffusion coefficient case of graphitee and the energy dependent diffusion coefficient case of beryllium. The time behavior of the energy spectrum for a number of geometrical bucklings has been obtained~ From this study, we attempt to answer the following questions. 1. Is there any equilibrium spectrum at the time of measurement in the finite medium? 2 At, what time does this equilibrium set in? In all the previous studies it has been assumed that at sufficiently long times the equilibrium spectrum does existe We test the validity of this statement.at the time of measurement,

-11The following characteristics of the equilibrium spectrum for a given buckling shall be studied. 1. The shape of the equilibrium spectrum. 2. The average speed. 3. The decay constant. 4. The diffusion cooling coefficient. The numerical method developed here is used for the heavy gas model. However, by similar use of the multigroup technique, the integral equation for the crystalline medium can also be solved. A medium other than the heavy gas requires the scattering frequency of thermal neutrons before the integral equation can be solved.. The determination of the scattering frequency, for media other than the heavy gas, is in itself a laborious t-ask. Thus, it is not undertaken here. A general formulation for treating the multi-velocity problem in a pulsed, multiplying medium is given. Campbell and Stelson's experimental data for water and uranium..mixtures are analyzed, In the following chapters, we shall attempt to study the problems, and answer the questions, that have been outlined briefly here.

CHAPTER IIT THE MATHEMATICAL F0M1UJLATTION OF'THE PROBLEM In this chapter we present'the mathematical formrlation which is used in the subsequent chapters~ The behavior of neutrons in a medium -can be completely described by the classical Boltzmann Equation~ This equation is based upon the neutron conservation principle. We shall write down the basic equation for the ne- tron distribution:as a function of energy, space, angle and time variables. We restrict our consideration here to the non-multiplying medium1 I. The Boltzmann Equation This eqzuation has been discussed in such great detail by (52) (14) Weinberg.and Wigner(52 and also by Davison(4), that only a brief discussion need be given here, The physical behavior of neutrons can be represented by the following:equa'tion:,, Q. t ) _ C([(E) t 2S:(E) + rL V g ~(E., t) + S(E.t) L'J-.Vj~cb(Ig9K~t) -&t U- ~ fJ t(E'g t)(E, t f) t) 4(El2E; J1 E d' (201) We define the following quantities: cp(E, Y, t>, t) dEd.rdS is the angular flux, i oe., the number of neutrons at time t of energy in the interval E to E + dE, having their direction — 3 - -- of motion along 2S in the solid angle interval Q to 92 + dQ, in the volume element dr.around r, multiplied by their speed v. -12 -

-13F(E' - E; 2' - >) dEd2 is the probability that the neutrons of energy E' along the direction 6shall have the final energy E in the interval E to E + dE along the direction of motion Q2 in the solid angle interval 2 to!2 + d.2, after the collision with the atoms of the moderator. We shall assume the collisions to be scattering collisions only. The above quantity is called the scattering frequency. S(E, r, 2, t) d.ddEdEt is the numiber of neutrons of energy E between E and E + dE in the volume element dr around r, along the direction 2 in the solid angle interval 2 to 2 + ad produced in the time interval between t and t + d.t., due to external sources. Ea(E) and 5s(E) are the energy dependent macroscopic absorption and scattering cross sections, respectively. The cross sections are assumed to be independent of space and angle variables. The rigorous solution of Equation (2.1) is in general not possible. It has only been given for a Few special cases. We shall examine the equation under different approximations. The derivation of the Boltzmnann Equation (2.1) is based upon the following assumptions. a) Neutron- neutron interaction is neglected. Since the volume occupied by the neutrons is very small compared to the total volume, we can neglect the mutual interaction between neutrons. Even for the high intensity sources, the neutron density is of the order of 107 or 108 neutrons per cubic centimeter. b) The neutron moves with a constant velocity between two collisions, along a straight line, A change in velocity occurs only at the instant of collision.

c) The collision of neutrons are instantaneous, and not delayed. Only in the case of fissionable material do we consider the delay effects, d) The medium is isotropic and homogeneous. e) The neutron density in the medium is very large, so that the effect of fluctuations can be neglected. II. Angular Distribution The problems connected with the angular distribution of neutrons lie in the domain of the transport theory, If the size of the scattering medium is large compared to the scattering mean free path of the neutrons, then diffusion theory can be used, We shall derive the diffusion theory equation from the integral Equation (2.1). We expand the angular flux, the scattering frequency and the external source in series of spherical harmnonics. 00 C ( t) S E( t Et) t —0 rn:-,.m 4Rn' F(E E ].p' ) F ( E -, J 3, ) x X 4 (2 4) We make the followlng assumptions. The extermal source s isotropic.

The scattering'frequency is a function'of the angle between the initial.and final directions,,o = Q' Q'. This assumption holds good in an isotropic, homogeneous medium. 1* = Cos 8. Cos i' + S SenCe c(m-') (2.6) For the sperical harmonics, we make use of the following relations. y i,n) tM) t (f~ Q) { L i =:omplex Conjugate of. -T: (2.7) L.). mr(QJ)? T, CJ) {Addition $2.8) Theorem {Orthogonality Relation. F,(E' i ) _= F(E'-E; ) d 2.1o) We:also write down the scalar flux and current. 2(r tC E(E,,i't) = Flux (2. 1) (E K J f(E1,J trI Current (2.12) We integrate Equation (2.1) over all angles. Making use of the above relations, we get _ TI t - j E,rtt)- SZ(E) tZd(E f(Ertr (-2-.13) + f FoQ(E) > U(E') +(E t)dE + S(E t)

-16In order to determine the neutron current J(E,r,t), we multiply Equation (2.1) by Q and integrate over all angles. Then we get $ UJ(f~ E, t~~fi) = ~ je) 2C(E) + ZS(E)i J(E) rt)+7~. L++4Ag + tJ (E') /(E4EC)'(E't)CE'( The external source term is absent, since the source is isotropic. We further assume 6 j(E, ) - O (2.15) 5* J2 Jm * EY S i) cs4 clL - 1 H (ER) (2.16) Equation (2.16) is justified because we propose to use only the diffusion theory approximation, i.e., we will consider only the PO and P1 components of the angular flux. Then Equation (2.16) follows from the: fact that J Q c J 4 =- $, (2.17) and $PJs< +1( J ) fR 0) (2.18) Waen Equations (2.15), (2.16) and (2.L7-18) are substituted into Equation (2.14) we obtain: 5 F (E E ) s (E') (E, rt)E -4 t) (.19 - fT.(E) +s(E)E J-E, t) - n

-17If we further assume that s'(E')J'(E',r,t) - ZsJ(E,rt), a basic assumption in the Selengut-Gortzel(22) and Fermi Age slowing down models, we obtain: 4-* T(E, t) L=-Y [ cE)+ s (E)+)- E)]N7(Er) (2.20) T(E t) - D(E) V ( E, 6t) (2.21), the average cosine of the scattering angle, is given by fFj(E'+ E)dE Equations (2.20) and(2.21) are sometimes referred to as Fick's law. D(E) is the diffusion coefficient. It depends upon the energy (E) and the scattering kernel of the medium F(E' e- E;c) through the mean cosine of the scattering angle po. We shall restrict our considerations to only the diffusion theory approximation. Equation (2,14) reduces to the following expression when the expression (2.21) for the neutron current is substituted into it. ej. (E,(t) - - [ y(E)h Z 5(E) - V()J00(Et) (2.22) Ati'L ~.;, ) 4'+ (6 t) It should be remarked that the above expression is exact, except for the expression of the neutron current, which contains the assumption of diffusion theory. III. Space Distribution It is possible to expand cp(E,,rt.) into a set of orthonormal functions of space variable. We have to assume that the neutron distribution

goes to zero at an extrapolated distance from the boundary, The extras polated distance is assumed to be constant''with en'ergy and.given'by'an average value'. We expand cp(Er,- t) and S(Er,Ft.)) as follows: )(E, re) =(E t) r):o (2.23) S(E rt) =r V %j(Y) =- 4) {B2 - Geometrical Buckling I? N YI l n (2.24) F(QY) = at the extrapolated boundary (r). where n represents the number of the spatial mode. Substituting expansion (2.23) into Equation (2.22) we obtain`5- ~ ir~ C a ~~E,t ) =- q( E_ t E5(E),tD(E) 3Q (at) E/ For the fundamental spatial mode n equals zero. If we consider,only the fundamental spatial mode, we have the following equation. L ~) Eat() - 4 - [5(E) ts)t) B i(2t) (2.26) t a E+ J~ z4(E') Fc(E4E)C(-E)dE'+ CSt) Equation (2.26) is also true for any n-th spatial mode, and can be used to obtain the neutron distribution for that spatial mode, provided we use the 2corresponding Bn corresponding Bn.

-1 9We shall later use the above equation for studying the time behavior of the neutron energy spectra in finite media characterized by a fixed geometrical buckling. IV. Time Dependent Problem The time dependent problem can be treated in principle by taking the Laplace transform with respect to the time variable. The Laplace transform of Equation (2.1) is- as follows* [,, ( )4S) f (E,,) =:~. _-)X9(; I~t+i @S d, (2.2y) -' —.j ~,J. r,, E 4 4rO,2 / d P. E'-, -' -T t ~E, ro) where: cp(E,r,Q,s) = Laplace transform of angular flux. 00 - -st = f cp(E, r,,t) est dt S(:,r,R) = Laplace transform.i o.l. sl.Lrce. 00 = f e-st S(E,r,,t)5(t) dt (E,r,2,t = 0) = cp(E,r,,O) = Angular flux at t 0 O =0 The solution of Equation (2.27) can be reduced to the equivalent stationary problem by replacing the l/v absorption term by an effective absorption term which includes the transform variable's'. The effective absorption becomes a sum of two terms. iOi = YQ+ R U-5~ ~r

-20-. Though in principle the'solution can be obtained, yet the problem of inversion of the transform is not a trivial one. This problem shall be discussed in the next chapter, for the case of hydrogen. If the source term is represented by a Dirac delta function -.r-f in energy by equating s(E,r,Q) to s(r,Q) 8(E-Eo), then for an energy smaller than Eo, Equation (2~27) becomes a homogeneous equation. The problem reduces to an eigen value problem with s as a paramneter. This gives the solution for the Iaplace transform of the angular flux and the angular flux itself, as follows: () r,, o e ~ 4J7$ E (s,~"m) (' 2') F.(Eirrr ) = oE (I.2.29) km is the eigen value associated with the eigen function m(r4,E,Q)o We shall apply the eigen value method to Equation (2.26) in the fourth -and fifth chapters. We shall assume the Am corresponds to an operator half bounded from above, i,e,, set- of eigen values for m = 0,1..., the zeroth eigen value is algebraically larger than any other eigen value, so that at long times it predominates, In the case of non-multiplying or subcritical systems, the neutron distribution decays in time. Am is negative for these cases. Equation (2~29) can be rewritten as follows: (9e (A>Qt) - tr ~Atv r- " E'(2.29) We shall make use of the above results in the fourth chapter.

V. Energy Ind.ependent Solution Upon integrating Equation (2.26) over -all energies, it reduces to: 00 00 00 S L a~e( t~t~de = - hEI SCI1_ + JXE) 8(9 ) -fj)ES)c/& + Jf 2?(E') I(Et) F(+E) E /o(2.30) oo 0o F apE)ul (Edt)ro so5ce, E) S (Ett) F(E=(E/)dE E o ~C/...(231) We further relate the flux cp(E,t) with the neutron density n(Et). ~ ( E t) C E nCE~t) (2.32) Substituting Equations (2( 31) and (2.32) into (2.30) Fr =a -ple nero.(E) S (Ej Ct) (2 33. The solution of the Equation (2e33) is given by l(t) = rd e (2.34) The decay constant is given by: A = To) t D(E) L 8 (.2-35) --- b S ntQt) 2 t E) Lo c hf heEqtion (2.)ie

-22n(t) = J (Et) d E (2.37) -- J (Et)de LT(t) - _-.....(2.38) J q(E,(t)d For constant diffusion coefficient and 1/v absorption, the decay ~constant X reduces to;A = z. <Vcl + D' 5"' (2.40) The determination of the diffusion coefficient D and the absorption cross section.ao from the measurements of the decay constant X in the pulsed assembly requires the knowledge of the average speed V, We shall generate cp(E,t) for a few geometrical bucklings, B2, and obtain 7. We note that since -V is a function of time, there will not be pure exponential decay until n(E,t) [and thus V(t)] attains its asymptotic value. VI. The Scattering Integral We shall discuss the scattering integral in the next chapter. It must be pointed out that the complexity of the integral Equation (2.1) arises from the complexity of the scattering frequency F(E'-> E;'; I ). The physics of the problem is determined by this kernel. A detailed theoretical or experimental knowledge of this kernel is essential for

-23 - the rigorous solution of the Equation (2.1). For the numerical solution of the problem, we decompose the scattering integral into two parts. f EI ( E') (E, t) F(E-E)E (2. 41) E:= sE 0 Tbt/(E.))ie (T )4E ofJ C+ E/ ET o where: ET = thermal cut off energy F (E'-E) dE = thermal neutron scattering frequency thermal F (E'-'E)dE = fast neutron scattering frequency. fast We shall evaluate the fast neutron scattering integral at the slowing down time (ts) in the third chapter. We also discuss the thermal neutron scattering integral in the next..chapter,

CHAPTER III. FAST AND THERMAL.NEUTRON SCATTERING We shall discuss fast and thermal neutron scattering, in detail, in this chapter, For -the study of the slowing down of fast neutrons, we can neglect the low energy -effects and assume the atoms of the moderator -to be free and at rest. The lower limit for the validity of the slowing down model is taken to be about 10 KTE, which corresponds to 0..25 ev. Below this energy, one muast take into account the thermal motion of the atoms of the moderator, the crystalline effects and the chemical binding of the moderator, In the low energy region below 10 KT, the neutrons lose their energy due to elastic as well as inelastic scattering with the atoms of the moderator and with the medium as a whole. We shall speak of fast neautron scattering as the slowing down process and of thermal neutron scattering as thermalization, These are two distinct processes valid in two different energy intervals. In the slowing down process, the neutron can lose its energy only by elastic scattering with the atoms of the moderator, Inelastic nuclear scattering is only important.for high energy (>IMev) neutrons, therefore it will be neglected. In the thermalization process, the neutron can lose as well as gain energy due to elastic or inelastic scattering with the atoms of the medium. The physical state of the medium plays an important role in neutron thermalization. I. The Slowing Dow-n Process The slowing down of neutrons in heavy -element and hydrogenous moderators constitute two separate problems. A neutron can lose all of -24

-25 - its energy in a single collision with an atom of hydrogen, but requires a large number of collisions wFith an atom of a heavy moderator, before it can lose all its energy.. No attempt is made to review slowing down., theory here. We discuss two specific problems. 1. Determination of'the heavy moderator fast neutron scattering integral. 2, Determination of the space-and time-dependent slowing down.solution for hydrogen. The fast scattering integral will be used in the following way. A burst of fast neutrons is inserted into the assembly -at t = 0, then at some later time of the order of the slowing.down time ts, these neutrons will begin to enter the thermal group. Those last collision slowing down neutrons arriving at energy E < ET -will serve as a souxce for the thermalization problem. Clearly those neutrons will arrive from the collisions at some energy between ET and ET/l, so we will consider the fast integral between those limits. Furthermore, this "source" for the thermalization problem must vanish below E = -cET, The scattering frequency for the elastic scattering of neutrons is: F(E YE)clE = clEE. (5.1) This scattering frequency is obtained under the following assumptions: 1. Neutrons lose energy by elastic scattering which is isotropic in the -center of mass system. 2.. Atoms of the moderator are free and.at rest, Low energy effects are ignored,

-26A. The Heavy Element Scattering Integral (A >> 1) We determine the following integral encountered in Equation (2.41). (3.2) S (Et) = j s(E') c( E,t) F ( E)E (5.2) The ranges of energies for El and. E as indicated in the discussion above, are as follows: E7 E E E >A X E T Substitut.ing -the scattering frequency(3,1) into Equation (3.2), we obtain S, ( t) Et F (3) E n E'( I- o<) The determination of S(E,t) involves the knowledge of cp(E',t) as a function.of energy E' at time t. The time dependent slowing down (31) ((48) (46) problem has been studied by Marshak, von Dardel, Svartholm( (25) (50) (16) Kazarnovsky, Waller, Eriksson, and others. All these authors have attempted to obtain the velocity distribution of neutrons as a function of the dimensionless quantity vt7Es in a non-absorbing, homogeneous, infinite median of heavy mass. An excellent review of most of these (1) studies has been made by Amaldi( The basic physics of'these studies is the same, but they differ in the methods employed to obtain final solutions. We give in Figure 3o1 the velocity distribution curve for a moderator of mass equal to 12, as given by Eriksson(l6) based upon Waller's(50) exact solution. For a very large mass (A > 1) the velocity distribution given by Kazarnovsky(3) reduces to the following distribution.

40.0... ERIKSSON'S CURVE 30.0 z a: GA US SIAfA 1 20.0 0 5.0 10.0 15.0 20.0 Figure 3.2 Comparison of Time Dependent Distribution Function cp(x) for Mass A = 12.

-28T(-t) =j "1E eF A (z-l) Y (3.4) where: lr(Z) = ~(yt) (LkL) ~= f-id z = _s C EL u = lethargy -corresponding to energy E. ts(E) = slowing down time corresponding to energy E. When converted int.o energy -units, the above equation reduces to the following expression0 We -co:mpare the Gaussian dist~ribation with Eriksson's exact -curve in Figu.re 31ol In conrtrast2 with the Ga.lassian distribution, the Fermi age slowing down model gives a line distribution in timeo According to this model cp(E',t) can be given as follows: 4 CE,t) /.i (t-ts(EZ (i36) We determine cp(E' ) by the time independent Fermi age model. The elementary calculations yield q(El ) as follows.,.,.. / -'~ -/ ( E_.'J..,L(o) p e t E[' ~ E j (3.7) In Equation (357) the following terms need to be defined: S(Eo) the source strength-of neutrons of energy Eo. p' = the resonance escape probability for;the energy interval Eo to E',

-29e = the non-escape probability factor = logarithmic decrement = the Fermi age for the energy interval Eo to E. The Dirac delta function 5(t - ts'(E')) in Equation (3.6) is characteristic of the Fermi age model. We evaluate S(E,t) based upon Equations (3.5) and (3.6) for t equal to ts, which is defined as the mean slowing down time for neutrons having a final energy between ET and.ET/l. 1. Determination of S(E, t) by the Fermi Age Model We obtain S(E,ts) using the approximations given by the Fermi age model as' mentioned above. S~(E, F) = _s~ Lt - ts) J zc )d (3.8) In writing Equation (3.8) we have made one more assumption. The slowing down time of neutrons having final energy in the interval between ET:and ET/a does not change with energy. It can be represented by ts, This is valid, provided ET/a does not differ appreciably from Et. It is true for heavy elements, for which a is close to unity, We also assume that the resonance escape and non-escape prob2, abilities, p' and e do not change in the above energy interval but this is an excellent assumption. We integrate the right hand side of Equation (3.8)........... CC-...... E(3.9) a CI-( E:T

-302, Determination of S(Et ) using Gaussian Distribution Substituting the Gaussian distribution, as given by Equation (3.5), into Equation (3.3) we obtain S(Elts) = ~r v i 3ex1- 0 tE' -I) ts( )E( ) ET T5 In order to integrate Equation (3.10) we make the following substitutions. (3.11) ( EZ yi 1 ) = X (5.12) x( Eif) S (Et) ir ( ~) etxe& ~ > x (5.15) x ( ET) We use the following formulae to integrate Equation (3513).. (s-.l1) S'- X5e (5.16) j~ X\ a= - -J (1Q- (15) StE. j — =:: 9 2 9' q ( T. - t3 tSi -/Q() ts) M O'C-~tsCE/ct t)S'ts`f; t5~~~~~~~~~~~~,~

-311.0 FERMI AGE THEORY --- GAUSSIAN DISTRIBUTION / Graphite (A = 12). z 0 0 0 0.5 0.75 1.0 Graphite (A = 12).

-32Values of S(E,ts) are listed in Table 3.1, for a few values of energy, obtained by using Equations (3.9) and (3.16). In these calculations, A was taken equal to 12 and ts equal to the mean slowing down time between ET and ET/o, TABLE 3.1 THE SOURCE DISTRIBUTION S(E,ts ). E S(E,ts)Gaussian S(Ets)Fermi age Et 1.000 1.000 0o 9Et 0.818 0.718 o,86Et 0.510 0.588 0o8Et 0.362 0.370 Et o0 O In Figure 3.2 S(E,tS)Gaussian and S(E,ts)Fermi age are p for the sake of comparison. The basic features of the shape of the energy distribution are almost the same in both curves so that the simpler source, determined by the Fermi age model, and given by Equation (3.9), is used for the generation of the time dependent spectrum. It must be pointed out here, that the behavior of neutrons for times longer than ts is governed by the physics of the problem and is independent of the details of the initial source. The neutrons lose memory of details of the source after a few collisions. B. Discussion of the Space and Time Dependent Slowing Down Solution for Eydrogen The time dependent slowing down solution in a non-absorbing, infinite homogeneous mediutm of hydrogen was given exactly by Ornstein and

-33Uhlenbeck(36) We modify their solution for the case of a finite medium having 1/v absorption, using the diffusion theory approximation. The neutron transport equation, as given in the second chapter, is as follows. L B_4E~r,-O) - / < + s, ~Ltt) D7 ( 2,,. -~ ~ (2.22) ~ o E/ c s I EEj F E E) QS(EI') +S(2t) Assunming the source to be a.delta function in space and time variables, we get (E,, t jr) t) 5 (E) 7 We take the Fourier -transform of Equation (2.22) with respect to the space variable, obtaining t - ) ( K,t) tie,K + F. t) ~ b~+t o / )( ) (35118) + S(t),S() L+ i2 <;xo F EE (E.; K.t) _ |ur(ert)Q dr =Fouiier transform (3.19) < -= Transform variable On taking the Laplace transform of Equation (3.18) -vAth respect to the time variable, we get 1 %(4' ~()l ) -; DKPM (E.I, K S) { "r.10' I' K3 4)E ( S.20E) L -Lr ~!_ E'~S ~'~)c +

Rearranging Equation (3.20) g _ + _C~.QL) +'(.D)IK< C K( 4,s) 1 ~ E) ~,)+ (-S) (3.21)... J -, ) ~& + (KE) On differentiating Equation (3.21,) we obtain a differential equation. X E) =S'E) (2_ 2).2 where X = le + t ) The solution of Equation (3.22) can be given in a simple manner.,/ We assume the source S(E) equal to 6(E -Eo). Transforming X to' (E,k,s) we get E(IK,.S).......... (Irt >OtL s3t i s E +lbg I- zLt(Eo)t eL'S Eo) + ok l [-T-/' J- t ~ &' J,) (3.24) E S/lr/ +'EFn+ V+' KV rs I s., "i,("\)+.,:~. _., _. For energies E less than Eo, Equation (3.24) reduces to the following equation..!...s..;r'3.25)., - + z, _ F.......' t. K,. _,

-35The following comments about Equation (3.25) be noted: 1. For Za = 0 and K2 = 0 (infinite, homogeneous medium containing a uniform distribution of source): Equation (3.25) reduces to the expression given by Ornstein and Uhlenbeck.(36) 2. For the time independent case (S = 0); this equation corresponds to the solution given by the Selengut-Gortzel model as discussed by Hurwitz and Zweifel (22)and also by Simon.(4) 3. If the spatial dependent part can be separated from cp(E,r,t) then the Fourier transform variable can be shown to be equal to the geometrical buckling of the medium. The Fourier Laplace transform of the slowing down kernel is related with the neutron flux as follows: /L IS*~dE 4 la,(E )+$;gD l C tE (2 A~p L 9 E 8 L S+ (E) Ts(Ff i+ tY (e,,S _ (3*26) The above equation can be written in words as follows: 2 q (E,K,s) = (First Flight Factor) x (Resonance Captive Escape Factor) x (Non-leakage Factor) x (Time Dependent Factor) (3.27)

-36II Thermal Neutron Scattering Neutron scattering in the thermal energy region depends upon the low energ: effects mentioned earlier, namely, the thermal motion of the atoms of the moderatory the chemical binding effects, and the crystalline nature of the medium, In the thermal energy region, neutrons can gain as well as lose energy, This process is characteristic of the thermalization phenomenon, If we consider only thermal motion of the atom and neglect.chemical and crystalline effects, then t the thermalization is carried out by means.of elastic scattering only". When the chemical bond and crystalline effects are taken into account, then inelastic scattering also becomes important* The neutrons lose or gain energy through exchange of quantas in mblecules and phonons in crystals. The scattering of thermal neutrons is a big subject in itself and covers a wide range of problems. A number of Amali(l) Singw,Zemchanu(27be(5 excellent articles by Amaldi(l) Kothari and Singwi(27) Zemach and Glauber(5 and Nelkin and Cohen(34) and several othei authors exist. We shall consider the effects of low energy neutron scattering upon the neutron transport problem, For l/v absorption, the total scattering cross section and scattering frequency are the only parameters which,are dependent upon the low energy -effects. The thermal scattering term of the transport equation can be given as follows. J S(E) i(tE) ( E E)dEl S E5(E) q() 3.28) We express cp(E,t) as a product of two terms twhere M( M ) is the Maellian (3.29) where M(E) is the Trxelisn d9stnibzalion,

-37The total scattering cross section Es(E) can be given by the following integral T:(Ej -,r 2stL) beE-r E'd d(3.30) We -make use of'the Detailed Balance Theorem of statistical therrmodynamics. This theorem gives the balance between two small energy intervals for'the equilibrium distribution. Accordingly, the following -equality holds: M(E) Is.E) (E —,,r- - M(E) ZsCE) FFE.E. j (3.31) We substitute Equations (3.29), (3.30) and (3.31) into Equation (3.28), then the following result is obtained. f WIEl) F E E)4 (Et) dIE - S(E) E(Et) (3.32) ~E =' ~ J M(i E) _ CE) F(E.PE*a C[ Lt)-V&t)2K E' Expanding 4.r(E',t) in a Taylor,series about E, we get Cf~ C E el = tE,t, +(E Ej,) +( A E ) VZ - On substituting the above expansion into Equation (3532), we obtain the following result for the scattering tenm.,4sE"J k(Et) F/E-")i 4 - - sbl) q(C)t ~E/~~ K-QI i b 1< s}; (3.3 4) "'" "':'~ I ~,;, (. The K-th energy transfer moment is d.efined as follows A V__ __.._s _ (. 35) In order to determine Ak integrals, we need the energy exchange (partial) scattering -cross section. The physical approximations involved

-38in the determination of this quantity defines the basic physics of the problem. We shall discuss briefly the energy transfer moments. The main interest of the reactor physicist lie.s in the determination of the energy transfer moments represented by the integrals AK(E), It is these integrals, which appear in the neutron transport problem. Determination of these integrals involve the knowledge of the partial scattering cross sections. Zemach and Glauber(55) have given a general formalism for the determination of these cross sections for (27) a large class of scattering nuclei. Kothari and Singwi have extensively reviewed the interaction of thermal neutrons with solids. All of these formalisms are based upon the concept of Fermi pseudo-potential) introduced by Fermi(17) to explain the chemical binding effect of the hydrogen atom bolund as a harmonic oscillator. Pseudo-potential is characterized by the Dirac Delta function. In addition to Fermi, Vaughn and Cohen (47) have also determined the scattering cross sections for the atom bound as an harmonic oscillator. This model has been found to be adequate for solid zirconium hydride. The quantum of vibration is of 0.13 evo No attempt is made here to discuss the various models employed to determine the partial scattering cross sections for liquids, solids and polyatomie gases. In the case of a monoatomic gas3Wigner and Wilkins(53) derived the partial scattering cross section, assuming that the velocity disctributioh of the atoms of the moderator has a Maxwellian distribution. They also assumed the scattering cross section to be constant and independent of relative speed between the neutron and the atom.

-39(9) Brown and St. John(9) modified the Wigner-Wilins kernel by taking into account the effect of the relative speed upon the scattering crosssection. A general method for obtaining the thermal neutron scattering frequency has been given by Osborn,(37) It must be pointed out that, under the assumptions involved, the Wigner-Wilkins scattering kernel is exact. The energy transfer moments can be easily obtained, provided the heavy mass approximation is made and terms of order (1/A)2 are neglected. Under this approximation, the integral equation is reduced to the Wilkins(54) second order differential equation, In this thesis, we shall deal with the heavy gas differential equation. extensively. An excellent review article summarizing the neutron thermalization studies dealing with the monoatomic gas has been given by Cohen. (13) For moderators of greatest interest one integral quantity, defined by Nelkin(3 33) as the thermalization parameter. has been determined. This parameter is given by the following integral o( X Mi = j mJ ME) ) CE-E) FLE E)4EdE (336) o 0 It will be shown in the fourth -chapter that neutron thermalization is governed by M2* In the expression for the time constant of neutron thermalization, M2 appears as the basic physical parameter, which determines the rate of thermalization.

CHAPTER IV. EIGEN VALUE SOLUTION OF THE TIME DEPENDENT PROBLEM The use of the eigen value method to solve the time dependent problem was pointed out in the second chapter. A general formalism for determining eigen values governing the time behavior of neutrons, in any physical medium, is presented in this chapter. We shall determine the first two eigen values; the zeroth eigen value governs the decay of a neutron pulse in the diffusion period, while the first eigen value is associated with the last stage of neutron thermalization. I. General Formulation The neutron distribution given by Equation (2.26) for a single spatial mode shall be employed to study the time dependent problem. Rewriting Equation (2.26), we get 1 }t - [~C %( E)*A D(E)5B] ~(E;t) r at + S(E) (E t) F(E/- E) dE (2.26) - Js(E') p(t) (Et) We expand cp(E,t) as a product of two functions as follows. iE~t) = tlts ) ME) (4.1) M (E) clE - E ci E Maxwellian distribution(4.2) The neutron source S(E,t) is pulsed in time and is of the following form. S'CEt) = EB) 8(E-Eo) ift) (4.3) For energies less than the source neutron energy Eo, the source term disappears. On substituting Equations (4.1) and (4.2) into Equation

-41(2.26) we get: t - L (E) ~t (E) 2B (Et) M(E) ~ J 2(E) F-(E.'t)frE' F(E E' - CEL)'(E,t) WACE) Scattering terms can be represented in terms of the energy transfer moment integrals Ak(E), which have been defined in Equation (3535). Using the results of Equations (3.34) and (3.35), we represent the scattering terms as follows.' Z (E" J (E~t) MI(E') Fr (E ) E sE)(Et) M(i) 4 5 8 -1 a!| iMSE)E) IE) Recapitulating, the definition of Ak(E) is given by Equation (3-35) as follows. A E ) Br (CAI>!s)tK P) F(E E')dE I _AK- -s C() We substitute the result given by Equation (4.5) into Equation (4.1), getting 1L K(EA M(;) = - L (hE) +D(E.) I' tLEt) M(E) ir gt i-;5 t~L it 3r A d iSM(S) MC (E-E) Ak(F) We make the following ansatz: ~(EL) =: E (4.7) cp(O,t) = O Neutron distribution for energy zero is equal to zero at all times, cp(: Xt) = O Neutron distribution for energy Eo is zero at all times,

-42 - In the above ansatz, the eigen function corresponding to the i-th eigen value is given by the sum of the associated Laguerre polynomials of first order. L(1)(E) is the n-th Laguerre polynomial of first order. ki is the corresponding eigen value. Substituting the above ansatz into Equation (4.6), we get A h-Zt E, L (L(E) M(E) ('4.8) - - Z a-X:'tL C. ahc | (Zs(E)t-D Sz) L~ (E)MFl.) ~ L I (Ic f4) 1 L+ (E) K1(E)M(L)(E{(_. Multiplying by L )(E) and integrating over all energies, we get | Cla -- L (IE) M (E) C E? n.[r t-.., + g >Ac( (E;) L4) eeE) dE We define the following integrals. J'~ ('O I lcn2 (4.10) I Ln (E) Mf(E) Lwy (E) DBCE "- (ED) (4.11) Qo O* 00 EK _I J L L'(-) L (E) M(E) SA, (IE J4E (4.12),,( E) = — (4.13) Equation (4.13) is given by the 1/v a'bsorption law.

-43 - Lr -E (Jo (4.14) Substituting these integrals in to (4.9) we get XS~ [zT,3; & AUtD[Zl -4~ 4rnr -hr x2)m 4 J( 3j: a4.15) r) Ul. The above set of equations holds good for every coefficient ani. ni It will have a non-vanishing solution provided the secular determinant of its coefficients vanishes. The secular determinant is given as: li- a + Eqo) [2Xven +D5)r,,) - FM (4.16) In principle, the above determinant can be solved, provided we know all the matrix elements involved. Wn for all m and n are known. However, (DB )mn and F require a detailed knowledge of the scattering mn kernels.. For a'. heavy gas model, all the energy moments of the scattering kernel can be obtained by using -the Wigner-Wilkins kernel, therefore, F inn can be obtained for the heavy gas model. In the case of -crystals and liquids, extensive numerical calculations using digital computers are required. Even then, in these cases, the problem is a formidable one, because -of the lack of experimental data, and inadequacy of the theoretical models used. It is possible to solve the above determinant for the 2 x 2 case, since the F.n's involved are known. However, the above solution is also possible for a definite energy dependence of the diffusion coefficient, When D is proportional to v", then (DB2) can be determined. For the case of a constant diffusion coefficient, we shall solve the 2 x 2

-44determinant and obtain the first two eigen values Xo and X1. Hurwitz and Nelkin(21) determined X for the infinite heavy gas 0 medium using perturbation techniques. With the help of a variational principle, using a Maxwellian distribution as a trial function, Nelkin(32) also determined Xo0 This method, however, employed the concept of'neutron temperature.' The use of this concept has been found objectionable. Singwi(43) obtained Xo without using the idea of'neutron temperature.' Singwits expression is exactly the same as that of Nelkin's, and therefore, represents a considerable advance, Hafele and Dresner(l9) calculated o for a heavy gas model with a better approximation. As pointed out above, the determination of Xo and X1 involves the calculation of Wmn and Fmno The first few values of m and n have been calculated and are given below, A few associated Laguerre polynomials are also listed. We define the n-th associated Laguerre polynomial of the first order as follows:(30) Ln, ( ft = U L S ) L We list the first four Laguerre polynomial L. (E.) n Lo (E) = L1 (E) - j L- V () (4.17) La (E) - (-. (k -E + E ) 31 E, ~s,

-45The following matrix elements have been calculated. Wo = 0.500 JWol = 0 o.1768 wzz = 0. 4375 Wlo = 0.1768 4 Wo2 = o.o10834; wL2- o _0- 1914 W22 = 0.3984; W21 = 0. 1914 W30 Wo3 = 0.0781Cr W31 = W13 = 0o 1271at,W32 = W23 = 01928 JW33 = 0.37064 W33 also for CVI~M) Wn =T We determine the F integrals. 00~~~~~~~~~~ F = O o o~ Fol =O 0 ~ F = a J {;M(E)Z5() F( E e E") (E E ) E =lo - { s =0 By symmetry and detailed balance theorem.

-46Ill= -J { f- } { tv M(E) N5 (L) - E- E) E E'-E jE'} e, X - j J (E) pF(a) r((CE: LvE1 E ~/ / +J E J"(3 E ) Z5 (F ) FI) l E ) = - c ) ) ) F(E 0 o II. Determination of Eigen Values X and X We shall determine Xo and X1 from Equation (4.16) by considering m equal to zero and one, and n equal to zero and one, The secular determinant (4.16) reduces to the 2 x 2 determinant, stant with energyo 2 Ecr + L) C22) Li-(V C ( 1c (4.18) AIl = EUo- +! Wi6/oLN LD,' o9 o(4.20) \Ale k~ ) p lLLtL o al + 4+26't These expressions of Equation (4.18) give the upper limits twfor the two eigen valshall + 4S" l2 6c 1I

-47We make the following observations: 1. For the case of zero absorption and infinite medium. X = 0 0 1 4 ( kul~Viik -W o j2) This means that the neutron energy spectrum emitted from a pulsed neutron source reaches the final equilibrium spectrum with the decay constant x1, which is characteristic of the physical model of the medium. For this particular case, the final energy spectrum is Maxwellian and the approach to the Maxwellian distribution is governed by the first eigen value X1. 2. For the case of infinite medium and 1/v absorber. In this case the two eigen values are as follows. xM = d9 L)7S x = ace b t For small absorptions, X1 is much larger than Xo0, therefore -the final distribution is given by -the smallerof these two eigen values, i.e. k0o The final spectrum is Maxwellian, but decays exponentially with.a decay constant X, which is proportional to the amount of absorption inside the medium. For large absorptions, X1 is not large compared to Xo, therefore, both decay constants play important roles in the decay of the spectrum. The decaying spectrum is not Maxwellian in this case. 3. For the case of a finite, absorbing medieum, ko and X1 depend upon B2 and B4. From Equations (4.19) and

-48(4.20). 2 I o = ao + boB2 - CoB4 1 = a1 + + C1B4 If ao and X1 are plotted against B2, we find that the deviation of k and X1 from straight lines are exactly of the same magnitude but of.opposite signs. Coefficient Co is the well known'diffusion cooling coefficient.' If the departure from the straight line for ko is due to the shift of the zeroth eigen function associated with ko towards the low energy side, then for X1, it is due to the shift of the first eigen function associated with kI to the high energy side. The slopes in both cases are of positive signs but of d.ifferent magnitudes. III. Relation Between the Diffusion Cooling Coefficient and k. The expression for the diffusion cooling coefficient, for constant diffusion coefficient is as follows, CL = (WL ), (4. 22) It is possible to obtain a relation between the diffusion cooling coefficient Co and the first eigen value i1 in the infinite non-absorbing medium case. (4.23) Eliminating M2, we get 2-. 2

49) If we substitute the values for Woo, Wol1 and Wll we get 0 hZ 7T (4.25) If X1 is defined as the reciprocal of the thermalization time constant (tth) in the infinite, non-absorbing medium, then A1 - -s 2 t= (4.26) Nelkin(32) derived the same relation between the diffusion cooling coefficient and the time constant, with which'the'neutron temperature' approaches the moderator temperature in the infinite, non-absorbing medium. He assumed that the neutron temperature can be given by the following relation. Tl( t jo) 1910Jn1 p 1(4.27) Using the variational principle, with P' as the variational parameter and a modified Maxwellian trial function, he obtained the following result for the neutron temperature yt TV, C tt 7 1 +' e.(4 28) M cl A ya irc where y = " If y-1 is the time constant, then the relation between y1l and the diffusion cooling coefficient is as follows. - 77( C Nelkin(32) (4.29)

-50This is exactly the same relation as given for tth by Equation (4.26). The difference in the two results is due to the methods employed to obtain them. In obtaining tth, we have not assumed the concept of neutron temperature, which is the weak point in Nelkin's derivation. tth gives the rate with which the eigen function associated with k1 is decaying in time. The information about B will give the diffusion cooling coefficient Co and vice versa. Equation (4.21) gives XB for the case of an infinite, non-absorbing medium. k1 = 8.278 x 104 M2 (sec-1) (4.30) Thermalization time constant can be given as follows: tth = 12.08/M2 (mu sec) (4.31) Determination of ki or tth by Equation (4.30) and (4.31) involves the knowledge of M2, which has to be obtained from the energy exchange scattering moments. For a few interesting moderators, M2 has been determined, using extensive numerical calculations. Table 4.1 gives the list of values of M2 for a few moderators, TABLE 4,1 THERMALIZATION PARAMETER No. Moderator M_2 (cm'1) Reference Remark 1 beryllium 0.30 Singwi & Kiothari(4 ) GD = 10000K 0.o36 Nelkin(D = 9300K 2 beryllium oxide 0*20 Singwi & Kothari 3 graphite 0o068 Nelkin(32) 0D = 900~K Vibrations L to latice planes:* D = 25000K Vibrations |] to latice planes. 4 water 1.05 x 10-1l Nelkin(33) Effective mass = 18 0D= Debye temperature

For the heavy gas model, M2 = 4 5 zso where = 2/A so = free atom scattering so cross sections We shall list in Table 4.2 X1 (in kilocycles/second) and tth in t seconds, along with earlier determinations by various other'authors using different methods. TAB34E 4.2 THE TEHERMALIZATION TIME CONSTANT (tth) No. Moderator _X1 (kilocycles/sec) tth(U sec) Remarks 1 beryllium 24.8 40.266 Singwi & Kothari(27,) obtained 28 ji secs. 29.8 33556 based upon neutron temperature idea 2 beryllium oxide 16.56 60.40 R.. Bhandari et al.(8) obtained 67 f secs 3 graphite 5.629 177.65 K.H. Beckurt obtained experimentally 185 + 45 IJ sees. Antonov et alo., obtained thermalization time' 200 kL sec. 4 water 87000 11,94 von Dardel s9) result = 7 ~ sec. As indicated above, the thermalization time constant determined by the eigen value method agrees with the independent results of other authors. The case of the heavy gas moderator willnmt be discussed here. It is discussed in detail in the following chapter.

-52IV. Determination of x~ Since the experimentalists are interested in the ratio of ko/kl, in order-to estimate the contamination of Xo by X1 in their measurements, we shall calculate this ratio for a given example. We shall take the experiments of de SauLssure and Silver(15) for beryllium assemblies. We assume the following constants for beryllium. Eo vo = 0.288 x 103/sec M2 = 0536 cml.2 D~ = 1.25 x 105cm'sec. o = 0.288 x o103 + 1o25 x 105 B2 _ o.8655 x 105 B4 (4.32) 1 = 300 o88 x 103 + 1. 458 x 105 B2 + 0.8653 x 105 B4 (4.3355) We have the following results for various size assemblies listed in Table 4,30 TABLE 4.3 COWTIPARISON OF Xo ARD 1 IN BERYLLIUM No. B2 ( 02cm 2) X ( 103/sec) X1(103/sec) 1 0 0,288 30.088 9.57 x 10 3 2 1.05 1. 589 31.627 5.024 x 10-2 3 2,02 2.778 33.0o68 8.40 x 10-2 4 o 31 4.331 35 009 1.237 x 10'1 5 5o36 6.739 38,152 1.776 x l0'l 6 7,18 8o817 41.002 2.150 x 10'1 In Figure 4.1 we have plotted Xo/X1 against B obtained from (4,32) and (4o,33) equations for beryllium assemblies used by de Saussure and Silver(l5) in their pulsed neutron experiments. For the smallest

-530.22....... 0,20 0.18 0.16.. 0.14 Z O 12 0 iG 0.10 w I-. I o 0.08 0.06 0.04 0.02 0 1.0 2.0 3.0 4.0.5.0 6.0 7.0 GEOMETRICAL BUCKLING B2 (102 cm2) Figure 4.1 Variation of lo/X1 with B2 for Beryllium.

assembly of geometrical buckling, B2 equal to 7.18 x 102 cm 2, ko/Xl is equal to 1/4,65. This clearly demonstrates that the contribution by the first and.other higher eigen values are small compared to the zeroth eigen value and, thus, may be neglected, V, Determination of M2 The neutron thermalization parameter is characterized by the integral M2. Ma = JN g AE) 5site iE g J F(EES C) 4 H Nelkin(32) proposed the determination of its value from the measurement'of the diffusion cooling coefficient, We propose an absorption method for determining its val',eo We shall discuss the merits of the two methods. A. Diffusion Cooling,IMethod Equation (4.22) relates M2 with Co, the diffusion cooling coefficient, a 9 7. d Cc, 4 l UO PI) Wo) The determination of the diffusion cooling coefficient depends on the accuracy of the determination of the decay constant ko' given by Equation (4.19). A small error in o can cause a large error in the value of Co. The diffusion coefficient is assumed to be independent of energy in Equation (4.22). This is not rigorously true at low energies. The geometrical buckling, B2, depends upon the extrapolated distance, The lack of complete knowledge of the extrapolated distance introduces a small uncertainty in the determination of B2 Equation (4.22) has been obtained considering only the fundamental spatial modes The effects of higher modes have been entirely neglected.

-55 5 A rigorous determination of M2 should be based upon the.corrections for the finite medium effects introduced by the leakage and boundary0 The above method is handicapped due to lack of infornation for obtaining these:corrections, B. Infinite Medi]"m 1/v Absorption Method For a very large assembly, which nmay be considered as infinite, two eigen values ko and Xl were given-earlier, At= e pso f; a ob t'4tOhe(WlO.ay In the presence of a small I/v absorber, X1 >> X Then the final decay.of the spectrmn is governed by the zeroth eigen value ko, If we increase the 1/v absorption by using boric acid in water, for example, then X1 becomes comparable to Xoo The absorption term in the expression of X1 becomes.larger thanem the thermalization term0 In that -case the final decay of'the spectrm is governed. not only by Xo but also by X1l It is possible to determine0 Eand Xl bY analyzing'the decay curve.e Since'the'absorption term is known, 2 can. be determined rigorously0 This:method avoids the finite medii'm corrections but requires intense neutron -sources, as a large amounrt:of ab'sorption will be needed before the decay can be governed by X0 and al In the next chapter we shall give an estimate of the absorption required for determining M2 for'the heavy gas'model1

CHAPTER V. NEUTRON THERMALIZATION IN THE HEAVY GAS TMEDIUM In the fourth chapter, a general method to obtain eigen values for any physical model was discussed. The explicit expressions for the zeroth and the first eigen values were obtained, considering only the first three energy transfer moments. In this chapter, the problem of the heavy gas model will be discussed in detail. Wigner and Wilkins(53) derived the scattering kernel based upon the Maxwellian distribution of velocities for the atoms of the moderator. Since they ignored the chemical and crystalline effects, their kernel can be used only for the monoatomic gases. As will be shown in this chapter, energy transfer moments numbering more than three involve terms of the order of (1/A2), where A is the mass of the scattering nucleus. For the case of a heavy gas, therefore, the energy transfer moments higher than three can be neglected. This leads to the reduction of the integral equation into a differential equation. The differential operator of the heavy gas model has the associated Laguerre polynomials of order one as its eigen function. In the case of a medium having 1/v absorption and of known energy dependence of the diffusion coefficient, the results can be obtained to any degree of accuracy, and all the eigen values can be determined. We shall determine the first three eigen values Xo, X1, and X2. The zeroth eigen value was determined by Hgfele and Dresner.(19) The first eigen value, X1, and the associated eigen function, c1(E), are obtained here rigorously. I. Transformation of the Integral Equation into the Differential Equation We have already outlined the method for the transformation of the integral equation into the differential equation in Chapter TV.

'57 Rewriting Equation (4.6) Vr @ q(E,t) IvI(E) X L(E) + 6S J t) MCE) -1,_ K- E 5 ( -E) M(E) Awcy) Where Ak(E) has been defined previously as the k-th energy transfer moment. For the heavy gas model, we have the following values for the various energy transfer moments, based upon the Wigner-Wilkins kernel. ~S E) Ao = Z CE) Is ~E) Al = ( 0 -E) 25 (LE) 8~ = q 8 E: 250(5.1) CaslE) POK = ( Z1o k1'-e) 5 = ( R ) Substituting these moments into Equation (4.6), we get Wilkin's(54) differential equation. \ a2?(Evt) t" [E a) tt,) M(>-/ ) E - 5 {iJ'E E)+ P(E) ( e ) 2 X&(5'.2) The heavy gas differential operator has as its eigen values the associated Laguerre polynomial of order one. This result has been used by Hafele and Dresner(l9) in their derivation of the diffusion cooling coefficient. Kazarnovsky et al~ (24) coefficient. Kazarnovsky et al(24) also made use of this result in their studies. The energy portion of the differential operator in Equation (5.2) satisfies the following equation, according to Murphy and

-58 - Margenau 30) 9, (t) (5) L bet t#2 E ~ Lp(E) = -K L~E (E). where Ln(1) (E) is the associated Laguerre polynomial of first order. In' is the eigen value associated with the n-th eigen function. The eigen value takes the following values' zero, one two.... We make the following ansatz: -Ajt (,);(E t) =aCo (E) M(E) I~(4t) i,'qn(54).~(ot) - o, ) (.t) - o for t > 0 Substituting the above ansatz into Equation (5.2), multiplying by Lm(1) (E) and integrating over all energies, we get We have used the no Mtations introduced in Chapter IV A new We have used the notations introduced in Chapter IV. A new symbol, 6mn is introduced in Equation (5.5). It is the Kronecker delta function, 0o 0 ~ (I) J L) (E) Ee LE (E)4E = W n (5.6) The secular determinant given by Equation (5.5) is as follows Inm m4-L -sA )Vtl) o|5T)

-59II. Determination of Eigen Values For the case of a known diffusion coefficient, all-the eigen values can be obtained, since all the matrix elements are known. Hafele and Dresner(19) obtained the zeroth eigen value and the diffusion cooling coefficient'c' by expanding the secular determinant in the minors of the first row. This is possible because the off-diagonal terms contain B2 or its higher powers. It is not possible, for the case of higher eigen values, to do higher order calculations except by tedious numerical calculations. For the infinite medium having zero absorption, ko equals zero, therefore, the next higher eigen value, X1, becomes important. It governs the rate of thermalization in this case. If we use the assumption of zero absorption and infinite medium, then Equation (5~7) reduces to the following simple determinant. h n V1e, - -o (5.8) The two by two determinant is obtained from determinant (5.8), by putting m equal to zero and one, and n equal to zero and one. i _ ^oo A WoV Io W t -SO @ Lode (O 9) |I tssoo ( (5.-9) This gives two values for X ~o-=0 v /\ l — W14koo =.-t (5.10)

660If we substitute the value of M2 equal to 4e4o into Equation (4.21), it then becomes identical to Equation (5.10) for the heavy gas model. By expressing the secular determinant in the minors of the first row, we get matrices of the k-th order, for n equal to k. AW TT l ____ __ kP(okQ(ko Te0W LWoo 7T (nla - To(:/..o% where A rnn = A Wrmn ann n + A w,, ~/< 7T in = o~ 11. k2.^.K This leads to the following equation for X. J- l The above equation can be reduced to the following two equations. A- =0 (5-13) -Woo = Z - + 0>nn (5.14) X equals zero gives the zeroth eigen value, therefore, we have the final equilibrium distribution constant in time. Equation (5.14) gives the first eigen value, which determines the rate of thermalization. We have calculated the following eigen values in various orders. They are listed in Table 5.1.

TABLE 5.1 EIGEN VALUES IN HEAVY GAS MODEL Determinant X X1 (svoV) X2 (~oVo) 2 x 2 0 1.506 3 x 3 0 1.329 4.088 4 x 4 0 1.274 We have calculated the first eigen value, X1, in various approximations. It is guessed that the next value of X1 will be within two percent of the value for the four by four case. X2 has also been calculated to estimate its relative importance compared to X1.(38) III. Determination of Eigen Functions The eigen functions associated with the eigen values Xo and X1 will be determined. With the help of Equation (5.5) all the coefficients ani for the i-th eigen value can be determined. Rewriting Equation (5.5) E) = _ i l ioh J'(5.8) Zeroth eigen function cpo(E). The zeroth eigen value, Xo, is equal to zero. This leads to all ano equaln eual to zero. This leads to all ano equal to zero, except for n equal to zero. Qoo3o -t Fo ec O for, ~ L

-62The zeroth eigen function is ) CE) = ao E e This eigen function is the Maxwellian distribution. First eigen function (P1(E) The first eigen function associated with the first eigen value k1 is obtained by determining the various coefficients, anl. Using the 4 x 4 case, the ratios of various anl's to aOl's have been determined. The coefficients all/aOl, a21/a0l and a31/a01 were evaluated. With the help of these coefficients we construct the first eigen function in the following manner.'(E)... EQ Ed L, L(EE) ~:LA4 ] (5.~5) + al L-:oE) -Ja: (5.15)The following set of equations, obtained from Equation (5.5'), was used to obtain the above coefficients. - tOlt VJoo Q + 11 l qi Jo ai Wo3 - 9~l h9 OW1 t l(W- 1, - a 2 h t iI tr -QohIL iot%2(1. 11..@l+ 4k(,- i)- iLq4A31 i. W- 3J5i6) - "'4 T ToNIL + Tlo of __30 -_ qiL 1- aP1 k l.- + Q (-o ~~~'~'-~- ~ - ~ ~ —~- A qF7,1 W3a +? 3 1 - A-~o-~ The last three equations were used to obtain all/aOl, a21/aOl and a31/ao1. The values of these coefficients are listed in Table 5.2.

-63TABLE 5.2 COEFFICIENTS OF THE FIRST EIGENFUNCTION Coefficients Values all/aOl -2.127 a21/ao1 -0.7516 a31/a01 -0.3515 The first equation may be used to determine the contribution of these coefficients to the eigen function. The first equation of Bet (5.16) is as follows: Q Wo a~i 40o Qj AJOj W l a VJot - WOO (5.17) Substituting the values for the coefficients from Table 5.2 and the values of the matrix elements Wrn from Chapter IV, we get a h}<} Wo + Clqi oaR -i- A. hams>; F- - o v4r/7- (5 8) ol. abi C' —y and W-, oo The difference is about three percent. Substituting the values of the coefficients into Equation (5.15) leads to E(S) - Es al - 1,o13, 3,86oE - o s6 (5.19) +~ oo.2'3&]3

The adjoint of the above eigen function cpj(E) is 4t(E~) _ ao [_- 4 13 ~ 386o$5E -aO'61SE (5.20) I e < O 93 Eo, The eigen functions CPo(E) and cpi(E) have been numerically determined for different values of energy. In Table 5.3, these values are listed. TABLE 5.3 EIGENFUNCTIONS IN HEAVY GAS MODEL Energy (KT) po (E) - Ee aOO Cpl(E) Ee-Ef(E)a0l x a01 0.1 0.09 -o0.327 0.2 0.164 -0.535 0.3 0.222 -o.646 o. 4 0.268 -o0.686 0.5 0.303 -0.674 o.6 0.329 -0.624 0.7 0.348 -0.549 0.8 0.357 -0.455 0.9 0.366 -0.358 1.0 0.368 -0.255 1.5 0.335 0.1998 2.0 0.271 0.452 2.5 0.205 0.522 3.0 0.149 0.486 4.0 0.073 0.308 5.0 0.034 0.158

The zeroth eigen function is a pure Maxwellian distribution with a maximum at 1.0 KT. On the other hand, when one examines cpl(E), one discovers that it has a maximum at about 2.5 KT and a minimum at about 0.4 KT. The zero of the curve occurs at 1.3 KT. The eigen value method discussed in the fourth and fifth chapters does not give the ratios between the amplitudes of the various eigen functions. It does give however, the eigen values and also the shape of the eigen functions, accurately. In order to obtain the ratios between the amplitudes of the eigen functions, one must have the complete knowledge of the function cp(E,t) at some time. In the case of pulsed neutron source assemblies, this information is not available unless extensive numerical calculations are made to generate cp(E,t) at different times. In the next few chapters the results of extensive numerical calculations, carried out to generate cp(E,t), will be given0

CHAPTER VI. THE TIME DEPENDENT ENERGY SPECTRUM IN THE INFINITE MEDIUM The rigprous solution of the time dependent neutron energy spectrum can only be given by numerical methods. The analytical methods have not proved to be sufficient to answer all questions, due to convergence difficulties. The fast neutrons exhibit 1/E energy behavior, and the slow neutrons have Maxwellian distribution represented by Ee-E It has not been found possible to represent the complete behavior of neutrons in all the energy intervals by a single function. In the fourth and fifth chapters we used the eigen value method to study neutron thermalization. Though elegant to use, the eigen value method has its own limitations, It can not give the complete time behavior of the neutron energy spectrum, unless one has the following information, a) A complete set of eigen values XO k1.... kn...n b) A complete set of the corresponding eigen functions 0P' P'"' l CPn c) A boundary condition to determine the coefficients of the eigen functions. We used the associated Laguerre polynomials of first order to determine the eigen values. In the case of the general model, we limited ourselves to the first two Laguerre polynomials, due to the lack of information about the higher energy transfer moments. In the case of a heavy gas, however, we used the first four associated Laguerre polynomials of the first order, It is not easy to get higher eigen

values than X1, for crystals and liquids, unless the problem of the scattering kernel is rigorously solved. For the heavy gas model, it is possible to determine all the eigen values and eigen functions analytically. This, however, would require the expansion of the eigen functions by -a large number of Laguerre polynomials. The eigen value method is useful in determining the lower eigen values and eigen functions but is inadequate to determine the coefficients of the eigen functions due to laqk of accurate knowledge of the boundary conditions, We generate the neutron distribution function c(Eg,t) numerically in the following pages, for a-few interesting cases. With the help of this distribution, we shall determine the thermalization time, the rate of thermalization, the equilibrium spectra and their decay rates. The diffusion cooling phenomenon in finite media will also be studied. We propose, in addition, to test the validity of the statement, used by earlier workers in the analysis of pulsed neutron experiments, that "after sufficiently long -time the asymptotic neutron energy spectrum can be given as follows:" I. Details of the Numerical Method. We shall discuss, in this section, the details of the numerical method used to generate the neutron distribution function cp(Et). The electrical analog computer was used to obtain these distributions, We shall use the heavy gas- model for our problem. In the funaemental spatial

-68mode, the time behavior of the neutrons in the above model is given by Equation (5.3), Rewriting it, we obtain; X + E AE + - + D 6) I+E t) -a IWe solve Equation (6,1) using the multigroup method. For one group of problems, the energy interval between E = 10KT to E -= 0 is divided into fourteen groups, For the other group of problems we used twenty groups, in the above energy interval. In the latter case, the energy interval was divided in the following manner, Between E = 0 and E = KT, ten groups of 0,1 KT interval, and between E = KT and E = 10 KT, ten groups of 10 KT interval, were made, The functions, 2 E,~) and ~ (gJ for each group are expressed by the finite difference method in the following way, Let AEm be the energy interval between the (n + m) th and n th groups; let APm, be the energy interval between the n th and (n - m) th groups. If cPnm't; CPm and Pn+m represent the distributions for n-m, n and n+m groups respectively, then we have for ) t @-L MIa- rA. +i n (6.2) The above expression has been employed due to different group intervals between the (n+m) th and n th groups and the (n-m) th and n th groups

For ( ), we have the following expression: (%)EEJ Fh #71E ] E'bh-h (6.3) The above expressions are substituted into Equation (6.1). This gives a set of coupled differential equations. 46 En (4Em,4.M/) L E,) 3/2 _- 3. + EE | ~+ — ~' L E L E - EEI c Et tZrat EEto 2 ECe a N t The above equation is the representative of the n th groupg which is coupled to the (n+m) th higher group and the (n-m) th lower group. We have expressed the energy in KT units. The symbols introduced stand for the following quantities.'(i~ /~:rso t(6.5) E LralJ Ear rscD E - I

-70Equation (6.4) is the basic equation for generating cp(E,t). We solve the set with the help of'the analog computer, using time as the continuous variable.. We use the following boundary conditions for generating the solution. 1. Neutrons of zero energy have zero distribution at all times. cp(o,t) = o for t > 2. The neutron distribution function at time t = ts is given. s (ts is the mean slowing down time for the source neutrons to have final energy in the interval ET to ET/ ). In the third chapter, we calculated this distribution. Let this given distribution be represented by S(Et) _ oThe neutron distribution function using the above source condition S(E,t) will give the complete time behavior of the neutron energy spectrum, for times greater than the mean slowing down time ts. We have generated the distribution function cp(E,t) for beryllium and graphite in the following cases. A. Infinite medium, zero absorption B. Infinite medium, l/v absorption C. Finite medium, zero absorpiton a) Graphite - constant diffusion coefficient. b) Beryllium - energy dependent diffusion coefficient.

The asymptotic behavior can be used to obtain general results for any other heavy medium in the above cases, except for the special case of beryllium. This is possible since the initial source condition does not influence the asympototic behavior, which is given by the physics of the medium, In order to transfer the distribution function from the known case to the unknown one, we have to carry out the following transformation for the time and geometrical buckling variables. Let to, t, -.,). and Do be the time, geometrical buckling, average logarithmic energy decrement, free atom scattering cross section, and diffusion coefficient respectively, for the known medium. For the unknown medium) the corresponding quantities are tX _ h), Ax;S(q)ssand D X' The relation between the time scales for the two cases is: - tO (6X7) -When the diffusion coefficient is constant, we get a relationship for the geometrical bucklings by equating the leakage of neuttrons in the two cases, e )r LF f~ b (6,8) II. The Infinite Medium, Zero Absorption Case. In Figure 6.1 the neutron energy spectrum at different timesin the interval of t = ts to t = 350 + ts micro seconds has been plotted for beryllium. In determination of these spectra, twenty group calculations were used, For graphite we obtained the spectra in the time interval of t = ts to t = 510 + ts micro seconds, using fourteen groups,(39) This is plotted in Figure 6.2~

-72100 5 = f,+20/~sec 50 ff + 350p/sec- 20 -s / Xf. ~ + 100 iisec -i~(> t = tI + 50,+sec DI C 2 z i 0.5 - -- MAXWELLIAN SPECTRUM f =14 I Isec 0. 2 0.1 0.1 0.2 0.5 1 2 5 10 ENERGY (KT UNITS) (E-) Figure 6.1 Time Behavior of Neutron Energy Spectrum in Beryllium for Infinite Medium and Zero Absorption.

100.0 o 50.0 / -f - Ts + 510/,sec'. t - + OLsec / J, -f:~ + 50usec' -/ XW,, /_ _ I 20.0 10.0 o 5.0 m I0'r 2.0 4.0 0.5 0.2 0.4 0 2.0 4.0 6.0 8.0 40.0 ENERGY (KT UNITS) Figure 6.2 Time Behavior of Neutron Energy Spectrum in Graphite for Infinite Medium and Zero Absorption.

-74The energy spectrum at earlier times approximates a Gaussian distribution, but converges finally into the Maxwellian distribution. The degree of convergence is very good. After attaining the Maxwellian distribution, the neutrons retain this spectrum at all further times. In the non-absorbing infinite medium, the equilibrium spectrum is the Maxwellian spectrum. This is a basic result of the H Theorem of Statistical Mechanics and holds good for any scattering model. It is also consistent with the detailed balance theorem. A verification of this important result is furnished by the distributions plotted in Figures 6.1 and 6.2. A. Asymptotic Energy Spectrum The time behavior of a few energy groups is given in Figures 6.3 and 6.4. We have plotted the behavior of the individual groups in order to demonstrate the equilibrium state. Parallel curves indicate the existance of the equilibrium. The slopes of these curves indicate that the spectrum is constant in time after a certain time known as the thermalization time. This thernalization time will be defined and discussed subsequently. The energy groups having energy equal to, or lower than one KT, always increase in time till they reach the equilibrium state. The groups of energy higher than one KT increase first, reaching a maximum value before decaying to the asymptotic value. It is interesting to know the rate with which the energy groups are reaching their asymptotic distribution. The curves in Figures 6.3 and 6.4 provide excellent data for determining the thermalization time constant. In the fifth chapter we have determined X1 for the non-absorbing, infinite heavy gas model. If the contribution of the eigen values higher

-7550 1.0 KT 0.5 0.3 20 0.2 3.0 0 40. i n, 4.0 0: 5.0 0 30 60 90 120 150 180 210 TIME (ase) Figure 6.3 Time Behavior of Energy Groups in Beryllium for Infinite Medium and Zero Absorption.

50 E= 4.0 KT,,,- Io ->0 o 2.0 _o 10 5.0 50 110 170 230 290 350 410 470 TIME (\sec) Figure 6.4 Time Behavior of Energy Groups in Graphite for Infinite Medium and Zero Absorption.

-77than B3 is neglected then the distribution function will be represented as follows. _Alt b E,,t) = () iE)- (6.9) The first eigen value, X1, gives the thermalization time constant. We shall fit our data to the type of distribution given by Equation (6:9). The exponential fit of the distribution function is discussed in detail in a following section of this chapter. B. The Thermalization Time. The thermalization time may be defined as the total time required by the source neutrons to have slowed down and to have been completely thermalized, so as to attain the equilibrium distribution and to retain it for all times thereafter. It can be determined with the help of Figures 6.3 and 6.4 for beryllium and graphite, respectively. It may be noted that the precise determination of the thermalization time is not possible. A good estimate, however, based upon the slopes of the cuirves mentioned can be given. The uncertainty involved is between five to ten percent. For beryllium, the thermalization time can be taken as equal to 100 micro seconds plus 14 micro seconds, the mean slowing down time defined previously, making the total equal to 114 micro seconds, Kothari and Singwi(27) report the thermalization time in beryllium to be about 145 micro seconds, The time required by the neutrons to acquire the moderator temperature 310K from 2000'K, was calculated by these authors to be about 130 mirrcro seconds.

-78We estimate the thermalization time to be equal to about 240 micro seconds in the case of graphite, It takes about 200 micro seconds for the thermalization process and 38.6 micro seconds for the slowing (2) down process, Antonov et al. determined the thermalization time to be equal to 200 micro seconds, Beckurts(4) gave the estimate of the thermalization time constant in graphite equal to 180 t 40 micro seconds, If one were to take the thermalization time constant to be of the same order as the thermalization time, then the results of Beckurts agree with Antonov's results, Both values were obtained experimentally. The beryllium curves were obtained using twenty energy groups, with the help of the Oak Ridge National Laboratory analog computer. In the case of graphite only fourteen groups were employed, This problem was carried out on the analog computer of the Nuclear Engineering Department of the University of Michigan. It is interesting to note that the estimate for the thermalization times for beryllium and graphite given by the neutron distributions cp(Eat) are of the same order as given by the other authors, It may be suspected that the thermalization time is independent of the nature of the scattering kernel; and depends upon the mean velocities of the initial and final distributions. This statement is of a speculative nature, The agreement between the results based upon the heavy gas model and the results obtained for the crystalline medium, forces one to conclude the above statement. C. The Exponential Fit of the Distribution Function p(E{,t) We fit the distribution function cp(Et) as a sum of two exponentials'. We use the beryllium data in the interval of 60 miearQ seconds, 3tE5t) 44g)2-A=tt *,(cE) A (6A10)

Since, X0 is very small (around 0,l x 10 sec ), the above expression reduces to Equation (6.9). The fit was made using Cornell's method. (12) It involved the determination of cpo(E), cp1(E) and Xl, which correspond to the zeroth eigen function, first eigen function and first eigen value of the analytical method, We record these quantities for each energy group in Table 6.1. The arithmetic mean of the decay constant 1 is equal to 4.75 x 15. sec-1 in beryllium. When expressed in terms of k;US units, X1 is found to be 1.176 A e$Ui. The reciprocal X1, defined as the thermalization time constant, is 21.05 micro seconds in beryllium. For graphite tth is 56.8 L secs. On the heavy gas model, tth for beryllium given by Kotlhari and Singwi (27) is 21 1t secse The examination of cpl(E) indicates that it has maximum between 2 and 3 KT, at around 2.5 KT, and a minimum at 0o4KT, The zero is at about 1.3 KT, We plot cpo(E) and cp (E) in Figure 6.5 4 From the values of cp0(E); cpl(E) and ki listed in Table 6.1 we can construct neutron distribution function cp(Egt) for times greater than 60 [ sees by the following equation. Al(t-6o) ~(L-t~) = siboa(lt) - c (6., D)

TABLE 6,1 NUMERICAL VALUES FOR cp and.1 IN A HEAVY GAS MEDIUM OF BERYLLIUM) E(KT) cpo(E) c(E) X1 x 10-5 Sec0.1 10.13 -3.09 0.448 0.2 18.42 -4.79 o.446 0.3 25.14 -5.22 0.456 0,4 30.59 -5-39 0.472 0,5 34.29 -5.03 0.479 o0,6 37.10 -4,495 o.476 0O7 38.74 -3.86 0o490.08 40.18 -3.12 0.496. 9 40.74 -2.37 0,518 1,0 41.04 -1.62 0a503 1.5 37.89 1.635 0,439 2.0 29.07 3.122 0.461 3,o0.13.44 2.953 0O476 4,0 6.19 1.78 0.483 5.0 2*53 0.903 0.477 6.0o 0.957 0,438 0.477 Arithmetic mean of x1 equals 0.475 x 105 sec.-1

50 40/~~4 - 30 s 2 0.._... _, -e.-1001 -20 0 1 2 3 4 5 ENERGY (KT UNITS) Figure 6.5 Exponential Fit of the Infinite Medium Data for Beryllium.

8 _ _ _ _ _ _ _ _ -— _41...... ANALYTICAL 2 NUMERICAL (ANALOG COMPUTER) - -2 -4 X, = 1.I76,0so V (NUMERICAL) x, =1.2738 EZsovo (ANALYTICAL) -6~(E) =fEe-E o [-4.0t3 + -6 X NORMALIZED AT 0.4 KT 38)05 - 0.5685[- + 3.8605E - 0.5685fE + 0.0293E3] -8 0 4 2 3 4 5 ENERGY (KT UNITS) Figure 6.6 Comparison of Numerical and Analytical Eigenfunction (p1(E).

50.0 40.0 20.0 ~~~~~~~~~~20.01~~~ (ANALYTICAL) MAXWELLIAN I 0 *O. 10.0I oe- 10.0= 4NUMERICAL \! z 8.0 4.0 2.0 I 1. 0 1.0 2.0 3.0 4.0 E ENERGY (KT UNITS) -- Figure 6.7 Comparison of Numerical and Analytical cpo(E).

-84D. Comparison Between the Analytical and Numerical Results. The analytical results obtained in the fifth chapter for the heavy gas model, by the eigen value method, are compared with the numerical results determined in the previous section. 1,176 5.:S, r(. Numerical result fot- X' by.the ahalog computer (6.12) 1.2738 ~ So. Analytical result for X1 by the eigen value method The difference between these two results is about eight percent. The analytical value will decrease when more than four Laguerre polynomials are used, thus bringing the two values closer together. In Figure 6.6, we compare cpl(E) obtained by the two different methods. The agreement is very good on the low energy side, but is only fair on the high energy side, The analytical curve will converge to the numerical curve when a large number of Laguerre polynomials are used to describe the first eigen function..We normalize the analytical and numerical distributions', 1(E), at 0,4 KT, The numerical curve has a minimum at 0,4 KT, a maximum at about 2,5 KT and zero at 1.23 KT, whereas the analytical curve has a minimum at 0.4 KT, a maximum at about 2.5 KI and zero at 1...27 KT. In addition, we also compare (cp(E) given by the analytical and numerical methods, in Figure 6.7. The analytical method gives cpo(E) equal to a pure Maxwellian distribution. We use the values of co(E) as listed in Tables 5,3 and 6.1 for purposes of comparison,

III. The Infinite Medium 1/v Absorber It was pointed out in the fourth chapter that the asymptotic energy spectrum is Maxwellian only in the presence of a small l/v absorber, and not for a large absorption. We give the ratio of k0 and X1 for the various values of absorption, for the heavy gas model, in Table 6o,2 k0 and X1 are given according to the following equations. AO-= 2Iaoo = /- a LrS (62L3) Al-.(' 1 7 To (6.14) TABLE 6.2 THE VALUES OF - FOR DIFFERENT ABSORPTIONS X0 A XO 0.0' - 0.1 13.74 0,4 4.185 1,0 2.27 5.0 1.2548 10.0 1,12738 We conclude from Table 6.2., that B1 is of the same order as k0 for A larger than 1.0 The importance of this result has been discussed in Chapter IV.

-86The neutron distribution function cp(E,t) has been generated for various values of the absorption. It was not possible to obtain the distribution for A larger than 0,4. This was due to the rapid decay of the spectrum. In Figure 6.8 we plot the energy spectrum at 300 micro seconds, for the different values of absorption. The spectra were normalized for E = KT. All the cases have the same spectrum except for the case A = 0.4. The agreement is very good on the low energy side, but not on the high energy side. The neutron distribution for high energy is so small that the uncertainty of -the computer data becomes of the order of 15 to 25 percent of the actual value. The deviation on the high energy side may, therefore, be neglected for A equals 0.4. The decay constants calculated from the distribution function for the various values of A are given in Table 6.3. TABLE 6.3 THE DECAY CONSTANTS FOR DIFFERENT ABSORPTIONS- INI' T:HE''INFINITEiI'.MEDIUM A 1 (sec-1) X0 (sec-1) 0,08 2.850 x 103 3.56 x 104 0.10 3.572 x 103 3.57 x 104 0.20 7.185 x 103 3.59 x 104 0.40 14.275 x 103 3.57 x 104 The constancy of k0 indicates the decay rate to be proportional to the amount of absorption.

-8730.0 20.0 Z-ao A = SYMBOL 0 10.0 6.0 0.20 ~& 0.08 0.10 "' 20.2 6.0 2.0.4 ~4.0 0 z 0.2 0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 E-ENERGY (KT UNITS)- — & Figure 6.8 Neutron Energy Spectrum at 300 Ci sees in Infinite Medium of Graphite for Different Amounts of Absorptions.

CHEAPTER VII. TIHE. TME DEPENDENT EENERGY SPECTRUM IN FINITE MEDTUM The studdy of the'time dependent energy spectra in the finite media assimes special.significaance in the light of pulsed neutron experiments. In order to understand and interpret the results obtained from these experiments, it is essential to study the characteristics of the energy spectra existing in the medium at different times. The experimentalist measures the decaying spectrum (which is integrated over -all energies) as a fuxnction of tiime. We attempt to simulate the physics of the experiment by determining the energy spectrum as a function of time. With the help of the multigroup formalism introduced in the last chapter, we generate the distribution function cp(E,t). We deal with two cases. One is the constant diffusion coefficient, graphite case and the second is the enrergy dependent diffusion coefficient beryllitm case. We -ase the heavy gas fmodel for the'treatment of these cases. I. Constant Diffusion Coefficient, Graphite Case The assumption of the constant diffusion coefficient (D) for graphite is very good. It follows from the constancy of the scattering cross section (4.8 + 0.2 barns) in the energy interval between 0Q005 ev. and 1 Mev. (56) The Bragg cut off for graphite occurs at 0.0015 ev. We may further assume that the absorption cross section is very small, since it is about 0.5 + 0,2:millibarns, It may therefore be neglected. The energy spectra, as pointed out earlier, are based upon the determination of the energy transfer in the heavy gas model. The spatial dlistribultion of neutrons is assumed to be given by'the fundamental spatial -88

-89modes The higher spatial modes have been as-sumed to have small contributions, compared to the fundamental one, The time behavior of the neutron energy spectra for a few geometrical bucklings have been Obtained, In Figure 7*1, we give the spectra for the case where the geometrical buckling, B2, is equal to 18.5 x 10-3 cm The time behavior of the individual energy groups has been plotted in Figure 7,2. We find that the energy spectrum decays with a single decay constant after 250 to 300 micro seconds after the introduction of the source pulse. The existence of'the same decay -constant, within two to three percent, for each -energy group, indicates the presence of an equilibrium spectrum. The following characteristics of'the equilibrium spectra have been studied 1/2 in detail: the equilibrium spectrum; the -.average speed, E; and the decay.constant. ~. A, The Equilibrium Spectrum The normalized equilibrium spectra for -the following geometrical bucklings, B2 equ.l to zero, 2.96 x.10-3, 4,62 x 103, 8.21 x 10-3 and 185 x.10"3, are recorded in Table 7.1. For B2 = -and B2 = 185.x l03 cm-2 equilibrium'specmtra are plotted in Figure 7..3o We note from'these;curves that the energy spectrum shifts to the low energy side, with increase in the geometrical buckling, B2, The maximum of the spectrum for B2 equal to zero occurs at 1.0 KT, while for the case of B2 equal to 18,5 x 10`3 cm-2, the maximum appears'at 0.9 KT. We have determined the' 1/2 average speed, E ford ifferent.geametrical'bucklings.

-9o100 _ 50 01,. ", I / f l | | t = tS + 700 F sec!I 05 - — 1 t A I I I I I I 11= t$+4oo IsI 0.1 0.2 0.5 1 2 5 10 ENERGY (KT UNITS) (Er) Figure 7.1 Time Behatvior of Neutron Energy Spectrum in Finite Medium of Geometrical Buckling B2 = 18.5 x 1003cm-2 in Graphite.

50 20 3 40 \ = 0.9 KT 0 300 KT400 500 600 700 m TIME (tEsec) 2 0 4.0 KT 0 100 200 300 400 500 600 700 TIME (/Ftsec) Figure 7.2 Time Behavior of Each Energy Group in the Finite Medium of Graphite of Geometrical Buckling B2 = 18.5 x 10-3cm-2.

TAB2LE 7,1 EQUILIBRITUM ENERGY DISTRIBUTION (GRAPHITE) E.2 = O B2 = 2,.96 x 1053 B2 4,62 x 10-3 B2 = 8.21 x 103 10 (KT) ( cmn-2) (cm-2 ) (cm'2)' ( cm"2) 0,1 9.80 10,38 11,00 11.65 11.86 0.2 17.94 18.50 18.60 20.67 21.03 0.3 24, 63 29,39 26.00 27.26 27.63 o04 30.00 30.93 32.00 32.31 32.75 0.5 33.83 34.88 34.95 35.86 36.30 0.6 36.70 37-31 37.50 38.42 38.64 0,7 38.65 39.20 39.40 39-97 40.20 0.8 39.85 40,19 40.40 4o.80 40.85 o,9 40.60 40.69 40.90 41.04 41.11 1.0 40.91 491 4o91 40.91 40.91 40.91 1o5 36.68 36.13 35.97 34.65 34.83 2,0 29.38 28,70 28.52 27.62 26.66 3+0 14.49 14.16 13.99 13.58 12.57 4.0 6.35 6.26 6.17 5.96 5.16 5.0 2.62 2.57 2.56 2.49 2.08 6.o 1,00 1.00 1.00 o.98 0.76 7.0 0.35 0.35 0.35 0.30 0.26 8.o 0.10 0.10 0.10 0.10 0.10 9.0 0.03 0.03 0.03 0.03 0.03 10.0

-9345 40 35 2 8=0 30 2 -3 -2 I I8.5X 40 cm 25 (II 0 z 40 5 0 0 0.5 4.0 4.5 2.0 2.5 3.0 3.5 E, ENERGY (KT UNITS) Figure 7.3 Equilibrium Neutron Energy Spectra in Graphite for Geometrical Bucklings. B = 0 and B2 = 18.5 x 10'3cm2.

-94B. Average Speed, E1/2 The average speed in a medium having the equilibriumn energy spectrum cpo(E) can be given as follows. E Y. - - (7.1) E na(E) den $CE) dE (7.2) In Table 7.2 the average speed is given for ta few geometrical bucklings, TAKB~ 7.2 AVERAGE SPEED FOR VARIOUS VALUES OF B2 (GRAPHITE) B2 (cmj2) Average Speed* 0 1.o914 2.96 x 10-3 1.0820 4.62 x lO3 1.0789 8,21 x 10-3 1.0655 18.50 x 10-3 1.0513 * in units of vo = speed for most probable energy (KT). For the Maxwellian distribution, E1/2 can be given exactly: -E 1.129 (723) When compared to the case where B equals zero, in Table 7.2, we discover that the difference between the two values is about 3530 percent, which is due to numerical errors in the solution.

-95The average speed given in Table 7.2 characterizes the equilibrium spectrum. With an increase in geometrical buckling, the average speed decreases. The data given for E in Table 7.2 has been fitted by the least square method in the following manner. EY2 i's 163f5 _3*=2 (7.4) Beckurts(5)also obtained the expression for v on the heavy gas V0 model by the iteritive process assuming that.the equilibrium energy spectrum exists. His results are as follows: YE 2 v I 12- 9 t- I a= Beckurts' equatiorn5) (7.5) ___:e "(7.6) Transforming'Equation (7.4) in an analogous manner, we obtain: E" -_. Igog { 82j (7.7) C]/' -m i, s.2 (7.8) The expression for V/vo given in this study has been obtained after establishing the existence of an equilibrium energy spectrum during the last stage of the decay of the neutron pulse. C. Average Energy The average energy of'the neutrons can be defined as follows E Eo (E),~~ i7(tI

The values for the average energy are given in Table 7~3. TABLE 7.3 AVERAGE ENERGY B2 (cml) Average Energy* (E) 0 1,4176 2.96 x 10-3 1.3727 4,62 x 10o-3 1.3659 8,21 x 10-3 1.3410 18,50 x 10-3 1.2960 * in units of KT, --- B2 The least square fit for E ersus 13 gives E = IH 45 [1- 5.'G15 (7.S10) For a Maxwellian, E = o15 KT, so again there is an error, this time 5% in the B2 = 0 value. However, it is felt that the relative values may be accurately expressed by Table 7.3. D). -Neutron Teynperature' The concept of neutron temperature, though useful in characterizing the spectrum, is not rigorous, except in one case. This case is the one of infinite medium and zero absorption, where the final energy distr'ibution of the neutrons;is Maxwelliano For this case, the moderator temperature To may be taken equal to the neutron temperature Tn for the neutron gas. For the sake of convenience, we extend this concept to the cases of non-Maxwellian distributions, For the Maxwellian distribution,, E is equal to 3/2 KT- and is

-97also equal to 3/2 KTn for (B = 0). Extending this to the non-Maxwellian distribution, we get __________ T ve(f~~ O) ( l t) With the help of Table 7.3 we calculate the percentage change in the neutron temperature 5Tn/To for different values of B2. These values are listed in Table 7.4. TABLE 7.4 NEUTRON TEMPERATURE B2 (cm-2) (K?) 0 0 2.96 x 10-3 -3.17 % 4.62 x 10-3 -3.65 % 8,21 x 10-3 -5.40 % 18,50 x 10-3 -8.58 % E. Decay Constant We have calculated the decay constarts from the time dependent energy spectrum curves. After about 250 to 300 micro seconds, there exists an equilibriumwn spectrum for all cases of the geometri.cal bucklings. The neutron energy spectrum decays thereafter. In Table 7.5 we list the decay constants.

-98TABLE 7.5 DECAY CONSTANTS B2 (cm-2) (sec-1) 0 O 2.,05 x 10-3 4,84 x 102 2.96 x 10-3 7.08 x 102 3.65 x 10-3 9,00 x 102 4.62 x 10-3 10.16 x 102 6,0o4 x 10-3 13.83 x 102 8.21 x 10-3 17,94 x lo2 18.5 x 10-3 34,40 x 102 The uncertainty involved in the determination of the decay constants, for all cases except the last one, was within 3 percent. In the last case, the decay was very fast due to the small size of the assembly (about 40 cm3). The uncertainty in the determination of X was more than 5 percent. This was due to the time lag between starting the computer and the stop watch. The time scale was 10 micro seconds equals about 2.6 seconds. It must be pointed out that the error involved in the determination of X has to be minimized by obtaining the time measurements accurately* In the case of.small assemblies, the decay is so fast that the uncertainty in the value of X becomes very large. In Beckurtst experiments on graphite, the smallest size assembly was of B2 equal to 8.21 x 10-3 cm-2,

-99We have made a least square fit using the first six values from Table 7.5 in the range of bucklings < 8.21 x 10-3 cm 3. x = 2..29 x 105 B2 [1 - 2.068 B2] (7.12) II. Energy Dependent Diffusion Coefficient (Beryllium) Case. In the case of beryllium, Bhandari (7) calculated the transport mean free path as a function of energy. The diffusion coefficient for beryllium is taken to be energy dependent as given by Bhandari. The time dependent neutron energy spectra are generated for a few geometrical bucklings. In Figure 7.4, the time behavior of the neutron energy spectrum is given for B2 equal to 7.18 x 10-'2 cm-2, The behavior of individual energy groups is plotted in Figure 7*5. The parallel lines indicate the same decay constant for all energy groups. After about 100 micro seconds plus the slowing down time, there does exist an equilibrium spectrum. The neutronenergy spectrum maintains its shape but -decays with a single decay constant B equal to 8.39 x 103 sec.. We have studied the characteristics of the equilibrium spectrum for a few geometrical bucklings. A. Equilibrium Spectrum For the sake of comparison, we plot the equilibri umn spectra for B2 equal to zero and B equal to 7.18 x 10-2 cm,2 in Figure 7 6. The energy spectrum shifts toward the low energy side, The shift is small compared to the case of graphite, We calculate the average speed associated with the above spectra.

NEUTRON DISTRIBUTION, (E, f) 0 0~~~~~~~~~'' 01.0 I U N -a __, (D 23' N~~~~~~~ cD Pi P. CD 1~~~~~~~~~~0 w* 0-I(~ ~ ~ \1 I i" I u~ 0 V7 11I~ P1 ~~ ~ ~ ~ N CDO3 Z O,:,:,i +~ o ~\ + UI H _ _ _~c __ 0 r~ (F~5 0 M U)~~~~~~~~~~~~~~~~~~~~~~ (D 0 + 0 0 rY G) v, II ~ \ \ (D 0 o _ _ _+0 _ _ _ _ _ CD _ _ _ _" __ F~~ o'Q ooo~ 0,/0 o i0 fl)~~~~~~~~Cf r~'C~~~~~~~~~~~~~~~~~~~~~~~~ O - 1 —ol I -0Y )) ~ -Il.00 I tdo _, H d ~~~~~~~~~~iIC 0J'cd ~~~~~~~-Tj~~~ 3 — qC COO

-101100 80 60 40 20 6.2 0 1.0....8.6.4.2 0 50 100 150 200 250 300 350 TIME (EL Sec) Figure 7.5 Time Behavior of Each Energy Group in the Finite Medium of Beryllium of Geometrical Buckling B2 = 7.18 x 10O2cm-2.

-10242.5 40.0 \X 30.0 -.2= O 0 ~e tlkF —-g = n = 7.18 16 C 20.0 zI z 10.0 0 0 1.0 2.0 3.0 ENERGY, KT. UNIT Figure 7.6 Equilibrium Neutron Energy Spectr~ in Beryllium for Geometrical Buckling B = 0 and B = 7.18 x 10-2cm2.

-103B. Average Speed v/vo or 1/2 We have already defined the average speed by Equation (7.1), In Table 7.6 we record the values for the average speed. TABLE 7.6 B cm Average Speed 0 1,106 2.94 x 10-2 1.0941 5-36 x 10-2 1.0919 71.8 x 10-2 1.0586 The least square fit gives ZE / - 109 l 5 i -o l7 (7.13) The difference between the average speed for the Maxwellian distribution and for the case of B2 equal to zero, from the Equation (7.13) is about 1.7%. The coefficient of B- for graphite from the Equation (7f4) is -1.852 compared to -0.5179 for beryllium from Equation (7.11). We can, therefore, conclude that the shift in the energy spectrum toward the low energy side is larger for graphite than for beryllium. Co Neutron Temperature The neutron temperature can be estimated from Equation (7.13)* 1/2 The average speed E is proportional to the square root of the neutron temperature, if one assumes the Maxwellian distribution. With this

-104assumption, we have QI-) = 1.1099 [1 - 0.5179 B2] (7*14) where To is the moderator temperature. For the smallest assembly of B2 equal to 7.18 x 10-2cm used by de Saussure, the change in neutron temperature is about -7.27%, D. Decay Constant X The value of the decay constant for a few geometrical bucklings has been obtained, The.se values are given in Table 7.7TABLE 7,7 DECAY CONSTANT (BERYLLIUM) 2 (-2 ) B (cm 2) X (sec ) 2..94 x 10-2 3,64 x 103 3-96 x 10-2 4.53 x 103 5.36 x 10-2 6.86 x 103 7,18 x 102 8,39 x 103 The least square fit of X versus B gives the following expression for X. x = 1.267 x 105 [1 - 0890 B2 ] B (7.15) III, Diffusion Cooling Coefficient According to Equation (2.35) the decay of a pulse of neutrons in a non-absorbing finite medium, can be given by the following decay constant. - = pv B2 (7.16)

-105(Dv) in the above equation is averaged over the energy distribution of the decaying pulse in the medium. We expand Dv as follows: Dv = ( 1)o [1-CB2] (7.17) where (f))o is averaged over the Maxwellian distribution, C determines the deviation from the Maxwellian spectrum If the diffusion coefficient is constant, then (v)o is equal to ()o( V)o. We shall determine the diffusion cooling coefficient C using the following two equations for v acid X. x = (Dv)0 [1 - CB2J B2 (7.18) v = (V)o [ - CB2] (7*19) In the case of a constant diffusion coefficient) the two values of C must be equal. We shall discuss in brief the results for graphite and beryllium. A, Graphite (Constant Diffusion Coefficient) From two methods, we get the following values of C for graphite, from Equations (7.4) and (7.12). C = 1.852 from v (7.20) C = 2.068 from X These two values differ by almost 10 percent. We observe that for (v)o and (fv)o we get 1.088 and 2.29 x 105 cm2/sec. instead of 1.129 and 2.13 x 105 cm2/sec,, respectively. If we make corrections of these orders to the diffusion cooling coefficients given in Equation(7*20), we obtain C as equal to 1.924 (from X), and C equal to 1.922 (from v). We shall consider C equal to 1.922 for graphite.

-io6B: Beryllium (.Energy Dependent Diffusion Coeffient). For beryllium we get the following values of C from Equations (7?13) and (7eoJ5)9 C = 0o518 from 7 (7,21) C = 0.890 from X We shall compare the values of C,obtained by Hurwitz and Nelkin(21) and by Hafele and Dresner( to our value,. We use the following values for the constants~ So3e) = 8= 28 x 10-1 cm 2Sa Lj) = 4 l x 1o1 cm~ - (Be) = 0.222 Based -po.n heavy gas model S ( ) = 0*16667 Dr (le)= 1026 x 105/cm2 sec (-)S 0LC) = 2l13 x 105 cm2/sec We use the following formulae, obtained under the conastant diffusion coefficient assumptio~n. CHN =<1 ~DouL1s) (7Y22) H N; z,~eso (oH D2- 7 ) ~(7,23) Values obtained for graphite and beryllium are as follows: C (Be) = o046 XHN (7.24) CH.D(Be) = 0.459

-107CH (graphite) = 2.102 CH*-D(graphite) = 2.092 We conclude from the above study, that the di~tffusion cooling coefficient, as calculated in this study, differ.s by 9 percent from that of other authors, This observation is true for the -case of c:onstant diffusion coefficient, The value of the diffusion cooling coefficient obtained from X for beryllium is considerably higher than from o. This can be explained by the energy dependence of the transport mean free path, The experimental values of C for beryllium, according to de Saussure (15) (4) asd Silver(5) and for graphite according to Beckurts are as follows: C (Be) = 1.1 + 0.9 (7.26) C (graphite) = 7.0 + 1 The experimental values for graphite *are considerably higher than the heavy gas model results;. But for beryllium the calculat;ed value for the energy dependent D(E) is closer to the experimental value than the constant.diffusion coefficient value.

CHA\PTER VIII. THE TIME DEPENDENT DISTRIBUTION IN THE PULSED MULTIPLYING MEDIUM In this chapter, we discuss the time behavior of the neutron energy spectrum, in the multiplying medium. In addition to the decrease of neutron density, there is a regeneration of neutrons due to the presence of fissionable material. The process of regeneration is characteristic only of the multiplying medium. I. General Formulation, The distribution of the neutrons according to the diffusion theory approximation can be given as follows, where cp(E,r,t) is the thermal flux. The energy E lies in thermal region. r ~VE,y(t) - - { 1E ) + KJE) - T vr- $(E) C t) t r zf, E, Y t t (8.1) c(Eik(E} tY) F(Ei-E) cE' + S,, (E, rt) The external source term Sext(E,r,t) is represented by the delta function in time. S tCE,(F r,t) 5 2AE, ) 8(t) Equation (8.1) contains the extra source of neutrons due to the presence of the fissionable material. We shall assume that the fission of the neutrons takes place at energy E. The slowing down of the neutrons -108

-109born with fission energy Eo is given by the slowing down kernel q(r'-~r; t't-t; Eo-~E). We have neglected the effect of the delayed neutrons, which is about 0.75%. In order to include this effect, the extra source term should be multiplied by a factor of (1-P), where 3 equals 0.0075 delayed neutron fraction. We define the kernel as follows. The slowing down kernel gives the probability that the fission neutrons of energy Eo, born at r' and at 4A will slow down to energy dE about E at dr about r and at a time dt at t. Let us make the following ansatz.,v ~r) = - 3Bp ( K}) 3 ~p(K) -o (8.2) cpp(r) vanishes at a fixed boundary r. This boundary is given by the actual geometrical dimension, plus the extrapolated distance equal to 0.71 Xtr' The symbol Xtr stands for the transport mean free path. 2 Bp in Equation (8.2) stands for the geometrical buckling for the p-th spatial mode. We have assumed in Equation (8.2) that the neutron distribution goes to zero at the fixed boundary r for all energies. That means that the extrapolated distance is independent.of energy. Substituting the set (8.2) into Equation (8.1) we get E_ ~ Sg21t a L ~3p(r)_()X (E) +Z E(E) tI DRE) b (E + KF'2 (E') (Et)F(EE) d (8-3) FE o

-110We shall use the results of the fundamental theorems of the asymptotic reactor theory to simplify the slowing down kernel. These theorems first given by Weinberg(5l) are discussed in detail by Weinberg and Wigner.(52) We give these theorems in brief here. First theorem: This theorem has two postulates. 1. The finite medium slowing down kernel can be replaced by a displacement kernel. The slowing down kernel is a displacement kernel in the infinite medium. In the finite medium, away from the boundary, the kernel can be taken to be the displacement kernel. 2. The effect of the finite medium is taken into account by making the distribution zero at a fixed boundary r. Mathematically this is expressed as 5) R(r. r v dr = ~)4P~" d Second theorem: This theorem states that the non escape probability during moderation in a uniform bare assembly is given by the Fourier transform of the slowing down kernel. In essence the second theorem replaces the infinite medium integral by the Fourier transform. - On P

-111iUsing the results of both theorems, the slowing down integral can be represented as follows. P t' r/ (8.6) (r ) ) C( E, ) ( te r p t-t E o E) E t = P i pt' where tJ(it13; t-toj E(+E) is equal to the Fourier transform of the slowing down kernel in the infinite medium with Bp as the Fourier transform variable. Substituting Equation (8.6) into Equation (8.3) for the p-th spatial mode, we get 2t~ i +~T2IE() ) J pCEQt~cj4 (4; &.t;Et' EoJ4E t_ -(-8.7) + J ES (Eit i +( (Et) 4~0~st dF, We make the following expansion for the time dependent function. ~ (p ( E t) L E ~9 At) U At t (8.8) p ( o, t= O0( t O Xi is the eigen value for the i-th energy eigen function, for the p-th spatial mode. Our assembly is subcritical, therefore, the time behavior is given by the decaying exponential.

-112Substituting the expansion for cp (E,t) into (8.7) E e- F e t E-r + zE.), e t (&-9) + E I = T-.f p~ jeE ) a ) a We simplify the slowing down integral. I = 9vea C Let t-t' equal to. We substitute this into the integral (8.10). We can identify the above integral as the Laplace-Fourier transform of the slowing down kernel in the infinite medium with s= -Xpi being the Laplace transform variable. — 7 ~( - Sto + % 9.>(3jA~+)()+^SFEE (C~oC) t3~iE We substitute (8.11) into Equation (8.9). Then we get for the p-th spatial mode and the i-th energy mode: - b(E,'),'EE We have removed the subscripts p and i for the sake of brevity.

-113II. The Monoenergetic Case We specialize Equation (8.12) for the monoenergetic thermal case. We assume that the experimentalist is measuring the density of monoenergetic neutrons. Then the decay constant X is given by the following expression. h ( C Sa +; 1) Lr _ SiT~r 2 i( S,-A, UTEo- (8.13) In order to use the above equation, we must determine the Laplace-Fourier transform of the slowing down kernel. The determination of this transform requires a detailed knowledge of the slowing down process. We shall treat two cases, namely, the heavy element moderated and hydrogen moderated systems. In the case of heavy element moderation, we can use the continuous slowing down model, i.e., the Fermi age model. For hydrogen moderated systems, this model is not applicable, because the neutron can lose all its energy in a single collision with the hydrogen nucleus. We require a different model for this case. A. Heavy Element Moderation. (A >~ 1) On the basis of the Fermi age slowing down model, the LaplaceFourier transform is here calculated. We rewrite Equation (8.10). 00 - A; ~~-~'E2 9 Ato 0 The Fourier transform of the slowing down kernel, in the infinite medium, can be given in a simple manner as follows: 5m ( 5 Lo {o E, ) c 4 _ B

-114The slowing down kernel according to the Fermi age model is as follows. 00 (C EE) - 77 - (8.15) where T equals the Fermi age corresponding to the energy interval from Eo to E, and p equals the resonance escape probability. Substituting Equation (8.15) into Equation (8.14) we get - 9 - g2 q, i -, Eo E) p e8.163 When the neutrons are slowing down, according to the Fermi age model, the time taken by neutrons of speed vo to have slowed down to speed v is as follows. ts ( L) [= Ur5 & Lo] It is characteristic of the Fermi age model that all neutrons having a speed vo at time zero, have the same speed v at the time ts. Mathematically, it can be expressed as follows. qR (b 3v to) *E)= e. p Cto-ts) (8.17) where b(to-ts) is the Dirac delta function. Using Equation (8.15) in Equation (8.10), we get 2-~~~ 9 -B'Z Ato eoBo (B3, - j uj E) -s Ato-t5) dto - p~ BiAt-5 -8.18)

i ll5 The Laplace-Fourier transform of the slowing down kernel in the infinite medium, is given by Equation (8.19). We substitute this into the expression for the decay constant X in Equation (8.13) -A = W t:G+ D 529 C)) - 2-z-a A ts (8.19) If we take into account the delayed neutrons and expres-% VZf in the conventional form, Equation (8.19) reduces to A = ~ ( f+ D BQ ) - Y ea (lF-) I e (8.19) U XJ where u/ Equation (8.19) should be used to analyse pulse neutron data in a multiplying medium containing the heavy element moderator, i.e., graphite or beryllium. It must be pointed out that it is not possible Xts to omit the exponential factor e, unless Xts << 1. This is not true in the case of small assemblies, which have large decay constants. B. Hydrogen Moderation (A = 1) As pointed out earlier, the hydrogen moderated systems can not be treated on the basis of the Fermi age model. The study of slowing down neutrons in these systems is still actively pursued. A great deal of effort is being spent to study this problem theoretically and experimentally. It is well known that the age of the fission neutrons to the indium resonance energy (1.44 e.v.) in water is 26 cm2, as determined by theoretical methods.(22) The earlier experimental values are higher, as obtained by Hill et a1o(20) (30 cm2) and by Barkov and Mukhin(3)

(29.4 + 1.5 cm2). Recent experiments at Argonne National Laboratory give a lower value of the age, about 27 cm2. It must be pointed out, however, that the Selengut-Gortzel model gives a theoretical value of 30 cm2, in agreement with the earlier experimental values. A discussion of the various models has been given by Hurwitz and Zweifel.(22) Due to lack of agreement about the use of a particular model, we represent the Fourier transform of the slowing down kernel in an arbitrary manner for hydrogeneous systems.., (E -; -E) --: ) LB r K](8.20o) We define the 2K-th spatial moment as follows. r _ _ C )., (8.21) $00:( (, E:,o E) rer The spatial moments are given by the type of model used to describe the slowing down kernel. In the case of water and other hydrogeneous systems Equation (8.20) should be employed. A choice must be made between an experimental or a theoretical model. For analyzing pulse neutron experiments in a multiplying medium moderated by hydrogeneous material, we need the Laplace-Fourier transform of the slowing down kernel. We have given the solution for this case in Chapter III. The results obtained were given in Equation ( 3.26) and are as follows. Emo(~,s~ly/p'l ) D EI - r_ EQ(E - 5-Er (8.22) -h)!'

-117We write the above equation as the product of three factors. + R fs(,-&(E) _I E' E[oz E A ~ U'(DE)~s X:,:, h =E A =-t )] We simplify the above equation. 2 z -62 Ats CCCE Bj-A) - p eE (8.24) Equation (8.24) is true only if it is possible to separate the tim and spae and space dependet soltions We assume that in the exponentials B2 involving resonance escape probability p and non-escape probability e, T the decay constant X is small compared to zaO + v(DB2 + Zs). This leads to a separation of Equation (8.23) into two parts. 4Eo EE _ qo E' ~~4~;t Zao)E t I -o +-)L3 = L_)'1/ c~',t- >,) -.e -- )( Eo iZ5c~oi+ g(~p)B"-t 8-)5)

Equation (8.25) can be expressed as the product of two functions, as follows. The function P(B2;- k) may be regarded as a function of X alone, since the leakage term vDB2 may be neglected in comparison to a~ + vZ - X. This leads to a simplification of Equation (8.26). %XtE; By->) = (B) p (-A) M8.27) / / Z T_ (..)=.. "; O ~v,.+'_ ) ].',.., (8.28)., +::: + +Es)]E assuming a +s ) A A e. 4 M__ (8.29) assuming oa, ZS T -' A >> L 2) I; E () + s )+(Eo)>) Do B The decay constant in pulse neutron assemblies consisting of hydrogenous and fissionable material can be given by the following equation.:; X am a L-r B- L T1 -aU,(R)PeA)(I-A.83o)

It must be pointed out that Equation (8.28) corresponds to the expression for the Selengut-Gortzel slowing down model in hydrogen. It is essentially the same expression as derived by Simon(42) for the case of constant cross section using the Selengut-Gortzel model for the timeindependent infinite medium problem. In the analysis of pulse neutron experiments on subcritical systems, Campbell and Stelson(ll) used the following equation for the decay constant X. A, -_ Lr 2: - T D ) (8.31) When the expressions for X given by Equations (8.30) and (8.31) are compared, we find that the factor P(- X) is missing in Equation (8.31). The effect of this factor is to decrease p(B2). Campbell and Stelson(ll) expressed P(B2) as equal to 1/1 + B2T and obtained T as equal to 21.2 cm2. If, however, we take into account the correction factor P(- X), then the value of T increases. We shall determine this correction subsequently. III. The Multi-Energy Case We shall use Equations (8.9) and (8.11) to study this case. Combining these two equations, we get Z -AI ('t A I #LjE) (8.32) ~+ J s(E) jE') F(E E) E' E'

-120A detailed knowledge of the variation of the kernel with the final energy E is required, in the multi-energy case. This variation depends upon the physical model used to describe the slowing down phenomenon. Even if the Fermi age slowing down model is used, the problem does not become a simple one unless a few drastic assumptions are made. In the final stage of the slowing down process, we assume that the slowing down time can be an average, ts. In the energy interval between 10 KT and 1 KT, the slowing down time ts(E) varies by almost a factor of five. The slowing down length, however, does not change more than a few percent in the above interval, thus it can justifiably be given by an average value. We shall use the above assumptions for simplifying the LaplaceFourier transform of the slowing down kernel. Let us make the following ansatz: ( E) - T_ an L (E) iM(E) % (o0) -o (&) - o833) For the scattering integral, we shall use the expression given previously, - +(E)+ Z s(E/ tE ) F(E )E = 1 M(g) AkE E' We substitute ansatz (8.33) into Equation (8.32), then multiply by Lm( (E) and finally integrate over all energies. Using the notation of the fourth chapter, we have + Fd +.tG mn,

-121Qo N 2 (I) G \,} = j 2r sf $x (B, A +,^- ) L ~46(M L (E)c4E(8.35) We shall simplify the GmS integral for a special case. The secular determinant for (8.34) is as follows: (e- 2 ) mn v ~ 2=),, -F~- =G -o (8.36) The solution of the secular determinant (8.36) will give the values for X. It will, however, require a detailed knowledge of Fn and Gmn integrals. For the Fermi age slowing down model, we have::oo A4(EZ _B ) z e L (E) (837) As discussed in detail above, we assume: 1. The neutron age TEoE does not change appreciably in the energy region (ET > E > O). ET equals thermal cut off energy. 2. The slowing down time ts can be given by an average value ts weighted over the Maxwellian distribution. J t(E) EeX (8.8) 08-38) 0 -leQ-d

On substituting these values into Equation (8.37), we get o tv crs setio l, M(m ) Lti (8.39a For 1/v cross section law, Gmn integrals can be evaluated. IV. Analysis of Experiments We shall analyze the data given by the experiments of Campbell and Stelson.(ll) The value of the decay constant X for different geometrical bucklings B2 at a concentration of 26.5 grams U235 per liter are given in Table 8.1. TABLE 8 1 2 (lo) THE VALUE OF THE DECAY CONSTANT FOR DIFFERENT B3 VALUES B2 cm2 10-4/sec0o0436 0.550 0~0483 0.630 0.0562 0.704 0.072 0.883 0.135 1.386 0.192 1.768 0.296 2.358 0.506 3.194 We calculate the experimental non-leakage probability P (B2) exp using the following equation. - u-a + a X: d -,- u- (, - "l:'-) (8.40)

-123We use the following values for the constants. VZa = 0.48 x 105 sec-1 for water VZU = 0.965 x 105 sec-l for U235 (26.5 gms 3235/liter) a r =- 2.07 4 2 -1 vD = 3.49 x 10 cm sec - 29.48 x 104 sec'l v = 2.2 x 105 cms/sec Thermal neutron speed We express % P (B2) as the product of two f-unctions, T(B2). P(-X) as given by Equation (8.27). P(-X) is calculated according to Equation (8.29), In Table 8.2 we give Pexp(B2); P(B2) and P(-X) for different values of the geometrical buckling, B2. TABLE.8 2 NON -ESCAPE PROBABILITIES B2 xp (B2) a) ( 0..o436 0,524 10o3697 0.5053 0,0483 0.490 1.0o4285 o0 4699 0,0562 o. 469 1,0480 0,4475 0,0o720 0o.407 10 o646 0.o3823 0.135 o,.265 1,0986 002412 0.192 0.174 1.1156 0o1560 0,296 0.060 1o.1676 0.0514 0.506 o0.009 1.2318 0.0073 In the region where B2 is less than 0.o14 cm-2, Campbell and Stelson fitted Pexp(B2) by 1/1 + B2r and obtained s as equal to 2L12 c2. (11)

-124If we take P(B2) instead of P xp(B2) we get equal to 22.82 cm2, It is exp important to note that the P(-X) correction increases r in the right direction. In small finite systems, the experimental value of F(B2), rather than the neutron -age, should be compared with the theoretical calculations. For.small systems it is no longer possible t.o express (B2) as equal to 1/1 + B2; Equation (8.20) should be used. In the case of heavy element moderators, P(-X) is equal to eXtS. The experimental value of Pexp (B2) can be used to obtain the slowing down time ts. Kloverstrom and Komoto (26) carried out the experiments on graphite enriched uranium systems and analyzed their data taking into account the eXt s factor.

APIW1ER IX. DISCUSSION OF RESULTS. We have studied the time behavior of the low energy spectrum, for- a pulse of fast neutrons, in a homogeneous medium, Analytical and numerical methods have been employed to undertake this study. Characteristics of neutron thermalization and diffusion of a pulse of neutrons in a medium have been determined. The analytical method is based upon the expansion of the neutron distribution function (p(E,t) for one spatial mode, in the following manner. In the above expansions, iP(E) is the eigen function associated with the i th eigen value Xi. Each eigen function is expanded by a sum of an infinite number of associated Laguerre polynomials of first order, weighted by the Maxwellian energy distribution. The choice of the i>guerre polynomial is arbitrary. I. First Eigen Value We have given a general formulation of the method for determining eigen values which characterize the decay of a pulse of neutrons in any medium. The zeroth eigen value was determined by Singwi.(43) The expression for the first eigen value, according to Equation (4,20), is as follows: A4 ( Wl woo - tol) (2 ~ D +8z u~o 2o q -F0 lo I J1' -1 Wl25 Woo -M wlAI -125

-126In a non-absorbing, infinite medium, the zeroth eigen value is zero, but the first eigen value is not zero and is equal to the following expression. F 14Woo Lr 1 4=q ( >~l oo - woi). = -'i78 XIo9 Mi ( 5ec-i) X1 determines the rate of thermalization. The reciprocal of B1 is the thermalization time constant, with which the neutron energy distribution approaches the Maxwellian distribution, We have derived the following relation between the thermalization time constant and the diffusion cooling coefficient, as given by Equation (4.25). tt- Co'T (9.4) Nelkin(32) derived this equation based upon the concept of neutron temperature, using the variational method, We have rigorously derived this expression without using the concept of neutron temperature, which has been criticized on physical grounds. The thermalization time constants for a number of moderators have been determined and are found to be in good agreement with the results of other authors. II, Determination of M2 (32) It was pointed out by Nelkin that the experimental value of the diffusion cooling coefficient can be used to determine the thermalization parameter M2. This parameter is related to the second

-127energy transfer moment, as shown by the following equation. j = —- JJ EZS(E) M(E) F E E) ( -E) EdE (9.5) Finite medium corrections, coupled with the dependence of the diffusion coefficient upon energy, limit the usefulness of this method. We have pointed out the importance of the determination of M2 by the absorbtion method. The eigen values for an infinite medium, containing a 1/v absorber, are given as follows. /\o = 1ao SoT Ad I W 0 0 t9.6) For large absorptions is comparable to k0 and the f in al decay of the For large absorptions, Xi is comparable to X0 and the final decay of the energy spectrum is governed by both, not by X0 alone, which is the case for small absorptions. The advantage of the absorption method lies in the determination of M2 by direct measurement of decay constants, which can be measured accurately. III. Thermalization in a Heavy Gas Medium For the heavy gas, we have determined k1 and X2 for a nonabsorbing and infinite medium case. The fact that the eigen functions of the heavy gas differential operator are the associated Laguerre polynomials of first order, simplifies the problem. The two eigen values for this case are as follows: Z1 - 1 2738 EsoL) Four Laguerre polynomials (9?7) h 4' 0 Y~ 8 ~ s-.pro (0Three Laguerre polynomials

-128The eigen function corresponding to the first eigen value was also obtained. The expression for cp1(E) is as follows. () - E C aol C- L oi~ t3,86l -o'569 a (9.8) + o 29l E 33 e e IV. The Time Dependent Energy Spectrum in the Infinite Medium A unique method of studying the time behavior of the energy spectrum of a pulse of neutrons has been developed. With its help, we have generated the neutron energy spectra for times greater than the slowing down time. The source at the slowing down time was calculated on the basis of the Fermi age slowing down model, It must, however, be pointed out that the results obtained here are independent of the details of the initial source at the slowing down time, In the infinite medium, we have studied two cases: zero absorption and 1/v absorption. For the case of zero absorption, we determined the total thermalization time required by the neutrons to reach the Maxwellian distribution, The thermalization time for beryllium and graphite was found to be 114 micro seconds and 240 micro seconds, respectively. Kothari and Singwi's calculations for beryllium gave 143 micro seconds, and Antonov et alts. experiments gave for graphite 200 micro seconds. It is speculated that the total thermalization time depends upon the mean speeds of the initial and final distributions, and is insensitive to the scattering kernel, The time constant with which the neutron energy

-129C distribution approaches the Maxwellian distribution is found to be equal to 1.176 se U~$o compared to 1.2738 f 5o 9 as given by the eigen value method. By means of an exponential fit, we have determined cp1(E) and have compared it with the first eigen function cpl(E). The agreement is fairly good. V. The Time Dependent Energy Spectrum in a Finite Medium The finiteness of the medium introduces the effects of leakage of neutrons, In all the previous theoretical studies, it was assumed that the asymptotic energy spectrum can be described by a single decay constant. In this thesis, however, it is established that the asymptotic energy spectrum decays with a single decay constant. The energy spectrum attains an equilibrium spectrum at about the thermalization time, The shape of the spectrum remains constant in time, but the amplitude deQreasfeso A shift of the equilibrium spectrum toward low energy was observed with the increase in geometrical buckling. This is an evidence of thetcooling of the spectrum', The following characteristics of the equilibrium specta were studied, 1. Average Speed 2 Average Energy E and I Neutron Temperaturet 3. Decay Constant 4. Diffusion Cooling Coefficient Two special cases were studied: firsts the constant diffusion coefficient case of graphite, and secondly, the energy dependent diffusion coefficient case of beryllium. The average speed v and the decay constant X, when fitted by the least square fit give the'diffusion cooling coefficient T.

-130The value for the diffusion cooling coefficient for graphite obtained here was 1.922. This differs by about 9 percent from the results of other authors, namely a value of 2.10 for the diffusion cooling coefficient. For beryllium, we obtained the following values of the diffusion cooling coefficient, O. 52 from v 0.89 from X Since we used the energy dependent diffusion coefficient, the two values are not expected to be identical. The results of other authors, based upon a constant diffusion coefficient, are equal to 0.46. The Multiplying Medium We have analyzed the data of Campbell and Stelson's experiments on U235 and water mixtures. The Laplace- Fourier transform should be used to analyze the pulse neutron experiments instead of the Fourier transform of the slowing down kernel. Taking into account this correction, 2 2 the age is found to be 22.82 sq cm2 instead of 21.2 cm as calculated by Campbell and Stelson. (ll) The pulse neutron experiments on multiplying media can be used to determine a large number of reactor physics parameteres: the slowing down time, the slowing down length, and r. However, the determination of these quantities requires a careful analysis of the experiments, as pointed out above.

-131VI. Conclusion The eigen value method is useful in determining the upper limits for lower eigen values, which characterize the last stage of neutron thermalization and diffusion processes. This method fails to give the neutron distribution function at all times. In order to obtain results of higher order, a large number of Laguerre polynomials and the higher energy transfer moments are required. However, the usefulness of this method lies in its exactness and in its general applicability. For the heavy gas model, this method can give exact results of all orders, though the calculations involved are tedious and require large machine calculations. In addition to the eigen value method, a numerical method to generate the behavior of the neutron energy spectrum at times greater than the slowing down time, has been developed'This method has been applied to the heavy gas model, and the' time dependent energy spectra were obtained for the various cases mentioned earlier. This method is unique, since it simulates the experiment and follows a pulse of neutrons during thermalization and diffusion periods. The advantage of this method is that it describes the neutron distribution function cp(E,t) completely, provided the initial source distribution is given exactly. It also gives the correct asymptotic behavior, which describes the last stage of neutron thermalization and diffusions and is independent of the details of the initial source. Results obtained by -this method differ by 5 to 10 percent when compared to the analytical, results. The limitations of this method are determined by the machine used to generate these distributions. It might be possible to generate these distributions for the crystalline model by a similar multigroup technique,

-132It should be pointed out that the two methods described are complimentary to one another. We have established in this study the existence of an asymptotic energy spectrumm characterized by a single decay constant X0. Evidence of the'diffusion cooling', based uoaq the time dependent energy spectra for the finite media, has been given. The energy spectrum deviates from the Maxwellian distribution and shifts toward the low energy side in the finite media, In the case of an infinite medium, for small absorbtions, the energy spectrum is Maxwellian and is characterized by a single decay constant XkO However, for large 1/v absorption, the energy spectrum is not Maxwellian, and its decay is governed by two decay constants, one of which depends upon the energy transfer properties of the medium. It should also be pointed out here, that for small absorptions, the pulse neutron experiments give the absorption cross section and diffusion coefficient averaged over the Maxwellian distribution. The deviation from the Maxwellian distribution is contained in the'diffusion cooling' term. It is hoped that this study would further the understanding of neutron thermalization and diffusion, and would stimulate and aid in the study of therrmalization and diffusion by pulse neutron techniques.

REFERENCES 1. Amaldi, E. Handbuch der Physik. 38, 1 (t959). 2. Antonov, A.V., Isakoff, A.I., Murin, IoD., Neupocoyev, B.A., Frank, I.M., Shapiro, F.L. and Shtranich, I.V. Proc. Int. Conf. on Peaceful Uses of Atomic Energy, Geneva. 5, 3, 82 (195 5). 3. Barkov, L.M. and Mukhin, K.N. Jour. Nuc. Energy, 4, 91 (1957). 4. Beckurts, K.H. Nucl. Sci. Engng, 2, 516 (1957). 5. Beckurts, K.H. Zeitschr. f. Naturforsch, 12a, 956 (1957). 6. Bengston, J., et al. UCRL-5159 (1958). 7. Bhandari, R.C. J. Nuclear Energy, 6, 104 (1958). 8. Bhandari, R.C., Kothari, L.S. and Singwi, K.S. J. of Nuclear Energy, 7, 45 (1958). 9. Brown, H.D. and St. John, D.S. du Pont Report DP-33 (1954). 10. Campbell, E.C. Private Communication. 11. Campbell, E.C. and Stelson, P.H. ORNL-2204 (1956). 12. Cornell, R.G. ORNL-2120 (1958). 13. Cohen, E. R. Proc. lst Intern. Conf on Peaceful Uses of Atomic Energy Genera (1955) p/611. 14. Davison, B. Neutron Transport Theory. London: Oxford Univ. Press (1956). 15. de Saussure, G. and Silver, E.G. ORNL-2609, 59 (1958). 16. Eriksson, K.E. Arkiv for Fysik, 16, 1 (1959). 17. Fermi, E. Ricera Scientifica, 2, 13 (1936); English translation: NP-2385 (1951). 18. Fultz, S.C. UCRL-4874 (1958). 19. IHffele, W. and Dresner, L. ORNL-2842, 124 (1959). 20. Hill, J.E., Roberts, L., and Fitch, T.E. Jour. Appl. Phys., 26, 1018 (1955). 21. Hurwitz, H., Jr. and Nelkin, M.S. Nuclear Sci. Engng., 3 (Jan. 1958). -133

-13422. Hurwitz, H., Jr. and Zweifel, P.F. Jour. Appl. Phys., 26, 923 (1955). 23. Kazarnovsky, M.V. Jr. Nucl. Energy II, 9, 293 (1959). 24. Kazarnovsky, M.V., Stepanov, A.V., and Shapiro, F.L. Proc. 2nd International Conf. on Peaceful Uses of Atomic Energy, 2148 (1958). 25. Kirschbaum, A.J. UCRL-4983-T (1957). 26. Kloverstrom, F.A. and Komoto, T.T. UCRL-5477 (1959). 27. Kothari, L.S. and Singwi, K.S. Solid State Physics, 8, (1959) New York: Academic Press. 28. Krieger, T.J. and Nelkin, M.S. KAPL-1597 (1956). 29. Krieger, T.J. and Zweifel, P.F. Nuc. Sc. and Engng., 1 (Jan. 1959). 30. Margenau, H. and Murphy, G.M. Mathematics for Physics and Chemistry, (1951) New York: D. Van Nostrand Company, Inc. 31. Marshak, R.E. Rev. Mod. Physics, 19, 218 (1947). 32. Nelkin, M.S. Journal of Nuclear Energy, 8, 48 (1958). 33. Nelkin, M.S. GA-746 (1959). 34. Nelkin, M.S. and Cohen, E.R. Proc. 2nd Int. Conf. on Peaceful Uses of Atomic Energy, Geneva. P 1839 (1958). 35. Nelkin, M.S. GA-816 (1959). 36. Ornstein, L.S. and Uhlenbeck, G.E. Physica, 4, 478 (1937). 37. Osborne, R.K. Nuc. Sc. and Engng., 3, (Jan. 1958). 38. Purohit, S.N. ORNL-2842 (1959). 39. Purohit, S.N. and Zweifel, P.F. Trans. of the Am. Nuc. Soc., 2, No. 2, 100 (1959). 40. Reynolds, H.L. UCRL-5175 (1958). 41. Sachs, R.G. and Teller, E. Phys. Rev., 60, 18 (1941). 42. Simon, A. ORNL-2098 (1956). 43. Singwi, K.S. Submitted to Arkiv for Fysik for publication. 44. Singwi, K.S. and Kothari, L.S. J. of Nuclear Energy, 8, 59 (1958).

-13545. Sjostrand, N.G. Arkiv for Fysik, 15, 12, 147 (1959). 46. Svartholm, N. and Tek, Chalmers. Hogsk Handl., 164 (1955). 47. Vaughan, E.V. and Cohen, E.R. Proc. Neut. Therm. Conf., ORNL-2739 (1958). 48. von Dardel, G.F. Trans. Roy. Inst. Technol., Stockholm, 75 (1954). 49. von Dardel, G.F. and Sjostrand, N.G. Progress in Nuclear Energy, Series I, Physics and Math., 2, Pergamon Press (1958). 50. Waller, I. Proc. 2nd Intern, Conf. on Peaceful Uses of Atomic Energy, Geneva, 153 (1958) 51. Weinberg, A.M. Am. Journal of Physics, 20, 401 (1952). 52. Weinberg, A.M. and Wigner, E.P. The Physical Theory of Neutron Chain Reactors, Chicago: The University of Chicago Press, (1958). 53. Wigner, E.P. and Wilkins, J.E., Jr. AECD-2275 (1944). 54. Wilkins, J.E., Jr. CP-2481 (1944) o 55. Zemach, A.C. and Glauber, R.J. Phy. Rev., 101, 118 (1956). 56. Neutron Cross Section Tables, BNL-325 (1955).

APPENDIX Io PROPERTIES OF ASSOCIATED LAGUERRE POLYNOMIALS We list the properties of the associated Laguerre polynomial of (X-th order, as given by Copson(l) and Sneddon. (2) Ln(Z) _ (+l+n) lF1 (-n; d+l; Z) (I.1) (-n) Zr Hypergeometric 1F1 (-n, I+, Z) ( ) r=O r! (Q+l)r Function 7 a+r' ()r +r (1.3) Zt % tn Ln (Z) (l-t) e () e-Z Z- dn (eZ Zn+) (I.5) Ln (Z) n! dZn (e If (a-n) -— is positive integer or zero, then (-1) eZ dn+a -Z n Ln (Z) - (n) e dn e-Z (I.6) n. dZn+ae f ZeZ L (z) Ln(Z) dZ = 5 bn m (I.7) 0 nO -136

-137If f(Z) = Z Cn Ln (Z) and term by term integration is valid, 0 then o GCn = eZ ZC f(Z) LC (Z) dZ (I.8) Ln n! r=O r! (2)r (I.9) L(l) (z) n (_1)n n! Zr n =(n+l) r0 (n-r)! (r!)(r+l)! (1.10) 00 (1) (Z) =n L) () ) (Z) Ze dZ nEl mn (I.11) The normalizing factor is T he normalized associated ()n (1)n Zr Y (1) (Z)= _( n+l ) (1)nZr n+ r=O (n-r)' (r!)(r+l) (I.12) n Ln (Z) is the solution of the following differential equaZ Y + (2-Z) dy+ ny = 0 (1.13) dZ2 dz The differential operator of the heavy gas model has the associated Laguerre polynomial of first order as its eigen function.

-138REFERENCES 1. Copson, E.T. Functions of a Complex Variable. Oxford University Press (1935). 2. Sneddon, I.N. Special Functions of Mathematical Physics and Chemistry. New York: Inter-Science Publishers, Inc. (1956).

APPENDIX II. TIME DEPENDENT ENERGY SPECTRA FOR DIFFERENT B2 We give the additional time dependent energy spectra for beryllium assemblies in Figures II.1, 2 and 3. All these spectra were obtained using twenty energy groups with the, ORNL analog computer. In Figure II.4 we give the neutron energy spectrum at 300 micro seconds in graphite for various geometrical bucklings, B2. All these curves were obtained using fourteen energy groups with the help of the analog computer of the Nuclear Engineering Department of the University of Michigan. In Chapter VII, a complete set of curves were given for the case of B2 = 18.5 x 10-2 cm2. These curves were obtained using twenty energy groups by the ORNL analog computer -139

NEUTRON DISTRIBUTION, (E, f) 0 I ~~~~~~~~~~~~~~~0 00 _ CD 0 r) 0 ND0 o'5 H H o bD H' b) 01 -jFjj Z_ _ _ HO-' P1O 0 m 0,n CD c 0 P02 0 z +u CD 1:i H' II i Y II H' 0 0 01~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~000 0 ~~~~~~~~~~~~~~~~~~~on CY,\ ~\ o v i0 k-~~, O1 lil CI I* ~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ i'o~.~~~~~~~~~~ dct- iN ( 1 -

-141100.0 80.0I I I I 1-I ttsQ 60.0 l I. 40.0 | < -;-| | | / t' (50 + ts,u sec XH 40.0 E0.0 E G__ 2Finite0 + )f G 10. 1 S.0 6.0 2 __ rJ 0.6 0.4 0.2 0 0.2 0.4 0. 1.0 2.0 4.0 6.0 10 E ENERGY (KT UNITS) —Finite Medium o! Geometrical Buckling B2 = 3+96 x

NEUTRON DISTRIBUTION,.(E, f) ooq, 09a~~~~~~~~~ H (D H,0 CD L CD O tm _ _ __o HO~~~~~~~~ H-, C)C~ O ~m0 u, + ~ ~ ~ ~ ~ (D c iao +' 0 ]:: 1-b t cn 0 Tz — T P~ 9~ Q~ U2 i o )o + (D~~~~~f cH 0 I' as ~~~~~~~~~sa i-iota0~~~~~~~~~~~~~ i' mCC td O QFj COH o\ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~~~w~ a~~~~~~~~ bN 01\~~1 o x

c(El) N -~~~~~~~~~~~~~~ O~~~~~~~~~~~) ~~~~)N N (D OD 0 N p0, _ _ cD~~~c *0 $t 0 PD W~ /1 d —~- rn uI \ n (D r I ~ an( z~~~~~~~ N ~ ~ ~~~~~~~~~~N (DCD x x H-c (D-.( D M FDas 0 0 / PY / mb o o *i Cf) M o-d r:: 3 3 3 CD O P\- 01 0 x {::bi FO -J- Mx III C) N ~~~~~~N Hd - c 0 0 0 -', 00 FI~ O k5

APPENDIX III. HIGHER SPATIAL MODE CONTRIBUTION The neutron distribution as a function of energy, space and time variables can be represented as follows: cp(E,rt) Z m cpm(r) Pm(Et) (III.1) m where m represents the m-th mode, Am is the amplitude of m-th spatial mode, cPm(r) spatial distribution of m-th mode, cpm(E,t) energy-time distribution of m-th mode. In order to determine cp(E,r,t) completely, we need cpm(E,t) and Am for each m-th mode. We have shown in Chapter VII that the neutron distribution cp(E,t) at times greater than the thermalization time can be represented by a single exponential decay. We also note that for the Fermi Age slowing down model, the slowing down time is constant and does not vary with the spatial mode or the geometrical buckling. This does not hold good for the other slowing down model. For times greater than the thermalization time III.1 can be expressed in the following manner for a rectangular block. -B(r) e (III cp.(E,r,t) = Z en(r) e" X mnt (111.2) where: Y,,m,n represent the integers which characterize the spatial mode in a rectangular assembly. -z144

-B2T e amplitude factor for i,m,n spatial mode, kX,m,n -= decay constant for O,m,n spatial mode. 2 -2 -2 We take the specific example of Bl = 2.96 x 10 cm This is one of the assemblies used by de Saussure and Silver (ORNL-2641) in their experiment on beryllium. We take T = 80 cm2 and calculate Xe,m,n for different spatial modes using the decay constant formulas given in Chapter VII. We rewrite it: X = 1.267 x 105 [1-0.89B2] B2 For t = 114. sees, we have: p(E,r,t) = [p1lll(r) + 6.37 x 10-2 cP112(r) + 4.17 x 10-3 22(r) + 6.97 x 1044 p13(r) + 2.77 x 104 p222(r) +.] (II"'3) For t = 200 p secs, we have: cp(E,r,t) = [111(r) + 4.77 x 10-2 112(r) + 2.37 x 10-3 (r) + 3.30 x 10-4 13(r) + 1.21 x 10 4 p222(r)] (III,4) From the above analysis, it is apparent that all the higher spatial modes except (112) contribute less than 1/2%, The contribution by (112) diminishes with time, leaving the fundamental as the dominant mode. Some of the next higher spatial modes are suppressed by the judicious placement of the source and the detector. In de Saussure and Silver's experiment next higher mode to the fundamental is (113). The contribution to this mode is less than 0.1% of the fundamental.

-146It may, therefore, be concluded that at times greater than the thermalization time, only the fundamental spatial mode is important. During the interval between the slowing down time and the thermalization time, one has to use the exact PQ,m,n(E,t) for each (e,m,n) mode. We have limited this study to the generation of c(E,t) for a few geometrical bucklings used in the experiments. In order to undertake a complete modal analysis, (P,m,n(E,t) will have to be generated by the additional machine calculations or by the interpolation method using the curves given in this study. It should, however, be pointed out that for times slightly greater than the slowing down time, the distribution function will be strongly dependent upon the initial source distribution. The initial source distribution given by the Fermi Age model will also have to be replaced by Eriksson's distribution described in Chapter III. No attempt is made to undertake this study, here' it is a complete problem in itself.

APPENDIX IV COMPARISON OF THE STEADY AND TIME-DEPENDENT ENERGY SPECTRUM IN THE FINITE MEDIUM The energy spectrum of a decaying pulse of neutrons was found to be softened in the finite medium compared to the Maxwellian distribution. However, the energy spectrum due to a steady source is hardened spectrum. In the heavy gas model, under the diffusion theory approximation, the neutron energy distribution for the time independent case can be given as follows: a2m(E)+ E ~( + DB2 E d + EdcE + cp(E) (p + p(E) - S(E) (IV.1) dE2 dE s() (v.l) where S(E) = External Source = b(E - Eo) This equation has been studied extensively by Hurwitz et al.(1) and also by Cohen(2) in some special cases. We shall deal with two cases. Case A: Absorption cross section a = constant Diffusion coefficient D(E) = constant Case B: Absorption cross section follows 1/V law Diffusion coefficient is constant., For Case A, the differential Equation (IV.l) can be transformed into the confluent hypergeometric equation form by using cp(E) = r(E)Ee-E. The solution of the transformed equation is given by -147

Jahrnke and Emde. (3) It is as follows~ cp(E) A Ee-E 1F1 [o,2,E] (IV.2) where 2 ) (E)B2 + B a(E) constant (e a(a+l) E2 lFl (X,2.,E) -1 + - E +) +.. For the Case B, we'make the following additional transformation E1/2 = x Transformed equation is as follows: x d.2(E) + (3 - 2x2) d'(E) - 4(x cB2 - O)(X) = O (IVo3) d.x2. dx where: D(E) B2 B2 =y a nd.o =a x We seek the series solution for the differential Equation (IJv.3). r(x) = Z an xn (IV.4) On substituting the above series into the differential -Equation (IV-3), and equating the coefficients of xn+l, we get the following recurrence relation between the coefficients. {2(n-2) + 4~B2) an_2 + 4t aan-l (Iv.5) n(n+2)

-149On evaluating the coefficients, we get the following solution for cp(E) cp(E)= ao Ee-E [1 + -Ca %/22+{ 1 E 3 8 B 3 E3/2 40{a 4CTa (4oc:)3 (IV.6) 5{ 3 (4a g2+2) + 8 (N{aB2 + ~15 3 8 3 E2 (_+_g2) (' _ ta 4L 4aa (4a)2 (o2 (42 + (I 15 ( 3 (4aB2+2)+ 8 (++B2 (4 )+e] For the case of zero absorption, the solutions given by the Equations (IV.6) and (IV.2) become identical. We note the following features of the solution for cp(E) as given by the Equation (IV.6). 1. For zero absorption and infinite medium, the neutron egergy distribution is Maxwellian. 2, For small 1/v absorption and infinite medium, the spectrum is Maxwellian. 3. For large 1/v absorption and for B2 p 0, the spectrum is hardened. For a special case of B2 = 2.96 x lo3 cm-2 for graphite, we calculated (1) according to the pulsed and steady state sources according v to the following definition. E = 10 KT 1 (v E= 10 7) f cp(E) dE 0

-150We have the following results for the pulsed,steady state and Maxwellian cases, (1) 936 X 10-1 v)Steady state = 8,o6 x 10-1 (!4Maxwellian = 8.84 x l0The degree of hardness of the spectrum is not equal to the degree of softness. Shift of the energy spectrum occurs in two different directions for these two cases. The deviation is not symmetrical about the Maxwellian spectrum. REFERENCES 1. Hurwitz, H, H, Jr., Nelkin, H, S., and Habetler, G. J. (1956) Nuclear Scio Engng. 1, 280. 2, Cohen, E, R,, (1957) Nuclear Sci. Engng. 2, 227. 3. Jahnke, E., and Emde, F., "Tables of Functions," (1945) Dover Publications, New York.

APPENDIX V. DIFFUSION PARAMETERS OF BERYLLIUM AND GRAPHITE We list the diffusion parameters of beryllium, as compiled by.( Saussure and Silier~l) CO (2) de Saussure and Silver(l) and of graphite as given by Dardel and Sjostrand.(2) TABLE V. 1 BERYLLIUM (DENSITY 1.85 gms/cc)(l) Diffusion Coefficient Diffusion Cooling (C) Authors (Do) x 10-5 (cm2se) Coefficient Antonov et al. 1.17 + 0.05 2.5 + 0.9 Campbell et al. 1.25 0 Kloverst,rom et al. 1.24 + 0.04 3.15 + 0.65 de Saussure and Silver 1.25 + 0.0o6 1.1 + 0.9 TABLE V.2 GRAPHITE (DENSITY 1.60 gms/cc) (2) Authors Do(cm2/sec) x 10-5 C Do(cm4/sec) x 10-5 Beckurts 2.13 + 0.017 16.3 + 2.5 Antonov et al. 2.07 + 0.03 12.5 + 2.0 Decay- constant X = ZaoVo + DoB2(1 CB2) REFERENCES 1 de Saussure, G. and Silver, E.G. ORNL-2641 (1959). 2 von Dardel, G. and Sjostrand, N.G. Progress in Nuclear Energy Series I, Vol. II, (1958) Pergamon Press, New York. -~L5~L

UNIVERSITY OF MICHIGAN 3 9015 03695 1641