THE UNIVERSITY OF MICHIGAN INDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING REACTIVITY EFFECTS OF MODERATOR EXPULSION IN AN ENRICHED, LIGHT WATER REACTOR Jerome L. Shapiro A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in The University of Michigan 1960 September, 1960 IP-459

Doctoral Committee: Professor William Kerr, Chairman Professor Henry J. Gomberg Assistant Professor Geza Gyorey Associate Professor Terry Kammash Professor Gordon J. Van Wylen

ACKNOWLEDGEMENTS The author wishes to gratefully acknowledge the assistance of the staff of the Phoenix Laboratory. In particular, W. R. Dunbar, P. Herman and R. H. White gave valuable aid in the construction and operation of the experimental equipment. A. Veliko, of The University of Michigan Computing Center, and B. Smith of the General Motors Research Laboratory gave much expert advice on running the computer problems. The encouragement of the members of the doctoral committee is also acknowledged with thanks. ii

TABLE OF CONTENTS Page ACKNOWLEDGMENTSD+. G N.........................................o. ii LIST OF TABLES*.....................0...........0.... vi LIST OF ILLUSTRATIONS,................................v.. vii I. INTRODUCTION............................. 1 II. EXPERIMENTAL PROGRAM............................... 5 A. Gamma Ray Densitometer.............................. 5 1. Principle..........5.................... 5 2. Densitometer Design.................... B. Air Flow Measurement Technique................1....... C. Density Measurement...........,................. 17 1. Calibration.................................... 17 2. Procedure................... 22 D. Air Supply...........2*....................... 24 E. Density Determination Results..................24 F. Reactor Experiment.............................. 28 1. Core and Instrumentation........................ 28 2. Control Rod Calibration............... 31 35 Procedure*.....................................* 31 G. Reactivity Measurement *........................*..... 35 H. Error Analysis......,,......................... 39 1. Void Determination..................... 39 2. Reactivity Determination,.................... 4o I. Results of Early Experiments................... 41 III. THEORETICAL ANALYSIS................................. 44 A. Introduction.........................44 B. Methods........................................ 45 iii

TABLE OF C ONTENTS CONT'D Page 1. Standard Group-Diffusion.,.................. 45 2. Group Perturbation.................................... 46 5. Buckling Calculation......................... 47 C. FNR Reactor Physical Properties...................... 49 1. Standard Fuel Element............................ 49 2. Control Elements..................... 49 3. Reflector Elements.................................. 49 4. Grid Plate....................................... 49 D. Nuclear Constants......................... 50 1. Basic Data.......................................... 50 2. Infinite Multiplication Factor.........54 E. Second Order Reactivity Effects......................... 55 1. Thermal Disadvantage Factor................... 55 2. Neutron Temperature.......... 58 5. Neutron Streaming E ffect63.... 65 4. Transverse Buckling,..... 65 F. Results of Analysis..................................... 67 1. Comparison of WANDA Perturbation, WANDA Criticality and PDQ(XY) Criticality Methods........... 67 2. Transverse Buckling Calculation................... 75 IV. SUMMARY OF RESULTS............................... 82 A. Localized Voids.......................................... 82 B. Distributed Voidss....................................... 82 V. CONCLUSIONS................................................ 84 A. Experimental Technique................................. 84 B. Analysis *.............................................. 84 RE7FEtENCES..........................................88 iv

TABLE OF CONTENTS CONT'D Page APPENDIX A. 7EFECT OF NON-HOMOGENEITY ON VOID MEASTJREMET WITH THE GAMMA RAY DENSITOMETER...,.......................,. 92 B. AGE CALCTJIATIONS IN MIXTURES OF ALUMTINUM, WATER, AND AIR.... 95 C. LITERATURE SURVEY*.....*............ *O.. O........ 98 D. SAMPLE DATA AND CALCULATIONSo o........................... 106 E. REACTOR PHYSICS NOMENCLATURE............................... 112 V

LIST OF TABLES Page 1. Relative Experimental Void Coefficients of First Runo..,,... 37 2. Parameter Substitution for Adjoint Calculationo.............o 47 3. Components of Standard FNR Fuel Element o...,.,,o...,,,.....o.. 50 4, FNR Core Nuclear Properties56o. ~ oo. o o. oooo o.o.. 56 5. FNR Reflector Nuclear Properties 00...000.....0...0..00.... 57 6 Relative Effect of Variation of F on Reactivity.......o... 58 7. Neutron Streaming Effect on the Void Coefficient............. 65 8. The Transverse Buckling Variation and Its Effect on Void Coefficient, by One Group Approximation............... 67 9, Nuclear Data for Core Used in Comparison of Analytic Methods.....6.0............0 0.00.0a0...o..Uaaoao a...o.o. 0 0 0 o 0 0o 68 10. Nuclear Data for Reflectors Used in Comparison of Analytic Methods............................ o 70 11. Transverse Bucklings Calculated by the Buckling Iteration Procedureoo O O oO.....................................................o 77 12. Results of WANDA Series 2 for Constant and Variable Bucklings o0 0 OoO o000o 0oD 78 135 Calculated Void Coefficients Assuming Constant Bucklings..... 79 B. Measured Age in AI-H20 Mixtures.,..,.. o,...,o.. o.. o,o,...o 96 C. Characteristics of SPERT-I and FNR Coreso. o,,,oo...ooOo...oo 103 vi

LIST OF ILLUSTRATIONS Page 1. Schematic representation of the densitometer.................. 7 2. Light output spectrum of a scintillation crystal due to a mono-energetic source................................... 7 3. Pulse height spectrum of a NaI crystal and multiplier phototube for a mono-energetic gamma ray source*. **...****... 7 4. Detector circuit schematic*............. *+.*........ 9 5. Gamma ray densitometer and traversing stand.......1...... 11 6. Traversing mechanism - Plan view,.*........*........ 12 7. Standard FNR fuel element......1................ 14 8. Air flow system - schematic*.*m*s.*.**c**.*a.***..**..***.* 15 9. Calibration of the window of the differential analyzer,,,,,,,, 18 10. Pulse height spectra for two values of AE using the Cr51 source,4...*.............*.*..........**.9 *.*... *9**.. 19 11. Measured distribution of void in element,................ 25 12. Flowmeter calibration in terms of average air density in element and relative pressure* *........ *... *.......* 27 153 Variation of average bubble rise velocity with average air volume fraction,,.......*.......-.***...*.......*.*..* 29 14. Reactor core configuration*........................... 530 15. Pulse channel (modified) of the FNR...................... 32 16. Differential control rod calibration..........,............ 533 17, Integral control rod worth...*................* *........... 34 18, Experimental effect of voids on reactivity,.,.3.**........ 36 19. Experimental effect of voids on reactivity (leaky tubes, unplugged element)..*..... *........*.... *...*... ** 38 vii

LIST OF ILLUSTRATIONS CONT'D Page 20. Experimental effect of voids on reactivity (from MMPP-110-l) for lattice position 355..............~....................... 43 21. Variation of thermal average cross sections with U/H (SOFOCATE data)................................52 22. Variation of thermal average cross sections with water density...........*.... o 60 23*. a) Configuration used in one-dimensional calculations.....*, 69 (b) Configuration used in two-dimensional calculations...,,,, 69 24. One-dimensional, two-group fluxes and adjoints,............... 71 25. Calculated reactivity effect of void assuming B2 =.0034106,, 72 26. Variation of transverse bucklings with void concentration,.... 76 27. Variation of void coefficient during buckling iteration subroutine of MAGNUM,........................................ 81 A. Maximum possible error in void measurement due to temporal fluctuations in density...................................... 94 viii

I. INTRODUCTION The study of the effect of moderator expulsion from the core of a thermal nuclear reactor is generally stimlated by concern over the inherent safety of the system. Given a situation where all external shut down mechanisms are removed and a positive reactivity applied, it is desirable to be able to predict the course of the energy release. In the case of water moderated reactors the characteristics involved are the neutron chain reaction kinetics, heat transfer, hydrodynamics and the coupling between the material properties and the reactivity. A thorough understanding of these four subjects is sufficient for solution of the problem. Another situation where these factors interact is in a reactor having a boiling moderator. As the boiling rate increases in reactors of this type the statistical fluctuations of power level increase. At a particular power level, varying from one reactor to another, periodic oscillations of an unstable nature begin. This boiling instability is not well understood at this time. There is no concrete evidence of what factors control the initiation of instability. However, it is perfectly clear that the key factor permitting power oscillations to exist is the interdependence of the power and the reactivity. This relationship may be expressed as: k(t) = k(o) + Cp(P) x [P - Po] where k(t) is the reactivity at a time t, k(o) is the reactivity at a point of reference t = 0 when the power P = Po and Cp(P) is the power coefficient of reactivity. -1

-2Ulnder some important circumstances, namely during a rapid transient from an initial subcooled state or for any transient from an initial boiling state, the major component of the reactivity coefficient is due to the expulsion of moderator by the formation of vapor. The reactivity at any time will be considered to be a function of the physical properties of the reactor at that time. In other words the neutron distribution may always be assumed to be the equilibrium distribution~ This is true for reactors small in units of the migration length of neutrons (a common circumstance), ThIerefore the study of the effects of physical changes on the reactivity, performed on a steady state reactor can be successfully applied to reactor kinetics0 This feature, in reality an approximation, is true for all but extremely large (larger than ordinarily encountered, even in hazards analyses) variations of reactivity0 The major purpose of this investigation is to develop and test techniquaes for accurate measuzrement and calculation of the effect of the presence of vapor bubbles on the reactivity of a water moderated reactor. Extensive use is made of the Ford Nuclear Reactor for experimental purposes and an IBM 704 digital computer with several codes programmed by others for reactor physics analyses. In order to evaluate the importance of each of the various effects which contribute to the void coefficient it is desirable to calculate each component separately. This is better from a research point of view, then lumping all the calculations into one program which, given all the input data will calculate the reactivity. Therefore, the procedure of separate calculations is followed.

There is some advantage to the use of the FJNR for experimental verification of the theory. In general it is more difficult to calculate the nuclear parameters of reactors of small dimensions rather than large, with reflectors rather than bare, with rectangular geometry rather than spherical or cylindrical, or having hydrogenous moderator rather than a heavier moderator. The FNR core is small, rectangular, completely reflected (with water and graphite) and water moderated. The successful method then, should be applicable to larger, less complicated reactors. The major feature absent in this core, which is present in many other reactors is resonance absorption, Modification to include this effect in the calculations must then be made for calculating void coefficients in those reactors. The experimental technique consists of bubbling air through the core, one element at a time, measuring the reactivity changes by means of a calibrated control rod. Common methods of producing voids in reactors (for void coefficient measurements) in the past have been to insert expanded plastics containing closed air bubbles, hollow tubes or boxes, aluminum, or magnesium. The use of heat, to produce boiling is possible but the consequent temperature effects are difficult to evaluate, The expanded plastics have been shown to be the best, requiring a very small correction for the difference between plastic and water, It is difficult however, to handle in sheets thin enough (l/10") to place between the plates of our fuel elements. There is some uncertainty in the sizeable correction for aluminum absorption and scattering when using that material. Magnesium requires a smaller correction than aluminum and is therefore the better metal for use in simulating voids. However, results of some SPERT investigations(28'50) show a serious discrepancy

-4between magnesium and STYROFOAM (plastic) coefficientso A rough comparison of the effect of magnesium sheet and air bubbles is made in the FNR during the course of the experiment. The prime difficulty associated with the use of air bubbles is the measurements of the void volume~ This was accomplished, by means of a gamma ray penetration measurement~ Although this principle has been used in the past, never has the requirement for accuracy been as stringent as in this application0 The unique features of the problem here are the large test section size (3" x 3") and the small air volumes (1-10o) being measured. The effect of the voids, or vapor (vapor is essentially void as far as neutrons and gammas are concerned) is stated in this study in terms of the reactivity coefficient of voids, called simply the Void Coefficient and is defined by Void Coefficient a / i k v where Lk is the change in multiplication factor from the original value k, due to a change Av in void volume from the original water volume vo This volume (v) is the water which is affected by boiling0 The water in the control elements and reflector does not ordinarily boil and is therefore not included in vO

II. EXPERIMENTAL PROGRAMS A, Gamma Ray Densitometer 1. Principle Given a well collimated beam of mono-energetic rays of intensity Ion the intensity I after passing through matter is given by I = Io e (1) In this formula p. is the interaction probability per unit length of the material and t is the total distance traversed by the rays through the material, There are, in general three different types of interactions that are possible between gamma radiation and matter, They are: (1) Photoelectric effect - Absorption of the photon by the electrons of an atom. (2) Compton effect - Scattering of the photon by an electron. The scattered photon has a different direction and a lower energy than the original photon. (3) Pair production - Conversion of the photon into an electron and positron pair. In all three types of interaction, although the primary photon disappears, secondary particles and/or photons are produced. The simple concept of formula (1) holds then only when I is the primary particle intensity, Therefore p is the total interaction probability per unit length, The product pt, the number of "mean free paths" between the source and the detector is directly proportional to the density of the material. -5

a6Schematically, the problem presented in the densitometer as applied here is shown in Figure 1o n aluminum plates are separated by n-l water channels. The intensity I as a function of the density of water pw is given by: I = Io exp - [kAln tA1 + w(n-l)twp/p] o (2) where Al = interaction coefficient per unit length of aluminum, tAl = thickness of aluminum plate, W 3= interaction coefficient per unit length of water of density po. If the counting rate of the detector when p = p0 is Co then the count rate C as a function of p will be LwCn-l)tw-(l-p/po) C = Co e () If;o = 1 and X 1-p, (4 then, rearranging (3), X = -_ 1 (5) pw^(n-l)t Co Since any attenuating material between source and detector other than water maintains constant density, it will not affect the ratio C/Co Thus the effect of structural material is only to reduce the magnitudes of C and CO but not their ratio. Of course, any background counts must be subtracted from both, before calculating the ratio, The detector used in this experiment was a NaI(Tl) scintillation crystal0 This material emits photons in the visible light spectrum when subjected

-7WATER SOURCE -7 a - W W | DETECTOR SHIELD - - AND COLLIMATOR - Figure 1. Schematic Representation of the Densitometer U. 0 a:r/ rPHOTOELECTRIC a Z COMPTON INTERACTIONS Cz COMPTON 3P Q INTERACTIONS z \ OUTPUT LIGHT INTENSITY PER INTERACTION Figure 2. Light Output Spectrum of a Scintillation Crystal Due to a Monoenergetic Source U. W. \ y^\ PHOTOPEAK aD Z Cm Zn / I GAUSSIAN-N w PULSE AMPLITUDE Figure 5. Pulse Height Spectrum of a Scintillation Crystal and Multiplier Phototube, For a Monoenergetic Source

to radiation. When an interaction occurs between a gamma photon and the crystal the number of photons emitted is directly proportional to the energy absorbed during the interaction. Therefore when a beam of monoenergetic gamma rays impinges on the scintillation crystal the relative number of interactions plotted as a function of the intensity of output light pulses per interaction would appear as in Figure 2. The vertical line represents the photoelectric type of interaction in which all the energy of the gamma photon goes into the light pulseo The continuous spectrum is due to compton interactions in which only a portion of the energy is re-emitted in the light pulse It is assumed for simplicity that the gamma ray energy is below the threshold energy for pair production (a 1 Mev), This output spectrum permits the counting of only primary gamma rays. Connecting the crystal to an electron multiplier phototube, a fraction of the light, pulses are collected and transformed into an electronic pulse whose magnitude is proportional to the input light intensity. Fluctuations in the multiplication of the phototube spread the single line of Figure 2 into a Guassian curve. The phototube output appears as in Figure 3, The pulses from the phototube are amplified and fed to a set of discriminator circuits, called a differential analyzer, which allows the counting of pulses within a particular amplitude band. The pulses falling within the desired band are then counted with an electronic scaler. A schematic diagram of the electronics is shown in Figure 4. 2. Densitometer Design To reproduce the hydraulic effect of the realtor pool water above the core the original design called for the densitometer to be submerged near

HIGH VOLTAGE POWER SUPPLY 40 v. D. C. CRYSTAL LINEAR DIFFERENTIAL SCALER 4HOT-OTUBE — AMPLIFIER ANALYZER PHOTOTUBE FPRE-AMPLIFIER 4 D Figure 4. Detector Circuit Schematic

-10the bottom of the pool itself. Therefore, the detecting crystal and phototube were packaged in a plastic watertight case along with a transistorized preamplifier. The traversing mechanism was also designed to operate remotely. It turned out that the radioactivity in the pool water was very high for the sensitive detector and furthermore was not a sufficiently constant background. The equipment was then modified to operate with its own water column outside the pool. The schematic drawing of the equipment is shown in Figure 5, The traversing stand is designed to allow vertical movement of the plate to which the source and detector are fixed. The water chamber, built of welded aluminum plate will house one fuel elements A standpipe, twenty feet high can be bolted to the water chamber. The chamber rests on a movable table which permits lateral motion between the source and detector of about four inches. Figure 6 shows a plan view of the traversing mechanism. The choice of source was determined by consideration of several factors. a) Energy - The lower the energy of the emitter the greater its interaction cross-section and the greater the sensitivity of the measurement as shown by Equation (5), Of course, as the energy is reduced, the source intensity required increases but this is not a problem in this case. b) Number of Gamma Ray Energies Emitted - For most accurate analysis, the gamma ray selected should be the highest energy present in the decay scheme, c) Half-life - The lifetime of the emitter should be as long as possible, at least several days.

ER"-TICAVERTICAL TRAVERSE LIFT RODS If I I STANDPIPE PRESSURE-7 r~r""" ~ -* GUAGE^~..!- I _ I I-' ~40 v. - I -TEST ELEMENT l BATTERY HIGH VOLTAGE eow'= <'.' ll~ POWER SUPPLY ARMETER,,/ -WATER * 2 — CHAMBER, o AMPLIFIER AIR, o PRESSURE o I.L. 111 1. ACCUMULATOR I I DETECTOR AND PRE-AMPUFIER ~ -' |..............,,,'- I ETIADIFFERENTIAL AMMARA I'.; ANALYZER SOURCE i HORIZONTAL TRAVERSE SCALER 1I MECHANISMCAL INLET PRESSURE i-TRAVERSING STAND REGULATING VALVE Figure 5. Gamma Ray Densitometer and Traversing Stand

T WATER DA Figure 6. TraCHAMBER 0 SOURCE CRYS MULTIPLIER PRE ^ PHOTOTUBE AMPLIFIER I~ A \^HORIZONTAL {c 3 JTTRAVERSE DIAL Figure 6. Traversing Mechanism-Plan View

-13d) Availability - A common material, preferably something that could be produced simply in the Ford Nuclear Reactor is desirable. For these reasons Cr-51 was chosen as the source. It emits a single gamma ray of energy 0*32 Mev with a half-life of 27 days. The source consisted of several pieces of chromium, that had been irradiated for about two weeks in the reactor. Analysis by induced radioactivity and x-ray flourescence techniques indicated that the metal consisted of 98% Cr, 2% Fe and a trace of W. Irradiation of this mixture produces predominately Cr-51 which has a half-life of 27 days, some Mn-56 having a halflife of 2.5 hours and some W-181, half-life 24 hours. The chromium pieces were inserted into a 3/8" diameter plastic test tube about one inch long. This tube was then placed in a lead container having a wall thickness of 1-1/2". A 1/4" square hole was left in the front cover. The detector, a 1" diameter by 1-1/2" long crystal was shielded in front by 3/4" lead having a 1/4" square hole in the center. Three 1/8" sheets of lead were wrapped around the detector and phototube with a grounded copper shield underneath to lower the background as much as possible, A dummy fuel element, identical to those used in the core except for the absence of uranium was used. This element is an eighteen curved plate assembly of the MTR type, as shown in Figure 7. B. Air Flow Measurement Technique A schematic diagram of the air flow system is shown in Figure 8. In out -pool measement the rchamber and stack are used to substitute for the reactor pool,

-14U -AI ALLOY,_ ~-~I I __ FUEL PLATE SECTION I I ~~~I ~ 2S ALIMINUM END VIEW I I I I I BOX I I I I I v I I I I I I I I I I SIDE VIEW Figure 7. Standard F Fuel Element SIDE VIEW

~-\ fFLOW ADJUSTING ROTAMETER CCUMULATOR o TANK 60PSI AIR SUPPLY PRESSURE o REGULATOR 20 0 REACTOR L I POOL 0 o 0 0 Figure 8. Air Flow System —Schematic | I.. 2I 1::i!! CORE 1

As the flow meter is a carefully calibrated device an attempt was made to predict the average density of air-water mixture in the fuel elements Assuming the air absolute temperature change from the flowmeter to the pool water is negligible, the average air volume fraction at a given height in the element is P1 F1 X =- _, (6) P2 V2 A where P = absolute pressure at the flowmeter, P2 = absolute pressure at the particular height in the element, F1 = air flow rate through the meter, VP = bubble rise velocity, A fuel element free area. This equation was used to calculate X in the original calibration of the Ford Nuclear Reactor (MMPP-llO-l)o There the bubble velocity V2 was measured by abruptly starting flow, timing the first bubble as it appeared at the top of the element. The velocity obtained was 1.1 feet per second, which is within 15% of results obtained for the rise of single bubbles in water (with no wall effect) by Miyagi(6) A review of the literature(28) concerning the rate of rise of air bubbles through water indicates that the rise velocity is a function of the density of the mixture. It is also shown that the distance between the walls strongly affects the rise velocity, This feature of the experiment makes the radiation penetration method the best way of measuring the density.

-17The gamma ray density measurement, combined with the calibrated flowmeter, can be used to calculate the average rise velocity of bubbles, The results of this calculation are shown in Figure 13. C. Density Measurement 1. Calibration Various tests were performed with the electronic equipment to check its characteristics and determine the best operating conditions. A precision test pulse generator and several known sources were used. The calibration of the E and AE dials of the differential analyzer is shown in Figure 9. The E dial determines the lowest energy pulse counted. The highest energy pulse counted will then be E + AE. The calibration determines the dimensions of the AE setting in terms of the E setting. This relationship is not a function of E setting. Therefore, for a particular gain and voltage, when the E dial is calibrated in energy units, the AE dial calibration is calibrated at the same time. As we were interested in counting all the pulses in the chromium photopeak it was necessary to adjust the gain and voltage of the instrument such that the entire photopeak appeared within E settings of no more than 55 divisions since as shown in Figure 9, the maximum window width of 1000 corresponds to 57 pulse height units on the E scale. An example is shown in Figure 10. For an amplifier gain of 3504 and a phototube voltage of 775 volts, counts were taken as a function of pulse height selector setting (E) for two values of window width (AE). It is demonstrated that with the widest window, almost all the pulses in the photopeak are counted at one time (PH. = 55). The advantages of using this setting

70 60 50 40 - ------- 0 0 I 60 L = 700 40: 4 —-200 400 600 800 1000 AE, WINDOW WIDTH DIAL SETTING Figure 9. Calibration of the Window of the Differential Analyzer i 200 400 6Q0 800 I000

-19E =1000 9 8 7 6 Iz 5 8 4 AE=40 8 6 5 -- 50 60 70 80 90 100 110 120 130 140 150 160 E, PULSE HEIGHT SELECTOR SETTING Figure 10. Pulse Height Spectra for Two Values of the LE using the Chromium-51 Source V = 775 v. Gain = 3.04

-20over the smaller window width are obvious here. The total number of counts is higher giving better counting statistics. Furthermore, since the upper curve is flatter than the lower one, small fluctuations in the pulse height setting do not cause as great an error, In the measurements the settings used were: V = 775 volts Amplifier gain = 3o04 E = 50 tE 900 Having checked and calibrated the instrument it remained to check the principle in operation before proceeding with the measurements, A method was evolved which not only checked the principle, but also corroborated the value of [w as given in the literature. With the equipment set up in operating condition four sets of counts were made: Co = count rate with water chamber empty, C1 count rate with water chamber full of water, C2 = count rate with fuel element inserted into water, C3 = count rate with fuel element inserted and water drained, It is simple to show, using Equation (1) for the four different cases that 1 Co, C2 w=.ln (9) n tAl C1 C The above method was carried out in two different ways. First, with AE = 20, a spectrum of count rate versus pulse height was made for each of the four conditions noted above, for pulse heights varying from 50 to 100o

-21The area under the Gaussian portion of each curve was numerically integrated, These areas were then substituted into Equation (9) and the result was 4w = 0.117cm'1, This coincides with the value given by Goldstein7) The reason that only the Gaussian portion of each curve was selected is that it is proportional to the primary photon intensity. The difference between the rest of the curve and the Gaussian (see Figure 5 for example) is due to two factors - first, to Compton type interactions of the primary photons in the crystal and second, to photons which are scattered (Compton) in the attenuating medium, to such small angles that they still enter the detector with energies close to the photopeak. But for this latter factor, the area under the entire curve would give exact results. Since the method used in the final measurements does, in fact, integrate over the entire photopeak, the extent of this buildup factor was measured by repeating the four sets of measurements with an instrument setting of E = 50, AE = 900. This gave the result tw = 0.107 +.002cm-1. The apparent reduction in [i is due to the buildup contribution. The effect of this factor on the measurements may be determined from this experiment. For this purpose, Equation (1) can be modified to (see Goldstein47, pg. 186): I = Io e t [1+ But] (1) The buildup factor B as applied here is in very much a function of the geometry used and the window width of the instrument and its analytic determination would be extremely involved. Therefore B was determined from the two measurements just described,

-22Using Equation (12) instead of (1) we find that for the measurement including buildup, =_ 1 S n ^~ C2 1 + B[}w(n-l)t + IAlntAl] In - in w n tpAl \ C C i [l+B(pw(n-l)tw + twntAl)][1 + BAlntAl] Knowing the true value of l., from the Gaussian analysis and using pA1 = 0.273 (from Goldstein, 7 pg. 152) thie factor B can be determined. Equation (3) is then rewritten: pw(n-l)twX I1 + Bl(n-l)tw(1-X)c 0e [l + B-(n-l)tw for B <.25, the numerator in the brackets of Equation (3') can be expanded as an exponential producing the following correction to Equation (5): X(l-B) + B = 1 fin C + in [1 + Btw(n-l)tw]} (5' ) vtw(n l-t1 l Co If B is small enough (< 0.15) the approximation ln[l + Btw(n-l)tw] ^ B4w(n-l)tw may be applied to Equation (5') with the result 1 1 C x =lB -(n-l)tw (5")CO 2. Procedure With the twenty foot header attached to the water chamber, the voiding tubes inserted in the fuel element and the element inside the chamber, the apparatus was filled to the overflow spout located one foot from the top of the header, As air is introduced, water will flow out the spout keeping

the hydraulic head constant. Care was taken that void volume was always increasing from run to run to keep the water level constant (or, if necessary, water was added), The source-detector platform was adjusted to the desired height. Then a background count was taken with the source removed. When the source was replaced it was not disturbed during the entire run since movement gave rise to noticeable changes in count rate. Although the source and detector each have 1/4" dimensions, the rack and pinion device which moves the water chamber in the horizontal plane translated 0.191" for every quarter turn of the pinion. Using the quarter turn as the basis then, counts were taken for twelve different positions. Therefore, the total length of the element surveyed was 12 x 0.191 + 0.250 - 0.191 = 2.35" The actual distance between side plates of the element is 2.62". Therefore, approximately 0.13" adjacent to each side plate was not surveyed. As will be shown, this had a negligible effect on the results. With the source in place one minute counts were made (approximately 12000 - 20000 counts) at each of the twelve positions of the chamber with no air flow. Then the flow was turned on, steadied and two minute counts were taken. The flow was then turned off and another set of one minute counts were taken, In this manner, each set of two minute counts at flow was bracketed by one minute "no-flow" readings. This was done to compensate for small, slow changes in the gain of the instrument and gave much more consistent results than the method of taking a set of no-flow readings at the beginning and end of the entire run (each run lasted several hours).

-24The flow meter is equipped with a pyrex float and a stainless steel float, to give extended range to the instlmento T1his was not needed though and only the stainless steel float position was used as a reading. Measurements were taken near the bottom and near the top of the element, D. Air Supply The compressed air supply line of the Phoenix building was used as a source. The 55 psig line air was passed through a pressure regulator into a 40 gallon dnrum The puzpose of the regulator and drum was to prevent severe pressure surgeso The air then passed through a flowmeter (Brooks Rotameter) monitored by a pr ssuez gauge, A needle valve on the outlet of the floTwmeter allowed adjustment of the Plow rate into a fifty foot length of high pressure rabber hose (cloth-lined) At the end of the hose a set of eight babbling t bes'was a attached, 2hese consisted of three foot lengths of 2S aluminum tibes 1/16" oDo x 0,015" wall thickness, One end of each tube was passe:d. through a cylindrical aluminum block and. sealed with sealing wax, The ribber connecting tube was then clamped to the aluminum block, The set of tubes was then inserted, between the plates of the fuel element, down to the bottom, protruding into the cylindrical end boxo E. Dens ity Determination Results The void measurement p-rocedure was followed thirty one times for varying flow rates, The result of each traverse were normalized by dividing the void concentration at each point by the void concentration integrated over the entire traverse., This gives the distribution profile shown in Figure 11, The numbers on the abscissa are the setting numbers of the pinion, which detezrmine th.e horizontal source-detector position~ The

-251.6 1.4 1.2 (I)v,~~~~~~~~~~~~~~~~~~~~~~~~~~~~~( z 1.0 w Iw cr Z 0 Z.6 w 0 o) 4 5 6 7 8 9 l0 11 12 13 14 15 1 >I 4'AVERAGE -1.62 HORIZONTAL TRAVERSE POINTS, —----- -MEASUREMENT LIMITS i~- --- -ELEMENT SIDE PLATES Figure 11. Measured Distribution of Void in Element

-26ordinate is in arbitrary unitso Extrapolation of the curve to the end plates gives a fair estimate of the error produced by not measuring the small space at either end. The ratio of the cross-hatched area to the rest of the area is the correction factor amounting to approximately 0.5%. The peak-to-average ratio of this distribution is about 1o6. Traverses were taken near the top and near the bottom of the fuel element. The resulting integrated densities were calculated using Equation (5")~ Figure 12 is a plot of X(P2/f4P1) as a function of the position of the stainless steel ball in the flowmeter. There does not appear to be any detectable difference between the measurements at the top and bottom of the element. This is understandable, since the pressure difference is relatively small (about 3%). It would be difficult to detect this small effect, There is another factor which would tend to minimize any differences, As will be shown, bubble velocities increase as the flow rate increases. Thus, as the pressure decreases in the course of the rise of a stream of bubbles the density does not decrease proportionally to the bubble volume increase. Another factor which was examined was the possible effect on the density of water circulation in the test stand, A set of measurements were made with a rubber stopper placed at the base of the fuel element, removing the possibility of water circulation through the test section. These data appear in Figure 12 also, There is no discernable effect. This indicates the lack of water circulation even without the plug in place. As will be shown later this is not true in the reactor pool,

70 60 --- - -A BOTTOM OF ELEMENT - -- -- __ -0- TOP OF ELEMENT -- TOP OF ELEMENT (PLUGGED) I- I UL 50 -I 0 LaJ w W I / C,) U) 40 w z L&J -J - CD 30L _ _ __ _ _z w w 20 w 0 LL 10 A I0 10 ------ -- -- -- ~ —- -- -- -- -- --- -- -- --— ~ — 0 0.02 0.04 0.06 0.08 0.10 0.1O2 0.14 0.16 0.18 020 0.22 024 0.26 0.28 0.30 0.32 0.34 W XI //p-, (PSI)"12 Figure 12. Flowmeter Calibration in Terms of Average Air Density in an Element and Relative Pressure

-28To give an idea of how these results compare with the experiments of others, Equation (6) was used to calculate the average velocity of the bubbles. The velocities, calculated directly from the data shown in Figure 12, are plotted in Figure 13. Although the actual velocities are strongly affected by the wall spacing, the trend of the curve agrees qualitatively with the results of Siemes(5) and Verschoor(3) and especially well with those of O'Brien and Gosling (4) F. Reactor Experiment 1. Core and Instrumentation A complete description of the reactor is given in the literature (54) It is of interest here to describe the details of the core and the specific instrumentation used during the measurements. The core is built up of twenty full fuel elements (140 gms U235) and four partial elements (39 gms U235) in a five-by-five array with one corner missing (see Figure 14)o A single row of graphite reflectors was on two sides of the core and a double row of reflectors was on the other -two sides. The core configuration was chosen to be as symmetrical as possible, simplifying experiment and calculation. With this configuration it was necessary to sample the effect of only one quadrant of the core. The test element contained a rubber stopper in its lower end box to mauke the hydraulic condition similar to that of the out-of-core runs used for void calibration. The instrumentation used to observe the reactor condition was the pulse channel of the FNR. Due to the restriction to operation at relatively low powers neither the Log N nor linear level channels were on scale during the experiment, Therefore the existing pulse channel of the FNR was modified

-29) 4w3 -- co _- - - _ _ _ __ co 4- o_1 o /o 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.10 0.11 X, AVERAGE VOLUME FRACTION OF AIR Figure 15. Variation of Average Bubble Rise Velocity with Average Air Volume Fraction

-30Figure 14. Reactor STANDARD FUEL ELEMENT Core Configuration CONTROL (PARTIAL) ELEMENT m~ GRAPHITE REFLECTOR

-501by the addition of a linear count rate branch. This permitted more accurate reactivity readings since small level changes were more readily observable. A schematic diagram of this channel is shown in Figure 15, During the course of the experiment the automatic control unit was disconnected and the control rod was adjusted manually, To measure changes in reactivity for a small change in the reactor condition, the position of the control rod was adjusted until the flux indication on the linear level recorder leveled off, The position of the control rod was then read on a scale attached directly to the rod drive extension, The positions of the shim-safety rods, which were inserted approximately one inch into the core, were not altered during measurements, 2. Control Rod Calibration The control rod was calibrated by the positive period method, at the conclusion of the experiment. Five points were taken which covered the region on the control rod in which the reactivity measurements were made, The differential reactivity curve of Figure 16, obtained from the data was graphically integrated to obtain the integral rod worth curve of Figure 17, In the calculations a value of peff equals 0o755% was assumed, This parameter, never measured in the FNR is within 3% of the value calculated by Reynolds(48) for a similar core, 3. Procedure With the bubbling tubes inserted in one element and the power level steady at about one watt the position of all the rods was noted. The air flow was then turned on and adjusted to the lowest flow rate, After a

-32FISSION CHAMBER LINEAR PRE-AMPLIFIER LINEAR T-THIS EQUIPMENT ADDED AMPLIFIER FOR THIS EXPERIMENT ONI LOG. COUNT RATE SCALER LINEAR METER AMPLIFIER LINEAR COUNT LOG. COUNT RATEl RATEMETER RECORDER (AZAR) LINEAR COUNT RATE RECORDER Figure 15. Pulse Channel of the FNR (modified)

-33U 0.02 -- Ia. i-) 0.01 L U. - LU 0 7 i 0 1 2 3 4 5 6 7 8 9 10 II 12 CONTROL ROD POSITION, INCHES WITHDRAWN Figure 16. Differential Control rod Calibration

0.20 0.15 u0I.II 0' 1 w cL 0.05 / 0 I 2 3 4 5 6 7 8 9 10 II 12 13 14 CONTROL ROD POSITION, INCHES WITHDRAWN Figure 17. Integral Control Rod Worth

-35steady state condition appeared, as indicated by constant flow rate and pressure readings, the operator was advised to readjust the control rod to maintain criticality. When this was completed, the operator (in the control room) notified an observer who read the control rod position. The flow was then increased and the measurement repeated. After the maximum flow rate condition had been measured the flow was turned off and the original (no-flow) criticality point was re-measured. This was done to show that no reactivity drifts (such as due to xenon decay) of significant magnitude had occurred during the measurement. Then the reactor was shut down, the test fuel element moved to a new position and the procedure repeated. The measurements at four positions in the core were done over a period of sixteen hours. At the completion, the test element was returned to its original position and another series of measurements was taken with the rubber stopper removed from the lower end box of the element. G. Reactivity Measurement The data and tabulated results are reproduced in the Appendix. The graphs of reactivity versus air volume content per fuel element are shown in Figure 18. The void coefficient, defined as the slope of these lines divided by the fraction of the core perturbed, varies from - 0.316 at the center of the core to - 0.120 at the corner position. The effect of the withdrawal of a magnesium (2% Zn) sheet from a water channel at the edge of the test element in position 35 is also shown on Figure 18.

0.20 I 3747573656 0.18 POSITION OF MAGNESIUM SHEET 0.16 CONFIGURATION/ A / - 0.10 I 0 08 -— 47, 56 0.08 0 M AGIU V / P LATTICE POSITIONS% AIR VOLUME CONCENTRATION PER ELEMENT, % Av/v Figure 18. Experimental Effect of Voids on Reactivity

-37The air bubbling experiment was actually done twice. The first time, all eight positions in the quadrant were examined (see Figure 19). It was found, however, that an air leak had developed. The leak was found by retesting the void apparatus in the densitometer. It appears, however, that some information can be gained from the results. That the leak was consistent was shown by measuring the void effect in one position (lattice position 35) at the beginning of the experiment and again at the end (on the following day). The results were the same, It seems also, from observation of Figure 19, that the effect of the leak on the air volume was linear with increasing volume. The reactivity data were therefore considered useful as to the ratios of the void effects in the various channels. Since the data in Figure 18 were found for four elements only, the void coefficients for lattice positions 36, 45, 47, and 56 were calculated by multiplying the measured void coefficient of position 35 by the ratios determined in the "leaky run." Table I shows the void coefficients of the "leaky run" relative to the coefficient in position 355 TABLE I RELATIVE EXPERIMENTAL VOID COEFFICIENTS OF.FIRST RUN Lattice Position 35 1 36, 45.778 55.610 47, 56, 37.500 57.343 Using this information to fill in the gaps in the data and assuming the reactivity effect of voids in the control element fuel plates are equal to that of the central element, the void coefficient due to uniformly distributed voids is - 0.190.

-380.20'- - 1374757 0.18 0.16 I — 1 CONFIGURATION 0.14 > 0.12 -O — -- >< I I I J 0.08 +_ 0.06 / 0.04/ _ 0.04 / & * LATTICE POSITIONS G -35 x - 36 0.02 + - 45 rE — 47 00A — 56....... -571 2 4 6 8 10 12 14 AIR VOLUME CONCENTRATION PER ELEMENT,% Av/v Figure 19. Experimental Effect of Voids on Reactivity (leaky tubes, unplugged test element)

-539 H. Error Analysis 1. Void Determination To calculate the integrated volume fraction of air from the scintillation counter data each total count with void is divided by the total count without void. Then the logarithm of the ratio is taken, The summation of these ratios for each step in the traverse is directly proportional to the net result. The experimental error must then be found for a function n _ _ v = In (I i) = in (Ii). (10) <- Io-i i Toi i 1 In reference (49) it is shown that 2 2 2 F..lii 1 (11) X[ v adj X (11) where t1 and.o are the standard deviations of 1 and Io respectively, and ois the standard deviation of ~o Using the usual assumption of a Gaussian distribution function for large values of I, Terslatmauedeatino - I he gen(12) The resultant measured deviation of v is then given by 1 1 c_ - ( + () 3) v (oi In the present case, the lowest value of I1 and Io is about 4 3.2 x 10 and n is twelve, for which a9 = 2.7%. As indicated in Appendix A, the error due to non-linear averaging of discontinuous voids by the gamma ray penetration is certainly less than 1% for air volume concentrations less than 10%o Therefore, this factor is neglected,

-40o There are obviously other factors which produce errors, as seen from the scatter of the data in Figure 12o These can be due to electronic fluctuations to some extent, but more than likely are due primarily to changes in the air stream, path* The approximate standard error, in units of W is 1.22o Assuming a value of P2/ JP1= 5 the error in void concentration is 0,4% Av/v. If this error was due to faults in the measurement technique, it would be meaningless as far as interpretation of the reactivity results is concerned. However, since it is probably due mainly to actual fluctuations of flow, this is directly applicable to the in-core measurements, 2o Reactivity Determination To determine that the reactor is just critical, the state of the linear count rate recorder is examinedo Since it was only observed for a finite time, it is possible that the reactor was on a period, Assume the chart was observed for a time At and the power level remained constant within AP/P. Then for large T At < g (14) T P where T is the shortest possible reactor period that will go undetected, This can then be related to the largest possible reactivity by the approximation, valid for small reactivities T = M (15) P combining these two equations we see that p -- (16) tP

The recorder was observed for 100 seconds or more, the average power level varying no more than 1%. Therefore, the error associated with each reading due to this factor is approximately + 0,001% Ak/ks The control rod scale is divided into 1/16 inch divisions and was read to the nearest 1/52 of an inch. The error associated with each reading due to rod reading is + 1/64", At the center of the rod the differential rod worth is about 0,024% Ak/k per inch, Therefore, this contribution to the error in each reading is + 0o0004% Ak/k. The largest contribution to the error is due to the uncertainty in the rod calibration and the uncertainty in the value of peff' This is estimated to be about 5% of the reactivity change being measured. I, Results of Early Experiments It is of interest to go back to the early void coefficient experiments on the FNR(55) which sparked the present investigation, to compare with the results of Figure 180 Two types of experiments were performed at that time, In the first all the water in a fuel element was pumped out and the reactivity change measured. In the second. air was bubbled through the core (unplugged) just as was done here. The difference was that a Wet Test Meter was used to measure the air flow and the bubble rise velocity was assumed to be constant at 1,1 feet per second. The resultant reactivity curve marked "uncorrected" for lattice position 35 is shown in Figure 204 Using bubble velocity versus flow rate data obtained. in this experiment this data was corrected producing the second curve, The dotted line is taken from the (unplugged) curve of Figure 18, reproduced here for comparison and corrected

.42for the difference in rfel loading, It would seem that the Wet Test Meter produces the remaining discrepancy,

-43Configuration (normalized for 0.14 configuration change 0.2 X/ by the relative masses of fuel 0.12.- I 0.10 ____- / _ w 0.08 Velocity _ _ —-- Wc Correctido I D4/ ( Uncorrected 0.06 0.02 -- Figure 20. Experimental Effect of Voids on Reactivity (From MMPP- 110-1) for Lattice Position 35

III. THEORETICAL ANALYSIS A. Introduction The purpose of the theoretical analysis is the establish the simplest means of accurately calculating the effects of localized or distributed voids. It is important that the method selected be readily separable into components, each of which considers a particular effect, such as neutron streaming, or variation of microscope cross sections. This allows modification of one part independent of the others, whenever necessary. It shows which features of the effect are small enough to be neglected in further calculations on the same system. It also gives a clear picture of the phenomena which make up the total effect, suggesting design modifications in some cases. For example, if the void coefficient were strongly affected by the axial reflector characteristics, it is conceivable that in a boiling reactor the disturbance of the reflector over the core due to the rising steam bubbles would strongly affect its stability. Redesigning might eliminate this effect. In order to limit the area of study somewhat, an attempt was made to employ two-group neutron diffusion theory as much as possible. This theory most readily permits the use of measured material constants. Naturally, the resulting implications of the theoretical approach will be limited to reactors similar to the FNR. For example, reactors with large resonance capture cross sections may require more energy groups for adequate description. Except for this limitation however, a theory which works well for this small, hydrogenous system should be adequate for many other reactors. -44

_45B. Methods 1. Standard Group-Diffusion The group-diffusion method, as applied to two neutron groups (42) involves the solution of the following equations (2) 2 Xlv^fl' l + V''f2O2 -V[DlV1] + [D1BZ1 + l + ZR ]0 =.1.. _ V..2.2 -V[D2V2] + [D2BZ2 + a2102 = ERll + X2 l + 2 (17) where the nuclear parameters may be discontinuously varying functions, but must be constant over some region. The subscript 1 refers to the epithermal, or fast group. The subscript 2 refers to the thermal groupo When more than two regions are present, the analytical solution of the above equations by manual computation grows tedious, involving the solution of a determinant of order 2n where n is the number of regions. The solution of these equations by difference techniques has been programmed on the IBM 704 computer. Three of these programs have been used here. The brief descriptions of these programs is included here to indicate their general properties and advantages. Complete descriptions of these codes are available in the literature. WANDA(43) is a one-dimensional code. It will solve the neutron diffusion problem in an infinite cylinder, sphere, or semi-infinite slab. The eigenvalue X is determined for a finite cylinder by specify2 2 ing a finite value of B1 and B2 o These are the transverse, or axial bucklings. The terms in which they appear define the axial leakage. As will be shown later the exact value of the transverse bucklings is a crucial point,

-46PDQ(42) is a two-dimensional version of WANDA. It can solve the problem in X-Y geometry given one value of the transverse buckling (only one value is permitted) Thus any cross sectional configuration which can be defined in rectangular coordinates is allowable. This of course is a closer approximation to the actual physical shape of the reactor than the cylindrical modelo It can also operate in R-Z geometry, assuming the shape of the core is a finite cylinder. In this model, the detail of the radial. flux is lost but it is not necessary to assume transverse bucklingo MAGNUM(56) is a 32 group code. It contains its own cross section. library. Coupling between energy groups uses the SelengutGoertzel approximation. This code has several different features but is used here mainly for its buckling iteration feature(57) which is discussed in more detail. in a later sectiono This code is one dimensionalo 2o Group Perturbation Perturbation theory,(58) valid for small changes, gives the reactivity effect in terms of the changes in neutron balance caused by the nuclear perturbation. Each effect - fast leakage, thermal leakage, absorption, etco - is calculated from the unperturbed fluxes, weighted by an importance function and integrated over the perturbed regiono This method is particularly valuable in that it separates the various factors which make up the total effect. It also allows estimation, by simple computation, of many types of changes using only one set of fluxes and importance functions. The importance functions are the

-47adjoint functions defined by the transpose of the flux operator. The adjoint functions can then be found by applying normal procedure as followed in the flux calculation by using the substitutions of Table II. TABLE II PARAMETER SUBSTITUTIONS FOR ADJOINT CALCULATION Flux Adjoint 02 f D1 D2 v R y R D2 D1 x1i v2f/(v2j + vZj) X2 VZ^ /(v Uf + V F f) vE, X-2 (v 2 + V ) va v a R Ia + aZ R 2 2 1 2 2 2 B22 B1 3. Buckling Calculation In the use of one or two dimensional methods of solving three dimensional reactor physics problems it is usually necessary to supply a term which fixes the leakage in the otherwise neglected dimension.

-48This term is the transverse buckling B2 included in Equations (1). The relationship between the bucklings and criticality is given in two-group theory by: k eff - k (18) eff 2 2 ( 2 [1 + L2(Bz + B2)] + ( B1) where the B-r are the group bucklings in the plane in which the fluxes are calculated by Equations (17)~ One way to arrive at a value of the transverse buckling is by an iterative process. Assume a first value of Bz2. Put it into Equations (17) solving for X, which for values close to unity is equal to keff. From Equation (18), calculate the Br2. Now turn the geometry over such that the original transverse direction is the one examined in detail by Equations (17) and the new values of the transverse bucklings are the Br2 just calculated This then, gives a new value of X, or keffo This procedure may be repeated until convergence of keff is found to be within tolerable limits. This method has been used(33) in the calculation of void coefficients, but with only one pass, i eo no iterating procedureo The Buckling Iteration Program which is part of the MAGNUM code follows a procedure similar to this. The difference is that it 2 does not use Equation (18) to calculate the Bri It integrates the leakage flux over the surface to determine these factors. As part of the output, this code prints the results of each iteration, showing how rapidly the void coefficient convergeso

-49C. FNR Reactor Physical Properties 1. Standard Fuel Element The standard fuel element is made up of 18 curved plates o.o6o inches thick, 24-5/8 inches long brazed into two side plates, 0.188 inches thick to form a box as shown in Figure 7. Each curved plate is actually a sandwich of two 0.020 inch thick 2S aluminum plates around one 0.020 inch thick section of U-A1 alloy. This latter section is usually referred to as the "meat" while the outer slabs are termed the cladding. The Uranium in the meat section does not extend the full length of the plate, but only to within about 1/2 inch of each edge. Each full fuel element contains 140 grams of U-235 and 16 grams of U-238, 2, Control Elements The four elements through which the control and safety rods pass are similar to the standard element but have 13 interior plates removed. Two 1/2 inch thick aluminum plates are substituted for structural purposes 3. Reflector Elements The graphite reflector consists of aluminum cans, of the same size and shape as fuel elements, filled with graphite and sealed. 4. Grid Plate The grid plate into which the fuel elements are placed is a 5-1/2 inch thick aluminum slab. The center to center spacing of the grid plate holes is 3.189" in a direction parallel to the side plates of the elements and 35035" in the other direction,

-50Table III lists the average contents of each fuel element (between the two planes defined by the top and bottom of the fuel region) TABLE III COMPONENTS OF STANDARD FNR FUEL ELEMENT Component Volume (cm3) Per Cent Water 2190 48.45 Aluminum 1549 41 33 U-235 7.407 0.1977 U-238 o08466 0o02259 TOTAL 3747o 2556 100o00029 D. Nuclear Constants 1. Basic Data The two-group nuclear parameters are generated first as a function of densityo Later, the effects of changes in neutron temperature, thermal disadvantage factor and streaming due to the finite bubble size are calculatedo Two major sources of information were used for the thermal microscopic cross sections used here Thermal cross sections averaged over a Wigner-Wilkins spectrum (SOFOCATE code) have been calculated for absorption and fission of U-235, scattering of hydrogen and absorption of a unit 1/v - type absorber~o' It is assumed here that aluminum * A unit 1/v - type absorber is a material whose absorption cross section, varying inversely with neutron velocity, is one barn at v = 2200 m/seco

-51and water are 1/v type absorbers. A list of 2200 m/sec cross sections tabulated by Argonne National Laboratory(44) is used to fill out the remaining necessary data. A plot of the SOFOCATE data as a function of the ratio of U-235 to hydrogen is given in Figure 21. This data assumes that the number of 1/v type barns (absorption) per hydrogen atom is very small. In the case of the FNR, this ratio is about 0.15 barns per hydrogen atom. The correction to be applied due to the amount of absorber is about 1% and is neglected. The scattering cross section for water is calculated using both SOFOCATE and ANL tabulated data. The SOFOCATE code gives as(H) which varies strongly with the neutron spectrum. It is assumed that the scattering cross section of oxygen and the chemical binding effect do not change significantly with energy in the thermal range. Thus, from the ANL tabulation, at 2200 m/sec as(H) = 45.28 b as(o) = 4.2 b oa(H20) = 103 b assuming that the molecular cross section can be represented by as(H2O) = [2as(H) + aS(0)] 7, (19) where 7 is a constant accounting for molecular binding, we find that 7 = 1.087. Now using os(H) given by SOFOCATE, we get the effective thermal cross section sg(H20) = 96.63 b. The thermal data for U-238 is taken from Weinberg and Wigner.(46

.1:.,~3 0.9 48 0.8 620 - 46 (I/V) 0.5 -- 5 00 -- 401'.0.4 A 460 -— 38 -- -_ 0.3 4 —20 -— 36 2 — 380 -— 34 O0.1 -- 340 32 o 01 300' 30 0.002.004.006.008.010.012.014.016.018.020 25 H N /N Figure 21. Variation of Thermal-Averaged Cross Sections with Ratio of U-2355 to H Atoms (from WAPD-185)

-53The epithermal parameter D1 for the A1-H20 mixture is given by Murray,(45) based on calculations at ORNL. The reported results are: D1 = 1.143 + 0.153 o (20) where (X = volume of aluminum/volume of water. The relationship used for T is based on results of the measurement of T in aluminum-water mixtures by Hill et al.(5) The following equation was derived from theoretical considerations (see Appendix B) and adjusted to fit the experiments: 100(1 + 2 (21) T = - (21).0207 a2 + 1.27 a + 5.72 When some of the water is displaced by air (assumed to have the properties of a void), the age of the mixture (for small voids) will vary as: Vtotal 2 T ix T [ total ] (22) mix 0Vtotal Vair where To is the age without the air present but for the same volume ratio of materials. Thus the previous formula for T is modified to: 1.00 (1 + )2 + 2 T ~~= (^P- ( —-- (25).0207() + 1.27 - + 35.72 a + p P P where p = water volume/(water + air) volumes. A similar correction to D1 gives: D] = [1.143 + 0.153 5 ] o (24) From these two numbers, Z1, the removal cross sections is found by = D1/7

-54The data for the long-lived fission product poison average thermal cross section is taken from Deutsch.(51) The FNR had been operated for approximately 15 Mw-days before the experiments were performed. The burnup (at the rate of 1.07 gm/Mw-days) is 0.55%. The P -1 average poison cross section ZP is 0.00107 cm using Deutsch's approximation. Calculations have been made of the thermal flux depression in fuel elements identical to thos in the FNR, by Reynolds(48) using the P5 approximation to the one-velocity transport equation. The average flux in the fuel region is 0.9877 and the average flux in the aluminum is 0.9887 (relative to the average flux in the water channels). Table IV shows the calculated properties of the core as a function of p. Table V shows the properties of the reflectors. 2. Infinite Multiplication Factor Calculations by Trubey and Lessig(52) have shown that for fully enriched, water moderated reactors the assumption that ep = 1 is excellent for a uranium density Nu < 2 x 1020 atoms/cc. In the present case as Nu = 0.957 x 1020 this assumption is valid. The cancelling of the resonance escape probability and the fast fission factor is valid in this two group approach since the gain and loss of neutrons occur in the same groupo The variation of k,(= v4f/Z ) with water density is shown in Table IV. The world consistent v = 2.47 is used as given by Weinberg and Wigner,(46)

-55E. Second Order Reactivity Effects The only density effects included in the cross section table are due to the changes in the macroscopic absorption and scattering cross sections of the moderator. It would be convenient indeed if these were the only effects to be considered, since these are easily calculated. Therefore the other, more subtle effects, are investigated individually, below, to determine whether they can be treated by approximation or neglected entirely. 1. Thermal Disadvantage Factor The effective multiplication factor keff is directly proportional to the thermal utilization f. This latter term is the ratio of thermal neutron absorption in fuel to total thermal neutron absorption. In a homogeneous reactor this is simply the ratio of fuel to total cross sections. In a heterogeneous reactor, each cross section must be weighted by the average thermal neutron flux in the region. Thus: f a. (25) where U is the fuel region M is the moderator region. In terms of the thermal disadvantage factor F =, (26) afJ F' a(27) a a

Now, given a change in F and Za the change in reactivity will be given by: k =df = (f-l) (dF + (28) k f F TABLE IV FNR CORE NUCLEAR PROPERTIES (a = 0.711) U -1 2 cm 0.04982 0.04982 0.04982 0.04982 0.04982 Zl,cm1 0.00451 0.00451 0.00451 0.00451 0.00451', cm-l 0.01036 0.00984 0.00953 0.00932 0.00881 Z, cm 1 0.00106 0.00106 0.00106 0.00106 0.00106 a,otalcm-l 0.06575 0.06523 0.06497 0.06471 0.06420 Zf, cm 1 0.04249 0.04249 0.04249 0.04249 0.04249 T, cm2 63.19 69.28 73.36 76.29 84.47 D1, cm 1.247 1.292 1.519 1.338 1.390 ZR, cm-1 0.01973 0.01865 0.01798 0.01754 0.01646 2r-,cm- 0.03489 0.03489 0.03489 0.03489 0.03489 Y,cmnl 1.2770 1.2132 1.1748 1.1493 1.0855 tr' otalcm-l 1.5119 1.2481 1.2097 1.1842 1.1204 D2. cm 0.2541 0.2671 0.2755 0.2815 0.2975 L2, cm2 3.865 4.095 4.244 4.350 4.634 Zf/Ia 0.6462 0.6514 0.6545 0.6566 0.6618 v4f/1a 1.5961 1.6090 1.6166 1.6218 1.6346

-57TABLE V FNR REFLECTOR NUCLEAR PROPERTIES Material Graphite (p=l.6) Water (p-l) 2, cm'1 0.00026 0.022 D2, cm 0.778 0.164 L2, cm2 2959 7.453 Dl, cm 1.470 1.145 2R' cm-1 0.00320.0423 T, cm2 460 27 As moderator is removed, the reactivity increases due to the last term. The removal of absorber in the moderator region also reduces F, also increasing the reactivity. The magnitude of this latter effect was calculated by the author using the P3 approximation to the one-velocity transport equation. Reynold's calculation by P3 and P5 show very little is to be gained by use of the higher approximation. To determine the order of magnitude of the effect the aluminum region was removed, simplifying the calculations further. The results in Table VI are for a two region assembly of U-A1 alloy and water having the same ratio H/U and the same plate spacing. Clearly, due to the close spacing and small fuel concentrations, the effect of changes in F for this reactor are always negligible. The void coefficient effect of this term is pF (f-l) - (V.C.)F = (f- 0.0017 (29) X

-58TABLE VI RELATIVE EFFECT OF VARIATION IN F ON REACTIVITY p 1 0.85 LaM o0.0195 0.0166 ~F ~ ~ o 0.618 o.618 F 1.020 1.019 A ZM/M -o.149 A F/F -0.001 2. Neutron Temperature The removal of moderating material tends to shift the thermal neutron energy distribution to a higher average kinetic energy ("hardening"). This shift has several effects. Since the uranium fission and absorption cross sections are not exactly "one-over-v" type and do not even change in the same proportions, there is a slight change in'. Since all the microscopic absorption cross sections vary with energy the thermal utilization is also affected. Finally, the change in the microscopic scattering cross section of hydrogen and the change in the microscopic absorption cross sections combine to increase the thermal diffusion length. All these effects are in addition to the direct, or macroscopic effects of the density change as previously discussed. The epithermal microscopic parameters are affected in this manner only by the slight shift of the "cut-off" energy, the energy at which the 1/E portion of the spectrum is assumed to intersect the Maxwell-Boltzmann portion of the spectrum.

-59From the SOFOCATE data of Figure 21 the average thermal microscopic cross sections are given as a function of the atomic ratio U/H. For the initial or unperturbed condition of water density p - 1 this ratio is 0.00245. Using this information, Figure 22 was constructed from Figure 21, to give the cross section variation as a function of p. Also included is the ratio af25/aa25 which is of course, T'/v. Examination of these curves shows that for p > 0.75 the assumption of linearity is quite good. It should be remembered that since these cross sections are calculated for an infinite medium, these curves are not to be trusted for low densities in a small, reflected reactor. The actual cross sections will be higher at low densities. In any case it must also be kept in mind that the use of the reported measurements of T limits us to p > 0.7 as this was the limit of measurement. The three temperature effects will now be examined separately. (a) Temperature effect on r. The linear approximation to the curve of n/v in Figure 22 is I/v = 0.85378 + 0.00147 p. (30) Therefore, d,/ = [ o.oo147 (51) 0.85378 P --..... + 0.00147 The bracketed number is the void coefficient effect and is given very closely by (V.C.) = -0.00172 p.

0.9 -b 0S85 0.8 -620 - 44 "' 0.853 0.7 -580 42 _0. 851 s _ _. _ _ _ _-__,_ _ f a 0.6 - 540 - 40 _ 1 0.849 0.5 -— 500 38 0.4 460 36 IL- 0.3 4 20 - 34 -H —0.2 - 380 32 \ k r 0.1 -340 30 0 300 -- 28...., I 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 WATER DENSITY, p Figure 22. Variation of thermal-averaged cross sections in FNR with water density (from Figure 21)

-61(b) Temperature effect of f. Looking again at Figure 22 we see that the linear approximations aa(v) = 0.725 + 0.071 p (32) and a5 = 469 + 61 p (33) are reasonable over the same range as before. Using these relationships and assuming water and aluminum to be 1/v materials, it turns out that df,-0.00004 (l-f) dx (34) f x or that (V.C.)f ) -0.00001 (c) Temperature effect on thermal leakage. The thermal non-leakage probability Pth is approximated (in the "single buckling" method) by Pth-1LB r(35) Pth - 1 + L2B2 where L2 = D/4. The reactivity effect of variations in Pth is given by dPth L2B2 dD da _ _[ _ - (36) Pth 1 + L2B2 D F Using the same linear approximations as before, we see that dZa 0.00132 p + 0.00584 dp ( Za 0.0148 p + 0.0507 p

-62Assuming that only the hydrogen scattering cross section varies, as shown in Figure 22: as = - 39.15 + 3.20 p (38) we find that r, = 0,034 + (1.184 + 0.0918 pt) Pd (39) where pt is meant as the variable and Pd is the macroscopic contribution to the cross section and is fixed in this portion of the analysis. From this it is simple to show that dD 0 _ 0918 p dp (40) D 0.034 P - + (1.184 + 0.0918 p) making the overall reactivity effect dPth L2B2 0.00132 p + 0.00584 Pth 1 + L2B2 0.0148 p + 0.0507 + [ i.1 7 dp (4_1) 1 + 129 + 371 p P P when p = 1, this becomes dPth L2B2 (0.1794) d (42) Pth 1 + L2B2 p For L2 = 3.865 and B2 = -0.002. (v.c.)L2 o0.00144 The value of the thermal buckling used here was a result of a buckling iteration calculation using MAGNUM. Had the geometric buckling been used here this term would have been five times larger and of the opposite sign.

-63The sum of these three neutron temperature - microscopic cross section effects is then (v.c.)T = -.00029 and is seen to be of the opposite sign from the previously discussed disadvantage factor effect. 3. Neutron Streaming Effect When voids are distributed throughout a diffusing medium, one of the major effects, as we have seen, is the increase in the mean free path of neutrons in the system. If the holes are small and closely spaced, the effective mean free path is inversely proportional to the average density of the system. If the holes are large compared to a mean free path, the arithmetic mean free path is inversely proportional to the average density, but the root-mean-square free path is a function of the size and shape of the voids. It is this latter quantity which is of interest in describing the random interactions which occur in the reactor. The description of this "streaming" effect was originated and subsequently amplified by Behrens(39) in connection with the coolant channels of gas-cooled reactors. For a uniform distribution of spherical holes in a homogeneous medium the effect on the mean square slowing down length is(39) T 1 + 2 + i [coth r + 0.125] (43) where 0 = volume ratio of holes to material of the reactor, X = mean free path in the material of the reactor, r = hydraulic radius of the holes (= 2 radius), To= the age in the material without the holeso

-64This correction can be applied to the diffusion area as well by using the thermal value of k. However, the major effect in this FNR reactor is the fast leakage. Due to the relative uncertainty in this correction factor, the thermal effect is neglected. It is observed that as r 0/x diminishes the correction becomes (1 + ~0) which is simply the inverse density effect already accounted for. Therefore, the only correction to be applied here is the ratio T/TQ from the above equation to (1 + )2. In heterogeneous reactors such as the FNR the bubbles, or voids are in the water portion of the moderator only. Therefore it is not correct to modify the total age directly. The total age can be separated into the water and aluminum components by rewriting the equation in the Appendix as 1 1 [2 + p2 4.579 (44) [-+ -_+ 45. (44) T (1+ca)2 TA1 Tw a1 w a2 P2 The effect of the Behrens correction to 1/T is then found by multiplying the value Tw by the ratio 1 + 20 + r [coth + 0.125] RB =........k -....... —' (45) (1 + 0)2 The relationship between 0 and p is: = 1/p - 1 (46) It is extremely difficult to determine the average bubble radius in the core. However, it is certainly less than half of the plate spacing. Thus the parameter r/k is less than 0.2. The effect

on the void coefficient is shown for two different values of r/X, as a function of p in Table VII. The fact that the effect decreases with increasing void content (at constant bubble size) indicates that the ratio RB does not increase as fast as (1 - p). In fact it reaches a plateau. TABLE VII NEUTRON STREAMING EFFECT ON THE VOID COEFFICIENT X r/X RB AT/T () (V.C.) stream.02 0.1 1.0018 0.136 -0.026 0.2 1.0049 o.469 -0.088.05 0.1 1.0031 0.286 -0.021 0.2 1.0081 0.701 -0.053.08 0.1 1.0036 0.365 -0,017 0.2 1.0105 0.892 -0o042.10 0.1 1.0057 0.390 -0.015 0.2 1.0113 1.013 -00038.15 0.1 1.0039 0.349 -0.009 0.2 1.0124 1.108 -0.028 4. Transverse Buckling Considering the shape of this reactor it is apparent that the spacial flux variations do not lead one to a convenient, simplified model such as a sphere or slab. This then must be treated as a decidedly

-66three-dimensional reactor. To solve the diffusion equations in this geometry it would be necessary to use a three dimensional code. Although this type of code exists, it requires. a larger computer than is presently available here. Even on machines of sufficient capacity it is limited to a relatively small number of space points and takes several hours to complete. Therefore, it was necessary to make use of the pseudo-one or two dimensional codes. Several methods were studied. In the one and two dimensional codes, the leakage in the dimensions not explicitly examined is taken into account by inserting a transverse buckling. In most applications this parameter is a quantity measured by extrapolation of neutron flux traverses. In the calculation of reactivity coefficients however, a most important feature is the variation of the transverse buckling as a specific parameter change is made. An estimate of this variation can be made by application of one group theory. For an infinitely reflected reactor the reflector savings 5 is given by Dc Bc Lr tan B6 = c Lr (47) c Dr where the subscripts refer to the core and reflector properties. Where the core is large enough so that tan Bc8 a BcS (48) we see that 8 - Lr. (49) Dr

-67To apply this to the fast group, in which we are primarily interested, the value fTr is used for Lro The buckling is then determined as that of a bare slab with dimensions (H + 2 6) where H is the actual core length: Bz 2 = ( 2 (50) The resulting effect on the reactivity is approximately: k +TB2 B2 (51) k 1 + TB2 B2 The results of this method are shown in Table VIII as a function of p. TABLE VIII THE TRANSVERSE BUCKLING VARIATION AND ITS EFFECT ON VOID COEFFICIENT, BY ONE GROUP APPROXIMATION P Do Bz2 k/k(%) V.C. 1 1.247 o001873 - o0164.97 1.273 oo001861 00326 o0163.95 1.292 o001850 0.0805 oO161.92 1.319 o001839 0.128 o0160.90 1.338.001831 0o158,0158.85 1.390.001802 0,229 o0153 F. Results of Analysis 1. Comparison of WANDA Perturbation, WANDA Criticality and PDQ(XY) Criticality Methods In order to establish the relative accuracy of three approaches to the calculation of the void coefficient, one problem was solved by three methods, using the same nuclear datao The data used in this portion of the study varied somewhat from those given in Tables IV and V

-68as they were made before those Tables were completed. The nuclear data used here are reproduced in Table IX and X. The transverse buckling was assumed in each case to be Bz2 = 0.0034106. This value was chosen to make the initial value of keff for the WANDA Criticality Perturbation methods unity. The geometry assumed for the WANDA problems is shown in Figure 23 (a). The water reflector thickness, shown as 8 cm was changed in later calculations to 16 cm. The dimensions were chosen to give roughly the same geometric buckling as the square configuration used in PDQ, shown in Figure 23 (b). The flux is set to zero on the outer boundaries and is symmetrical about the center lines. TABLE IX NUCLEAR DATA FOR CORE USED IN COMPARISON OF ANALYTIC METHODS P 1.97.95.92.90.85 Dl 1.247 1.273 1.292 1.319 1.338 1.390 Y^1.001.001.001.001.001.001 ZR.02028.01971.01931.01871.01850.01755 vfl.00227.00227.00227.00227.00227.00227 D2.2222.2284.2335.2411.2474.2606 2 o.0811.0807.0805.0801.0798.0792 V4f2.1243.1343.1343.41343.1343

-69F =FUEL G = GRAPHITE W= WATER /~)~ \ (ALL DIMENSIONS ARE IN Cm.) A. sb F\ G \w _ V _, _,w Figure 23a. Configuration Used in One Dimensional Calculations F = FUEL 30 G =GRAPHITE W W=WATER (ALL DIMENSIONS ARE IN Cm) 1524 F --—.o~.. W 5.72.1 9 7621 - 430 kI —19.O-5 162 -- 30 Figure 23b. Configuration Used in Two Dimensional Calculations

-70TABLE X NUCLEAR DATA FOR REFLECTOR USED IN COMPARISON OF ANALYTIC METHODS Water Graphite D1 1.143 1.470 Zai.000001.00001 >R...o 0341.0032 D2.e164.778,a2.021.00036 The one-dimensional, two-group fluxes and adjoints, as determined by the WANDA program are shown in Figure 24 (for p = 1). Using these functions the reactivity changes due to the water density changes were determined by first-order perturbation theory. The portions of variable density in the one dimensional approximation are first, the inner cylinder of radius 6 cm., second the outer annular cylinder between r = 16 and R 20 and thirdly, the entire fuel bearing section of the core. In the PDQ problems the central and corner elements as well as the entire fuel bearing section were voided in turn. The results of these calculations are plotted, in Figure 25 in terms of reactivity for a void fraction assuming, for the localized void that it is inserted throughout a fraction of the core equal to 1/21 (approximately one fuel element)~ The values for the entire core are reduced by that fraction to compare relative magnitudes.

r10 I I 9 #j|th THERMAL FLUX l.l EPITHERMAL FLUX 8 / | th -THERMAL ADJOINT s7 1 IFUEL E I,f, EPITHERMAL ADJOINT 7_ ________ ^ ____tFUEL Q 6 > 4 82 RADIUS, Cm.GRAPHITEWATER Figure 2. One Dimensional, Two-Group Fluxes and Adjoints 2 5 12 —- — 0 24-2 —- --- -- --- -- 4 8 12 16 20 24 28 32 36 40 44 RADIUS, Cm. Figure 24. One Dimensional, Two-Group Fluxes and Adjoints

-72- --- BY ONE DIMENSIONAL PERTURBATION -A —.. N CRITICALITY (WANDA) 0.4 El " TWO. " (PDO) (ALL ASSUMING Bg.0.0034106) EXPERIMENTAL RESULTS 0.3/ / >O.2 / 0. / I I E /E T LEMENT 0 _______ 2 4 6 8 10 12 14 16 AIR VOLUME CONCENTRATION, % AV/V Figure 25. Calculated Reactivity Effects of Void

The comparison shows that first order perturbation theory is excellent for this type of perturbation, especially if the voids are inserted uniformly over the coreo The error as compared to the criticality method varies from +12,5% in the center to +10% at the edge of the core, to -2% for the entire core (for 10% void). Obviously the good result for the entire core is due to the slowly changing fluxes and adjoints. The larger errors in the cases of localized voiding are due to the considerably local changes in flux and adjointo A more surprising result is found in the comparison of one and two dimensional calculations. The difference between the results is less than 6% in every case studied, for both localized and uniformly distributed voidso In comparing the overall value of one method over another the effort and expense involved in each should be included as factors. The WANDA perturbation scheme requires the running of two problem sets, one each for fluxes and adjointso These each take about 30 to 45 seconds for the number of space points used here (45). The resulting curves must then be manipulated and integrated to obtain the integrals of the various products of the fluxes and adjoints and their derivatives over the regions of interest. Manually this is a tedious and time consuming job - of the order of twenty five hours. If this is to be done repeatedly, the entire procedure can be programmed on the digital computero For any one reactor however, this need only be done once. The integrals can then be used for any combination of perturbations with relative easeo To generate the 13 points included in Figure 25, about thirty hours of

-74manual calculating (assuming no errors are made) plus two minutes of IBM 704 time and about fifteen minutes of data card punching time are required. The WANDA criticality results require no additional analyses. The fifteen points included in Figure 25 required the running of fifteen problems. The total time is about six minutes on the IBM 704 and thirty minutes of data card punching. The PDQ program has some odd features. The length of time a problem takes is a function the number of points in the mesh (roughly the square of the number of WANDA points to get the same resolution) and the amount of final deviation from criticality. A strongly subcritical core can iterate for hours without approaching the final result. The seven problems run in this study each required between thirty minutes and one hour to complete on the IBM 704. It is possible to do a two-dimensional perturbation calculation, using PDQ in a manner similar to the WANDA perturbation technique. This introduces the problem of calculating surface integrals of fluxadjoint products and gradients. For best results this should be computerprogrammed. Based on the results of this comparison, however, it is doubtful that sufficient accuracy would be gained by this complication. The obvious conclusion is that the WANDA criticality technique is best unless data on a particular interior, off center element are desired. If no digital computer were available the two-group perturbation technique using manually calculated fluxes and adjoints (as in Spinrad and Kurath(53)) could be employed.

-75Two more notes may be added at this point. First, the fluxes and adjoints in the water channels have been shown(59) to be extremely inadequate when determined by the method used here. Therefore the prediction of void coefficients in these channels is not expected to be accurate. Second, although the three techniques compare very well with each other, they obviously do not agree well with the experimental data, which are also plotted on Figure 25 for reference. The error decreases with the radial position, varying from +58% at the center to +18% at the edge. 2. Transverse Buckling Calculation An examination of the effect of the transverse buckling value on the void coefficient was made. The problem is not only to show the effects of various assumed bucklings on the reactivity but to determine the best value of the buckling. The method selected was the Buckling Iteration Program. The output of this program is a series of 32 transverse bucklings, one for each group. Also tabulated by group are the values of D, A, Zt and NH. These results were used to determine the two group transverse bucklings by the following relationships: 52 (1 + B2 Ti) = (1 + Bi AT) i=4 AT = UiDi (52) i Zc + EH B2 was found by trial and error in each case. The thermal buckling is given directly in the program output. These bucklings were then included as input to WANDA calculations. The results of the BIP calculations are shown in Table XI and Figure 26. These were calculated from the final pass results.

-76E' 0.0025 VOID IN REFLECTOR AS WELL AS CORE E 0.0024.. I 0.0023 -- -0.0036 E -0.0038 a: -0.0044 10.0044 I-0.0046 0 5 10 15 20 AIR VOLUME CONCENTRATION,% Av/v Figure 26. Variation of Transverse Bucklings with Void Concentration (as calculated by the buckling interaction subroutine of MAGNUM, assuming the void to be uniformly distributed).

-77TABLE XI TRANSVERSE BUCKLINGS CALCULATED BY THE BUCKLING ITERATION PROCEDURE P B2 (cm2) B2(cm-2) epi th 1.002491 -.003599.95*.002449 -.003676.90.002418 -.003995.85*.002366 -.003825.80.002346 -.004464 There are two limitations in this code which restricted its use somewhat, for this application. First, the number of regions is limited to five. Second, the buckling iteration can be done only on the inner region. To calculate the effect of voiding the entire core, it was therefore necessary to eliminate the water annulus. The geometry was identical to the WANDA geometry in other respects. The effect of voids in the axial reflector was examined by calculating the bucklings with and without voids present during the axial iteration. The resultant effect on bucklings, shown in Figure 26 is, due to symmetry conditions, for voids above and below the core. At best this is only a guess of the effect of the actual voids in the core since the water density in the reflector (in the air bubble experiments) is less than in the core. The bucklings determined by this method were then used as input to a second series of WANDA calculations, in two forms. At first, * Void assumed to extend into the axial reflector region.

-78the bucklings were assumed constant at the values of p = 1. Then the bucklings in the perturbed regions were assumed to vary as shown in Table XI, while the bucklings in the other regions were assumed to remain constant. The results are given in Table XII. TABLE XII RESULTS OF WANDA SERIES 2 FOR CONSTANT AND VARIABLE BUCKLINGS Procedure p Void Region Buckling Assumption Void Coefficient.95 Central Constant - 442.90 Central Constant -.456.85 Central Constant -o466.95 Outer Constant -.108.90 Outer Constant -.109.85 Outer Constant -o113.95 Central (& Reflector) Variable -o312.90 Central Variable - 353.85 Central (& Reflector) Variable -o381.95 Outer (& Reflector) Variable -.0735.90 Outer Variable -o0766.85 Outer (& Reflector) Variable -.0870.90 Entire Core Variable - 190.85 Entire Core (& Reflo) Variable -.185 At the same time, calculations were made with constant buckling approximation assuming several values of buckling; with the further 2 2 assumption that Bth Bpi. These results are tabulated, along with the similar number from WANDA series 1, in Table XIII.

-79TABLE XIII CALCULATED VOID COEFFICIENTS ASSUMING CONSTANT BUCKLINGS (p = O.9, CENTRAL VOID) 2 2 Case B B Void Coefficient epi th 1.002000.002000 -.476 2.002491 -.003599 -.456 5.002800.002800 -.567 4.0034106.0034106 -.538 2 2 It is likely that had Bth Bepi in case 2, the resulting void coefficient would have been increased to somewhere between those of cases 1 and 3, giving a somewhat linear relationship with buckling. 2 This would indicate that the result of employing the calculated B rather than assuming it to be equal to B is to reduce the void coepi efficient by approximately 10%. One more interesting result was obtained from the Buckling Iteration Procedure. The program which calculates the bucklings also calculates keff during each pass. It uses as criteria for convergence, both the bucklings and k. The actual values of the void coefficient eff' as calculated are not too reliable due to the limitations previously mentioned which require slight liberties to be taken with the geometry. However, it is of interest to observe how the void coefficient as calculated from the keff values in each iteration vary. In the five problems run, the keff values in each iteration were used with the values in corresponding iterations to find the resultant void coefficient. In other words, this gives the void coefficient as calculated if the machine

-80had stopped after that iteration in each problem. The normalized results in each case were the sameo These are shown in Figure 27, The indication is that at least two complete sets of Radial-Axial iterations were necessary to give accurate results,

-811.8,- -- -- -- j - 1.7 l - — INITIAL BUCKLING GUESSES: 2 2 \1.6L l l l | l [ Bth = Bepi =0.00238 1.6 -o- RADIAL SOLUTION 1.5 1.'85 t ~-a-D AXIAL SOLUTION 1,4 Wt 1. L-.J in F.3 1.2'Li 1.1 -I I.0 9.0 0.8 0.7 I 2 0 2 3 4 5 6 7 8 9 10 ITERATION NUMBER Figure 27. Variation of Void Coefficient During Buckling Interation Subroutine of MAGNUM. (Normalized to zeroth iteration). Final Bucklings as shown in Figure 26.

IV. SUMMARY OF RESULTS A. Localized Voids The experimental results show linearity of the reactivity effect with void, up tO the 11 void introduced, indicating a constant, negative value of the void coefficient. The magnitude of the coefficient varies monotonically with radial distance from -.316 at the center of the core to -.120 at the corner. The calculated coefficients vary from -.45 to -.11 with constant buckling and -.35 to -.08 assuming a spacially uniform buckling change. The calculations indicate a coefficient of increasing magnitude as void volume increases. The experimental coefficient due to magnesium insertion (one point) is -.284. Calculations using perturbation theory indicate that the magnitude of the (uniformly distributed) magnesium coefficient is about 6% lower than the void coefficient in this reactor. The data for all but the central element show a tendency for an increasing void coefficient. However more data points would be'required to show whether these are true effects rather than the result of errors. B. Distributed Voids No actual measurement of the effect of voids uniformly distributed over the core was undertaken. Using air bubbles this would be an extremely complex as well as unsafe experiment. Due to the small mean free path in the core, with the low void concentrations used here, it is reasonable to assume that interaction effects between elements is small. Therefore, the addition of the measured localized results is used as the -82

-835 distributed coefficient. The resultant coefficient is -.190. The calculated void coefficient using the spacially uniform buckling variation is -.190.

V. CONCLUSIONS A. Experimental Technique The feasibility of the use of a stream of air bubbles as a void simulator has been demonstratedo Furthermore, in pool type reactors having plate spacingssimilar to the FNR, flow data generated in this study can be used, In such reactors this technique would actually be the easiest to carry outo If distributed void coefficients need to be measured directly some other method should be adopted, Either expanded plastic or magnesium sheets should be employed, If the plate spacing differs significantly it is still possible to use this data to generate a new bubble velocity versus density curve without requiring the use of a gamma ray densitometero The method is simply to compare (as a function of density) total void volumes in a section having the new spacing but the same water fraction, This can be done in an open water tank by measuring differences in total water height. The Wet Test Meter should obviously not be used for this application. This is mainly due to the high pressure required to drive reasonable flow through the small bubbling orifices. The efficiency of the Wet Test Meter used in the early experiments on the FNR began to drop rapidly when the bubble volume was increased over 2%, B, Analysis The combination of the buckling iteration process with two-group, one dimensional diffusion theory appears to be quite accurate for predicting the distributed void coefficient. Perturbation and eigenvalue -84

-85variation methods give equally accurate results. The choice of which to use should depend on the number of problems to be investigated and the availability of a program for evaluating the perturbation integrals. For void fractions up to 25% it is not necessary to take into account cross section variations due to changes in neutron temperature or thermal disadvantage factor. A most significant parameter in this reactor, however, is the transverse buckling. It is important not only to choose the correct bucklings initially, but to include the buckling changes as the water density varies. The buckling change is not well approximated by the simple relationships of reflector savings calculations using one group concepts. For calculation of localized void coefficients, the same technique does not give sufficiently accurate results. The coefficient is overestimated in the central position and underestimated on the edge, The error in both cases is approximately 25%. In each of these calculations it was assumed that the transverse buckling variation with density was the same as the variation of the entire core buckling under a distributed void condition. This assumption is equivalent to the separability of flux in the axial and radial direction. Separability implies that the solution of the flux equation V 2 (r,z) + Bg2 (rz) = 0 (53) can be found by writing 0 (r,z) = e (r) Z (z). (54) In a small, completely reflected core it is not surprising to find that this assumption is poor. A method, suggested by the results of

-86this investigation, with which to take into account the non-separability, ioe,, the dependence of transverse bucklings on radial position, is to extend the buckling iteration scheme so as to allow a, Iteration of any region in the core (rather than the central region only) and b, iteration on more than one region at a time (three region iteration would be sufficient for most problems)o Modification (a) is presently being inserted into the MAGNUM Buckling Iteration subroutine(57) This is important not only from the standpoint of non-separability but for consideration of reactors having a central non-fuel zoneo Modification (b) has not been attempted to date. The extension would involve separate treatment for each region (similar to the present one region treatment)o Regions having interior as well as exterior surfaces would include the gradient of the interior surface flux in the buckling calculationo Each pass would include one iteration for each region, The convergence criteria would then monitor the bucklings by region and the eigenvalues by region, The value of this type of program is that it would allow accurate calculation of reactivity coefficients for small, two dimensional coreso It should be less time con.suming than existing two dimensional codes since it does not require a two-dimensional, difference equation solution of flux. For example, the PDQ code has an R-Z geometry option which should give the results required here. However, attempts to use this code were unsuccessful due to the long machine running times involved, One R-Z

-87problem of 2500 mesh points iterated for over two hours before converging.* These were two-group problems. The MAGNUM 32 group code, using the (one region) Buckling Iteration subroutine required only 15 minutes to converge With the equivalent of 1440 mesh points (36 radial, 40 axial). It is not likely that the modification of multiple region buckling iteration would multiply t4e running time by much more than the number of regions undergoing iteration. A less important feature of the comparison of the analysis with the experiment concerns the neutron streaming effect. It is observed that the void coefficient as calculated neglecting streaming, increases with increasing void volume. In the experiment no consistent trend of increasing void coefficient was observed, This feature has been observed by others. Conceivably it is the neutron streaming effect which tends to flatten the curve. Using Behren's analysis it is shown (Table VII) that the effect of neutron streaming decreases, for constant bubble size, as void content increases —at least over the region examined, While this experiment could not be deemed a confirmation of Behren's analysis, the indication from this limited information is that this analysis has the correct properties, * Other, similar problems took so long that they had to be removed before convergence was complete.

REFERENCES 1. Newman, P. C. and Whelan, P. F. "Relation of Volume of Bubble to the Diameter of the Orifice at Which It Is Formed,: Nature, 169, (letter), (February 23, 1952), 326. 2. Ladyzhensky, R. M. "Study of the Movement of Air Bubbles in Water at High Re Values," Zhur. Priklad, Khim., 21, (1954), 22-32. 3. Verschoor, H. "Some Aspects of the Motion of a Swarm of Gas Bubbles Rising through a Vertical Liquid Column," Trans. Ins, of Chem. Eng. (London), 28, (1950), 52. 4. O'Brien, M. P. and Gosling, J. E. "Velocity of Large Bubbles in Vertical Tubes," Industrial and Engineering Chemistry, 27, (1935), 1436. 5. Siemes, W. "Gasblasen in flussigkeiten," Chemie-Ing-Tech, (January 26, 1954/Nr. 11), 614. 6. Miyagi, 0. "The Motion of an Air Bubble Rising Through Water," Technical Reports, Tohoku Imp. Univ., 4, n. 2, (1924), 1. 7. Barr, G. "The Air-Bubble Viscosometer," Philosophical Magazine, 7, n-l, (1925), 3950 8. Allen, H. S. "On the Motion of a Sphere in a Viscous Fluid," Philosophical Magazine, 5, n. 50, (1900), 323. 9. Greenfield, M. L, et al. Studies on Density Transients in Volume-Heated Boiling Systems, AECU-2950 (USAEC), October 1954. 10. Hooker, Ho H. and Popper, G, Fo A Gamma Ray Attenuation Method for Void Fraction Determinations in Experimental Boiling Heat Transfer Test Facilities, ANL(5766 USAEC), November 1958. 11. Egen, R, A., Dingee, D. A, and Chastain, J. W. Vapor Formation and Behavior in Boiling Heat Transfer) BMI-1163 (Navy Department, Bu Ships) February 1957. 12. Cook, W. H, Boiling Density in Vertical Rectangular Multichannel Sections with Natural Circulation, ANL-5621 (USAEC) November 1956. 13. Petrick, M. Two-Phase Air-Water Flow Phenomena, ANL-5787 (USAEC) March 1958. 14. Yanez, J. Determination of the Density of a Mixture Water-Air by Using Gamma Rays Interaction with Matter. Master's Thesis, University of Michigan, February 1959. -88.

-89&15. De Shong, J. A., Jr. Styrofoam Simulation of Boiling and Temperature Effects in the EBWR Cold Critical Experiments, AN -5697 (USAECj March 1957. 16. Marlow, K. Wo, et al. Effects of Air Voids and Fuel-Free Regions (Naval Research Laboratory Research Reactor), NRL Report 5142 (U.S. Naval Research Laboratory), June 26, 1958. 17. Bogart, D. and Hallman, T, M. Reactivity Measurements with the Bulk Shielding Reactor: Worth of Voids within the Core, NASA MEMO 12-10-58E (NASA Reactor Facility Hazards Summary, Volume II, Iewis Research Center), (November 15, 1956), 109. 18. Beckjord, E. et al. Operation of a High Performance Light Water Boiling Reactor, A/Conf. 15/P/1923 (Geneva Conference, USA) June 1958. 19. Kaufmann, S. Go Void Coefficient Measurements, ZPR-VII Memo 5 (Argonne National Laboratory), September 11, 1957. 20. Thie, J. A. and Redman, W. C. Properties of Exponential and Critical Systems of Thoria-Urania and Heavy Water and Their Application to Reactor Design, A/Conf. 15/P/600 (Geneva Conference, USA), June 1958. 21. Mattson, R. A. Preliminary Dynamic Void Measurements, ZPR-VII Memo 6, (Argonne National Laboratory), November 4 1957. 22. Eriksen, V. 0. Personal Communication. 23. Williams, R. 0. et al. Reactor Analysis of the Organic Moderated Reactor Experiment and Comparison of Experimental Results, A/Conf. 15/P/630, (Geneva Conference, USA), June 1958. 24. French, N. E., Ruane, R. Fo and Storm, Mo L. Distributed Reactivity Coefficients in a PMA Slab Core: Measurement and Analysis, KAPL-MMLS-12, (USAEC), February 9, 1960. 25. Myer, Wo E., et al. Experimental Investigations of Reactor Transients, IDO 16285, (USAEC) (April 20, 1956), 9. 26. Grund, Jo and Neal, W, J. Void Coefficients (Quarterly Progress Report), IDO 16416 (USAEC) (October 1, 1956), 30. 27. Watanabe, T. Calculated Void and Temperature Effects, Ibid., 46. 28. Grund, J. and Neal, W. J, Void Coefficient (Quarterly Progress Report), IDO 16437 (USAEC), (lMarch 21, 1958), 21. 29. Watanabe, To Effect of Voids on Reactivity, Ibid., 58.

-9,030. Grund, J. E. Average Void Coefficient (Quarterly Progress Report), IDO 16452, (USAEC), September 10, 1958. 31. Quigley, To et al. "B" Core Experimental Data (Quarterly Progress Report), IDO 16512, (USAEC), (May 6 1959), 8, 32. Stephan, L. A. Average Void Coefficient (Quarterly Progress Report), IDO 16537, (USAEC), (September 1, 1959), 13. 33. Spano, A. H. and Metcalf, D, R. Calculations of Void and Temperature Coefficients for the B-12/64 Core (Quarterly Progress Report}, IDO 16539, (USAC) (November 20, 1959), 21. 34. Lewis Research Center Staff, Reactor Physics Calculations, NASA-Memo 12-8-58E (NASA Reactor Facility Hazards Summary, Volume I, Lewis Research Center), (October 15, 1956), 149. 35. Pigford, To H., et al. "Neutron Lifetimes and Void Coefficients for Research Reactors," Preprint 154 Nuclear Engineering and Science Congress, December 1955, (American Institute of Chemical Engineers). 36. Silver, E. Go Comparison of Pool-Type Reactor Critical Experiments with Two-Group and Thirty Group Calculations, ORNL-2499, (USAEC), June 9, 1958. 37. Thie, J. Ao Theoretical Reactor Statics and Kinetics of Boiling Reactors, A/Conf. 15/P/638, (Geneva Conference, USA), June 19580 38. Webster, J. W. An Analysis of the Accuracy of Perturbation Theory, ID-16173, (USAEC) June 17, 19540 39. Behrens, Do J. "The Effect of Holes in a Reacting Material on the Passage of Neutrons," Proc. Phys. Soc. (London), 62, 358 A, (October 1949), 607. 40. Deutsch, R. W. "Computing 3-Group Constants for Neutron Diffusion," Nucleonics, 15, no 1, (January 1957), 47o 41. Amster, H. Jo A Compendium of Thermal Neutrons Cross Sections Averaged over the Spectra of Wigner and Wilkins, WAPD-185, (USAEC), January 1958. 42. Bilodeau, G. G. et al. PDQ —An IBM-704 Code to Solve the Two-Dimensional Few-Group Neutron Diffusion Equations, WAPD-TM-70, (USAEC), August 1957. 43. Marlowe, 0. J. WANDA —A One-Dimensional Few Group Diffusion Equation Code for the IBM-704, WAPD-TM-28, (USAEC), November 1956. 44. Reactor Physics Constants Center, Reactor Physics Constants, ANL-5800, (USAEC), no date.

-9145. Murray, R. L. Nuclear Reactor Physics, Prentice Hall, Inc., 1957. 46. Weinberg, A. M. and Wigner, E. P. The Physical Theory of Neutron Chain Reactors, University of Chicago Press, 1958. 47. Goldstein, H. The Attenuation of Gamma Rays and Neutrons in Reactor Shields, Division of Reactor Development (USAEC), May 1, 1957. 48. Reynolds, A. B. Reactivity Effects of Large Voids in the Reflector of a Light Water Moderated and Reflected Reactor, AECU-4391, (USAEC), 1959. 49. Worthing, A. G. and Geffner, J. Treattnent of Experimental Data, New York: John Wiley and Sons, (1943), 209. 50. Roberts, L. D. et al. The Slowing Down Distribution to Indium Resonance of Fission Neutrons from a Point Fission Source in Two Mixtures of Light Water 1:1 and 1:2 by Volume, AECD-34o8, (USAEC), February 23, 1949. 51. Deutsch, R. W.. "Fission Product Buildup in Enriched Thermal Reactors, Nucleonics, 14, n. 9, (September 1956), 89. 52. Trubey, P. K. and Lessig, R. H. A Calculation of the Infinite Medium Neutron Multiplication Factor in a U35-H20 Solution, ORNL-2609, (USAEC), (October 1958), 73. 53. Spinrad, B. I. and Kurath, D. Computation Forms for Solution of Critical Problems by Two-Group Diffusion Theory, ANL-4352, (USAEC), March 1952. 54. Shapiro, J. L. and Gomberg, H. J. "Operating a Research Reactor," Mechanical Engineering, 80, n. 10, (October 1948), 68. 55. Shapiro, J. L. et al. Initial Calibration of the Ford Nuclear Reactor, MMPP-110-1, (University of Michigan), June 1958. 56. Bookston, J. M., Davis, C. L. and Smith, B. E. GNU, A Multigroup One-Dimensional-Diffusion Program for the IBM-704, GMR-101, (General Motors Research Laboratory), -April 8, 1957. 57. Foderero, A. "An Iteration Method for the Specification of Multigroup Bucklings," Nuclear Science and Engineering, 6, n. 6, (December 1959), 514. 58. Glasstone, S. and Edlund, M. C. The Elements of Nuctar Reactor Theory, D. Van Nostrand Co., Inc., (1952), 372. 59. Calame, G. P. A Few Group Theory of Water Gap Peaking, KAPL-2000-8, (USAEC), (December 1959), A-6. 60. Lombard, D. B. and Blanchard, C. H. "Fission to Indium Age in Water," Nuclear Science and Engineering, 7, n. 5, (May 1960), 448.

APPENDIX A EFFECT OF NON-HOMOGENEITY ON VOID MEASUREMENTS WITH THE GAMMA RAY DENSITOMETER It is extremely difficult to accurately determine the error due to heterogeneity of the air-water mixture. The error is due to the fact that the density variations about the mean are weighted exponentially by the attenuated gamma ray population, but are arithmetically averaged by the counting device. To be truly accurate, the instantaneous count rate should be logarithmically averaged. The variation from mean density appears to be a statistically random function with some unknown properties. However, as will be shown, the maximum possible error due to this effect is simple to evaluate and turns out to be small enough to allow its use as an upper bound. The assumption is made that the density of the mixture varies in the worst possible way —it is either 1 or 0 (of course, the probability that this occurs is small when the volume is divided into many small channels as is the case here). Let Io be the count rate when the mixture density p = 1 and Io exp (p2o) be the count rate when p = O0. ^o is the total voidable length of the test section and i is the total attenuation coefficient. Under the conditions stated, the average count rate will be I = (l-p)Ioe0~ + PIo (A-l) where p is the actual (time) average density of the mixture. The pex measured by the attenuation method is calculated by Pe =- 1 - 1 n (A-2) -9ex I -92

-93which, expressed in terms of Equation (A-l) is; Pex = 1 - gn [p + (1-P)e e~. (A-3) I-lo Thus, the deviation from the true density is Pex-P = 1_ n [p e4- o(l1-) + (l-p)ekoP] (A-4) P PAtop This function is negative for all values of p and increases as the product (Qo increases for a given value of p. For this experiment 4Ao = 0.59. Figure A shows the variation of the error with p.

-94-26 --- -24 -22 -20 0 \ -18 -6 -4 -2 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 p, AVERAGE WATER DENSITY Figure A. Maximum Possible Error in Void Measurement Due to Temporal Fluxtuations in Density

APPENDIX B AGE CALCULATIONS IN MIXTURES OF ALUMINUM, WATER AND AIR The age of neutrons at energy Ef due to a monoenergetic source of energy Eo in a homogeneous medium may be defined as Eo D dE T C-Fs -- (B-l) Ef Zs E (B) Under the same approximation (Infinite medium, no capture) for which Equation (B-l) is accurate we can substitute 1/3 Ztr for the diffusion coefficient D. Now, if "suitably averaged" cross sections for the medium can be found, we can write T. n Eo/Ef (B-2) 3 Z s Ztr If we mix two materials, each of which having, in its pure state, cross sections ZSi and Stri. logarithmic decrement.i and an age Ti, the age may be written TIX n Eo/Ef (B-3) [~zr + E* ] 3[1 s1 + 52 S2] tr1 + r where +-V2'+ *- (B-4) and the vi are the volumes of each material. in the mixture, Expressing the age in the mixture in terms of the ages in the pure materials we find the result: 7 1 = >... =1 (1+) + + 1L T1 2 + ( 2) + 2(l l) TMIX 2 (1+JT2) 21(1-I2) (B-5) -95

-96where a v- tri = Zs(l-"i) v2 ri This derivation includes the implicit assumption that the microscopic material cross sections do not chang as other materials are added. Given two or more measurements of TMI.X as a function of a, it is possible to fit the constants to the data to derive a working formula. Measurements of the age in pure water and in aluminum-water mixtures have been made (50). Table B-l shows the results of these measurements. TABLE B-l Measured Age in A1.-H20 Mixtures a rT(cm2) 0 31o8 0.5 51o6 1.0 79.8 The age in pure water is difficult to measure because of the low age. Much theoretical work has been done to show that the age is closer to 27 cm2. New experiments have been reported() which give a result in closer agreement with theory. Due to the present uncertainty in this parameter and the fact that we are interested in the age with 0.5 < a < 1.0, the data for a = 0 is neglected. The resultant fit of Equation (B-5) to the other two points gives MX 100(,+c)2-) 2(B-6) 0.0207a2 + 1.27c + 3572 (a =v v/vW) which is equivalent to saying that TW = 26.9 cm2, (vA = 0)

-97and TA1 = 4830 cm2, (VW 0). The data SW = 0.948 1 - SW = 0.676 A1l = 0.0723 1 - TA1 = 0-9754 was obtained from reference (4). As seen in Equation (B-5) if one of the components has an infinite age (e.g., air) the correction on TMiX is simply the square of the relative density. Therefore, Equation (B-6) can be corrected for the density change of the water as follows: 100 + + TMIX =..I j (B-7).0207 - + 1.27 ) + 3.72 This correction consists of correcting the volume ratio of metal to water from a to a/p. Then, the effect of the change in density of the entire medium is taken into account by the last term in the numerator.

APPENDIX C LITERATURE SURVEY (1) Characteristics of bubbles rising in water. Newman and Whelan ) discuss the ratio of bubble volume to the area of the orifice. For "very slow" rates of formation this ratio varies between 0.2 cm. and 0.3 cm. As the bubble rate is increased, the ratio decreases then increases until a random size distribution sets in. Experimental studies of the rate of rise of air bubbles in various liquids, including water have been made over a period of fifty years. Of those listed here (Reference (2) - (8)) several(2''7' ) discuss single bubbles only whereas the others study the rise of bubbles in swarms. The results are difficult to compare due to the different experimental apparatus used. The rise velocity is a function of the distance between the walls(3). In general though a velocity with a lower limit of about 1 f.p.s,, increasing with increasing air volumes is indicated. The rate of change of velocity with air volume is.sharpest at low volumes, decreasing as volume increases. (2) Measurement of vapor or gas densities in liquids using radiation attenuation. As radiation detection instrumentation has become more reliable, attenuation methods of measuring density have become more pupular. Use is made primarily of x-rays(9) and gamma rays(11-14) In all the cases referred to here, the radiation penetration is the prime measurement of the experiment (rather than as a calibration device as used in this paper). This increases the importance of the error analysis, in terms of -98~~

-99the possibility of localized flow patterns. It is shown in these studies that large errors may arise when separated flow occurso In all but one case, no traverses of the test section were performed. Petrick(l3), in comparing the "one-shot" to traverse techniques shows the superiority of the latter, especially when localized density variation occur. In all the experiments referred to here the scintillation detector was employed. (14) Only Yanez( used a single channel differential analyzer to eliminate counting of scattered radiation. The others depended on collimation only. (3) Measurements of void effects on reactivityo Voids have been introduced or simulated in reactors using several techniques. De Shong(15) used STYROFOAM, a rigid low density expanded plastic containing small air bubbles. The average density was varied by making holes in the plastic allowing water to fill the hole, There was some seepage of water into the plastic as time passed. The results of the measurement, on the Experimental Boiling Water Reactor were felt to be adequate, falling within 10 per cent of a "modified two-group" calculation (for distributed case)o Several groups have employed air as the void using varying means (1.6-25) (16) of insertion( o3) At NRL ) fuel elements were removed and empty sealed cans were inserted, yielding positive void coefficients. This is associated with the effedt of insertion of experiments into the reactor, rather than boiling effects, The same type of air can insertion was done (17) on the BSR 7) In addition, a special fuel element was constructed in which the water in any water channel between fuel plates could be voided, Almost an entire element could be voided, one channel at a time, The

-100effect of interactions between channels was measured showing an increase of about 10 per cent in the (negative) coefficient with maximum interaction, over the measurement with minimum interaction. The latter number is -0.57 at the center of the BeO and water reflected core. No void coefficient distributions or uniform void coefficients were measured. In the Vallecitos BWR a thin, air can was used to replace water between the plates. This was moved around the strongly heterogeneous core to measure void coefficients. Calculations by two-group perturbation theory accurately predict the distributed coefficient. The local void coefficients, similarly calculated, show an error of 30 per cent to 50 per cent too high in the central position. It is also shown that the central coefficient is positive if the core diameter is 2-6 feet but negative if the diameter is reduced to 1.6 feet. Void coefficients in heavy water reactors have also been studied(19 22) In one case(9 20), hollow aluminum tubes were spaced throughout the core to simulate distributed voids. The effect of the aluminum walls was considerably larger than the effect of the D20 removal. Therefore, a dynamic type of reactivity measurement, oscillating the water in and out of one tube by means of a pluner, was attempted(21) This gave relative void coefficients with "good precision." In another type of dynamic void insertion, air bubbles were (22) introduced into a D20 reactor. The bubbling was turned on and off by means of a solenoid operated valve. The power oscillation was monitored to give the reactivity effect. The main reason given for this type of measurement was that the reactor contained xenon fiom previous high power

-101operation. Measurements of the void concentration were not made directly. In one oscillation, the total air was collected at the top, its volume measured and the void concentration calculated from this data (details of this calculation omitted). In an organic moderated reactor (23) a fuel element was closed off to the coolant, thus presenting a void condition. Calculations were attempted but were not in agreement with the measurements. At KAPL distributed reactivity coefficients for various materials were examined in the Plastic Mockup Assembly(24)) Since the moderator is solid (polyethylene), voids are made by removing strips of moderator, The measurements were compared with calculations using 3 group (2 dimension) perturbation theory with two different methods of cross section calculations. The axial buckling, assumed to be constant, was determined from measured flux distributions. The calculations show excellent agreement when the MUFT-SOFOCATE codes are used to generate cross sections, but poor agreement with other techniques for cross section generation. (4) Experimental and analytical void coefficient program of Philips Petroleum Co,, Atomic Energy Division, Arco, Idaho. An extensive program of investigation of void coefficients has been under way for several years by personnel of the BORAX and SPERT projects(25-33)o (25) Nyer(25) reports on experimental results of the SPERT-I reactor. The method of inserting voids is not given. Linearity of the distributed void condition with total void content was noted (up to 6 per cent voided), One central element was voided up to 90 per cent with a quite linear reactivity loss,

-102(26,28,30) (32) Grund ) and Stephan3) have used aluminum, lead, magnesium, inflated air bags and STYROFOAM for void simulation. The results show conclusively that lead and aluminum are poor materials due to difficulties in calculating the correction factors for absorption and scattering. Inflated metal bags developed leaks quite frequently. Magnesium required a much smaller correction than other metals but in one core (B-12/64) a large discrepency was found between STYROFOAM and magnesium. It is not clear whether it was the technique or the material at fault. The authors apparently feel that STYROFOAM is the best material for void simulation, but magnesium is used for small plate spacing cores. Stephan(32) reports that the average void coefficient as measured by moving a single void element from place to place in the core agreed with the measurement of the uniformly distributed core. Quigley summarizes the (measured) characteristics of the SPERT I cdoes. Table C gives some of these results along with the corresponding features of the FNR core. Both the SPERT A and B type fuel elements are similar to the FNR element from a nuclear standpoint. In Table C, the first number in the core designation is the number of fuel plates per element. The second number is the number of elements in the core. (27,29) Calculations by Watanabe ) were consistently high by a factor of 80-100 per cent. Although the details of the method are not supplied the calculations were evidently two-group, one-dimensional, using a constant axial buckling. Spano and Metcalf(33) report the results of a much more detailed method of calculation. The procedure is outlined as follows:

TABLE C Characteristics of SPERT-I(31) and FNR Cores B-24/32 B-16/40 B-12/64 A-17/28 FNR-18/21 Moderator Thickness (mils) 65 190/65 190 117 117 Hydrogen/Uranium Ratio 269 537 759 325 408 Metal/Water Ratio 1.14 0.63 0.46 0.79 0.711 Temperature coefficient -1.1xl0-2 -1.7x10-2 -1.8xlO -2 -067xlO 2 -l.01x102 at 20~C ($/~C) Total Moderator Volume (cm3) 5.45x104 9.05x104 1.82x105 5.40x104 4.5xl04 Central Void Coefficient -17xl0-4 -4.7x104 0.8x10-4 -9.3xl0-4 -9.3x104 $/cm3 \Ak /v -0,700 -0.321 0.110 -0.379 -0.315 k v Average Void Coefficient -7.3x10-4 -2.9x10-4 -0.93x10-4 -4.6x10-4 -5.6x10-4 $/cm3 i_/-Av -0 300 -0.198 -0.127 -0.188 -0.190 k v

-104a. Use MUFT and SOFOCATE to generate the two group cross sections. b. Apply the cross sections to a two-dimensional PDQ code cal, culation using a guess for the axial buckling. c. With the resultant X, calculate the material buckling. d. Let B2 -B2 radial B2material axial e. Perform a one-dimensional, axial WANDA calculation using the B2radial as the transverse buckling, to find X. This last value of X is assumed to be the correct keff. The entire procedure is followed once for the unvoided condition and once for every voided condition. Results are given for only two points. One, a distributed void calculation was within 3 per cent of the measured value. The other, a centrally localized void had no experimental measurement to compare. (5) Other calculation techniques. At the Lewis Research Center(34) two-group calculations were made for a core similar to that of the Bulk Shielding Reactor (BeO reflected, a = 0.733) A ore-dimensional (cylindrical) analog simulator was used for the calculations. The resultant uniformly distributed void coefficient was -.18. The transverse buckling was assumed constant and was determined by measurements on the MTR (a similar core). Only density effects were accounted for. Pigford et al(35) have calculated distributed void coefficients by the two-group perturbation method, accounting only for density changes. For a core similar to that of the FNR, but without graphite reflector the V.C. = -0.277. With a one foot thick graphite reflector the V.C, = -0.110.

-105This calculation and comparison was somewhat unrealistic since it assumed the same core radius in each case. However, the results do bracket the results of this study. A comparison of two-group and thirty group calculations with (36) experiment is reported by Silver Unfortunately, only the predicted excess reactivities were compared. No clear cut advantage of one method over another was found. It is likely that a distributed reactivity coefficient calculation may have shown the superiority of one technique or another. One method reported here was the iteration of axial buckling, as described in Section III-B (3), Results of calculations of the EBWR void coefficient are given by Thie(37) Here resonance escape probability variation is a most important feature. The modified two-group calculations used resulted in a distributed void coefficient within 10 per cent of the measured value. Calculations and experiments involving localized voids disagreed by as much as 30 per cent.

APPENDIX D SAMPLE DATA AND CALCULATIONS (1) Measurement of effective buildup factor in densitometer (Section II-C (1)) Data: Date: 12/7/59 E =50 AE = 900 Gain = 4:86 Condition Count Time (min.) Corrected Count Background 7101 10 Empty chamber 165446 5 161896 (=Co) Water inserted 72077 5 68527 (=C1) Element inserted 50906 5 47556 (=C) Water drained 87316 5 83766 (=C3) n = 18 plates tAl = 0.060 inches = 0.1524 cm. ntA = 2.7437 Calculations: Trial and error calculation is made for B, using Equation (9') with W = 0.117 cm-1, A = 0.273 cm-1 0.117 = Ln 1.3356-en 1+B[0.117(5.83)+0.273(2.744)] l 2.744' [1+B(o.117) (8. 574) ] [l+B(O. 273) (2.744)]J B = o.o880 -106

-107(2) Measurement of air volume fraction (Section II-E) Data: 9/15/59 E = 50 Barometpr: 75.5 cm. AE = 900 Water height: 21'7-1/2" (Bottom) Gain = 4:86 Air Pressure: 50 psig. Background: 1097: 2 min. 1128 i' 2 min. Densitometer Counts Pinion (1 min.) (2 min.) (1 min,) Setting F=0 F=40 F=0 4 20779 42291 20558 5 20682 43645 20864 6 21075 44093 20954 7 20993 44882 20707 8 21202 45427 21045 9 21022 45425 21019 10 20989 45635 21105 11 20969 45792 21083 12 21261 44676 21267 13 21233 44029 21020 14 21063 42846 20741 15 21198 42226 21031

-108Calculation: Pinion Corrected Counts Setting F -0 F = 40 Ratio 4 40237 41191 1.0237 5 40446 42545 1.0519 6 42029 42993 1.0229 7 40600 43782 1.0784 8 41147 44327 1.0773 9 40941 44325 1.0827 10 40994 44535 1.0864 11 40952 44692 1.0913 12 41428 43676 1.0518 13 41153 42929 1.0432 14 40704 41746 1.0256 15 41129 41126 0.9999 i 1.8482 in i o.6 45 The average air volume fraction xi of each 0.191 inch wide section is found using Equation (5"). To find the overall fraction, x, the xi are averaged over the entire water area of one element. 1 1 C. Thus x n L T 1-0.079 (0.117)17(0.117)(2.54) Coi (0.191)2.54(0.117)18(2.54) z 35.6 i x = 0.134 in H Ci - 0.0823 i Co, 1

-109P2 = 14.62 + 0.4335 (21.625) = 23.55 psia. P1 = 50 + 14.62 = 64.62 psia.. P2/P1 = 2.93 and W = XP2/.P1 = 0.241 (3) Calibration of control rod (Section II-F(2)). Data: Date: 4/1/60 Shim Rod Positions Control Rod Positions Run A B C Initial Final Period(sec.) 1 15.98 15.89 15.77 36" 31-9/16" 303 2 16.04 16.00 15.85 32-23/32" 29-1/16" 140.5 5 16.17 16.16 15.99 30-5/16" 26-23/32" 99.7 4 16.58 16.38 16.21 27-7/8" 24-7/16" 89.2 5 ___ 16.46 16.49 16.28 26-3/4" 22-3/4" 73.9 Control rod in: 36" out: 12" Calculations: Average Rod Run Position A L A k(%) Ak/AL(%/in) 1 2.22" 4.44" 0.0291.00655 2 5.11" 3.66" 0.0568.0155 3 7.48" 3,5911 0.0747.0208 4 9.84" 3.44" 0.0815.0237 5 11.25" 4.00"o 0.0939.0235

-110(4) Reactivity measurements (Section II-6) Data: Date: 4/1/60 Barometer: 75.6 cm. Pool height: 1" below trough Test Position: 35 Element condition: plugged Pressure: 50 psig. Run F Control Rod Position 0 0 31-25/32 1 10 28-29/32 2 20 27-3/4 3 50 26-3/4 4 40 25-27/32 5 50 24-31/32 Shim Rod Positions 6 60 24-7/16 A = 16.10 7 70 23-7/8 B = 15.93 8 55 24-11/16 C = 15.87 9 45 25-1/2 10 35 26-1/16 11 25 27-3/8 12 15 28-1/4 13 0 31-13/16

-111Calculations: P2/P1 = 2.96 The following table is then filled in using the densitometer calibration of Figure 12 and the control rod calibration of Figure 17. Run XP2,/1 lX(%) k(%) A k(%) 0 0 0.0285 -.0001 1.113 3.82.0773 -.0489 2.146 4.94.1019 -o0735 3.191 6.45.1246 -.0962 4.237 8.01.1487 -.1183 5.266 8.99.1680 -.1395 6.295 9.97.1802 -.1518 7.324 10.95.1931 -o1648 8.282 9.53.1745 -.1462 9.250 8.45,1550 -,1266 10.214 7.23.1410 -,1126 11.163 5.51.1100 -o0816 12.127 4.29.0910 - 0626 13 0 0.0283 +.0001

APPENDIX E REACTOR PHYSICS NOMENCLATURE The following nomenclature refers to the reactor physics notation only. The remaining symbols are defined as they appear. 0, neutron flux, adjoint function v, average number of neutrons released per fission q average number of neutrons released per absorption in fuel f, thermal utilization e- fast fission factor p, resonance escape probability k,, infinite multiplication factor (= nfEp) keff, effective multiplication factor D, diffusion coefficient L, diffusion length v, Fermi age B- geometric buckling a, microscopic atomic cross section N, atomic density Z, macroscopic cross section (= Na) X, neutron mean free path (= 1/Z), logarithmic energy decrement X, eigenvalue of diffusion equations Xi, fraction of neutrons born in energy group i -112

-113P, power level F, thermal disadvantage factor b, reflector savings Peff' effective delayed neutron fraction Subscripts a, thermal absorption R, removal tr, transport s, scattering 1, epithermal neutrons 2, thermal neutrons

UNIVERSIY OF MICHIGAN 3 9015 03525 0177